Actual source code: ex1.c


  2: static char help[] = "VecTagger interface routines.\n\n";

  4: #include <petscvec.h>

  6: static PetscErrorCode ISGetBlockGlobalIS(IS is, Vec vec, PetscInt bs, IS *isBlockGlobal)
  7: {
  8:   const PetscInt *idxin;
  9:   PetscInt       *idxout, i, n, rstart;
 10:   PetscLayout     map;

 12:   PetscFunctionBegin;

 14:   PetscCall(VecGetLayout(vec, &map));
 15:   rstart = map->rstart / bs;
 16:   PetscCall(ISGetLocalSize(is, &n));
 17:   PetscCall(PetscMalloc1(n, &idxout));
 18:   PetscCall(ISGetIndices(is, &idxin));
 19:   for (i = 0; i < n; i++) idxout[i] = rstart + idxin[i];
 20:   PetscCall(ISRestoreIndices(is, &idxin));
 21:   PetscCall(ISCreateBlock(PetscObjectComm((PetscObject)vec), bs, n, idxout, PETSC_OWN_POINTER, isBlockGlobal));
 22:   PetscFunctionReturn(PETSC_SUCCESS);
 23: }

 25: int main(int argc, char **argv)
 26: {
 27:   Vec           vec, tagged, untagged;
 28:   VecScatter    taggedScatter, untaggedScatter;
 29:   PetscInt      bs;
 30:   PetscInt      n, nloc, nint, i, j, k, localStart, localEnd, ntagged, nuntagged;
 31:   MPI_Comm      comm;
 32:   VecTagger     tagger;
 33:   PetscScalar  *array;
 34:   PetscRandom   rand;
 35:   VecTaggerBox *defaultBox;
 36:   VecTaggerBox *boxes;
 37:   IS            is, isBlockGlobal, isComp;
 38:   PetscBool     listed;

 40:   PetscFunctionBeginUser;
 41:   PetscCall(PetscInitialize(&argc, &argv, NULL, help));
 42:   n    = 10.;
 43:   bs   = 1;
 44:   comm = PETSC_COMM_WORLD;
 45:   PetscOptionsBegin(comm, "", "VecTagger Test Options", "Vec");
 46:   PetscCall(PetscOptionsInt("-bs", "The block size of the vector", "ex1.c", bs, &bs, NULL));
 47:   PetscCall(PetscOptionsInt("-n", "The size of the vector (in blocks)", "ex1.c", n, &n, NULL));
 48:   PetscOptionsEnd();

 50:   PetscCall(PetscRandomCreate(comm, &rand));
 51:   PetscCall(PetscRandomSetFromOptions(rand));

 53:   PetscCall(VecCreate(comm, &vec));
 54:   PetscCall(PetscObjectSetName((PetscObject)vec, "Vec to Tag"));
 55:   PetscCall(VecSetBlockSize(vec, bs));
 56:   PetscCall(VecSetSizes(vec, PETSC_DECIDE, n));
 57:   PetscCall(VecSetUp(vec));
 58:   PetscCall(VecGetLocalSize(vec, &nloc));
 59:   PetscCall(VecGetArray(vec, &array));
 60:   for (i = 0; i < nloc; i++) {
 61:     PetscScalar val;

 63:     PetscCall(PetscRandomGetValue(rand, &val));
 64:     array[i] = val;
 65:   }
 66:   PetscCall(VecRestoreArray(vec, &array));
 67:   PetscCall(PetscRandomDestroy(&rand));
 68:   PetscCall(VecViewFromOptions(vec, NULL, "-vec_view"));

 70:   PetscCall(VecTaggerCreate(comm, &tagger));
 71:   PetscCall(VecTaggerSetBlockSize(tagger, bs));
 72:   PetscCall(VecTaggerSetType(tagger, VECTAGGERABSOLUTE));
 73:   PetscCall(PetscMalloc1(bs, &defaultBox));
 74:   for (i = 0; i < bs; i++) {
 75: #if !defined(PETSC_USE_COMPLEX)
 76:     defaultBox[i].min = 0.1;
 77:     defaultBox[i].max = 1.5;
 78: #else
 79:     defaultBox[i].min = PetscCMPLX(0.1, 0.1);
 80:     defaultBox[i].max = PetscCMPLX(1.5, 1.5);
 81: #endif
 82:   }
 83:   PetscCall(VecTaggerAbsoluteSetBox(tagger, defaultBox));
 84:   PetscCall(PetscFree(defaultBox));
 85:   PetscCall(VecTaggerSetFromOptions(tagger));
 86:   PetscCall(VecTaggerSetUp(tagger));
 87:   PetscCall(PetscObjectViewFromOptions((PetscObject)tagger, NULL, "-vec_tagger_view"));
 88:   PetscCall(VecTaggerGetBlockSize(tagger, &bs));

 90:   PetscCall(VecTaggerComputeBoxes(tagger, vec, &nint, &boxes, &listed));
 91:   if (listed) {
 92:     PetscViewer viewer = NULL;

 94:     PetscCall(PetscOptionsGetViewer(comm, NULL, NULL, "-vec_tagger_boxes_view", &viewer, NULL, NULL));
 95:     if (viewer) {
 96:       PetscBool iascii;

 98:       PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
 99:       if (iascii) {
100:         PetscCall(PetscViewerASCIIPrintf(viewer, "Num boxes: %" PetscInt_FMT "\n", nint));
101:         PetscCall(PetscViewerASCIIPushTab(viewer));
102:         for (i = 0, k = 0; i < nint; i++) {
103:           PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT ": ", i));
104:           for (j = 0; j < bs; j++, k++) {
105:             if (j) PetscCall(PetscViewerASCIIPrintf(viewer, " x "));
106: #if !defined(PETSC_USE_COMPLEX)
107:             PetscCall(PetscViewerASCIIPrintf(viewer, "[%g,%g]", (double)boxes[k].min, (double)boxes[k].max));
108: #else
109:             PetscCall(PetscViewerASCIIPrintf(viewer, "[%g+%gi,%g+%gi]", (double)PetscRealPart(boxes[k].min), (double)PetscImaginaryPart(boxes[k].min), (double)PetscRealPart(boxes[k].max), (double)PetscImaginaryPart(boxes[k].max)));
110: #endif
111:           }
112:           PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
113:         }
114:         PetscCall(PetscViewerASCIIPopTab(viewer));
115:       }
116:     }
117:     PetscCall(PetscViewerDestroy(&viewer));
118:     PetscCall(PetscFree(boxes));
119:   }

121:   PetscCall(VecTaggerComputeIS(tagger, vec, &is, &listed));
122:   PetscCall(ISGetBlockGlobalIS(is, vec, bs, &isBlockGlobal));
123:   PetscCall(PetscObjectSetName((PetscObject)isBlockGlobal, "Tagged IS (block global)"));
124:   PetscCall(ISViewFromOptions(isBlockGlobal, NULL, "-tagged_is_view"));

126:   PetscCall(VecGetOwnershipRange(vec, &localStart, &localEnd));
127:   PetscCall(ISComplement(isBlockGlobal, localStart, localEnd, &isComp));
128:   PetscCall(PetscObjectSetName((PetscObject)isComp, "Untagged IS (global)"));
129:   PetscCall(ISViewFromOptions(isComp, NULL, "-untagged_is_view"));

131:   PetscCall(ISGetLocalSize(isBlockGlobal, &ntagged));
132:   PetscCall(ISGetLocalSize(isComp, &nuntagged));

134:   PetscCall(VecCreate(comm, &tagged));
135:   PetscCall(PetscObjectSetName((PetscObject)tagged, "Tagged selection"));
136:   PetscCall(VecSetSizes(tagged, ntagged, PETSC_DETERMINE));
137:   PetscCall(VecSetUp(tagged));

139:   PetscCall(VecCreate(comm, &untagged));
140:   PetscCall(PetscObjectSetName((PetscObject)untagged, "Untagged selection"));
141:   PetscCall(VecSetSizes(untagged, nuntagged, PETSC_DETERMINE));
142:   PetscCall(VecSetUp(untagged));

144:   PetscCall(VecScatterCreate(vec, isBlockGlobal, tagged, NULL, &taggedScatter));
145:   PetscCall(VecScatterCreate(vec, isComp, untagged, NULL, &untaggedScatter));

147:   PetscCall(VecScatterBegin(taggedScatter, vec, tagged, INSERT_VALUES, SCATTER_FORWARD));
148:   PetscCall(VecScatterEnd(taggedScatter, vec, tagged, INSERT_VALUES, SCATTER_FORWARD));
149:   PetscCall(VecScatterBegin(untaggedScatter, vec, untagged, INSERT_VALUES, SCATTER_FORWARD));
150:   PetscCall(VecScatterEnd(untaggedScatter, vec, untagged, INSERT_VALUES, SCATTER_FORWARD));

152:   PetscCall(VecViewFromOptions(tagged, NULL, "-tagged_vec_view"));
153:   PetscCall(VecViewFromOptions(untagged, NULL, "-untagged_vec_view"));

155:   PetscCall(VecScatterDestroy(&untaggedScatter));
156:   PetscCall(VecScatterDestroy(&taggedScatter));

158:   PetscCall(VecDestroy(&untagged));
159:   PetscCall(VecDestroy(&tagged));
160:   PetscCall(ISDestroy(&isComp));
161:   PetscCall(ISDestroy(&isBlockGlobal));
162:   PetscCall(ISDestroy(&is));

164:   PetscCall(VecTaggerDestroy(&tagger));
165:   PetscCall(VecDestroy(&vec));
166:   PetscCall(PetscFinalize());
167:   return 0;
168: }

170: /*TEST

172:   test:
173:     suffix: 0
174:     requires: !complex
175:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view

177:   test:
178:     suffix: 1
179:     requires: !complex
180:     nsize: 3
181:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view

183:   test:
184:     suffix: 2
185:     requires: !complex
186:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -bs 2

188:   test:
189:     suffix: 3
190:     requires: !complex
191:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -vec_tagger_block_size 2 -vec_tagger_box 0.1,1.5,0.1,1.5

193:   test:
194:     suffix: 4
195:     requires: !complex
196:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -vec_tagger_invert

198:   test:
199:     suffix: 5
200:     requires: !complex
201:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -vec_tagger_type relative -vec_tagger_box 0.25,0.75

203:   test:
204:     suffix: 6
205:     requires: !complex
206:     nsize: 3
207:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -vec_tagger_type cdf -vec_tagger_box 0.25,0.75

209:   test:
210:     suffix: 7
211:     requires: !complex
212:     nsize: 3
213:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -vec_tagger_type cdf -vec_tagger_box 0.25,0.75 -vec_tagger_cdf_method iterative -vec_tagger_cdf_max_it 10

215:   test:
216:     suffix: 8
217:     requires: !complex
218:     nsize: 3
219:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -vec_tagger_type or -vec_tagger_num_subs 2 -sub_0_vec_tagger_type absolute -sub_0_vec_tagger_box 0.0,0.25 -sub_1_vec_tagger_type relative -sub_1_vec_tagger_box 0.75,inf
220:     filter: sed -e s~Inf~inf~g

222:   test:
223:     suffix: 9
224:     requires: !complex
225:     nsize: 3
226:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -vec_tagger_type and -vec_tagger_num_subs 2 -sub_0_vec_tagger_type absolute -sub_0_vec_tagger_box -inf,0.5 -sub_1_vec_tagger_type relative -sub_1_vec_tagger_box 0.25,0.75
227:     filter: sed -e s~Inf~inf~g

229:   test:
230:     suffix: 10
231:     requires: complex
232:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view

234:   test:
235:     suffix: 11
236:     requires: complex
237:     nsize: 3
238:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view

240:   test:
241:     suffix: 12
242:     requires: complex
243:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -bs 2

245:   test:
246:     suffix: 13
247:     requires: complex
248:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -vec_tagger_block_size 2 -vec_tagger_box 0.1+0.1i,1.5+1.5i,0.1+0.1i,1.5+1.5i

250:   test:
251:     suffix: 14
252:     requires: complex
253:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -vec_tagger_invert

255:   test:
256:     suffix: 15
257:     requires: complex
258:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -vec_tagger_type relative -vec_tagger_box 0.25+0.25i,0.75+0.75i

260:   test:
261:     suffix: 16
262:     requires: complex
263:     nsize: 3
264:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -vec_tagger_type cdf -vec_tagger_box 0.25+0.25i,0.75+0.75i

266:   test:
267:     suffix: 17
268:     requires: complex
269:     nsize: 3
270:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -vec_tagger_type cdf -vec_tagger_box 0.25+0.25i,0.75+0.75i -vec_tagger_cdf_method iterative -vec_tagger_cdf_max_it 10

272:   test:
273:     suffix: 18
274:     requires: complex
275:     nsize: 3
276:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -vec_tagger_type or -vec_tagger_num_subs 2 -sub_0_vec_tagger_type absolute -sub_0_vec_tagger_box 0.0+0.0i,0.25+0.25i -sub_1_vec_tagger_type relative -sub_1_vec_tagger_box 0.75+0.75i,inf+infi
277:     filter: sed -e s~Inf~inf~g

279:   test:
280:     suffix: 19
281:     requires: complex
282:     nsize: 3
283:     args: -n 12 -vec_view -vec_tagger_view -vec_tagger_boxes_view -tagged_is_view -untagged_is_view -tagged_vec_view -untagged_vec_view -vec_tagger_type and -vec_tagger_num_subs 2 -sub_0_vec_tagger_type absolute -sub_0_vec_tagger_box -inf-infi,0.75+0.75i -sub_1_vec_tagger_type relative -sub_1_vec_tagger_box 0.25+0.25i,1.+1.i
284:     filter: sed -e s~Inf~inf~g

286: TEST*/