Actual source code: part1d.c

  1: #ifdef PETSC_RCS_HEADER
  2: static char vcid[] = "$Id: part1d.c,v 1.38 2000/10/17 13:48:55 knepley Exp $";
  3: #endif

  5: /* Implements 1d unstructured grids */
  6: #include "src/mesh/impls/triangular/1d/1dimpl.h"         /*I "partition.h" I*/

  8: #undef  __FUNCT__
 10: static int PartitionView_Triangular_1D_File(Partition p, PetscViewer viewer) {
 11:   Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
 12:   FILE                    *fd;
 13:   int                      numLocElements = p->numLocElements;
 14:   int                      numLocNodes    = q->numLocNodes;
 15:   int                      i;
 16:   int                      ierr;

 19:   PetscViewerASCIIPrintf(viewer, "Partition Object:n");
 20:   PetscViewerASCIIPrintf(viewer, "  Partition of triangular 1D grid with %d elements and %d nodesn", p->numElements, q->numNodes);
 21:   PetscViewerASCIIGetPointer(viewer, &fd);
 22:   PetscSynchronizedFPrintf(p->comm, fd, "    Proc %d: %d elements %d nodesn", p->rank, numLocElements, numLocNodes);
 23:   PetscSynchronizedFlush(p->comm);
 24:   if (p->ordering != PETSC_NULL) {
 25:     PetscViewerASCIIPrintf(viewer, "  Global element renumbering:n");
 26:     AOView(p->ordering, viewer);
 27:   }
 28:   if (q->nodeOrdering != PETSC_NULL) {
 29:     PetscViewerASCIIPrintf(viewer, "  Global node renumbering:n");
 30:     AOView(q->nodeOrdering, viewer);
 31:   }
 32:   PetscSynchronizedFPrintf(p->comm, fd, "  %d ghost elements on proc %dn", p->numOverlapElements - numLocElements, p->rank);
 33:   for(i = 0; i < p->numOverlapElements - numLocElements; i++)
 34:     PetscSynchronizedFPrintf(p->comm, fd, "  %d %d %dn", i, p->ghostElements[i], p->ghostElementProcs[i]);
 35:   PetscSynchronizedFlush(p->comm);
 36:   PetscSynchronizedFPrintf(p->comm, fd, "  %d ghost nodes on proc %dn", q->numOverlapNodes - numLocNodes, p->rank);
 37:   for(i = 0; i < q->numOverlapNodes - numLocNodes; i++)
 38:     PetscSynchronizedFPrintf(p->comm, fd, "  %d %dn", i, q->ghostNodes[i]);
 39:   PetscSynchronizedFlush(p->comm);

 41:   return(0);
 42: }

 44: #undef  __FUNCT__
 46: static int PartitionView_Triangular_1D(Partition p, PetscViewer viewer) {
 47:   PetscTruth isascii;
 48:   int        ierr;

 51:   PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_ASCII, &isascii);
 52:   if (isascii == PETSC_TRUE) {
 53:     PartitionView_Triangular_1D_File(p, viewer);
 54:   }
 55:   return(0);
 56: }

 58: #undef  __FUNCT__
 60: static int PartitionDestroy_Triangular_1D(Partition p) {
 61:   Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
 62:   int                      ierr;

 65:   PetscFree(p->firstElement);
 66:   if (p->ordering != PETSC_NULL)
 67:     {AODestroy(p->ordering);                                                                      }
 68:   if (p->ghostElements != PETSC_NULL)
 69:     {PetscFree(p->ghostElements);                                                                 }
 70:   if (p->ghostElementProcs != PETSC_NULL)
 71:     {PetscFree(p->ghostElementProcs);                                                             }
 72:   PetscFree(q->firstNode);
 73:   PetscFree(q->firstBdNode);
 74:   if (q->nodeOrdering != PETSC_NULL)
 75:     {AODestroy(q->nodeOrdering);                                                                  }
 76:   if (q->ghostNodes != PETSC_NULL)
 77:     {PetscFree(q->ghostNodes);                                                                    }
 78:   if (q->ghostNodeProcs != PETSC_NULL)
 79:     {PetscFree(q->ghostNodeProcs);                                                                }
 80:   if (q->ghostBdNodes != PETSC_NULL)
 81:     {PetscFree(q->ghostBdNodes);                                                                  }
 82:   PetscFree(q);

 84:   return(0);
 85: }

 87: #undef  __FUNCT__
 89: static int PartitionGetTotalNodes_Triangular_1D(Partition p, int *size) {
 90:   Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;

 93:   *size = q->numNodes;
 94:   return(0);
 95: }

 97: #undef  __FUNCT__
 99: static int PartitionGetStartNode_Triangular_1D(Partition p, int *node) {
100:   Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;

103:   *node = q->firstNode[p->rank];
104:   return(0);
105: }

107: #undef  __FUNCT__
109: static int PartitionGetEndNode_Triangular_1D(Partition p, int *node) {
110:   Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;

113:   *node = q->firstNode[p->rank+1];
114:   return(0);
115: }

117: #undef  __FUNCT__
119: static int PartitionGetNumNodes_Triangular_1D(Partition p, int *size) {
120:   Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;

123:   *size = q->numLocNodes;
124:   return(0);
125: }

127: #undef  __FUNCT__
129: static int PartitionGetNumOverlapNodes_Triangular_1D(Partition p, int *size) {
130:   Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;

133:   *size = q->numOverlapNodes;
134:   return(0);
135: }

137: #undef  __FUNCT__
139: static int PartitionGhostNodeIndex_Private(Partition p, int node, int *gNode) {
140:   Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
141:   int                      low, high, mid;

144:   /* Use bisection since the array is assumed to be sorted */
145:   low  = 0;
146:   high = q->numOverlapNodes - (q->firstNode[p->rank+1] - q->firstNode[p->rank]) - 1;
147:   while (low <= high) {
148:     mid = (low + high)/2;
149:     if (node == q->ghostNodes[mid]) {
150:       *gNode = mid;
151:       return(0);
152:     } else if (node < q->ghostNodes[mid]) {
153:       high = mid - 1;
154:     } else {
155:       low  = mid + 1;
156:     }
157:   }
158:   *gNode = -1;
159:   /* Flag for ghost node not present */
160:   PetscFunctionReturn(1);
161: }

163: #undef  __FUNCT__
165: static int PartitionGlobalToLocalNodeIndex_Triangular_1D(Partition p, int node, int *locNode) {
166:   Partition_Triangular_1D *q           = (Partition_Triangular_1D *) p->data;
167:   int                      numLocNodes = q->numLocNodes;
168:   int                      gNode; /* Local ghost node number */
169:   int                      ierr;

172:   if (node < 0) {
173:     *locNode = node;
174:     return(0);
175:   }
176:   /* Check for ghost node */
177:   if ((node < q->firstNode[p->rank]) || (node >= q->firstNode[p->rank+1])) {
178:     /* Search for canonical number */
179:     PartitionGhostNodeIndex_Private(p, node, &gNode);
180:     *locNode = numLocNodes + gNode;
181:   } else {
182:     *locNode = node - q->firstNode[p->rank];
183:   }
184:   return(0);
185: }

187: #undef  __FUNCT__
189: static int PartitionLocalToGlobalNodeIndex_Triangular_1D(Partition p, int locNode, int *node) {
190:   Partition_Triangular_1D *q           = (Partition_Triangular_1D *) p->data;
191:   int                      numLocNodes = q->numLocNodes;

194:   if (locNode < 0) {
195:     *node = locNode;
196:     return(0);
197:   }
198:   /* Check for ghost node */
199:   if (locNode >= numLocNodes) {
200:     *node = q->ghostNodes[locNode - numLocNodes];
201:   } else {
202:     *node = locNode + q->firstNode[p->rank];
203:   }
204:   return(0);
205: }

207: #undef  __FUNCT__
209: static int PartitionGlobalToGhostNodeIndex_Triangular_1D(Partition p, int node, int *ghostNode, int *ghostProc) {
210:   Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
211:   int                      ierr;

214:   if (node < 0) {
215:     *ghostNode = node;
216:     *ghostProc = -1;
217:     return(0);
218:   }
219:   /* Check for ghost node */
220:   if ((node < q->firstNode[p->rank]) || (node >= q->firstNode[p->rank+1])) {
221:     PartitionGhostNodeIndex_Private(p, node, ghostNode);
222:     *ghostProc = q->ghostNodeProcs[*ghostNode];
223:   } else {
224:     SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE, "Global node %d is not a ghost node", node);
225:   }
226:   return(0);
227: }

229: #undef  __FUNCT__
231: static int PartitionGhostToGlobalNodeIndex_Triangular_1D(Partition p, int ghostNode, int *node, int *ghostProc) {
232:   Partition_Triangular_1D *q             = (Partition_Triangular_1D *) p->data;
233:   int                      numGhostNodes = q->numOverlapNodes - q->numLocNodes;

236:   if (ghostNode < 0) {
237:     *node      = ghostNode;
238:     *ghostProc = -1;
239:     return(0);
240:   }
241:   /* Check for ghost node */
242:   if (ghostNode < numGhostNodes) {
243:     *node      = q->ghostNodes[ghostNode];
244:     *ghostProc = q->ghostNodeProcs[ghostNode];
245:   } else {
246:     SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE, "Ghost node %d does not exist", ghostNode);
247:   }
248:   return(0);
249: }

251: #undef  __FUNCT__
253: static int PartitionGetNodeOrdering_Triangular_1D(Partition p, AO *order) {
254:   Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;

257:   *order = q->nodeOrdering;
258:   return(0);
259: }

261: static struct _PartitionOps POps = {/* Generic Operations */
262:                                     PETSC_NULL/* PartitionSetup */,
263:                                     PETSC_NULL/* PartitionSetFromOptions */,
264:                                     PartitionView_Triangular_1D,
265:                                     PETSC_NULL/* PartitionCopy */,
266:                                     PETSC_NULL/* PartitionDuplicate */,
267:                                     PartitionDestroy_Triangular_1D,
268:                                     /* Partition-Specific Operations */
269:                                     PETSC_NULL/* PartitionGhostNodeExchange */,
270:                                     /* Node Query Functions */
271:                                     PartitionGetTotalNodes_Triangular_1D,
272:                                     PartitionGetStartNode_Triangular_1D,
273:                                     PartitionGetEndNode_Triangular_1D,
274:                                     PartitionGetNumNodes_Triangular_1D,
275:                                     PartitionGetNumOverlapNodes_Triangular_1D,
276:                                     PartitionGlobalToLocalNodeIndex_Triangular_1D,
277:                                     PartitionLocalToGlobalNodeIndex_Triangular_1D,
278:                                     PartitionGlobalToGhostNodeIndex_Triangular_1D,
279:                                     PartitionGhostToGlobalNodeIndex_Triangular_1D,
280:                                     PartitionGetNodeOrdering_Triangular_1D,
281:                                     /* Face Query Functions */
282:                                     PETSC_NULL/* PartitionGetTotalFaces */,
283:                                     PETSC_NULL/* PartitionGetStartFace */,
284:                                     PETSC_NULL/* PartitionGetEndFace */,
285:                                     PETSC_NULL/* PartitionGetNumFaces */,
286:                                     PETSC_NULL/* PartitionGetNumOverlapFaces */,
287:                                     PETSC_NULL/* PartitionGlobalToLocalFaceIndex */,
288:                                     PETSC_NULL/* PartitionLocalToGlobalFaceIndex */,
289:                                     PETSC_NULL/* PartitionGetFaceOrdering */,
290:                                     /* Edge Query Functions */
291:                                     PartitionGetTotalElements,
292:                                     PartitionGetStartElement,
293:                                     PartitionGetEndElement,
294:                                     PartitionGetNumElements,
295:                                     PartitionGetNumOverlapElements,
296:                                     PartitionGlobalToLocalElementIndex,
297:                                     PartitionLocalToGlobalElementIndex,
298:                                     PartitionGetElementOrdering};

300: EXTERN_C_BEGIN
301: #undef  __FUNCT__
303: int PartitionCreate_Triangular_1D(Partition p) {
304:   Partition_Triangular_1D *q;
305:   Mesh                     mesh = p->mesh;
306:   int                      numProcs, rank, rem;
307:   int                      proc;
308:   int                      ierr;

311:   PetscNew(Partition_Triangular_1D, &q);
312:   PetscLogObjectMemory(p, sizeof(Partition_Triangular_1D));
313:   PetscMemcpy(p->ops, &POps, sizeof(struct _PartitionOps));
314:   p->data = (void *) q;
315:   PetscStrallocpy(PARTITION_SER_TRIANGULAR_1D_BINARY, &p->serialize_name);
316:   PetscLogObjectParent(mesh, p);

318:   /* Initialize structure */
319:   MPI_Comm_size(p->comm, &numProcs);
320:   MPI_Comm_rank(p->comm, &rank);
321:   p->numProcs             = numProcs;
322:   p->rank                 = rank;
323:   p->isElementPartitioned = PETSC_FALSE;
324:   p->ordering             = PETSC_NULL;
325:   p->ghostElements        = PETSC_NULL;
326:   p->ghostElementProcs    = PETSC_NULL;
327:   q->isNodePartitioned    = PETSC_FALSE;
328:   q->nodeOrdering         = PETSC_NULL;
329:   q->ghostNodes           = PETSC_NULL;
330:   q->ghostNodeProcs       = PETSC_NULL;
331:   q->ghostBdNodes         = PETSC_NULL;
332:   PetscMalloc((numProcs+1) * sizeof(int), &p->firstElement);
333:   PetscMalloc((numProcs+1) * sizeof(int), &q->firstNode);
334:   PetscMalloc((numProcs+1) * sizeof(int), &q->firstBdNode);
335:   PetscLogObjectMemory(p, (numProcs+1)*4 * sizeof(int));
336:   PetscMemzero(p->firstElement, (numProcs+1) * sizeof(int));
337:   PetscMemzero(q->firstNode,    (numProcs+1) * sizeof(int));
338:   PetscMemzero(q->firstBdNode,  (numProcs+1) * sizeof(int));

340:   /* Setup crude preliminary partition */
341:   for(proc = 0; proc < numProcs; proc++) {
342:     rem                   = (mesh->numEdges%numProcs);
343:     p->firstElement[proc] = (mesh->numEdges/numProcs)*proc + PetscMin(rem, proc);
344:     rem                   = (mesh->numNodes%numProcs);
345:     q->firstNode[proc]    = (mesh->numNodes/numProcs)*proc + PetscMin(rem, proc);
346:     rem                   = (mesh->numBdNodes%numProcs);
347:     q->firstBdNode[proc]  = (mesh->numBdNodes/numProcs)*proc + PetscMin(rem, proc);
348:   }
349:   p->firstElement[numProcs] = mesh->numEdges;
350:   q->firstNode[numProcs]    = mesh->numNodes;
351:   q->firstBdNode[numProcs]  = mesh->numBdNodes;

353:   p->numLocElements         = p->firstElement[rank+1] - p->firstElement[rank];
354:   p->numElements            = p->firstElement[numProcs];
355:   p->numOverlapElements     = p->numLocElements;
356:   q->numLocNodes            = q->firstNode[rank+1] - q->firstNode[rank];
357:   q->numNodes               = q->firstNode[numProcs];
358:   q->numOverlapNodes        = q->numLocNodes;
359:   q->numLocBdNodes          = q->firstBdNode[rank+1] - q->firstBdNode[rank];
360:   q->numBdNodes             = q->firstBdNode[numProcs];
361:   q->numOverlapBdNodes      = q->numLocBdNodes;
362:   return(0);
363: }
364: EXTERN_C_END