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: PetscWeakFormGetResidual(wf, key.label, key.value, key.field, key.part, &n0, &f0, &n1, &f1);
35: return 0;
36: }
38: static PetscErrorCode TestSetIndex(PetscWeakForm wf)
39: {
40: PetscPointFunc f[4] = {f0, f1, f2, f3};
41: DMLabel label;
42: const PetscInt value = 3, field = 1, part = 2;
43: PetscFormKey key;
44: PetscInt i, j, k, l;
46: DMLabelCreate(PETSC_COMM_SELF, "Test", &label);
47: key.label = label;
48: key.value = value;
49: key.field = field;
50: key.part = part;
51: /* Check f0 */
52: for (i = 0; i < 4; ++i) {
53: for (j = 0; j < 4; ++j) {
54: if (j == i) continue;
55: for (k = 0; k < 4; ++k) {
56: if ((k == i) || (k == j)) continue;
57: for (l = 0; l < 4; ++l) {
58: if ((l == i) || (l == j) || (l == k)) continue;
59: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, i, f[i], 0, NULL);
60: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, j, f[j], 0, NULL);
61: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, k, f[k], 0, NULL);
62: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, l, f[l], 0, NULL);
63: CheckResidual(wf, key, 4, f, 0, NULL);
64: PetscWeakFormClear(wf);
65: }
66: }
67: }
68: }
69: /* Check f1 */
70: for (i = 0; i < 4; ++i) {
71: for (j = 0; j < 4; ++j) {
72: if (j == i) continue;
73: for (k = 0; k < 4; ++k) {
74: if ((k == i) || (k == j)) continue;
75: for (l = 0; l < 4; ++l) {
76: if ((l == i) || (l == j) || (l == k)) continue;
77: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, i, f[i]);
78: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, j, f[j]);
79: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, k, f[k]);
80: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, l, f[l]);
81: CheckResidual(wf, key, 0, NULL, 4, f);
82: PetscWeakFormClear(wf);
83: }
84: }
85: }
86: }
87: /* Check f0 and f1 */
88: for (i = 0; i < 4; ++i) {
89: for (j = 0; j < 4; ++j) {
90: if (j == i) continue;
91: for (k = 0; k < 4; ++k) {
92: if ((k == i) || (k == j)) continue;
93: for (l = 0; l < 4; ++l) {
94: if ((l == i) || (l == j) || (l == k)) continue;
95: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, i, f[i], i, f[i]);
96: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, j, f[j], j, f[j]);
97: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, k, f[k], k, f[k]);
98: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, l, f[l], l, f[l]);
99: CheckResidual(wf, key, 4, f, 4, f);
100: PetscWeakFormClear(wf);
101: }
102: }
103: }
104: }
105: /* Check f0 and f1 in different orders */
106: for (i = 0; i < 4; ++i) {
107: for (j = 0; j < 4; ++j) {
108: if (j == i) continue;
109: for (k = 0; k < 4; ++k) {
110: if ((k == i) || (k == j)) continue;
111: for (l = 0; l < 4; ++l) {
112: if ((l == i) || (l == j) || (l == k)) continue;
113: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, l, f[l], i, f[i]);
114: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, k, f[k], j, f[j]);
115: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, j, f[j], k, f[k]);
116: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, i, f[i], l, f[l]);
117: CheckResidual(wf, key, 4, f, 4, f);
118: PetscWeakFormClear(wf);
119: }
120: }
121: }
122: }
123: DMLabelDestroy(&label);
124: return 0;
125: }
127: static PetscErrorCode TestAdd(PetscWeakForm wf)
128: {
129: PetscPointFunc f[4] = {f0, f1, f2, f3}, fp[4];
130: DMLabel label;
131: const PetscInt value = 3, field = 1, part = 2;
132: PetscFormKey key;
133: PetscInt i, j, k, l;
135: DMLabelCreate(PETSC_COMM_SELF, "Test", &label);
136: key.label = label;
137: key.value = value;
138: key.field = field;
139: key.part = part;
140: /* Check f0 */
141: for (i = 0; i < 4; ++i) {
142: for (j = 0; j < 4; ++j) {
143: if (j == i) continue;
144: for (k = 0; k < 4; ++k) {
145: if ((k == i) || (k == j)) continue;
146: for (l = 0; l < 4; ++l) {
147: if ((l == i) || (l == j) || (l == k)) continue;
148: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[i], NULL);
149: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[j], NULL);
150: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[k], NULL);
151: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[l], NULL);
152: fp[0] = f[i];
153: fp[1] = f[j];
154: fp[2] = f[k];
155: fp[3] = f[l];
156: CheckResidual(wf, key, 4, fp, 0, NULL);
157: PetscWeakFormClear(wf);
158: }
159: }
160: }
161: }
162: /* Check f1 */
163: for (i = 0; i < 4; ++i) {
164: for (j = 0; j < 4; ++j) {
165: if (j == i) continue;
166: for (k = 0; k < 4; ++k) {
167: if ((k == i) || (k == j)) continue;
168: for (l = 0; l < 4; ++l) {
169: if ((l == i) || (l == j) || (l == k)) continue;
170: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[i]);
171: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[j]);
172: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[k]);
173: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[l]);
174: fp[0] = f[i];
175: fp[1] = f[j];
176: fp[2] = f[k];
177: fp[3] = f[l];
178: CheckResidual(wf, key, 0, NULL, 4, fp);
179: PetscWeakFormClear(wf);
180: }
181: }
182: }
183: }
184: /* Check f0 and f1 */
185: for (i = 0; i < 4; ++i) {
186: for (j = 0; j < 4; ++j) {
187: if (j == i) continue;
188: for (k = 0; k < 4; ++k) {
189: if ((k == i) || (k == j)) continue;
190: for (l = 0; l < 4; ++l) {
191: if ((l == i) || (l == j) || (l == k)) continue;
192: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[i], f[i]);
193: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[j], f[j]);
194: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[k], f[k]);
195: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[l], f[l]);
196: fp[0] = f[i];
197: fp[1] = f[j];
198: fp[2] = f[k];
199: fp[3] = f[l];
200: CheckResidual(wf, key, 4, fp, 4, fp);
201: PetscWeakFormClear(wf);
202: }
203: }
204: }
205: }
206: DMLabelDestroy(&label);
207: return 0;
208: }
210: static PetscErrorCode TestSetIndexAdd(PetscWeakForm wf)
211: {
212: PetscPointFunc f[4] = {f0, f1, f2, f3};
213: DMLabel label;
214: const PetscInt value = 3, field = 1, part = 2;
215: PetscFormKey key;
217: DMLabelCreate(PETSC_COMM_SELF, "Test", &label);
218: key.label = label;
219: key.value = value;
220: key.field = field;
221: key.part = part;
222: /* Check f0 */
223: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, f[0], 0, NULL);
224: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 1, f[1], 0, NULL);
225: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[2], NULL);
226: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[3], NULL);
227: CheckResidual(wf, key, 4, f, 0, NULL);
228: PetscWeakFormClear(wf);
229: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, f[0], 0, NULL);
230: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[1], NULL);
231: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 2, f[2], 0, NULL);
232: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[3], NULL);
233: CheckResidual(wf, key, 4, f, 0, NULL);
234: PetscWeakFormClear(wf);
235: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, f[0], 0, NULL);
236: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[1], NULL);
237: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[2], NULL);
238: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 3, f[3], 0, NULL);
239: CheckResidual(wf, key, 4, f, 0, NULL);
240: PetscWeakFormClear(wf);
241: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[0], NULL);
242: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 1, f[1], 0, NULL);
243: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[2], NULL);
244: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 3, f[3], 0, NULL);
245: CheckResidual(wf, key, 4, f, 0, NULL);
246: PetscWeakFormClear(wf);
247: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[0], NULL);
248: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, f[1], NULL);
249: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 2, f[2], 0, NULL);
250: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 3, f[3], 0, NULL);
251: CheckResidual(wf, key, 4, f, 0, NULL);
252: PetscWeakFormClear(wf);
253: /* Check f1 */
254: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 0, f[0]);
255: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 1, f[1]);
256: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[2]);
257: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[3]);
258: CheckResidual(wf, key, 0, NULL, 4, f);
259: PetscWeakFormClear(wf);
260: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 0, f[0]);
261: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[1]);
262: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 2, f[2]);
263: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[3]);
264: CheckResidual(wf, key, 0, NULL, 4, f);
265: PetscWeakFormClear(wf);
266: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 0, f[0]);
267: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[1]);
268: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[2]);
269: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 3, f[3]);
270: CheckResidual(wf, key, 0, NULL, 4, f);
271: PetscWeakFormClear(wf);
272: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[0]);
273: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 1, f[1]);
274: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[2]);
275: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 3, f[3]);
276: CheckResidual(wf, key, 0, NULL, 4, f);
277: PetscWeakFormClear(wf);
278: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[0]);
279: PetscWeakFormAddResidual(wf, key.label, key.value, key.field, key.part, NULL, f[1]);
280: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 2, f[2]);
281: PetscWeakFormSetIndexResidual(wf, key.label, key.value, key.field, key.part, 0, NULL, 3, f[3]);
282: CheckResidual(wf, key, 0, NULL, 4, f);
283: PetscWeakFormClear(wf);
285: DMLabelDestroy(&label);
286: return 0;
287: }
289: int main(int argc, char **argv)
290: {
291: PetscWeakForm wf;
294: PetscInitialize(&argc, &argv, NULL, help);
295: PetscWeakFormCreate(PETSC_COMM_SELF, &wf);
296: TestSetIndex(wf);
297: TestAdd(wf);
298: TestSetIndexAdd(wf);
299: PetscWeakFormDestroy(&wf);
300: PetscFinalize();
301: return 0;
302: }
304: /*TEST
306: test:
307: suffix: 0
309: TEST*/