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