Actual source code: spacesum.c
1: #include <petsc/private/petscfeimpl.h>
3: /*@
4: PetscSpaceSumGetNumSubspaces - Get the number of spaces in the sum space
6: Input Parameter:
7: . sp - the function space object
9: Output Parameter:
10: . numSumSpaces - the number of spaces
12: Level: intermediate
14: Note:
15: The name NumSubspaces is slightly misleading because it is actually getting the number of defining spaces of the sum, not a number of Subspaces of it
17: .seealso: `PETSCSPACESUM`, `PetscSpace`, `PetscSpaceSumSetNumSubspaces()`, `PetscSpaceSetDegree()`, `PetscSpaceSetNumVariables()`
18: @*/
19: PetscErrorCode PetscSpaceSumGetNumSubspaces(PetscSpace sp, PetscInt *numSumSpaces)
20: {
21: PetscFunctionBegin;
23: PetscAssertPointer(numSumSpaces, 2);
24: PetscTryMethod(sp, "PetscSpaceSumGetNumSubspaces_C", (PetscSpace, PetscInt *), (sp, numSumSpaces));
25: PetscFunctionReturn(PETSC_SUCCESS);
26: }
28: /*@
29: PetscSpaceSumSetNumSubspaces - Set the number of spaces in the sum space
31: Input Parameters:
32: + sp - the function space object
33: - numSumSpaces - the number of spaces
35: Level: intermediate
37: Note:
38: The name NumSubspaces is slightly misleading because it is actually setting the number of defining spaces of the sum, not a number of Subspaces of it
40: .seealso: `PETSCSPACESUM`, `PetscSpace`, `PetscSpaceSumGetNumSubspaces()`, `PetscSpaceSetDegree()`, `PetscSpaceSetNumVariables()`
41: @*/
42: PetscErrorCode PetscSpaceSumSetNumSubspaces(PetscSpace sp, PetscInt numSumSpaces)
43: {
44: PetscFunctionBegin;
46: PetscTryMethod(sp, "PetscSpaceSumSetNumSubspaces_C", (PetscSpace, PetscInt), (sp, numSumSpaces));
47: PetscFunctionReturn(PETSC_SUCCESS);
48: }
50: /*@
51: PetscSpaceSumGetConcatenate - Get the concatenate flag for this space.
53: Input Parameter:
54: . sp - the function space object
56: Output Parameter:
57: . concatenate - flag indicating whether subspaces are concatenated.
59: Level: intermediate
61: Notes:
62: A concatenated sum space will have the number of components equal to the sum of the number of
63: components of all subspaces. A non-concatenated, or direct sum space will have the same
64: number of components as its subspaces.
66: .seealso: `PETSCSPACESUM`, `PetscSpace`, `PetscSpaceSumSetConcatenate()`
67: @*/
68: PetscErrorCode PetscSpaceSumGetConcatenate(PetscSpace sp, PetscBool *concatenate)
69: {
70: PetscFunctionBegin;
72: PetscTryMethod(sp, "PetscSpaceSumGetConcatenate_C", (PetscSpace, PetscBool *), (sp, concatenate));
73: PetscFunctionReturn(PETSC_SUCCESS);
74: }
76: /*@
77: PetscSpaceSumSetConcatenate - Sets the concatenate flag for this space.
79: Input Parameters:
80: + sp - the function space object
81: - concatenate - are subspaces concatenated components (true) or direct summands (false)
83: Level: intermediate
85: Notes:
86: A concatenated sum space will have the number of components equal to the sum of the number of
87: components of all subspaces. A non-concatenated, or direct sum space will have the same
88: number of components as its subspaces .
90: .seealso: `PETSCSPACESUM`, `PetscSpace`, `PetscSpaceSumGetConcatenate()`
91: @*/
92: PetscErrorCode PetscSpaceSumSetConcatenate(PetscSpace sp, PetscBool concatenate)
93: {
94: PetscFunctionBegin;
96: PetscTryMethod(sp, "PetscSpaceSumSetConcatenate_C", (PetscSpace, PetscBool), (sp, concatenate));
97: PetscFunctionReturn(PETSC_SUCCESS);
98: }
100: /*@
101: PetscSpaceSumGetSubspace - Get a space in the sum space
103: Input Parameters:
104: + sp - the function space object
105: - s - The space number
107: Output Parameter:
108: . subsp - the `PetscSpace`
110: Level: intermediate
112: Note:
113: The name GetSubspace is slightly misleading because it is actually getting one of the defining spaces of the sum, not a Subspace of it
115: .seealso: `PETSCSPACESUM`, `PetscSpace`, `PetscSpaceSumSetSubspace()`, `PetscSpaceSetDegree()`, `PetscSpaceSetNumVariables()`
116: @*/
117: PetscErrorCode PetscSpaceSumGetSubspace(PetscSpace sp, PetscInt s, PetscSpace *subsp)
118: {
119: PetscFunctionBegin;
121: PetscAssertPointer(subsp, 3);
122: PetscTryMethod(sp, "PetscSpaceSumGetSubspace_C", (PetscSpace, PetscInt, PetscSpace *), (sp, s, subsp));
123: PetscFunctionReturn(PETSC_SUCCESS);
124: }
126: /*@
127: PetscSpaceSumSetSubspace - Set a space in the sum space
129: Input Parameters:
130: + sp - the function space object
131: . s - The space number
132: - subsp - the number of spaces
134: Level: intermediate
136: Note:
137: The name SetSubspace is slightly misleading because it is actually setting one of the defining spaces of the sum, not a Subspace of it
139: .seealso: `PETSCSPACESUM`, `PetscSpace`, `PetscSpaceSumGetSubspace()`, `PetscSpaceSetDegree()`, `PetscSpaceSetNumVariables()`
140: @*/
141: PetscErrorCode PetscSpaceSumSetSubspace(PetscSpace sp, PetscInt s, PetscSpace subsp)
142: {
143: PetscFunctionBegin;
146: PetscTryMethod(sp, "PetscSpaceSumSetSubspace_C", (PetscSpace, PetscInt, PetscSpace), (sp, s, subsp));
147: PetscFunctionReturn(PETSC_SUCCESS);
148: }
150: static PetscErrorCode PetscSpaceSumGetNumSubspaces_Sum(PetscSpace space, PetscInt *numSumSpaces)
151: {
152: PetscSpace_Sum *sum = (PetscSpace_Sum *)space->data;
154: PetscFunctionBegin;
155: *numSumSpaces = sum->numSumSpaces;
156: PetscFunctionReturn(PETSC_SUCCESS);
157: }
159: static PetscErrorCode PetscSpaceSumSetNumSubspaces_Sum(PetscSpace space, PetscInt numSumSpaces)
160: {
161: PetscSpace_Sum *sum = (PetscSpace_Sum *)space->data;
162: PetscInt Ns = sum->numSumSpaces;
164: PetscFunctionBegin;
165: PetscCheck(!sum->setupcalled, PetscObjectComm((PetscObject)space), PETSC_ERR_ARG_WRONGSTATE, "Cannot change number of subspaces after setup called");
166: if (numSumSpaces == Ns) PetscFunctionReturn(PETSC_SUCCESS);
167: if (Ns >= 0) {
168: PetscInt s;
169: for (s = 0; s < Ns; ++s) PetscCall(PetscSpaceDestroy(&sum->sumspaces[s]));
170: PetscCall(PetscFree(sum->sumspaces));
171: }
173: Ns = sum->numSumSpaces = numSumSpaces;
174: PetscCall(PetscCalloc1(Ns, &sum->sumspaces));
175: PetscFunctionReturn(PETSC_SUCCESS);
176: }
178: static PetscErrorCode PetscSpaceSumGetConcatenate_Sum(PetscSpace sp, PetscBool *concatenate)
179: {
180: PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
182: PetscFunctionBegin;
183: *concatenate = sum->concatenate;
184: PetscFunctionReturn(PETSC_SUCCESS);
185: }
187: static PetscErrorCode PetscSpaceSumSetConcatenate_Sum(PetscSpace sp, PetscBool concatenate)
188: {
189: PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
191: PetscFunctionBegin;
192: PetscCheck(!sum->setupcalled, PetscObjectComm((PetscObject)sp), PETSC_ERR_ARG_WRONGSTATE, "Cannot change space concatenation after setup called.");
194: sum->concatenate = concatenate;
195: PetscFunctionReturn(PETSC_SUCCESS);
196: }
198: static PetscErrorCode PetscSpaceSumGetSubspace_Sum(PetscSpace space, PetscInt s, PetscSpace *subspace)
199: {
200: PetscSpace_Sum *sum = (PetscSpace_Sum *)space->data;
201: PetscInt Ns = sum->numSumSpaces;
203: PetscFunctionBegin;
204: PetscCheck(Ns >= 0, PetscObjectComm((PetscObject)space), PETSC_ERR_ARG_WRONGSTATE, "Must call PetscSpaceSumSetNumSubspaces() first");
205: PetscCheck(s >= 0 && s < Ns, PetscObjectComm((PetscObject)space), PETSC_ERR_ARG_OUTOFRANGE, "Invalid subspace number %" PetscInt_FMT, s);
207: *subspace = sum->sumspaces[s];
208: PetscFunctionReturn(PETSC_SUCCESS);
209: }
211: static PetscErrorCode PetscSpaceSumSetSubspace_Sum(PetscSpace space, PetscInt s, PetscSpace subspace)
212: {
213: PetscSpace_Sum *sum = (PetscSpace_Sum *)space->data;
214: PetscInt Ns = sum->numSumSpaces;
216: PetscFunctionBegin;
217: PetscCheck(!sum->setupcalled, PetscObjectComm((PetscObject)space), PETSC_ERR_ARG_WRONGSTATE, "Cannot change subspace after setup called");
218: PetscCheck(Ns >= 0, PetscObjectComm((PetscObject)space), PETSC_ERR_ARG_WRONGSTATE, "Must call PetscSpaceSumSetNumSubspaces() first");
219: PetscCheck(s >= 0 && s < Ns, PetscObjectComm((PetscObject)space), PETSC_ERR_ARG_OUTOFRANGE, "Invalid subspace number %" PetscInt_FMT, s);
221: PetscCall(PetscObjectReference((PetscObject)subspace));
222: PetscCall(PetscSpaceDestroy(&sum->sumspaces[s]));
223: sum->sumspaces[s] = subspace;
224: PetscFunctionReturn(PETSC_SUCCESS);
225: }
227: static PetscErrorCode PetscSpaceSetFromOptions_Sum(PetscSpace sp, PetscOptionItems PetscOptionsObject)
228: {
229: PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
230: PetscInt Ns, Nc, Nv, deg, i;
231: PetscBool concatenate = PETSC_TRUE;
232: const char *prefix;
234: PetscFunctionBegin;
235: PetscCall(PetscSpaceGetNumVariables(sp, &Nv));
236: if (!Nv) PetscFunctionReturn(PETSC_SUCCESS);
237: PetscCall(PetscSpaceGetNumComponents(sp, &Nc));
238: PetscCall(PetscSpaceSumGetNumSubspaces(sp, &Ns));
239: PetscCall(PetscSpaceGetDegree(sp, °, NULL));
240: Ns = (Ns == PETSC_DEFAULT) ? 1 : Ns;
242: PetscOptionsHeadBegin(PetscOptionsObject, "PetscSpace sum options");
243: PetscCall(PetscOptionsBoundedInt("-petscspace_sum_spaces", "The number of subspaces", "PetscSpaceSumSetNumSubspaces", Ns, &Ns, NULL, 0));
244: PetscCall(PetscOptionsBool("-petscspace_sum_concatenate", "Subspaces are concatenated components of the final space", "PetscSpaceSumSetFromOptions", concatenate, &concatenate, NULL));
245: PetscOptionsHeadEnd();
247: PetscCheck(Ns >= 0 && (Nv <= 0 || Ns != 0), PetscObjectComm((PetscObject)sp), PETSC_ERR_ARG_OUTOFRANGE, "Cannot have a sum space of %" PetscInt_FMT " spaces", Ns);
248: if (Ns != sum->numSumSpaces) PetscCall(PetscSpaceSumSetNumSubspaces(sp, Ns));
249: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)sp, &prefix));
250: for (i = 0; i < Ns; ++i) {
251: PetscInt sNv;
252: PetscSpace subspace;
254: PetscCall(PetscSpaceSumGetSubspace(sp, i, &subspace));
255: if (!subspace) {
256: char subspacePrefix[256];
258: PetscCall(PetscSpaceCreate(PetscObjectComm((PetscObject)sp), &subspace));
259: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)subspace, prefix));
260: PetscCall(PetscSNPrintf(subspacePrefix, 256, "sumcomp_%" PetscInt_FMT "_", i));
261: PetscCall(PetscObjectAppendOptionsPrefix((PetscObject)subspace, subspacePrefix));
262: } else PetscCall(PetscObjectReference((PetscObject)subspace));
263: PetscCall(PetscSpaceSetFromOptions(subspace));
264: PetscCall(PetscSpaceGetNumVariables(subspace, &sNv));
265: PetscCheck(sNv, PetscObjectComm((PetscObject)sp), PETSC_ERR_ARG_WRONGSTATE, "Subspace %" PetscInt_FMT " has not been set properly, number of variables is 0.", i);
266: PetscCall(PetscSpaceSumSetSubspace(sp, i, subspace));
267: PetscCall(PetscSpaceDestroy(&subspace));
268: }
269: PetscFunctionReturn(PETSC_SUCCESS);
270: }
272: static PetscErrorCode PetscSpaceSetUp_Sum(PetscSpace sp)
273: {
274: PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
275: PetscBool concatenate = PETSC_TRUE;
276: PetscBool uniform;
277: PetscInt Nv, Ns, Nc, i, sum_Nc = 0, deg = PETSC_INT_MAX, maxDeg = PETSC_INT_MIN;
278: PetscInt minNc, maxNc;
280: PetscFunctionBegin;
281: if (sum->setupcalled) PetscFunctionReturn(PETSC_SUCCESS);
283: PetscCall(PetscSpaceGetNumVariables(sp, &Nv));
284: PetscCall(PetscSpaceGetNumComponents(sp, &Nc));
285: PetscCall(PetscSpaceSumGetNumSubspaces(sp, &Ns));
286: if (Ns == PETSC_DEFAULT) {
287: Ns = 1;
288: PetscCall(PetscSpaceSumSetNumSubspaces(sp, Ns));
289: }
290: PetscCheck(Ns >= 0, PetscObjectComm((PetscObject)sp), PETSC_ERR_ARG_OUTOFRANGE, "Cannot have %" PetscInt_FMT " subspaces", Ns);
291: uniform = PETSC_TRUE;
292: if (Ns) {
293: PetscSpace s0;
295: PetscCall(PetscSpaceSumGetSubspace(sp, 0, &s0));
296: for (PetscInt i = 1; i < Ns; i++) {
297: PetscSpace si;
299: PetscCall(PetscSpaceSumGetSubspace(sp, i, &si));
300: if (si != s0) {
301: uniform = PETSC_FALSE;
302: break;
303: }
304: }
305: }
307: minNc = Nc;
308: maxNc = Nc;
309: for (i = 0; i < Ns; ++i) {
310: PetscInt sNv, sNc, iDeg, iMaxDeg;
311: PetscSpace si;
313: PetscCall(PetscSpaceSumGetSubspace(sp, i, &si));
314: PetscCall(PetscSpaceSetUp(si));
315: PetscCall(PetscSpaceGetNumVariables(si, &sNv));
316: PetscCheck(sNv == Nv, PetscObjectComm((PetscObject)sp), PETSC_ERR_ARG_WRONGSTATE, "Subspace %" PetscInt_FMT " has %" PetscInt_FMT " variables, space has %" PetscInt_FMT ".", i, sNv, Nv);
317: PetscCall(PetscSpaceGetNumComponents(si, &sNc));
318: if (i == 0 && sNc == Nc) concatenate = PETSC_FALSE;
319: minNc = PetscMin(minNc, sNc);
320: maxNc = PetscMax(maxNc, sNc);
321: sum_Nc += sNc;
322: PetscCall(PetscSpaceSumGetSubspace(sp, i, &si));
323: PetscCall(PetscSpaceGetDegree(si, &iDeg, &iMaxDeg));
324: deg = PetscMin(deg, iDeg);
325: maxDeg = PetscMax(maxDeg, iMaxDeg);
326: }
328: if (concatenate) PetscCheck(sum_Nc == Nc, PetscObjectComm((PetscObject)sp), PETSC_ERR_ARG_OUTOFRANGE, "Total number of subspace components (%" PetscInt_FMT ") does not match number of target space components (%" PetscInt_FMT ").", sum_Nc, Nc);
329: else PetscCheck(minNc == Nc && maxNc == Nc, PetscObjectComm((PetscObject)sp), PETSC_ERR_ARG_OUTOFRANGE, "Subspaces must have same number of components as the target space.");
331: sp->degree = deg;
332: sp->maxDegree = maxDeg;
333: sum->concatenate = concatenate;
334: sum->uniform = uniform;
335: sum->setupcalled = PETSC_TRUE;
336: PetscFunctionReturn(PETSC_SUCCESS);
337: }
339: static PetscErrorCode PetscSpaceSumView_Ascii(PetscSpace sp, PetscViewer v)
340: {
341: PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
342: PetscBool concatenate = sum->concatenate;
343: PetscInt i, Ns = sum->numSumSpaces;
345: PetscFunctionBegin;
346: if (concatenate) PetscCall(PetscViewerASCIIPrintf(v, "Sum space of %" PetscInt_FMT " concatenated subspaces%s\n", Ns, sum->uniform ? " (all identical)" : ""));
347: else PetscCall(PetscViewerASCIIPrintf(v, "Sum space of %" PetscInt_FMT " subspaces%s\n", Ns, sum->uniform ? " (all identical)" : ""));
348: for (i = 0; i < (sum->uniform ? (Ns > 0 ? 1 : 0) : Ns); ++i) {
349: PetscCall(PetscViewerASCIIPushTab(v));
350: PetscCall(PetscSpaceView(sum->sumspaces[i], v));
351: PetscCall(PetscViewerASCIIPopTab(v));
352: }
353: PetscFunctionReturn(PETSC_SUCCESS);
354: }
356: static PetscErrorCode PetscSpaceView_Sum(PetscSpace sp, PetscViewer viewer)
357: {
358: PetscBool isascii;
360: PetscFunctionBegin;
361: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
362: if (isascii) PetscCall(PetscSpaceSumView_Ascii(sp, viewer));
363: PetscFunctionReturn(PETSC_SUCCESS);
364: }
366: static PetscErrorCode PetscSpaceDestroy_Sum(PetscSpace sp)
367: {
368: PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
369: PetscInt i, Ns = sum->numSumSpaces;
371: PetscFunctionBegin;
372: for (i = 0; i < Ns; ++i) PetscCall(PetscSpaceDestroy(&sum->sumspaces[i]));
373: PetscCall(PetscFree(sum->sumspaces));
374: if (sum->heightsubspaces) {
375: /* sp->Nv is the spatial dimension, so it is equal to the number
376: * of subspaces on higher co-dimension points */
377: for (PetscInt d = 0; d < sp->Nv; ++d) PetscCall(PetscSpaceDestroy(&sum->heightsubspaces[d]));
378: }
379: PetscCall(PetscFree(sum->heightsubspaces));
380: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumSetSubspace_C", NULL));
381: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumGetSubspace_C", NULL));
382: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumSetNumSubspaces_C", NULL));
383: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumGetNumSubspaces_C", NULL));
384: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumGetConcatenate_C", NULL));
385: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumSetConcatenate_C", NULL));
386: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumGetInterleave_C", NULL));
387: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumSetInterleave_C", NULL));
388: PetscCall(PetscFree(sum));
389: PetscFunctionReturn(PETSC_SUCCESS);
390: }
392: static PetscErrorCode PetscSpaceGetDimension_Sum(PetscSpace sp, PetscInt *dim)
393: {
394: PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
395: PetscInt i, d = 0, Ns = sum->numSumSpaces;
397: PetscFunctionBegin;
398: if (!sum->setupcalled) {
399: PetscCall(PetscSpaceSetUp(sp));
400: PetscCall(PetscSpaceGetDimension(sp, dim));
401: PetscFunctionReturn(PETSC_SUCCESS);
402: }
404: for (i = 0; i < Ns; ++i) {
405: PetscInt id;
407: PetscCall(PetscSpaceGetDimension(sum->sumspaces[i], &id));
408: d += id;
409: }
411: *dim = d;
412: PetscFunctionReturn(PETSC_SUCCESS);
413: }
415: static PetscErrorCode PetscSpaceEvaluate_Sum(PetscSpace sp, PetscInt npoints, const PetscReal points[], PetscReal B[], PetscReal D[], PetscReal H[])
416: {
417: PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
418: PetscBool concatenate = sum->concatenate;
419: DM dm = sp->dm;
420: PetscInt Nc = sp->Nc, Nv = sp->Nv, Ns = sum->numSumSpaces;
421: PetscInt i, s, offset, ncoffset, pdimfull, numelB, numelD, numelH;
422: PetscReal *sB = NULL, *sD = NULL, *sH = NULL;
424: PetscFunctionBegin;
425: if (!sum->setupcalled) {
426: PetscCall(PetscSpaceSetUp(sp));
427: PetscCall(PetscSpaceEvaluate(sp, npoints, points, B, D, H));
428: PetscFunctionReturn(PETSC_SUCCESS);
429: }
430: PetscCall(PetscSpaceGetDimension(sp, &pdimfull));
431: numelB = npoints * pdimfull * Nc;
432: numelD = numelB * Nv;
433: numelH = numelD * Nv;
434: if (B || D || H) PetscCall(DMGetWorkArray(dm, numelB, MPIU_REAL, &sB));
435: if (D || H) PetscCall(DMGetWorkArray(dm, numelD, MPIU_REAL, &sD));
436: if (H) PetscCall(DMGetWorkArray(dm, numelH, MPIU_REAL, &sH));
437: if (B)
438: for (i = 0; i < numelB; ++i) B[i] = 0.;
439: if (D)
440: for (i = 0; i < numelD; ++i) D[i] = 0.;
441: if (H)
442: for (i = 0; i < numelH; ++i) H[i] = 0.;
444: for (s = 0, offset = 0, ncoffset = 0; s < Ns; ++s) {
445: PetscInt sNv, spdim, sNc, p;
447: PetscCall(PetscSpaceGetNumVariables(sum->sumspaces[s], &sNv));
448: PetscCall(PetscSpaceGetNumComponents(sum->sumspaces[s], &sNc));
449: PetscCall(PetscSpaceGetDimension(sum->sumspaces[s], &spdim));
450: PetscCheck(offset + spdim <= pdimfull, PetscObjectComm((PetscObject)sp), PETSC_ERR_ARG_OUTOFRANGE, "Subspace dimensions exceed target space dimension.");
451: if (s == 0 || !sum->uniform) PetscCall(PetscSpaceEvaluate(sum->sumspaces[s], npoints, points, sB, sD, sH));
452: if (B || D || H) {
453: for (p = 0; p < npoints; ++p) {
454: for (PetscInt j = 0; j < spdim; ++j) {
455: PetscInt b = sum->interleave_basis ? (j * Ns + s) : (j + offset);
457: for (PetscInt c = 0; c < sNc; ++c) {
458: PetscInt compoffset, BInd, sBInd;
460: compoffset = concatenate ? (sum->interleave_components ? (c * Ns + s) : (c + ncoffset)) : c;
461: BInd = (p * pdimfull + b) * Nc + compoffset;
462: sBInd = (p * spdim + j) * sNc + c;
463: if (B) B[BInd] = sB[sBInd];
464: if (D || H) {
465: for (PetscInt v = 0; v < Nv; ++v) {
466: PetscInt DInd, sDInd;
468: DInd = BInd * Nv + v;
469: sDInd = sBInd * Nv + v;
470: if (D) D[DInd] = sD[sDInd];
471: if (H) {
472: for (PetscInt v2 = 0; v2 < Nv; ++v2) {
473: PetscInt HInd, sHInd;
475: HInd = DInd * Nv + v2;
476: sHInd = sDInd * Nv + v2;
477: H[HInd] = sH[sHInd];
478: }
479: }
480: }
481: }
482: }
483: }
484: }
485: }
486: offset += spdim;
487: ncoffset += sNc;
488: }
490: if (H) PetscCall(DMRestoreWorkArray(dm, numelH, MPIU_REAL, &sH));
491: if (D || H) PetscCall(DMRestoreWorkArray(dm, numelD, MPIU_REAL, &sD));
492: if (B || D || H) PetscCall(DMRestoreWorkArray(dm, numelB, MPIU_REAL, &sB));
493: PetscFunctionReturn(PETSC_SUCCESS);
494: }
496: static PetscErrorCode PetscSpaceGetHeightSubspace_Sum(PetscSpace sp, PetscInt height, PetscSpace *subsp)
497: {
498: PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
499: PetscInt Nc, dim, order;
500: PetscBool tensor;
502: PetscFunctionBegin;
503: PetscCall(PetscSpaceGetNumComponents(sp, &Nc));
504: PetscCall(PetscSpaceGetNumVariables(sp, &dim));
505: PetscCall(PetscSpaceGetDegree(sp, &order, NULL));
506: PetscCall(PetscSpacePolynomialGetTensor(sp, &tensor));
507: PetscCheck(height <= dim && height >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Asked for space at height %" PetscInt_FMT " for dimension %" PetscInt_FMT " space", height, dim);
508: if (!sum->heightsubspaces) PetscCall(PetscCalloc1(dim, &sum->heightsubspaces));
509: if (height <= dim) {
510: if (!sum->heightsubspaces[height - 1]) {
511: PetscSpace sub;
512: const char *name;
514: PetscCall(PetscSpaceCreate(PetscObjectComm((PetscObject)sp), &sub));
515: PetscCall(PetscObjectGetName((PetscObject)sp, &name));
516: PetscCall(PetscObjectSetName((PetscObject)sub, name));
517: PetscCall(PetscSpaceSetType(sub, PETSCSPACESUM));
518: PetscCall(PetscSpaceSumSetNumSubspaces(sub, sum->numSumSpaces));
519: PetscCall(PetscSpaceSumSetConcatenate(sub, sum->concatenate));
520: PetscCall(PetscSpaceSetNumComponents(sub, Nc));
521: PetscCall(PetscSpaceSetNumVariables(sub, dim - height));
522: for (PetscInt i = 0; i < sum->numSumSpaces; i++) {
523: PetscSpace subh;
525: PetscCall(PetscSpaceGetHeightSubspace(sum->sumspaces[i], height, &subh));
526: PetscCall(PetscSpaceSumSetSubspace(sub, i, subh));
527: }
528: PetscCall(PetscSpaceSetUp(sub));
529: sum->heightsubspaces[height - 1] = sub;
530: }
531: *subsp = sum->heightsubspaces[height - 1];
532: } else {
533: *subsp = NULL;
534: }
535: PetscFunctionReturn(PETSC_SUCCESS);
536: }
538: /*@
539: PetscSpaceSumSetInterleave - Set whether the basis functions and components of a uniform sum are interleaved
541: Logically collective
543: Input Parameters:
544: + sp - a `PetscSpace` of type `PETSCSPACESUM`
545: . interleave_basis - if `PETSC_TRUE`, the basis vectors of the subspaces are interleaved
546: - interleave_components - if `PETSC_TRUE` and the space concatenates components (`PetscSpaceSumGetConcatenate()`),
547: interleave the concatenated components
549: Level: developer
551: .seealso: `PetscSpace`, `PETSCSPACESUM`, `PETSCFEVECTOR`, `PetscSpaceSumGetInterleave()`
552: @*/
553: PetscErrorCode PetscSpaceSumSetInterleave(PetscSpace sp, PetscBool interleave_basis, PetscBool interleave_components)
554: {
555: PetscFunctionBegin;
557: PetscTryMethod(sp, "PetscSpaceSumSetInterleave_C", (PetscSpace, PetscBool, PetscBool), (sp, interleave_basis, interleave_components));
558: PetscFunctionReturn(PETSC_SUCCESS);
559: }
561: static PetscErrorCode PetscSpaceSumSetInterleave_Sum(PetscSpace sp, PetscBool interleave_basis, PetscBool interleave_components)
562: {
563: PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
565: PetscFunctionBegin;
566: sum->interleave_basis = interleave_basis;
567: sum->interleave_components = interleave_components;
568: PetscFunctionReturn(PETSC_SUCCESS);
569: }
571: /*@
572: PetscSpaceSumGetInterleave - Get whether the basis functions and components of a uniform sum are interleaved
574: Logically collective
576: Input Parameter:
577: . sp - a `PetscSpace` of type `PETSCSPACESUM`
579: Output Parameters:
580: + interleave_basis - if `PETSC_TRUE`, the basis vectors of the subspaces are interleaved
581: - interleave_components - if `PETSC_TRUE` and the space concatenates components (`PetscSpaceSumGetConcatenate()`),
582: interleave the concatenated components
584: Level: developer
586: .seealso: `PetscSpace`, `PETSCSPACESUM`, `PETSCFEVECTOR`, `PetscSpaceSumSetInterleave()`
587: @*/
588: PetscErrorCode PetscSpaceSumGetInterleave(PetscSpace sp, PeOp PetscBool *interleave_basis, PeOp PetscBool *interleave_components)
589: {
590: PetscFunctionBegin;
592: if (interleave_basis) PetscAssertPointer(interleave_basis, 2);
593: if (interleave_components) PetscAssertPointer(interleave_components, 3);
594: PetscTryMethod(sp, "PetscSpaceSumGetInterleave_C", (PetscSpace, PetscBool *, PetscBool *), (sp, interleave_basis, interleave_components));
595: PetscFunctionReturn(PETSC_SUCCESS);
596: }
598: static PetscErrorCode PetscSpaceSumGetInterleave_Sum(PetscSpace sp, PetscBool *interleave_basis, PetscBool *interleave_components)
599: {
600: PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
602: PetscFunctionBegin;
603: if (interleave_basis) *interleave_basis = sum->interleave_basis;
604: if (interleave_components) *interleave_components = sum->interleave_components;
605: PetscFunctionReturn(PETSC_SUCCESS);
606: }
608: static PetscErrorCode PetscSpaceInitialize_Sum(PetscSpace sp)
609: {
610: PetscFunctionBegin;
611: sp->ops->setfromoptions = PetscSpaceSetFromOptions_Sum;
612: sp->ops->setup = PetscSpaceSetUp_Sum;
613: sp->ops->view = PetscSpaceView_Sum;
614: sp->ops->destroy = PetscSpaceDestroy_Sum;
615: sp->ops->getdimension = PetscSpaceGetDimension_Sum;
616: sp->ops->evaluate = PetscSpaceEvaluate_Sum;
617: sp->ops->getheightsubspace = PetscSpaceGetHeightSubspace_Sum;
619: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumGetNumSubspaces_C", PetscSpaceSumGetNumSubspaces_Sum));
620: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumSetNumSubspaces_C", PetscSpaceSumSetNumSubspaces_Sum));
621: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumGetSubspace_C", PetscSpaceSumGetSubspace_Sum));
622: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumSetSubspace_C", PetscSpaceSumSetSubspace_Sum));
623: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumGetConcatenate_C", PetscSpaceSumGetConcatenate_Sum));
624: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumSetConcatenate_C", PetscSpaceSumSetConcatenate_Sum));
625: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumGetInterleave_C", PetscSpaceSumGetInterleave_Sum));
626: PetscCall(PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumSetInterleave_C", PetscSpaceSumSetInterleave_Sum));
627: PetscFunctionReturn(PETSC_SUCCESS);
628: }
630: /*MC
631: PETSCSPACESUM = "sum" - A `PetscSpace` object that encapsulates a sum of subspaces.
633: Level: intermediate
635: Note:
636: That sum can either be direct or a concatenation. For example if A and B are spaces each with 2 components,
637: the direct sum of A and B will also have 2 components while the concatenated sum will have 4 components. In both cases A and B must be defined over the
638: same number of variables.
640: .seealso: `PetscSpace`, `PetscSpaceType`, `PetscSpaceCreate()`, `PetscSpaceSetType()`, `PetscSpaceSumGetNumSubspaces()`, `PetscSpaceSumSetNumSubspaces()`,
641: `PetscSpaceSumGetConcatenate()`, `PetscSpaceSumSetConcatenate()`
642: M*/
643: PETSC_EXTERN PetscErrorCode PetscSpaceCreate_Sum(PetscSpace sp)
644: {
645: PetscSpace_Sum *sum;
647: PetscFunctionBegin;
649: PetscCall(PetscNew(&sum));
650: sum->numSumSpaces = PETSC_DEFAULT;
651: sp->data = sum;
652: PetscCall(PetscSpaceInitialize_Sum(sp));
653: PetscFunctionReturn(PETSC_SUCCESS);
654: }
656: PETSC_EXTERN PetscErrorCode PetscSpaceCreateSum(PetscInt numSubspaces, const PetscSpace subspaces[], PetscBool concatenate, PetscSpace *sumSpace)
657: {
658: PetscInt i, Nv, Nc = 0;
660: PetscFunctionBegin;
661: PetscCall(PetscSpaceCreate(PetscObjectComm((PetscObject)subspaces[0]), sumSpace));
662: PetscCall(PetscSpaceSetType(*sumSpace, PETSCSPACESUM));
663: PetscCall(PetscSpaceSumSetNumSubspaces(*sumSpace, numSubspaces));
664: PetscCall(PetscSpaceSumSetConcatenate(*sumSpace, concatenate));
665: for (i = 0; i < numSubspaces; ++i) {
666: PetscInt sNc;
668: PetscCall(PetscSpaceSumSetSubspace(*sumSpace, i, subspaces[i]));
669: PetscCall(PetscSpaceGetNumComponents(subspaces[i], &sNc));
670: if (concatenate) Nc += sNc;
671: else Nc = sNc;
672: }
673: PetscCall(PetscSpaceGetNumVariables(subspaces[0], &Nv));
674: PetscCall(PetscSpaceSetNumComponents(*sumSpace, Nc));
675: PetscCall(PetscSpaceSetNumVariables(*sumSpace, Nv));
676: PetscCall(PetscSpaceSetUp(*sumSpace));
677: PetscFunctionReturn(PETSC_SUCCESS);
678: }