Actual source code: discretization.c

  1: #ifdef PETSC_RCS_HEADER
  2: static char vcid[] = "$Id: disc.c,v 1.3 2000/01/10 03:54:16 knepley Exp $";
  3: #endif

  5: /*
  6:   This is the interface to discretization objects
  7: */

 9:  #include src/grid/discretization/discimpl.h

 11: /* Logging support */
 12: int DISCRETIZATION_COOKIE;
 13: int DISCRETIZATION_EvaluateOperatorGalerkin;

 15: /*------------------------------------------------- Generic Functions ------------------------------------------------*/
 16: #undef  __FUNCT__
 18: /*@
 19:   DiscretizationDestroy - Destroys a discretization object.

 21:   Collective on Discretization

 23:   Input Parameter:
 24: . disc - The discretization

 26:   Level: beginner

 28: .keywords: discretization, destroy
 29: .seealso: DiscretizationView()
 30: @*/
 31: int DiscretizationDestroy(Discretization disc) {
 32:   int op, arg;

 37:   if (--disc->refct > 0) return(0);
 38:   if (disc->bdDisc != PETSC_NULL) {
 39:     DiscretizationDestroy(disc->bdDisc);
 40:   }
 41:   PetscFree(disc->quadPoints);
 42:   PetscFree(disc->quadWeights);
 43:   PetscFree(disc->quadShapeFuncs);
 44:   PetscFree(disc->quadShapeFuncDers);
 45:   PetscFree(disc->funcVal);
 46:   for(op = 0; op < disc->numOps; op++) {
 47:     if (disc->operators[op]->precompInt) {
 48:       PetscFree(disc->operators[op]->precompInt);
 49:     }
 50:     PetscFree(disc->operators[op]);
 51:   }
 52:   PetscFree(disc->operators);
 53:   for(arg = 0; arg < 2; arg++) {
 54:     PetscFree(disc->fieldVal[arg]);
 55:   }
 56:   PetscFree(disc->fieldVal);
 57:   (*disc->ops->destroy)(disc);
 58:   PetscLogObjectDestroy(disc);
 59:   PetscHeaderDestroy(disc);
 60:   return(0);
 61: }

 63: #undef  __FUNCT__
 65: /*@
 66:   DiscretizationView - Views a discretization object.

 68:   Collective on Discretization

 70:   Input Parameters:
 71: + disc   - The disc context to distroy
 72: - viewer - The viewer

 74:   Level: beginner

 76: .keywords: discretization, view
 77: .seealso: DiscretizationDestroy()
 78: @*/
 79: int DiscretizationView(Discretization disc, PetscViewer viewer)
 80: {

 85:   if (viewer == PETSC_NULL) {
 86:     viewer = PETSC_VIEWER_STDOUT_SELF;
 87:   } else {
 89:   }
 90:   (*disc->ops->view)(disc, viewer);
 91:   return(0);
 92: }

 94: EXTERN_C_BEGIN
 95: #undef  __FUNCT__
 97: int DiscretizationSerialize_Generic(MPI_Comm comm, Discretization *disc, PetscViewer viewer, PetscTruth store) {
 98:   Discretization d;
 99:   int            fd;
100:   int            len;
101:   char          *type_name = PETSC_NULL;
102:   PetscTruth     match;
103:   int            ierr;


109:   PetscTypeCompare((PetscObject) viewer, PETSC_VIEWER_BINARY, &match);
110:   if (match == PETSC_FALSE) SETERRQ(PETSC_ERR_ARG_WRONG, "Must be binary viewer");
111:   PetscViewerBinaryGetDescriptor(viewer, &fd);

113:   if (store) {
114:     d = *disc;
115:     PetscStrlen(d->type_name, &len);
116:     PetscBinaryWrite(fd, &len,          1,   PETSC_INT,  0);
117:     PetscBinaryWrite(fd,  d->type_name, len, PETSC_CHAR, 0);
118:     PetscBinaryWrite(fd, &d->comp,      1,   PETSC_INT,  0);
119:     PetscBinaryWrite(fd, &d->field,     1,   PETSC_INT,  0);
120:   } else {
121:     DiscretizationCreate(comm, &d);
122:     PetscBinaryRead(fd, &len,       1,   PETSC_INT);
123:     if (len) {
124:       PetscMalloc((len+1) * sizeof(char), &type_name);
125:       type_name[len] = 0;
126:     }
127:     PetscBinaryRead(fd,  type_name, len, PETSC_CHAR);
128:     PetscBinaryRead(fd, &d->comp,   1,   PETSC_INT);
129:     PetscBinaryRead(fd, &d->field,  1,   PETSC_INT);
130:     DiscretizationSetType(d, type_name);
131:     if (len) {
132:       PetscFree(type_name);
133:     }
134:     *disc = d;
135:   }

137:   return(0);
138: }
139: EXTERN_C_END

141: #undef __FUNCT__  
143: /*
144:   DiscretizationSetTypeFromOptions - Sets the type of discretization from user options.

146:   Collective on Discretization

148:   Input Parameter:
149: . disc - The discretization

151:   Level: intermediate

153: .keywords: Discretization, set, options, database, type
154: .seealso: DiscretizationSetFromOptions(), DiscretizationSetType()
155: */
156: static int DiscretizationSetTypeFromOptions(Discretization disc) {
157:   PetscTruth opt;
158:   char      *defaultType;
159:   char       typeName[256];
160:   int        ierr;

163:   if (disc->type_name != PETSC_NULL) {
164:     defaultType = disc->type_name;
165:   } else {
166:     defaultType = DISCRETIZATION_TRIANGULAR_2D_LINEAR;
167:   }

169:   if (DiscretizationRegisterAllCalled == PETSC_FALSE) {
170:     DiscretizationRegisterAll(PETSC_NULL);
171:   }
172:   PetscOptionsList("-disc_type", "Discretization method"," DiscretizationSetType", DiscretizationList, defaultType, typeName, 256, &opt);
173:   if (opt == PETSC_TRUE) {
174:     DiscretizationSetType(disc, typeName);
175:   } else {
176:     DiscretizationSetType(disc, defaultType);
177:   }
178:   return(0);
179: }

181: #undef __FUNCT__  
183: /*
184:   DiscretizationSetSerializeTypeFromOptions - Sets the type of discretization serialization from user options.

186:   Collective on Discretization

188:   Input Parameter:
189: . disc - The discretization

191:   Level: intermediate

193: .keywords: Discretization, set, options, database, type
194: .seealso: DiscretizationSetFromOptions(), DiscretizationSetSerializeType()
195: */
196: static int DiscretizationSetSerializeTypeFromOptions(Discretization disc) {
197:   PetscTruth opt;
198:   char      *defaultType;
199:   char       typeName[256];
200:   int        ierr;

203:   if (disc->serialize_name != PETSC_NULL) {
204:     defaultType = disc->serialize_name;
205:   } else {
206:     defaultType = DISCRETIZATION_SER_GENERIC;
207:   }

209:   if (DiscretizationSerializeRegisterAllCalled == PETSC_FALSE) {
210:     DiscretizationSerializeRegisterAll(PETSC_NULL);
211:   }
212:   PetscOptionsList("-disc_serialize_type", "Discretization serialization method"," DiscretizationSetSerializeType",
213:                           DiscretizationSerializeList, defaultType, typeName, 256, &opt);
214:   if (opt == PETSC_TRUE) {
215:     DiscretizationSetSerializeType(disc, typeName);
216:   } else {
217:     DiscretizationSetSerializeType(disc, defaultType);
218:   }
219:   return(0);
220: }

222: #undef __FUNCT__  
224: /*@
225:   DiscretizationSetFromOptions - Sets various Discretization parameters from user options.

227:   Collective on Discretization

229:   Input Parameter:
230: . disc - The discretization

232:   Notes:  To see all options, run your program with the -help option, or consult the users manual.
233:           Must be called after DiscretizationCreate() but before the Discretization is used.

235:   Level: intermediate

237: .keywords: Discretization, set, options, database
238: .seealso: DiscretizationCreate(), DiscretizationPrintHelp(), DiscretizationSetOptionsPrefix()
239: @*/
240: int DiscretizationSetFromOptions(Discretization disc) {
241:   Discretization bdDisc;
242:   PetscTruth     opt;
243:   int            ierr;

247:   PetscOptionsBegin(disc->comm, disc->prefix, "Discretization options", "Discretization");

249:   /* Handle generic disc options */
250:   PetscOptionsHasName(PETSC_NULL, "-help", &opt);
251:   if (opt == PETSC_TRUE) {
252:     DiscretizationPrintHelp(disc);
253:   }

255:   /* Handle disc type options */
256:   DiscretizationSetTypeFromOptions(disc);

258:   /* Handle disc serialization options */
259:   DiscretizationSetSerializeTypeFromOptions(disc);

261:   /* Handle specific disc options */
262:   if (disc->ops->setfromoptions != PETSC_NULL) {
263:     (*disc->ops->setfromoptions)(disc);
264:   }

266:   PetscOptionsEnd();

268:   /* Handle subobject options */
269:   DiscretizationGetBoundaryDiscretization(disc, &bdDisc);
270:   if (bdDisc != PETSC_NULL) {
271:     DiscretizationSetFromOptions(bdDisc);
272:   }

274:   DiscretizationViewFromOptions(disc, disc->name);
275:   return(0);
276: }

278: #undef  __FUNCT__
280: /*@
281:   DiscretizationViewFromOptions - This function visualizes the discretization based upon user options.

283:   Collective on Discretization

285:   Input Parameter:
286: . disc - The discretization

288:   Level: intermediate

290: .keywords: Discretization, view, options, database
291: .seealso: DiscretizationSetFromOptions(), DiscretizationView()
292: @*/
293: int DiscretizationViewFromOptions(Discretization disc, char *title) {
294:   PetscViewer viewer;
295:   PetscDraw   draw;
296:   PetscTruth  opt;
297:   char       *titleStr;
298:   char        typeName[1024];
299:   char        fileName[PETSC_MAX_PATH_LEN];
300:   int         len;
301:   int         ierr;

304:   PetscOptionsHasName(disc->prefix, "-disc_view", &opt);
305:   if (opt == PETSC_TRUE) {
306:     PetscOptionsGetString(disc->prefix, "-disc_view", typeName, 1024, &opt);
307:     PetscStrlen(typeName, &len);
308:     if (len > 0) {
309:       PetscViewerCreate(disc->comm, &viewer);
310:       PetscViewerSetType(viewer, typeName);
311:       PetscOptionsGetString(disc->prefix, "-disc_view_file", fileName, 1024, &opt);
312:       if (opt == PETSC_TRUE) {
313:         PetscViewerSetFilename(viewer, fileName);
314:       } else {
315:         PetscViewerSetFilename(viewer, disc->name);
316:       }
317:       DiscretizationView(disc, viewer);
318:       PetscViewerFlush(viewer);
319:       PetscViewerDestroy(viewer);
320:     } else {
321:       DiscretizationView(disc, PETSC_NULL);
322:     }
323:   }
324:   PetscOptionsHasName(disc->prefix, "-disc_view_draw", &opt);
325:   if (opt == PETSC_TRUE) {
326:     PetscViewerDrawOpen(disc->comm, 0, 0, 0, 0, 300, 300, &viewer);
327:     PetscViewerDrawGetDraw(viewer, 0, &draw);
328:     if (title != PETSC_NULL) {
329:       titleStr = title;
330:     } else {
331:       PetscObjectName((PetscObject) disc);                                                         CHKERRQ(ierr) ;
332:       titleStr = disc->name;
333:     }
334:     PetscDrawSetTitle(draw, titleStr);
335:     DiscretizationView(disc, viewer);
336:     PetscViewerFlush(viewer);
337:     PetscDrawPause(draw);
338:     PetscViewerDestroy(viewer);
339:   }
340:   return(0);
341: }

343: #undef __FUNCT__  
345: /*@
346:   DiscretizationPrintHelp - Prints all options for the discretization.

348:   Input Parameter:
349: . disc - The discretization

351:   Options Database Keys:
352: $  -help, -h

354:   Level: intermediate

356: .keywords: Discretization, help
357: .seealso: DiscretizationSetFromOptions()
358: @*/
359: int DiscretizationPrintHelp(Discretization disc) {
360:   char p[64];
361:   int  ierr;


366:   PetscStrcpy(p, "-");
367:   if (disc->prefix != PETSC_NULL) {
368:     PetscStrcat(p, disc->prefix);
369:   }

371:   (*PetscHelpPrintf)(disc->comm, "Discretization options ------------------------------------------------n");
372:   (*PetscHelpPrintf)(disc->comm,"   %sdisc_type <typename> : Sets the discretization typen", p);
373:   return(0);
374: }

376: #undef __FUNCT__  
378: /*@C
379:   DiscretizationSetOptionsPrefix - Sets the prefix used for searching for all discretization options in the database.

381:   Input Parameters:
382: + disc   - The discretization
383: - prefix - The prefix to prepend to all option names

385:   Notes:
386:   A hyphen (-) must NOT be given at the beginning of the prefix name.
387:   The first character of all runtime options is AUTOMATICALLY the hyphen.

389:   Level: intermediate

391: .keywords: Discretization, set, options, prefix, database
392: .seealso: DiscretizationSetFromOptions()
393: @*/
394: int DiscretizationSetOptionsPrefix(Discretization disc, char *prefix) {

399:   PetscObjectSetOptionsPrefix((PetscObject) disc, prefix);
400:   return(0);
401: }

403: #undef __FUNCT__  
405: /*@C
406:   DiscretizationAppendOptionsPrefix - Appends to the prefix used for searching for all discretization options in the database.

408:   Collective on Discretization

410:   Input Parameters:
411: + disc   - The discretization
412: - prefix - The prefix to prepend to all option names

414:   Notes:
415:   A hyphen (-) must NOT be given at the beginning of the prefix name.
416:   The first character of all runtime options is AUTOMATICALLY the hyphen.

418:   Level: intermediate

420: .keywords: Discretization, append, options, prefix, database
421: .seealso: DiscretizationGetOptionsPrefix()
422: @*/
423: int DiscretizationAppendOptionsPrefix(Discretization disc, char *prefix) {
425: 
428:   PetscObjectAppendOptionsPrefix((PetscObject) disc, prefix);
429:   return(0);
430: }

432: #undef __FUNCT__  
434: /*@C
435:   DiscretizationGetOptionsPrefix - Sets the prefix used for searching for all discretization options in the database.

437:   Input Parameter:
438: . disc   - The discretization

440:   Output Parameter:
441: . prefix - A pointer to the prefix string used

443:   Level: intermediate

445: .keywords: Discretization, get, options, prefix, database
446: .seealso: DiscretizationAppendOptionsPrefix()
447: @*/
448: int DiscretizationGetOptionsPrefix(Discretization disc, char **prefix) {

453:   PetscObjectGetOptionsPrefix((PetscObject) disc, prefix);
454:   return(0);
455: }

457: #undef __FUNCT__  
459: /*@
460:   DiscretizationSetupDefaultOperators - Adds the default set of operators to the discretization if they are
461:   not already present.

463:   Input Parameter:
464: . disc   - The discretization

466:   Level: intermediate

468: .keywords: Discretization, default, operator
469: .seealso: DiscretizationSetup(), DiscretizationSerialize()
470: @*/
471: int DiscretizationSetupDefaultOperators(Discretization disc) {

476:   if (disc->numOps > 0) return(0);
477:   (*disc->ops->setupoperators)(disc);
478:   return(0);
479: }

481: /*-------------------------------------------------- Query Functions -------------------------------------------------*/
482: #undef  __FUNCT__
484: /*@C
485:   DiscretizationSetNumComponents - This function sets the number of components in the field being discretized.

487:   Not collective

489:   Input Parameters:
490: . disc    - The discretization
491: . numComp - The number of components in the field

493:   Level: intermediate

495: .keywords discretization, component, field
496: .seealso DiscretizationCreate(), DiscretizationRegisterOperator()
497: @*/
498: int DiscretizationSetNumComponents(Discretization disc, int numComp) {
501:   if (numComp <= 0) SETERRQ(PETSC_ERR_ARG_WRONG, "Discretization msut have at least one component");
502:   disc->comp = numComp;
503:   return(0);
504: }

506: #undef  __FUNCT__
508: /*@C
509:   DiscretizationGetNumComponents - This function returns the number of components in the field being discretized.

511:   Not collective

513:   Input Parameter:
514: . disc    - The discretization

516:   Output Parameter:
517: . numComp - The number of components in the field

519:   Level: intermediate

521: .keywords discretization, component, field
522: .seealso DiscretizationCreate(), DiscretizationRegisterOperator()
523: @*/
524: int DiscretizationGetNumComponents(Discretization disc, int *numComp) {
528:   *numComp = disc->comp;
529:   return(0);
530: }

532: #undef  __FUNCT__
534: /*@C
535:   DiscretizationGetNumFunctions - This function returns the number of functions per element in the Discretization.

537:   Not collective

539:   Input Parameter:
540: . disc    - The discretization

542:   Output Parameter:
543: . numFunc - The number of functions per element

545:   Level: intermediate

547: .keywords discretization, functions
548: .seealso DiscretizationCreate(), DiscretizationRegisterOperator()
549: @*/
550: int DiscretizationGetNumFunctions(Discretization disc, int *numFunc) {
554:   *numFunc = disc->funcs;
555:   return(0);
556: }

558: #undef  __FUNCT__
560: /*@C
561:   DiscretizationGetBoundaryDiscretization - This function returns boundary discretization, or PETSC_NULL if none exists.

563:   Not collective

565:   Input Parameter:
566: . disc   - The discretization

568:   Output Parameter:
569: . bdDisc - The boundary discretization, or PETSC_NULL

571:   Level: intermediate

573: .keywords discretization, functions
574: .seealso DiscretizationCreate(), DiscretizationRegisterOperator()
575: @*/
576: int DiscretizationGetBoundaryDiscretization(Discretization disc, Discretization *bdDisc) {
580:   *bdDisc = disc->bdDisc;
581:   return(0);
582: }

584: #undef  __FUNCT__
586: /*@C
587:   DiscretizationSetField - This function sets the corresponding field in the Grid for this Discretization.

589:   Not collective

591:   Input Parameters:
592: . disc  - The discretization
593: . field - The corresponding field in the Grid

595:   Level: intermediate

597: .keywords Discretization, component, field
598: .seealso DiscretizationGetField(), DiscretizationCreate(), DiscretizationRegisterOperator()
599: @*/
600: int DiscretizationSetField(Discretization disc, int field) {
603:   /* Could possibly validate with an associated Grid */
604:   disc->field = field;
605:   return(0);
606: }

608: #undef  __FUNCT__
610: /*@C
611:   DiscretizationGetField - This function gets the corresponding field in the Grid for this Discretization.

613:   Not collective

615:   Input Parameter:
616: . disc  - The discretization

618:   Output Parameter:
619: . field - The corresponding field in the Grid

621:   Level: intermediate

623: .keywords discretization, component, field
624: .seealso DiscretizationSetField(), DiscretizationCreate(), DiscretizationRegisterOperator()
625: @*/
626: int DiscretizationGetField(Discretization disc, int *field) {
630:   *field = disc->field;
631:   return(0);
632: }

634: #undef  __FUNCT__
636: /*@C
637:   DiscretizationGetFieldValues - This function gets the values for the corresponding field in the Grid for this Discretization.

639:   Not collective

641:   Input Parameter:
642: . disc  - The discretization

644:   Output Parameter:
645: . fieldValues - The values for the corresponding field in the Grid

647:   Level: intermediate

649: .keywords discretization, component, field
650: .seealso DiscretizationSetField(), DiscretizationCreate(), DiscretizationRegisterOperator()
651: @*/
652: int DiscretizationGetFieldValues(Discretization disc, PetscScalar **fieldValues) {
656:   *fieldValues = disc->fieldVal[0];
657:   return(0);
658: }

660: /*--------------------------------------------- Quadrature Query Functions -------------------------------------------*/
661: #undef  __FUNCT__
663: /*@C
664:   DiscretizationGetNumQuadraturePoints - This function returns the number of points used in the element quadrature.

666:   Not collective

668:   Input Parameter:
669: . disc      - The discretization

671:   Output Parameter:
672: . numPoints - The number of quadrature points

674:   Level: intermediate

676: .keywords discretization, quadrature, points
677: .seealso DiscretizationCreate(), DiscretizationRegisterOperator()
678: @*/
679: int DiscretizationGetNumQuadraturePoints(Discretization disc, int *numPoints) {
683:   *numPoints = disc->numQuadPoints;
684:   return(0);
685: }

687: #undef  __FUNCT__
689: /*@C
690:   DiscretizationGetQuadraturePoints - This function returns the coordinates of each quadrature point.

692:   Not collective

694:   Input Parameter:
695: . disc   - The discretization

697:   Output Parameter:
698: . coords - The quadrature point coordinates

700:   Level: intermediate

702: .keywords discretization, quadrature, coordinates, points
703: .seealso DiscretizationCreate(), DiscretizationRegisterOperator()
704: @*/
705: int DiscretizationGetQuadraturePoints(Discretization disc, double **coords) {
709:   *coords = disc->quadPoints;
710:   return(0);
711: }

713: #undef  __FUNCT__
715: /*@C
716:   DiscretizationGetQuadratureWeights - This function returns the weight associated with each quadrature point.

718:   Not collective

720:   Input Parameter:
721: . disc    - The discretization

723:   Output Parameter:
724: . weights - The quadrature weights

726:   Level: intermediate

728: .keywords discretization, quadrature, weights
729: .seealso DiscretizationCreate(), DiscretizationRegisterOperator()
730: @*/
731: int DiscretizationGetQuadratureWeights(Discretization disc, double **weights) {
735:   *weights = disc->quadWeights;
736:   return(0);
737: }

739: #undef  __FUNCT__
741: /*@C
742:   DiscretizationGetQuadratureFunctions - This function returns each shape function evaluated at the quadrature points.

744:   Not collective

746:   Input Parameter:
747: . disc  - The discretization

749:   Output Parameter:
750: . funcs - The shape functions evaluated at the quadrature points

752:   Level: intermediate

754: .keywords discretization, quadrature, functions
755: .seealso DiscretizationCreate(), DiscretizationRegisterOperator()
756: @*/
757: int DiscretizationGetQuadratureFunctions(Discretization disc, double **funcs) {
761:   *funcs = disc->quadShapeFuncs;
762:   return(0);
763: }

765: #undef  __FUNCT__
767: /*@C
768:   DiscretizationGetQuadratureDerivatives - This function returns the derivatives of each shape function evaluated
769:   at the quadrature points.

771:   Not collective

773:   Input Parameter:
774: . disc - The discretization

776:   Output Parameter:
777: . ders - The shape function derivatives evaluated at the quadrature points

779:   Level: intermediate

781: .keywords discretization, quadrature, derivatives
782: .seealso DiscretizationCreate(), DiscretizationRegisterOperator()
783: @*/
784: int DiscretizationGetQuadratureDerivatives(Discretization disc, double **ders) {
788:   *ders = disc->quadShapeFuncDers;
789:   return(0);
790: }

792: /*----------------------------------------------- Registration Functions ---------------------------------------------*/
793: #undef  __FUNCT__
795: static int DiscretizationRegisterOperator_Private(Discretization disc, PetscScalar *precompInt, OperatorFunction func,
796:                                                   ALEOperatorFunction ALEFunc, int *newOp)
797: {
798:   Operator *tempOperators;
799:   Operator  op;
800:   int       ierr;

803:   while (disc->numOps >= disc->maxOps) {
804:     PetscMalloc(disc->maxOps*2 * sizeof(Operator), &tempOperators);
805:     PetscLogObjectMemory(disc, disc->maxOps * sizeof(Operator));
806:     PetscMemcpy(tempOperators, disc->operators, disc->maxOps * sizeof(Operator));
807:     PetscFree(disc->operators);
808:     disc->operators = tempOperators;
809:     disc->maxOps   *= 2;
810:   }
811:   /* Create the new operator */
812:   PetscNew(struct _Operator, &op);
813:   PetscLogObjectMemory(disc, sizeof(struct _Operator));
814:   op->precompInt = precompInt;
815:   op->test       = PETSC_NULL;
816:   op->opFunc     = func;
817:   op->ALEOpFunc  = ALEFunc;
818:   disc->operators[disc->numOps] = op;
819:   *newOp = disc->numOps;
820:   disc->numOps++;
821:   return(0);
822: }

824: #undef  __FUNCT__
826: /*@C
827:   DiscretizationRegisterOperator - This function defines a new operator

829:   Collective on Discretization

831:   Input Parameters:
832: + disc  - The discretization
833: - func  - The function defining the operator

835:   Output Parameter:
836: . newOp - The index of the new operator

838:   Level: advanced

840: .keywords Discretization, operator, register
841: .seealso DiscretizationRegisterALEOperator(), GridAddMatOperator(), GridAddRhsOperator()
842: @*/
843: int DiscretizationRegisterOperator(Discretization disc, OperatorFunction func, int *newOp) {

849:   DiscretizationRegisterOperator_Private(disc, PETSC_NULL, func, PETSC_NULL, newOp);
850:   return(0);
851: }

853: #undef  __FUNCT__
855: /*@C
856:   DiscretizationRegisterALEOperator - This function defines a new ALE operator

858:   Collective on Discretization

860:   Input Parameters:
861: + disc  - The discretization
862: - func  - The function deinfing the operator

864:   Output Parameter:
865: . newOp - The index of the new operator

867:   Level: advanced

869: .keywords Discretization, operator, register
870: .seealso DiscretizationRegisterOperator(), GridAddMatOperator(), GridAddRhsOperator()
871: @*/
872: int DiscretizationRegisterALEOperator(Discretization disc, ALEOperatorFunction func, int *newOp) {

878:   DiscretizationRegisterOperator_Private(disc, PETSC_NULL, PETSC_NULL, func, newOp);
879:   return(0);
880: }

882: #undef  __FUNCT__
884: /*@C
885:   DiscretizationRegisterPrecomputedOperator - This function defines a new operator for which the integral can be precomputed .
886:   This requires knowledge of the discretization itself and the storage format, and thus is only useful for developers.

888:   Collective on Discretization

890:   Input Parameters:
891: + disc       - The discretization
892: - precompInt - The array of precomputed values

894:   Output Parameter:
895: . newOp - The index of the new operator

897:   Level: developer

899: .keywords Discretization, operator, register
900: .seealso DiscretizationRegisterALEOperator(), GridAddMatOperator(), GridAddRhsOperator()
901: @*/
902: int DiscretizationRegisterPrecomputedOperator(Discretization disc, PetscScalar *precompInt, int *newOp) {

909:   DiscretizationRegisterOperator_Private(disc, precompInt, PETSC_NULL, PETSC_NULL, newOp);
910:   return(0);
911: }

913: /*------------------------------------------------ Evaluation Functions ----------------------------------------------*/
914: #undef  __FUNCT__
916: /*@
917:   DiscretizationEvaluateFunctionGalerkin - This function calculates the
918:   weak form of a function over a single element.

920:   Collective on Discretization

922:   Input Parameters:
923: + disc  - The discretization
924: . mesh  - The associated mesh
925: . f     - The PointFunction of which we want the weak form
926: . alpha - A scalar multiplier
927: . elem  - The local element number
928: - ctx   - The user-supplied context

930:   Output Parameter:
931: . v     - The element vector for the given element

933:   Level: beginner

935: .keywords: discretization, function, weak form
936: .seealso: DiscretizationEvaluateOperatorGalerkin(), DiscretizationEvaluateNonlinearOperatorGalerkin()
937: @*/
938: int DiscretizationEvaluateFunctionGalerkin(Discretization disc, Mesh mesh, PointFunction f, PetscScalar alpha, int elem, PetscScalar *v, void *ctx)
939: {

946:   (*(disc->ops->evaluatefunctiongalerkin))(disc, mesh, f, alpha, elem, v, ctx);
947:   return(0);
948: }

950: #undef  __FUNCT__
952: /*@
953:   DiscretizationEvaluateOperatorGalerkin - This function calculates the
954:   weak form of an operator over a single element.

956:   Collective on Discretization

958:   Input Parameters:
959: + disc     - The discretization for the basis functions
960: . mesh     - The associated mesh
961: . elemSize - The size of the element matrix
962: . rowStart - The starting row index in the element matrix
963: . colStart - The starting column index in the element matrix
964: . op       - The operator index (of registered operators)
965: . alpha    - The scalar multiple of the operator
966: . elem     - The local element number
967: . field    - The field values
968: - ctx      - The user-supplied context

970:   Output Parameter:
971: . mat      - The element matrix for the given element

973:   Level: beginner

975: .keywords: discretization, operator, weak form
976: .seealso: DiscretizationEvaluateFunctionGalerkin(), DiscretizationEvaluateALEOperatorGalerkin(), DiscretizationEvaluateOperatorGalerkinMF(),
977:           DiscretizationEvaluateNonlinearOperatorGalerkin()
978: @*/
979: int DiscretizationEvaluateOperatorGalerkin(Discretization disc, Mesh mesh, int elemSize, int rowStart, int colStart,
980:                                            int op, PetscScalar alpha, int elem, PetscScalar *field, PetscScalar *mat, void *ctx)
981: {

987:   /* field may be PETSC_NULL */
989:   (*(disc->ops->evaluateoperatorgalerkin))(disc, mesh, elemSize, rowStart, colStart, op, alpha, elem,
990:                                                       field, mat, ctx);
991: 
992:   return(0);
993: }

995: #undef  __FUNCT__
997: /*@
998:   DiscretizationEvaluateALEOperatorGalerkin - This function calculates the
999:   weak form of an operator over a single element.

1001:   Collective on Discretization

1003:   Input Parameters:
1004: + disc     - The discretization for the basis functions
1005: . mesh     - The associated mesh
1006: . elemSize - The size of the element matrix
1007: . rowStart - The starting row index in the element matrix
1008: . colStart - The starting column index in the element matrix
1009: . op       - The operator index (of registered operators)
1010: . alpha    - The scalar multiple of the operator
1011: . elem     - The local element number
1012: . field    - The field values
1013: . ALEfield - The mesh velocity field values
1014: - ctx      - The user-supplied context

1016:   Output Parameter:
1017: . mat      - The element matrix for the given element

1019:   Level: beginner

1021: .keywords: discretization, operator, ALE, weak form
1022: .seealso: DiscretizationEvaluateOperatorGalerkin()
1023: @*/
1024: int DiscretizationEvaluateALEOperatorGalerkin(Discretization disc, Mesh mesh, int elemSize, int rowStart, int colStart,
1025:                                               int op, PetscScalar alpha, int elem, PetscScalar *field, PetscScalar *ALEfield, PetscScalar *mat, void *ctx)
1026: {

1032:   /* field may be PETSC_NULL */
1035:   (*(disc->ops->evaluatealeoperatorgalerkin))(disc, mesh, elemSize, rowStart, colStart, op, alpha, elem,
1036:                                                          field, ALEfield, mat, ctx);
1037: 
1038:   return(0);
1039: }

1041: #undef  __FUNCT__
1043: /*@
1044:   DiscretizationEvaluateOperatorGalerkinMF - This function calculates the
1045:   weak form of the action of an operator over a single element.

1047:   Collective on Discretization

1049:   Input Parameters:
1050: + disc     - The discretization
1051: . mesh     - The associated mesh
1052: . elemSize - The size of the element matrix
1053: . rowStart - The starting row index in the element matrix
1054: . colStart - The starting column index in the element matrix
1055: . op       - The operator index (of registered operators)
1056: . alpha    - The scalar multiple of the operator
1057: . elem     - The local element number
1058: . field    - The element vector for the grid vector acted upon
1059: . app      - The element vector to which to apply the matrix, usually identical to field
1060: . mat      - A temporary element matrix for work space
1061: - ctx      - The user-supplied context

1063:   Output Parameter:
1064: . vec      - The element vector for the given element

1066:   Level: beginner

1068: .keywords: discretization, operator, MF, weak form
1069: .seealso: DiscretizationEvaluateOperatorGalerkin()
1070: @*/
1071: int DiscretizationEvaluateOperatorGalerkinMF(Discretization disc, Mesh mesh, int elemSize, int rowStart, int colStart, int op,
1072:                                              PetscScalar alpha, int elem, PetscScalar *field, PetscScalar *app, PetscScalar *vec, PetscScalar *mat, void *ctx)
1073: {
1074:   Discretization test;    /* The space of test functions */
1075:   int            rowSize; /* The number of shape functions per element */
1076:   int            colSize; /* The number of test  functions per element */
1077:   int            i, j;
1078:   int            ierr;

1081:   /* Get discretization info */
1084:   /* field may be PETSC_NULL */
1087:   test    = disc->operators[op]->test;
1089:   rowSize = test->size;
1090:   colSize = disc->size;

1092:   DiscretizationEvaluateOperatorGalerkin(disc, mesh, elemSize, rowStart, colStart, op, alpha, elem, field, mat, ctx);
1093: 
1094:   /* Perform the local element matvec */
1095:   for(i = 0; i < rowSize; i++)
1096:     for(j = 0; j < colSize; j++)
1097:       vec[i+rowStart] += mat[(i+rowStart)*elemSize+j+colStart]*app[j];
1098:   PetscLogFlops(2*rowSize*colSize);
1099:   return(0);
1100: }

1102: #undef  __FUNCT__
1104: /*@
1105:   DiscretizationEvaluateALEOperatorGalerkinMF - This function calculates the
1106:   weak form of the action of an operator over a single element.

1108:   Collective on Discretization

1110:   Input Parameters:
1111: + disc     - The discretization
1112: . mesh     - The associated mesh
1113: . elemSize - The size of the element matrix
1114: . rowStart - The starting row index in the element matrix
1115: . colStart - The starting column index in the element matrix
1116: . op       - The operator index (of registered operators)
1117: . alpha    - The scalar multiple of the operator
1118: . elem     - The local element number
1119: . field    - The element vector for the grid vector acted upon
1120: . app      - The element vector to which to apply the matrix, usually identical to field
1121: . ALEfield - The element vector for the ALE velocity field
1122: . mat      - A temporary element matrix for work space
1123: - ctx      - The user-supplied context

1125:   Output Parameter:
1126: . vec      - The element vector for the given element

1128:   Level: beginner

1130: .keywords: discretization, operator, ALE, MF, weak form
1131: .seealso: DiscretizationEvaluateOperatorGalerkin()
1132: @*/
1133: int DiscretizationEvaluateALEOperatorGalerkinMF(Discretization disc, Mesh mesh, int elemSize, int rowStart, int colStart, int op,
1134:                                                 PetscScalar alpha, int elem, PetscScalar *field, PetscScalar *app, PetscScalar *ALEfield, PetscScalar *vec,
1135:                                                 PetscScalar *mat, void *ctx)
1136: {
1137:   Discretization test;    /* The space of test functions */
1138:   int            rowSize; /* The number of shape functions per element */
1139:   int            colSize; /* The number of test  functions per element */
1140:   int            i, j;
1141:   int            ierr;

1144:   /* Get discretization info */
1147:   /* field may be PETSC_NULL */
1151:   test    = disc->operators[op]->test;
1153:   rowSize = test->size;
1154:   colSize = disc->size;

1156:   DiscretizationEvaluateALEOperatorGalerkin(disc, mesh, elemSize, rowStart, colStart, op, alpha, elem, field, ALEfield, mat, ctx);
1157: 
1158:   /* Perform the local element matvec */
1159:   for(i = 0; i < rowSize; i++)
1160:     for(j = 0; j < colSize; j++)
1161:       vec[i+rowStart] += mat[(i+rowStart)*elemSize+j+colStart]*app[j];
1162:   PetscLogFlops(2*rowSize*colSize);
1163:   return(0);
1164: }

1166: #undef  __FUNCT__
1168: /*@
1169:   DiscretizationEvaluateNonlinearOperatorGalerkin - This function calculates the
1170:   weak form of a nonlinear operator over a single element.

1172:   Collective on Discretization

1174:   Input Parameters:
1175: + disc    - The discretization
1176: . mesh    - The associated mesh
1177: . op      - The function defining the nonlinear operator
1178: . alpha   - A scalar multiplier
1179: . elem    - The local element number
1180: . numArgs - THe number of input element vectors
1181: . field   - The element vectors for the grid vectors acted upon
1182: - ctx     - The user-supplied context

1184:   Output Parameter:
1185: . vec     - The global element vector for the given element

1187:   Level: beginner

1189: .keywords: discretization, operator, nonlinear, weak form
1190: .seealso: DiscretizationEvaluateFunctionGalerkin(), DiscretizationEvaluateOperatorGalerkin(), DiscretizationEvaluateNonlinearALEOperatorGalerkin()
1191: @*/
1192: int DiscretizationEvaluateNonlinearOperatorGalerkin(Discretization disc, Mesh mesh, NonlinearOperator op, PetscScalar alpha, int elem,
1193:                                                     int numArgs, PetscScalar **field, PetscScalar *vec, void *ctx)
1194: {

1202:   (*(disc->ops->evaluatenonlinearoperatorgalerkin))(disc, mesh, op, alpha, elem, numArgs, field, vec, ctx);
1203: 
1204:   return(0);
1205: }

1207: #undef  __FUNCT__
1209: /*@
1210:   DiscretizationEvaluateNonlinearALEOperatorGalerkin - This function calculates the
1211:   weak form of a nonlinear operator over a single element.

1213:   Input Parameters:
1214: + disc     - The discretization
1215: . mesh     - The associated mesh
1216: . op       - The function defining the nonlinear operator
1217: . alpha    - A scalar multiplier
1218: . elem     - The local element number
1219: . numArgs  - THe number of input element vectors
1220: . field    - The element vectors for the grid vectors acted upon
1221: . ALEfield - The element vector for the mesh velocity field
1222: - ctx      - The user-supplied context

1224:   Output Parameter:
1225: . vec      - The global element vector for the given element

1227:   Level: beginner

1229: .keywords: discretization, operator, nonlinear, ALE, weak form
1230: .seealso: DiscretizationEvaluateFunctionGalerkin(), DiscretizationEvaluateOperatorGalerkin(), DiscretizationEvaluateNonlinearOperatorGalerkin()
1231: @*/
1232: int DiscretizationEvaluateNonlinearALEOperatorGalerkin(Discretization disc, Mesh mesh, NonlinearOperator op, PetscScalar alpha, int elem,
1233:                                                        int numArgs, PetscScalar **field, PetscScalar *ALEfield, PetscScalar *vec, void *ctx)
1234: {

1243:   (*(disc->ops->evaluatenonlinearaleoperatorgalerkin))(disc, mesh, op, alpha, elem, numArgs, field, ALEfield, vec, ctx);
1244: 
1245:   return(0);
1246: }

1248: /*---------------------------------------------- Interpolation Functions ---------------------------------------------*/
1249: #undef  __FUNCT__
1251: /*@
1252:   DiscretizationInterpolateField - This function interpolates a field at a given point.

1254:   Collective on Discretization

1256:   Input Parameters:
1257: + disc        - The discretization
1258: . mesh        - The associated mesh
1259: . elem        - The element containing the point
1260: . x,y,z       - The interpolation point
1261: . oldFieldVal - The element vector for 'elem'
1262: - type        - The method of interpolation, e.g. INTERPOLATION_LOCAL

1264:   Output Parameter:
1265: . newFieldVal - The field components at the given point

1267:   Level: beginner

1269: .keywords: discretization, interpolation
1270: .seealso: DiscretizationInterpolateElementVec()
1271: @*/
1272: int DiscretizationInterpolateField(Discretization disc, Mesh mesh, int elem, double x, double y, double z,
1273:                                    PetscScalar *oldFieldVal, PetscScalar *newFieldVal, InterpolationType type)
1274: {

1282:   (*disc->ops->interpolatefield)(disc, mesh, elem, x, y, z, oldFieldVal, newFieldVal, type);
1283:   return(0);
1284: }

1286: #undef  __FUNCT__
1288: /*@
1289:   DiscretizationInterpolateElementVec - Interpolates a given element vector from one discretization to another.

1291:   Input Parameters:
1292: + disc    - The original discretization
1293: . vec     - The original element vector
1294: - newDisc - The discretization defining the new element vector

1296:   Output Parameter:
1297: . newVec  - The interpolated element vector

1299:   Level: beginner

1301: .keywords: discretization, interpolation
1302: .seealso: DiscretizationInterpolateField()
1303: @*/
1304: int DiscretizationInterpolateElementVec(Discretization disc, ElementVec vec, Discretization newDisc, ElementVec newVec) {
1305:   int size, newSize;

1313:   if (disc->comp != newDisc->comp) SETERRQ(PETSC_ERR_ARG_INCOMP, "Fields must have the same number of components");
1314:   ElementVecGetSize(vec,    &size);
1315:   ElementVecGetSize(newVec, &newSize);
1316:   if ((disc->size > size) || (newDisc->size > newSize)) SETERRQ(PETSC_ERR_ARG_WRONG, "Element vector is too small to contain field");
1317:   (*disc->ops->interpolateelemvec)(disc, vec, newDisc, newVec);
1318:   return(0);
1319: }