1 #include <petsc/private/petscdsimpl.h> /*I "petscds.h" I*/ 2 3 PetscClassId PETSCWEAKFORM_CLASSID = 0; 4 5 const char *const PetscWeakFormKinds[] = {"objective", "residual_f0", "residual_f1", "jacobian_g0", "jacobian_g1", "jacobian_g2", "jacobian_g3", "jacobian_preconditioner_g0", "jacobian_preconditioner_g1", "jacobian_preconditioner_g2", "jacobian_preconditioner_g3", "dynamic_jacobian_g0", "dynamic_jacobian_g1", "dynamic_jacobian_g2", "dynamic_jacobian_g3", "boundary_residual_f0", "boundary_residual_f1", "boundary_jacobian_g0", "boundary_jacobian_g1", "boundary_jacobian_g2", "boundary_jacobian_g3", "boundary_jacobian_preconditioner_g0", "boundary_jacobian_preconditioner_g1", "boundary_jacobian_preconditioner_g2", "boundary_jacobian_preconditioner_g3", "riemann_solver", "PetscWeakFormKind", "PETSC_WF_", NULL}; 6 7 static PetscErrorCode PetscChunkBufferCreate(size_t unitbytes, size_t expected, PetscChunkBuffer **buffer) 8 { 9 PetscFunctionBegin; 10 PetscCall(PetscNew(buffer)); 11 PetscCall(PetscCalloc1(expected * unitbytes, &(*buffer)->array)); 12 (*buffer)->size = expected; 13 (*buffer)->unitbytes = unitbytes; 14 (*buffer)->alloc = expected * unitbytes; 15 PetscFunctionReturn(0); 16 } 17 18 static PetscErrorCode PetscChunkBufferDuplicate(PetscChunkBuffer *buffer, PetscChunkBuffer **bufferNew) 19 { 20 PetscFunctionBegin; 21 PetscCall(PetscNew(bufferNew)); 22 PetscCall(PetscCalloc1(buffer->size * buffer->unitbytes, &(*bufferNew)->array)); 23 PetscCall(PetscMemcpy((*bufferNew)->array, buffer->array, buffer->size * buffer->unitbytes)); 24 (*bufferNew)->size = buffer->size; 25 (*bufferNew)->unitbytes = buffer->unitbytes; 26 (*bufferNew)->alloc = buffer->size * buffer->unitbytes; 27 PetscFunctionReturn(0); 28 } 29 30 static PetscErrorCode PetscChunkBufferDestroy(PetscChunkBuffer **buffer) 31 { 32 PetscFunctionBegin; 33 PetscCall(PetscFree((*buffer)->array)); 34 PetscCall(PetscFree(*buffer)); 35 PetscFunctionReturn(0); 36 } 37 38 static PetscErrorCode PetscChunkBufferCreateChunk(PetscChunkBuffer *buffer, PetscInt size, PetscChunk *chunk) 39 { 40 PetscFunctionBegin; 41 if ((buffer->size + size) * buffer->unitbytes > buffer->alloc) { 42 char *tmp; 43 44 if (!buffer->alloc) buffer->alloc = (buffer->size + size) * buffer->unitbytes; 45 while ((buffer->size + size) * buffer->unitbytes > buffer->alloc) buffer->alloc *= 2; 46 PetscCall(PetscMalloc(buffer->alloc, &tmp)); 47 PetscCall(PetscMemcpy(tmp, buffer->array, buffer->size * buffer->unitbytes)); 48 PetscCall(PetscFree(buffer->array)); 49 buffer->array = tmp; 50 } 51 chunk->start = buffer->size * buffer->unitbytes; 52 chunk->size = size; 53 chunk->reserved = size; 54 buffer->size += size; 55 PetscFunctionReturn(0); 56 } 57 58 static PetscErrorCode PetscChunkBufferEnlargeChunk(PetscChunkBuffer *buffer, PetscInt size, PetscChunk *chunk) 59 { 60 size_t siz = size; 61 62 PetscFunctionBegin; 63 if (chunk->size + size > chunk->reserved) { 64 PetscChunk newchunk; 65 PetscInt reserved = chunk->size; 66 67 /* TODO Here if we had a chunk list, we could update them all to reclaim unused space */ 68 while (reserved < chunk->size + size) reserved *= 2; 69 PetscCall(PetscChunkBufferCreateChunk(buffer, (size_t)reserved, &newchunk)); 70 newchunk.size = chunk->size + size; 71 PetscCall(PetscMemcpy(&buffer->array[newchunk.start], &buffer->array[chunk->start], chunk->size * buffer->unitbytes)); 72 *chunk = newchunk; 73 } else { 74 chunk->size += siz; 75 } 76 PetscFunctionReturn(0); 77 } 78 79 /*@C 80 PetscFormKeySort - Sorts an array of PetscFormKey in place in increasing order. 81 82 Not Collective 83 84 Input Parameters: 85 + n - number of values 86 - X - array of PetscFormKey 87 88 Level: intermediate 89 90 .seealso: `PetscIntSortSemiOrdered()`, `PetscSortInt()` 91 @*/ 92 PetscErrorCode PetscFormKeySort(PetscInt n, PetscFormKey arr[]) 93 { 94 PetscFunctionBegin; 95 if (n <= 1) PetscFunctionReturn(0); 96 PetscValidPointer(arr, 2); 97 PetscCall(PetscTimSort(n, arr, sizeof(PetscFormKey), Compare_PetscFormKey_Private, NULL)); 98 PetscFunctionReturn(0); 99 } 100 101 PetscErrorCode PetscWeakFormGetFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt *n, void (***func)()) 102 { 103 PetscFormKey key; 104 PetscChunk chunk; 105 106 PetscFunctionBegin; 107 key.label = label; 108 key.value = value; 109 key.field = f; 110 key.part = part; 111 PetscCall(PetscHMapFormGet(ht, key, &chunk)); 112 if (chunk.size < 0) { 113 *n = 0; 114 *func = NULL; 115 } else { 116 *n = chunk.size; 117 *func = (void (**)()) & wf->funcs->array[chunk.start]; 118 } 119 PetscFunctionReturn(0); 120 } 121 122 /* A NULL argument for func causes this to clear the key */ 123 PetscErrorCode PetscWeakFormSetFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt n, void (**func)()) 124 { 125 PetscFormKey key; 126 PetscChunk chunk; 127 PetscInt i; 128 129 PetscFunctionBegin; 130 key.label = label; 131 key.value = value; 132 key.field = f; 133 key.part = part; 134 if (!func) { 135 PetscCall(PetscHMapFormDel(ht, key)); 136 PetscFunctionReturn(0); 137 } else PetscCall(PetscHMapFormGet(ht, key, &chunk)); 138 if (chunk.size < 0) { 139 PetscCall(PetscChunkBufferCreateChunk(wf->funcs, n, &chunk)); 140 PetscCall(PetscHMapFormSet(ht, key, chunk)); 141 } else if (chunk.size <= n) { 142 PetscCall(PetscChunkBufferEnlargeChunk(wf->funcs, n - chunk.size, &chunk)); 143 PetscCall(PetscHMapFormSet(ht, key, chunk)); 144 } 145 for (i = 0; i < n; ++i) ((void (**)()) & wf->funcs->array[chunk.start])[i] = func[i]; 146 PetscFunctionReturn(0); 147 } 148 149 PetscErrorCode PetscWeakFormAddFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, void (*func)()) 150 { 151 PetscFormKey key; 152 PetscChunk chunk; 153 154 PetscFunctionBegin; 155 if (!func) PetscFunctionReturn(0); 156 key.label = label; 157 key.value = value; 158 key.field = f; 159 key.part = part; 160 PetscCall(PetscHMapFormGet(ht, key, &chunk)); 161 if (chunk.size < 0) { 162 PetscCall(PetscChunkBufferCreateChunk(wf->funcs, 1, &chunk)); 163 PetscCall(PetscHMapFormSet(ht, key, chunk)); 164 ((void (**)()) & wf->funcs->array[chunk.start])[0] = func; 165 } else { 166 PetscCall(PetscChunkBufferEnlargeChunk(wf->funcs, 1, &chunk)); 167 PetscCall(PetscHMapFormSet(ht, key, chunk)); 168 ((void (**)()) & wf->funcs->array[chunk.start])[chunk.size - 1] = func; 169 } 170 PetscFunctionReturn(0); 171 } 172 173 PetscErrorCode PetscWeakFormGetIndexFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt ind, void (**func)()) 174 { 175 PetscFormKey key; 176 PetscChunk chunk; 177 178 PetscFunctionBegin; 179 key.label = label; 180 key.value = value; 181 key.field = f; 182 key.part = part; 183 PetscCall(PetscHMapFormGet(ht, key, &chunk)); 184 if (chunk.size < 0) { 185 *func = NULL; 186 } else { 187 PetscCheck(ind < chunk.size, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Index %" PetscInt_FMT " not in [0, %" PetscInt_FMT ")", ind, chunk.size); 188 *func = ((void (**)()) & wf->funcs->array[chunk.start])[ind]; 189 } 190 PetscFunctionReturn(0); 191 } 192 193 /* Ignore a NULL func */ 194 PetscErrorCode PetscWeakFormSetIndexFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt ind, void (*func)()) 195 { 196 PetscFormKey key; 197 PetscChunk chunk; 198 199 PetscFunctionBegin; 200 if (!func) PetscFunctionReturn(0); 201 key.label = label; 202 key.value = value; 203 key.field = f; 204 key.part = part; 205 PetscCall(PetscHMapFormGet(ht, key, &chunk)); 206 if (chunk.size < 0) { 207 PetscCall(PetscChunkBufferCreateChunk(wf->funcs, ind + 1, &chunk)); 208 PetscCall(PetscHMapFormSet(ht, key, chunk)); 209 } else if (chunk.size <= ind) { 210 PetscCall(PetscChunkBufferEnlargeChunk(wf->funcs, ind - chunk.size + 1, &chunk)); 211 PetscCall(PetscHMapFormSet(ht, key, chunk)); 212 } 213 ((void (**)()) & wf->funcs->array[chunk.start])[ind] = func; 214 PetscFunctionReturn(0); 215 } 216 217 PetscErrorCode PetscWeakFormClearIndexFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt ind) 218 { 219 PetscFormKey key; 220 PetscChunk chunk; 221 222 PetscFunctionBegin; 223 key.label = label; 224 key.value = value; 225 key.field = f; 226 key.part = part; 227 PetscCall(PetscHMapFormGet(ht, key, &chunk)); 228 if (chunk.size < 0) { 229 PetscFunctionReturn(0); 230 } else if (!ind && chunk.size == 1) { 231 PetscCall(PetscHMapFormDel(ht, key)); 232 PetscFunctionReturn(0); 233 } else if (chunk.size <= ind) { 234 PetscFunctionReturn(0); 235 } 236 ((void (**)()) & wf->funcs->array[chunk.start])[ind] = NULL; 237 PetscFunctionReturn(0); 238 } 239 240 /*@ 241 PetscWeakFormCopy - Copy the pointwise functions to another PetscWeakForm 242 243 Not Collective 244 245 Input Parameter: 246 . wf - The original PetscWeakForm 247 248 Output Parameter: 249 . wfNew - The copy PetscWeakForm 250 251 Level: intermediate 252 253 .seealso: `PetscWeakFormCreate()`, `PetscWeakFormDestroy()` 254 @*/ 255 PetscErrorCode PetscWeakFormCopy(PetscWeakForm wf, PetscWeakForm wfNew) 256 { 257 PetscInt f; 258 259 PetscFunctionBegin; 260 wfNew->Nf = wf->Nf; 261 PetscCall(PetscChunkBufferDestroy(&wfNew->funcs)); 262 PetscCall(PetscChunkBufferDuplicate(wf->funcs, &wfNew->funcs)); 263 for (f = 0; f < PETSC_NUM_WF; ++f) { 264 PetscCall(PetscHMapFormDestroy(&wfNew->form[f])); 265 PetscCall(PetscHMapFormDuplicate(wf->form[f], &wfNew->form[f])); 266 } 267 PetscFunctionReturn(0); 268 } 269 270 /*@ 271 PetscWeakFormClear - Clear all functions from the PetscWeakForm 272 273 Not Collective 274 275 Input Parameter: 276 . wf - The original PetscWeakForm 277 278 Level: intermediate 279 280 .seealso: `PetscWeakFormCopy()`, `PetscWeakFormCreate()`, `PetscWeakFormDestroy()` 281 @*/ 282 PetscErrorCode PetscWeakFormClear(PetscWeakForm wf) 283 { 284 PetscInt f; 285 286 PetscFunctionBegin; 287 for (f = 0; f < PETSC_NUM_WF; ++f) PetscCall(PetscHMapFormClear(wf->form[f])); 288 PetscFunctionReturn(0); 289 } 290 291 static PetscErrorCode PetscWeakFormRewriteKeys_Internal(PetscWeakForm wf, PetscHMapForm hmap, DMLabel label, PetscInt Nv, const PetscInt values[]) 292 { 293 PetscFormKey *keys; 294 PetscInt n, i, v, off = 0; 295 296 PetscFunctionBegin; 297 PetscCall(PetscHMapFormGetSize(hmap, &n)); 298 PetscCall(PetscMalloc1(n, &keys)); 299 PetscCall(PetscHMapFormGetKeys(hmap, &off, keys)); 300 for (i = 0; i < n; ++i) { 301 if (keys[i].label == label) { 302 PetscBool clear = PETSC_TRUE; 303 void (**funcs)(); 304 PetscInt Nf; 305 306 PetscCall(PetscWeakFormGetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &Nf, &funcs)); 307 for (v = 0; v < Nv; ++v) { 308 PetscCall(PetscWeakFormSetFunction_Private(wf, hmap, keys[i].label, values[v], keys[i].field, keys[i].part, Nf, funcs)); 309 if (values[v] == keys[i].value) clear = PETSC_FALSE; 310 } 311 if (clear) PetscCall(PetscWeakFormSetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, 0, NULL)); 312 } 313 } 314 PetscCall(PetscFree(keys)); 315 PetscFunctionReturn(0); 316 } 317 318 /*@C 319 PetscWeakFormRewriteKeys - Change any key on the given label to use the new set of label values 320 321 Not Collective 322 323 Input Parameters: 324 + wf - The original PetscWeakForm 325 . label - The label to change keys for 326 . Nv - The number of new label values 327 - values - The set of new values to relabel keys with 328 329 Note: This is used internally when boundary label values are specified from the command line. 330 331 Level: intermediate 332 333 .seealso: `PetscWeakFormReplaceLabel()`, `PetscWeakFormCreate()`, `PetscWeakFormDestroy()` 334 @*/ 335 PetscErrorCode PetscWeakFormRewriteKeys(PetscWeakForm wf, DMLabel label, PetscInt Nv, const PetscInt values[]) 336 { 337 PetscInt f; 338 339 PetscFunctionBegin; 340 for (f = 0; f < PETSC_NUM_WF; ++f) PetscCall(PetscWeakFormRewriteKeys_Internal(wf, wf->form[f], label, Nv, values)); 341 PetscFunctionReturn(0); 342 } 343 344 static PetscErrorCode PetscWeakFormReplaceLabel_Internal(PetscWeakForm wf, PetscHMapForm hmap, DMLabel label) 345 { 346 PetscFormKey *keys; 347 PetscInt n, i, off = 0, maxFuncs = 0; 348 void (**tmpf)(); 349 const char *name = NULL; 350 351 PetscFunctionBegin; 352 if (label) PetscCall(PetscObjectGetName((PetscObject)label, &name)); 353 PetscCall(PetscHMapFormGetSize(hmap, &n)); 354 PetscCall(PetscMalloc1(n, &keys)); 355 PetscCall(PetscHMapFormGetKeys(hmap, &off, keys)); 356 for (i = 0; i < n; ++i) { 357 PetscBool match = PETSC_FALSE; 358 const char *lname = NULL; 359 360 if (label == keys[i].label) continue; 361 if (keys[i].label) PetscCall(PetscObjectGetName((PetscObject)keys[i].label, &lname)); 362 PetscCall(PetscStrcmp(name, lname, &match)); 363 if ((!name && !lname) || match) { 364 void (**funcs)(); 365 PetscInt Nf; 366 367 PetscCall(PetscWeakFormGetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &Nf, &funcs)); 368 maxFuncs = PetscMax(maxFuncs, Nf); 369 } 370 } 371 /* Need temp space because chunk buffer can be reallocated in SetFunction() call */ 372 PetscCall(PetscMalloc1(maxFuncs, &tmpf)); 373 for (i = 0; i < n; ++i) { 374 PetscBool match = PETSC_FALSE; 375 const char *lname = NULL; 376 377 if (label == keys[i].label) continue; 378 if (keys[i].label) PetscCall(PetscObjectGetName((PetscObject)keys[i].label, &lname)); 379 PetscCall(PetscStrcmp(name, lname, &match)); 380 if ((!name && !lname) || match) { 381 void (**funcs)(); 382 PetscInt Nf, j; 383 384 PetscCall(PetscWeakFormGetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &Nf, &funcs)); 385 for (j = 0; j < Nf; ++j) tmpf[j] = funcs[j]; 386 PetscCall(PetscWeakFormSetFunction_Private(wf, hmap, label, keys[i].value, keys[i].field, keys[i].part, Nf, tmpf)); 387 PetscCall(PetscWeakFormSetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, 0, NULL)); 388 } 389 } 390 PetscCall(PetscFree(tmpf)); 391 PetscCall(PetscFree(keys)); 392 PetscFunctionReturn(0); 393 } 394 395 /*@C 396 PetscWeakFormReplaceLabel - Change any key on a label of the same name to use the new label 397 398 Not Collective 399 400 Input Parameters: 401 + wf - The original PetscWeakForm 402 - label - The label to change keys for 403 404 Note: This is used internally when meshes are modified 405 406 Level: intermediate 407 408 .seealso: `PetscWeakFormRewriteKeys()`, `PetscWeakFormCreate()`, `PetscWeakFormDestroy()` 409 @*/ 410 PetscErrorCode PetscWeakFormReplaceLabel(PetscWeakForm wf, DMLabel label) 411 { 412 PetscInt f; 413 414 PetscFunctionBegin; 415 for (f = 0; f < PETSC_NUM_WF; ++f) PetscCall(PetscWeakFormReplaceLabel_Internal(wf, wf->form[f], label)); 416 PetscFunctionReturn(0); 417 } 418 419 PetscErrorCode PetscWeakFormClearIndex(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscWeakFormKind kind, PetscInt ind) 420 { 421 PetscFunctionBegin; 422 PetscCall(PetscWeakFormClearIndexFunction_Private(wf, wf->form[kind], label, val, f, part, ind)); 423 PetscFunctionReturn(0); 424 } 425 426 PetscErrorCode PetscWeakFormGetObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt *n, void (***obj)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 427 { 428 PetscFunctionBegin; 429 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, n, (void (***)(void))obj)); 430 PetscFunctionReturn(0); 431 } 432 433 PetscErrorCode PetscWeakFormSetObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt n, void (**obj)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 434 { 435 PetscFunctionBegin; 436 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, n, (void (**)(void))obj)); 437 PetscFunctionReturn(0); 438 } 439 440 PetscErrorCode PetscWeakFormAddObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, void (*obj)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 441 { 442 PetscFunctionBegin; 443 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, (void (*)(void))obj)); 444 PetscFunctionReturn(0); 445 } 446 447 PetscErrorCode PetscWeakFormGetIndexObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt ind, void (**obj)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 448 { 449 PetscFunctionBegin; 450 PetscCall(PetscWeakFormGetIndexFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, ind, (void (**)(void))obj)); 451 PetscFunctionReturn(0); 452 } 453 454 PetscErrorCode PetscWeakFormSetIndexObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt ind, void (*obj)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 455 { 456 PetscFunctionBegin; 457 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, ind, (void (*)(void))obj)); 458 PetscFunctionReturn(0); 459 } 460 461 PetscErrorCode PetscWeakFormGetResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt *n0, void (***f0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n1, void (***f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 462 { 463 PetscFunctionBegin; 464 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, n0, (void (***)(void))f0)); 465 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, n1, (void (***)(void))f1)); 466 PetscFunctionReturn(0); 467 } 468 469 PetscErrorCode PetscWeakFormAddResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, void (*f0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 470 { 471 PetscFunctionBegin; 472 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, (void (*)(void))f0)); 473 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, (void (*)(void))f1)); 474 PetscFunctionReturn(0); 475 } 476 477 PetscErrorCode PetscWeakFormSetResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt n0, void (**f0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n1, void (**f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 478 { 479 PetscFunctionBegin; 480 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, n0, (void (**)(void))f0)); 481 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, n1, (void (**)(void))f1)); 482 PetscFunctionReturn(0); 483 } 484 485 PetscErrorCode PetscWeakFormSetIndexResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt i0, void (*f0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i1, void (*f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 486 { 487 PetscFunctionBegin; 488 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, i0, (void (*)(void))f0)); 489 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, i1, (void (*)(void))f1)); 490 PetscFunctionReturn(0); 491 } 492 493 PetscErrorCode PetscWeakFormGetBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt *n0, void (***f0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n1, void (***f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 494 { 495 PetscFunctionBegin; 496 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, n0, (void (***)(void))f0)); 497 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, n1, (void (***)(void))f1)); 498 PetscFunctionReturn(0); 499 } 500 501 PetscErrorCode PetscWeakFormAddBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, void (*f0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 502 { 503 PetscFunctionBegin; 504 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, (void (*)(void))f0)); 505 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, (void (*)(void))f1)); 506 PetscFunctionReturn(0); 507 } 508 509 PetscErrorCode PetscWeakFormSetBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt n0, void (**f0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n1, void (**f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 510 { 511 PetscFunctionBegin; 512 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, n0, (void (**)(void))f0)); 513 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, n1, (void (**)(void))f1)); 514 PetscFunctionReturn(0); 515 } 516 517 PetscErrorCode PetscWeakFormSetIndexBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt i0, void (*f0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i1, void (*f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 518 { 519 PetscFunctionBegin; 520 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, i0, (void (*)(void))f0)); 521 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, i1, (void (*)(void))f1)); 522 PetscFunctionReturn(0); 523 } 524 525 PetscErrorCode PetscWeakFormHasJacobian(PetscWeakForm wf, PetscBool *hasJac) 526 { 527 PetscInt n0, n1, n2, n3; 528 529 PetscFunctionBegin; 530 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 531 PetscValidBoolPointer(hasJac, 2); 532 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_G0], &n0)); 533 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_G1], &n1)); 534 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_G2], &n2)); 535 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_G3], &n3)); 536 *hasJac = n0 + n1 + n2 + n3 ? PETSC_TRUE : PETSC_FALSE; 537 PetscFunctionReturn(0); 538 } 539 540 PetscErrorCode PetscWeakFormGetJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt *n0, void (***g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n1, void (***g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n2, void (***g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n3, void (***g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 541 { 542 PetscInt find = f * wf->Nf + g; 543 544 PetscFunctionBegin; 545 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, n0, (void (***)(void))g0)); 546 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, n1, (void (***)(void))g1)); 547 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, n2, (void (***)(void))g2)); 548 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, n3, (void (***)(void))g3)); 549 PetscFunctionReturn(0); 550 } 551 552 PetscErrorCode PetscWeakFormAddJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, void (*g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 553 { 554 PetscInt find = f * wf->Nf + g; 555 556 PetscFunctionBegin; 557 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, (void (*)(void))g0)); 558 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, (void (*)(void))g1)); 559 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, (void (*)(void))g2)); 560 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, (void (*)(void))g3)); 561 PetscFunctionReturn(0); 562 } 563 564 PetscErrorCode PetscWeakFormSetJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt n0, void (**g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n1, void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n2, void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n3, void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 565 { 566 PetscInt find = f * wf->Nf + g; 567 568 PetscFunctionBegin; 569 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, n0, (void (**)(void))g0)); 570 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, n1, (void (**)(void))g1)); 571 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, n2, (void (**)(void))g2)); 572 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, n3, (void (**)(void))g3)); 573 PetscFunctionReturn(0); 574 } 575 576 PetscErrorCode PetscWeakFormSetIndexJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt i0, void (*g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i1, void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i2, void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i3, void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 577 { 578 PetscInt find = f * wf->Nf + g; 579 580 PetscFunctionBegin; 581 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, i0, (void (*)(void))g0)); 582 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, i1, (void (*)(void))g1)); 583 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, i2, (void (*)(void))g2)); 584 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, i3, (void (*)(void))g3)); 585 PetscFunctionReturn(0); 586 } 587 588 PetscErrorCode PetscWeakFormHasJacobianPreconditioner(PetscWeakForm wf, PetscBool *hasJacPre) 589 { 590 PetscInt n0, n1, n2, n3; 591 592 PetscFunctionBegin; 593 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 594 PetscValidBoolPointer(hasJacPre, 2); 595 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GP0], &n0)); 596 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GP1], &n1)); 597 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GP2], &n2)); 598 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GP3], &n3)); 599 *hasJacPre = n0 + n1 + n2 + n3 ? PETSC_TRUE : PETSC_FALSE; 600 PetscFunctionReturn(0); 601 } 602 603 PetscErrorCode PetscWeakFormGetJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt *n0, void (***g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n1, void (***g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n2, void (***g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n3, void (***g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 604 { 605 PetscInt find = f * wf->Nf + g; 606 607 PetscFunctionBegin; 608 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, n0, (void (***)(void))g0)); 609 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, n1, (void (***)(void))g1)); 610 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, n2, (void (***)(void))g2)); 611 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, n3, (void (***)(void))g3)); 612 PetscFunctionReturn(0); 613 } 614 615 PetscErrorCode PetscWeakFormAddJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, void (*g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 616 { 617 PetscInt find = f * wf->Nf + g; 618 619 PetscFunctionBegin; 620 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, (void (*)(void))g0)); 621 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, (void (*)(void))g1)); 622 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, (void (*)(void))g2)); 623 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, (void (*)(void))g3)); 624 PetscFunctionReturn(0); 625 } 626 627 PetscErrorCode PetscWeakFormSetJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt n0, void (**g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n1, void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n2, void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n3, void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 628 { 629 PetscInt find = f * wf->Nf + g; 630 631 PetscFunctionBegin; 632 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, n0, (void (**)(void))g0)); 633 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, n1, (void (**)(void))g1)); 634 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, n2, (void (**)(void))g2)); 635 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, n3, (void (**)(void))g3)); 636 PetscFunctionReturn(0); 637 } 638 639 PetscErrorCode PetscWeakFormSetIndexJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt i0, void (*g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i1, void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i2, void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i3, void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 640 { 641 PetscInt find = f * wf->Nf + g; 642 643 PetscFunctionBegin; 644 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, i0, (void (*)(void))g0)); 645 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, i1, (void (*)(void))g1)); 646 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, i2, (void (*)(void))g2)); 647 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, i3, (void (*)(void))g3)); 648 PetscFunctionReturn(0); 649 } 650 651 PetscErrorCode PetscWeakFormHasBdJacobian(PetscWeakForm wf, PetscBool *hasJac) 652 { 653 PetscInt n0, n1, n2, n3; 654 655 PetscFunctionBegin; 656 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 657 PetscValidBoolPointer(hasJac, 2); 658 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDG0], &n0)); 659 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDG1], &n1)); 660 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDG2], &n2)); 661 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDG3], &n3)); 662 *hasJac = n0 + n1 + n2 + n3 ? PETSC_TRUE : PETSC_FALSE; 663 PetscFunctionReturn(0); 664 } 665 666 PetscErrorCode PetscWeakFormGetBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt *n0, void (***g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n1, void (***g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n2, void (***g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n3, void (***g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 667 { 668 PetscInt find = f * wf->Nf + g; 669 670 PetscFunctionBegin; 671 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, n0, (void (***)(void))g0)); 672 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, n1, (void (***)(void))g1)); 673 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, n2, (void (***)(void))g2)); 674 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, n3, (void (***)(void))g3)); 675 PetscFunctionReturn(0); 676 } 677 678 PetscErrorCode PetscWeakFormAddBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, void (*g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 679 { 680 PetscInt find = f * wf->Nf + g; 681 682 PetscFunctionBegin; 683 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, (void (*)(void))g0)); 684 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, (void (*)(void))g1)); 685 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, (void (*)(void))g2)); 686 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, (void (*)(void))g3)); 687 PetscFunctionReturn(0); 688 } 689 690 PetscErrorCode PetscWeakFormSetBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt n0, void (**g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n1, void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n2, void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n3, void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 691 { 692 PetscInt find = f * wf->Nf + g; 693 694 PetscFunctionBegin; 695 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, n0, (void (**)(void))g0)); 696 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, n1, (void (**)(void))g1)); 697 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, n2, (void (**)(void))g2)); 698 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, n3, (void (**)(void))g3)); 699 PetscFunctionReturn(0); 700 } 701 702 PetscErrorCode PetscWeakFormSetIndexBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt i0, void (*g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i1, void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i2, void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i3, void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 703 { 704 PetscInt find = f * wf->Nf + g; 705 706 PetscFunctionBegin; 707 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, i0, (void (*)(void))g0)); 708 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, i1, (void (*)(void))g1)); 709 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, i2, (void (*)(void))g2)); 710 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, i3, (void (*)(void))g3)); 711 PetscFunctionReturn(0); 712 } 713 714 PetscErrorCode PetscWeakFormHasBdJacobianPreconditioner(PetscWeakForm wf, PetscBool *hasJacPre) 715 { 716 PetscInt n0, n1, n2, n3; 717 718 PetscFunctionBegin; 719 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 720 PetscValidBoolPointer(hasJacPre, 2); 721 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP0], &n0)); 722 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP1], &n1)); 723 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP2], &n2)); 724 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP3], &n3)); 725 *hasJacPre = n0 + n1 + n2 + n3 ? PETSC_TRUE : PETSC_FALSE; 726 PetscFunctionReturn(0); 727 } 728 729 PetscErrorCode PetscWeakFormGetBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt *n0, void (***g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n1, void (***g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n2, void (***g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n3, void (***g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 730 { 731 PetscInt find = f * wf->Nf + g; 732 733 PetscFunctionBegin; 734 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, n0, (void (***)(void))g0)); 735 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, n1, (void (***)(void))g1)); 736 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, n2, (void (***)(void))g2)); 737 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, n3, (void (***)(void))g3)); 738 PetscFunctionReturn(0); 739 } 740 741 PetscErrorCode PetscWeakFormAddBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, void (*g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 742 { 743 PetscInt find = f * wf->Nf + g; 744 745 PetscFunctionBegin; 746 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, (void (*)(void))g0)); 747 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, (void (*)(void))g1)); 748 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, (void (*)(void))g2)); 749 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, (void (*)(void))g3)); 750 PetscFunctionReturn(0); 751 } 752 753 PetscErrorCode PetscWeakFormSetBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt n0, void (**g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n1, void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n2, void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n3, void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 754 { 755 PetscInt find = f * wf->Nf + g; 756 757 PetscFunctionBegin; 758 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, n0, (void (**)(void))g0)); 759 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, n1, (void (**)(void))g1)); 760 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, n2, (void (**)(void))g2)); 761 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, n3, (void (**)(void))g3)); 762 PetscFunctionReturn(0); 763 } 764 765 PetscErrorCode PetscWeakFormSetIndexBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt i0, void (*g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i1, void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i2, void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i3, void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 766 { 767 PetscInt find = f * wf->Nf + g; 768 769 PetscFunctionBegin; 770 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, i0, (void (*)(void))g0)); 771 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, i1, (void (*)(void))g1)); 772 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, i2, (void (*)(void))g2)); 773 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, i3, (void (*)(void))g3)); 774 PetscFunctionReturn(0); 775 } 776 777 PetscErrorCode PetscWeakFormHasDynamicJacobian(PetscWeakForm wf, PetscBool *hasDynJac) 778 { 779 PetscInt n0, n1, n2, n3; 780 781 PetscFunctionBegin; 782 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 783 PetscValidBoolPointer(hasDynJac, 2); 784 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GT0], &n0)); 785 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GT1], &n1)); 786 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GT2], &n2)); 787 PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GT3], &n3)); 788 *hasDynJac = n0 + n1 + n2 + n3 ? PETSC_TRUE : PETSC_FALSE; 789 PetscFunctionReturn(0); 790 } 791 792 PetscErrorCode PetscWeakFormGetDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt *n0, void (***g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n1, void (***g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n2, void (***g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt *n3, void (***g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 793 { 794 PetscInt find = f * wf->Nf + g; 795 796 PetscFunctionBegin; 797 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, n0, (void (***)(void))g0)); 798 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, n1, (void (***)(void))g1)); 799 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, n2, (void (***)(void))g2)); 800 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, n3, (void (***)(void))g3)); 801 PetscFunctionReturn(0); 802 } 803 804 PetscErrorCode PetscWeakFormAddDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, void (*g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 805 { 806 PetscInt find = f * wf->Nf + g; 807 808 PetscFunctionBegin; 809 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, (void (*)(void))g0)); 810 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, (void (*)(void))g1)); 811 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, (void (*)(void))g2)); 812 PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, (void (*)(void))g3)); 813 PetscFunctionReturn(0); 814 } 815 816 PetscErrorCode PetscWeakFormSetDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt n0, void (**g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n1, void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n2, void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt n3, void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 817 { 818 PetscInt find = f * wf->Nf + g; 819 820 PetscFunctionBegin; 821 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, n0, (void (**)(void))g0)); 822 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, n1, (void (**)(void))g1)); 823 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, n2, (void (**)(void))g2)); 824 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, n3, (void (**)(void))g3)); 825 PetscFunctionReturn(0); 826 } 827 828 PetscErrorCode PetscWeakFormSetIndexDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, PetscInt i0, void (*g0)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i1, void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i2, void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), PetscInt i3, void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 829 { 830 PetscInt find = f * wf->Nf + g; 831 832 PetscFunctionBegin; 833 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, i0, (void (*)(void))g0)); 834 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, i1, (void (*)(void))g1)); 835 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, i2, (void (*)(void))g2)); 836 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, i3, (void (*)(void))g3)); 837 PetscFunctionReturn(0); 838 } 839 840 PetscErrorCode PetscWeakFormGetRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt *n, void (***r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *)) 841 { 842 PetscFunctionBegin; 843 PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, n, (void (***)(void))r)); 844 PetscFunctionReturn(0); 845 } 846 847 PetscErrorCode PetscWeakFormSetRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt n, void (**r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *)) 848 { 849 PetscFunctionBegin; 850 PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, n, (void (**)(void))r)); 851 PetscFunctionReturn(0); 852 } 853 854 PetscErrorCode PetscWeakFormSetIndexRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt i, void (*r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *)) 855 { 856 PetscFunctionBegin; 857 PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, i, (void (*)(void))r)); 858 PetscFunctionReturn(0); 859 } 860 861 /*@ 862 PetscWeakFormGetNumFields - Returns the number of fields 863 864 Not collective 865 866 Input Parameter: 867 . wf - The PetscWeakForm object 868 869 Output Parameter: 870 . Nf - The number of fields 871 872 Level: beginner 873 874 .seealso: `PetscWeakFormSetNumFields()`, `PetscWeakFormCreate()` 875 @*/ 876 PetscErrorCode PetscWeakFormGetNumFields(PetscWeakForm wf, PetscInt *Nf) 877 { 878 PetscFunctionBegin; 879 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 880 PetscValidIntPointer(Nf, 2); 881 *Nf = wf->Nf; 882 PetscFunctionReturn(0); 883 } 884 885 /*@ 886 PetscWeakFormSetNumFields - Sets the number of fields 887 888 Not collective 889 890 Input Parameters: 891 + wf - The PetscWeakForm object 892 - Nf - The number of fields 893 894 Level: beginner 895 896 .seealso: `PetscWeakFormGetNumFields()`, `PetscWeakFormCreate()` 897 @*/ 898 PetscErrorCode PetscWeakFormSetNumFields(PetscWeakForm wf, PetscInt Nf) 899 { 900 PetscFunctionBegin; 901 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 902 wf->Nf = Nf; 903 PetscFunctionReturn(0); 904 } 905 906 /*@ 907 PetscWeakFormDestroy - Destroys a PetscWeakForm object 908 909 Collective on wf 910 911 Input Parameter: 912 . wf - the PetscWeakForm object to destroy 913 914 Level: developer 915 916 .seealso `PetscWeakFormCreate()`, `PetscWeakFormView()` 917 @*/ 918 PetscErrorCode PetscWeakFormDestroy(PetscWeakForm *wf) 919 { 920 PetscInt f; 921 922 PetscFunctionBegin; 923 if (!*wf) PetscFunctionReturn(0); 924 PetscValidHeaderSpecific((*wf), PETSCWEAKFORM_CLASSID, 1); 925 926 if (--((PetscObject)(*wf))->refct > 0) { 927 *wf = NULL; 928 PetscFunctionReturn(0); 929 } 930 ((PetscObject)(*wf))->refct = 0; 931 PetscCall(PetscChunkBufferDestroy(&(*wf)->funcs)); 932 for (f = 0; f < PETSC_NUM_WF; ++f) PetscCall(PetscHMapFormDestroy(&(*wf)->form[f])); 933 PetscCall(PetscFree((*wf)->form)); 934 PetscCall(PetscHeaderDestroy(wf)); 935 PetscFunctionReturn(0); 936 } 937 938 static PetscErrorCode PetscWeakFormViewTable_Ascii(PetscWeakForm wf, PetscViewer viewer, PetscBool splitField, const char tableName[], PetscHMapForm map) 939 { 940 PetscInt Nf = wf->Nf, Nk, k; 941 942 PetscFunctionBegin; 943 PetscCall(PetscHMapFormGetSize(map, &Nk)); 944 if (Nk) { 945 PetscFormKey *keys; 946 void (**funcs)(void) = NULL; 947 const char **names; 948 PetscInt *values, *idx1, *idx2, *idx; 949 PetscBool showPart = PETSC_FALSE, showPointer = PETSC_FALSE; 950 PetscInt off = 0; 951 952 PetscCall(PetscMalloc6(Nk, &keys, Nk, &names, Nk, &values, Nk, &idx1, Nk, &idx2, Nk, &idx)); 953 PetscCall(PetscHMapFormGetKeys(map, &off, keys)); 954 /* Sort keys by label name and value */ 955 { 956 /* First sort values */ 957 for (k = 0; k < Nk; ++k) { 958 values[k] = keys[k].value; 959 idx1[k] = k; 960 } 961 PetscCall(PetscSortIntWithPermutation(Nk, values, idx1)); 962 /* If the string sort is stable, it will be sorted correctly overall */ 963 for (k = 0; k < Nk; ++k) { 964 if (keys[idx1[k]].label) PetscCall(PetscObjectGetName((PetscObject)keys[idx1[k]].label, &names[k])); 965 else names[k] = ""; 966 idx2[k] = k; 967 } 968 PetscCall(PetscSortStrWithPermutation(Nk, names, idx2)); 969 for (k = 0; k < Nk; ++k) { 970 if (keys[k].label) PetscCall(PetscObjectGetName((PetscObject)keys[k].label, &names[k])); 971 else names[k] = ""; 972 idx[k] = idx1[idx2[k]]; 973 } 974 } 975 PetscCall(PetscViewerASCIIPrintf(viewer, "%s\n", tableName)); 976 PetscCall(PetscViewerASCIIPushTab(viewer)); 977 for (k = 0; k < Nk; ++k) { 978 if (keys[k].part != 0) showPart = PETSC_TRUE; 979 } 980 for (k = 0; k < Nk; ++k) { 981 const PetscInt i = idx[k]; 982 PetscInt n, f; 983 984 if (keys[i].label) { 985 if (showPointer) PetscCall(PetscViewerASCIIPrintf(viewer, "(%s:%p, %" PetscInt_FMT ") ", names[i], keys[i].label, keys[i].value)); 986 else PetscCall(PetscViewerASCIIPrintf(viewer, "(%s, %" PetscInt_FMT ") ", names[i], keys[i].value)); 987 } 988 PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE)); 989 if (splitField) PetscCall(PetscViewerASCIIPrintf(viewer, "(%" PetscInt_FMT ", %" PetscInt_FMT ") ", keys[i].field / Nf, keys[i].field % Nf)); 990 else PetscCall(PetscViewerASCIIPrintf(viewer, "(%" PetscInt_FMT ") ", keys[i].field)); 991 if (showPart) PetscCall(PetscViewerASCIIPrintf(viewer, "(%" PetscInt_FMT ") ", keys[i].part)); 992 PetscCall(PetscWeakFormGetFunction_Private(wf, map, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &n, &funcs)); 993 for (f = 0; f < n; ++f) { 994 char *fname; 995 size_t len, l; 996 997 if (f > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", ")); 998 PetscCall(PetscDLAddr(funcs[f], &fname)); 999 if (fname) { 1000 /* Eliminate argument types */ 1001 PetscCall(PetscStrlen(fname, &len)); 1002 for (l = 0; l < len; ++l) 1003 if (fname[l] == '(') { 1004 fname[l] = '\0'; 1005 break; 1006 } 1007 PetscCall(PetscViewerASCIIPrintf(viewer, "%s", fname)); 1008 } else if (showPointer) { 1009 PetscCall(PetscViewerASCIIPrintf(viewer, "%p", funcs[f])); 1010 } 1011 PetscCall(PetscFree(fname)); 1012 } 1013 PetscCall(PetscViewerASCIIPrintf(viewer, "\n")); 1014 PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE)); 1015 } 1016 PetscCall(PetscViewerASCIIPopTab(viewer)); 1017 PetscCall(PetscFree6(keys, names, values, idx1, idx2, idx)); 1018 } 1019 PetscFunctionReturn(0); 1020 } 1021 1022 static PetscErrorCode PetscWeakFormView_Ascii(PetscWeakForm wf, PetscViewer viewer) 1023 { 1024 PetscViewerFormat format; 1025 PetscInt f; 1026 1027 PetscFunctionBegin; 1028 PetscCall(PetscViewerGetFormat(viewer, &format)); 1029 PetscCall(PetscViewerASCIIPrintf(viewer, "Weak Form System with %" PetscInt_FMT " fields\n", wf->Nf)); 1030 PetscCall(PetscViewerASCIIPushTab(viewer)); 1031 for (f = 0; f < PETSC_NUM_WF; ++f) PetscCall(PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE, PetscWeakFormKinds[f], wf->form[f])); 1032 PetscCall(PetscViewerASCIIPopTab(viewer)); 1033 PetscFunctionReturn(0); 1034 } 1035 1036 /*@C 1037 PetscWeakFormView - Views a PetscWeakForm 1038 1039 Collective on wf 1040 1041 Input Parameters: 1042 + wf - the PetscWeakForm object to view 1043 - v - the viewer 1044 1045 Level: developer 1046 1047 .seealso `PetscWeakFormDestroy()`, `PetscWeakFormCreate()` 1048 @*/ 1049 PetscErrorCode PetscWeakFormView(PetscWeakForm wf, PetscViewer v) 1050 { 1051 PetscBool iascii; 1052 1053 PetscFunctionBegin; 1054 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 1055 if (!v) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)wf), &v)); 1056 else PetscValidHeaderSpecific(v, PETSC_VIEWER_CLASSID, 2); 1057 PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &iascii)); 1058 if (iascii) PetscCall(PetscWeakFormView_Ascii(wf, v)); 1059 PetscTryTypeMethod(wf, view, v); 1060 PetscFunctionReturn(0); 1061 } 1062 1063 /*@ 1064 PetscWeakFormCreate - Creates an empty PetscWeakForm object. 1065 1066 Collective 1067 1068 Input Parameter: 1069 . comm - The communicator for the PetscWeakForm object 1070 1071 Output Parameter: 1072 . wf - The PetscWeakForm object 1073 1074 Level: beginner 1075 1076 .seealso: `PetscDS`, `PetscWeakFormDestroy()` 1077 @*/ 1078 PetscErrorCode PetscWeakFormCreate(MPI_Comm comm, PetscWeakForm *wf) 1079 { 1080 PetscWeakForm p; 1081 PetscInt f; 1082 1083 PetscFunctionBegin; 1084 PetscValidPointer(wf, 2); 1085 *wf = NULL; 1086 PetscCall(PetscDSInitializePackage()); 1087 1088 PetscCall(PetscHeaderCreate(p, PETSCWEAKFORM_CLASSID, "PetscWeakForm", "Weak Form System", "PetscWeakForm", comm, PetscWeakFormDestroy, PetscWeakFormView)); 1089 1090 p->Nf = 0; 1091 PetscCall(PetscChunkBufferCreate(sizeof(&PetscWeakFormCreate), 2, &p->funcs)); 1092 PetscCall(PetscMalloc1(PETSC_NUM_WF, &p->form)); 1093 for (f = 0; f < PETSC_NUM_WF; ++f) PetscCall(PetscHMapFormCreate(&p->form[f])); 1094 *wf = p; 1095 PetscFunctionReturn(0); 1096 } 1097