1c4762a1bSJed Brown static char help[] = "Tests MAIJ matrix for large DOF\n\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscdm.h> 4c4762a1bSJed Brown #include <petscdmda.h> 5c4762a1bSJed Brown 6d71ae5a4SJacob Faibussowitsch int main(int argc, char *argv[]) 7d71ae5a4SJacob Faibussowitsch { 8c4762a1bSJed Brown Mat M; 9c4762a1bSJed Brown Vec x, y; 10c4762a1bSJed Brown DM da, daf; 11c4762a1bSJed Brown 12327415f7SBarry Smith PetscFunctionBeginUser; 139566063dSJacob Faibussowitsch PetscCall(PetscInitialize(&argc, &argv, 0, help)); 149566063dSJacob Faibussowitsch PetscCall(DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_STAR, 4, 5, PETSC_DECIDE, PETSC_DECIDE, 41, 1, 0, 0, &da)); 159566063dSJacob Faibussowitsch PetscCall(DMSetFromOptions(da)); 169566063dSJacob Faibussowitsch PetscCall(DMSetUp(da)); 179566063dSJacob Faibussowitsch PetscCall(DMRefine(da, PETSC_COMM_WORLD, &daf)); 189566063dSJacob Faibussowitsch PetscCall(DMCreateInterpolation(da, daf, &M, NULL)); 199566063dSJacob Faibussowitsch PetscCall(DMCreateGlobalVector(da, &x)); 209566063dSJacob Faibussowitsch PetscCall(DMCreateGlobalVector(daf, &y)); 21c4762a1bSJed Brown 229566063dSJacob Faibussowitsch PetscCall(MatMult(M, x, y)); 239566063dSJacob Faibussowitsch PetscCall(MatMultTranspose(M, y, x)); 249566063dSJacob Faibussowitsch PetscCall(DMDestroy(&da)); 259566063dSJacob Faibussowitsch PetscCall(DMDestroy(&daf)); 269566063dSJacob Faibussowitsch PetscCall(VecDestroy(&x)); 279566063dSJacob Faibussowitsch PetscCall(VecDestroy(&y)); 289566063dSJacob Faibussowitsch PetscCall(MatDestroy(&M)); 299566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 30b122ec5aSJacob Faibussowitsch return 0; 31c4762a1bSJed Brown } 32c4762a1bSJed Brown 33c4762a1bSJed Brown /*TEST 34c4762a1bSJed Brown 35c4762a1bSJed Brown test: 36c4762a1bSJed Brown nsize: 2 37*3886731fSPierre Jolivet output_file: output/empty.out 38c4762a1bSJed Brown 39c4762a1bSJed Brown TEST*/ 40