Actual source code: ex4.c
petsc-3.13.4 2020-08-01
1: static char help[]= "Test PetscSFFCompose when the ilocal array is not the identity\n\n";
3: #include <petsc.h>
4: #include <petscsf.h>
6: int main(int argc, char **argv)
7: {
8: PetscErrorCode ierr;
9: PetscSF sfA, sfB, sfBA;
10: PetscInt nrootsA, nleavesA, nrootsB, nleavesB;
11: PetscInt *ilocalA, *ilocalB;
12: PetscSFNode *iremoteA, *iremoteB;
13: Vec a, b, ba;
14: const PetscScalar *arrayR;
15: PetscScalar *arrayW;
16: PetscMPIInt size;
17: PetscInt i;
18: PetscInt maxleafB;
19: PetscBool flag = PETSC_FALSE;
21: PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr;
22: MPI_Comm_size(PETSC_COMM_WORLD,&size);
24: if (size != 1) SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_USER, "Only coded for one MPI process");
26: PetscSFCreate(PETSC_COMM_WORLD, &sfA);
27: PetscSFCreate(PETSC_COMM_WORLD, &sfB);
28: PetscSFSetFromOptions(sfA);
29: PetscSFSetFromOptions(sfB);
31: PetscOptionsGetBool(NULL,NULL,"-sparse_sfB",&flag,NULL);
33: if (flag) {
34: /* sfA permutes indices, sfB has sparse leaf space. */
35: nrootsA = 3;
36: nleavesA = 3;
37: nrootsB = 3;
38: nleavesB = 2;
39: } else {
40: /* sfA reverses indices, sfB is identity */
41: nrootsA = nrootsB = nleavesA = nleavesB = 4;
42: }
43: PetscMalloc1(nleavesA, &ilocalA);
44: PetscMalloc1(nleavesA, &iremoteA);
45: PetscMalloc1(nleavesB, &ilocalB);
46: PetscMalloc1(nleavesB, &iremoteB);
48: for (i = 0; i < nleavesA; i++) {
49: iremoteA[i].rank = 0;
50: iremoteA[i].index = i;
51: if (flag) {
52: ilocalA[i] = (i + 1) % nleavesA;
53: } else {
54: ilocalA[i] = nleavesA - i - 1;
55: }
56: }
58: for (i = 0; i < nleavesB; i++) {
59: iremoteB[i].rank = 0;
60: if (flag) {
61: ilocalB[i] = nleavesB - i;
62: iremoteB[i].index = nleavesB - i - 1;
63: } else {
64: ilocalB[i] = i;
65: iremoteB[i].index = i;
66: }
67: }
69: PetscSFSetGraph(sfA, nrootsA, nleavesA, ilocalA, PETSC_OWN_POINTER, iremoteA, PETSC_OWN_POINTER);
70: PetscSFSetGraph(sfB, nrootsB, nleavesB, ilocalB, PETSC_OWN_POINTER, iremoteB, PETSC_OWN_POINTER);
71: PetscSFSetUp(sfA);
72: PetscSFSetUp(sfB);
73: PetscObjectSetName((PetscObject)sfA, "sfA");
74: PetscObjectSetName((PetscObject)sfB, "sfB");
76: VecCreateSeq(PETSC_COMM_WORLD, nrootsA, &a);
77: VecCreateSeq(PETSC_COMM_WORLD, nleavesA, &b);
78: PetscSFGetLeafRange(sfB, NULL, &maxleafB);
79: VecCreateSeq(PETSC_COMM_WORLD, maxleafB+1, &ba);
80: VecGetArray(a, &arrayW);
81: for ( i = 0; i < nrootsA; i++ ) {
82: arrayW[i] = (PetscScalar)i;
83: }
84: VecRestoreArray(a, &arrayW);
86: PetscPrintf(PETSC_COMM_WORLD, "Initial Vec A\n");
87: VecView(a, NULL);
88: VecGetArrayRead(a, &arrayR);
89: VecGetArray(b, &arrayW);
91: PetscSFBcastBegin(sfA, MPIU_SCALAR, arrayR, arrayW);
92: PetscSFBcastEnd(sfA, MPIU_SCALAR, arrayR, arrayW);
93: VecRestoreArray(b, &arrayW);
94: VecRestoreArrayRead(a, &arrayR);
95: PetscPrintf(PETSC_COMM_WORLD, "\nBroadcast A->B over sfA\n");
96: VecView(b, NULL);
98: VecGetArrayRead(b, &arrayR);
99: VecGetArray(ba, &arrayW);
100: arrayW[0] = 10.0; /* Not touched by bcast */
101: PetscSFBcastBegin(sfB, MPIU_SCALAR, arrayR, arrayW);
102: PetscSFBcastEnd(sfB, MPIU_SCALAR, arrayR, arrayW);
103: VecRestoreArrayRead(b, &arrayR);
104: VecRestoreArray(ba, &arrayW);
106: PetscPrintf(PETSC_COMM_WORLD, "\nBroadcast B->BA over sfB\n");
107: VecView(ba, NULL);
109: PetscSFCompose(sfA, sfB, &sfBA);
110: PetscSFSetFromOptions(sfBA);
111: PetscObjectSetName((PetscObject)sfBA, "(sfB o sfA)");
112: VecGetArrayRead(a, &arrayR);
113: VecGetArray(ba, &arrayW);
114: arrayW[0] = 11.0; /* Not touched by bcast */
115: PetscSFBcastBegin(sfBA, MPIU_SCALAR, arrayR, arrayW);
116: PetscSFBcastEnd(sfBA, MPIU_SCALAR, arrayR, arrayW);
117: VecRestoreArray(ba, &arrayW);
118: VecRestoreArrayRead(a, &arrayR);
119: PetscPrintf(PETSC_COMM_WORLD, "\nBroadcast A->BA over sfBA (sfB o sfA)\n");
120: VecView(ba, NULL);
122: VecDestroy(&ba);
123: VecDestroy(&b);
124: VecDestroy(&a);
126: PetscSFView(sfA, NULL);
127: PetscSFView(sfB, NULL);
128: PetscSFView(sfBA, NULL);
129: PetscSFDestroy(&sfA);
130: PetscSFDestroy(&sfB);
131: PetscSFDestroy(&sfBA);
133: PetscFinalize();
135: return ierr;
136: }
138: /*TEST
140: test:
141: suffix: 1
143: test:
144: suffix: 2
145: filter: grep -v "type" | grep -v "sort"
146: args: -sparse_sfB
148: test:
149: suffix: 2_window
150: filter: grep -v "type" | grep -v "sort"
151: output_file: output/ex4_2.out
152: args: -sparse_sfB -sf_type window -sf_window_sync {{fence active lock}} -sf_window_flavor {{create dynamic allocate}}
153: requires: define(PETSC_HAVE_MPI_ONE_SIDED)
155: # The nightly test suite with MPICH uses ch3:sock, which is broken when winsize == 0 in some of the processes
156: test:
157: suffix: 2_window_shared
158: filter: grep -v "type" | grep -v "sort"
159: output_file: output/ex4_2.out
160: args: -sparse_sfB -sf_type window -sf_window_sync {{fence active lock}} -sf_window_flavor shared
161: requires: define(PETSC_HAVE_MPI_PROCESS_SHARED_MEMORY) !define(PETSC_HAVE_MPICH_NUMVERSION) define(PETSC_HAVE_MPI_ONE_SIDED)
163: TEST*/