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: PetscFunctionBegin;
12: PetscCall(PetscViewerGetFormat(viewer, &format));
13: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERGLVIS, &isglvis));
14: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
15: if (isglvis) {
16: DM dm;
18: PetscCall(VecGetDM(v, &dm));
19: /* DMView() cannot be tested, as DMView_Shell defaults to VecView */
20: if (!dm) PetscFunctionReturn(PETSC_SUCCESS);
21: PetscCall(VecView_GLVis(v, viewer));
22: } else if (isascii) {
23: const char *name;
24: PetscInt n;
26: PetscCall(VecGetLocalSize(v, &n));
27: PetscCall(PetscObjectGetName((PetscObject)v, &name));
28: if (!PetscGlobalRank) PetscCall(PetscViewerASCIIPrintf(viewer, "Hello from rank 0 -> vector name %s, size %" PetscInt_FMT "\n", name, n));
29: }
30: PetscFunctionReturn(PETSC_SUCCESS);
31: }
33: PetscErrorCode DMSetUpGLVisViewer_Shell(PetscObject odm, PetscViewer viewer)
34: {
35: DM dm = (DM)odm;
36: Vec V;
37: PetscInt dim = 2;
38: const char *fec_type = {"testme"};
40: PetscFunctionBegin;
41: PetscCall(DMCreateGlobalVector(dm, &V));
42: PetscCall(PetscObjectSetName((PetscObject)V, "sample"));
43: PetscCall(PetscViewerGLVisSetFields(viewer, 1, &fec_type, &dim, NULL, (PetscObject *)&V, NULL, NULL));
44: PetscCall(VecDestroy(&V));
45: PetscFunctionReturn(PETSC_SUCCESS);
46: }
48: int main(int argc, char **argv)
49: {
50: DM dm;
51: Vec v;
53: PetscFunctionBeginUser;
54: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
55: PetscCall(DMShellCreate(PETSC_COMM_WORLD, &dm));
56: PetscCall(PetscObjectComposeFunction((PetscObject)dm, "DMSetUpGLVisViewer_C", DMSetUpGLVisViewer_Shell));
57: PetscCall(VecCreateMPI(PETSC_COMM_WORLD, 1, PETSC_DECIDE, &v));
58: PetscCall(PetscObjectSetName((PetscObject)v, "seed"));
59: PetscCall(VecSetOperation(v, VECOP_VIEW, (void (*)(void))VecView_Shell));
60: PetscCall(DMShellSetGlobalVector(dm, v));
61: PetscCall(VecDestroy(&v));
62: PetscCall(DMViewFromOptions(dm, NULL, "-dm_view"));
63: PetscCall(DMGetGlobalVector(dm, &v));
64: PetscCall(VecViewFromOptions(v, NULL, "-vec_view"));
65: PetscCall(DMRestoreGlobalVector(dm, &v));
66: PetscCall(PetscObjectComposeFunction((PetscObject)dm, "DMSetUpGLVisViewer_C", NULL));
67: PetscCall(DMDestroy(&dm));
68: PetscCall(PetscFinalize());
69: return 0;
70: }
72: /*TEST
74: test:
75: suffix: glvis_par
76: nsize: {{1 2}}
77: args: -dm_view glvis: -vec_view glvis:
78: output_file: output/ex8_glvis.out
80: TEST*/