1 #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/
2
TSRHSSplitGetRHSSplit(TS ts,const char splitname[],TS_RHSSplitLink * isplit)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 @*/
TSRHSSplitSetIS(TS ts,const char splitname[],IS is)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
55 ierr = PetscObjectIncrementTabLevel((PetscObject)newsplit->ts,(PetscObject)ts,1);CHKERRQ(ierr);
56 ierr = PetscLogObjectParent((PetscObject)ts,(PetscObject)newsplit->ts);CHKERRQ(ierr);
57 ierr = PetscSNPrintf(prefix,sizeof(prefix),"%srhsplit_%s_",((PetscObject)ts)->prefix ? ((PetscObject)ts)->prefix : "",newsplit->splitname);CHKERRQ(ierr);
58 ierr = TSSetOptionsPrefix(newsplit->ts,prefix);CHKERRQ(ierr);
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 PetscFunctionReturn(0);
66 }
67
68 /*@C
69 TSRHSSplitGetIS - Retrieves the elements for a split as an IS
70
71 Logically Collective on TS
72
73 Input Parameters:
74 + ts - the TS context obtained from TSCreate()
75 - splitname - name of this split
76
77 Output Parameters:
78 - is - the index set for part of the solution vector
79
80 Level: intermediate
81
82 .seealso: TSRHSSplitSetIS()
83
84 @*/
TSRHSSplitGetIS(TS ts,const char splitname[],IS * is)85 PetscErrorCode TSRHSSplitGetIS(TS ts,const char splitname[],IS *is)
86 {
87 TS_RHSSplitLink isplit;
88 PetscErrorCode ierr;
89
90 PetscFunctionBegin;
91 PetscValidHeaderSpecific(ts,TS_CLASSID,1);
92 *is = NULL;
93 /* look up the split */
94 ierr = TSRHSSplitGetRHSSplit(ts,splitname,&isplit);CHKERRQ(ierr);
95 if (isplit) *is = isplit->is;
96 PetscFunctionReturn(0);
97 }
98
99 /*@C
100 TSRHSSplitSetRHSFunction - Set the split right-hand-side functions.
101
102 Logically Collective on TS
103
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)
110
111 Calling sequence of fun:
112 $ rhsfunc(TS ts,PetscReal t,Vec u,Vec f,ctx);
113
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)
118
119 Level: beginner
120
121 @*/
TSRHSSplitSetRHSFunction(TS ts,const char splitname[],Vec r,TSRHSFunction rhsfunc,void * ctx)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;
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 ierr = TSSetRHSFunction(isplit->ts,r,rhsfunc,ctx);CHKERRQ(ierr);
143 ierr = VecDestroy(&ralloc);CHKERRQ(ierr);
144 PetscFunctionReturn(0);
145 }
146
147 /*@C
148 TSRHSSplitGetSubTS - Get the sub-TS by split name.
149
150 Logically Collective on TS
151
152 Output Parameters:
153 + splitname - the number of the split
154 - subts - the array of TS contexts
155
156 Level: advanced
157
158 .seealso: TSGetRHSSplitFunction()
159 @*/
TSRHSSplitGetSubTS(TS ts,const char splitname[],TS * subts)160 PetscErrorCode TSRHSSplitGetSubTS(TS ts,const char splitname[],TS *subts)
161 {
162 TS_RHSSplitLink isplit;
163 PetscErrorCode ierr;
164
165 PetscFunctionBegin;
166 PetscValidHeaderSpecific(ts,TS_CLASSID,1);
167 PetscValidPointer(subts,3);
168 *subts = NULL;
169 /* look up the split */
170 ierr = TSRHSSplitGetRHSSplit(ts,splitname,&isplit);CHKERRQ(ierr);
171 if (isplit) *subts = isplit->ts;
172 PetscFunctionReturn(0);
173 }
174
175 /*@C
176 TSRHSSplitGetSubTSs - Get an array of all sub-TS contexts.
177
178 Logically Collective on TS
179
180 Output Parameters:
181 + n - the number of splits
182 - subksp - the array of TS contexts
183
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).
187
188 Level: advanced
189
190 .seealso: TSGetRHSSplitFunction()
191 @*/
TSRHSSplitGetSubTSs(TS ts,PetscInt * n,TS * subts[])192 PetscErrorCode TSRHSSplitGetSubTSs(TS ts,PetscInt *n,TS *subts[])
193 {
194 TS_RHSSplitLink ilink = ts->tsrhssplit;
195 PetscInt i = 0;
196 PetscErrorCode ierr;
197
198 PetscFunctionBegin;
199 PetscValidHeaderSpecific(ts,TS_CLASSID,1);
200 if (subts) {
201 ierr = PetscMalloc1(ts->num_rhs_splits,subts);CHKERRQ(ierr);
202 while (ilink) {
203 (*subts)[i++] = ilink->ts;
204 ilink = ilink->next;
205 }
206 }
207 if (n) *n = ts->num_rhs_splits;
208 PetscFunctionReturn(0);
209 }
210