Actual source code: ex8.c

  1: static char help[] = "Test parallel ruotines for GLVis\n\n";

  3: #include <petscdmshell.h>
  4: #include <petsc/private/glvisvecimpl.h>

  6: PetscErrorCode VecView_Shell(Vec v, PetscViewer viewer)
  7: {
  8:   PetscViewerFormat format;
  9:   PetscBool         isglvis, isascii;

 11:   PetscViewerGetFormat(viewer, &format);
 12:   PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERGLVIS, &isglvis);
 13:   PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii);
 14:   if (isglvis) {
 15:     DM dm;

 17:     VecGetDM(v, &dm);
 18:     /* DMView() cannot be tested, as DMView_Shell defaults to VecView */
 19:     if (!dm) return 0;
 20:     VecView_GLVis(v, viewer);
 21:   } else if (isascii) {
 22:     const char *name;
 23:     PetscInt    n;

 25:     VecGetLocalSize(v, &n);
 26:     PetscObjectGetName((PetscObject)v, &name);
 27:     if (!PetscGlobalRank) PetscViewerASCIIPrintf(viewer, "Hello from rank 0 -> vector name %s, size %" PetscInt_FMT "\n", name, n);
 28:   }
 29:   return 0;
 30: }

 32: PetscErrorCode DMSetUpGLVisViewer_Shell(PetscObject odm, PetscViewer viewer)
 33: {
 34:   DM          dm = (DM)odm;
 35:   Vec         V;
 36:   PetscInt    dim      = 2;
 37:   const char *fec_type = {"testme"};

 39:   DMCreateGlobalVector(dm, &V);
 40:   PetscObjectSetName((PetscObject)V, "sample");
 41:   PetscViewerGLVisSetFields(viewer, 1, &fec_type, &dim, NULL, (PetscObject *)&V, NULL, NULL);
 42:   VecDestroy(&V);
 43:   return 0;
 44: }

 46: int main(int argc, char **argv)
 47: {
 48:   DM  dm;
 49:   Vec v;

 52:   PetscInitialize(&argc, &argv, NULL, help);
 53:   DMShellCreate(PETSC_COMM_WORLD, &dm);
 54:   PetscObjectComposeFunction((PetscObject)dm, "DMSetUpGLVisViewer_C", DMSetUpGLVisViewer_Shell);
 55:   VecCreateMPI(PETSC_COMM_WORLD, 1, PETSC_DECIDE, &v);
 56:   PetscObjectSetName((PetscObject)v, "seed");
 57:   VecSetOperation(v, VECOP_VIEW, (void (*)(void))VecView_Shell);
 58:   DMShellSetGlobalVector(dm, v);
 59:   VecDestroy(&v);
 60:   DMViewFromOptions(dm, NULL, "-dm_view");
 61:   DMGetGlobalVector(dm, &v);
 62:   VecViewFromOptions(v, NULL, "-vec_view");
 63:   DMRestoreGlobalVector(dm, &v);
 64:   PetscObjectComposeFunction((PetscObject)dm, "DMSetUpGLVisViewer_C", NULL);
 65:   DMDestroy(&dm);
 66:   PetscFinalize();
 67:   return 0;
 68: }

 70: /*TEST

 72:   test:
 73:     suffix: glvis_par
 74:     nsize: {{1 2}}
 75:     args: -dm_view glvis: -vec_view glvis:
 76:     output_file: output/ex8_glvis.out

 78: TEST*/