Actual source code: spacesum.c

  1: #include <petsc/private/petscfeimpl.h>
  2: /*@
  3:   PetscSpaceSumGetNumSubspaces - Get the number of spaces in the sum

  5:   Input Parameter:
  6: . sp  - the function space object

  8:   Output Parameter:
  9: . numSumSpaces - the number of spaces

 11: Level: intermediate

 13: .seealso: `PETSCSPACESUM`, `PetscSpace`, `PetscSpaceSumSetNumSubspaces()`, `PetscSpaceSetDegree()`, `PetscSpaceSetNumVariables()`
 14: @*/
 15: PetscErrorCode PetscSpaceSumGetNumSubspaces(PetscSpace sp, PetscInt *numSumSpaces)
 16: {
 19:   PetscTryMethod(sp, "PetscSpaceSumGetNumSubspaces_C", (PetscSpace, PetscInt *), (sp, numSumSpaces));
 20:   return 0;
 21: }

 23: /*@
 24:   PetscSpaceSumSetNumSubspaces - Set the number of spaces in the sum

 26:   Input Parameters:
 27: + sp  - the function space object
 28: - numSumSpaces - the number of spaces

 30: Level: intermediate

 32: .seealso: `PETSCSPACESUM`, `PetscSpace`, `PetscSpaceSumGetNumSubspaces()`, `PetscSpaceSetDegree()`, `PetscSpaceSetNumVariables()`
 33: @*/
 34: PetscErrorCode PetscSpaceSumSetNumSubspaces(PetscSpace sp, PetscInt numSumSpaces)
 35: {
 37:   PetscTryMethod(sp, "PetscSpaceSumSetNumSubspaces_C", (PetscSpace, PetscInt), (sp, numSumSpaces));
 38:   return 0;
 39: }

 41: /*@
 42:  PetscSpaceSumGetConcatenate - Get the concatenate flag for this space.
 43:  A concatenated sum space will have number of components equal to the sum of the number of components of all subspaces. A non-concatenated,
 44:  or direct sum space will have the same number of components as its subspaces .

 46:  Input Parameters:
 47: . sp - the function space object

 49:  Output Parameters:
 50: . concatenate - flag indicating whether subspaces are concatenated.

 52: Level: intermediate

 54: .seealso: `PETSCSPACESUM`, `PetscSpace`, `PetscSpaceSumSetConcatenate()`
 55: @*/
 56: PetscErrorCode PetscSpaceSumGetConcatenate(PetscSpace sp, PetscBool *concatenate)
 57: {
 59:   PetscTryMethod(sp, "PetscSpaceSumGetConcatenate_C", (PetscSpace, PetscBool *), (sp, concatenate));
 60:   return 0;
 61: }

 63: /*@
 64:   PetscSpaceSumSetConcatenate - Sets the concatenate flag for this space.
 65:   A concatenated sum space will have number of components equal to the sum of the number of components of all subspaces. A non-concatenated,
 66:   or direct sum space will have the same number of components as its subspaces .

 68:  Input Parameters:
 69: + sp - the function space object
 70: - concatenate - are subspaces concatenated components (true) or direct summands (false)

 72:   Level: intermediate

 74: .seealso: `PETSCSPACESUM`, `PetscSpace`, `PetscSpaceSumGetConcatenate()`
 75: @*/
 76: PetscErrorCode PetscSpaceSumSetConcatenate(PetscSpace sp, PetscBool concatenate)
 77: {
 79:   PetscTryMethod(sp, "PetscSpaceSumSetConcatenate_C", (PetscSpace, PetscBool), (sp, concatenate));
 80:   return 0;
 81: }

 83: /*@
 84:   PetscSpaceSumGetSubspace - Get a space in the sum

 86:   Input Parameters:
 87: + sp - the function space object
 88: - s  - The space number

 90:   Output Parameter:
 91: . subsp - the PetscSpace

 93:   Level: intermediate

 95: .seealso: `PETSCSPACESUM`, `PetscSpace`, `PetscSpaceSumSetSubspace()`, `PetscSpaceSetDegree()`, `PetscSpaceSetNumVariables()`
 96: @*/
 97: PetscErrorCode PetscSpaceSumGetSubspace(PetscSpace sp, PetscInt s, PetscSpace *subsp)
 98: {
101:   PetscTryMethod(sp, "PetscSpaceSumGetSubspace_C", (PetscSpace, PetscInt, PetscSpace *), (sp, s, subsp));
102:   return 0;
103: }

105: /*@
106:   PetscSpaceSumSetSubspace - Set a space in the sum

108:   Input Parameters:
109: + sp    - the function space object
110: . s     - The space number
111: - subsp - the number of spaces

113:   Level: intermediate

115: .seealso: `PETSCSPACESUM`, `PetscSpace`, `PetscSpaceSumGetSubspace()`, `PetscSpaceSetDegree()`, `PetscSpaceSetNumVariables()`
116: @*/
117: PetscErrorCode PetscSpaceSumSetSubspace(PetscSpace sp, PetscInt s, PetscSpace subsp)
118: {
121:   PetscTryMethod(sp, "PetscSpaceSumSetSubspace_C", (PetscSpace, PetscInt, PetscSpace), (sp, s, subsp));
122:   return 0;
123: }

125: static PetscErrorCode PetscSpaceSumGetNumSubspaces_Sum(PetscSpace space, PetscInt *numSumSpaces)
126: {
127:   PetscSpace_Sum *sum = (PetscSpace_Sum *)space->data;

129:   *numSumSpaces = sum->numSumSpaces;
130:   return 0;
131: }

133: static PetscErrorCode PetscSpaceSumSetNumSubspaces_Sum(PetscSpace space, PetscInt numSumSpaces)
134: {
135:   PetscSpace_Sum *sum = (PetscSpace_Sum *)space->data;
136:   PetscInt        Ns  = sum->numSumSpaces;

139:   if (numSumSpaces == Ns) return 0;
140:   if (Ns >= 0) {
141:     PetscInt s;
142:     for (s = 0; s < Ns; ++s) PetscSpaceDestroy(&sum->sumspaces[s]);
143:     PetscFree(sum->sumspaces);
144:   }

146:   Ns = sum->numSumSpaces = numSumSpaces;
147:   PetscCalloc1(Ns, &sum->sumspaces);
148:   return 0;
149: }

151: static PetscErrorCode PetscSpaceSumGetConcatenate_Sum(PetscSpace sp, PetscBool *concatenate)
152: {
153:   PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;

155:   *concatenate = sum->concatenate;
156:   return 0;
157: }

159: static PetscErrorCode PetscSpaceSumSetConcatenate_Sum(PetscSpace sp, PetscBool concatenate)
160: {
161:   PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;


165:   sum->concatenate = concatenate;
166:   return 0;
167: }

169: static PetscErrorCode PetscSpaceSumGetSubspace_Sum(PetscSpace space, PetscInt s, PetscSpace *subspace)
170: {
171:   PetscSpace_Sum *sum = (PetscSpace_Sum *)space->data;
172:   PetscInt        Ns  = sum->numSumSpaces;


177:   *subspace = sum->sumspaces[s];
178:   return 0;
179: }

181: static PetscErrorCode PetscSpaceSumSetSubspace_Sum(PetscSpace space, PetscInt s, PetscSpace subspace)
182: {
183:   PetscSpace_Sum *sum = (PetscSpace_Sum *)space->data;
184:   PetscInt        Ns  = sum->numSumSpaces;


190:   PetscObjectReference((PetscObject)subspace);
191:   PetscSpaceDestroy(&sum->sumspaces[s]);
192:   sum->sumspaces[s] = subspace;
193:   return 0;
194: }

196: static PetscErrorCode PetscSpaceSetFromOptions_Sum(PetscSpace sp, PetscOptionItems *PetscOptionsObject)
197: {
198:   PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
199:   PetscInt        Ns, Nc, Nv, deg, i;
200:   PetscBool       concatenate = PETSC_TRUE;
201:   const char     *prefix;

203:   PetscSpaceGetNumVariables(sp, &Nv);
204:   if (!Nv) return 0;
205:   PetscSpaceGetNumComponents(sp, &Nc);
206:   PetscSpaceSumGetNumSubspaces(sp, &Ns);
207:   PetscSpaceGetDegree(sp, &deg, NULL);
208:   Ns = (Ns == PETSC_DEFAULT) ? 1 : Ns;

210:   PetscOptionsHeadBegin(PetscOptionsObject, "PetscSpace sum options");
211:   PetscOptionsBoundedInt("-petscspace_sum_spaces", "The number of subspaces", "PetscSpaceSumSetNumSubspaces", Ns, &Ns, NULL, 0);
212:   PetscOptionsBool("-petscspace_sum_concatenate", "Subspaces are concatenated components of the final space", "PetscSpaceSumSetFromOptions", concatenate, &concatenate, NULL);
213:   PetscOptionsHeadEnd();

216:   if (Ns != sum->numSumSpaces) PetscSpaceSumSetNumSubspaces(sp, Ns);
217:   PetscObjectGetOptionsPrefix((PetscObject)sp, &prefix);
218:   for (i = 0; i < Ns; ++i) {
219:     PetscInt   sNv;
220:     PetscSpace subspace;

222:     PetscSpaceSumGetSubspace(sp, i, &subspace);
223:     if (!subspace) {
224:       char subspacePrefix[256];

226:       PetscSpaceCreate(PetscObjectComm((PetscObject)sp), &subspace);
227:       PetscObjectSetOptionsPrefix((PetscObject)subspace, prefix);
228:       PetscSNPrintf(subspacePrefix, 256, "sumcomp_%" PetscInt_FMT "_", i);
229:       PetscObjectAppendOptionsPrefix((PetscObject)subspace, subspacePrefix);
230:     } else PetscObjectReference((PetscObject)subspace);
231:     PetscSpaceSetFromOptions(subspace);
232:     PetscSpaceGetNumVariables(subspace, &sNv);
234:     PetscSpaceSumSetSubspace(sp, i, subspace);
235:     PetscSpaceDestroy(&subspace);
236:   }
237:   return 0;
238: }

240: static PetscErrorCode PetscSpaceSetUp_Sum(PetscSpace sp)
241: {
242:   PetscSpace_Sum *sum         = (PetscSpace_Sum *)sp->data;
243:   PetscBool       concatenate = PETSC_TRUE;
244:   PetscBool       uniform;
245:   PetscInt        Nv, Ns, Nc, i, sum_Nc = 0, deg = PETSC_MAX_INT, maxDeg = PETSC_MIN_INT;
246:   PetscInt        minNc, maxNc;

248:   if (sum->setupCalled) return 0;

250:   PetscSpaceGetNumVariables(sp, &Nv);
251:   PetscSpaceGetNumComponents(sp, &Nc);
252:   PetscSpaceSumGetNumSubspaces(sp, &Ns);
253:   if (Ns == PETSC_DEFAULT) {
254:     Ns = 1;
255:     PetscSpaceSumSetNumSubspaces(sp, Ns);
256:   }
258:   uniform = PETSC_TRUE;
259:   if (Ns) {
260:     PetscSpace s0;

262:     PetscSpaceSumGetSubspace(sp, 0, &s0);
263:     for (PetscInt i = 1; i < Ns; i++) {
264:       PetscSpace si;

266:       PetscSpaceSumGetSubspace(sp, i, &si);
267:       if (si != s0) {
268:         uniform = PETSC_FALSE;
269:         break;
270:       }
271:     }
272:   }

274:   minNc = Nc;
275:   maxNc = Nc;
276:   for (i = 0; i < Ns; ++i) {
277:     PetscInt   sNv, sNc, iDeg, iMaxDeg;
278:     PetscSpace si;

280:     PetscSpaceSumGetSubspace(sp, i, &si);
281:     PetscSpaceSetUp(si);
282:     PetscSpaceGetNumVariables(si, &sNv);
284:     PetscSpaceGetNumComponents(si, &sNc);
285:     if (i == 0 && sNc == Nc) concatenate = PETSC_FALSE;
286:     minNc = PetscMin(minNc, sNc);
287:     maxNc = PetscMax(maxNc, sNc);
288:     sum_Nc += sNc;
289:     PetscSpaceSumGetSubspace(sp, i, &si);
290:     PetscSpaceGetDegree(si, &iDeg, &iMaxDeg);
291:     deg    = PetscMin(deg, iDeg);
292:     maxDeg = PetscMax(maxDeg, iMaxDeg);
293:   }


298:   sp->degree       = deg;
299:   sp->maxDegree    = maxDeg;
300:   sum->concatenate = concatenate;
301:   sum->uniform     = uniform;
302:   sum->setupCalled = PETSC_TRUE;
303:   return 0;
304: }

306: static PetscErrorCode PetscSpaceSumView_Ascii(PetscSpace sp, PetscViewer v)
307: {
308:   PetscSpace_Sum *sum         = (PetscSpace_Sum *)sp->data;
309:   PetscBool       concatenate = sum->concatenate;
310:   PetscInt        i, Ns = sum->numSumSpaces;

312:   if (concatenate) PetscViewerASCIIPrintf(v, "Sum space of %" PetscInt_FMT " concatenated subspaces%s\n", Ns, sum->uniform ? " (all identical)" : "");
313:   else PetscViewerASCIIPrintf(v, "Sum space of %" PetscInt_FMT " subspaces%s\n", Ns, sum->uniform ? " (all identical)" : "");
314:   for (i = 0; i < (sum->uniform ? (Ns > 0 ? 1 : 0) : Ns); ++i) {
315:     PetscViewerASCIIPushTab(v);
316:     PetscSpaceView(sum->sumspaces[i], v);
317:     PetscViewerASCIIPopTab(v);
318:   }
319:   return 0;
320: }

322: static PetscErrorCode PetscSpaceView_Sum(PetscSpace sp, PetscViewer viewer)
323: {
324:   PetscBool iascii;

326:   PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii);
327:   if (iascii) PetscSpaceSumView_Ascii(sp, viewer);
328:   return 0;
329: }

331: static PetscErrorCode PetscSpaceDestroy_Sum(PetscSpace sp)
332: {
333:   PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
334:   PetscInt        i, Ns = sum->numSumSpaces;

336:   for (i = 0; i < Ns; ++i) PetscSpaceDestroy(&sum->sumspaces[i]);
337:   PetscFree(sum->sumspaces);
338:   if (sum->heightsubspaces) {
339:     PetscInt d;

341:     /* sp->Nv is the spatial dimension, so it is equal to the number
342:      * of subspaces on higher co-dimension points */
343:     for (d = 0; d < sp->Nv; ++d) PetscSpaceDestroy(&sum->heightsubspaces[d]);
344:   }
345:   PetscFree(sum->heightsubspaces);
346:   PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumSetSubspace_C", NULL);
347:   PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumGetSubspace_C", NULL);
348:   PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumSetNumSubspaces_C", NULL);
349:   PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumGetNumSubspaces_C", NULL);
350:   PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumGetConcatenate_C", NULL);
351:   PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumSetConcatenate_C", NULL);
352:   PetscFree(sum);
353:   return 0;
354: }

356: static PetscErrorCode PetscSpaceGetDimension_Sum(PetscSpace sp, PetscInt *dim)
357: {
358:   PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
359:   PetscInt        i, d = 0, Ns = sum->numSumSpaces;

361:   if (!sum->setupCalled) {
362:     PetscSpaceSetUp(sp);
363:     PetscSpaceGetDimension(sp, dim);
364:     return 0;
365:   }

367:   for (i = 0; i < Ns; ++i) {
368:     PetscInt id;

370:     PetscSpaceGetDimension(sum->sumspaces[i], &id);
371:     d += id;
372:   }

374:   *dim = d;
375:   return 0;
376: }

378: static PetscErrorCode PetscSpaceEvaluate_Sum(PetscSpace sp, PetscInt npoints, const PetscReal points[], PetscReal B[], PetscReal D[], PetscReal H[])
379: {
380:   PetscSpace_Sum *sum         = (PetscSpace_Sum *)sp->data;
381:   PetscBool       concatenate = sum->concatenate;
382:   DM              dm          = sp->dm;
383:   PetscInt        Nc = sp->Nc, Nv = sp->Nv, Ns = sum->numSumSpaces;
384:   PetscInt        i, s, offset, ncoffset, pdimfull, numelB, numelD, numelH;
385:   PetscReal      *sB = NULL, *sD = NULL, *sH = NULL;

387:   if (!sum->setupCalled) {
388:     PetscSpaceSetUp(sp);
389:     PetscSpaceEvaluate(sp, npoints, points, B, D, H);
390:     return 0;
391:   }
392:   PetscSpaceGetDimension(sp, &pdimfull);
393:   numelB = npoints * pdimfull * Nc;
394:   numelD = numelB * Nv;
395:   numelH = numelD * Nv;
396:   if (B || D || H) DMGetWorkArray(dm, numelB, MPIU_REAL, &sB);
397:   if (D || H) DMGetWorkArray(dm, numelD, MPIU_REAL, &sD);
398:   if (H) DMGetWorkArray(dm, numelH, MPIU_REAL, &sH);
399:   if (B)
400:     for (i = 0; i < numelB; ++i) B[i] = 0.;
401:   if (D)
402:     for (i = 0; i < numelD; ++i) D[i] = 0.;
403:   if (H)
404:     for (i = 0; i < numelH; ++i) H[i] = 0.;

406:   for (s = 0, offset = 0, ncoffset = 0; s < Ns; ++s) {
407:     PetscInt sNv, spdim, sNc, p;

409:     PetscSpaceGetNumVariables(sum->sumspaces[s], &sNv);
410:     PetscSpaceGetNumComponents(sum->sumspaces[s], &sNc);
411:     PetscSpaceGetDimension(sum->sumspaces[s], &spdim);
413:     if (s == 0 || !sum->uniform) PetscSpaceEvaluate(sum->sumspaces[s], npoints, points, sB, sD, sH);
414:     if (B || D || H) {
415:       for (p = 0; p < npoints; ++p) {
416:         PetscInt j;

418:         for (j = 0; j < spdim; ++j) {
419:           PetscInt c;

421:           for (c = 0; c < sNc; ++c) {
422:             PetscInt compoffset, BInd, sBInd;

424:             compoffset = concatenate ? c + ncoffset : c;
425:             BInd       = (p * pdimfull + j + offset) * Nc + compoffset;
426:             sBInd      = (p * spdim + j) * sNc + c;
427:             if (B) B[BInd] = sB[sBInd];
428:             if (D || H) {
429:               PetscInt v;

431:               for (v = 0; v < Nv; ++v) {
432:                 PetscInt DInd, sDInd;

434:                 DInd  = BInd * Nv + v;
435:                 sDInd = sBInd * Nv + v;
436:                 if (D) D[DInd] = sD[sDInd];
437:                 if (H) {
438:                   PetscInt v2;

440:                   for (v2 = 0; v2 < Nv; ++v2) {
441:                     PetscInt HInd, sHInd;

443:                     HInd    = DInd * Nv + v2;
444:                     sHInd   = sDInd * Nv + v2;
445:                     H[HInd] = sH[sHInd];
446:                   }
447:                 }
448:               }
449:             }
450:           }
451:         }
452:       }
453:     }
454:     offset += spdim;
455:     ncoffset += sNc;
456:   }

458:   if (H) DMRestoreWorkArray(dm, numelH, MPIU_REAL, &sH);
459:   if (D || H) DMRestoreWorkArray(dm, numelD, MPIU_REAL, &sD);
460:   if (B || D || H) DMRestoreWorkArray(dm, numelB, MPIU_REAL, &sB);
461:   return 0;
462: }

464: static PetscErrorCode PetscSpaceGetHeightSubspace_Sum(PetscSpace sp, PetscInt height, PetscSpace *subsp)
465: {
466:   PetscSpace_Sum *sum = (PetscSpace_Sum *)sp->data;
467:   PetscInt        Nc, dim, order;
468:   PetscBool       tensor;

470:   PetscSpaceGetNumComponents(sp, &Nc);
471:   PetscSpaceGetNumVariables(sp, &dim);
472:   PetscSpaceGetDegree(sp, &order, NULL);
473:   PetscSpacePolynomialGetTensor(sp, &tensor);
475:   if (!sum->heightsubspaces) PetscCalloc1(dim, &sum->heightsubspaces);
476:   if (height <= dim) {
477:     if (!sum->heightsubspaces[height - 1]) {
478:       PetscSpace  sub;
479:       const char *name;

481:       PetscSpaceCreate(PetscObjectComm((PetscObject)sp), &sub);
482:       PetscObjectGetName((PetscObject)sp, &name);
483:       PetscObjectSetName((PetscObject)sub, name);
484:       PetscSpaceSetType(sub, PETSCSPACESUM);
485:       PetscSpaceSumSetNumSubspaces(sub, sum->numSumSpaces);
486:       PetscSpaceSumSetConcatenate(sub, sum->concatenate);
487:       PetscSpaceSetNumComponents(sub, Nc);
488:       PetscSpaceSetNumVariables(sub, dim - height);
489:       for (PetscInt i = 0; i < sum->numSumSpaces; i++) {
490:         PetscSpace subh;

492:         PetscSpaceGetHeightSubspace(sum->sumspaces[i], height, &subh);
493:         PetscSpaceSumSetSubspace(sub, i, subh);
494:       }
495:       PetscSpaceSetUp(sub);
496:       sum->heightsubspaces[height - 1] = sub;
497:     }
498:     *subsp = sum->heightsubspaces[height - 1];
499:   } else {
500:     *subsp = NULL;
501:   }
502:   return 0;
503: }

505: static PetscErrorCode PetscSpaceInitialize_Sum(PetscSpace sp)
506: {
507:   sp->ops->setfromoptions    = PetscSpaceSetFromOptions_Sum;
508:   sp->ops->setup             = PetscSpaceSetUp_Sum;
509:   sp->ops->view              = PetscSpaceView_Sum;
510:   sp->ops->destroy           = PetscSpaceDestroy_Sum;
511:   sp->ops->getdimension      = PetscSpaceGetDimension_Sum;
512:   sp->ops->evaluate          = PetscSpaceEvaluate_Sum;
513:   sp->ops->getheightsubspace = PetscSpaceGetHeightSubspace_Sum;

515:   PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumGetNumSubspaces_C", PetscSpaceSumGetNumSubspaces_Sum);
516:   PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumSetNumSubspaces_C", PetscSpaceSumSetNumSubspaces_Sum);
517:   PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumGetSubspace_C", PetscSpaceSumGetSubspace_Sum);
518:   PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumSetSubspace_C", PetscSpaceSumSetSubspace_Sum);
519:   PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumGetConcatenate_C", PetscSpaceSumGetConcatenate_Sum);
520:   PetscObjectComposeFunction((PetscObject)sp, "PetscSpaceSumSetConcatenate_C", PetscSpaceSumSetConcatenate_Sum);
521:   return 0;
522: }

524: /*MC
525:   PETSCSPACESUM = "sum" - A `PetscSpace` object that encapsulates a sum of subspaces.
526:   That sum can either be direct or concatenate a concatenation. For example if A and B are spaces each with 2 components,
527:   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
528:   same number of variables.

530:   Level: intermediate

532: .seealso: `PetscSpace`, `PetscSpaceType`, `PetscSpaceCreate()`, `PetscSpaceSetType()`
533: M*/
534: PETSC_EXTERN PetscErrorCode PetscSpaceCreate_Sum(PetscSpace sp)
535: {
536:   PetscSpace_Sum *sum;

539:   PetscNew(&sum);
540:   sum->numSumSpaces = PETSC_DEFAULT;
541:   sp->data          = sum;
542:   PetscSpaceInitialize_Sum(sp);
543:   return 0;
544: }

546: PETSC_EXTERN PetscErrorCode PetscSpaceCreateSum(PetscInt numSubspaces, const PetscSpace subspaces[], PetscBool concatenate, PetscSpace *sumSpace)
547: {
548:   PetscInt i, Nv, Nc = 0;

550:   if (sumSpace) PetscSpaceDestroy(sumSpace);
551:   PetscSpaceCreate(PetscObjectComm((PetscObject)subspaces[0]), sumSpace);
552:   PetscSpaceSetType(*sumSpace, PETSCSPACESUM);
553:   PetscSpaceSumSetNumSubspaces(*sumSpace, numSubspaces);
554:   PetscSpaceSumSetConcatenate(*sumSpace, concatenate);
555:   for (i = 0; i < numSubspaces; ++i) {
556:     PetscInt sNc;

558:     PetscSpaceSumSetSubspace(*sumSpace, i, subspaces[i]);
559:     PetscSpaceGetNumComponents(subspaces[i], &sNc);
560:     if (concatenate) Nc += sNc;
561:     else Nc = sNc;
562:   }
563:   PetscSpaceGetNumVariables(subspaces[0], &Nv);
564:   PetscSpaceSetNumComponents(*sumSpace, Nc);
565:   PetscSpaceSetNumVariables(*sumSpace, Nv);
566:   PetscSpaceSetUp(*sumSpace);

568:   return 0;
569: }