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: }