Actual source code: varorder.c
1: #ifdef PETSC_RCS_HEADER
2: static char vcid[] = "$Id: varorder.c,v 1.9 2000/01/10 03:54:18 knepley Exp $";
3: #endif
5: #include src/grid/gridimpl.h
7: /* Loggging support */
8: int VAR_ORDER_COOKIE;
10: /*------------------------------------------------- Variable Ordering ------------------------------------------------*/
11: /*@C
12: VarOrderingCreate - This function creates a global variable ordering.
14: Collective on Grid
16: Input Parameter:
17: . grid - The underlying grid for the ordering
19: Output Parameter:
20: . order - The ordering
22: Level: beginner
24: .keywords: variable ordering, create
25: .seealso: VarOrderingCreateConstrained(), VarOrderingDestroy()
26: @*/
27: int VarOrderingCreate(Grid grid, VarOrdering *order)
28: {
34: GridSetUp(grid);
35: (*grid->ops->gridcreatevarordering)(grid, grid->cm, order);
36: return(0);
37: }
39: /*@C
40: VarOrderingCreateConstrained - This function creates a global variable ordering on the constraint space.
42: Collective on Grid
44: Input Parameter:
45: . grid - The underlying grid for the ordering
47: Output Parameter:
48: . order - The ordering
50: Level: beginner
52: .keywords: variable ordering, create, constraint
53: .seealso: VarOrderingCreate(), VarOrderingDestroy()
54: @*/
55: int VarOrderingCreateConstrained(Grid grid, VarOrdering *order)
56: {
57: FieldClassMap cm;
58: int ierr;
63: GridSetUp(grid);
64: GridGetClassMap(grid, &cm);
65: (*grid->ops->gridcreatevarordering)(grid, cm, order);
66: return(0);
67: }
69: /*@C
70: VarOrderingCreateGeneral - This function creates a global variable ordering on the specified fields.
72: Collective on Grid
74: Input Parameters:
75: + grid - The underlying grid for the ordering
76: . numFields - The number of fields
77: - fields - The fields
79: Output Parameter:
80: . order - The ordering
82: Level: beginner
84: .keywords: variable ordering, create, constraint
85: .seealso: VarOrderingCreate(), VarOrderingDestroy()
86: @*/
87: int VarOrderingCreateGeneral(Grid grid, int numFields, int *fields, VarOrdering *order)
88: {
89: FieldClassMap map, tempMap;
90: PetscTruth isConstrained, reduceSystem;
91: int numTotalFields, f;
92: int ierr;
97: GridGetNumFields(grid, &numTotalFields);
98: if (numFields > numTotalFields) {
99: SETERRQ2(PETSC_ERR_ARG_WRONG, "Invalid number %d of fields, grid only has %d", numFields, numTotalFields);
100: }
101: for(f = 0; f < numFields; f++) {
102: if ((fields[f] < 0) || (fields[f] >= numTotalFields)) {
103: SETERRQ1(PETSC_ERR_ARG_WRONG, "Invalid field number %d", fields[f]);
104: }
105: }
106: GridSetUp(grid);
107: FieldClassMapCreateTriangular2D(grid, numFields, fields, &tempMap);
108: GridIsConstrained(grid, &isConstrained);
109: GridGetReduceSystem(grid, &reduceSystem);
110: if (reduceSystem == PETSC_TRUE) {
111: FieldClassMapConstrain(tempMap, grid, reduceSystem, isConstrained, &map);
112: FieldClassMapDestroy(tempMap);
113: } else {
114: map = tempMap;
115: }
116: (*grid->ops->gridcreatevarordering)(grid, map, order);
117: FieldClassMapDestroy(map);
118: return(0);
119: }
121: /*@C
122: VarOrderingCreateSubset - This function creates a new global variable ordering from an existing
123: one on a subset of the fields.
125: Collective on Grid
127: Input Parameters:
128: + order - The original ordering
129: . numFields - The number of fields in the new ordering
130: . fields - The fields in the new ordering
131: - contract - The flag for contracting the indices
133: Output Parameter:
134: . newOrder - The new ordering
136: Note:
137: If contract is PETSC_FALSE, this function returns a copy of the original
138: ordering, but with some fields no longer active.
140: Level: beginner
142: .keywords: variable ordering, create
143: .seealso: VarOrderingCreate(), VarOrderingDestroy()
144: @*/
145: int VarOrderingCreateSubset(VarOrdering order, int numFields, int *fields, PetscTruth contract, VarOrdering *newOrder)
146: {
147: FieldClassMap map, cm;
148: int f, g;
149: int ierr;
154: VarOrderingGetClassMap(order, &map);
155: /* Check for consistency */
156: for(f = 0; f < numFields; f++) {
157: for(g = 0; g < map->numFields; g++) {
158: if (map->fields[g] == fields[f])
159: break;
160: }
161: if (g == map->numFields) {
162: SETERRQ1(PETSC_ERR_ARG_INCOMP, "Fields not a subset: field %d not in the existing order", fields[f]);
163: }
164: }
165: if (contract == PETSC_TRUE) {
166: SETERRQ(PETSC_ERR_SUP, "Coming soon");
167: } else {
168: FieldClassMapCreateSubset(map, numFields, fields, &cm);
169: VarOrderingDuplicate(order, newOrder);
170: PetscObjectCompose((PetscObject) *newOrder, "ClassMap", (PetscObject) cm);
171: FieldClassMapDestroy(cm);
172: /* Free extra localStart entries */
173: for(g = 0; g < map->numFields; g++) {
174: for(f = 0; f < numFields; f++) {
175: if (map->fields[g] == fields[f])
176: break;
177: }
178: if (f == numFields) {
179: PetscFree((*newOrder)->localStart[map->fields[g]]);
180: }
181: }
182: }
183: return(0);
184: }
186: /*@C
187: VarOrderingConstrain - This function constrains a previous variable ordering using the constraints
188: already defined in the grid.
190: Collective on VarOrdering
192: Input Parameters:
193: + grid - The underlying grid for the ordering
194: - oldOrder - The original ordering
196: Output Parameter:
197: . order - The constrained ordering
199: Level: beginner
201: .keywords variable ordering, finite element
202: .seealso VarOrderingDestroy()
203: @*/
204: int VarOrderingConstrain(Grid grid, VarOrdering oldOrder, VarOrdering *order)
205: {
206: PetscTruth isConstrained;
207: FieldClassMap cm;
208: int numFields;
209: int field;
210: int ierr;
215: GridSetUp(grid);
216: GridIsConstrained(grid, &isConstrained);
217: if (isConstrained == PETSC_TRUE) {
218: GridGetNumFields(grid, &numFields);
219: GridGetClassMap(grid, &cm);
220: for(field = 0; field < numFields; field++) {
221: if (grid->fields[field].isConstrained == PETSC_TRUE) break;
222: }
223: if (field == numFields) SETERRQ(PETSC_ERR_ARG_WRONG, "No constrained fields");
224: (*grid->ops->gridvarorderingconstrain)(grid, field, grid->constraintCtx, cm, oldOrder, order);
225: } else {
226: VarOrderingDuplicate(oldOrder, order);
227: }
228: return(0);
229: }
231: /*@C
232: VarOrderingCreateReduce - This function creates a global variable ordering on the reduction space.
234: Collective on Grid
236: Input Parameter:
237: . grid - The underlying grid for the ordering
239: Output Parameter:
240: . order - The reduction ordering
242: Level: beginner
244: .keywords: variable ordering, create, reduction
245: .seealso: VarOrderingCreate(), VarOrderingDestroy()
246: @*/
247: int VarOrderingCreateReduce(Grid grid, VarOrdering *order)
248: {
254: GridSetUp(grid);
255: if (grid->reduceSystem == PETSC_TRUE) {
256: (*grid->ops->gridcreatevarordering)(grid, grid->reductionCM, order);
257: } else {
258: (*grid->ops->gridcreatevarordering)(grid, grid->cm, order);
259: }
260: return(0);
261: }
263: /*@C
264: VarOrderingCreateReduceGeneral - This function creates a global variable ordering on the reduction space.
266: Collective on Grid
268: Input Parameters:
269: + grid - The underlying grid for the ordering
270: . cm - The class map
271: - reductionCM - The class map for the reduction space
273: Output Parameter:
274: . order - The reduction ordering
276: Level: beginner
278: .keywords: variable ordering, create, reduction
279: .seealso: VarOrderingCreate(), VarOrderingDestroy()
280: @*/
281: int VarOrderingCreateReduceGeneral(Grid grid, FieldClassMap cm, FieldClassMap reductionCM, VarOrdering *order)
282: {
290: GridSetUp(grid);
291: if (grid->reduceSystem == PETSC_TRUE) {
292: (*grid->ops->gridcreatevarordering)(grid, reductionCM, order);
293: } else {
294: (*grid->ops->gridcreatevarordering)(grid, cm, order);
295: }
296: return(0);
297: }
299: /*@C
300: VarOrderingDestroy - This function destroys a global variable ordering.
302: Not collective
304: Input Parameter:
305: . order - The ordering
307: Level: beginner
309: .keywords: variable ordering
310: .seealso: VarOrderingCreate()
311: @*/
312: int VarOrderingDestroy(VarOrdering order)
313: {
314: FieldClassMap map;
315: int f;
316: int ierr;
320: if (--order->refct > 0)
321: return(0);
322: PetscObjectQuery((PetscObject) order, "ClassMap", (PetscObject *) &map);
323: PetscFree(order->firstVar);
324: PetscFree(order->offsets);
325: if (order->localOffsets) {
326: PetscFree(order->localOffsets);
327: }
328: for(f = 0; f < map->numFields; f++) {
329: PetscFree(order->localStart[map->fields[f]]);
330: }
331: PetscFree(order->localStart);
332: PetscLogObjectDestroy(order);
333: PetscHeaderDestroy(order);
334: return(0);
335: }
337: /*@
338: VarOrderingGetClassMap - This function returns the class map for the variable ordering.
340: Not collective
342: Input Parameter:
343: . order - The variable ordering
345: Output Parameter:
346: . map - The field class map
348: Level: intermediate
350: .keywords: variable ordering, class map
351: .seealso: VarOrderingGetLocalSize()
352: @*/
353: int VarOrderingGetClassMap(VarOrdering order, FieldClassMap *map)
354: {
360: PetscObjectQuery((PetscObject) order, "ClassMap", (PetscObject *) map);
362: return(0);
363: }
365: /*@
366: VarOrderingGetNumTotalFields - Gets the total number of fields in the associated grid.
368: Not collective
370: Input Parameter:
371: . order - The variable ordering
373: Output Parameter:
374: . numFields - The total number fields
376: Level: intermediate
378: .keywords: variable ordering, grid, field
379: .seealso: VarOrderingGetSize(), VarOrderingGetLocalSize()
380: @*/
381: int VarOrderingGetNumTotalFields(VarOrdering order, int *numFields)
382: {
386: *numFields = order->numTotalFields;
387: return(0);
388: }
390: /*@
391: VarOrderingGetSize - Gets the number of global variables in a variable ordering.
393: Not collective
395: Input Parameter:
396: . order - The variable ordering
398: Output Parameter:
399: . size - The number of global variables
401: Level: intermediate
403: .keywords: variable ordering, grid
404: .seealso: VarOrderingGetLocalSize()
405: @*/
406: int VarOrderingGetSize(VarOrdering order, int *size)
407: {
411: *size = order->numVars;
412: return(0);
413: }
415: /*@
416: VarOrderingGetLocalSize - Gets the number of local variables in a variable ordering.
418: Not collective
420: Input Parameter:
421: . order - The variable ordering
423: Output Parameter:
424: . size - The number of local variables
426: Level: intermediate
428: .keywords: variable ordering, grid
429: .seealso: VarOrderingGetSize()
430: @*/
431: int VarOrderingGetLocalSize(VarOrdering order, int *size)
432: {
436: *size = order->numLocVars;
437: return(0);
438: }
440: /*@
441: VarOrderingGetFirst - Gets the layout of variables across domains.
443: Not collective
445: Input Parameter:
446: . order - The variable ordering
448: Output Parameter:
449: . first - An array of the first variable in each domain, and the
450: total number of variables at the end
452: Level: intermediate
454: .keywords: variable ordering, grid
455: .seealso: VarOrderingGetSize()
456: @*/
457: int VarOrderingGetFirst(VarOrdering order, int **first)
458: {
462: *first = order->firstVar;
463: return(0);
464: }
466: /*@
467: VarOrderingSerialize - This function stores or recreates a variable ordering using a viewer for a binary file.
469: Collective on Grid
471: Input Parameter:
472: + map - The class map
473: . viewer - The viewer context
474: - store - This flag is PETSC_TRUE is data is being written, otherwise it will be read
476: Output Parameter:
477: . order - The ordering
479: Level: beginner
481: .keywords: variable ordering, serialize
482: .seealso: GridSerialize()
483: @*/
484: int VarOrderingSerialize(FieldClassMap map, VarOrdering *order, PetscViewer viewer, PetscTruth store)
485: {
486: VarOrdering o;
487: int fd;
488: int cookie;
489: int numOverlapNodes = map->numOverlapNodes;
490: int numGhostNodes = map->numGhostNodes;
491: int numClasses = map->numClasses;
492: int zero = 0;
493: int one = 1;
494: int numProcs, hasLocOffsets;
495: int f, field;
496: PetscTruth match;
497: int ierr;
503: PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_BINARY, &match);
504: if (match == PETSC_FALSE) SETERRQ(PETSC_ERR_ARG_WRONG, "Must be binary viewer");
505: PetscViewerBinaryGetDescriptor(viewer, &fd);
506: MPI_Comm_size(map->comm, &numProcs);
508: if (store) {
509: o = *order;
511: PetscBinaryWrite(fd, &o->cookie, 1, PETSC_INT, 0);
512: PetscBinaryWrite(fd, &o->numVars, 1, PETSC_INT, 0);
513: PetscBinaryWrite(fd, &o->numLocVars, 1, PETSC_INT, 0);
514: PetscBinaryWrite(fd, &o->numOverlapVars, 1, PETSC_INT, 0);
515: PetscBinaryWrite(fd, &o->numNewVars, 1, PETSC_INT, 0);
516: PetscBinaryWrite(fd, &o->numLocNewVars, 1, PETSC_INT, 0);
517: PetscBinaryWrite(fd, &o->numTotalFields, 1, PETSC_INT, 0);
518: PetscBinaryWrite(fd, o->firstVar, numProcs+1, PETSC_INT, 0);
519: PetscBinaryWrite(fd, o->offsets, numOverlapNodes, PETSC_INT, 0);
520: if (o->localOffsets == PETSC_NULL) {
521: PetscBinaryWrite(fd, &zero, 1, PETSC_INT, 0);
522: } else {
523: PetscBinaryWrite(fd, &one, 1, PETSC_INT, 0);
524: PetscBinaryWrite(fd, o->localOffsets, numGhostNodes, PETSC_INT, 0);
525: }
526: for(f = 0; f < map->numFields; f++) {
527: field = map->fields[f];
528: PetscBinaryWrite(fd, o->localStart[field], numClasses, PETSC_INT, 0);
529: }
530: } else {
531: PetscBinaryRead(fd, &cookie, 1, PETSC_INT);
532: if (cookie != VAR_ORDER_COOKIE) SETERRQ(PETSC_ERR_ARG_WRONG, "Non-order object");
534: PetscHeaderCreate(o, _VarOrdering, int, VAR_ORDER_COOKIE, 0, "VarOrdering", map->comm, VarOrderingDestroy, 0);
535: PetscLogObjectCreate(o);
536: PetscObjectCompose((PetscObject) o, "ClassMap", (PetscObject) map);
538: PetscBinaryRead(fd, &o->numVars, 1, PETSC_INT);
539: PetscBinaryRead(fd, &o->numLocVars, 1, PETSC_INT);
540: PetscBinaryRead(fd, &o->numOverlapVars, 1, PETSC_INT);
541: PetscBinaryRead(fd, &o->numNewVars, 1, PETSC_INT);
542: PetscBinaryRead(fd, &o->numLocNewVars, 1, PETSC_INT);
543: PetscBinaryRead(fd, &o->numTotalFields, 1, PETSC_INT);
544: PetscMalloc((numProcs+1) * sizeof(int), &o->firstVar);
545: PetscBinaryRead(fd, o->firstVar, numProcs+1, PETSC_INT);
546: PetscMalloc(numOverlapNodes * sizeof(int), &o->offsets);
547: PetscBinaryRead(fd, o->offsets, numOverlapNodes, PETSC_INT);
548: PetscBinaryRead(fd, &hasLocOffsets, 1, PETSC_INT);
549: if (hasLocOffsets) {
550: PetscMalloc(numGhostNodes * sizeof(int), &o->localOffsets);
551: PetscBinaryRead(fd, o->localOffsets, numGhostNodes, PETSC_INT);
552: }
553: PetscMalloc(o->numTotalFields * sizeof(int *), &o->localStart);
554: for(f = 0; f < map->numFields; f++) {
555: field = map->fields[f];
556: PetscMalloc(numClasses * sizeof(int), &o->localStart[field]);
557: PetscBinaryRead(fd, o->localStart[field], numClasses, PETSC_INT);
558: }
559: PetscLogObjectMemory(o, (numProcs+1 + numOverlapNodes + o->numTotalFields*numClasses) * sizeof(int)
560: + o->numTotalFields * sizeof(int *));
561: *order = o;
562: }
564: return(0);
565: }
567: /*@C VarOrderingDuplicate
568: This function duplicates a global variable ordering.
570: Collective on VarOrdering
572: Input Parameter:
573: . order - The ordering
575: Output Parameter:
576: . neworder - The new ordering
578: Level: beginner
580: .keywords variable ordering, finite element
581: .seealso VarOrderingDuplicateIndices()
582: @*/
583: int VarOrderingDuplicate(VarOrdering order, VarOrdering *neworder)
584: {
585: VarOrdering o;
586: FieldClassMap map;
587: int numFields;
588: int *fields;
589: int numOverlapNodes;
590: int numGhostNodes;
591: int numClasses;
592: int numProcs;
593: int field, newField;
594: int ierr;
600: MPI_Comm_size(order->comm, &numProcs);
602: /* Create the ordering */
603: PetscHeaderCreate(o, _VarOrdering, int, VAR_ORDER_COOKIE, 0, "VarOrdering", order->comm, VarOrderingDestroy, 0);
604: PetscLogObjectCreate(o);
605: PetscObjectQuery((PetscObject) order, "ClassMap", (PetscObject *) &map);
606: PetscObjectCompose((PetscObject) o, "ClassMap", (PetscObject) map);
607: numFields = map->numFields;
608: fields = map->fields;
609: numOverlapNodes = map->numOverlapNodes;
610: numGhostNodes = map->numGhostNodes;
611: numClasses = map->numClasses;
613: /* Set sizes */
614: o->numVars = order->numVars;
615: o->numLocVars = order->numLocVars;
616: o->numOverlapVars = order->numOverlapVars;
617: o->numNewVars = order->numNewVars;
618: o->numLocNewVars = order->numLocNewVars;
619: o->numOverlapNewVars = order->numOverlapNewVars;
620: o->numTotalFields = order->numTotalFields;
621: /* Allocate memory */
622: PetscMalloc((numProcs+1) * sizeof(int), &o->firstVar);
623: PetscMalloc(numOverlapNodes * sizeof(int), &o->offsets);
624: /* Copy data */
625: PetscMemcpy(o->firstVar, order->firstVar, (numProcs+1) * sizeof(int));
626: PetscMemcpy(o->offsets, order->offsets, numOverlapNodes * sizeof(int));
627: PetscLogObjectMemory(o, (numProcs+1 + numOverlapNodes) * sizeof(int));
629: if (order->localOffsets != PETSC_NULL) {
630: PetscMalloc(numGhostNodes * sizeof(int), &o->localOffsets);
631: PetscMemcpy(o->localOffsets, order->localOffsets, numGhostNodes * sizeof(int));
632: PetscLogObjectMemory(o, numGhostNodes * sizeof(int));
633: }
635: PetscMalloc(o->numTotalFields * sizeof(int *), &o->localStart);
636: for(field = 0; field < numFields; field++) {
637: newField = fields[field];
638: PetscMalloc(numClasses * sizeof(int), &o->localStart[newField]);
639: PetscMemcpy(o->localStart[newField], order->localStart[newField], numClasses * sizeof(int));
640: }
641: PetscLogObjectMemory(o, o->numTotalFields * sizeof(int *) + numClasses*numFields * sizeof(int));
643: *neworder = o;
644: return(0);
645: }
647: /*@C VarOrderingDuplicateIndices
648: This function copies the global variable ordering indices to another ordering.
650: Collective on VarOrdering
652: Input Parameter:
653: . order - The original ordering
655: Output Parameter:
656: . neworder - The ordering with updated indices
658: Level: intermediate
660: .keywords variable ordering, finite element
661: .seealso VarOrderingDuplicate()
662: @*/
663: int VarOrderingDuplicateIndices(VarOrdering order, VarOrdering neworder)
664: {
666: SETERRQ(PETSC_ERR_SUP, "Not yet implemented");
667: #if 0
668: if ((mat->rowSize != newmat->rowSize) || (mat->colSize != newmat->colSize)) SETERRQ(1, "Incompatible matrix sizes");
669: PetscMemcpy(newmat->rowIndices, mat->rowIndices, mat->rowSize * sizeof(int));
670: PetscMemcpy(newmat->colIndices, mat->colIndices, mat->colSize * sizeof(int));
671: return(0);
672: #endif
673: }
675: /*@
676: VarOrderingCompatible - This function determines whether the two orderings are shape compatible.
678: Collective on VarOrdering
680: Input Parameters:
681: + order - The variable ordering
682: - anOrder - Another ordering
684: Output Parameter:
685: . compatible - The compatibility flag
687: Level: intermediate
689: .keywords: variable ordering, compatible
690: .seealso: VarOrderingGetLocalSize()
691: @*/
692: int VarOrderingCompatible(VarOrdering order, VarOrdering anOrder, PetscTruth *compatible)
693: {
694: int locCompat = 1;
695: int compat;
702: if ((order->numVars != anOrder->numVars) ||
703: (order->numLocVars != anOrder->numLocVars) ||
704: (order->numOverlapVars != anOrder->numOverlapVars)) {
705: locCompat = 0;
706: }
707: MPI_Allreduce(&locCompat, &compat, 1, MPI_INT, MPI_LAND, order->comm);
708: if (compat) {
709: *compatible = PETSC_TRUE;
710: } else {
711: *compatible = PETSC_FALSE;
712: }
713: return(0);
714: }
716: /*---------------------------------------------- Local Variable Ordering ---------------------------------------------*/
717: /*@C LocalVarOrderingCreate
718: This function creates a local variable ordering for a finite element calculation.
720: Collective on Grid
722: Input Parameter:
723: + grid - The underlying grid for the ordering
724: . numFields - The number of fields in the ordering
725: - fields - The fields in the ordering
727: Output Parameter :
728: . order - The ordering
730: Level: beginner
732: .keywords variable ordering, finite element
733: .seealso LocalVarOrderingDestroy
734: @*/
735: int LocalVarOrderingCreate(Grid grid, int numFields, int *fields, LocalVarOrdering *order)
736: {
743: GridSetUp(grid);
744: (*grid->ops->gridcreatelocalvarordering)(grid, numFields, fields, order);
745: return(0);
746: }
748: /*@C LocalVarOrderingDestroy
749: This function destroys a local variable ordering.
751: Not collective
753: Input Parameter:
754: . order - The ordering
756: Level: beginner
758: .keywords variable ordering, finite element
759: .seealso LocalVarOrderingCreate
760: @*/
761: int LocalVarOrderingDestroy(LocalVarOrdering order)
762: {
767: PetscFree(order->elemStart);
768: PetscFree(order->fields);
769: PetscLogObjectDestroy(order);
770: PetscHeaderDestroy(order);
771: return(0);
772: }
774: /*@
775: LocalVarOrderingSerialize - This function stores or recreates a local variable ordering using a viewer for a binary file.
777: Collective on Grid
779: Input Parameters:
780: + grid - The grid for the ordering
781: . viewer - The viewer context
782: - store - This flag is PETSC_TRUE is data is being written, otherwise it will be read
784: Output Parameter:
785: . order - The ordering
787: Level: beginner
789: .keywords: grid vector, serialize
790: .seealso: GridSerialize()
791: @*/
792: int LocalVarOrderingSerialize(Grid grid, LocalVarOrdering *order, PetscViewer viewer, PetscTruth store)
793: {
794: LocalVarOrdering o;
795: MPI_Comm comm;
796: int fd;
797: int cookie, numFields;
798: PetscTruth match;
799: int ierr;
805: PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_BINARY, &match);
806: if (match == PETSC_FALSE) SETERRQ(PETSC_ERR_ARG_WRONG, "Must be binary viewer");
807: PetscViewerBinaryGetDescriptor(viewer, &fd);
808: GridGetNumFields(grid, &numFields);
810: if (store) {
811: o = *order;
813: PetscBinaryWrite(fd, &o->cookie, 1, PETSC_INT, 0);
815: PetscBinaryWrite(fd, &o->numFields, 1, PETSC_INT, 0);
816: PetscBinaryWrite(fd, o->fields, o->numFields, PETSC_INT, 0);
817: PetscBinaryWrite(fd, &o->elemSize, 1, PETSC_INT, 0);
818: PetscBinaryWrite(fd, o->elemStart, numFields, PETSC_INT, 0);
819: } else {
820: PetscBinaryRead(fd, &cookie, 1, PETSC_INT);
821: if (cookie != VAR_ORDER_COOKIE) SETERRQ(PETSC_ERR_ARG_WRONG, "Non-order object");
823: PetscObjectGetComm((PetscObject) grid, &comm);
824: PetscHeaderCreate(o, _LocalVarOrdering, int, VAR_ORDER_COOKIE, 0, "VarOrdering", comm, LocalVarOrderingDestroy, PETSC_NULL);
825: PetscLogObjectCreate(o);
826: PetscBinaryRead(fd, &o->numFields, 1, PETSC_INT);
828: PetscMalloc(o->numFields * sizeof(int), &o->fields);
829: PetscMalloc(numFields * sizeof(int), &o->elemStart);
830: PetscLogObjectMemory(o, (o->numFields + numFields) * sizeof(int));
831: PetscBinaryRead(fd, o->fields, o->numFields, PETSC_INT);
832: PetscBinaryRead(fd, &o->elemSize, 1, PETSC_INT);
833: PetscBinaryRead(fd, o->elemStart, numFields, PETSC_INT);
834: *order = o;
835: }
837: return(0);
838: }
840: /*@C LocalVarOrderingDuplicate
841: This function duplicates a local variable ordering.
843: Collective on LocalVarOrdering
845: Input Parameter:
846: . order - The ordering
848: Output Parameter:
849: . neworder - The new ordering
851: Level: beginner
853: .keywords variable ordering, finite element
854: .seealso LocalVarOrderingDuplicateIndices()
855: @*/
856: int LocalVarOrderingDuplicate(LocalVarOrdering order, LocalVarOrdering *neworder)
857: {
858: SETERRQ(PETSC_ERR_SUP, "Not yet implemented");
859: #if 0
861: return(0);
862: #endif
863: }
865: /*@C LocalVarOrderingDuplicateIndices
866: This function copies the local variable ordering indices to another ordering.
868: Collective on LocalVarOrdering
870: Input Parameter:
871: . order - The original ordering
873: Output Parameter:
874: . neworder - The ordering with updated indices
876: Level: intermediate
878: .keywords variable ordering, finite element
879: .seealso LocalVarOrderingDuplicate()
880: @*/
881: int LocalVarOrderingDuplicateIndices(LocalVarOrdering order, LocalVarOrdering neworder)
882: {
883: SETERRQ(PETSC_ERR_SUP, "Not yet implemented");
884: #if 0
886: if ((mat->rowSize != newmat->rowSize) || (mat->colSize != newmat->colSize)) SETERRQ(1, "Incompatible matrix sizes");
887: PetscMemcpy(newmat->rowIndices, mat->rowIndices, mat->rowSize * sizeof(int));
888: PetscMemcpy(newmat->colIndices, mat->colIndices, mat->colSize * sizeof(int));
889: return(0);
890: #endif
891: }
893: /*@C
894: LocalVarOrderingGetSize - This function returns the greatest number of variables on any node,
895: or equivalently the number of shape functions in the element matrix.
897: Not collective
899: Input Parameter:
900: . order - The ordering
902: Output Parameter:
903: . size - The size
905: Level: intermediate
907: .keywords local, order, size
908: .seealso LocalVarOrderingGetFieldStart()
909: @*/
910: int LocalVarOrderingGetSize(LocalVarOrdering order, int *size) {
913: *size = order->elemSize;
914: return(0);
915: }
917: /*@C
918: LocalVarOrderingGetField - This function returns the field at a given index in this ordering.
920: Not collective
922: Input Parameters:
923: + order - The ordering
924: - f - The field index
926: Output Parameter:
927: . field - The field
929: Level: intermediate
931: .keywords local, order, field
932: .seealso LocalVarOrderingGetFieldStart()
933: @*/
934: int LocalVarOrderingGetField(LocalVarOrdering order, int f, int *field) {
937: if ((f < 0) || (f >= order->numFields)) {
938: SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE, "Invalid field index %d must be in [0,%d)", f, order->numFields);
939: }
940: *field = order->fields[f];
941: return(0);
942: }
944: /*@C
945: LocalVarOrderingGetFieldIndex - This function returns the field index in this ordering of the given field.
947: Not collective
949: Input Parameters:
950: + order - The ordering
951: - field - The field
953: Output Parameter:
954: . f - The field index
956: Level: intermediate
958: .keywords local, order, field, index
959: .seealso LocalVarOrderingGetFieldStart()
960: @*/
961: int LocalVarOrderingGetFieldIndex(LocalVarOrdering order, int field, int *f) {
962: int index;
966: for(index = 0; index < order->numFields; index++) {
967: if (order->fields[index] == field) break;
968: }
969: if (index >= order->numFields) SETERRQ1(PETSC_ERR_ARG_WRONG, "Invalid field %d not present in order", field);
970: *f = index;
971: return(0);
972: }
974: /*@C
975: LocalVarOrderingGetFieldStart - This function returns the first index of the given field in the element vector or matrix.
977: Not collective
979: Input Parameters:
980: + order - The ordering
981: - field - The field
983: Output Parameter:
984: . start - The element vector or matrix index, or -1 if the field is not active
986: Level: intermediate
988: .keywords local, order, start, field
989: .seealso LocalVarGetSize()
990: @*/
991: int LocalVarOrderingGetFieldStart(LocalVarOrdering order, int field, int *start) {
994: *start = order->elemStart[field];
995: return(0);
996: }