Actual source code: ex12.c
1: static char help[] = "Tests for PetscWeakForm.\n\n";
3: #include <petscds.h>
5: static void f0(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
6: {
7: f0[0] = 0.0;
8: }
10: static void f1(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
11: {
12: f0[0] = 0.0;
13: }
15: static void f2(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
16: {
17: f0[0] = 0.0;
18: }
20: static void f3(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[])
21: {
22: f0[0] = 0.0;
23: }
25: static PetscErrorCode CheckResidual(PetscWeakForm wf, PetscFormKey key, PetscInt in0, PetscPointFunc if0[], PetscInt in1, PetscPointFunc if1[])
26: {
27: PetscPointFunc *f0, *f1;
28: PetscInt n0, n1, i;
30: PetscFunctionBegin;
31: PetscCall(PetscWeakFormGetResidual(wf, key.label, key.value, key.field, key.part, &n0, &f0, &n1, &f1));
32: PetscCheck(n0 == in0, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Found %" PetscInt_FMT " f0 functions != %" PetscInt_FMT " functions input", n0, in0);
33: for (i = 0; i < n0; ++i) PetscCheck(f0[i] == if0[i], PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "f0[%" PetscInt_FMT "] != input f0[%" PetscInt_FMT "]", i, i);
34: PetscCheck(n1 == in1, PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "Found %" PetscInt_FMT " f1 functions != %" PetscInt_FMT " functions input", n0, in0);
35: for (i = 0; i < n1; ++i) PetscCheck(f1[i] == if1[i], PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "f1[%" PetscInt_FMT "] != input f1[%" PetscInt_FMT "]", i, i);
36: PetscFunctionReturn(PETSC_SUCCESS);
37: }
39: static PetscErrorCode TestSetIndex(PetscWeakForm wf)
40: {
41: PetscPointFunc f[4] = {f0, f1, f2, f3};
42: DMLabel label;
43: const PetscInt value = 3, field = 1, part = 2;
44: PetscFormKey key;
45: PetscInt i, j, k, l;
47: PetscFunctionBegin;
48: PetscCall(DMLabelCreate(PETSC_COMM_SELF, "Test", &label));
49: key.label = label;
50: key.value = value;
51: key.field = field;
52: key.part = part;
53: /* Check f0 */
54: for (i = 0; i < 4; ++i) {
55: for (j = 0; j < 4; ++j) {
56: if (j == i) continue;
57: for (k = 0; k < 4; ++k) {
58: if ((k == i) || (k == j)) continue;
59: for (l = 0; l < 4; ++l) {
60: if ((l == i) || (l == j) || (l == k)) continue;
61: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, i, f[i], 0, NULL));
62: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, j, f[j], 0, NULL));
63: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, k, f[k], 0, NULL));
64: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, l, f[l], 0, NULL));
65: PetscCall(CheckResidual(wf, key, 4, f, 0, NULL));
66: PetscCall(PetscWeakFormClear(wf));
67: }
68: }
69: }
70: }
71: /* Check f1 */
72: for (i = 0; i < 4; ++i) {
73: for (j = 0; j < 4; ++j) {
74: if (j == i) continue;
75: for (k = 0; k < 4; ++k) {
76: if ((k == i) || (k == j)) continue;
77: for (l = 0; l < 4; ++l) {
78: if ((l == i) || (l == j) || (l == k)) continue;
79: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, i, f[i]));
80: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, j, f[j]));
81: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, k, f[k]));
82: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, l, f[l]));
83: PetscCall(CheckResidual(wf, key, 0, NULL, 4, f));
84: PetscCall(PetscWeakFormClear(wf));
85: }
86: }
87: }
88: }
89: /* Check f0 and f1 */
90: for (i = 0; i < 4; ++i) {
91: for (j = 0; j < 4; ++j) {
92: if (j == i) continue;
93: for (k = 0; k < 4; ++k) {
94: if ((k == i) || (k == j)) continue;
95: for (l = 0; l < 4; ++l) {
96: if ((l == i) || (l == j) || (l == k)) continue;
97: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, i, f[i], i, f[i]));
98: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, j, f[j], j, f[j]));
99: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, k, f[k], k, f[k]));
100: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, l, f[l], l, f[l]));
101: PetscCall(CheckResidual(wf, key, 4, f, 4, f));
102: PetscCall(PetscWeakFormClear(wf));
103: }
104: }
105: }
106: }
107: /* Check f0 and f1 in different orders */
108: for (i = 0; i < 4; ++i) {
109: for (j = 0; j < 4; ++j) {
110: if (j == i) continue;
111: for (k = 0; k < 4; ++k) {
112: if ((k == i) || (k == j)) continue;
113: for (l = 0; l < 4; ++l) {
114: if ((l == i) || (l == j) || (l == k)) continue;
115: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, l, f[l], i, f[i]));
116: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, k, f[k], j, f[j]));
117: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, j, f[j], k, f[k]));
118: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, i, f[i], l, f[l]));
119: PetscCall(CheckResidual(wf, key, 4, f, 4, f));
120: PetscCall(PetscWeakFormClear(wf));
121: }
122: }
123: }
124: }
125: PetscCall(DMLabelDestroy(&label));
126: PetscFunctionReturn(PETSC_SUCCESS);
127: }
129: static PetscErrorCode TestAdd(PetscWeakForm wf)
130: {
131: PetscPointFunc f[4] = {f0, f1, f2, f3}, fp[4];
132: DMLabel label;
133: const PetscInt value = 3, field = 1, part = 2;
134: PetscFormKey key;
135: PetscInt i, j, k, l;
137: PetscFunctionBegin;
138: PetscCall(DMLabelCreate(PETSC_COMM_SELF, "Test", &label));
139: key.label = label;
140: key.value = value;
141: key.field = field;
142: key.part = part;
143: /* Check f0 */
144: for (i = 0; i < 4; ++i) {
145: for (j = 0; j < 4; ++j) {
146: if (j == i) continue;
147: for (k = 0; k < 4; ++k) {
148: if ((k == i) || (k == j)) continue;
149: for (l = 0; l < 4; ++l) {
150: if ((l == i) || (l == j) || (l == k)) continue;
151: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[i], NULL));
152: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[j], NULL));
153: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[k], NULL));
154: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[l], NULL));
155: fp[0] = f[i];
156: fp[1] = f[j];
157: fp[2] = f[k];
158: fp[3] = f[l];
159: PetscCall(CheckResidual(wf, key, 4, fp, 0, NULL));
160: PetscCall(PetscWeakFormClear(wf));
161: }
162: }
163: }
164: }
165: /* Check f1 */
166: for (i = 0; i < 4; ++i) {
167: for (j = 0; j < 4; ++j) {
168: if (j == i) continue;
169: for (k = 0; k < 4; ++k) {
170: if ((k == i) || (k == j)) continue;
171: for (l = 0; l < 4; ++l) {
172: if ((l == i) || (l == j) || (l == k)) continue;
173: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[i]));
174: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[j]));
175: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[k]));
176: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[l]));
177: fp[0] = f[i];
178: fp[1] = f[j];
179: fp[2] = f[k];
180: fp[3] = f[l];
181: PetscCall(CheckResidual(wf, key, 0, NULL, 4, fp));
182: PetscCall(PetscWeakFormClear(wf));
183: }
184: }
185: }
186: }
187: /* Check f0 and f1 */
188: for (i = 0; i < 4; ++i) {
189: for (j = 0; j < 4; ++j) {
190: if (j == i) continue;
191: for (k = 0; k < 4; ++k) {
192: if ((k == i) || (k == j)) continue;
193: for (l = 0; l < 4; ++l) {
194: if ((l == i) || (l == j) || (l == k)) continue;
195: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[i], f[i]));
196: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[j], f[j]));
197: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[k], f[k]));
198: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[l], f[l]));
199: fp[0] = f[i];
200: fp[1] = f[j];
201: fp[2] = f[k];
202: fp[3] = f[l];
203: PetscCall(CheckResidual(wf, key, 4, fp, 4, fp));
204: PetscCall(PetscWeakFormClear(wf));
205: }
206: }
207: }
208: }
209: PetscCall(DMLabelDestroy(&label));
210: PetscFunctionReturn(PETSC_SUCCESS);
211: }
213: static PetscErrorCode TestSetIndexAdd(PetscWeakForm wf)
214: {
215: PetscPointFunc f[4] = {f0, f1, f2, f3};
216: DMLabel label;
217: const PetscInt value = 3, field = 1, part = 2;
218: PetscFormKey key;
220: PetscFunctionBegin;
221: PetscCall(DMLabelCreate(PETSC_COMM_SELF, "Test", &label));
222: key.label = label;
223: key.value = value;
224: key.field = field;
225: key.part = part;
226: /* Check f0 */
227: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, f[0], 0, NULL));
228: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 1, f[1], 0, NULL));
229: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[2], NULL));
230: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[3], NULL));
231: PetscCall(CheckResidual(wf, key, 4, f, 0, NULL));
232: PetscCall(PetscWeakFormClear(wf));
233: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, f[0], 0, NULL));
234: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[1], NULL));
235: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 2, f[2], 0, NULL));
236: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[3], NULL));
237: PetscCall(CheckResidual(wf, key, 4, f, 0, NULL));
238: PetscCall(PetscWeakFormClear(wf));
239: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, f[0], 0, NULL));
240: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[1], NULL));
241: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[2], NULL));
242: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 3, f[3], 0, NULL));
243: PetscCall(CheckResidual(wf, key, 4, f, 0, NULL));
244: PetscCall(PetscWeakFormClear(wf));
245: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[0], NULL));
246: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 1, f[1], 0, NULL));
247: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[2], NULL));
248: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 3, f[3], 0, NULL));
249: PetscCall(CheckResidual(wf, key, 4, f, 0, NULL));
250: PetscCall(PetscWeakFormClear(wf));
251: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[0], NULL));
252: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[1], NULL));
253: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 2, f[2], 0, NULL));
254: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 3, f[3], 0, NULL));
255: PetscCall(CheckResidual(wf, key, 4, f, 0, NULL));
256: PetscCall(PetscWeakFormClear(wf));
257: /* Check f1 */
258: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 0, f[0]));
259: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 1, f[1]));
260: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[2]));
261: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[3]));
262: PetscCall(CheckResidual(wf, key, 0, NULL, 4, f));
263: PetscCall(PetscWeakFormClear(wf));
264: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 0, f[0]));
265: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[1]));
266: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 2, f[2]));
267: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[3]));
268: PetscCall(CheckResidual(wf, key, 0, NULL, 4, f));
269: PetscCall(PetscWeakFormClear(wf));
270: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 0, f[0]));
271: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[1]));
272: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[2]));
273: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 3, f[3]));
274: PetscCall(CheckResidual(wf, key, 0, NULL, 4, f));
275: PetscCall(PetscWeakFormClear(wf));
276: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[0]));
277: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 1, f[1]));
278: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[2]));
279: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 3, f[3]));
280: PetscCall(CheckResidual(wf, key, 0, NULL, 4, f));
281: PetscCall(PetscWeakFormClear(wf));
282: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[0]));
283: PetscCall(PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[1]));
284: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 2, f[2]));
285: PetscCall(PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 3, f[3]));
286: PetscCall(CheckResidual(wf, key, 0, NULL, 4, f));
287: PetscCall(PetscWeakFormClear(wf));
289: PetscCall(DMLabelDestroy(&label));
290: PetscFunctionReturn(PETSC_SUCCESS);
291: }
293: int main(int argc, char **argv)
294: {
295: PetscWeakForm wf;
297: PetscFunctionBeginUser;
298: PetscCall(PetscInitialize(&argc, &argv, NULL, help));
299: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &wf));
300: PetscCall(TestSetIndex(wf));
301: PetscCall(TestAdd(wf));
302: PetscCall(TestSetIndexAdd(wf));
303: PetscCall(PetscWeakFormDestroy(&wf));
304: PetscCall(PetscFinalize());
305: return 0;
306: }
308: /*TEST
310: test:
311: suffix: 0
313: TEST*/