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 PetscBool found = PETSC_FALSE; 5 6 PetscFunctionBegin; 7 *isplit = ts->tsrhssplit; 8 /* look up the split */ 9 while (*isplit) { 10 PetscCall(PetscStrcmp((*isplit)->splitname, splitname, &found)); 11 if (found) break; 12 *isplit = (*isplit)->next; 13 } 14 PetscFunctionReturn(0); 15 } 16 17 /*@C 18 TSRHSSplitSetIS - Set the index set for the specified split 19 20 Logically Collective on TS 21 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 26 27 Level: intermediate 28 29 .seealso: `TSRHSSplitGetIS()` 30 31 @*/ 32 PetscErrorCode TSRHSSplitSetIS(TS ts, const char splitname[], IS is) { 33 TS_RHSSplitLink newsplit, next = ts->tsrhssplit; 34 char prefix[128]; 35 36 PetscFunctionBegin; 37 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 38 PetscValidHeaderSpecific(is, IS_CLASSID, 3); 39 40 PetscCall(PetscNew(&newsplit)); 41 if (splitname) { 42 PetscCall(PetscStrallocpy(splitname, &newsplit->splitname)); 43 } else { 44 PetscCall(PetscMalloc1(8, &newsplit->splitname)); 45 PetscCall(PetscSNPrintf(newsplit->splitname, 7, "%" PetscInt_FMT, ts->num_rhs_splits)); 46 } 47 PetscCall(PetscObjectReference((PetscObject)is)); 48 newsplit->is = is; 49 PetscCall(TSCreate(PetscObjectComm((PetscObject)ts), &newsplit->ts)); 50 51 PetscCall(PetscObjectIncrementTabLevel((PetscObject)newsplit->ts, (PetscObject)ts, 1)); 52 PetscCall(PetscSNPrintf(prefix, sizeof(prefix), "%srhsplit_%s_", ((PetscObject)ts)->prefix ? ((PetscObject)ts)->prefix : "", newsplit->splitname)); 53 PetscCall(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 PetscFunctionReturn(0); 61 } 62 63 /*@C 64 TSRHSSplitGetIS - Retrieves the elements for a split as an IS 65 66 Logically Collective on TS 67 68 Input Parameters: 69 + ts - the TS context obtained from TSCreate() 70 - splitname - name of this split 71 72 Output Parameters: 73 - is - the index set for part of the solution vector 74 75 Level: intermediate 76 77 .seealso: `TSRHSSplitSetIS()` 78 79 @*/ 80 PetscErrorCode TSRHSSplitGetIS(TS ts, const char splitname[], IS *is) { 81 TS_RHSSplitLink isplit; 82 83 PetscFunctionBegin; 84 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 85 *is = NULL; 86 /* look up the split */ 87 PetscCall(TSRHSSplitGetRHSSplit(ts, splitname, &isplit)); 88 if (isplit) *is = isplit->is; 89 PetscFunctionReturn(0); 90 } 91 92 /*@C 93 TSRHSSplitSetRHSFunction - Set the split right-hand-side functions. 94 95 Logically Collective on TS 96 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) 103 104 Calling sequence of fun: 105 $ rhsfunc(TS ts,PetscReal t,Vec u,Vec f,ctx); 106 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) 111 112 Level: beginner 113 114 @*/ 115 PetscErrorCode TSRHSSplitSetRHSFunction(TS ts, const char splitname[], Vec r, TSRHSFunction rhsfunc, void *ctx) { 116 TS_RHSSplitLink isplit; 117 DM dmc; 118 Vec subvec, ralloc = NULL; 119 120 PetscFunctionBegin; 121 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 122 if (r) PetscValidHeaderSpecific(r, VEC_CLASSID, 3); 123 124 /* look up the split */ 125 PetscCall(TSRHSSplitGetRHSSplit(ts, splitname, &isplit)); 126 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); 127 128 if (!r && ts->vec_sol) { 129 PetscCall(VecGetSubVector(ts->vec_sol, isplit->is, &subvec)); 130 PetscCall(VecDuplicate(subvec, &ralloc)); 131 r = ralloc; 132 PetscCall(VecRestoreSubVector(ts->vec_sol, isplit->is, &subvec)); 133 } 134 135 if (ts->dm) { 136 PetscInt dim; 137 138 PetscCall(DMGetDimension(ts->dm, &dim)); 139 if (dim != -1) { 140 PetscCall(DMClone(ts->dm, &dmc)); 141 PetscCall(TSSetDM(isplit->ts, dmc)); 142 PetscCall(DMDestroy(&dmc)); 143 } 144 } 145 146 PetscCall(TSSetRHSFunction(isplit->ts, r, rhsfunc, ctx)); 147 PetscCall(VecDestroy(&ralloc)); 148 PetscFunctionReturn(0); 149 } 150 151 /*@C 152 TSRHSSplitGetSubTS - Get the sub-TS by split name. 153 154 Logically Collective on TS 155 156 Input Parameter: 157 . ts - the TS context obtained from TSCreate() 158 159 Output Parameters: 160 + splitname - the number of the split 161 - subts - the array of TS contexts 162 163 Level: advanced 164 165 .seealso: `TSGetRHSSplitFunction()` 166 @*/ 167 PetscErrorCode TSRHSSplitGetSubTS(TS ts, const char splitname[], TS *subts) { 168 TS_RHSSplitLink isplit; 169 170 PetscFunctionBegin; 171 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 172 PetscValidPointer(subts, 3); 173 *subts = NULL; 174 /* look up the split */ 175 PetscCall(TSRHSSplitGetRHSSplit(ts, splitname, &isplit)); 176 if (isplit) *subts = isplit->ts; 177 PetscFunctionReturn(0); 178 } 179 180 /*@C 181 TSRHSSplitGetSubTSs - Get an array of all sub-TS contexts. 182 183 Logically Collective on TS 184 185 Input Parameter: 186 . ts - the TS context obtained from TSCreate() 187 188 Output Parameters: 189 + n - the number of splits 190 - subksp - the array of TS contexts 191 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). 195 196 Level: advanced 197 198 .seealso: `TSGetRHSSplitFunction()` 199 @*/ 200 PetscErrorCode TSRHSSplitGetSubTSs(TS ts, PetscInt *n, TS *subts[]) { 201 TS_RHSSplitLink ilink = ts->tsrhssplit; 202 PetscInt i = 0; 203 204 PetscFunctionBegin; 205 PetscValidHeaderSpecific(ts, TS_CLASSID, 1); 206 if (subts) { 207 PetscCall(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 PetscFunctionReturn(0); 215 } 216