1c4762a1bSJed Brown static char help[] = "Test MATMFFD for the rectangular case\n\n";
2c4762a1bSJed Brown
3c4762a1bSJed Brown #include <petscmat.h>
4c4762a1bSJed Brown
myF(PetscCtx ctx,Vec x,Vec y)5*2a8381b2SBarry Smith static PetscErrorCode myF(PetscCtx ctx, Vec x, Vec y)
6d71ae5a4SJacob Faibussowitsch {
7c4762a1bSJed Brown const PetscScalar *ax;
8c4762a1bSJed Brown PetscScalar *ay;
9c4762a1bSJed Brown PetscInt i, j, m, n;
10c4762a1bSJed Brown
11c4762a1bSJed Brown PetscFunctionBegin;
129566063dSJacob Faibussowitsch PetscCall(VecGetArrayRead(x, &ax));
139566063dSJacob Faibussowitsch PetscCall(VecGetArray(y, &ay));
149566063dSJacob Faibussowitsch PetscCall(VecGetLocalSize(y, &m));
159566063dSJacob Faibussowitsch PetscCall(VecGetLocalSize(x, &n));
16c4762a1bSJed Brown for (i = 0; i < m; i++) {
17c4762a1bSJed Brown PetscScalar xx, yy;
18c4762a1bSJed Brown
19c4762a1bSJed Brown yy = 0.0;
20c4762a1bSJed Brown for (j = 0; j < n; j++) {
21c4762a1bSJed Brown xx = PetscPowScalarInt(ax[j], i + 1);
22c4762a1bSJed Brown yy += xx;
23c4762a1bSJed Brown }
24c4762a1bSJed Brown ay[i] = yy;
25c4762a1bSJed Brown }
269566063dSJacob Faibussowitsch PetscCall(VecRestoreArray(y, &ay));
279566063dSJacob Faibussowitsch PetscCall(VecRestoreArrayRead(x, &ax));
283ba16761SJacob Faibussowitsch PetscFunctionReturn(PETSC_SUCCESS);
29c4762a1bSJed Brown }
30c4762a1bSJed Brown
main(int argc,char ** args)31d71ae5a4SJacob Faibussowitsch int main(int argc, char **args)
32d71ae5a4SJacob Faibussowitsch {
33c4762a1bSJed Brown Mat A, B;
34c4762a1bSJed Brown Vec base;
35c4762a1bSJed Brown PetscInt m = 3, n = 2;
36c4762a1bSJed Brown
37327415f7SBarry Smith PetscFunctionBeginUser;
38c8025a54SPierre Jolivet PetscCall(PetscInitialize(&argc, &args, NULL, help));
399566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL));
409566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL));
419566063dSJacob Faibussowitsch PetscCall(MatCreateMFFD(PETSC_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE, m, n, &A));
429566063dSJacob Faibussowitsch PetscCall(MatCreateVecs(A, &base, NULL));
439566063dSJacob Faibussowitsch PetscCall(VecSet(base, 2.0));
449566063dSJacob Faibussowitsch PetscCall(MatMFFDSetFunction(A, myF, NULL));
459566063dSJacob Faibussowitsch PetscCall(MatMFFDSetBase(A, base, NULL));
469566063dSJacob Faibussowitsch PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
479566063dSJacob Faibussowitsch PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
489566063dSJacob Faibussowitsch PetscCall(MatComputeOperator(A, NULL, &B));
499566063dSJacob Faibussowitsch PetscCall(VecDestroy(&base));
509566063dSJacob Faibussowitsch PetscCall(MatDestroy(&A));
519566063dSJacob Faibussowitsch PetscCall(MatDestroy(&B));
529566063dSJacob Faibussowitsch PetscCall(PetscFinalize());
53b122ec5aSJacob Faibussowitsch return 0;
54c4762a1bSJed Brown }
55c4762a1bSJed Brown
56c4762a1bSJed Brown /*TEST
57c4762a1bSJed Brown
58c4762a1bSJed Brown test:
59c4762a1bSJed Brown nsize: {{1 2 3 4}}
603886731fSPierre Jolivet output_file: output/empty.out
61c4762a1bSJed Brown
62c4762a1bSJed Brown TEST*/
63