1 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/ 2 #include <petscdm.h> 3 static PetscErrorCode TSRHSSplitGetRHSSplit(TS ts,const char splitname[],TS_RHSSplitLink *isplit) 4 { 5 PetscBool found = PETSC_FALSE; 6 7 PetscFunctionBegin; 8 *isplit = ts->tsrhssplit; 9 /* look up the split */ 10 while (*isplit) { 11 PetscCall(PetscStrcmp((*isplit)->splitname,splitname,&found)); 12 if (found) break; 13 *isplit = (*isplit)->next; 14 } 15 PetscFunctionReturn(0); 16 } 17 18 /*@C 19 TSRHSSplitSetIS - Set the index set for the specified split 20 21 Logically Collective on TS 22 23 Input Parameters: 24 + ts - the TS context obtained from TSCreate() 25 . splitname - name of this split, if NULL the number of the split is used 26 - is - the index set for part of the solution vector 27 28 Level: intermediate 29 30 .seealso: TSRHSSplitGetIS() 31 32 @*/ 33 PetscErrorCode TSRHSSplitSetIS(TS ts,const char splitname[],IS is) 34 { 35 TS_RHSSplitLink newsplit,next = ts->tsrhssplit; 36 char prefix[128]; 37 38 PetscFunctionBegin; 39 PetscValidHeaderSpecific(ts,TS_CLASSID,1); 40 PetscValidHeaderSpecific(is,IS_CLASSID,3); 41 42 PetscCall(PetscNew(&newsplit)); 43 if (splitname) { 44 PetscCall(PetscStrallocpy(splitname,&newsplit->splitname)); 45 } else { 46 PetscCall(PetscMalloc1(8,&newsplit->splitname)); 47 PetscCall(PetscSNPrintf(newsplit->splitname,7,"%D",ts->num_rhs_splits)); 48 } 49 PetscCall(PetscObjectReference((PetscObject)is)); 50 newsplit->is = is; 51 PetscCall(TSCreate(PetscObjectComm((PetscObject)ts),&newsplit->ts)); 52 53 PetscCall(PetscObjectIncrementTabLevel((PetscObject)newsplit->ts,(PetscObject)ts,1)); 54 PetscCall(PetscLogObjectParent((PetscObject)ts,(PetscObject)newsplit->ts)); 55 PetscCall(PetscSNPrintf(prefix,sizeof(prefix),"%srhsplit_%s_",((PetscObject)ts)->prefix ? ((PetscObject)ts)->prefix : "",newsplit->splitname)); 56 PetscCall(TSSetOptionsPrefix(newsplit->ts,prefix)); 57 if (!next) ts->tsrhssplit = newsplit; 58 else { 59 while (next->next) next = next->next; 60 next->next = newsplit; 61 } 62 ts->num_rhs_splits++; 63 PetscFunctionReturn(0); 64 } 65 66 /*@C 67 TSRHSSplitGetIS - Retrieves the elements for a split as an IS 68 69 Logically Collective on TS 70 71 Input Parameters: 72 + ts - the TS context obtained from TSCreate() 73 - splitname - name of this split 74 75 Output Parameters: 76 - is - the index set for part of the solution vector 77 78 Level: intermediate 79 80 .seealso: TSRHSSplitSetIS() 81 82 @*/ 83 PetscErrorCode TSRHSSplitGetIS(TS ts,const char splitname[],IS *is) 84 { 85 TS_RHSSplitLink isplit; 86 87 PetscFunctionBegin; 88 PetscValidHeaderSpecific(ts,TS_CLASSID,1); 89 *is = NULL; 90 /* look up the split */ 91 PetscCall(TSRHSSplitGetRHSSplit(ts,splitname,&isplit)); 92 if (isplit) *is = isplit->is; 93 PetscFunctionReturn(0); 94 } 95 96 /*@C 97 TSRHSSplitSetRHSFunction - Set the split right-hand-side functions. 98 99 Logically Collective on TS 100 101 Input Parameters: 102 + ts - the TS context obtained from TSCreate() 103 . splitname - name of this split 104 . r - vector to hold the residual (or NULL to have it created internally) 105 . rhsfunc - the RHS function evaluation routine 106 - ctx - user-defined context for private data for the split function evaluation routine (may be NULL) 107 108 Calling sequence of fun: 109 $ rhsfunc(TS ts,PetscReal t,Vec u,Vec f,ctx); 110 111 + t - time at step/stage being solved 112 . u - state vector 113 . f - function vector 114 - ctx - [optional] user-defined context for matrix evaluation routine (may be NULL) 115 116 Level: beginner 117 118 @*/ 119 PetscErrorCode TSRHSSplitSetRHSFunction(TS ts,const char splitname[],Vec r,TSRHSFunction rhsfunc,void *ctx) 120 { 121 TS_RHSSplitLink isplit; 122 DM dmc; 123 Vec subvec,ralloc = NULL; 124 125 PetscFunctionBegin; 126 PetscValidHeaderSpecific(ts,TS_CLASSID,1); 127 if (r) PetscValidHeaderSpecific(r,VEC_CLASSID,3); 128 129 /* look up the split */ 130 PetscCall(TSRHSSplitGetRHSSplit(ts,splitname,&isplit)); 131 PetscCheck(isplit,PETSC_COMM_SELF,PETSC_ERR_USER,"The split %s is not created, check the split name or call TSRHSSplitSetIS() to create one",splitname); 132 133 if (!r && ts->vec_sol) { 134 PetscCall(VecGetSubVector(ts->vec_sol,isplit->is,&subvec)); 135 PetscCall(VecDuplicate(subvec,&ralloc)); 136 r = ralloc; 137 PetscCall(VecRestoreSubVector(ts->vec_sol,isplit->is,&subvec)); 138 } 139 140 if (ts->dm) { 141 PetscInt dim; 142 143 PetscCall(DMGetDimension(ts->dm, &dim)); 144 if (dim != -1) { 145 PetscCall(DMClone(ts->dm, &dmc)); 146 PetscCall(TSSetDM(isplit->ts, dmc)); 147 PetscCall(DMDestroy(&dmc)); 148 } 149 } 150 151 PetscCall(TSSetRHSFunction(isplit->ts,r,rhsfunc,ctx)); 152 PetscCall(VecDestroy(&ralloc)); 153 PetscFunctionReturn(0); 154 } 155 156 /*@C 157 TSRHSSplitGetSubTS - Get the sub-TS by split name. 158 159 Logically Collective on TS 160 161 Input Parameter: 162 . ts - the TS context obtained from TSCreate() 163 164 Output Parameters: 165 + splitname - the number of the split 166 - subts - the array of TS contexts 167 168 Level: advanced 169 170 .seealso: TSGetRHSSplitFunction() 171 @*/ 172 PetscErrorCode TSRHSSplitGetSubTS(TS ts,const char splitname[],TS *subts) 173 { 174 TS_RHSSplitLink isplit; 175 176 PetscFunctionBegin; 177 PetscValidHeaderSpecific(ts,TS_CLASSID,1); 178 PetscValidPointer(subts,3); 179 *subts = NULL; 180 /* look up the split */ 181 PetscCall(TSRHSSplitGetRHSSplit(ts,splitname,&isplit)); 182 if (isplit) *subts = isplit->ts; 183 PetscFunctionReturn(0); 184 } 185 186 /*@C 187 TSRHSSplitGetSubTSs - Get an array of all sub-TS contexts. 188 189 Logically Collective on TS 190 191 Input Parameter: 192 . ts - the TS context obtained from TSCreate() 193 194 Output Parameters: 195 + n - the number of splits 196 - subksp - the array of TS contexts 197 198 Note: 199 After TSRHSSplitGetSubTS() the array of TSs is to be freed by the user with PetscFree() 200 (not the TS just the array that contains them). 201 202 Level: advanced 203 204 .seealso: TSGetRHSSplitFunction() 205 @*/ 206 PetscErrorCode TSRHSSplitGetSubTSs(TS ts,PetscInt *n,TS *subts[]) 207 { 208 TS_RHSSplitLink ilink = ts->tsrhssplit; 209 PetscInt i = 0; 210 211 PetscFunctionBegin; 212 PetscValidHeaderSpecific(ts,TS_CLASSID,1); 213 if (subts) { 214 PetscCall(PetscMalloc1(ts->num_rhs_splits,subts)); 215 while (ilink) { 216 (*subts)[i++] = ilink->ts; 217 ilink = ilink->next; 218 } 219 } 220 if (n) *n = ts->num_rhs_splits; 221 PetscFunctionReturn(0); 222 } 223