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, "%" PetscInt_FMT, 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:   PetscSNPrintf(prefix, sizeof(prefix), "%srhsplit_%s_", ((PetscObject)ts)->prefix ? ((PetscObject)ts)->prefix : "", newsplit->splitname);
 53:   TSSetOptionsPrefix(newsplit->ts, prefix);
 54:   if (!next) ts->tsrhssplit = newsplit;
 55:   else {
 56:     while (next->next) next = next->next;
 57:     next->next = newsplit;
 58:   }
 59:   ts->num_rhs_splits++;
 60:   return 0;
 61: }

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

 66:    Logically Collective on TS

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

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

 75:    Level: intermediate

 77: .seealso: `TSRHSSplitSetIS()`

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

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

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

 95:    Logically Collective on TS

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

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

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

112:  Level: beginner

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


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

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

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

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

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

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

154:    Logically Collective on TS

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

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

163:    Level: advanced

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

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

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

183:    Logically Collective on TS

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

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

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

196:    Level: advanced

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

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