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