xref: /petsc/src/dm/dt/tests/ex12.c (revision 609caa7c8c030312b00807b4f015fd827bb80932)
1 static char help[] = "Tests for PetscWeakForm.\n\n";
2 
3 #include <petscds.h>
4 
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[])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 }
9 
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[])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 }
14 
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[])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 }
19 
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[])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 }
24 
CheckResidual(PetscWeakForm wf,PetscFormKey key,PetscInt in0,PetscPointFn * if0[],PetscInt in1,PetscPointFn * if1[])25 static PetscErrorCode CheckResidual(PetscWeakForm wf, PetscFormKey key, PetscInt in0, PetscPointFn *if0[], PetscInt in1, PetscPointFn *if1[])
26 {
27   PetscPointFn **f0, **f1;
28   PetscInt       n0, n1, i;
29 
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 }
38 
TestSetIndex(PetscWeakForm wf)39 static PetscErrorCode TestSetIndex(PetscWeakForm wf)
40 {
41   PetscPointFn  *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;
46 
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 }
128 
TestAdd(PetscWeakForm wf)129 static PetscErrorCode TestAdd(PetscWeakForm wf)
130 {
131   PetscPointFn  *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;
136 
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 }
212 
TestSetIndexAdd(PetscWeakForm wf)213 static PetscErrorCode TestSetIndexAdd(PetscWeakForm wf)
214 {
215   PetscPointFn  *f[4] = {f0, f1, f2, f3};
216   DMLabel        label;
217   const PetscInt value = 3, field = 1, part = 2;
218   PetscFormKey   key;
219 
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));
288 
289   PetscCall(DMLabelDestroy(&label));
290   PetscFunctionReturn(PETSC_SUCCESS);
291 }
292 
main(int argc,char ** argv)293 int main(int argc, char **argv)
294 {
295   PetscWeakForm wf;
296 
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 }
307 
308 /*TEST
309 
310   test:
311     suffix: 0
312     output_file: output/empty.out
313 
314 TEST*/
315