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: static int PartitionView_Triangular_1D_File(Partition p, PetscViewer viewer) {
9: Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
10: FILE *fd;
11: int numLocElements = p->numLocElements;
12: int numLocNodes = q->numLocNodes;
13: int i;
14: int ierr;
17: PetscViewerASCIIPrintf(viewer, "Partition Object:n");
18: PetscViewerASCIIPrintf(viewer, " Partition of triangular 1D grid with %d elements and %d nodesn", p->numElements, q->numNodes);
19: PetscViewerASCIIGetPointer(viewer, &fd);
20: PetscSynchronizedFPrintf(p->comm, fd, " Proc %d: %d elements %d nodesn", p->rank, numLocElements, numLocNodes);
21: PetscSynchronizedFlush(p->comm);
22: if (p->ordering != PETSC_NULL) {
23: PetscViewerASCIIPrintf(viewer, " Global element renumbering:n");
24: AOView(p->ordering, viewer);
25: }
26: if (q->nodeOrdering != PETSC_NULL) {
27: PetscViewerASCIIPrintf(viewer, " Global node renumbering:n");
28: AOView(q->nodeOrdering, viewer);
29: }
30: PetscSynchronizedFPrintf(p->comm, fd, " %d ghost elements on proc %dn", p->numOverlapElements - numLocElements, p->rank);
31: for(i = 0; i < p->numOverlapElements - numLocElements; i++)
32: PetscSynchronizedFPrintf(p->comm, fd, " %d %d %dn", i, p->ghostElements[i], p->ghostElementProcs[i]);
33: PetscSynchronizedFlush(p->comm);
34: PetscSynchronizedFPrintf(p->comm, fd, " %d ghost nodes on proc %dn", q->numOverlapNodes - numLocNodes, p->rank);
35: for(i = 0; i < q->numOverlapNodes - numLocNodes; i++)
36: PetscSynchronizedFPrintf(p->comm, fd, " %d %dn", i, q->ghostNodes[i]);
37: PetscSynchronizedFlush(p->comm);
39: return(0);
40: }
42: static int PartitionView_Triangular_1D(Partition p, PetscViewer viewer) {
43: PetscTruth isascii;
44: int ierr;
47: PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_ASCII, &isascii);
48: if (isascii == PETSC_TRUE) {
49: PartitionView_Triangular_1D_File(p, viewer);
50: }
51: return(0);
52: }
54: static int PartitionDestroy_Triangular_1D(Partition p) {
55: Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
56: int ierr;
59: PetscFree(p->firstElement);
60: if (p->ordering != PETSC_NULL)
61: {AODestroy(p->ordering); }
62: if (p->ghostElements != PETSC_NULL)
63: {PetscFree(p->ghostElements); }
64: if (p->ghostElementProcs != PETSC_NULL)
65: {PetscFree(p->ghostElementProcs); }
66: PetscFree(q->firstNode);
67: PetscFree(q->firstBdNode);
68: if (q->nodeOrdering != PETSC_NULL)
69: {AODestroy(q->nodeOrdering); }
70: if (q->ghostNodes != PETSC_NULL)
71: {PetscFree(q->ghostNodes); }
72: if (q->ghostNodeProcs != PETSC_NULL)
73: {PetscFree(q->ghostNodeProcs); }
74: if (q->ghostBdNodes != PETSC_NULL)
75: {PetscFree(q->ghostBdNodes); }
76: PetscFree(q);
78: return(0);
79: }
81: static int PartitionGetTotalNodes_Triangular_1D(Partition p, int *size) {
82: Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
85: *size = q->numNodes;
86: return(0);
87: }
89: static int PartitionGetStartNode_Triangular_1D(Partition p, int *node) {
90: Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
93: *node = q->firstNode[p->rank];
94: return(0);
95: }
97: static int PartitionGetEndNode_Triangular_1D(Partition p, int *node) {
98: Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
101: *node = q->firstNode[p->rank+1];
102: return(0);
103: }
105: static int PartitionGetNumNodes_Triangular_1D(Partition p, int *size) {
106: Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
109: *size = q->numLocNodes;
110: return(0);
111: }
113: static int PartitionGetNumOverlapNodes_Triangular_1D(Partition p, int *size) {
114: Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
117: *size = q->numOverlapNodes;
118: return(0);
119: }
121: static int PartitionGhostNodeIndex_Private(Partition p, int node, int *gNode) {
122: Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
123: int low, high, mid;
126: /* Use bisection since the array is assumed to be sorted */
127: low = 0;
128: high = q->numOverlapNodes - (q->firstNode[p->rank+1] - q->firstNode[p->rank]) - 1;
129: while (low <= high) {
130: mid = (low + high)/2;
131: if (node == q->ghostNodes[mid]) {
132: *gNode = mid;
133: return(0);
134: } else if (node < q->ghostNodes[mid]) {
135: high = mid - 1;
136: } else {
137: low = mid + 1;
138: }
139: }
140: *gNode = -1;
141: /* Flag for ghost node not present */
142: PetscFunctionReturn(1);
143: }
145: static int PartitionGlobalToLocalNodeIndex_Triangular_1D(Partition p, int node, int *locNode) {
146: Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
147: int numLocNodes = q->numLocNodes;
148: int gNode; /* Local ghost node number */
149: int ierr;
152: if (node < 0) {
153: *locNode = node;
154: return(0);
155: }
156: /* Check for ghost node */
157: if ((node < q->firstNode[p->rank]) || (node >= q->firstNode[p->rank+1])) {
158: /* Search for canonical number */
159: PartitionGhostNodeIndex_Private(p, node, &gNode);
160: *locNode = numLocNodes + gNode;
161: } else {
162: *locNode = node - q->firstNode[p->rank];
163: }
164: return(0);
165: }
167: static int PartitionLocalToGlobalNodeIndex_Triangular_1D(Partition p, int locNode, int *node) {
168: Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
169: int numLocNodes = q->numLocNodes;
172: if (locNode < 0) {
173: *node = locNode;
174: return(0);
175: }
176: /* Check for ghost node */
177: if (locNode >= numLocNodes) {
178: *node = q->ghostNodes[locNode - numLocNodes];
179: } else {
180: *node = locNode + q->firstNode[p->rank];
181: }
182: return(0);
183: }
185: static int PartitionGlobalToGhostNodeIndex_Triangular_1D(Partition p, int node, int *ghostNode, int *ghostProc) {
186: Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
187: int ierr;
190: if (node < 0) {
191: *ghostNode = node;
192: *ghostProc = -1;
193: return(0);
194: }
195: /* Check for ghost node */
196: if ((node < q->firstNode[p->rank]) || (node >= q->firstNode[p->rank+1])) {
197: PartitionGhostNodeIndex_Private(p, node, ghostNode);
198: *ghostProc = q->ghostNodeProcs[*ghostNode];
199: } else {
200: SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE, "Global node %d is not a ghost node", node);
201: }
202: return(0);
203: }
205: static int PartitionGhostToGlobalNodeIndex_Triangular_1D(Partition p, int ghostNode, int *node, int *ghostProc) {
206: Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
207: int numGhostNodes = q->numOverlapNodes - q->numLocNodes;
210: if (ghostNode < 0) {
211: *node = ghostNode;
212: *ghostProc = -1;
213: return(0);
214: }
215: /* Check for ghost node */
216: if (ghostNode < numGhostNodes) {
217: *node = q->ghostNodes[ghostNode];
218: *ghostProc = q->ghostNodeProcs[ghostNode];
219: } else {
220: SETERRQ1(PETSC_ERR_ARG_OUTOFRANGE, "Ghost node %d does not exist", ghostNode);
221: }
222: return(0);
223: }
225: static int PartitionGetNodeOrdering_Triangular_1D(Partition p, AO *order) {
226: Partition_Triangular_1D *q = (Partition_Triangular_1D *) p->data;
229: *order = q->nodeOrdering;
230: return(0);
231: }
233: static struct _PartitionOps POps = {/* Generic Operations */
234: PETSC_NULL/* PartitionSetup */,
235: PETSC_NULL/* PartitionSetFromOptions */,
236: PartitionView_Triangular_1D,
237: PETSC_NULL/* PartitionCopy */,
238: PETSC_NULL/* PartitionDuplicate */,
239: PartitionDestroy_Triangular_1D,
240: /* Partition-Specific Operations */
241: PETSC_NULL/* PartitionGhostNodeExchange */,
242: /* Node Query Functions */
243: PartitionGetTotalNodes_Triangular_1D,
244: PartitionGetStartNode_Triangular_1D,
245: PartitionGetEndNode_Triangular_1D,
246: PartitionGetNumNodes_Triangular_1D,
247: PartitionGetNumOverlapNodes_Triangular_1D,
248: PartitionGlobalToLocalNodeIndex_Triangular_1D,
249: PartitionLocalToGlobalNodeIndex_Triangular_1D,
250: PartitionGlobalToGhostNodeIndex_Triangular_1D,
251: PartitionGhostToGlobalNodeIndex_Triangular_1D,
252: PartitionGetNodeOrdering_Triangular_1D,
253: /* Face Query Functions */
254: PETSC_NULL/* PartitionGetTotalFaces */,
255: PETSC_NULL/* PartitionGetStartFace */,
256: PETSC_NULL/* PartitionGetEndFace */,
257: PETSC_NULL/* PartitionGetNumFaces */,
258: PETSC_NULL/* PartitionGetNumOverlapFaces */,
259: PETSC_NULL/* PartitionGlobalToLocalFaceIndex */,
260: PETSC_NULL/* PartitionLocalToGlobalFaceIndex */,
261: PETSC_NULL/* PartitionGetFaceOrdering */,
262: /* Edge Query Functions */
263: PartitionGetTotalElements,
264: PartitionGetStartElement,
265: PartitionGetEndElement,
266: PartitionGetNumElements,
267: PartitionGetNumOverlapElements,
268: PartitionGlobalToLocalElementIndex,
269: PartitionLocalToGlobalElementIndex,
270: PartitionGetElementOrdering};
272: EXTERN_C_BEGIN
273: int PartitionCreate_Triangular_1D(Partition p) {
274: Partition_Triangular_1D *q;
275: Mesh mesh = p->mesh;
276: int numProcs, rank, rem;
277: int proc;
278: int ierr;
281: PetscNew(Partition_Triangular_1D, &q);
282: PetscLogObjectMemory(p, sizeof(Partition_Triangular_1D));
283: PetscMemcpy(p->ops, &POps, sizeof(struct _PartitionOps));
284: p->data = (void *) q;
285: PetscStrallocpy(PARTITION_SER_TRIANGULAR_1D_BINARY, &p->serialize_name);
286: PetscLogObjectParent(mesh, p);
288: /* Initialize structure */
289: MPI_Comm_size(p->comm, &numProcs);
290: MPI_Comm_rank(p->comm, &rank);
291: p->numProcs = numProcs;
292: p->rank = rank;
293: p->isElementPartitioned = PETSC_FALSE;
294: p->ordering = PETSC_NULL;
295: p->ghostElements = PETSC_NULL;
296: p->ghostElementProcs = PETSC_NULL;
297: q->isNodePartitioned = PETSC_FALSE;
298: q->nodeOrdering = PETSC_NULL;
299: q->ghostNodes = PETSC_NULL;
300: q->ghostNodeProcs = PETSC_NULL;
301: q->ghostBdNodes = PETSC_NULL;
302: PetscMalloc((numProcs+1) * sizeof(int), &p->firstElement);
303: PetscMalloc((numProcs+1) * sizeof(int), &q->firstNode);
304: PetscMalloc((numProcs+1) * sizeof(int), &q->firstBdNode);
305: PetscLogObjectMemory(p, (numProcs+1)*4 * sizeof(int));
306: PetscMemzero(p->firstElement, (numProcs+1) * sizeof(int));
307: PetscMemzero(q->firstNode, (numProcs+1) * sizeof(int));
308: PetscMemzero(q->firstBdNode, (numProcs+1) * sizeof(int));
310: /* Setup crude preliminary partition */
311: for(proc = 0; proc < numProcs; proc++) {
312: rem = (mesh->numEdges%numProcs);
313: p->firstElement[proc] = (mesh->numEdges/numProcs)*proc + PetscMin(rem, proc);
314: rem = (mesh->numNodes%numProcs);
315: q->firstNode[proc] = (mesh->numNodes/numProcs)*proc + PetscMin(rem, proc);
316: rem = (mesh->numBdNodes%numProcs);
317: q->firstBdNode[proc] = (mesh->numBdNodes/numProcs)*proc + PetscMin(rem, proc);
318: }
319: p->firstElement[numProcs] = mesh->numEdges;
320: q->firstNode[numProcs] = mesh->numNodes;
321: q->firstBdNode[numProcs] = mesh->numBdNodes;
323: p->numLocElements = p->firstElement[rank+1] - p->firstElement[rank];
324: p->numElements = p->firstElement[numProcs];
325: p->numOverlapElements = p->numLocElements;
326: q->numLocNodes = q->firstNode[rank+1] - q->firstNode[rank];
327: q->numNodes = q->firstNode[numProcs];
328: q->numOverlapNodes = q->numLocNodes;
329: q->numLocBdNodes = q->firstBdNode[rank+1] - q->firstBdNode[rank];
330: q->numBdNodes = q->firstBdNode[numProcs];
331: q->numOverlapBdNodes = q->numLocBdNodes;
332: return(0);
333: }
334: EXTERN_C_END