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) {
 28:       PetscViewerASCIIPrintf(viewer,"Hello from rank 0 -> vector name %s, size %D\n",name,n);
 29:     }
 30:   }
 31:   return 0;
 32: }

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

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

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

 53:   PetscInitialize(&argc,&argv,NULL,help);
 54:   DMShellCreate(PETSC_COMM_WORLD,&dm);
 55:   PetscObjectComposeFunction((PetscObject)dm,"DMSetUpGLVisViewer_C",DMSetUpGLVisViewer_Shell);
 56:   VecCreateMPI(PETSC_COMM_WORLD,1,PETSC_DECIDE,&v);
 57:   PetscObjectSetName((PetscObject)v,"seed");
 58:   VecSetOperation(v,VECOP_VIEW,(void (*)(void))VecView_Shell);
 59:   DMShellSetGlobalVector(dm,v);
 60:   VecDestroy(&v);
 61:   DMViewFromOptions(dm,NULL,"-dm_view");
 62:   DMGetGlobalVector(dm,&v);
 63:   VecViewFromOptions(v,NULL,"-vec_view");
 64:   DMRestoreGlobalVector(dm,&v);
 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*/