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: #undef  __FUNCT__
 13: /*@C
 14:   VarOrderingCreate - This function creates a global variable ordering.

 16:   Collective on Grid

 18:   Input Parameter:
 19: . grid  - The underlying grid for the ordering

 21:   Output Parameter:
 22: . order - The ordering

 24:   Level: beginner

 26: .keywords: variable ordering, create
 27: .seealso: VarOrderingCreateConstrained(), VarOrderingDestroy()
 28: @*/
 29: int VarOrderingCreate(Grid grid, VarOrdering *order)
 30: {

 36:   GridSetUp(grid);
 37:   (*grid->ops->gridcreatevarordering)(grid, grid->cm, order);
 38:   return(0);
 39: }

 41: #undef  __FUNCT__
 43: /*@C
 44:   VarOrderingCreateConstrained - This function creates a global variable ordering on the constraint space.

 46:   Collective on Grid

 48:   Input Parameter:
 49: . grid  - The underlying grid for the ordering

 51:   Output Parameter:
 52: . order - The ordering

 54:   Level: beginner

 56: .keywords: variable ordering, create, constraint
 57: .seealso: VarOrderingCreate(), VarOrderingDestroy()
 58: @*/
 59: int VarOrderingCreateConstrained(Grid grid, VarOrdering *order)
 60: {
 61:   FieldClassMap cm;
 62:   int           ierr;

 67:   GridSetUp(grid);
 68:   GridGetClassMap(grid, &cm);
 69:   (*grid->ops->gridcreatevarordering)(grid, cm, order);
 70:   return(0);
 71: }

 73: #undef  __FUNCT__
 75: /*@C
 76:   VarOrderingCreateGeneral - This function creates a global variable ordering on the specified fields.

 78:   Collective on Grid

 80:   Input Parameters:
 81: + grid      - The underlying grid for the ordering
 82: . numFields - The number of fields
 83: - fields    - The fields

 85:   Output Parameter:
 86: . order     - The ordering

 88:   Level: beginner

 90: .keywords: variable ordering, create, constraint
 91: .seealso: VarOrderingCreate(), VarOrderingDestroy()
 92: @*/
 93: int VarOrderingCreateGeneral(Grid grid, int numFields, int *fields, VarOrdering *order)
 94: {
 95:   FieldClassMap map, tempMap;
 96:   PetscTruth    isConstrained, reduceSystem;
 97:   int           numTotalFields, f;
 98:   int           ierr;

103:   GridGetNumFields(grid, &numTotalFields);
104:   if (numFields > numTotalFields) {
105:     SETERRQ2(PETSC_ERR_ARG_WRONG, "Invalid number  %d of fields, grid only has %d", numFields, numTotalFields);
106:   }
107:   for(f = 0; f < numFields; f++) {
108:     if ((fields[f] < 0) || (fields[f] >= numTotalFields)) {
109:       SETERRQ1(PETSC_ERR_ARG_WRONG, "Invalid field number %d", fields[f]);
110:     }
111:   }
112:   GridSetUp(grid);
113:   FieldClassMapCreateTriangular2D(grid, numFields, fields, &tempMap);
114:   GridIsConstrained(grid, &isConstrained);
115:   GridGetReduceSystem(grid, &reduceSystem);
116:   if (reduceSystem == PETSC_TRUE) {
117:     FieldClassMapConstrain(tempMap, grid, reduceSystem, isConstrained, &map);
118:     FieldClassMapDestroy(tempMap);
119:   } else {
120:     map  = tempMap;
121:   }
122:   (*grid->ops->gridcreatevarordering)(grid, map, order);
123:   FieldClassMapDestroy(map);
124:   return(0);
125: }

127: #undef  __FUNCT__
129: /*@C
130:   VarOrderingCreateSubset - This function creates a new global variable ordering from an existing
131:   one on a subset of the fields.

133:   Collective on Grid

135:   Input Parameters:
136: + order     - The original ordering
137: . numFields - The number of fields in the new ordering
138: . fields    - The fields in the new ordering
139: - contract  - The flag for contracting the indices

141:   Output Parameter:
142: . newOrder  - The new ordering

144:   Note:
145:   If contract is PETSC_FALSE, this function returns a copy of the original
146:   ordering, but with some fields no longer active.

148:   Level: beginner

150: .keywords: variable ordering, create
151: .seealso: VarOrderingCreate(), VarOrderingDestroy()
152: @*/
153: int VarOrderingCreateSubset(VarOrdering order, int numFields, int *fields, PetscTruth contract, VarOrdering *newOrder)
154: {
155:   FieldClassMap map, cm;
156:   int           f, g;
157:   int           ierr;

162:   VarOrderingGetClassMap(order, &map);
163:   /* Check for consistency */
164:   for(f = 0; f < numFields; f++) {
165:     for(g = 0; g < map->numFields; g++) {
166:       if (map->fields[g] == fields[f])
167:         break;
168:     }
169:     if (g == map->numFields) {
170:       SETERRQ1(PETSC_ERR_ARG_INCOMP, "Fields not a subset: field %d not in the existing order", fields[f]);
171:     }
172:   }
173:   if (contract == PETSC_TRUE) {
174:     SETERRQ(PETSC_ERR_SUP, "Coming soon");
175:   } else {
176:     FieldClassMapCreateSubset(map, numFields, fields, &cm);
177:     VarOrderingDuplicate(order, newOrder);
178:     PetscObjectCompose((PetscObject) *newOrder, "ClassMap", (PetscObject) cm);
179:     FieldClassMapDestroy(cm);
180:     /* Free extra localStart entries */
181:     for(g = 0; g < map->numFields; g++) {
182:       for(f = 0; f < numFields; f++) {
183:         if (map->fields[g] == fields[f])
184:           break;
185:       }
186:       if (f == numFields) {
187:         PetscFree((*newOrder)->localStart[map->fields[g]]);
188:       }
189:     }
190:   }
191:   return(0);
192: }

194: #undef  __FUNCT__
196: /*@C
197:   VarOrderingConstrain - This function constrains a previous variable ordering using the constraints
198:   already defined in the grid.

200:   Collective on VarOrdering

202:   Input Parameters:
203: + grid     - The underlying grid for the ordering
204: - oldOrder - The original ordering

206:   Output Parameter:
207: . order    - The constrained ordering

209:   Level: beginner

211: .keywords variable ordering, finite element
212: .seealso VarOrderingDestroy()
213: @*/
214: int VarOrderingConstrain(Grid grid, VarOrdering oldOrder, VarOrdering *order)
215: {
216:   PetscTruth    isConstrained;
217:   FieldClassMap cm;
218:   int           numFields;
219:   int           field;
220:   int           ierr;

225:   GridSetUp(grid);
226:   GridIsConstrained(grid, &isConstrained);
227:   if (isConstrained == PETSC_TRUE) {
228:     GridGetNumFields(grid, &numFields);
229:     GridGetClassMap(grid, &cm);
230:     for(field = 0; field < numFields; field++) {
231:       if (grid->fields[field].isConstrained == PETSC_TRUE) break;
232:     }
233:     if (field == numFields) SETERRQ(PETSC_ERR_ARG_WRONG, "No constrained fields");
234:     (*grid->ops->gridvarorderingconstrain)(grid, field, grid->constraintCtx, cm, oldOrder, order);
235:   } else {
236:     VarOrderingDuplicate(oldOrder, order);
237:   }
238:   return(0);
239: }

241: #undef  __FUNCT__
243: /*@C
244:   VarOrderingCreateReduce - This function creates a global variable ordering on the reduction space.

246:   Collective on Grid

248:   Input Parameter:
249: . grid  - The underlying grid for the ordering

251:   Output Parameter:
252: . order - The reduction ordering

254:   Level: beginner

256: .keywords: variable ordering, create, reduction
257: .seealso: VarOrderingCreate(), VarOrderingDestroy()
258: @*/
259: int VarOrderingCreateReduce(Grid grid, VarOrdering *order)
260: {

266:   GridSetUp(grid);
267:   if (grid->reduceSystem == PETSC_TRUE) {
268:     (*grid->ops->gridcreatevarordering)(grid, grid->reductionCM, order);
269:   } else {
270:     (*grid->ops->gridcreatevarordering)(grid, grid->cm,          order);
271:   }
272:   return(0);
273: }

275: #undef  __FUNCT__
277: /*@C
278:   VarOrderingCreateReduceGeneral - This function creates a global variable ordering on the reduction space.

280:   Collective on Grid

282:   Input Parameters:
283: + grid        - The underlying grid for the ordering
284: . cm          - The class map
285: - reductionCM - The class map for the reduction space

287:   Output Parameter:
288: . order       - The reduction ordering

290:   Level: beginner

292: .keywords: variable ordering, create, reduction
293: .seealso: VarOrderingCreate(), VarOrderingDestroy()
294: @*/
295: int VarOrderingCreateReduceGeneral(Grid grid, FieldClassMap cm, FieldClassMap reductionCM, VarOrdering *order)
296: {

304:   GridSetUp(grid);
305:   if (grid->reduceSystem == PETSC_TRUE) {
306:     (*grid->ops->gridcreatevarordering)(grid, reductionCM, order);
307:   } else {
308:     (*grid->ops->gridcreatevarordering)(grid, cm,          order);
309:   }
310:   return(0);
311: }

313: #undef  __FUNCT__
315: /*@C
316:   VarOrderingDestroy - This function destroys a global variable ordering.

318:   Not collective

320:   Input Parameter:
321: . order - The ordering

323:   Level: beginner

325: .keywords: variable ordering
326: .seealso: VarOrderingCreate()
327: @*/
328: int VarOrderingDestroy(VarOrdering order)
329: {
330:   FieldClassMap map;
331:   int           f;
332:   int           ierr;

336:   if (--order->refct > 0)
337:     return(0);
338:   PetscObjectQuery((PetscObject) order, "ClassMap", (PetscObject *) &map);
339:   PetscFree(order->firstVar);
340:   PetscFree(order->offsets);
341:   if (order->localOffsets) {
342:     PetscFree(order->localOffsets);
343:   }
344:   for(f = 0; f < map->numFields; f++) {
345:     PetscFree(order->localStart[map->fields[f]]);
346:   }
347:   PetscFree(order->localStart);
348:   PetscLogObjectDestroy(order);
349:   PetscHeaderDestroy(order);
350:   return(0);
351: }

353: #undef  __FUNCT__
355: /*@
356:   VarOrderingGetClassMap - This function returns the class map for the variable ordering.

358:   Not collective

360:   Input Parameter:
361: . order - The variable ordering

363:   Output Parameter:
364: . map   - The field class map

366:   Level: intermediate

368: .keywords: variable ordering, class map
369: .seealso: VarOrderingGetLocalSize()
370: @*/
371: int VarOrderingGetClassMap(VarOrdering order, FieldClassMap *map)
372: {

378:   PetscObjectQuery((PetscObject) order, "ClassMap", (PetscObject *) map);
380:   return(0);
381: }

383: #undef  __FUNCT__
385: /*@
386:   VarOrderingGetNumTotalFields - Gets the total number of fields in the associated grid.

388:   Not collective

390:   Input Parameter:
391: . order     - The variable ordering

393:   Output Parameter:
394: . numFields - The total number fields

396:   Level: intermediate

398: .keywords: variable ordering, grid, field
399: .seealso: VarOrderingGetSize(), VarOrderingGetLocalSize()
400: @*/
401: int VarOrderingGetNumTotalFields(VarOrdering order, int *numFields)
402: {
406:   *numFields = order->numTotalFields;
407:   return(0);
408: }

410: #undef  __FUNCT__
412: /*@
413:   VarOrderingGetSize - Gets the number of global variables in a variable ordering.

415:   Not collective

417:   Input Parameter:
418: . order - The variable ordering

420:   Output Parameter:
421: . size  - The number of global variables

423:   Level: intermediate

425: .keywords: variable ordering, grid
426: .seealso: VarOrderingGetLocalSize()
427: @*/
428: int VarOrderingGetSize(VarOrdering order, int *size)
429: {
433:   *size = order->numVars;
434:   return(0);
435: }

437: #undef  __FUNCT__
439: /*@
440:   VarOrderingGetLocalSize - Gets the number of local variables in a variable ordering.

442:   Not collective

444:   Input Parameter:
445: . order - The variable ordering

447:   Output Parameter:
448: . size  - The number of local variables

450:   Level: intermediate

452: .keywords: variable ordering, grid
453: .seealso: VarOrderingGetSize()
454: @*/
455: int VarOrderingGetLocalSize(VarOrdering order, int *size)
456: {
460:   *size = order->numLocVars;
461:   return(0);
462: }

464: #undef  __FUNCT__
466: /*@
467:   VarOrderingGetFirst - Gets the layout of variables across domains.

469:   Not collective

471:   Input Parameter:
472: . order - The variable ordering

474:   Output Parameter:
475: . first - An array of the first variable in each domain, and the
476:           total number of variables at the end

478:   Level: intermediate

480: .keywords: variable ordering, grid
481: .seealso: VarOrderingGetSize()
482: @*/
483: int VarOrderingGetFirst(VarOrdering order, int **first)
484: {
488:   *first = order->firstVar;
489:   return(0);
490: }

492: #undef  __FUNCT__
494: /*@ 
495:   VarOrderingSerialize - This function stores or recreates a variable ordering using a viewer for a binary file.

497:   Collective on Grid

499:   Input Parameter:
500: + map    - The class map
501: . viewer - The viewer context
502: - store  - This flag is PETSC_TRUE is data is being written, otherwise it will be read

504:   Output Parameter:
505: . order  - The ordering

507:   Level: beginner

509: .keywords: variable ordering, serialize
510: .seealso: GridSerialize()
511: @*/
512: int VarOrderingSerialize(FieldClassMap map, VarOrdering *order, PetscViewer viewer, PetscTruth store)
513: {
514:   VarOrdering o;
515:   int         fd;
516:   int         cookie;
517:   int         numOverlapNodes = map->numOverlapNodes;
518:   int         numGhostNodes   = map->numGhostNodes;
519:   int         numClasses      = map->numClasses;
520:   int         zero            = 0;
521:   int         one             = 1;
522:   int         numProcs, hasLocOffsets;
523:   int         f, field;
524:   PetscTruth  match;
525:   int         ierr;


531:   PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_BINARY, &match);
532:   if (match == PETSC_FALSE) SETERRQ(PETSC_ERR_ARG_WRONG, "Must be binary viewer");
533:   PetscViewerBinaryGetDescriptor(viewer, &fd);
534:   MPI_Comm_size(map->comm, &numProcs);

536:   if (store) {
537:     o = *order;
539:     PetscBinaryWrite(fd, &o->cookie,             1,               PETSC_INT, 0);
540:     PetscBinaryWrite(fd, &o->numVars,            1,               PETSC_INT, 0);
541:     PetscBinaryWrite(fd, &o->numLocVars,         1,               PETSC_INT, 0);
542:     PetscBinaryWrite(fd, &o->numOverlapVars,     1,               PETSC_INT, 0);
543:     PetscBinaryWrite(fd, &o->numNewVars,         1,               PETSC_INT, 0);
544:     PetscBinaryWrite(fd, &o->numLocNewVars,      1,               PETSC_INT, 0);
545:     PetscBinaryWrite(fd, &o->numTotalFields,     1,               PETSC_INT, 0);
546:     PetscBinaryWrite(fd,  o->firstVar,           numProcs+1,      PETSC_INT, 0);
547:     PetscBinaryWrite(fd,  o->offsets,            numOverlapNodes, PETSC_INT, 0);
548:     if (o->localOffsets == PETSC_NULL) {
549:       PetscBinaryWrite(fd, &zero,                1,               PETSC_INT, 0);
550:     } else {
551:       PetscBinaryWrite(fd, &one,                 1,               PETSC_INT, 0);
552:       PetscBinaryWrite(fd,  o->localOffsets,     numGhostNodes,   PETSC_INT, 0);
553:     }
554:     for(f = 0; f < map->numFields; f++) {
555:       field = map->fields[f];
556:       PetscBinaryWrite(fd, o->localStart[field], numClasses,      PETSC_INT, 0);
557:     }
558:   } else {
559:     PetscBinaryRead(fd, &cookie,                1,               PETSC_INT);
560:     if (cookie != VAR_ORDER_COOKIE) SETERRQ(PETSC_ERR_ARG_WRONG, "Non-order object");

562:     PetscHeaderCreate(o, _VarOrdering, int, VAR_ORDER_COOKIE, 0, "VarOrdering", map->comm, VarOrderingDestroy, 0);
563:     PetscLogObjectCreate(o);
564:     PetscObjectCompose((PetscObject) o, "ClassMap", (PetscObject) map);

566:     PetscBinaryRead(fd, &o->numVars,            1,               PETSC_INT);
567:     PetscBinaryRead(fd, &o->numLocVars,         1,               PETSC_INT);
568:     PetscBinaryRead(fd, &o->numOverlapVars,     1,               PETSC_INT);
569:     PetscBinaryRead(fd, &o->numNewVars,         1,               PETSC_INT);
570:     PetscBinaryRead(fd, &o->numLocNewVars,      1,               PETSC_INT);
571:     PetscBinaryRead(fd, &o->numTotalFields,     1,               PETSC_INT);
572:     PetscMalloc((numProcs+1)      * sizeof(int),   &o->firstVar);
573:     PetscBinaryRead(fd,  o->firstVar,           numProcs+1,      PETSC_INT);
574:     PetscMalloc(numOverlapNodes   * sizeof(int),   &o->offsets);
575:     PetscBinaryRead(fd,  o->offsets,            numOverlapNodes, PETSC_INT);
576:     PetscBinaryRead(fd, &hasLocOffsets,         1,               PETSC_INT);
577:     if (hasLocOffsets) {
578:       PetscMalloc(numGhostNodes   * sizeof(int),   &o->localOffsets);
579:       PetscBinaryRead(fd,  o->localOffsets,     numGhostNodes,   PETSC_INT);
580:     }
581:     PetscMalloc(o->numTotalFields * sizeof(int *), &o->localStart);
582:     for(f = 0; f < map->numFields; f++) {
583:       field = map->fields[f];
584:       PetscMalloc(numClasses      * sizeof(int),   &o->localStart[field]);
585:       PetscBinaryRead(fd, o->localStart[field], numClasses,      PETSC_INT);
586:     }
587:     PetscLogObjectMemory(o, (numProcs+1 + numOverlapNodes + o->numTotalFields*numClasses) * sizeof(int)
588:                      + o->numTotalFields * sizeof(int *));
589:     *order = o;
590:   }

592:   return(0);
593: }

595: #undef  __FUNCT__
597: /*@C VarOrderingDuplicate
598:   This function duplicates a global variable ordering.

600:   Collective on VarOrdering

602:   Input Parameter:
603: . order    - The ordering

605:   Output Parameter:
606: . neworder - The new ordering

608:   Level: beginner

610: .keywords variable ordering, finite element
611: .seealso VarOrderingDuplicateIndices()
612: @*/
613: int VarOrderingDuplicate(VarOrdering order, VarOrdering *neworder)
614: {
615:   VarOrdering   o;
616:   FieldClassMap map;
617:   int           numFields;
618:   int          *fields;
619:   int           numOverlapNodes;
620:   int           numGhostNodes;
621:   int           numClasses;
622:   int           numProcs;
623:   int           field, newField;
624:   int           ierr;


630:   MPI_Comm_size(order->comm, &numProcs);

632:   /* Create the ordering */
633:   PetscHeaderCreate(o, _VarOrdering, int, VAR_ORDER_COOKIE, 0, "VarOrdering", order->comm, VarOrderingDestroy, 0);
634:   PetscLogObjectCreate(o);
635:   PetscObjectQuery((PetscObject) order, "ClassMap", (PetscObject *) &map);
636:   PetscObjectCompose((PetscObject) o, "ClassMap", (PetscObject) map);
637:   numFields       = map->numFields;
638:   fields          = map->fields;
639:   numOverlapNodes = map->numOverlapNodes;
640:   numGhostNodes   = map->numGhostNodes;
641:   numClasses      = map->numClasses;

643:   /* Set sizes */
644:   o->numVars           = order->numVars;
645:   o->numLocVars        = order->numLocVars;
646:   o->numOverlapVars    = order->numOverlapVars;
647:   o->numNewVars        = order->numNewVars;
648:   o->numLocNewVars     = order->numLocNewVars;
649:   o->numOverlapNewVars = order->numOverlapNewVars;
650:   o->numTotalFields    = order->numTotalFields;
651:   /* Allocate memory */
652:   PetscMalloc((numProcs+1)    * sizeof(int), &o->firstVar);
653:   PetscMalloc(numOverlapNodes * sizeof(int), &o->offsets);
654:   /* Copy data */
655:   PetscMemcpy(o->firstVar, order->firstVar, (numProcs+1)    * sizeof(int));
656:   PetscMemcpy(o->offsets,  order->offsets,  numOverlapNodes * sizeof(int));
657:   PetscLogObjectMemory(o, (numProcs+1 + numOverlapNodes) * sizeof(int));

659:   if (order->localOffsets != PETSC_NULL) {
660:     PetscMalloc(numGhostNodes * sizeof(int), &o->localOffsets);
661:     PetscMemcpy(o->localOffsets, order->localOffsets, numGhostNodes * sizeof(int));
662:     PetscLogObjectMemory(o, numGhostNodes * sizeof(int));
663:   }

665:   PetscMalloc(o->numTotalFields * sizeof(int *), &o->localStart);
666:   for(field = 0; field < numFields; field++) {
667:     newField = fields[field];
668:     PetscMalloc(numClasses * sizeof(int), &o->localStart[newField]);
669:     PetscMemcpy(o->localStart[newField], order->localStart[newField], numClasses * sizeof(int));
670:   }
671:   PetscLogObjectMemory(o, o->numTotalFields * sizeof(int *) + numClasses*numFields * sizeof(int));

673:   *neworder = o;
674:   return(0);
675: }

677: #undef  __FUNCT__
679: /*@C VarOrderingDuplicateIndices
680:   This function copies the global variable ordering indices to another ordering.

682:   Collective on VarOrdering

684:   Input Parameter:
685: . order    - The original ordering

687:   Output Parameter:
688: . neworder - The ordering with updated indices

690:   Level: intermediate

692: .keywords variable ordering, finite element
693: .seealso VarOrderingDuplicate()
694: @*/
695: int VarOrderingDuplicateIndices(VarOrdering order, VarOrdering neworder)
696: {
698:   SETERRQ(PETSC_ERR_SUP, "Not yet implemented");
699: #if 0
700:   if ((mat->rowSize != newmat->rowSize) || (mat->colSize != newmat->colSize)) SETERRQ(1, "Incompatible matrix sizes");
701:   PetscMemcpy(newmat->rowIndices, mat->rowIndices, mat->rowSize * sizeof(int));
702:   PetscMemcpy(newmat->colIndices, mat->colIndices, mat->colSize * sizeof(int));
703:   return(0);
704: #endif
705: }

707: #undef  __FUNCT__
709: /*@
710:   VarOrderingCompatible - This function determines whether the two orderings are shape compatible.

712:   Collective on VarOrdering

714:   Input Parameters:
715: + order      - The variable ordering
716: - anOrder    - Another ordering

718:   Output Parameter:
719: . compatible - The compatibility flag

721:   Level: intermediate

723: .keywords: variable ordering, compatible
724: .seealso: VarOrderingGetLocalSize()
725: @*/
726: int VarOrderingCompatible(VarOrdering order, VarOrdering anOrder, PetscTruth *compatible)
727: {
728:   int locCompat = 1;
729:   int compat;

736:   if ((order->numVars        != anOrder->numVars)    ||
737:       (order->numLocVars     != anOrder->numLocVars) ||
738:       (order->numOverlapVars != anOrder->numOverlapVars)) {
739:     locCompat = 0;
740:   }
741:   MPI_Allreduce(&locCompat, &compat, 1, MPI_INT, MPI_LAND, order->comm);
742:   if (compat) {
743:     *compatible = PETSC_TRUE;
744:   } else {
745:     *compatible = PETSC_FALSE;
746:   }
747:   return(0);
748: }

750: /*---------------------------------------------- Local Variable Ordering ---------------------------------------------*/
751: #undef  __FUNCT__
753: /*@C LocalVarOrderingCreate
754:   This function creates a local variable ordering for a finite element calculation.

756:   Collective on Grid

758:   Input Parameter:
759: + grid      - The underlying grid for the ordering
760: . numFields - The number of fields in the ordering
761: - fields    - The fields in the ordering

763:   Output Parameter :
764: . order     - The ordering

766:   Level: beginner

768: .keywords variable ordering, finite element
769: .seealso LocalVarOrderingDestroy
770: @*/
771: int LocalVarOrderingCreate(Grid grid, int numFields, int *fields, LocalVarOrdering *order)
772: {

779:   GridSetUp(grid);
780:   (*grid->ops->gridcreatelocalvarordering)(grid, numFields, fields, order);
781:   return(0);
782: }

784: #undef  __FUNCT__
786: /*@C LocalVarOrderingDestroy
787:   This function destroys a local variable ordering.

789:   Not collective

791:   Input Parameter:
792: . order - The ordering

794:   Level: beginner

796: .keywords variable ordering, finite element
797: .seealso LocalVarOrderingCreate
798: @*/
799: int LocalVarOrderingDestroy(LocalVarOrdering order)
800: {

805:   PetscFree(order->elemStart);
806:   PetscFree(order->fields);
807:   PetscLogObjectDestroy(order);
808:   PetscHeaderDestroy(order);
809:   return(0);
810: }

812: #undef  __FUNCT__
814: /*@ 
815:   LocalVarOrderingSerialize - This function stores or recreates a local variable ordering using a viewer for a binary file.

817:   Collective on Grid

819:   Input Parameters:
820: + grid   - The grid for the ordering
821: . viewer - The viewer context
822: - store  - This flag is PETSC_TRUE is data is being written, otherwise it will be read

824:   Output Parameter:
825: . order  - The ordering

827:   Level: beginner

829: .keywords: grid vector, serialize
830: .seealso: GridSerialize()
831: @*/
832: int LocalVarOrderingSerialize(Grid grid, LocalVarOrdering *order, PetscViewer viewer, PetscTruth store)
833: {
834:   LocalVarOrdering o;
835:   MPI_Comm         comm;
836:   int              fd;
837:   int              cookie, numFields;
838:   PetscTruth       match;
839:   int              ierr;


845:   PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_BINARY, &match);
846:   if (match == PETSC_FALSE) SETERRQ(PETSC_ERR_ARG_WRONG, "Must be binary viewer");
847:   PetscViewerBinaryGetDescriptor(viewer, &fd);
848:   GridGetNumFields(grid, &numFields);

850:   if (store) {
851:     o = *order;
853:     PetscBinaryWrite(fd, &o->cookie,    1,             PETSC_INT,  0);

855:     PetscBinaryWrite(fd, &o->numFields, 1,             PETSC_INT,  0);
856:     PetscBinaryWrite(fd,  o->fields,    o->numFields,  PETSC_INT,  0);
857:     PetscBinaryWrite(fd, &o->elemSize,  1,             PETSC_INT,  0);
858:     PetscBinaryWrite(fd,  o->elemStart, numFields,     PETSC_INT,  0);
859:   } else {
860:     PetscBinaryRead(fd, &cookie,       1,             PETSC_INT);
861:     if (cookie != VAR_ORDER_COOKIE) SETERRQ(PETSC_ERR_ARG_WRONG, "Non-order object");

863:     PetscObjectGetComm((PetscObject) grid, &comm);
864:     PetscHeaderCreate(o, _LocalVarOrdering, int, VAR_ORDER_COOKIE, 0, "VarOrdering", comm, LocalVarOrderingDestroy, PETSC_NULL);
865:     PetscLogObjectCreate(o);
866:     PetscBinaryRead(fd, &o->numFields, 1,             PETSC_INT);

868:     PetscMalloc(o->numFields * sizeof(int), &o->fields);
869:     PetscMalloc(numFields    * sizeof(int), &o->elemStart);
870:     PetscLogObjectMemory(o, (o->numFields + numFields) * sizeof(int));
871:     PetscBinaryRead(fd,  o->fields,    o->numFields,  PETSC_INT);
872:     PetscBinaryRead(fd, &o->elemSize,  1,             PETSC_INT);
873:     PetscBinaryRead(fd,  o->elemStart, numFields,     PETSC_INT);
874:     *order = o;
875:   }

877:   return(0);
878: }

880: #undef  __FUNCT__
882: /*@C LocalVarOrderingDuplicate
883:   This function duplicates a local variable ordering.

885:   Collective on LocalVarOrdering

887:   Input Parameter:
888: . order    - The ordering

890:   Output Parameter:
891: . neworder - The new ordering

893:   Level: beginner

895: .keywords variable ordering, finite element
896: .seealso LocalVarOrderingDuplicateIndices()
897: @*/
898: int LocalVarOrderingDuplicate(LocalVarOrdering order, LocalVarOrdering *neworder)
899: {
900:   SETERRQ(PETSC_ERR_SUP, "Not yet implemented");
901: #if 0
903:   return(0);
904: #endif
905: }

907: #undef  __FUNCT__
909: /*@C LocalVarOrderingDuplicateIndices
910:   This function copies the local variable ordering indices to another ordering.

912:   Collective on LocalVarOrdering

914:   Input Parameter:
915: . order    - The original ordering

917:   Output Parameter:
918: . neworder - The ordering with updated indices

920:   Level: intermediate

922: .keywords variable ordering, finite element
923: .seealso LocalVarOrderingDuplicate()
924: @*/
925: int LocalVarOrderingDuplicateIndices(LocalVarOrdering order, LocalVarOrdering neworder)
926: {
927:   SETERRQ(PETSC_ERR_SUP, "Not yet implemented");
928: #if 0
930:   if ((mat->rowSize != newmat->rowSize) || (mat->colSize != newmat->colSize)) SETERRQ(1, "Incompatible matrix sizes");
931:   PetscMemcpy(newmat->rowIndices, mat->rowIndices, mat->rowSize * sizeof(int));
932:   PetscMemcpy(newmat->colIndices, mat->colIndices, mat->colSize * sizeof(int));
933:   return(0);
934: #endif
935: }

937: #undef  __FUNCT__
939: /*@C
940:   LocalVarOrderingGetSize - This function returns the greatest number of variables on any node,
941:   or equivalently the number of shape functions in the element matrix.

943:   Not collective

945:   Input Parameter:
946: . order - The ordering

948:   Output Parameter:
949: . size  - The size

951:   Level: intermediate

953: .keywords local, order, size
954: .seealso LocalVarOrderingGetFieldStart()
955: @*/
956: int LocalVarOrderingGetSize(LocalVarOrdering order, int *size) {
959:   *size = order->elemSize;
960:   return(0);
961: }

963: #undef  __FUNCT__
965: /*@C
966:   LocalVarOrderingGetField - This function returns the field at a given index in this ordering.

968:   Not collective

970:   Input Parameters:
971: + order - The ordering
972: - f     - The field index

974:   Output Parameter:
975: . field - The field

977:   Level: intermediate

979: .keywords local, order, field
980: .seealso LocalVarOrderingGetFieldStart()
981: @*/
982: int LocalVarOrderingGetField(LocalVarOrdering order, int f, int *field) {
985:   if ((f < 0) || (f >= order->numFields)) {
986:     SETERRQ2(PETSC_ERR_ARG_OUTOFRANGE, "Invalid field index %d must be in [0,%d)", f, order->numFields);
987:   }
988:   *field = order->fields[f];
989:   return(0);
990: }

992: #undef  __FUNCT__
994: /*@C
995:   LocalVarOrderingGetFieldIndex - This function returns the field index in this ordering of the given field.

997:   Not collective

999:   Input Parameters:
1000: + order - The ordering
1001: - field - The field

1003:   Output Parameter:
1004: . f     - The field index

1006:   Level: intermediate

1008: .keywords local, order, field, index
1009: .seealso LocalVarOrderingGetFieldStart()
1010: @*/
1011: int LocalVarOrderingGetFieldIndex(LocalVarOrdering order, int field, int *f) {
1012:   int index;

1016:   for(index = 0; index < order->numFields; index++) {
1017:     if (order->fields[index] == field) break;
1018:   }
1019:   if (index >= order->numFields) SETERRQ1(PETSC_ERR_ARG_WRONG, "Invalid field %d not present in order", field);
1020:   *f = index;
1021:   return(0);
1022: }

1024: #undef  __FUNCT__
1026: /*@C
1027:   LocalVarOrderingGetFieldStart - This function returns the first index of the given field in the element vector or matrix.

1029:   Not collective

1031:   Input Parameters:
1032: + order - The ordering
1033: - field - The field

1035:   Output Parameter:
1036: . start - The element vector or matrix index, or -1 if the field is not active

1038:   Level: intermediate

1040: .keywords local, order, start, field
1041: .seealso LocalVarGetSize()
1042: @*/
1043: int LocalVarOrderingGetFieldStart(LocalVarOrdering order, int field, int *start) {
1046:   *start = order->elemStart[field];
1047:   return(0);
1048: }