Actual source code: ex229.c
1: static char help[] = "Test MATMFFD for the rectangular case\n\n";
3: #include <petscmat.h>
5: static PetscErrorCode myF(void *ctx, Vec x, Vec y)
6: {
7: const PetscScalar *ax;
8: PetscScalar *ay;
9: PetscInt i, j, m, n;
11: PetscFunctionBegin;
12: PetscCall(VecGetArrayRead(x, &ax));
13: PetscCall(VecGetArray(y, &ay));
14: PetscCall(VecGetLocalSize(y, &m));
15: PetscCall(VecGetLocalSize(x, &n));
16: for (i = 0; i < m; i++) {
17: PetscScalar xx, yy;
19: yy = 0.0;
20: for (j = 0; j < n; j++) {
21: xx = PetscPowScalarInt(ax[j], i + 1);
22: yy += xx;
23: }
24: ay[i] = yy;
25: }
26: PetscCall(VecRestoreArray(y, &ay));
27: PetscCall(VecRestoreArrayRead(x, &ax));
28: PetscFunctionReturn(PETSC_SUCCESS);
29: }
31: int main(int argc, char **args)
32: {
33: Mat A, B;
34: Vec base;
35: PetscInt m = 3, n = 2;
37: PetscFunctionBeginUser;
38: PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
39: PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL));
40: PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL));
41: PetscCall(MatCreateMFFD(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, m, n, &A));
42: PetscCall(MatCreateVecs(A, &base, NULL));
43: PetscCall(VecSet(base, 2.0));
44: PetscCall(MatMFFDSetFunction(A, myF, NULL));
45: PetscCall(MatMFFDSetBase(A, base, NULL));
46: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
47: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
48: PetscCall(MatComputeOperator(A, NULL, &B));
49: PetscCall(VecDestroy(&base));
50: PetscCall(MatDestroy(&A));
51: PetscCall(MatDestroy(&B));
52: PetscCall(PetscFinalize());
53: return 0;
54: }
56: /*TEST
58: test:
59: nsize: {{1 2 3 4}}
61: TEST*/