Actual source code: ex15.c


  2: static char help[] = "Tests DMDA interpolation.\n\n";

  4: #include <petscdm.h>
  5: #include <petscdmda.h>

  7: int main(int argc, char **argv)
  8: {
  9:   PetscInt       M1 = 3, M2, dof = 1, s = 1, ratio = 2, dim = 1;
 10:   DM             da_c, da_f;
 11:   Vec            v_c, v_f;
 12:   Mat            Interp;
 13:   PetscScalar    one = 1.0;
 14:   PetscBool      pt;
 15:   DMBoundaryType bx = DM_BOUNDARY_NONE, by = DM_BOUNDARY_NONE, bz = DM_BOUNDARY_NONE;

 18:   PetscInitialize(&argc, &argv, (char *)0, help);
 19:   PetscOptionsGetInt(NULL, NULL, "-dim", &dim, NULL);
 20:   PetscOptionsGetInt(NULL, NULL, "-M", &M1, NULL);
 21:   PetscOptionsGetInt(NULL, NULL, "-stencil_width", &s, NULL);
 22:   PetscOptionsGetInt(NULL, NULL, "-ratio", &ratio, NULL);
 23:   PetscOptionsGetInt(NULL, NULL, "-dof", &dof, NULL);
 24:   PetscOptionsGetBool(NULL, NULL, "-periodic", (PetscBool *)&pt, NULL);

 26:   if (pt) {
 27:     if (dim > 0) bx = DM_BOUNDARY_PERIODIC;
 28:     if (dim > 1) by = DM_BOUNDARY_PERIODIC;
 29:     if (dim > 2) bz = DM_BOUNDARY_PERIODIC;
 30:   }
 31:   if (bx == DM_BOUNDARY_NONE) {
 32:     M2 = ratio * (M1 - 1) + 1;
 33:   } else {
 34:     M2 = ratio * M1;
 35:   }

 37:   /* Set up the array */
 38:   if (dim == 1) {
 39:     DMDACreate1d(PETSC_COMM_WORLD, bx, M1, dof, s, NULL, &da_c);
 40:     DMDACreate1d(PETSC_COMM_WORLD, bx, M2, dof, s, NULL, &da_f);
 41:   } else if (dim == 2) {
 42:     DMDACreate2d(PETSC_COMM_WORLD, bx, by, DMDA_STENCIL_BOX, M1, M1, PETSC_DECIDE, PETSC_DECIDE, dof, s, NULL, NULL, &da_c);
 43:     DMDACreate2d(PETSC_COMM_WORLD, bx, by, DMDA_STENCIL_BOX, M2, M2, PETSC_DECIDE, PETSC_DECIDE, dof, s, NULL, NULL, &da_f);
 44:   } else if (dim == 3) {
 45:     DMDACreate3d(PETSC_COMM_WORLD, bx, by, bz, DMDA_STENCIL_BOX, M1, M1, M1, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, dof, s, NULL, NULL, NULL, &da_c);
 46:     DMDACreate3d(PETSC_COMM_WORLD, bx, by, bz, DMDA_STENCIL_BOX, M2, M2, M2, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, dof, s, NULL, NULL, NULL, &da_f);
 47:   } else SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "dim must be 1,2, or 3");
 48:   DMSetFromOptions(da_c);
 49:   DMSetUp(da_c);
 50:   DMSetFromOptions(da_f);
 51:   DMSetUp(da_f);

 53:   DMCreateGlobalVector(da_c, &v_c);
 54:   DMCreateGlobalVector(da_f, &v_f);

 56:   VecSet(v_c, one);
 57:   DMCreateInterpolation(da_c, da_f, &Interp, NULL);
 58:   MatMult(Interp, v_c, v_f);
 59:   VecView(v_f, PETSC_VIEWER_STDOUT_WORLD);
 60:   MatMultTranspose(Interp, v_f, v_c);
 61:   VecView(v_c, PETSC_VIEWER_STDOUT_WORLD);

 63:   MatDestroy(&Interp);
 64:   VecDestroy(&v_c);
 65:   DMDestroy(&da_c);
 66:   VecDestroy(&v_f);
 67:   DMDestroy(&da_f);
 68:   PetscFinalize();
 69:   return 0;
 70: }