Actual source code: tsrhssplit.c
1: #include <petsc/private/tsimpl.h>
2: #include <petscdm.h>
3: static PetscErrorCode TSRHSSplitGetRHSSplit(TS ts, const char splitname[], TS_RHSSplitLink *isplit)
4: {
5: PetscBool found = PETSC_FALSE;
7: *isplit = ts->tsrhssplit;
8: /* look up the split */
9: while (*isplit) {
10: PetscStrcmp((*isplit)->splitname, splitname, &found);
11: if (found) break;
12: *isplit = (*isplit)->next;
13: }
14: return 0;
15: }
17: /*@C
18: TSRHSSplitSetIS - Set the index set for the specified split
20: Logically Collective
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
27: Level: intermediate
29: .seealso: [](chapter_ts), `TS`, `IS`, `TSRHSSplitGetIS()`
30: @*/
31: PetscErrorCode TSRHSSplitSetIS(TS ts, const char splitname[], IS is)
32: {
33: TS_RHSSplitLink newsplit, next = ts->tsrhssplit;
34: char prefix[128];
39: PetscNew(&newsplit);
40: if (splitname) {
41: PetscStrallocpy(splitname, &newsplit->splitname);
42: } else {
43: PetscMalloc1(8, &newsplit->splitname);
44: PetscSNPrintf(newsplit->splitname, 7, "%" PetscInt_FMT, ts->num_rhs_splits);
45: }
46: PetscObjectReference((PetscObject)is);
47: newsplit->is = is;
48: TSCreate(PetscObjectComm((PetscObject)ts), &newsplit->ts);
50: PetscObjectIncrementTabLevel((PetscObject)newsplit->ts, (PetscObject)ts, 1);
51: PetscSNPrintf(prefix, sizeof(prefix), "%srhsplit_%s_", ((PetscObject)ts)->prefix ? ((PetscObject)ts)->prefix : "", newsplit->splitname);
52: TSSetOptionsPrefix(newsplit->ts, prefix);
53: if (!next) ts->tsrhssplit = newsplit;
54: else {
55: while (next->next) next = next->next;
56: next->next = newsplit;
57: }
58: ts->num_rhs_splits++;
59: return 0;
60: }
62: /*@C
63: TSRHSSplitGetIS - Retrieves the elements for a split as an `IS`
65: Logically Collective
67: Input Parameters:
68: + ts - the `TS` context obtained from `TSCreate()`
69: - splitname - name of this split
71: Output Parameters:
72: - is - the index set for part of the solution vector
74: Level: intermediate
76: .seealso: [](chapter_ts), `TS`, `IS`, `TSRHSSplitSetIS()`
77: @*/
78: PetscErrorCode TSRHSSplitGetIS(TS ts, const char splitname[], IS *is)
79: {
80: TS_RHSSplitLink isplit;
83: *is = NULL;
84: /* look up the split */
85: TSRHSSplitGetRHSSplit(ts, splitname, &isplit);
86: if (isplit) *is = isplit->is;
87: return 0;
88: }
90: /*@C
91: TSRHSSplitSetRHSFunction - Set the split right-hand-side functions.
93: Logically Collective
95: Input Parameters:
96: + ts - the `TS` context obtained from `TSCreate()`
97: . splitname - name of this split
98: . r - vector to hold the residual (or NULL to have it created internally)
99: . rhsfunc - the RHS function evaluation routine
100: - ctx - user-defined context for private data for the split function evaluation routine (may be NULL)
102: Calling sequence of fun:
103: $ rhsfunc(TS ts,PetscReal t,Vec u,Vec f,ctx);
105: + t - time at step/stage being solved
106: . u - state vector
107: . f - function vector
108: - ctx - [optional] user-defined context for matrix evaluation routine (may be NULL)
110: Level: beginner
112: .seealso: [](chapter_ts), `TS`, `IS`, `TSRHSSplitSetIS()`
113: @*/
114: PetscErrorCode TSRHSSplitSetRHSFunction(TS ts, const char splitname[], Vec r, TSRHSFunction rhsfunc, void *ctx)
115: {
116: TS_RHSSplitLink isplit;
117: DM dmc;
118: Vec subvec, ralloc = NULL;
123: /* look up the split */
124: TSRHSSplitGetRHSSplit(ts, splitname, &isplit);
127: if (!r && ts->vec_sol) {
128: VecGetSubVector(ts->vec_sol, isplit->is, &subvec);
129: VecDuplicate(subvec, &ralloc);
130: r = ralloc;
131: VecRestoreSubVector(ts->vec_sol, isplit->is, &subvec);
132: }
134: if (ts->dm) {
135: PetscInt dim;
137: DMGetDimension(ts->dm, &dim);
138: if (dim != -1) {
139: DMClone(ts->dm, &dmc);
140: TSSetDM(isplit->ts, dmc);
141: DMDestroy(&dmc);
142: }
143: }
145: TSSetRHSFunction(isplit->ts, r, rhsfunc, ctx);
146: VecDestroy(&ralloc);
147: return 0;
148: }
150: /*@C
151: TSRHSSplitGetSubTS - Get the sub-`TS` by split name.
153: Logically Collective
155: Input Parameter:
156: . ts - the `TS` context obtained from `TSCreate()`
158: Output Parameters:
159: + splitname - the number of the split
160: - subts - the array of `TS` contexts
162: Level: advanced
164: .seealso: [](chapter_ts), `TS`, `IS`, `TSGetRHSSplitFunction()`
165: @*/
166: PetscErrorCode TSRHSSplitGetSubTS(TS ts, const char splitname[], TS *subts)
167: {
168: TS_RHSSplitLink isplit;
172: *subts = NULL;
173: /* look up the split */
174: TSRHSSplitGetRHSSplit(ts, splitname, &isplit);
175: if (isplit) *subts = isplit->ts;
176: return 0;
177: }
179: /*@C
180: TSRHSSplitGetSubTSs - Get an array of all sub-`TS` contexts.
182: Logically Collective
184: Input Parameter:
185: . ts - the `TS` context obtained from `TSCreate()`
187: Output Parameters:
188: + n - the number of splits
189: - subksp - the array of `TS` contexts
191: Level: advanced
193: Note:
194: After `TSRHSSplitGetSubTS()` the array of `TS`s is to be freed by the user with `PetscFree()`
195: (not the `TS` just the array that contains them).
197: .seealso: [](chapter_ts), `TS`, `IS`, `TSGetRHSSplitFunction()`
198: @*/
199: PetscErrorCode TSRHSSplitGetSubTSs(TS ts, PetscInt *n, TS *subts[])
200: {
201: TS_RHSSplitLink ilink = ts->tsrhssplit;
202: PetscInt i = 0;
205: if (subts) {
206: PetscMalloc1(ts->num_rhs_splits, subts);
207: while (ilink) {
208: (*subts)[i++] = ilink->ts;
209: ilink = ilink->next;
210: }
211: }
212: if (n) *n = ts->num_rhs_splits;
213: return 0;
214: }