Actual source code: ex20.c

  1: static char help[] = "Test DMStag transfer operators.\n\n";

  3: #include <petscdm.h>
  4: #include <petscdmstag.h>

  6: int main(int argc, char **argv)
  7: {
  8:   DM        dmc, dmf;
  9:   PetscInt  dim;
 10:   PetscBool flg;

 12:   PetscFunctionBeginUser;
 13:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 14:   PetscCall(PetscOptionsGetInt(NULL, NULL, "-dim", &dim, &flg));
 15:   PetscCheck(flg, PETSC_COMM_WORLD, PETSC_ERR_ARG_WRONG, "Supply -dim option");
 16:   if (dim == 1) PetscCall(DMStagCreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, 3, 1, 1, DMSTAG_STENCIL_BOX, 1, NULL, &dmc));
 17:   else if (dim == 2) PetscCall(DMStagCreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, 3, 3, PETSC_DECIDE, PETSC_DECIDE, 1, 1, 1, DMSTAG_STENCIL_BOX, 1, NULL, NULL, &dmc));
 18:   else if (dim == 3) PetscCall(DMStagCreate3d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, 3, 3, 3, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, 1, 1, 1, 1, DMSTAG_STENCIL_BOX, 1, NULL, NULL, NULL, &dmc));
 19:   else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "dim must be 1, 2, or 3");
 20:   PetscCall(DMSetFromOptions(dmc));
 21:   PetscCall(DMSetUp(dmc));

 23:   /* Directly create a coarsened DM and transfer operators */
 24:   PetscCall(DMRefine(dmc, MPI_COMM_NULL, &dmf));
 25:   {
 26:     Mat       Ai;
 27:     Vec       vc, vf;
 28:     PetscInt  size;
 29:     PetscReal norm;

 31:     PetscCall(DMCreateInterpolation(dmc, dmf, &Ai, NULL));
 32:     PetscCall(MatCreateVecs(Ai, &vc, &vf));
 33:     PetscCall(VecSet(vc, 1.0));
 34:     PetscCall(MatMult(Ai, vc, vf));
 35:     PetscCall(VecGetSize(vf, &size));
 36:     PetscCall(VecNorm(vf, NORM_1, &norm));
 37:     PetscCheck((norm - size) / (PetscReal)size <= PETSC_MACHINE_EPSILON * 10.0, PetscObjectComm((PetscObject)dmc), PETSC_ERR_PLIB, "Numerical test failed");
 38:     PetscCall(MatDestroy(&Ai));
 39:     PetscCall(VecDestroy(&vc));
 40:     PetscCall(VecDestroy(&vf));
 41:   }
 42:   {
 43:     Mat       Ar;
 44:     Vec       vf, vc;
 45:     PetscInt  size;
 46:     PetscReal norm;

 48:     PetscCall(DMCreateRestriction(dmc, dmf, &Ar));
 49:     PetscCall(MatCreateVecs(Ar, &vf, &vc));
 50:     PetscCall(VecSet(vf, 1.0));
 51:     PetscCall(MatMult(Ar, vf, vc));
 52:     PetscCall(VecGetSize(vc, &size));
 53:     PetscCall(VecNorm(vc, NORM_1, &norm));
 54:     PetscCheck((norm - size) / (PetscReal)size <= PETSC_MACHINE_EPSILON * 10.0, PetscObjectComm((PetscObject)dmc), PETSC_ERR_PLIB, "Numerical test failed");
 55:     PetscCall(MatDestroy(&Ar));
 56:     PetscCall(VecDestroy(&vf));
 57:     PetscCall(VecDestroy(&vc));
 58:   }
 59:   PetscCall(DMDestroy(&dmf));

 61:   PetscCall(DMDestroy(&dmc));
 62:   PetscCall(PetscFinalize());
 63:   return 0;
 64: }

 66: /*TEST

 68:    test:
 69:       suffix: 1d
 70:       nsize: 1
 71:       args: -dim 1

 73:    test:
 74:       suffix: 1d_ratio
 75:       nsize: 1
 76:       args: -dim 1 -stag_refine_x 3
 77:       output_file: output/ex20_1d.out

 79:    test:
 80:       suffix: 1d_par
 81:       nsize: 2
 82:       args: -dim 1 -stag_grid_x 6
 83:       output_file: output/ex20_1d.out

 85:    test:
 86:       suffix: 2d
 87:       nsize: 1
 88:       args: -dim 2
 89:       output_file: output/ex20_1d.out

 91:    test:
 92:       suffix: 2d_ratio
 93:       nsize: 1
 94:       args: -dim 2 -stag_refine_x 3 -stag_refine_y 4
 95:       output_file: output/ex20_1d.out

 97:    test:
 98:       suffix: 2d_par
 99:       nsize: 4
100:       args: -dim 2 -stag_grid_x 6 -stag_grid_y 7
101:       output_file: output/ex20_1d.out

103:    test:
104:       suffix: 3d
105:       nsize: 1
106:       args: -dim 3
107:       output_file: output/ex20_1d.out

109:    test:
110:       suffix: 3d_ratio
111:       nsize: 1
112:       args: -dim 3 -stag_refine_x 3 -stag_refine_y 4 -stag_refine_z 5
113:       output_file: output/ex20_1d.out

115:    test:
116:       suffix: 3d_par
117:       nsize: 8
118:       args: -dim 3 -stag_grid_x 6 -stag_grid_y 7 -stag_grid_z 8
119:       output_file: output/ex20_1d.out

121: TEST*/