Actual source code: ex4.c

  1: static char help[]= "Test PetscSFFCompose when the ilocal array is not the identity\n\n";

  3: #include <petscsf.h>

  5: int main(int argc, char **argv)
  6: {
  7:   PetscSF            sfA, sfB, sfBA;
  8:   PetscInt           nrootsA, nleavesA, nrootsB, nleavesB;
  9:   PetscInt          *ilocalA, *ilocalB;
 10:   PetscSFNode       *iremoteA, *iremoteB;
 11:   Vec                a, b, ba;
 12:   const PetscScalar *arrayR;
 13:   PetscScalar       *arrayW;
 14:   PetscMPIInt        size;
 15:   PetscInt           i;
 16:   PetscInt           maxleafB;
 17:   PetscBool          flag = PETSC_FALSE;

 19:   PetscInitialize(&argc,&argv,NULL,help);
 20:   MPI_Comm_size(PETSC_COMM_WORLD,&size);

 23:   PetscSFCreate(PETSC_COMM_WORLD, &sfA);
 24:   PetscSFCreate(PETSC_COMM_WORLD, &sfB);
 25:   PetscSFSetFromOptions(sfA);
 26:   PetscSFSetFromOptions(sfB);

 28:   PetscOptionsGetBool(NULL,NULL,"-sparse_sfB",&flag,NULL);

 30:   if (flag) {
 31:     /* sfA permutes indices, sfB has sparse leaf space. */
 32:     nrootsA = 3;
 33:     nleavesA = 3;
 34:     nrootsB = 3;
 35:     nleavesB = 2;
 36:   } else {
 37:     /* sfA reverses indices, sfB is identity */
 38:     nrootsA = nrootsB = nleavesA = nleavesB = 4;
 39:   }
 40:   PetscMalloc1(nleavesA, &ilocalA);
 41:   PetscMalloc1(nleavesA, &iremoteA);
 42:   PetscMalloc1(nleavesB, &ilocalB);
 43:   PetscMalloc1(nleavesB, &iremoteB);

 45:   for (i = 0; i < nleavesA; i++) {
 46:     iremoteA[i].rank = 0;
 47:     iremoteA[i].index = i;
 48:     if (flag) {
 49:       ilocalA[i] = (i + 1) % nleavesA;
 50:     } else {
 51:       ilocalA[i] = nleavesA - i - 1;
 52:     }
 53:   }

 55:   for (i = 0; i < nleavesB; i++) {
 56:     iremoteB[i].rank = 0;
 57:     if (flag) {
 58:       ilocalB[i] = nleavesB - i;
 59:       iremoteB[i].index = nleavesB - i - 1;
 60:     } else {
 61:       ilocalB[i] = i;
 62:       iremoteB[i].index = i;
 63:     }
 64:   }

 66:   PetscSFSetGraph(sfA, nrootsA, nleavesA, ilocalA, PETSC_OWN_POINTER, iremoteA, PETSC_OWN_POINTER);
 67:   PetscSFSetGraph(sfB, nrootsB, nleavesB, ilocalB, PETSC_OWN_POINTER, iremoteB, PETSC_OWN_POINTER);
 68:   PetscSFSetUp(sfA);
 69:   PetscSFSetUp(sfB);
 70:   PetscObjectSetName((PetscObject)sfA, "sfA");
 71:   PetscObjectSetName((PetscObject)sfB, "sfB");

 73:   VecCreateSeq(PETSC_COMM_WORLD, nrootsA, &a);
 74:   VecCreateSeq(PETSC_COMM_WORLD, nleavesA, &b);
 75:   PetscSFGetLeafRange(sfB, NULL, &maxleafB);
 76:   VecCreateSeq(PETSC_COMM_WORLD, maxleafB+1, &ba);
 77:   VecGetArray(a, &arrayW);
 78:   for (i = 0; i < nrootsA; i++) {
 79:     arrayW[i] = (PetscScalar)i;
 80:   }
 81:   VecRestoreArray(a, &arrayW);

 83:   PetscPrintf(PETSC_COMM_WORLD, "Initial Vec A\n");
 84:   VecView(a, NULL);
 85:   VecGetArrayRead(a, &arrayR);
 86:   VecGetArray(b, &arrayW);

 88:   PetscSFBcastBegin(sfA, MPIU_SCALAR, arrayR, arrayW,MPI_REPLACE);
 89:   PetscSFBcastEnd(sfA, MPIU_SCALAR, arrayR, arrayW,MPI_REPLACE);
 90:   VecRestoreArray(b, &arrayW);
 91:   VecRestoreArrayRead(a, &arrayR);
 92:   PetscPrintf(PETSC_COMM_WORLD, "\nBroadcast A->B over sfA\n");
 93:   VecView(b, NULL);

 95:   VecGetArrayRead(b, &arrayR);
 96:   VecGetArray(ba, &arrayW);
 97:   arrayW[0] = 10.0;             /* Not touched by bcast */
 98:   PetscSFBcastBegin(sfB, MPIU_SCALAR, arrayR, arrayW,MPI_REPLACE);
 99:   PetscSFBcastEnd(sfB, MPIU_SCALAR, arrayR, arrayW,MPI_REPLACE);
100:   VecRestoreArrayRead(b, &arrayR);
101:   VecRestoreArray(ba, &arrayW);

103:   PetscPrintf(PETSC_COMM_WORLD, "\nBroadcast B->BA over sfB\n");
104:   VecView(ba, NULL);

106:   PetscSFCompose(sfA, sfB, &sfBA);
107:   PetscSFSetFromOptions(sfBA);
108:   PetscObjectSetName((PetscObject)sfBA, "(sfB o sfA)");
109:   VecGetArrayRead(a, &arrayR);
110:   VecGetArray(ba, &arrayW);
111:   arrayW[0] = 11.0;             /* Not touched by bcast */
112:   PetscSFBcastBegin(sfBA, MPIU_SCALAR, arrayR, arrayW,MPI_REPLACE);
113:   PetscSFBcastEnd(sfBA, MPIU_SCALAR, arrayR, arrayW,MPI_REPLACE);
114:   VecRestoreArray(ba, &arrayW);
115:   VecRestoreArrayRead(a, &arrayR);
116:   PetscPrintf(PETSC_COMM_WORLD, "\nBroadcast A->BA over sfBA (sfB o sfA)\n");
117:   VecView(ba, NULL);

119:   VecDestroy(&ba);
120:   VecDestroy(&b);
121:   VecDestroy(&a);

123:   PetscSFView(sfA, NULL);
124:   PetscSFView(sfB, NULL);
125:   PetscSFView(sfBA, NULL);
126:   PetscSFDestroy(&sfA);
127:   PetscSFDestroy(&sfB);
128:   PetscSFDestroy(&sfBA);

130:   PetscFinalize();
131:   return 0;
132: }

134: /*TEST

136:    test:
137:      suffix: 1

139:    test:
140:      suffix: 2
141:      filter: grep -v "type" | grep -v "sort"
142:      args: -sparse_sfB

144:    test:
145:      suffix: 2_window
146:      filter: grep -v "type" | grep -v "sort"
147:      output_file: output/ex4_2.out
148:      args: -sparse_sfB -sf_type window -sf_window_sync {{fence active lock}} -sf_window_flavor {{create dynamic allocate}}
149:      requires: defined(PETSC_HAVE_MPI_ONE_SIDED) defined(PETSC_HAVE_MPI_FEATURE_DYNAMIC_WINDOW)

151:    # The nightly test suite with MPICH uses ch3:sock, which is broken when winsize == 0 in some of the processes
152:    test:
153:      suffix: 2_window_shared
154:      filter: grep -v "type" | grep -v "sort"
155:      output_file: output/ex4_2.out
156:      args: -sparse_sfB -sf_type window -sf_window_sync {{fence active lock}} -sf_window_flavor shared
157:      requires: defined(PETSC_HAVE_MPI_PROCESS_SHARED_MEMORY) !defined(PETSC_HAVE_MPICH_NUMVERSION) defined(PETSC_HAVE_MPI_ONE_SIDED)

159: TEST*/