Actual source code: tsrhssplit.c

  1: #include <petsc/private/tsimpl.h>
  2: #include <petscdm.h>
  3: static PetscErrorCode TSRHSSplitGetRHSSplit(TS ts,const char splitname[],TS_RHSSplitLink *isplit)
  4: {
  5:   PetscBool       found = PETSC_FALSE;

  7:   *isplit = ts->tsrhssplit;
  8:   /* look up the split */
  9:   while (*isplit) {
 10:     PetscStrcmp((*isplit)->splitname,splitname,&found);
 11:     if (found) break;
 12:     *isplit = (*isplit)->next;
 13:   }
 14:   return 0;
 15: }

 17: /*@C
 18:    TSRHSSplitSetIS - Set the index set for the specified split

 20:    Logically Collective on TS

 22:    Input Parameters:
 23: +  ts        - the TS context obtained from TSCreate()
 24: .  splitname - name of this split, if NULL the number of the split is used
 25: -  is        - the index set for part of the solution vector

 27:    Level: intermediate

 29: .seealso: TSRHSSplitGetIS()

 31: @*/
 32: PetscErrorCode TSRHSSplitSetIS(TS ts,const char splitname[],IS is)
 33: {
 34:   TS_RHSSplitLink newsplit,next = ts->tsrhssplit;
 35:   char            prefix[128];


 40:   PetscNew(&newsplit);
 41:   if (splitname) {
 42:     PetscStrallocpy(splitname,&newsplit->splitname);
 43:   } else {
 44:     PetscMalloc1(8,&newsplit->splitname);
 45:     PetscSNPrintf(newsplit->splitname,7,"%D",ts->num_rhs_splits);
 46:   }
 47:   PetscObjectReference((PetscObject)is);
 48:   newsplit->is = is;
 49:   TSCreate(PetscObjectComm((PetscObject)ts),&newsplit->ts);

 51:   PetscObjectIncrementTabLevel((PetscObject)newsplit->ts,(PetscObject)ts,1);
 52:   PetscLogObjectParent((PetscObject)ts,(PetscObject)newsplit->ts);
 53:   PetscSNPrintf(prefix,sizeof(prefix),"%srhsplit_%s_",((PetscObject)ts)->prefix ? ((PetscObject)ts)->prefix : "",newsplit->splitname);
 54:   TSSetOptionsPrefix(newsplit->ts,prefix);
 55:   if (!next) ts->tsrhssplit = newsplit;
 56:   else {
 57:     while (next->next) next = next->next;
 58:     next->next = newsplit;
 59:   }
 60:   ts->num_rhs_splits++;
 61:   return 0;
 62: }

 64: /*@C
 65:    TSRHSSplitGetIS - Retrieves the elements for a split as an IS

 67:    Logically Collective on TS

 69:    Input Parameters:
 70: +  ts        - the TS context obtained from TSCreate()
 71: -  splitname - name of this split

 73:    Output Parameters:
 74: -  is        - the index set for part of the solution vector

 76:    Level: intermediate

 78: .seealso: TSRHSSplitSetIS()

 80: @*/
 81: PetscErrorCode TSRHSSplitGetIS(TS ts,const char splitname[],IS *is)
 82: {
 83:   TS_RHSSplitLink isplit;

 86:   *is = NULL;
 87:   /* look up the split */
 88:   TSRHSSplitGetRHSSplit(ts,splitname,&isplit);
 89:   if (isplit) *is = isplit->is;
 90:   return 0;
 91: }

 93: /*@C
 94:    TSRHSSplitSetRHSFunction - Set the split right-hand-side functions.

 96:    Logically Collective on TS

 98:    Input Parameters:
 99: +  ts        - the TS context obtained from TSCreate()
100: .  splitname - name of this split
101: .  r         - vector to hold the residual (or NULL to have it created internally)
102: .  rhsfunc   - the RHS function evaluation routine
103: -  ctx       - user-defined context for private data for the split function evaluation routine (may be NULL)

105:  Calling sequence of fun:
106: $  rhsfunc(TS ts,PetscReal t,Vec u,Vec f,ctx);

108: +  t    - time at step/stage being solved
109: .  u    - state vector
110: .  f    - function vector
111: -  ctx  - [optional] user-defined context for matrix evaluation routine (may be NULL)

113:  Level: beginner

115: @*/
116: PetscErrorCode TSRHSSplitSetRHSFunction(TS ts,const char splitname[],Vec r,TSRHSFunction rhsfunc,void *ctx)
117: {
118:   TS_RHSSplitLink isplit;
119:   DM              dmc;
120:   Vec             subvec,ralloc = NULL;


125:   /* look up the split */
126:   TSRHSSplitGetRHSSplit(ts,splitname,&isplit);

129:   if (!r && ts->vec_sol) {
130:     VecGetSubVector(ts->vec_sol,isplit->is,&subvec);
131:     VecDuplicate(subvec,&ralloc);
132:     r    = ralloc;
133:     VecRestoreSubVector(ts->vec_sol,isplit->is,&subvec);
134:   }

136:   if (ts->dm) {
137:     PetscInt dim;

139:     DMGetDimension(ts->dm, &dim);
140:     if (dim != -1) {
141:       DMClone(ts->dm, &dmc);
142:       TSSetDM(isplit->ts, dmc);
143:       DMDestroy(&dmc);
144:     }
145:   }

147:   TSSetRHSFunction(isplit->ts,r,rhsfunc,ctx);
148:   VecDestroy(&ralloc);
149:   return 0;
150: }

152: /*@C
153:    TSRHSSplitGetSubTS - Get the sub-TS by split name.

155:    Logically Collective on TS

157:    Input Parameter:
158: .  ts - the TS context obtained from TSCreate()

160:    Output Parameters:
161: +  splitname - the number of the split
162: -  subts - the array of TS contexts

164:    Level: advanced

166: .seealso: TSGetRHSSplitFunction()
167: @*/
168: PetscErrorCode TSRHSSplitGetSubTS(TS ts,const char splitname[],TS *subts)
169: {
170:   TS_RHSSplitLink isplit;

174:   *subts = NULL;
175:   /* look up the split */
176:   TSRHSSplitGetRHSSplit(ts,splitname,&isplit);
177:   if (isplit) *subts = isplit->ts;
178:   return 0;
179: }

181: /*@C
182:    TSRHSSplitGetSubTSs - Get an array of all sub-TS contexts.

184:    Logically Collective on TS

186:    Input Parameter:
187: .  ts - the TS context obtained from TSCreate()

189:    Output Parameters:
190: +  n - the number of splits
191: -  subksp - the array of TS contexts

193:    Note:
194:    After TSRHSSplitGetSubTS() the array of TSs is to be freed by the user with PetscFree()
195:    (not the TS just the array that contains them).

197:    Level: advanced

199: .seealso: TSGetRHSSplitFunction()
200: @*/
201: PetscErrorCode TSRHSSplitGetSubTSs(TS ts,PetscInt *n,TS *subts[])
202: {
203:   TS_RHSSplitLink ilink = ts->tsrhssplit;
204:   PetscInt        i = 0;

207:   if (subts) {
208:     PetscMalloc1(ts->num_rhs_splits,subts);
209:     while (ilink) {
210:       (*subts)[i++] = ilink->ts;
211:       ilink = ilink->next;
212:     }
213:   }
214:   if (n) *n = ts->num_rhs_splits;
215:   return 0;
216: }