Actual source code: space.c
1: #include <petsc/private/petscfeimpl.h>
2: #include <petscdmshell.h>
4: PetscClassId PETSCSPACE_CLASSID = 0;
6: PetscFunctionList PetscSpaceList = NULL;
7: PetscBool PetscSpaceRegisterAllCalled = PETSC_FALSE;
9: /*@C
10: PetscSpaceRegister - Adds a new `PetscSpace` implementation
12: Not Collective
14: Input Parameters:
15: + name - The name of a new user-defined creation routine
16: - create_func - The creation routine for the implementation type
18: Sample usage:
19: .vb
20: PetscSpaceRegister("my_space", MyPetscSpaceCreate);
21: .ve
23: Then, your PetscSpace type can be chosen with the procedural interface via
24: .vb
25: PetscSpaceCreate(MPI_Comm, PetscSpace *);
26: PetscSpaceSetType(PetscSpace, "my_space");
27: .ve
28: or at runtime via the option
29: .vb
30: -petscspace_type my_space
31: .ve
33: Level: advanced
35: Note:
36: `PetscSpaceRegister()` may be called multiple times to add several user-defined types of `PetscSpace`. The creation function is called
37: when the type is set to 'name'.
39: .seealso: `PetscSpace`, `PetscSpaceRegisterAll()`, `PetscSpaceRegisterDestroy()`
40: @*/
41: PetscErrorCode PetscSpaceRegister(const char sname[], PetscErrorCode (*function)(PetscSpace))
42: {
43: PetscFunctionBegin;
44: PetscCall(PetscFunctionListAdd(&PetscSpaceList, sname, function));
45: PetscFunctionReturn(PETSC_SUCCESS);
46: }
48: /*@C
49: PetscSpaceSetType - Builds a particular `PetscSpace`
51: Collective
53: Input Parameters:
54: + sp - The `PetscSpace` object
55: - name - The kind of space
57: Options Database Key:
58: . -petscspace_type <type> - Sets the `PetscSpace` type; use -help for a list of available types
60: Level: intermediate
62: .seealso: `PetscSpace`, `PetscSpaceType`, `PetscSpaceGetType()`, `PetscSpaceCreate()`
63: @*/
64: PetscErrorCode PetscSpaceSetType(PetscSpace sp, PetscSpaceType name)
65: {
66: PetscErrorCode (*r)(PetscSpace);
67: PetscBool match;
69: PetscFunctionBegin;
71: PetscCall(PetscObjectTypeCompare((PetscObject)sp, name, &match));
72: if (match) PetscFunctionReturn(PETSC_SUCCESS);
74: PetscCall(PetscSpaceRegisterAll());
75: PetscCall(PetscFunctionListFind(PetscSpaceList, name, &r));
76: PetscCheck(r, PetscObjectComm((PetscObject)sp), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscSpace type: %s", name);
78: PetscTryTypeMethod(sp, destroy);
79: sp->ops->destroy = NULL;
81: sp->dim = PETSC_DETERMINE;
82: PetscCall((*r)(sp));
83: PetscCall(PetscObjectChangeTypeName((PetscObject)sp, name));
84: PetscFunctionReturn(PETSC_SUCCESS);
85: }
87: /*@C
88: PetscSpaceGetType - Gets the `PetscSpaceType` (as a string) from the object.
90: Not Collective
92: Input Parameter:
93: . sp - The `PetscSpace`
95: Output Parameter:
96: . name - The `PetscSpace` type name
98: Level: intermediate
100: .seealso: `PetscSpaceType`, `PetscSpace`, `PetscSpaceSetType()`, `PetscSpaceCreate()`
101: @*/
102: PetscErrorCode PetscSpaceGetType(PetscSpace sp, PetscSpaceType *name)
103: {
104: PetscFunctionBegin;
107: if (!PetscSpaceRegisterAllCalled) PetscCall(PetscSpaceRegisterAll());
108: *name = ((PetscObject)sp)->type_name;
109: PetscFunctionReturn(PETSC_SUCCESS);
110: }
112: /*@C
113: PetscSpaceViewFromOptions - View a `PetscSpace` based on values in the options database
115: Collective
117: Input Parameters:
118: + A - the `PetscSpace` object
119: . obj - Optional object that provides the options name prefix
120: - name - command line option name
122: Level: intermediate
124: .seealso: `PetscSpace`, `PetscSpaceView()`, `PetscObjectViewFromOptions()`, `PetscSpaceCreate()`
125: @*/
126: PetscErrorCode PetscSpaceViewFromOptions(PetscSpace A, PetscObject obj, const char name[])
127: {
128: PetscFunctionBegin;
130: PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
131: PetscFunctionReturn(PETSC_SUCCESS);
132: }
134: /*@C
135: PetscSpaceView - Views a `PetscSpace`
137: Collective
139: Input Parameters:
140: + sp - the `PetscSpace` object to view
141: - v - the viewer
143: Level: beginner
145: .seealso: `PetscSpace`, `PetscViewer`, `PetscSpaceViewFromOptions()`, `PetscSpaceDestroy()`
146: @*/
147: PetscErrorCode PetscSpaceView(PetscSpace sp, PetscViewer v)
148: {
149: PetscInt pdim;
150: PetscBool iascii;
152: PetscFunctionBegin;
155: if (!v) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)sp), &v));
156: PetscCall(PetscSpaceGetDimension(sp, &pdim));
157: PetscCall(PetscObjectPrintClassNamePrefixType((PetscObject)sp, v));
158: PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &iascii));
159: PetscCall(PetscViewerASCIIPushTab(v));
160: if (iascii) PetscCall(PetscViewerASCIIPrintf(v, "Space in %" PetscInt_FMT " variables with %" PetscInt_FMT " components, size %" PetscInt_FMT "\n", sp->Nv, sp->Nc, pdim));
161: PetscTryTypeMethod(sp, view, v);
162: PetscCall(PetscViewerASCIIPopTab(v));
163: PetscFunctionReturn(PETSC_SUCCESS);
164: }
166: /*@
167: PetscSpaceSetFromOptions - sets parameters in a `PetscSpace` from the options database
169: Collective
171: Input Parameter:
172: . sp - the `PetscSpace` object to set options for
174: Options Database Keys:
175: + -petscspace_degree <deg> - the approximation order of the space
176: . -petscspace_variables <n> - the number of different variables, e.g. x and y
177: - -petscspace_components <c> - the number of components, say d for a vector field
179: Level: intermediate
181: .seealso: `PetscSpace`, `PetscSpaceView()`
182: @*/
183: PetscErrorCode PetscSpaceSetFromOptions(PetscSpace sp)
184: {
185: const char *defaultType;
186: char name[256];
187: PetscBool flg;
189: PetscFunctionBegin;
191: if (!((PetscObject)sp)->type_name) {
192: defaultType = PETSCSPACEPOLYNOMIAL;
193: } else {
194: defaultType = ((PetscObject)sp)->type_name;
195: }
196: if (!PetscSpaceRegisterAllCalled) PetscCall(PetscSpaceRegisterAll());
198: PetscObjectOptionsBegin((PetscObject)sp);
199: PetscCall(PetscOptionsFList("-petscspace_type", "Linear space", "PetscSpaceSetType", PetscSpaceList, defaultType, name, 256, &flg));
200: if (flg) {
201: PetscCall(PetscSpaceSetType(sp, name));
202: } else if (!((PetscObject)sp)->type_name) {
203: PetscCall(PetscSpaceSetType(sp, defaultType));
204: }
205: {
206: PetscCall(PetscOptionsDeprecated("-petscspace_order", "-petscspace_degree", "3.11", NULL));
207: PetscCall(PetscOptionsBoundedInt("-petscspace_order", "DEPRECATED: The approximation order", "PetscSpaceSetDegree", sp->degree, &sp->degree, NULL, 0));
208: }
209: PetscCall(PetscOptionsBoundedInt("-petscspace_degree", "The (maximally included) polynomial degree", "PetscSpaceSetDegree", sp->degree, &sp->degree, NULL, 0));
210: PetscCall(PetscOptionsBoundedInt("-petscspace_variables", "The number of different variables, e.g. x and y", "PetscSpaceSetNumVariables", sp->Nv, &sp->Nv, NULL, 0));
211: PetscCall(PetscOptionsBoundedInt("-petscspace_components", "The number of components", "PetscSpaceSetNumComponents", sp->Nc, &sp->Nc, NULL, 0));
212: PetscTryTypeMethod(sp, setfromoptions, PetscOptionsObject);
213: /* process any options handlers added with PetscObjectAddOptionsHandler() */
214: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)sp, PetscOptionsObject));
215: PetscOptionsEnd();
216: PetscCall(PetscSpaceViewFromOptions(sp, NULL, "-petscspace_view"));
217: PetscFunctionReturn(PETSC_SUCCESS);
218: }
220: /*@C
221: PetscSpaceSetUp - Construct data structures for the `PetscSpace`
223: Collective
225: Input Parameter:
226: . sp - the `PetscSpace` object to setup
228: Level: intermediate
230: .seealso: `PetscSpace`, `PetscSpaceView()`, `PetscSpaceDestroy()`
231: @*/
232: PetscErrorCode PetscSpaceSetUp(PetscSpace sp)
233: {
234: PetscFunctionBegin;
236: PetscTryTypeMethod(sp, setup);
237: PetscFunctionReturn(PETSC_SUCCESS);
238: }
240: /*@
241: PetscSpaceDestroy - Destroys a `PetscSpace` object
243: Collective
245: Input Parameter:
246: . sp - the `PetscSpace` object to destroy
248: Level: beginner
250: .seealso: `PetscSpace`, `PetscSpaceCreate()`
251: @*/
252: PetscErrorCode PetscSpaceDestroy(PetscSpace *sp)
253: {
254: PetscFunctionBegin;
255: if (!*sp) PetscFunctionReturn(PETSC_SUCCESS);
258: if (--((PetscObject)(*sp))->refct > 0) {
259: *sp = NULL;
260: PetscFunctionReturn(PETSC_SUCCESS);
261: }
262: ((PetscObject)(*sp))->refct = 0;
263: PetscCall(DMDestroy(&(*sp)->dm));
265: PetscCall((*(*sp)->ops->destroy)(*sp));
266: PetscCall(PetscHeaderDestroy(sp));
267: PetscFunctionReturn(PETSC_SUCCESS);
268: }
270: /*@
271: PetscSpaceCreate - Creates an empty `PetscSpace` object. The type can then be set with `PetscSpaceSetType()`.
273: Collective
275: Input Parameter:
276: . comm - The communicator for the `PetscSpace` object
278: Output Parameter:
279: . sp - The `PetscSpace` object
281: Level: beginner
283: .seealso: `PetscSpace`, `PetscSpaceSetType()`, `PETSCSPACEPOLYNOMIAL`
284: @*/
285: PetscErrorCode PetscSpaceCreate(MPI_Comm comm, PetscSpace *sp)
286: {
287: PetscSpace s;
289: PetscFunctionBegin;
291: PetscCall(PetscCitationsRegister(FECitation, &FEcite));
292: *sp = NULL;
293: PetscCall(PetscFEInitializePackage());
295: PetscCall(PetscHeaderCreate(s, PETSCSPACE_CLASSID, "PetscSpace", "Linear Space", "PetscSpace", comm, PetscSpaceDestroy, PetscSpaceView));
297: s->degree = 0;
298: s->maxDegree = PETSC_DETERMINE;
299: s->Nc = 1;
300: s->Nv = 0;
301: s->dim = PETSC_DETERMINE;
302: PetscCall(DMShellCreate(comm, &s->dm));
303: PetscCall(PetscSpaceSetType(s, PETSCSPACEPOLYNOMIAL));
305: *sp = s;
306: PetscFunctionReturn(PETSC_SUCCESS);
307: }
309: /*@
310: PetscSpaceGetDimension - Return the dimension of this space, i.e. the number of basis vectors
312: Input Parameter:
313: . sp - The `PetscSpace`
315: Output Parameter:
316: . dim - The dimension
318: Level: intermediate
320: .seealso: `PetscSpace`, `PetscSpaceGetDegree()`, `PetscSpaceCreate()`, `PetscSpace`
321: @*/
322: PetscErrorCode PetscSpaceGetDimension(PetscSpace sp, PetscInt *dim)
323: {
324: PetscFunctionBegin;
327: if (sp->dim == PETSC_DETERMINE) PetscTryTypeMethod(sp, getdimension, &sp->dim);
328: *dim = sp->dim;
329: PetscFunctionReturn(PETSC_SUCCESS);
330: }
332: /*@
333: PetscSpaceGetDegree - Return the polynomial degrees that characterize this space
335: Input Parameter:
336: . sp - The `PetscSpace`
338: Output Parameters:
339: + minDegree - The degree of the largest polynomial space contained in the space
340: - maxDegree - The degree of the smallest polynomial space containing the space
342: Level: intermediate
344: .seealso: `PetscSpace`, `PetscSpaceSetDegree()`, `PetscSpaceGetDimension()`, `PetscSpaceCreate()`, `PetscSpace`
345: @*/
346: PetscErrorCode PetscSpaceGetDegree(PetscSpace sp, PetscInt *minDegree, PetscInt *maxDegree)
347: {
348: PetscFunctionBegin;
352: if (minDegree) *minDegree = sp->degree;
353: if (maxDegree) *maxDegree = sp->maxDegree;
354: PetscFunctionReturn(PETSC_SUCCESS);
355: }
357: /*@
358: PetscSpaceSetDegree - Set the degree of approximation for this space.
360: Input Parameters:
361: + sp - The `PetscSpace`
362: . degree - The degree of the largest polynomial space contained in the space
363: - maxDegree - The degree of the largest polynomial space containing the space. One of degree and maxDegree can be `PETSC_DETERMINE`.
365: Level: intermediate
367: .seealso: `PetscSpace`, `PetscSpaceGetDegree()`, `PetscSpaceCreate()`, `PetscSpace`
368: @*/
369: PetscErrorCode PetscSpaceSetDegree(PetscSpace sp, PetscInt degree, PetscInt maxDegree)
370: {
371: PetscFunctionBegin;
373: sp->degree = degree;
374: sp->maxDegree = maxDegree;
375: PetscFunctionReturn(PETSC_SUCCESS);
376: }
378: /*@
379: PetscSpaceGetNumComponents - Return the number of components for this space
381: Input Parameter:
382: . sp - The `PetscSpace`
384: Output Parameter:
385: . Nc - The number of components
387: Level: intermediate
389: Note:
390: A vector space, for example, will have d components, where d is the spatial dimension
392: .seealso: `PetscSpace`, `PetscSpaceSetNumComponents()`, `PetscSpaceGetNumVariables()`, `PetscSpaceGetDimension()`, `PetscSpaceCreate()`, `PetscSpace`
393: @*/
394: PetscErrorCode PetscSpaceGetNumComponents(PetscSpace sp, PetscInt *Nc)
395: {
396: PetscFunctionBegin;
399: *Nc = sp->Nc;
400: PetscFunctionReturn(PETSC_SUCCESS);
401: }
403: /*@
404: PetscSpaceSetNumComponents - Set the number of components for this space
406: Input Parameters:
407: + sp - The `PetscSpace`
408: - order - The number of components
410: Level: intermediate
412: .seealso: `PetscSpace`, `PetscSpaceGetNumComponents()`, `PetscSpaceSetNumVariables()`, `PetscSpaceCreate()`, `PetscSpace`
413: @*/
414: PetscErrorCode PetscSpaceSetNumComponents(PetscSpace sp, PetscInt Nc)
415: {
416: PetscFunctionBegin;
418: sp->Nc = Nc;
419: PetscFunctionReturn(PETSC_SUCCESS);
420: }
422: /*@
423: PetscSpaceSetNumVariables - Set the number of variables for this space
425: Input Parameters:
426: + sp - The `PetscSpace`
427: - n - The number of variables, e.g. x, y, z...
429: Level: intermediate
431: .seealso: `PetscSpace`, `PetscSpaceGetNumVariables()`, `PetscSpaceSetNumComponents()`, `PetscSpaceCreate()`, `PetscSpace`
432: @*/
433: PetscErrorCode PetscSpaceSetNumVariables(PetscSpace sp, PetscInt n)
434: {
435: PetscFunctionBegin;
437: sp->Nv = n;
438: PetscFunctionReturn(PETSC_SUCCESS);
439: }
441: /*@
442: PetscSpaceGetNumVariables - Return the number of variables for this space
444: Input Parameter:
445: . sp - The `PetscSpace`
447: Output Parameter:
448: . Nc - The number of variables, e.g. x, y, z...
450: Level: intermediate
452: .seealso: `PetscSpace`, `PetscSpaceSetNumVariables()`, `PetscSpaceGetNumComponents()`, `PetscSpaceGetDimension()`, `PetscSpaceCreate()`, `PetscSpace`
453: @*/
454: PetscErrorCode PetscSpaceGetNumVariables(PetscSpace sp, PetscInt *n)
455: {
456: PetscFunctionBegin;
459: *n = sp->Nv;
460: PetscFunctionReturn(PETSC_SUCCESS);
461: }
463: /*@C
464: PetscSpaceEvaluate - Evaluate the basis functions and their derivatives (jet) at each point
466: Input Parameters:
467: + sp - The `PetscSpace`
468: . npoints - The number of evaluation points, in reference coordinates
469: - points - The point coordinates
471: Output Parameters:
472: + B - The function evaluations in a npoints x nfuncs array
473: . D - The derivative evaluations in a npoints x nfuncs x dim array
474: - H - The second derivative evaluations in a npoints x nfuncs x dim x dim array
476: Level: beginner
478: Note:
479: Above nfuncs is the dimension of the space, and dim is the spatial dimension. The coordinates are given
480: on the reference cell, not in real space.
482: .seealso: `PetscSpace`, `PetscFECreateTabulation()`, `PetscFEGetCellTabulation()`, `PetscSpaceCreate()`
483: @*/
484: PetscErrorCode PetscSpaceEvaluate(PetscSpace sp, PetscInt npoints, const PetscReal points[], PetscReal B[], PetscReal D[], PetscReal H[])
485: {
486: PetscFunctionBegin;
487: if (!npoints) PetscFunctionReturn(PETSC_SUCCESS);
493: PetscTryTypeMethod(sp, evaluate, npoints, points, B, D, H);
494: PetscFunctionReturn(PETSC_SUCCESS);
495: }
497: /*@
498: PetscSpaceGetHeightSubspace - Get the subset of the primal space basis that is supported on a mesh point of a given height.
500: Not Collective
502: Input Parameters:
503: + sp - the `PetscSpace` object
504: - height - the height of the mesh point for which the subspace is desired
506: Output Parameter:
507: . subsp - the subspace
509: Level: advanced
511: Notes:
512: If the space is not defined on mesh points of the given height (e.g. if the space is discontinuous and
513: pointwise values are not defined on the element boundaries), or if the implementation of `PetscSpace` does not
514: support extracting subspaces, then NULL is returned.
516: This does not increment the reference count on the returned space, and the user should not destroy it.
518: .seealso: `PetscDualSpaceGetHeightSubspace()`, `PetscSpace`
519: @*/
520: PetscErrorCode PetscSpaceGetHeightSubspace(PetscSpace sp, PetscInt height, PetscSpace *subsp)
521: {
522: PetscFunctionBegin;
525: *subsp = NULL;
526: PetscTryTypeMethod(sp, getheightsubspace, height, subsp);
527: PetscFunctionReturn(PETSC_SUCCESS);
528: }