Actual source code: tsrhssplit.c

petsc-3.14.6 2021-03-30
Report Typos and Errors
  1: #include <petsc/private/tsimpl.h>

  3: static PetscErrorCode TSRHSSplitGetRHSSplit(TS ts,const char splitname[],TS_RHSSplitLink *isplit)
  4: {
  5:   PetscBool       found = PETSC_FALSE;
  6:   PetscErrorCode  ierr;

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

 19: /*@C
 20:    TSRHSSplitSetIS - Set the index set for the specified split

 22:    Logically Collective on TS

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

 29:    Level: intermediate

 31: .seealso: TSRHSSplitGetIS()

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


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

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

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

 71:    Logically Collective on TS

 73:    Input Parameters:
 74: +  ts        - the TS context obtained from TSCreate()
 75: -  splitname - name of this split

 77:    Output Parameters:
 78: -  is        - the index set for part of the solution vector

 80:    Level: intermediate

 82: .seealso: TSRHSSplitSetIS()

 84: @*/
 85: PetscErrorCode TSRHSSplitGetIS(TS ts,const char splitname[],IS *is)
 86: {
 87:   TS_RHSSplitLink isplit;
 88:   PetscErrorCode  ierr;

 92:   *is = NULL;
 93:   /* look up the split */
 94:   TSRHSSplitGetRHSSplit(ts,splitname,&isplit);
 95:   if (isplit) *is = isplit->is;
 96:   return(0);
 97: }

 99: /*@C
100:    TSRHSSplitSetRHSFunction - Set the split right-hand-side functions.

102:    Logically Collective on TS

104:    Input Parameters:
105: +  ts        - the TS context obtained from TSCreate()
106: .  splitname - name of this split
107: .  r         - vector to hold the residual (or NULL to have it created internally)
108: .  rhsfunc   - the RHS function evaluation routine
109: -  ctx       - user-defined context for private data for the split function evaluation routine (may be NULL)

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

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

119:  Level: beginner

121: @*/
122: PetscErrorCode TSRHSSplitSetRHSFunction(TS ts,const char splitname[],Vec r,TSRHSFunction rhsfunc,void *ctx)
123: {
124:   TS_RHSSplitLink isplit;
125:   Vec             subvec,ralloc = NULL;
126:   PetscErrorCode  ierr;


132:   /* look up the split */
133:   TSRHSSplitGetRHSSplit(ts,splitname,&isplit);
134:   if (!isplit) SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_USER,"The split %s is not created, check the split name or call TSRHSSplitSetIS() to create one",splitname);

136:   if (!r && ts->vec_sol) {
137:     VecGetSubVector(ts->vec_sol,isplit->is,&subvec);
138:     VecDuplicate(subvec,&ralloc);
139:     r    = ralloc;
140:     VecRestoreSubVector(ts->vec_sol,isplit->is,&subvec);
141:   }
142:   TSSetRHSFunction(isplit->ts,r,rhsfunc,ctx);
143:   VecDestroy(&ralloc);
144:   return(0);
145: }

147: /*@C
148:    TSRHSSplitGetSubTS - Get the sub-TS by split name.

150:    Logically Collective on TS

152:    Output Parameters:
153: +  splitname - the number of the split
154: -  subts - the array of TS contexts

156:    Level: advanced

158: .seealso: TSGetRHSSplitFunction()
159: @*/
160: PetscErrorCode TSRHSSplitGetSubTS(TS ts,const char splitname[],TS *subts)
161: {
162:   TS_RHSSplitLink isplit;
163:   PetscErrorCode  ierr;

168:   *subts = NULL;
169:   /* look up the split */
170:   TSRHSSplitGetRHSSplit(ts,splitname,&isplit);
171:   if (isplit) *subts = isplit->ts;
172:   return(0);
173: }

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

178:    Logically Collective on TS

180:    Output Parameters:
181: +  n - the number of splits
182: -  subksp - the array of TS contexts

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

188:    Level: advanced

190: .seealso: TSGetRHSSplitFunction()
191: @*/
192: PetscErrorCode TSRHSSplitGetSubTSs(TS ts,PetscInt *n,TS *subts[])
193: {
194:   TS_RHSSplitLink ilink = ts->tsrhssplit;
195:   PetscInt        i = 0;
196:   PetscErrorCode  ierr;

200:   if (subts) {
201:     PetscMalloc1(ts->num_rhs_splits,subts);
202:     while (ilink) {
203:       (*subts)[i++] = ilink->ts;
204:       ilink = ilink->next;
205:     }
206:   }
207:   if (n) *n = ts->num_rhs_splits;
208:   return(0);
209: }