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: VecGetArrayRead(x, &ax);
12: VecGetArray(y, &ay);
13: VecGetLocalSize(y, &m);
14: VecGetLocalSize(x, &n);
15: for (i = 0; i < m; i++) {
16: PetscScalar xx, yy;
18: yy = 0.0;
19: for (j = 0; j < n; j++) {
20: xx = PetscPowScalarInt(ax[j], i + 1);
21: yy += xx;
22: }
23: ay[i] = yy;
24: }
25: VecRestoreArray(y, &ay);
26: VecRestoreArrayRead(x, &ax);
27: return 0;
28: }
30: int main(int argc, char **args)
31: {
32: Mat A, B;
33: Vec base;
34: PetscInt m = 3, n = 2;
37: PetscInitialize(&argc, &args, (char *)0, help);
38: PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL);
39: PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL);
40: MatCreateMFFD(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, m, n, &A);
41: MatCreateVecs(A, &base, NULL);
42: VecSet(base, 2.0);
43: MatMFFDSetFunction(A, myF, NULL);
44: MatMFFDSetBase(A, base, NULL);
45: MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
46: MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
47: MatComputeOperator(A, NULL, &B);
48: VecDestroy(&base);
49: MatDestroy(&A);
50: MatDestroy(&B);
51: PetscFinalize();
52: return 0;
53: }
55: /*TEST
57: test:
58: nsize: {{1 2 3 4}}
60: TEST*/