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: PetscSpaceSumSetNumSubspaces(), PetscSpaceSetDegree(), PetscSpaceSetNumVariables()
     14: @*/
     15: PetscErrorCode PetscSpaceSumGetNumSubspaces(PetscSpace sp,PetscInt *numSumSpaces)
     16: {

     22:   PetscTryMethod(sp,"PetscSpaceSumGetNumSubspaces_C",(PetscSpace,PetscInt*),(sp,numSumSpaces));
     23:   return(0);
     24: }

     26: /*@
     27:   PetscSpaceSumSetNumSubspaces - Set the number of spaces in the sum

     29:   Input Parameters:
     30:   + sp  - the function space object
     31:   - numSumSpaces - the number of spaces

     33: Level: intermediate

     35: .seealso: PetscSpaceSumGetNumSubspaces(), PetscSpaceSetDegree(), PetscSpaceSetNumVariables()
     36: @*/
     37: PetscErrorCode PetscSpaceSumSetNumSubspaces(PetscSpace sp,PetscInt numSumSpaces)
     38: {

     43:   PetscTryMethod(sp,"PetscSpaceSumSetNumSubspaces_C",(PetscSpace,PetscInt),(sp,numSumSpaces));
     44:   return(0);
     45: }

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

     52:  Input Parameters:
     53:  . sp - the function space object

     55:  Output Parameters:
     56:  . concatenate - flag indicating whether subspaces are concatenated.

     58: Level: intermediate

     60: .seealso: PetscSpaceSumSetConcatenate()
     61: @*/
     62: PetscErrorCode PetscSpaceSumGetConcatenate(PetscSpace sp,PetscBool *concatenate)
     63: {

     68:   PetscTryMethod(sp,"PetscSpaceSumGetConcatenate_C",(PetscSpace,PetscBool*),(sp,concatenate));
     69:   return(0);
     70: }

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

     77:  Input Parameters:
     78:   + sp - the function space object
     79:   - concatenate - are subspaces concatenated components (true) or direct summands (false)

     81: Level: intermediate
     82: .seealso: PetscSpaceSumGetConcatenate()
     83: @*/
     84: PetscErrorCode PetscSpaceSumSetConcatenate(PetscSpace sp,PetscBool concatenate)
     85: {

     90:   PetscTryMethod(sp,"PetscSpaceSumSetConcatenate_C",(PetscSpace,PetscBool),(sp,concatenate));
     91:   return(0);
     92: }

     94: /*@
     95:   PetscSpaceSumGetSubspace - Get a space in the sum

     97:   Input Parameters:
     98:   + sp - the function space object
     99:   - s  - The space number

    101:   Output Parameter:
    102:   . subsp - the PetscSpace

    104: Level: intermediate

    106: .seealso: PetscSpaceSumSetSubspace(), PetscSpaceSetDegree(), PetscSpaceSetNumVariables()
    107: @*/
    108: PetscErrorCode PetscSpaceSumGetSubspace(PetscSpace sp,PetscInt s,PetscSpace *subsp)
    109: {

    115:   PetscTryMethod(sp,"PetscSpaceSumGetSubspace_C",(PetscSpace,PetscInt,PetscSpace*),(sp,s,subsp));
    116:   return(0);
    117: }

    119: /*@
    120:   PetscSpaceSumSetSubspace - Set a space in the sum

    122:   Input Parameters:
    123:   + sp    - the function space object
    124:   . s     - The space number
    125:   - subsp - the number of spaces

    127: Level: intermediate

    129: .seealso: PetscSpaceSumGetSubspace(), PetscSpaceSetDegree(), PetscSpaceSetNumVariables()
    130: @*/
    131: PetscErrorCode PetscSpaceSumSetSubspace(PetscSpace sp,PetscInt s,PetscSpace subsp)
    132: {

    138:   PetscTryMethod(sp,"PetscSpaceSumSetSubspace_C",(PetscSpace,PetscInt,PetscSpace),(sp,s,subsp));
    139:   return(0);
    140: }

    142: static PetscErrorCode PetscSpaceSumGetNumSubspaces_Sum(PetscSpace space,PetscInt *numSumSpaces)
    143: {
    144:   PetscSpace_Sum *sum = (PetscSpace_Sum*)space->data;

    147:   *numSumSpaces = sum->numSumSpaces;
    148:   return(0);
    149: }

    151: static PetscErrorCode PetscSpaceSumSetNumSubspaces_Sum(PetscSpace space,PetscInt numSumSpaces)
    152: {
    153:   PetscSpace_Sum *sum = (PetscSpace_Sum*)space->data;
    154:   PetscInt       Ns   = sum->numSumSpaces;

    158:   if (sum->setupCalled) SETERRQ(PetscObjectComm((PetscObject)space),PETSC_ERR_ARG_WRONGSTATE,
    159:                                 "Cannot change number of subspaces after setup called\n");
    160:   if (numSumSpaces == Ns) return(0);
    161:   if (Ns >= 0) {
    162:     PetscInt s;
    163:     for (s=0; s<Ns; ++s) {
    164:       PetscSpaceDestroy(&sum->sumspaces[s]);
    165:     }
    166:     PetscFree(sum->sumspaces);
    167:   }

    169:   Ns   = sum->numSumSpaces = numSumSpaces;
    170:   PetscCalloc1(Ns,&sum->sumspaces);
    171:   return(0);
    172: }

    174: static PetscErrorCode PetscSpaceSumGetConcatenate_Sum(PetscSpace sp,PetscBool *concatenate)
    175: {
    176:   PetscSpace_Sum *sum = (PetscSpace_Sum*)sp->data;

    179:   *concatenate = sum->concatenate;
    180:   return(0);
    181: }

    183: static PetscErrorCode PetscSpaceSumSetConcatenate_Sum(PetscSpace sp,PetscBool concatenate)
    184: {
    185:   PetscSpace_Sum *sum = (PetscSpace_Sum*)sp->data;

    188:   if (sum->setupCalled) SETERRQ(PetscObjectComm((PetscObject)sp),PETSC_ERR_ARG_WRONGSTATE,"Cannot change space concatenation after setup called.\n");

    190:   sum->concatenate = concatenate;
    191:   return(0);
    192: }

    194: static PetscErrorCode PetscSpaceSumGetSubspace_Sum(PetscSpace space,PetscInt s,PetscSpace *subspace)
    195: {
    196:   PetscSpace_Sum *sum = (PetscSpace_Sum*)space->data;
    197:   PetscInt       Ns   = sum->numSumSpaces;

    200:   if (Ns < 0) SETERRQ(PetscObjectComm((PetscObject)space),PETSC_ERR_ARG_WRONGSTATE,"Must call PetscSpaceSumSetNumSubspaces() first\n");
    201:   if (s<0 || s>=Ns) SETERRQ1 (PetscObjectComm((PetscObject)space),PETSC_ERR_ARG_OUTOFRANGE,"Invalid subspace number %D\n",subspace);

    203:   *subspace = sum->sumspaces[s];
    204:   return(0);
    205: }

    207: static PetscErrorCode PetscSpaceSumSetSubspace_Sum(PetscSpace space,PetscInt s,PetscSpace subspace)
    208: {
    209:   PetscSpace_Sum *sum = (PetscSpace_Sum*)space->data;
    210:   PetscInt       Ns   = sum->numSumSpaces;

    214:   if (sum->setupCalled) SETERRQ(PetscObjectComm((PetscObject)space),PETSC_ERR_ARG_WRONGSTATE,"Cannot change subspace after setup called\n");
    215:   if (Ns < 0) SETERRQ(PetscObjectComm((PetscObject)space),PETSC_ERR_ARG_WRONGSTATE,"Must call PetscSpaceSumSetNumSubspaces() first\n");
    216:   if (s < 0 || s >= Ns) SETERRQ1(PetscObjectComm((PetscObject)space),PETSC_ERR_ARG_OUTOFRANGE,"Invalid subspace number %D\n",subspace);

    218:   PetscObjectReference((PetscObject)subspace);
    219:   PetscSpaceDestroy(&sum->sumspaces[s]);
    220:   sum->sumspaces[s] = subspace;
    221:   return(0);
    222: }

    224: static PetscErrorCode PetscSpaceSetFromOptions_Sum(PetscOptionItems *PetscOptionsObject,PetscSpace sp)
    225: {
    226:   PetscSpace_Sum *sum = (PetscSpace_Sum*)sp->data;
    227:   PetscInt       Ns,Nc,Nv,deg,i;
    228:   PetscBool      concatenate = PETSC_TRUE;
    229:   const char     *prefix;

    233:   PetscSpaceGetNumVariables(sp,&Nv);
    234:   if (!Nv) return(0);
    235:   PetscSpaceGetNumComponents(sp,&Nc);
    236:   PetscSpaceSumGetNumSubspaces(sp,&Ns);
    237:   PetscSpaceGetDegree(sp,&deg,NULL);
    238:   Ns   = (Ns == PETSC_DEFAULT) ? 1 : Ns;

    240:   PetscOptionsHead(PetscOptionsObject,"PetscSpace sum options");
    241:   PetscOptionsBoundedInt("-petscspace_sum_spaces","The number of subspaces","PetscSpaceSumSetNumSubspaces",Ns,&Ns,NULL,0);
    242:   PetscOptionsBool("-petscspace_sum_concatenate","Subspaces are concatenated components of the final space","PetscSpaceSumSetFromOptions",
    243:                           concatenate,&concatenate,NULL);
    244:   PetscOptionsTail();

    246:   if (Ns < 0 || (Nv > 0 && Ns == 0)) SETERRQ1(PetscObjectComm((PetscObject)sp),PETSC_ERR_ARG_OUTOFRANGE,"Cannot have a sum space of %D spaces\n",Ns);
    247:   if (Ns != sum->numSumSpaces) {
    248:     PetscSpaceSumSetNumSubspaces(sp,Ns);
    249:   }
    250:   PetscObjectGetOptionsPrefix((PetscObject)sp,&prefix);
    251:   for (i=0; i<Ns; ++i) {
    252:     PetscInt   sNv;
    253:     PetscSpace subspace;

    255:     PetscSpaceSumGetSubspace(sp,i,&subspace);
    256:     if (!subspace) {
    257:       char subspacePrefix[256];

    259:       PetscSpaceCreate(PetscObjectComm((PetscObject)sp),&subspace);
    260:       PetscObjectSetOptionsPrefix((PetscObject)subspace,prefix);
    261:       PetscSNPrintf(subspacePrefix,256,"subspace%d_",i);
    262:       PetscObjectAppendOptionsPrefix((PetscObject)subspace,subspacePrefix);
    263:     } else {
    264:       PetscObjectReference((PetscObject)subspace);
    265:     }
    266:     PetscSpaceSetFromOptions(subspace);
    267:     PetscSpaceGetNumVariables(subspace,&sNv);
    268:     if (!sNv) SETERRQ1(PetscObjectComm(
    269:                          (PetscObject)sp),PETSC_ERR_ARG_WRONGSTATE,"Subspace %D has not been set properly, number of variables is 0.\n",i);
    270:     PetscSpaceSumSetSubspace(sp,i,subspace);
    271:     PetscSpaceDestroy(&subspace);
    272:   }
    273:   return(0);
    274: }

    276: static PetscErrorCode PetscSpaceSetUp_Sum(PetscSpace sp)
    277: {
    278:   PetscSpace_Sum *sum = (PetscSpace_Sum*)sp->data;
    279:   PetscBool      concatenate = PETSC_TRUE;
    280:   PetscInt       Nv,Ns,Nc,i,sum_Nc = 0,deg = PETSC_MIN_INT,maxDeg = PETSC_MIN_INT;

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

    286:   PetscSpaceGetNumVariables(sp,&Nv);
    287:   PetscSpaceGetNumComponents(sp,&Nc);
    288:   PetscSpaceSumGetNumSubspaces(sp,&Ns);
    289:   if (Ns == PETSC_DEFAULT) {
    290:     Ns   = 1;
    291:     PetscSpaceSumSetNumSubspaces(sp,Ns);
    292:   }
    293:   if (!Ns && Nv) SETERRQ(PetscObjectComm((PetscObject)sp),PETSC_ERR_ARG_OUTOFRANGE,"Cannot have zero subspaces\n");

    295:   for (i=0; i<Ns; ++i) {
    296:     PetscInt   sNv,sNc,iDeg,iMaxDeg;
    297:     PetscSpace si;

    299:     PetscSpaceSumGetSubspace(sp,i,&si);
    300:     PetscSpaceSetUp(si);
    301:     PetscSpaceGetNumVariables(si,&sNv);
    302:     if (sNv != Nv) SETERRQ3(PetscObjectComm((PetscObject)sp),PETSC_ERR_ARG_WRONGSTATE,"Subspace %D has %D variables, space has %D.\n",i,sNv,Nv);
    303:     PetscSpaceGetNumComponents(si,&sNc);
    304:     if (i == 0 && sNc == Nc) concatenate = PETSC_FALSE;
    305:     sum_Nc += sNc;
    306:     PetscSpaceSumGetSubspace(sp,i,&si);
    307:     PetscSpaceGetDegree(si,&iDeg,&iMaxDeg);
    308:     deg     = PetscMax(deg,iDeg);
    309:     maxDeg  = PetscMax(maxDeg,iMaxDeg);
    310:   }

    312:   if (concatenate) {
    313:     if (sum_Nc != Nc) {
    314:       SETERRQ2(PetscObjectComm((PetscObject)sp),PETSC_ERR_ARG_OUTOFRANGE,
    315:                "Total number of subspace components (%D) does not match number of target space components (%D).",sum_Nc,Nc);
    316:     }
    317:   } else {
    318:     if (sum_Nc != Ns*Nc) SETERRQ(PetscObjectComm(
    319:                                    (PetscObject)sp),PETSC_ERR_ARG_OUTOFRANGE,"Subspaces must have same number of components as the target space.");
    320:   }

    322:   sp->degree       = deg;
    323:   sp->maxDegree    = maxDeg;
    324:   sum->concatenate = concatenate;
    325:   sum->setupCalled = PETSC_TRUE;
    326:   return(0);
    327: }

    329: static PetscErrorCode PetscSpaceSumView_Ascii(PetscSpace sp,PetscViewer v)
    330: {
    331:   PetscSpace_Sum *sum = (PetscSpace_Sum*)sp->data;
    332:   PetscBool      concatenate = sum->concatenate;
    333:   PetscInt       i,Ns         = sum->numSumSpaces;

    337:   if (concatenate) {
    338:     PetscViewerASCIIPrintf(v,"Sum space of %D concatenated subspaces\n",Ns);
    339:   } else {
    340:     PetscViewerASCIIPrintf(v,"Sum space of %D subspaces\n",Ns);
    341:   }
    342:   for (i=0; i<Ns; ++i) {
    343:     PetscViewerASCIIPushTab(v);
    344:     PetscSpaceView(sum->sumspaces[i],v);
    345:     PetscViewerASCIIPopTab(v);
    346:   }
    347:   return(0);
    348: }

    350: static PetscErrorCode PetscSpaceView_Sum(PetscSpace sp,PetscViewer viewer)
    351: {
    352:   PetscBool      iascii;

    356:   PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);
    357:   if (iascii) {
    358:     PetscSpaceSumView_Ascii(sp,viewer);
    359:   }
    360:   return(0);
    361: }

    363: static PetscErrorCode PetscSpaceDestroy_Sum(PetscSpace sp)
    364: {
    365:   PetscSpace_Sum *sum = (PetscSpace_Sum*)sp->data;
    366:   PetscInt       i,Ns   = sum->numSumSpaces;

    370:   for (i=0; i<Ns; ++i) {
    371:     PetscSpaceDestroy(&sum->sumspaces[i]);
    372:   }
    373:   PetscFree(sum->sumspaces);
    374:   PetscObjectComposeFunction((PetscObject)sp,"PetscSpaceSumSetSubspace_C",NULL);
    375:   PetscObjectComposeFunction((PetscObject)sp,"PetscSpaceSumGetSubspace_C",NULL);
    376:   PetscObjectComposeFunction((PetscObject)sp,"PetscSpaceSumSetNumSubspaces_C",NULL);
    377:   PetscObjectComposeFunction((PetscObject)sp,"PetscSpaceSumGetNumSubspaces_C",NULL);
    378:   PetscFree(sum);
    379:   return(0);
    380: }

    382: static PetscErrorCode PetscSpaceGetDimension_Sum(PetscSpace sp,PetscInt *dim)
    383: {
    384:   PetscSpace_Sum *sum = (PetscSpace_Sum*)sp->data;
    385:   PetscInt       i,d = 0,Ns = sum->numSumSpaces;

    389:   PetscSpaceSetUp(sp);

    391:   for (i=0; i<Ns; ++i) {
    392:     PetscInt id;

    394:     PetscSpaceGetDimension(sum->sumspaces[i],&id);
    395:     d   += id;
    396:   }

    398:   *dim = d;
    399:   return(0);
    400: }

    402: static PetscErrorCode PetscSpaceEvaluate_Sum(PetscSpace sp,PetscInt npoints,const PetscReal points[],PetscReal B[],PetscReal D[],PetscReal H[])
    403: {
    404:   PetscSpace_Sum *sum = (PetscSpace_Sum*)sp->data;
    405:   PetscBool      concatenate = sum->concatenate;
    406:   DM             dm = sp->dm;
    407:   PetscInt       Nc = sp->Nc,Nv = sp->Nv,Ns = sum->numSumSpaces;
    408:   PetscInt       i,s,offset,ncoffset,pdimfull,numelB,numelD,numelH;
    409:   PetscReal      *sB = NULL,*sD = NULL,*sH = NULL;

    413:   if (!sum->setupCalled) {
    414:     PetscSpaceSetUp(sp);
    415:   }
    416:   PetscSpaceGetDimension(sp,&pdimfull);
    417:   numelB = npoints*pdimfull*Nc;
    418:   numelD = numelB*Nv;
    419:   numelH = numelD*Nv;
    420:   if (B || D || H) {
    421:     DMGetWorkArray(dm,numelB,MPIU_REAL,&sB);
    422:   }
    423:   if (D || H) {
    424:     DMGetWorkArray(dm,numelD,MPIU_REAL,&sD);
    425:   }
    426:   if (H) {
    427:     DMGetWorkArray(dm,numelH,MPIU_REAL,&sH);
    428:   }
    429:   if (B)
    430:     for (i=0; i<numelB; ++i) B[i] = 0.;
    431:   if (D)
    432:     for (i=0; i<numelD; ++i) D[i] = 0.;
    433:   if (H)
    434:     for (i=0; i<numelH; ++i) H[i] = 0.;

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

    439:     PetscSpaceGetNumVariables(sum->sumspaces[s],&sNv);
    440:     PetscSpaceGetNumComponents(sum->sumspaces[s],&sNc);
    441:     PetscSpaceGetDimension(sum->sumspaces[s],&spdim);
    442:     if (offset + spdim > pdimfull) SETERRQ(PetscObjectComm(
    443:                                              (PetscObject)sp),PETSC_ERR_ARG_OUTOFRANGE,"Subspace dimensions exceed target space dimension.\n");
    444:     PetscSpaceEvaluate(sum->sumspaces[s],npoints,points,sB,sD,sH);
    445:     if (B || D || H) {
    446:       for (p=0; p<npoints; ++p) {
    447:         PetscInt j;

    449:         for (j=0; j<spdim; ++j) {
    450:           PetscInt c;

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

    455:             compoffset = concatenate ? c+ncoffset : c;
    456:             BInd       = (p*pdimfull + j + offset)*Nc + compoffset;
    457:             sBInd      = (p*spdim + j)*sNc + c;
    458:             if (B) B[BInd] += sB[sBInd];
    459:             if (D || H) {
    460:               PetscInt v;

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

    465:                 DInd  = BInd*Nv + v;
    466:                 sDInd = sBInd*Nv + v;
    467:                 if (D) D[DInd] +=sD[sDInd];
    468:                 if (H) {
    469:                   PetscInt v2;

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

    474:                     HInd     = DInd*Nv + v2;
    475:                     sHInd    = sDInd*Nv + v2;
    476:                     H[HInd] += sH[sHInd];
    477:                   }
    478:                 }
    479:               }
    480:             }
    481:           }
    482:         }
    483:       }
    484:     }
    485:     offset   += spdim;
    486:     ncoffset += sNc;
    487:   }

    489:   if (H) {
    490:     DMRestoreWorkArray(dm,numelH,MPIU_REAL,&sH);
    491:   }
    492:   if (D || H) {
    493:     DMRestoreWorkArray(dm,numelD,MPIU_REAL,&sD);
    494:   }
    495:   if (B || D || H) {
    496:     DMRestoreWorkArray(dm,numelB,MPIU_REAL,&sB);
    497:   }
    498:   return(0);
    499: }

    501: static PetscErrorCode PetscSpaceInitialize_Sum(PetscSpace sp)
    502: {

    506:   sp->ops->setfromoptions = PetscSpaceSetFromOptions_Sum;
    507:   sp->ops->setup          = PetscSpaceSetUp_Sum;
    508:   sp->ops->view           = PetscSpaceView_Sum;
    509:   sp->ops->destroy        = PetscSpaceDestroy_Sum;
    510:   sp->ops->getdimension   = PetscSpaceGetDimension_Sum;
    511:   sp->ops->evaluate       = PetscSpaceEvaluate_Sum;

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

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

    528: Level: intermediate

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

    539:   PetscNewLog(sp,&sum);
    540:   sp->data = sum;
    541:   PetscSpaceInitialize_Sum(sp);
    542:   return(0);
    543: }

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

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

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

    571:   return(0);
    572: }