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 PetscErrorCode ierr; 10 11 PetscFunctionBegin; 12 ierr = PetscNew(buffer);CHKERRQ(ierr); 13 ierr = PetscCalloc1(expected*unitbytes, &(*buffer)->array);CHKERRQ(ierr); 14 (*buffer)->size = expected; 15 (*buffer)->unitbytes = unitbytes; 16 (*buffer)->alloc = expected*unitbytes; 17 PetscFunctionReturn(0); 18 } 19 20 static PetscErrorCode PetscChunkBufferDuplicate(PetscChunkBuffer *buffer, PetscChunkBuffer **bufferNew) 21 { 22 PetscErrorCode ierr; 23 24 PetscFunctionBegin; 25 ierr = PetscNew(bufferNew);CHKERRQ(ierr); 26 ierr = PetscCalloc1(buffer->size*buffer->unitbytes, &(*bufferNew)->array);CHKERRQ(ierr); 27 ierr = PetscMemcpy((*bufferNew)->array, buffer->array, buffer->size*buffer->unitbytes);CHKERRQ(ierr); 28 (*bufferNew)->size = buffer->size; 29 (*bufferNew)->unitbytes = buffer->unitbytes; 30 (*bufferNew)->alloc = buffer->size*buffer->unitbytes; 31 PetscFunctionReturn(0); 32 } 33 34 static PetscErrorCode PetscChunkBufferDestroy(PetscChunkBuffer **buffer) 35 { 36 PetscErrorCode ierr; 37 38 PetscFunctionBegin; 39 ierr = PetscFree((*buffer)->array);CHKERRQ(ierr); 40 ierr = PetscFree(*buffer);CHKERRQ(ierr); 41 PetscFunctionReturn(0); 42 } 43 44 static PetscErrorCode PetscChunkBufferCreateChunk(PetscChunkBuffer *buffer, PetscInt size, PetscChunk *chunk) 45 { 46 PetscErrorCode ierr; 47 48 PetscFunctionBegin; 49 if ((buffer->size + size)*buffer->unitbytes > buffer->alloc) { 50 char *tmp; 51 52 if (!buffer->alloc) buffer->alloc = (buffer->size + size)*buffer->unitbytes; 53 while ((buffer->size + size)*buffer->unitbytes > buffer->alloc) buffer->alloc *= 2; 54 ierr = PetscMalloc(buffer->alloc, &tmp);CHKERRQ(ierr); 55 ierr = PetscMemcpy(tmp, buffer->array, buffer->size*buffer->unitbytes);CHKERRQ(ierr); 56 ierr = PetscFree(buffer->array);CHKERRQ(ierr); 57 buffer->array = tmp; 58 } 59 chunk->start = buffer->size*buffer->unitbytes; 60 chunk->size = size; 61 chunk->reserved = size; 62 buffer->size += size; 63 PetscFunctionReturn(0); 64 } 65 66 static PetscErrorCode PetscChunkBufferEnlargeChunk(PetscChunkBuffer *buffer, PetscInt size, PetscChunk *chunk) 67 { 68 size_t siz = size; 69 PetscErrorCode ierr; 70 71 PetscFunctionBegin; 72 if (chunk->size + size > chunk->reserved) { 73 PetscChunk newchunk; 74 PetscInt reserved = chunk->size; 75 76 /* TODO Here if we had a chunk list, we could update them all to reclaim unused space */ 77 while (reserved < chunk->size+size) reserved *= 2; 78 ierr = PetscChunkBufferCreateChunk(buffer, (size_t) reserved, &newchunk);CHKERRQ(ierr); 79 newchunk.size = chunk->size+size; 80 ierr = PetscMemcpy(&buffer->array[newchunk.start], &buffer->array[chunk->start], chunk->size * buffer->unitbytes);CHKERRQ(ierr); 81 *chunk = newchunk; 82 } else { 83 chunk->size += siz; 84 } 85 PetscFunctionReturn(0); 86 } 87 88 /*@C 89 PetscFormKeySort - Sorts an array of PetscFormKey in place in increasing order. 90 91 Not Collective 92 93 Input Parameters: 94 + n - number of values 95 - X - array of PetscFormKey 96 97 Level: intermediate 98 99 .seealso: PetscIntSortSemiOrdered(), PetscSortInt() 100 @*/ 101 PetscErrorCode PetscFormKeySort(PetscInt n, PetscFormKey arr[]) 102 { 103 PetscErrorCode ierr; 104 105 PetscFunctionBegin; 106 if (n <= 1) PetscFunctionReturn(0); 107 PetscValidPointer(arr, 2); 108 ierr = PetscTimSort(n, arr, sizeof(PetscFormKey), Compare_PetscFormKey_Private, NULL);CHKERRQ(ierr); 109 PetscFunctionReturn(0); 110 } 111 112 PetscErrorCode PetscWeakFormGetFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt *n, void (***func)()) 113 { 114 PetscFormKey key; 115 PetscChunk chunk; 116 PetscErrorCode ierr; 117 118 PetscFunctionBegin; 119 key.label = label; key.value = value; key.field = f; key.part = part; 120 ierr = PetscHMapFormGet(ht, key, &chunk);CHKERRQ(ierr); 121 if (chunk.size < 0) {*n = 0; *func = NULL;} 122 else {*n = chunk.size; *func = (void (**)()) &wf->funcs->array[chunk.start];} 123 PetscFunctionReturn(0); 124 } 125 126 /* A NULL argument for func causes this to clear the key */ 127 PetscErrorCode PetscWeakFormSetFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt n, void (**func)()) 128 { 129 PetscFormKey key; 130 PetscChunk chunk; 131 PetscInt i; 132 PetscErrorCode ierr; 133 134 PetscFunctionBegin; 135 key.label = label; key.value = value; key.field = f; key.part = part; 136 if (!func) { 137 ierr = PetscHMapFormDel(ht, key);CHKERRQ(ierr); 138 PetscFunctionReturn(0); 139 } else { 140 ierr = PetscHMapFormGet(ht, key, &chunk);CHKERRQ(ierr); 141 } 142 if (chunk.size < 0) { 143 ierr = PetscChunkBufferCreateChunk(wf->funcs, n, &chunk);CHKERRQ(ierr); 144 ierr = PetscHMapFormSet(ht, key, chunk);CHKERRQ(ierr); 145 } else if (chunk.size <= n) { 146 ierr = PetscChunkBufferEnlargeChunk(wf->funcs, n - chunk.size, &chunk);CHKERRQ(ierr); 147 ierr = PetscHMapFormSet(ht, key, chunk);CHKERRQ(ierr); 148 } 149 for (i = 0; i < n; ++i) ((void (**)()) &wf->funcs->array[chunk.start])[i] = func[i]; 150 PetscFunctionReturn(0); 151 } 152 153 PetscErrorCode PetscWeakFormAddFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, void (*func)()) 154 { 155 PetscFormKey key; 156 PetscChunk chunk; 157 PetscErrorCode ierr; 158 159 PetscFunctionBegin; 160 if (!func) PetscFunctionReturn(0); 161 key.label = label; key.value = value; key.field = f; key.part = part; 162 ierr = PetscHMapFormGet(ht, key, &chunk);CHKERRQ(ierr); 163 if (chunk.size < 0) { 164 ierr = PetscChunkBufferCreateChunk(wf->funcs, 1, &chunk);CHKERRQ(ierr); 165 ierr = PetscHMapFormSet(ht, key, chunk);CHKERRQ(ierr); 166 ((void (**)()) &wf->funcs->array[chunk.start])[0] = func; 167 } else { 168 ierr = PetscChunkBufferEnlargeChunk(wf->funcs, 1, &chunk);CHKERRQ(ierr); 169 ierr = PetscHMapFormSet(ht, key, chunk);CHKERRQ(ierr); 170 ((void (**)()) &wf->funcs->array[chunk.start])[chunk.size-1] = func; 171 } 172 PetscFunctionReturn(0); 173 } 174 175 PetscErrorCode PetscWeakFormGetIndexFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt ind, void (**func)()) 176 { 177 PetscFormKey key; 178 PetscChunk chunk; 179 PetscErrorCode ierr; 180 181 PetscFunctionBegin; 182 key.label = label; key.value = value; key.field = f; key.part = part; 183 ierr = PetscHMapFormGet(ht, key, &chunk);CHKERRQ(ierr); 184 if (chunk.size < 0) {*func = NULL;} 185 else { 186 if (ind >= chunk.size) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Index %D not in [0, %D)", ind, chunk.size); 187 *func = ((void (**)()) &wf->funcs->array[chunk.start])[ind]; 188 } 189 PetscFunctionReturn(0); 190 } 191 192 /* Ignore a NULL func */ 193 PetscErrorCode PetscWeakFormSetIndexFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt ind, void (*func)()) 194 { 195 PetscFormKey key; 196 PetscChunk chunk; 197 PetscErrorCode ierr; 198 199 PetscFunctionBegin; 200 if (!func) PetscFunctionReturn(0); 201 key.label = label; key.value = value; key.field = f; key.part = part; 202 ierr = PetscHMapFormGet(ht, key, &chunk);CHKERRQ(ierr); 203 if (chunk.size < 0) { 204 ierr = PetscChunkBufferCreateChunk(wf->funcs, ind+1, &chunk);CHKERRQ(ierr); 205 ierr = PetscHMapFormSet(ht, key, chunk);CHKERRQ(ierr); 206 } else if (chunk.size <= ind) { 207 ierr = PetscChunkBufferEnlargeChunk(wf->funcs, ind - chunk.size + 1, &chunk);CHKERRQ(ierr); 208 ierr = PetscHMapFormSet(ht, key, chunk);CHKERRQ(ierr); 209 } 210 ((void (**)()) &wf->funcs->array[chunk.start])[ind] = func; 211 PetscFunctionReturn(0); 212 } 213 214 PetscErrorCode PetscWeakFormClearIndexFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt ind) 215 { 216 PetscFormKey key; 217 PetscChunk chunk; 218 PetscErrorCode ierr; 219 220 PetscFunctionBegin; 221 key.label = label; key.value = value; key.field = f; key.part = part; 222 ierr = PetscHMapFormGet(ht, key, &chunk);CHKERRQ(ierr); 223 if (chunk.size < 0) { 224 PetscFunctionReturn(0); 225 } else if (!ind && chunk.size == 1) { 226 ierr = PetscHMapFormDel(ht, key);CHKERRQ(ierr); 227 PetscFunctionReturn(0); 228 } else if (chunk.size <= ind) { 229 PetscFunctionReturn(0); 230 } 231 ((void (**)()) &wf->funcs->array[chunk.start])[ind] = NULL; 232 PetscFunctionReturn(0); 233 } 234 235 /*@ 236 PetscWeakFormCopy - Copy the pointwise functions to another PetscWeakForm 237 238 Not Collective 239 240 Input Parameter: 241 . wf - The original PetscWeakForm 242 243 Output Parameter: 244 . wfNew - The copy PetscWeakForm 245 246 Level: intermediate 247 248 .seealso: PetscWeakFormCreate(), PetscWeakFormDestroy() 249 @*/ 250 PetscErrorCode PetscWeakFormCopy(PetscWeakForm wf, PetscWeakForm wfNew) 251 { 252 PetscInt f; 253 PetscErrorCode ierr; 254 255 PetscFunctionBegin; 256 wfNew->Nf = wf->Nf; 257 ierr = PetscChunkBufferDestroy(&wfNew->funcs);CHKERRQ(ierr); 258 ierr = PetscChunkBufferDuplicate(wf->funcs, &wfNew->funcs);CHKERRQ(ierr); 259 for (f = 0; f < PETSC_NUM_WF; ++f) { 260 ierr = PetscHMapFormDestroy(&wfNew->form[f]);CHKERRQ(ierr); 261 ierr = PetscHMapFormDuplicate(wf->form[f], &wfNew->form[f]);CHKERRQ(ierr); 262 } 263 PetscFunctionReturn(0); 264 } 265 266 /*@ 267 PetscWeakFormClear - Clear all functions from the PetscWeakForm 268 269 Not Collective 270 271 Input Parameter: 272 . wf - The original PetscWeakForm 273 274 Level: intermediate 275 276 .seealso: PetscWeakFormCopy(), PetscWeakFormCreate(), PetscWeakFormDestroy() 277 @*/ 278 PetscErrorCode PetscWeakFormClear(PetscWeakForm wf) 279 { 280 PetscInt f; 281 PetscErrorCode ierr; 282 283 PetscFunctionBegin; 284 for (f = 0; f < PETSC_NUM_WF; ++f) {ierr = PetscHMapFormClear(wf->form[f]);CHKERRQ(ierr);} 285 PetscFunctionReturn(0); 286 } 287 288 static PetscErrorCode PetscWeakFormRewriteKeys_Internal(PetscWeakForm wf, PetscHMapForm hmap, DMLabel label, PetscInt Nv, const PetscInt values[]) 289 { 290 PetscFormKey *keys; 291 PetscInt n, i, v, off = 0; 292 PetscErrorCode ierr; 293 294 PetscFunctionBegin; 295 ierr = PetscHMapFormGetSize(hmap, &n);CHKERRQ(ierr); 296 ierr = PetscMalloc1(n, &keys);CHKERRQ(ierr); 297 ierr = PetscHMapFormGetKeys(hmap, &off, keys);CHKERRQ(ierr); 298 for (i = 0; i < n; ++i) { 299 if (keys[i].label == label) { 300 PetscBool clear = PETSC_TRUE; 301 void (**funcs)(); 302 PetscInt Nf; 303 304 ierr = PetscWeakFormGetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &Nf, &funcs);CHKERRQ(ierr); 305 for (v = 0; v < Nv; ++v) { 306 ierr = PetscWeakFormSetFunction_Private(wf, hmap, keys[i].label, values[v], keys[i].field, keys[i].part, Nf, funcs);CHKERRQ(ierr); 307 if (values[v] == keys[i].value) clear = PETSC_FALSE; 308 } 309 if (clear) {ierr = PetscWeakFormSetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, 0, NULL);CHKERRQ(ierr);} 310 } 311 } 312 ierr = PetscFree(keys);CHKERRQ(ierr); 313 PetscFunctionReturn(0); 314 } 315 316 /*@C 317 PetscWeakFormRewriteKeys - Change any key on the given label to use the new set of label values 318 319 Not Collective 320 321 Input Parameters: 322 + wf - The original PetscWeakForm 323 . label - The label to change keys for 324 . Nv - The number of new label values 325 - values - The set of new values to relabel keys with 326 327 Note: This is used internally when boundary label values are specified from the command line. 328 329 Level: intermediate 330 331 .seealso: PetscWeakFormReplaceLabel(), PetscWeakFormCreate(), PetscWeakFormDestroy() 332 @*/ 333 PetscErrorCode PetscWeakFormRewriteKeys(PetscWeakForm wf, DMLabel label, PetscInt Nv, const PetscInt values[]) 334 { 335 PetscInt f; 336 PetscErrorCode ierr; 337 338 PetscFunctionBegin; 339 for (f = 0; f < PETSC_NUM_WF; ++f) {ierr = PetscWeakFormRewriteKeys_Internal(wf, wf->form[f], label, Nv, values);CHKERRQ(ierr);} 340 PetscFunctionReturn(0); 341 } 342 343 static PetscErrorCode PetscWeakFormReplaceLabel_Internal(PetscWeakForm wf, PetscHMapForm hmap, DMLabel label) 344 { 345 PetscFormKey *keys; 346 PetscInt n, i, off = 0, maxFuncs = 0; 347 void (**tmpf)(); 348 const char *name = NULL; 349 PetscErrorCode ierr; 350 351 PetscFunctionBegin; 352 if (label) {ierr = PetscObjectGetName((PetscObject) label, &name);CHKERRQ(ierr);} 353 ierr = PetscHMapFormGetSize(hmap, &n);CHKERRQ(ierr); 354 ierr = PetscMalloc1(n, &keys);CHKERRQ(ierr); 355 ierr = PetscHMapFormGetKeys(hmap, &off, keys);CHKERRQ(ierr); 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) {ierr = PetscObjectGetName((PetscObject) keys[i].label, &lname);CHKERRQ(ierr);} 362 ierr = PetscStrcmp(name, lname, &match);CHKERRQ(ierr); 363 if ((!name && !lname) || match) { 364 void (**funcs)(); 365 PetscInt Nf; 366 367 ierr = PetscWeakFormGetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &Nf, &funcs);CHKERRQ(ierr); 368 maxFuncs = PetscMax(maxFuncs, Nf); 369 } 370 } 371 /* Need temp space because chunk buffer can be reallocated in SetFunction() call */ 372 ierr = PetscMalloc1(maxFuncs, &tmpf);CHKERRQ(ierr); 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) {ierr = PetscObjectGetName((PetscObject) keys[i].label, &lname);CHKERRQ(ierr);} 379 ierr = PetscStrcmp(name, lname, &match);CHKERRQ(ierr); 380 if ((!name && !lname) || match) { 381 void (**funcs)(); 382 PetscInt Nf, j; 383 384 ierr = PetscWeakFormGetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &Nf, &funcs);CHKERRQ(ierr); 385 for (j = 0; j < Nf; ++j) tmpf[j] = funcs[j]; 386 ierr = PetscWeakFormSetFunction_Private(wf, hmap, label, keys[i].value, keys[i].field, keys[i].part, Nf, tmpf);CHKERRQ(ierr); 387 ierr = PetscWeakFormSetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, 0, NULL);CHKERRQ(ierr); 388 } 389 } 390 ierr = PetscFree(tmpf);CHKERRQ(ierr); 391 ierr = PetscFree(keys);CHKERRQ(ierr); 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 PetscErrorCode ierr; 414 415 PetscFunctionBegin; 416 for (f = 0; f < PETSC_NUM_WF; ++f) {ierr = PetscWeakFormReplaceLabel_Internal(wf, wf->form[f], label);CHKERRQ(ierr);} 417 PetscFunctionReturn(0); 418 } 419 420 PetscErrorCode PetscWeakFormClearIndex(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscWeakFormKind kind, PetscInt ind) 421 { 422 PetscErrorCode ierr; 423 424 PetscFunctionBegin; 425 ierr = PetscWeakFormClearIndexFunction_Private(wf, wf->form[kind], label, val, f, part, ind);CHKERRQ(ierr); 426 PetscFunctionReturn(0); 427 } 428 429 PetscErrorCode PetscWeakFormGetObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt *n, 430 void (***obj)(PetscInt, PetscInt, PetscInt, 431 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 432 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 433 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 434 { 435 PetscErrorCode ierr; 436 437 PetscFunctionBegin; 438 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, n, (void (***)(void)) obj);CHKERRQ(ierr); 439 PetscFunctionReturn(0); 440 } 441 442 PetscErrorCode PetscWeakFormSetObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt n, 443 void (**obj)(PetscInt, PetscInt, PetscInt, 444 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 445 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 446 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 447 { 448 PetscErrorCode ierr; 449 450 PetscFunctionBegin; 451 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, n, (void (**)(void)) obj);CHKERRQ(ierr); 452 PetscFunctionReturn(0); 453 } 454 455 PetscErrorCode PetscWeakFormAddObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 456 void (*obj)(PetscInt, PetscInt, PetscInt, 457 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 458 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 459 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 460 { 461 PetscErrorCode ierr; 462 463 PetscFunctionBegin; 464 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, (void (*)(void)) obj);CHKERRQ(ierr); 465 PetscFunctionReturn(0); 466 } 467 468 PetscErrorCode PetscWeakFormGetIndexObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt ind, 469 void (**obj)(PetscInt, PetscInt, PetscInt, 470 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 471 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 472 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 473 { 474 PetscErrorCode ierr; 475 476 PetscFunctionBegin; 477 ierr = PetscWeakFormGetIndexFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, ind, (void (**)(void)) obj);CHKERRQ(ierr); 478 PetscFunctionReturn(0); 479 } 480 481 PetscErrorCode PetscWeakFormSetIndexObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt ind, 482 void (*obj)(PetscInt, PetscInt, PetscInt, 483 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 484 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 485 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 486 { 487 PetscErrorCode ierr; 488 489 PetscFunctionBegin; 490 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, ind, (void (*)(void)) obj);CHKERRQ(ierr); 491 PetscFunctionReturn(0); 492 } 493 494 PetscErrorCode PetscWeakFormGetResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 495 PetscInt *n0, 496 void (***f0)(PetscInt, PetscInt, PetscInt, 497 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 498 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 499 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 500 PetscInt *n1, 501 void (***f1)(PetscInt, PetscInt, PetscInt, 502 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 503 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 504 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 505 { 506 PetscErrorCode ierr; 507 508 PetscFunctionBegin; 509 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, n0, (void (***)(void)) f0);CHKERRQ(ierr); 510 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, n1, (void (***)(void)) f1);CHKERRQ(ierr); 511 PetscFunctionReturn(0); 512 } 513 514 PetscErrorCode PetscWeakFormAddResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 515 void (*f0)(PetscInt, PetscInt, PetscInt, 516 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 517 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 518 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 519 void (*f1)(PetscInt, PetscInt, PetscInt, 520 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 521 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 522 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 523 { 524 PetscErrorCode ierr; 525 526 PetscFunctionBegin; 527 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, (void (*)(void)) f0);CHKERRQ(ierr); 528 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, (void (*)(void)) f1);CHKERRQ(ierr); 529 PetscFunctionReturn(0); 530 } 531 532 PetscErrorCode PetscWeakFormSetResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 533 PetscInt n0, 534 void (**f0)(PetscInt, PetscInt, PetscInt, 535 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 536 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 537 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 538 PetscInt n1, 539 void (**f1)(PetscInt, PetscInt, PetscInt, 540 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 541 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 542 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 543 { 544 PetscErrorCode ierr; 545 546 PetscFunctionBegin; 547 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, n0, (void (**)(void)) f0);CHKERRQ(ierr); 548 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, n1, (void (**)(void)) f1);CHKERRQ(ierr); 549 PetscFunctionReturn(0); 550 } 551 552 PetscErrorCode PetscWeakFormSetIndexResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 553 PetscInt i0, 554 void (*f0)(PetscInt, PetscInt, PetscInt, 555 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 556 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 557 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 558 PetscInt i1, 559 void (*f1)(PetscInt, PetscInt, PetscInt, 560 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 561 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 562 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 563 { 564 PetscErrorCode ierr; 565 566 PetscFunctionBegin; 567 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, i0, (void (*)(void)) f0);CHKERRQ(ierr); 568 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, i1, (void (*)(void)) f1);CHKERRQ(ierr); 569 PetscFunctionReturn(0); 570 } 571 572 PetscErrorCode PetscWeakFormGetBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 573 PetscInt *n0, 574 void (***f0)(PetscInt, PetscInt, PetscInt, 575 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 576 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 577 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 578 PetscInt *n1, 579 void (***f1)(PetscInt, PetscInt, PetscInt, 580 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 581 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 582 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 583 { 584 PetscErrorCode ierr; 585 586 PetscFunctionBegin; 587 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, n0, (void (***)(void)) f0);CHKERRQ(ierr); 588 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, n1, (void (***)(void)) f1);CHKERRQ(ierr); 589 PetscFunctionReturn(0); 590 } 591 592 PetscErrorCode PetscWeakFormAddBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 593 void (*f0)(PetscInt, PetscInt, PetscInt, 594 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 595 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 596 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 597 void (*f1)(PetscInt, PetscInt, PetscInt, 598 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 599 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 600 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 601 { 602 PetscErrorCode ierr; 603 604 PetscFunctionBegin; 605 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, (void (*)(void)) f0);CHKERRQ(ierr); 606 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, (void (*)(void)) f1);CHKERRQ(ierr); 607 PetscFunctionReturn(0); 608 } 609 610 PetscErrorCode PetscWeakFormSetBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 611 PetscInt n0, 612 void (**f0)(PetscInt, PetscInt, PetscInt, 613 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 614 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 615 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 616 PetscInt n1, 617 void (**f1)(PetscInt, PetscInt, PetscInt, 618 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 619 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 620 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 621 { 622 PetscErrorCode ierr; 623 624 PetscFunctionBegin; 625 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, n0, (void (**)(void)) f0);CHKERRQ(ierr); 626 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, n1, (void (**)(void)) f1);CHKERRQ(ierr); 627 PetscFunctionReturn(0); 628 } 629 630 PetscErrorCode PetscWeakFormSetIndexBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 631 PetscInt i0, 632 void (*f0)(PetscInt, PetscInt, PetscInt, 633 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 634 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 635 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 636 PetscInt i1, 637 void (*f1)(PetscInt, PetscInt, PetscInt, 638 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 639 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 640 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 641 { 642 PetscErrorCode ierr; 643 644 PetscFunctionBegin; 645 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, i0, (void (*)(void)) f0);CHKERRQ(ierr); 646 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, i1, (void (*)(void)) f1);CHKERRQ(ierr); 647 PetscFunctionReturn(0); 648 } 649 650 PetscErrorCode PetscWeakFormHasJacobian(PetscWeakForm wf, PetscBool *hasJac) 651 { 652 PetscInt n0, n1, n2, n3; 653 PetscErrorCode ierr; 654 655 PetscFunctionBegin; 656 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 657 PetscValidBoolPointer(hasJac, 2); 658 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_G0], &n0);CHKERRQ(ierr); 659 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_G1], &n1);CHKERRQ(ierr); 660 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_G2], &n2);CHKERRQ(ierr); 661 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_G3], &n3);CHKERRQ(ierr); 662 *hasJac = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE; 663 PetscFunctionReturn(0); 664 } 665 666 PetscErrorCode PetscWeakFormGetJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 667 PetscInt *n0, 668 void (***g0)(PetscInt, PetscInt, PetscInt, 669 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 670 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 671 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 672 PetscInt *n1, 673 void (***g1)(PetscInt, PetscInt, PetscInt, 674 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 675 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 676 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 677 PetscInt *n2, 678 void (***g2)(PetscInt, PetscInt, PetscInt, 679 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 680 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 681 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 682 PetscInt *n3, 683 void (***g3)(PetscInt, PetscInt, PetscInt, 684 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 685 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 686 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 687 { 688 PetscInt find = f*wf->Nf + g; 689 PetscErrorCode ierr; 690 691 PetscFunctionBegin; 692 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, n0, (void (***)(void)) g0);CHKERRQ(ierr); 693 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, n1, (void (***)(void)) g1);CHKERRQ(ierr); 694 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, n2, (void (***)(void)) g2);CHKERRQ(ierr); 695 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, n3, (void (***)(void)) g3);CHKERRQ(ierr); 696 PetscFunctionReturn(0); 697 } 698 699 PetscErrorCode PetscWeakFormAddJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 700 void (*g0)(PetscInt, PetscInt, PetscInt, 701 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 702 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 703 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 704 void (*g1)(PetscInt, PetscInt, PetscInt, 705 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 706 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 707 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 708 void (*g2)(PetscInt, PetscInt, PetscInt, 709 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 710 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 711 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 712 void (*g3)(PetscInt, PetscInt, PetscInt, 713 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 714 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 715 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 716 { 717 PetscInt find = f*wf->Nf + g; 718 PetscErrorCode ierr; 719 720 PetscFunctionBegin; 721 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, (void (*)(void)) g0);CHKERRQ(ierr); 722 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, (void (*)(void)) g1);CHKERRQ(ierr); 723 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, (void (*)(void)) g2);CHKERRQ(ierr); 724 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, (void (*)(void)) g3);CHKERRQ(ierr); 725 PetscFunctionReturn(0); 726 } 727 728 PetscErrorCode PetscWeakFormSetJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 729 PetscInt n0, 730 void (**g0)(PetscInt, PetscInt, PetscInt, 731 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 732 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 733 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 734 PetscInt n1, 735 void (**g1)(PetscInt, PetscInt, PetscInt, 736 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 737 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 738 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 739 PetscInt n2, 740 void (**g2)(PetscInt, PetscInt, PetscInt, 741 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 742 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 743 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 744 PetscInt n3, 745 void (**g3)(PetscInt, PetscInt, PetscInt, 746 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 747 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 748 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 749 { 750 PetscInt find = f*wf->Nf + g; 751 PetscErrorCode ierr; 752 753 PetscFunctionBegin; 754 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, n0, (void (**)(void)) g0);CHKERRQ(ierr); 755 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, n1, (void (**)(void)) g1);CHKERRQ(ierr); 756 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, n2, (void (**)(void)) g2);CHKERRQ(ierr); 757 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, n3, (void (**)(void)) g3);CHKERRQ(ierr); 758 PetscFunctionReturn(0); 759 } 760 761 PetscErrorCode PetscWeakFormSetIndexJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 762 PetscInt i0, 763 void (*g0)(PetscInt, PetscInt, PetscInt, 764 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 765 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 766 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 767 PetscInt i1, 768 void (*g1)(PetscInt, PetscInt, PetscInt, 769 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 770 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 771 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 772 PetscInt i2, 773 void (*g2)(PetscInt, PetscInt, PetscInt, 774 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 775 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 776 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 777 PetscInt i3, 778 void (*g3)(PetscInt, PetscInt, PetscInt, 779 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 780 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 781 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 782 { 783 PetscInt find = f*wf->Nf + g; 784 PetscErrorCode ierr; 785 786 PetscFunctionBegin; 787 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, i0, (void (*)(void)) g0);CHKERRQ(ierr); 788 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, i1, (void (*)(void)) g1);CHKERRQ(ierr); 789 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, i2, (void (*)(void)) g2);CHKERRQ(ierr); 790 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, i3, (void (*)(void)) g3);CHKERRQ(ierr); 791 PetscFunctionReturn(0); 792 } 793 794 PetscErrorCode PetscWeakFormHasJacobianPreconditioner(PetscWeakForm wf, PetscBool *hasJacPre) 795 { 796 PetscInt n0, n1, n2, n3; 797 PetscErrorCode ierr; 798 799 PetscFunctionBegin; 800 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 801 PetscValidBoolPointer(hasJacPre, 2); 802 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GP0], &n0);CHKERRQ(ierr); 803 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GP1], &n1);CHKERRQ(ierr); 804 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GP2], &n2);CHKERRQ(ierr); 805 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GP3], &n3);CHKERRQ(ierr); 806 *hasJacPre = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE; 807 PetscFunctionReturn(0); 808 } 809 810 PetscErrorCode PetscWeakFormGetJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 811 PetscInt *n0, 812 void (***g0)(PetscInt, PetscInt, PetscInt, 813 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 814 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 815 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 816 PetscInt *n1, 817 void (***g1)(PetscInt, PetscInt, PetscInt, 818 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 819 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 820 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 821 PetscInt *n2, 822 void (***g2)(PetscInt, PetscInt, PetscInt, 823 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 824 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 825 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 826 PetscInt *n3, 827 void (***g3)(PetscInt, PetscInt, PetscInt, 828 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 829 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 830 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 831 { 832 PetscInt find = f*wf->Nf + g; 833 PetscErrorCode ierr; 834 835 PetscFunctionBegin; 836 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, n0, (void (***)(void)) g0);CHKERRQ(ierr); 837 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, n1, (void (***)(void)) g1);CHKERRQ(ierr); 838 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, n2, (void (***)(void)) g2);CHKERRQ(ierr); 839 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, n3, (void (***)(void)) g3);CHKERRQ(ierr); 840 PetscFunctionReturn(0); 841 } 842 843 PetscErrorCode PetscWeakFormAddJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 844 void (*g0)(PetscInt, PetscInt, PetscInt, 845 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 846 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 847 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 848 void (*g1)(PetscInt, PetscInt, PetscInt, 849 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 850 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 851 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 852 void (*g2)(PetscInt, PetscInt, PetscInt, 853 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 854 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 855 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 856 void (*g3)(PetscInt, PetscInt, PetscInt, 857 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 858 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 859 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 860 { 861 PetscInt find = f*wf->Nf + g; 862 PetscErrorCode ierr; 863 864 PetscFunctionBegin; 865 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, (void (*)(void)) g0);CHKERRQ(ierr); 866 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, (void (*)(void)) g1);CHKERRQ(ierr); 867 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, (void (*)(void)) g2);CHKERRQ(ierr); 868 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, (void (*)(void)) g3);CHKERRQ(ierr); 869 PetscFunctionReturn(0); 870 } 871 872 PetscErrorCode PetscWeakFormSetJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 873 PetscInt n0, 874 void (**g0)(PetscInt, PetscInt, PetscInt, 875 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 876 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 877 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 878 PetscInt n1, 879 void (**g1)(PetscInt, PetscInt, PetscInt, 880 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 881 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 882 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 883 PetscInt n2, 884 void (**g2)(PetscInt, PetscInt, PetscInt, 885 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 886 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 887 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 888 PetscInt n3, 889 void (**g3)(PetscInt, PetscInt, PetscInt, 890 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 891 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 892 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 893 { 894 PetscInt find = f*wf->Nf + g; 895 PetscErrorCode ierr; 896 897 PetscFunctionBegin; 898 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, n0, (void (**)(void)) g0);CHKERRQ(ierr); 899 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, n1, (void (**)(void)) g1);CHKERRQ(ierr); 900 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, n2, (void (**)(void)) g2);CHKERRQ(ierr); 901 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, n3, (void (**)(void)) g3);CHKERRQ(ierr); 902 PetscFunctionReturn(0); 903 } 904 905 PetscErrorCode PetscWeakFormSetIndexJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 906 PetscInt i0, 907 void (*g0)(PetscInt, PetscInt, PetscInt, 908 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 909 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 910 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 911 PetscInt i1, 912 void (*g1)(PetscInt, PetscInt, PetscInt, 913 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 914 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 915 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 916 PetscInt i2, 917 void (*g2)(PetscInt, PetscInt, PetscInt, 918 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 919 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 920 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 921 PetscInt i3, 922 void (*g3)(PetscInt, PetscInt, PetscInt, 923 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 924 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 925 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 926 { 927 PetscInt find = f*wf->Nf + g; 928 PetscErrorCode ierr; 929 930 PetscFunctionBegin; 931 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, i0, (void (*)(void)) g0);CHKERRQ(ierr); 932 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, i1, (void (*)(void)) g1);CHKERRQ(ierr); 933 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, i2, (void (*)(void)) g2);CHKERRQ(ierr); 934 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, i3, (void (*)(void)) g3);CHKERRQ(ierr); 935 PetscFunctionReturn(0); 936 } 937 938 PetscErrorCode PetscWeakFormHasBdJacobian(PetscWeakForm wf, PetscBool *hasJac) 939 { 940 PetscInt n0, n1, n2, n3; 941 PetscErrorCode ierr; 942 943 PetscFunctionBegin; 944 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 945 PetscValidBoolPointer(hasJac, 2); 946 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDG0], &n0);CHKERRQ(ierr); 947 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDG1], &n1);CHKERRQ(ierr); 948 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDG2], &n2);CHKERRQ(ierr); 949 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDG3], &n3);CHKERRQ(ierr); 950 *hasJac = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE; 951 PetscFunctionReturn(0); 952 } 953 954 PetscErrorCode PetscWeakFormGetBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 955 PetscInt *n0, 956 void (***g0)(PetscInt, PetscInt, PetscInt, 957 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 958 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 959 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 960 PetscInt *n1, 961 void (***g1)(PetscInt, PetscInt, PetscInt, 962 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 963 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 964 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 965 PetscInt *n2, 966 void (***g2)(PetscInt, PetscInt, PetscInt, 967 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 968 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 969 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 970 PetscInt *n3, 971 void (***g3)(PetscInt, PetscInt, PetscInt, 972 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 973 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 974 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 975 { 976 PetscInt find = f*wf->Nf + g; 977 PetscErrorCode ierr; 978 979 PetscFunctionBegin; 980 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, n0, (void (***)(void)) g0);CHKERRQ(ierr); 981 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, n1, (void (***)(void)) g1);CHKERRQ(ierr); 982 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, n2, (void (***)(void)) g2);CHKERRQ(ierr); 983 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, n3, (void (***)(void)) g3);CHKERRQ(ierr); 984 PetscFunctionReturn(0); 985 } 986 987 PetscErrorCode PetscWeakFormAddBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 988 void (*g0)(PetscInt, PetscInt, PetscInt, 989 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 990 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 991 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 992 void (*g1)(PetscInt, PetscInt, PetscInt, 993 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 994 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 995 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 996 void (*g2)(PetscInt, PetscInt, PetscInt, 997 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 998 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 999 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1000 void (*g3)(PetscInt, PetscInt, PetscInt, 1001 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1002 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1003 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1004 { 1005 PetscInt find = f*wf->Nf + g; 1006 PetscErrorCode ierr; 1007 1008 PetscFunctionBegin; 1009 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, (void (*)(void)) g0);CHKERRQ(ierr); 1010 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, (void (*)(void)) g1);CHKERRQ(ierr); 1011 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, (void (*)(void)) g2);CHKERRQ(ierr); 1012 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, (void (*)(void)) g3);CHKERRQ(ierr); 1013 PetscFunctionReturn(0); 1014 } 1015 1016 PetscErrorCode PetscWeakFormSetBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1017 PetscInt n0, 1018 void (**g0)(PetscInt, PetscInt, PetscInt, 1019 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1020 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1021 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1022 PetscInt n1, 1023 void (**g1)(PetscInt, PetscInt, PetscInt, 1024 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1025 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1026 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1027 PetscInt n2, 1028 void (**g2)(PetscInt, PetscInt, PetscInt, 1029 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1030 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1031 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1032 PetscInt n3, 1033 void (**g3)(PetscInt, PetscInt, PetscInt, 1034 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1035 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1036 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1037 { 1038 PetscInt find = f*wf->Nf + g; 1039 PetscErrorCode ierr; 1040 1041 PetscFunctionBegin; 1042 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, n0, (void (**)(void)) g0);CHKERRQ(ierr); 1043 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, n1, (void (**)(void)) g1);CHKERRQ(ierr); 1044 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, n2, (void (**)(void)) g2);CHKERRQ(ierr); 1045 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, n3, (void (**)(void)) g3);CHKERRQ(ierr); 1046 PetscFunctionReturn(0); 1047 } 1048 1049 PetscErrorCode PetscWeakFormSetIndexBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1050 PetscInt i0, 1051 void (*g0)(PetscInt, PetscInt, PetscInt, 1052 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1053 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1054 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1055 PetscInt i1, 1056 void (*g1)(PetscInt, PetscInt, PetscInt, 1057 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1058 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1059 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1060 PetscInt i2, 1061 void (*g2)(PetscInt, PetscInt, PetscInt, 1062 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1063 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1064 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1065 PetscInt i3, 1066 void (*g3)(PetscInt, PetscInt, PetscInt, 1067 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1068 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1069 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1070 { 1071 PetscInt find = f*wf->Nf + g; 1072 PetscErrorCode ierr; 1073 1074 PetscFunctionBegin; 1075 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, i0, (void (*)(void)) g0);CHKERRQ(ierr); 1076 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, i1, (void (*)(void)) g1);CHKERRQ(ierr); 1077 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, i2, (void (*)(void)) g2);CHKERRQ(ierr); 1078 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, i3, (void (*)(void)) g3);CHKERRQ(ierr); 1079 PetscFunctionReturn(0); 1080 } 1081 1082 PetscErrorCode PetscWeakFormHasBdJacobianPreconditioner(PetscWeakForm wf, PetscBool *hasJacPre) 1083 { 1084 PetscInt n0, n1, n2, n3; 1085 PetscErrorCode ierr; 1086 1087 PetscFunctionBegin; 1088 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 1089 PetscValidBoolPointer(hasJacPre, 2); 1090 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP0], &n0);CHKERRQ(ierr); 1091 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP1], &n1);CHKERRQ(ierr); 1092 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP2], &n2);CHKERRQ(ierr); 1093 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP3], &n3);CHKERRQ(ierr); 1094 *hasJacPre = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE; 1095 PetscFunctionReturn(0); 1096 } 1097 1098 PetscErrorCode PetscWeakFormGetBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1099 PetscInt *n0, 1100 void (***g0)(PetscInt, PetscInt, PetscInt, 1101 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1102 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1103 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1104 PetscInt *n1, 1105 void (***g1)(PetscInt, PetscInt, PetscInt, 1106 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1107 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1108 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1109 PetscInt *n2, 1110 void (***g2)(PetscInt, PetscInt, PetscInt, 1111 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1112 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1113 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1114 PetscInt *n3, 1115 void (***g3)(PetscInt, PetscInt, PetscInt, 1116 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1117 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1118 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1119 { 1120 PetscInt find = f*wf->Nf + g; 1121 PetscErrorCode ierr; 1122 1123 PetscFunctionBegin; 1124 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, n0, (void (***)(void)) g0);CHKERRQ(ierr); 1125 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, n1, (void (***)(void)) g1);CHKERRQ(ierr); 1126 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, n2, (void (***)(void)) g2);CHKERRQ(ierr); 1127 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, n3, (void (***)(void)) g3);CHKERRQ(ierr); 1128 PetscFunctionReturn(0); 1129 } 1130 1131 PetscErrorCode PetscWeakFormAddBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1132 void (*g0)(PetscInt, PetscInt, PetscInt, 1133 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1134 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1135 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1136 void (*g1)(PetscInt, PetscInt, PetscInt, 1137 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1138 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1139 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1140 void (*g2)(PetscInt, PetscInt, PetscInt, 1141 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1142 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1143 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1144 void (*g3)(PetscInt, PetscInt, PetscInt, 1145 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1146 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1147 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1148 { 1149 PetscInt find = f*wf->Nf + g; 1150 PetscErrorCode ierr; 1151 1152 PetscFunctionBegin; 1153 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, (void (*)(void)) g0);CHKERRQ(ierr); 1154 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, (void (*)(void)) g1);CHKERRQ(ierr); 1155 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, (void (*)(void)) g2);CHKERRQ(ierr); 1156 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, (void (*)(void)) g3);CHKERRQ(ierr); 1157 PetscFunctionReturn(0); 1158 } 1159 1160 PetscErrorCode PetscWeakFormSetBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1161 PetscInt n0, 1162 void (**g0)(PetscInt, PetscInt, PetscInt, 1163 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1164 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1165 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1166 PetscInt n1, 1167 void (**g1)(PetscInt, PetscInt, PetscInt, 1168 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1169 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1170 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1171 PetscInt n2, 1172 void (**g2)(PetscInt, PetscInt, PetscInt, 1173 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1174 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1175 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1176 PetscInt n3, 1177 void (**g3)(PetscInt, PetscInt, PetscInt, 1178 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1179 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1180 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1181 { 1182 PetscInt find = f*wf->Nf + g; 1183 PetscErrorCode ierr; 1184 1185 PetscFunctionBegin; 1186 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, n0, (void (**)(void)) g0);CHKERRQ(ierr); 1187 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, n1, (void (**)(void)) g1);CHKERRQ(ierr); 1188 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, n2, (void (**)(void)) g2);CHKERRQ(ierr); 1189 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, n3, (void (**)(void)) g3);CHKERRQ(ierr); 1190 PetscFunctionReturn(0); 1191 } 1192 1193 PetscErrorCode PetscWeakFormSetIndexBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1194 PetscInt i0, 1195 void (*g0)(PetscInt, PetscInt, PetscInt, 1196 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1197 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1198 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1199 PetscInt i1, 1200 void (*g1)(PetscInt, PetscInt, PetscInt, 1201 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1202 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1203 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1204 PetscInt i2, 1205 void (*g2)(PetscInt, PetscInt, PetscInt, 1206 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1207 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1208 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1209 PetscInt i3, 1210 void (*g3)(PetscInt, PetscInt, PetscInt, 1211 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1212 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1213 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1214 { 1215 PetscInt find = f*wf->Nf + g; 1216 PetscErrorCode ierr; 1217 1218 PetscFunctionBegin; 1219 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, i0, (void (*)(void)) g0);CHKERRQ(ierr); 1220 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, i1, (void (*)(void)) g1);CHKERRQ(ierr); 1221 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, i2, (void (*)(void)) g2);CHKERRQ(ierr); 1222 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, i3, (void (*)(void)) g3);CHKERRQ(ierr); 1223 PetscFunctionReturn(0); 1224 } 1225 1226 PetscErrorCode PetscWeakFormHasDynamicJacobian(PetscWeakForm wf, PetscBool *hasDynJac) 1227 { 1228 PetscInt n0, n1, n2, n3; 1229 PetscErrorCode ierr; 1230 1231 PetscFunctionBegin; 1232 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 1233 PetscValidBoolPointer(hasDynJac, 2); 1234 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GT0], &n0);CHKERRQ(ierr); 1235 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GT1], &n1);CHKERRQ(ierr); 1236 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GT2], &n2);CHKERRQ(ierr); 1237 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GT3], &n3);CHKERRQ(ierr); 1238 *hasDynJac = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE; 1239 PetscFunctionReturn(0); 1240 } 1241 1242 PetscErrorCode PetscWeakFormGetDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1243 PetscInt *n0, 1244 void (***g0)(PetscInt, PetscInt, PetscInt, 1245 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1246 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1247 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1248 PetscInt *n1, 1249 void (***g1)(PetscInt, PetscInt, PetscInt, 1250 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1251 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1252 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1253 PetscInt *n2, 1254 void (***g2)(PetscInt, PetscInt, PetscInt, 1255 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1256 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1257 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1258 PetscInt *n3, 1259 void (***g3)(PetscInt, PetscInt, PetscInt, 1260 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1261 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1262 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1263 { 1264 PetscInt find = f*wf->Nf + g; 1265 PetscErrorCode ierr; 1266 1267 PetscFunctionBegin; 1268 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, n0, (void (***)(void)) g0);CHKERRQ(ierr); 1269 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, n1, (void (***)(void)) g1);CHKERRQ(ierr); 1270 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, n2, (void (***)(void)) g2);CHKERRQ(ierr); 1271 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, n3, (void (***)(void)) g3);CHKERRQ(ierr); 1272 PetscFunctionReturn(0); 1273 } 1274 1275 PetscErrorCode PetscWeakFormAddDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1276 void (*g0)(PetscInt, PetscInt, PetscInt, 1277 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1278 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1279 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1280 void (*g1)(PetscInt, PetscInt, PetscInt, 1281 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1282 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1283 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1284 void (*g2)(PetscInt, PetscInt, PetscInt, 1285 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1286 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1287 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1288 void (*g3)(PetscInt, PetscInt, PetscInt, 1289 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1290 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1291 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1292 { 1293 PetscInt find = f*wf->Nf + g; 1294 PetscErrorCode ierr; 1295 1296 PetscFunctionBegin; 1297 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, (void (*)(void)) g0);CHKERRQ(ierr); 1298 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, (void (*)(void)) g1);CHKERRQ(ierr); 1299 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, (void (*)(void)) g2);CHKERRQ(ierr); 1300 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, (void (*)(void)) g3);CHKERRQ(ierr); 1301 PetscFunctionReturn(0); 1302 } 1303 1304 PetscErrorCode PetscWeakFormSetDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1305 PetscInt n0, 1306 void (**g0)(PetscInt, PetscInt, PetscInt, 1307 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1308 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1309 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1310 PetscInt n1, 1311 void (**g1)(PetscInt, PetscInt, PetscInt, 1312 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1313 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1314 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1315 PetscInt n2, 1316 void (**g2)(PetscInt, PetscInt, PetscInt, 1317 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1318 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1319 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1320 PetscInt n3, 1321 void (**g3)(PetscInt, PetscInt, PetscInt, 1322 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1323 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1324 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1325 { 1326 PetscInt find = f*wf->Nf + g; 1327 PetscErrorCode ierr; 1328 1329 PetscFunctionBegin; 1330 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, n0, (void (**)(void)) g0);CHKERRQ(ierr); 1331 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, n1, (void (**)(void)) g1);CHKERRQ(ierr); 1332 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, n2, (void (**)(void)) g2);CHKERRQ(ierr); 1333 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, n3, (void (**)(void)) g3);CHKERRQ(ierr); 1334 PetscFunctionReturn(0); 1335 } 1336 1337 PetscErrorCode PetscWeakFormSetIndexDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1338 PetscInt i0, 1339 void (*g0)(PetscInt, PetscInt, PetscInt, 1340 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1341 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1342 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1343 PetscInt i1, 1344 void (*g1)(PetscInt, PetscInt, PetscInt, 1345 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1346 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1347 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1348 PetscInt i2, 1349 void (*g2)(PetscInt, PetscInt, PetscInt, 1350 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1351 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1352 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1353 PetscInt i3, 1354 void (*g3)(PetscInt, PetscInt, PetscInt, 1355 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1356 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1357 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1358 { 1359 PetscInt find = f*wf->Nf + g; 1360 PetscErrorCode ierr; 1361 1362 PetscFunctionBegin; 1363 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, i0, (void (*)(void)) g0);CHKERRQ(ierr); 1364 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, i1, (void (*)(void)) g1);CHKERRQ(ierr); 1365 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, i2, (void (*)(void)) g2);CHKERRQ(ierr); 1366 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, i3, (void (*)(void)) g3);CHKERRQ(ierr); 1367 PetscFunctionReturn(0); 1368 } 1369 1370 PetscErrorCode PetscWeakFormGetRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt *n, 1371 void (***r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *)) 1372 { 1373 PetscErrorCode ierr; 1374 1375 PetscFunctionBegin; 1376 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, n, (void (***)(void)) r);CHKERRQ(ierr); 1377 PetscFunctionReturn(0); 1378 } 1379 1380 PetscErrorCode PetscWeakFormSetRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 1381 PetscInt n, 1382 void (**r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *)) 1383 { 1384 PetscErrorCode ierr; 1385 1386 PetscFunctionBegin; 1387 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, n, (void (**)(void)) r);CHKERRQ(ierr); 1388 PetscFunctionReturn(0); 1389 } 1390 1391 PetscErrorCode PetscWeakFormSetIndexRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 1392 PetscInt i, 1393 void (*r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *)) 1394 { 1395 PetscErrorCode ierr; 1396 1397 PetscFunctionBegin; 1398 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, i, (void (*)(void)) r);CHKERRQ(ierr); 1399 PetscFunctionReturn(0); 1400 } 1401 1402 /*@ 1403 PetscWeakFormGetNumFields - Returns the number of fields 1404 1405 Not collective 1406 1407 Input Parameter: 1408 . wf - The PetscWeakForm object 1409 1410 Output Parameter: 1411 . Nf - The number of fields 1412 1413 Level: beginner 1414 1415 .seealso: PetscWeakFormSetNumFields(), PetscWeakFormCreate() 1416 @*/ 1417 PetscErrorCode PetscWeakFormGetNumFields(PetscWeakForm wf, PetscInt *Nf) 1418 { 1419 PetscFunctionBegin; 1420 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 1421 PetscValidPointer(Nf, 2); 1422 *Nf = wf->Nf; 1423 PetscFunctionReturn(0); 1424 } 1425 1426 /*@ 1427 PetscWeakFormSetNumFields - Sets the number of fields 1428 1429 Not collective 1430 1431 Input Parameters: 1432 + wf - The PetscWeakForm object 1433 - Nf - The number of fields 1434 1435 Level: beginner 1436 1437 .seealso: PetscWeakFormGetNumFields(), PetscWeakFormCreate() 1438 @*/ 1439 PetscErrorCode PetscWeakFormSetNumFields(PetscWeakForm wf, PetscInt Nf) 1440 { 1441 PetscFunctionBegin; 1442 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 1443 wf->Nf = Nf; 1444 PetscFunctionReturn(0); 1445 } 1446 1447 /*@ 1448 PetscWeakFormDestroy - Destroys a PetscWeakForm object 1449 1450 Collective on wf 1451 1452 Input Parameter: 1453 . wf - the PetscWeakForm object to destroy 1454 1455 Level: developer 1456 1457 .seealso PetscWeakFormCreate(), PetscWeakFormView() 1458 @*/ 1459 PetscErrorCode PetscWeakFormDestroy(PetscWeakForm *wf) 1460 { 1461 PetscInt f; 1462 PetscErrorCode ierr; 1463 1464 PetscFunctionBegin; 1465 if (!*wf) PetscFunctionReturn(0); 1466 PetscValidHeaderSpecific((*wf), PETSCWEAKFORM_CLASSID, 1); 1467 1468 if (--((PetscObject)(*wf))->refct > 0) {*wf = NULL; PetscFunctionReturn(0);} 1469 ((PetscObject) (*wf))->refct = 0; 1470 ierr = PetscChunkBufferDestroy(&(*wf)->funcs);CHKERRQ(ierr); 1471 for (f = 0; f < PETSC_NUM_WF; ++f) {ierr = PetscHMapFormDestroy(&(*wf)->form[f]);CHKERRQ(ierr);} 1472 ierr = PetscFree((*wf)->form);CHKERRQ(ierr); 1473 ierr = PetscHeaderDestroy(wf);CHKERRQ(ierr); 1474 PetscFunctionReturn(0); 1475 } 1476 1477 static PetscErrorCode PetscWeakFormViewTable_Ascii(PetscWeakForm wf, PetscViewer viewer, PetscBool splitField, const char tableName[], PetscHMapForm map) 1478 { 1479 PetscInt Nf = wf->Nf, Nk, k; 1480 PetscErrorCode ierr; 1481 1482 PetscFunctionBegin; 1483 ierr = PetscHMapFormGetSize(map, &Nk);CHKERRQ(ierr); 1484 if (Nk) { 1485 PetscFormKey *keys; 1486 void (**funcs)(void); 1487 const char **names; 1488 PetscInt *values, *idx1, *idx2, *idx; 1489 PetscBool showPart = PETSC_FALSE, showPointer = PETSC_FALSE; 1490 PetscInt off = 0; 1491 1492 ierr = PetscMalloc6(Nk, &keys, Nk, &names, Nk, &values, Nk, &idx1, Nk, &idx2, Nk, &idx);CHKERRQ(ierr); 1493 ierr = PetscHMapFormGetKeys(map, &off, keys);CHKERRQ(ierr); 1494 /* Sort keys by label name and value */ 1495 { 1496 /* First sort values */ 1497 for (k = 0; k < Nk; ++k) {values[k] = keys[k].value; idx1[k] = k;} 1498 ierr = PetscSortIntWithPermutation(Nk, values, idx1);CHKERRQ(ierr); 1499 /* If the string sort is stable, it will be sorted correctly overall */ 1500 for (k = 0; k < Nk; ++k) { 1501 if (keys[idx1[k]].label) {ierr = PetscObjectGetName((PetscObject) keys[idx1[k]].label, &names[k]);CHKERRQ(ierr);} 1502 else {names[k] = "";} 1503 idx2[k] = k; 1504 } 1505 ierr = PetscSortStrWithPermutation(Nk, names, idx2);CHKERRQ(ierr); 1506 for (k = 0; k < Nk; ++k) { 1507 if (keys[k].label) {ierr = PetscObjectGetName((PetscObject) keys[k].label, &names[k]);CHKERRQ(ierr);} 1508 else {names[k] = "";} 1509 idx[k] = idx1[idx2[k]]; 1510 } 1511 } 1512 ierr = PetscViewerASCIIPrintf(viewer, "%s\n", tableName);CHKERRQ(ierr); 1513 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 1514 for (k = 0; k < Nk; ++k) { 1515 if (keys[k].part != 0) showPart = PETSC_TRUE; 1516 } 1517 for (k = 0; k < Nk; ++k) { 1518 const PetscInt i = idx[k]; 1519 PetscInt n, f; 1520 1521 if (keys[i].label) { 1522 if (showPointer) {ierr = PetscViewerASCIIPrintf(viewer, "(%s:%p, %D) ", names[i], keys[i].label, keys[i].value);CHKERRQ(ierr);} 1523 else {ierr = PetscViewerASCIIPrintf(viewer, "(%s, %D) ", names[i], keys[i].value);CHKERRQ(ierr);} 1524 } else {ierr = PetscViewerASCIIPrintf(viewer, "");CHKERRQ(ierr);} 1525 ierr = PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);CHKERRQ(ierr); 1526 if (splitField) {ierr = PetscViewerASCIIPrintf(viewer, "(%D, %D) ", keys[i].field/Nf, keys[i].field%Nf);CHKERRQ(ierr);} 1527 else {ierr = PetscViewerASCIIPrintf(viewer, "(%D) ", keys[i].field);CHKERRQ(ierr);} 1528 if (showPart) {ierr = PetscViewerASCIIPrintf(viewer, "(%D) ", keys[i].part);CHKERRQ(ierr);} 1529 ierr = PetscWeakFormGetFunction_Private(wf, map, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &n, &funcs);CHKERRQ(ierr); 1530 for (f = 0; f < n; ++f) { 1531 char *fname; 1532 size_t len, l; 1533 1534 if (f > 0) {ierr = PetscViewerASCIIPrintf(viewer, ", ");CHKERRQ(ierr);} 1535 ierr = PetscDLAddr(funcs[f], &fname);CHKERRQ(ierr); 1536 if (fname) { 1537 /* Eliminate argument types */ 1538 ierr = PetscStrlen(fname, &len);CHKERRQ(ierr); 1539 for (l = 0; l < len; ++l) if (fname[l] == '(') {fname[l] = '\0'; break;} 1540 ierr = PetscViewerASCIIPrintf(viewer, "%s", fname);CHKERRQ(ierr); 1541 } else if (showPointer) { 1542 ierr = PetscViewerASCIIPrintf(viewer, "%p", funcs[f]);CHKERRQ(ierr); 1543 } 1544 ierr = PetscFree(fname);CHKERRQ(ierr); 1545 } 1546 ierr = PetscViewerASCIIPrintf(viewer, "\n");CHKERRQ(ierr); 1547 ierr = PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);CHKERRQ(ierr); 1548 } 1549 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 1550 ierr = PetscFree6(keys, names, values, idx1, idx2, idx);CHKERRQ(ierr); 1551 } 1552 PetscFunctionReturn(0); 1553 } 1554 1555 static PetscErrorCode PetscWeakFormView_Ascii(PetscWeakForm wf, PetscViewer viewer) 1556 { 1557 PetscViewerFormat format; 1558 PetscInt f; 1559 PetscErrorCode ierr; 1560 1561 PetscFunctionBegin; 1562 ierr = PetscViewerGetFormat(viewer, &format);CHKERRQ(ierr); 1563 ierr = PetscViewerASCIIPrintf(viewer, "Weak Form System with %d fields\n", wf->Nf);CHKERRQ(ierr); 1564 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 1565 for (f = 0; f < PETSC_NUM_WF; ++f) { 1566 ierr = PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE, PetscWeakFormKinds[f], wf->form[f]);CHKERRQ(ierr); 1567 } 1568 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 1569 PetscFunctionReturn(0); 1570 } 1571 1572 /*@C 1573 PetscWeakFormView - Views a PetscWeakForm 1574 1575 Collective on wf 1576 1577 Input Parameters: 1578 + wf - the PetscWeakForm object to view 1579 - v - the viewer 1580 1581 Level: developer 1582 1583 .seealso PetscWeakFormDestroy(), PetscWeakFormCreate() 1584 @*/ 1585 PetscErrorCode PetscWeakFormView(PetscWeakForm wf, PetscViewer v) 1586 { 1587 PetscBool iascii; 1588 PetscErrorCode ierr; 1589 1590 PetscFunctionBegin; 1591 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 1592 if (!v) {ierr = PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject) wf), &v);CHKERRQ(ierr);} 1593 else {PetscValidHeaderSpecific(v, PETSC_VIEWER_CLASSID, 2);} 1594 ierr = PetscObjectTypeCompare((PetscObject) v, PETSCVIEWERASCII, &iascii);CHKERRQ(ierr); 1595 if (iascii) {ierr = PetscWeakFormView_Ascii(wf, v);CHKERRQ(ierr);} 1596 if (wf->ops->view) {ierr = (*wf->ops->view)(wf, v);CHKERRQ(ierr);} 1597 PetscFunctionReturn(0); 1598 } 1599 1600 /*@ 1601 PetscWeakFormCreate - Creates an empty PetscWeakForm object. 1602 1603 Collective 1604 1605 Input Parameter: 1606 . comm - The communicator for the PetscWeakForm object 1607 1608 Output Parameter: 1609 . wf - The PetscWeakForm object 1610 1611 Level: beginner 1612 1613 .seealso: PetscDS, PetscWeakFormDestroy() 1614 @*/ 1615 PetscErrorCode PetscWeakFormCreate(MPI_Comm comm, PetscWeakForm *wf) 1616 { 1617 PetscWeakForm p; 1618 PetscInt f; 1619 PetscErrorCode ierr; 1620 1621 PetscFunctionBegin; 1622 PetscValidPointer(wf, 2); 1623 *wf = NULL; 1624 ierr = PetscDSInitializePackage();CHKERRQ(ierr); 1625 1626 ierr = PetscHeaderCreate(p, PETSCWEAKFORM_CLASSID, "PetscWeakForm", "Weak Form System", "PetscWeakForm", comm, PetscWeakFormDestroy, PetscWeakFormView);CHKERRQ(ierr); 1627 1628 p->Nf = 0; 1629 ierr = PetscChunkBufferCreate(sizeof(&PetscWeakFormCreate), 2, &p->funcs);CHKERRQ(ierr); 1630 ierr = PetscMalloc1(PETSC_NUM_WF, &p->form);CHKERRQ(ierr); 1631 for (f = 0; f < PETSC_NUM_WF; ++f) {ierr = PetscHMapFormCreate(&p->form[f]);CHKERRQ(ierr);} 1632 *wf = p; 1633 PetscFunctionReturn(0); 1634 } 1635