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: 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 PetscErrorCode PetscWeakFormClearIndex(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscWeakFormKind kind, PetscInt ind) 344 { 345 PetscErrorCode ierr; 346 347 PetscFunctionBegin; 348 ierr = PetscWeakFormClearIndexFunction_Private(wf, wf->form[kind], label, val, f, part, ind);CHKERRQ(ierr); 349 PetscFunctionReturn(0); 350 } 351 352 PetscErrorCode PetscWeakFormGetObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt *n, 353 void (***obj)(PetscInt, PetscInt, PetscInt, 354 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 355 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 356 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 357 { 358 PetscErrorCode ierr; 359 360 PetscFunctionBegin; 361 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, n, (void (***)(void)) obj);CHKERRQ(ierr); 362 PetscFunctionReturn(0); 363 } 364 365 PetscErrorCode PetscWeakFormSetObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt n, 366 void (**obj)(PetscInt, PetscInt, PetscInt, 367 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 368 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 369 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 370 { 371 PetscErrorCode ierr; 372 373 PetscFunctionBegin; 374 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, n, (void (**)(void)) obj);CHKERRQ(ierr); 375 PetscFunctionReturn(0); 376 } 377 378 PetscErrorCode PetscWeakFormAddObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 379 void (*obj)(PetscInt, PetscInt, PetscInt, 380 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 381 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 382 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 383 { 384 PetscErrorCode ierr; 385 386 PetscFunctionBegin; 387 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, (void (*)(void)) obj);CHKERRQ(ierr); 388 PetscFunctionReturn(0); 389 } 390 391 PetscErrorCode PetscWeakFormGetIndexObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt ind, 392 void (**obj)(PetscInt, PetscInt, PetscInt, 393 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 394 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 395 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 396 { 397 PetscErrorCode ierr; 398 399 PetscFunctionBegin; 400 ierr = PetscWeakFormGetIndexFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, ind, (void (**)(void)) obj);CHKERRQ(ierr); 401 PetscFunctionReturn(0); 402 } 403 404 PetscErrorCode PetscWeakFormSetIndexObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt ind, 405 void (*obj)(PetscInt, PetscInt, PetscInt, 406 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 407 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 408 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 409 { 410 PetscErrorCode ierr; 411 412 PetscFunctionBegin; 413 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, ind, (void (*)(void)) obj);CHKERRQ(ierr); 414 PetscFunctionReturn(0); 415 } 416 417 PetscErrorCode PetscWeakFormGetResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 418 PetscInt *n0, 419 void (***f0)(PetscInt, PetscInt, PetscInt, 420 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 421 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 422 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 423 PetscInt *n1, 424 void (***f1)(PetscInt, PetscInt, PetscInt, 425 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 426 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 427 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 428 { 429 PetscErrorCode ierr; 430 431 PetscFunctionBegin; 432 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, n0, (void (***)(void)) f0);CHKERRQ(ierr); 433 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, n1, (void (***)(void)) f1);CHKERRQ(ierr); 434 PetscFunctionReturn(0); 435 } 436 437 PetscErrorCode PetscWeakFormAddResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 438 void (*f0)(PetscInt, PetscInt, PetscInt, 439 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 440 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 441 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 442 void (*f1)(PetscInt, PetscInt, PetscInt, 443 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 444 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 445 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 446 { 447 PetscErrorCode ierr; 448 449 PetscFunctionBegin; 450 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, (void (*)(void)) f0);CHKERRQ(ierr); 451 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, (void (*)(void)) f1);CHKERRQ(ierr); 452 PetscFunctionReturn(0); 453 } 454 455 PetscErrorCode PetscWeakFormSetResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 456 PetscInt n0, 457 void (**f0)(PetscInt, PetscInt, PetscInt, 458 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 459 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 460 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 461 PetscInt n1, 462 void (**f1)(PetscInt, PetscInt, PetscInt, 463 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 464 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 465 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 466 { 467 PetscErrorCode ierr; 468 469 PetscFunctionBegin; 470 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, n0, (void (**)(void)) f0);CHKERRQ(ierr); 471 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, n1, (void (**)(void)) f1);CHKERRQ(ierr); 472 PetscFunctionReturn(0); 473 } 474 475 PetscErrorCode PetscWeakFormSetIndexResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 476 PetscInt i0, 477 void (*f0)(PetscInt, PetscInt, PetscInt, 478 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 479 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 480 PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 481 PetscInt i1, 482 void (*f1)(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_F0], label, val, f, part, i0, (void (*)(void)) f0);CHKERRQ(ierr); 491 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, i1, (void (*)(void)) f1);CHKERRQ(ierr); 492 PetscFunctionReturn(0); 493 } 494 495 PetscErrorCode PetscWeakFormGetBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 496 PetscInt *n0, 497 void (***f0)(PetscInt, PetscInt, PetscInt, 498 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 499 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 500 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 501 PetscInt *n1, 502 void (***f1)(PetscInt, PetscInt, PetscInt, 503 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 504 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 505 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 506 { 507 PetscErrorCode ierr; 508 509 PetscFunctionBegin; 510 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, n0, (void (***)(void)) f0);CHKERRQ(ierr); 511 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, n1, (void (***)(void)) f1);CHKERRQ(ierr); 512 PetscFunctionReturn(0); 513 } 514 515 PetscErrorCode PetscWeakFormAddBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 516 void (*f0)(PetscInt, PetscInt, PetscInt, 517 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 518 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 519 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 520 void (*f1)(PetscInt, PetscInt, PetscInt, 521 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 522 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 523 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 524 { 525 PetscErrorCode ierr; 526 527 PetscFunctionBegin; 528 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, (void (*)(void)) f0);CHKERRQ(ierr); 529 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, (void (*)(void)) f1);CHKERRQ(ierr); 530 PetscFunctionReturn(0); 531 } 532 533 PetscErrorCode PetscWeakFormSetBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 534 PetscInt n0, 535 void (**f0)(PetscInt, PetscInt, PetscInt, 536 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 537 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 538 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 539 PetscInt n1, 540 void (**f1)(PetscInt, PetscInt, PetscInt, 541 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 542 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 543 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 544 { 545 PetscErrorCode ierr; 546 547 PetscFunctionBegin; 548 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, n0, (void (**)(void)) f0);CHKERRQ(ierr); 549 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, n1, (void (**)(void)) f1);CHKERRQ(ierr); 550 PetscFunctionReturn(0); 551 } 552 553 PetscErrorCode PetscWeakFormSetIndexBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 554 PetscInt i0, 555 void (*f0)(PetscInt, PetscInt, PetscInt, 556 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 557 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 558 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 559 PetscInt i1, 560 void (*f1)(PetscInt, PetscInt, PetscInt, 561 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 562 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 563 PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 564 { 565 PetscErrorCode ierr; 566 567 PetscFunctionBegin; 568 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, i0, (void (*)(void)) f0);CHKERRQ(ierr); 569 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, i1, (void (*)(void)) f1);CHKERRQ(ierr); 570 PetscFunctionReturn(0); 571 } 572 573 PetscErrorCode PetscWeakFormHasJacobian(PetscWeakForm wf, PetscBool *hasJac) 574 { 575 PetscInt n0, n1, n2, n3; 576 PetscErrorCode ierr; 577 578 PetscFunctionBegin; 579 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 580 PetscValidBoolPointer(hasJac, 2); 581 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_G0], &n0);CHKERRQ(ierr); 582 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_G1], &n1);CHKERRQ(ierr); 583 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_G2], &n2);CHKERRQ(ierr); 584 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_G3], &n3);CHKERRQ(ierr); 585 *hasJac = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE; 586 PetscFunctionReturn(0); 587 } 588 589 PetscErrorCode PetscWeakFormGetJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 590 PetscInt *n0, 591 void (***g0)(PetscInt, PetscInt, PetscInt, 592 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 593 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 594 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 595 PetscInt *n1, 596 void (***g1)(PetscInt, PetscInt, PetscInt, 597 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 598 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 599 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 600 PetscInt *n2, 601 void (***g2)(PetscInt, PetscInt, PetscInt, 602 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 603 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 604 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 605 PetscInt *n3, 606 void (***g3)(PetscInt, PetscInt, PetscInt, 607 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 608 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 609 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 610 { 611 PetscInt find = f*wf->Nf + g; 612 PetscErrorCode ierr; 613 614 PetscFunctionBegin; 615 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, n0, (void (***)(void)) g0);CHKERRQ(ierr); 616 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, n1, (void (***)(void)) g1);CHKERRQ(ierr); 617 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, n2, (void (***)(void)) g2);CHKERRQ(ierr); 618 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, n3, (void (***)(void)) g3);CHKERRQ(ierr); 619 PetscFunctionReturn(0); 620 } 621 622 PetscErrorCode PetscWeakFormAddJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 623 void (*g0)(PetscInt, PetscInt, PetscInt, 624 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 625 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 626 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 627 void (*g1)(PetscInt, PetscInt, PetscInt, 628 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 629 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 630 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 631 void (*g2)(PetscInt, PetscInt, PetscInt, 632 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 633 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 634 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 635 void (*g3)(PetscInt, PetscInt, PetscInt, 636 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 637 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 638 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 639 { 640 PetscInt find = f*wf->Nf + g; 641 PetscErrorCode ierr; 642 643 PetscFunctionBegin; 644 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, (void (*)(void)) g0);CHKERRQ(ierr); 645 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, (void (*)(void)) g1);CHKERRQ(ierr); 646 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, (void (*)(void)) g2);CHKERRQ(ierr); 647 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, (void (*)(void)) g3);CHKERRQ(ierr); 648 PetscFunctionReturn(0); 649 } 650 651 PetscErrorCode PetscWeakFormSetJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 652 PetscInt n0, 653 void (**g0)(PetscInt, PetscInt, PetscInt, 654 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 655 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 656 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 657 PetscInt n1, 658 void (**g1)(PetscInt, PetscInt, PetscInt, 659 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 660 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 661 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 662 PetscInt n2, 663 void (**g2)(PetscInt, PetscInt, PetscInt, 664 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 665 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 666 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 667 PetscInt n3, 668 void (**g3)(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 { 673 PetscInt find = f*wf->Nf + g; 674 PetscErrorCode ierr; 675 676 PetscFunctionBegin; 677 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, n0, (void (**)(void)) g0);CHKERRQ(ierr); 678 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, n1, (void (**)(void)) g1);CHKERRQ(ierr); 679 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, n2, (void (**)(void)) g2);CHKERRQ(ierr); 680 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, n3, (void (**)(void)) g3);CHKERRQ(ierr); 681 PetscFunctionReturn(0); 682 } 683 684 PetscErrorCode PetscWeakFormSetIndexJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 685 PetscInt i0, 686 void (*g0)(PetscInt, PetscInt, PetscInt, 687 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 688 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 689 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 690 PetscInt i1, 691 void (*g1)(PetscInt, PetscInt, PetscInt, 692 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 693 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 694 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 695 PetscInt i2, 696 void (*g2)(PetscInt, PetscInt, PetscInt, 697 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 698 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 699 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 700 PetscInt i3, 701 void (*g3)(PetscInt, PetscInt, PetscInt, 702 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 703 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 704 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 705 { 706 PetscInt find = f*wf->Nf + g; 707 PetscErrorCode ierr; 708 709 PetscFunctionBegin; 710 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, i0, (void (*)(void)) g0);CHKERRQ(ierr); 711 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, i1, (void (*)(void)) g1);CHKERRQ(ierr); 712 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, i2, (void (*)(void)) g2);CHKERRQ(ierr); 713 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, i3, (void (*)(void)) g3);CHKERRQ(ierr); 714 PetscFunctionReturn(0); 715 } 716 717 PetscErrorCode PetscWeakFormHasJacobianPreconditioner(PetscWeakForm wf, PetscBool *hasJacPre) 718 { 719 PetscInt n0, n1, n2, n3; 720 PetscErrorCode ierr; 721 722 PetscFunctionBegin; 723 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 724 PetscValidBoolPointer(hasJacPre, 2); 725 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GP0], &n0);CHKERRQ(ierr); 726 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GP1], &n1);CHKERRQ(ierr); 727 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GP2], &n2);CHKERRQ(ierr); 728 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GP3], &n3);CHKERRQ(ierr); 729 *hasJacPre = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE; 730 PetscFunctionReturn(0); 731 } 732 733 PetscErrorCode PetscWeakFormGetJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 734 PetscInt *n0, 735 void (***g0)(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 *n1, 740 void (***g1)(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 *n2, 745 void (***g2)(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 PetscInt *n3, 750 void (***g3)(PetscInt, PetscInt, PetscInt, 751 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 752 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 753 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 754 { 755 PetscInt find = f*wf->Nf + g; 756 PetscErrorCode ierr; 757 758 PetscFunctionBegin; 759 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, n0, (void (***)(void)) g0);CHKERRQ(ierr); 760 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, n1, (void (***)(void)) g1);CHKERRQ(ierr); 761 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, n2, (void (***)(void)) g2);CHKERRQ(ierr); 762 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, n3, (void (***)(void)) g3);CHKERRQ(ierr); 763 PetscFunctionReturn(0); 764 } 765 766 PetscErrorCode PetscWeakFormAddJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 767 void (*g0)(PetscInt, PetscInt, PetscInt, 768 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 769 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 770 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 771 void (*g1)(PetscInt, PetscInt, PetscInt, 772 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 773 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 774 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 775 void (*g2)(PetscInt, PetscInt, PetscInt, 776 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 777 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 778 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 779 void (*g3)(PetscInt, PetscInt, PetscInt, 780 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 781 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 782 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 783 { 784 PetscInt find = f*wf->Nf + g; 785 PetscErrorCode ierr; 786 787 PetscFunctionBegin; 788 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, (void (*)(void)) g0);CHKERRQ(ierr); 789 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, (void (*)(void)) g1);CHKERRQ(ierr); 790 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, (void (*)(void)) g2);CHKERRQ(ierr); 791 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, (void (*)(void)) g3);CHKERRQ(ierr); 792 PetscFunctionReturn(0); 793 } 794 795 PetscErrorCode PetscWeakFormSetJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 796 PetscInt n0, 797 void (**g0)(PetscInt, PetscInt, PetscInt, 798 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 799 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 800 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 801 PetscInt n1, 802 void (**g1)(PetscInt, PetscInt, PetscInt, 803 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 804 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 805 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 806 PetscInt n2, 807 void (**g2)(PetscInt, PetscInt, PetscInt, 808 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 809 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 810 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 811 PetscInt n3, 812 void (**g3)(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 { 817 PetscInt find = f*wf->Nf + g; 818 PetscErrorCode ierr; 819 820 PetscFunctionBegin; 821 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, n0, (void (**)(void)) g0);CHKERRQ(ierr); 822 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, n1, (void (**)(void)) g1);CHKERRQ(ierr); 823 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, n2, (void (**)(void)) g2);CHKERRQ(ierr); 824 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, n3, (void (**)(void)) g3);CHKERRQ(ierr); 825 PetscFunctionReturn(0); 826 } 827 828 PetscErrorCode PetscWeakFormSetIndexJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 829 PetscInt i0, 830 void (*g0)(PetscInt, PetscInt, PetscInt, 831 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 832 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 833 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 834 PetscInt i1, 835 void (*g1)(PetscInt, PetscInt, PetscInt, 836 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 837 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 838 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 839 PetscInt i2, 840 void (*g2)(PetscInt, PetscInt, PetscInt, 841 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 842 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 843 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 844 PetscInt i3, 845 void (*g3)(PetscInt, PetscInt, PetscInt, 846 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 847 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 848 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 849 { 850 PetscInt find = f*wf->Nf + g; 851 PetscErrorCode ierr; 852 853 PetscFunctionBegin; 854 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, i0, (void (*)(void)) g0);CHKERRQ(ierr); 855 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, i1, (void (*)(void)) g1);CHKERRQ(ierr); 856 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, i2, (void (*)(void)) g2);CHKERRQ(ierr); 857 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, i3, (void (*)(void)) g3);CHKERRQ(ierr); 858 PetscFunctionReturn(0); 859 } 860 861 PetscErrorCode PetscWeakFormHasBdJacobian(PetscWeakForm wf, PetscBool *hasJac) 862 { 863 PetscInt n0, n1, n2, n3; 864 PetscErrorCode ierr; 865 866 PetscFunctionBegin; 867 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 868 PetscValidBoolPointer(hasJac, 2); 869 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDG0], &n0);CHKERRQ(ierr); 870 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDG1], &n1);CHKERRQ(ierr); 871 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDG2], &n2);CHKERRQ(ierr); 872 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDG3], &n3);CHKERRQ(ierr); 873 *hasJac = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE; 874 PetscFunctionReturn(0); 875 } 876 877 PetscErrorCode PetscWeakFormGetBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 878 PetscInt *n0, 879 void (***g0)(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[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 883 PetscInt *n1, 884 void (***g1)(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[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 888 PetscInt *n2, 889 void (***g2)(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[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 893 PetscInt *n3, 894 void (***g3)(PetscInt, PetscInt, PetscInt, 895 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 896 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 897 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 898 { 899 PetscInt find = f*wf->Nf + g; 900 PetscErrorCode ierr; 901 902 PetscFunctionBegin; 903 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, n0, (void (***)(void)) g0);CHKERRQ(ierr); 904 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, n1, (void (***)(void)) g1);CHKERRQ(ierr); 905 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, n2, (void (***)(void)) g2);CHKERRQ(ierr); 906 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, n3, (void (***)(void)) g3);CHKERRQ(ierr); 907 PetscFunctionReturn(0); 908 } 909 910 PetscErrorCode PetscWeakFormAddBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 911 void (*g0)(PetscInt, PetscInt, PetscInt, 912 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 913 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 914 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 915 void (*g1)(PetscInt, PetscInt, PetscInt, 916 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 917 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 918 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 919 void (*g2)(PetscInt, PetscInt, PetscInt, 920 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 921 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 922 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 923 void (*g3)(PetscInt, PetscInt, PetscInt, 924 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 925 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 926 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 927 { 928 PetscInt find = f*wf->Nf + g; 929 PetscErrorCode ierr; 930 931 PetscFunctionBegin; 932 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, (void (*)(void)) g0);CHKERRQ(ierr); 933 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, (void (*)(void)) g1);CHKERRQ(ierr); 934 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, (void (*)(void)) g2);CHKERRQ(ierr); 935 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, (void (*)(void)) g3);CHKERRQ(ierr); 936 PetscFunctionReturn(0); 937 } 938 939 PetscErrorCode PetscWeakFormSetBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 940 PetscInt n0, 941 void (**g0)(PetscInt, PetscInt, PetscInt, 942 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 943 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 944 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 945 PetscInt n1, 946 void (**g1)(PetscInt, PetscInt, PetscInt, 947 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 948 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 949 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 950 PetscInt n2, 951 void (**g2)(PetscInt, PetscInt, PetscInt, 952 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 953 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 954 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 955 PetscInt n3, 956 void (**g3)(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 { 961 PetscInt find = f*wf->Nf + g; 962 PetscErrorCode ierr; 963 964 PetscFunctionBegin; 965 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, n0, (void (**)(void)) g0);CHKERRQ(ierr); 966 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, n1, (void (**)(void)) g1);CHKERRQ(ierr); 967 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, n2, (void (**)(void)) g2);CHKERRQ(ierr); 968 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, n3, (void (**)(void)) g3);CHKERRQ(ierr); 969 PetscFunctionReturn(0); 970 } 971 972 PetscErrorCode PetscWeakFormSetIndexBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 973 PetscInt i0, 974 void (*g0)(PetscInt, PetscInt, PetscInt, 975 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 976 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 977 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 978 PetscInt i1, 979 void (*g1)(PetscInt, PetscInt, PetscInt, 980 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 981 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 982 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 983 PetscInt i2, 984 void (*g2)(PetscInt, PetscInt, PetscInt, 985 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 986 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 987 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 988 PetscInt i3, 989 void (*g3)(PetscInt, PetscInt, PetscInt, 990 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 991 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 992 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 993 { 994 PetscInt find = f*wf->Nf + g; 995 PetscErrorCode ierr; 996 997 PetscFunctionBegin; 998 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, i0, (void (*)(void)) g0);CHKERRQ(ierr); 999 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, i1, (void (*)(void)) g1);CHKERRQ(ierr); 1000 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, i2, (void (*)(void)) g2);CHKERRQ(ierr); 1001 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, i3, (void (*)(void)) g3);CHKERRQ(ierr); 1002 PetscFunctionReturn(0); 1003 } 1004 1005 PetscErrorCode PetscWeakFormHasBdJacobianPreconditioner(PetscWeakForm wf, PetscBool *hasJacPre) 1006 { 1007 PetscInt n0, n1, n2, n3; 1008 PetscErrorCode ierr; 1009 1010 PetscFunctionBegin; 1011 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 1012 PetscValidBoolPointer(hasJacPre, 2); 1013 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP0], &n0);CHKERRQ(ierr); 1014 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP1], &n1);CHKERRQ(ierr); 1015 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP2], &n2);CHKERRQ(ierr); 1016 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP3], &n3);CHKERRQ(ierr); 1017 *hasJacPre = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE; 1018 PetscFunctionReturn(0); 1019 } 1020 1021 PetscErrorCode PetscWeakFormGetBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1022 PetscInt *n0, 1023 void (***g0)(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 *n1, 1028 void (***g1)(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 *n2, 1033 void (***g2)(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 PetscInt *n3, 1038 void (***g3)(PetscInt, PetscInt, PetscInt, 1039 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1040 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1041 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1042 { 1043 PetscInt find = f*wf->Nf + g; 1044 PetscErrorCode ierr; 1045 1046 PetscFunctionBegin; 1047 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, n0, (void (***)(void)) g0);CHKERRQ(ierr); 1048 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, n1, (void (***)(void)) g1);CHKERRQ(ierr); 1049 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, n2, (void (***)(void)) g2);CHKERRQ(ierr); 1050 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, n3, (void (***)(void)) g3);CHKERRQ(ierr); 1051 PetscFunctionReturn(0); 1052 } 1053 1054 PetscErrorCode PetscWeakFormAddBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1055 void (*g0)(PetscInt, PetscInt, PetscInt, 1056 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1057 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1058 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1059 void (*g1)(PetscInt, PetscInt, PetscInt, 1060 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1061 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1062 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1063 void (*g2)(PetscInt, PetscInt, PetscInt, 1064 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1065 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1066 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1067 void (*g3)(PetscInt, PetscInt, PetscInt, 1068 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1069 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1070 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1071 { 1072 PetscInt find = f*wf->Nf + g; 1073 PetscErrorCode ierr; 1074 1075 PetscFunctionBegin; 1076 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, (void (*)(void)) g0);CHKERRQ(ierr); 1077 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, (void (*)(void)) g1);CHKERRQ(ierr); 1078 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, (void (*)(void)) g2);CHKERRQ(ierr); 1079 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, (void (*)(void)) g3);CHKERRQ(ierr); 1080 PetscFunctionReturn(0); 1081 } 1082 1083 PetscErrorCode PetscWeakFormSetBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1084 PetscInt n0, 1085 void (**g0)(PetscInt, PetscInt, PetscInt, 1086 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1087 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1088 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1089 PetscInt n1, 1090 void (**g1)(PetscInt, PetscInt, PetscInt, 1091 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1092 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1093 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1094 PetscInt n2, 1095 void (**g2)(PetscInt, PetscInt, PetscInt, 1096 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1097 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1098 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1099 PetscInt n3, 1100 void (**g3)(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 { 1105 PetscInt find = f*wf->Nf + g; 1106 PetscErrorCode ierr; 1107 1108 PetscFunctionBegin; 1109 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, n0, (void (**)(void)) g0);CHKERRQ(ierr); 1110 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, n1, (void (**)(void)) g1);CHKERRQ(ierr); 1111 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, n2, (void (**)(void)) g2);CHKERRQ(ierr); 1112 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, n3, (void (**)(void)) g3);CHKERRQ(ierr); 1113 PetscFunctionReturn(0); 1114 } 1115 1116 PetscErrorCode PetscWeakFormSetIndexBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1117 PetscInt i0, 1118 void (*g0)(PetscInt, PetscInt, PetscInt, 1119 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1120 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1121 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1122 PetscInt i1, 1123 void (*g1)(PetscInt, PetscInt, PetscInt, 1124 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1125 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1126 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1127 PetscInt i2, 1128 void (*g2)(PetscInt, PetscInt, PetscInt, 1129 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1130 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1131 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1132 PetscInt i3, 1133 void (*g3)(PetscInt, PetscInt, PetscInt, 1134 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1135 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1136 PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1137 { 1138 PetscInt find = f*wf->Nf + g; 1139 PetscErrorCode ierr; 1140 1141 PetscFunctionBegin; 1142 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, i0, (void (*)(void)) g0);CHKERRQ(ierr); 1143 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, i1, (void (*)(void)) g1);CHKERRQ(ierr); 1144 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, i2, (void (*)(void)) g2);CHKERRQ(ierr); 1145 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, i3, (void (*)(void)) g3);CHKERRQ(ierr); 1146 PetscFunctionReturn(0); 1147 } 1148 1149 PetscErrorCode PetscWeakFormHasDynamicJacobian(PetscWeakForm wf, PetscBool *hasDynJac) 1150 { 1151 PetscInt n0, n1, n2, n3; 1152 PetscErrorCode ierr; 1153 1154 PetscFunctionBegin; 1155 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 1156 PetscValidBoolPointer(hasDynJac, 2); 1157 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GT0], &n0);CHKERRQ(ierr); 1158 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GT1], &n1);CHKERRQ(ierr); 1159 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GT2], &n2);CHKERRQ(ierr); 1160 ierr = PetscHMapFormGetSize(wf->form[PETSC_WF_GT3], &n3);CHKERRQ(ierr); 1161 *hasDynJac = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE; 1162 PetscFunctionReturn(0); 1163 } 1164 1165 PetscErrorCode PetscWeakFormGetDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1166 PetscInt *n0, 1167 void (***g0)(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[], PetscInt, const PetscScalar[], PetscScalar[]), 1171 PetscInt *n1, 1172 void (***g1)(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[], PetscInt, const PetscScalar[], PetscScalar[]), 1176 PetscInt *n2, 1177 void (***g2)(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[], PetscInt, const PetscScalar[], PetscScalar[]), 1181 PetscInt *n3, 1182 void (***g3)(PetscInt, PetscInt, PetscInt, 1183 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1184 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1185 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1186 { 1187 PetscInt find = f*wf->Nf + g; 1188 PetscErrorCode ierr; 1189 1190 PetscFunctionBegin; 1191 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, n0, (void (***)(void)) g0);CHKERRQ(ierr); 1192 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, n1, (void (***)(void)) g1);CHKERRQ(ierr); 1193 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, n2, (void (***)(void)) g2);CHKERRQ(ierr); 1194 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, n3, (void (***)(void)) g3);CHKERRQ(ierr); 1195 PetscFunctionReturn(0); 1196 } 1197 1198 PetscErrorCode PetscWeakFormAddDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1199 void (*g0)(PetscInt, PetscInt, PetscInt, 1200 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1201 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1202 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1203 void (*g1)(PetscInt, PetscInt, PetscInt, 1204 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1205 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1206 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1207 void (*g2)(PetscInt, PetscInt, PetscInt, 1208 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1209 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1210 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1211 void (*g3)(PetscInt, PetscInt, PetscInt, 1212 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1213 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1214 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1215 { 1216 PetscInt find = f*wf->Nf + g; 1217 PetscErrorCode ierr; 1218 1219 PetscFunctionBegin; 1220 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, (void (*)(void)) g0);CHKERRQ(ierr); 1221 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, (void (*)(void)) g1);CHKERRQ(ierr); 1222 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, (void (*)(void)) g2);CHKERRQ(ierr); 1223 ierr = PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, (void (*)(void)) g3);CHKERRQ(ierr); 1224 PetscFunctionReturn(0); 1225 } 1226 1227 PetscErrorCode PetscWeakFormSetDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1228 PetscInt n0, 1229 void (**g0)(PetscInt, PetscInt, PetscInt, 1230 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1231 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1232 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1233 PetscInt n1, 1234 void (**g1)(PetscInt, PetscInt, PetscInt, 1235 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1236 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1237 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1238 PetscInt n2, 1239 void (**g2)(PetscInt, PetscInt, PetscInt, 1240 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1241 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1242 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1243 PetscInt n3, 1244 void (**g3)(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 { 1249 PetscInt find = f*wf->Nf + g; 1250 PetscErrorCode ierr; 1251 1252 PetscFunctionBegin; 1253 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, n0, (void (**)(void)) g0);CHKERRQ(ierr); 1254 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, n1, (void (**)(void)) g1);CHKERRQ(ierr); 1255 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, n2, (void (**)(void)) g2);CHKERRQ(ierr); 1256 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, n3, (void (**)(void)) g3);CHKERRQ(ierr); 1257 PetscFunctionReturn(0); 1258 } 1259 1260 PetscErrorCode PetscWeakFormSetIndexDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part, 1261 PetscInt i0, 1262 void (*g0)(PetscInt, PetscInt, PetscInt, 1263 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1264 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1265 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1266 PetscInt i1, 1267 void (*g1)(PetscInt, PetscInt, PetscInt, 1268 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1269 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1270 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1271 PetscInt i2, 1272 void (*g2)(PetscInt, PetscInt, PetscInt, 1273 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1274 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1275 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), 1276 PetscInt i3, 1277 void (*g3)(PetscInt, PetscInt, PetscInt, 1278 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1279 const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], 1280 PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) 1281 { 1282 PetscInt find = f*wf->Nf + g; 1283 PetscErrorCode ierr; 1284 1285 PetscFunctionBegin; 1286 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, i0, (void (*)(void)) g0);CHKERRQ(ierr); 1287 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, i1, (void (*)(void)) g1);CHKERRQ(ierr); 1288 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, i2, (void (*)(void)) g2);CHKERRQ(ierr); 1289 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, i3, (void (*)(void)) g3);CHKERRQ(ierr); 1290 PetscFunctionReturn(0); 1291 } 1292 1293 PetscErrorCode PetscWeakFormGetRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt *n, 1294 void (***r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *)) 1295 { 1296 PetscErrorCode ierr; 1297 1298 PetscFunctionBegin; 1299 ierr = PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, n, (void (***)(void)) r);CHKERRQ(ierr); 1300 PetscFunctionReturn(0); 1301 } 1302 1303 PetscErrorCode PetscWeakFormSetRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 1304 PetscInt n, 1305 void (**r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *)) 1306 { 1307 PetscErrorCode ierr; 1308 1309 PetscFunctionBegin; 1310 ierr = PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, n, (void (**)(void)) r);CHKERRQ(ierr); 1311 PetscFunctionReturn(0); 1312 } 1313 1314 PetscErrorCode PetscWeakFormSetIndexRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, 1315 PetscInt i, 1316 void (*r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *)) 1317 { 1318 PetscErrorCode ierr; 1319 1320 PetscFunctionBegin; 1321 ierr = PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, i, (void (*)(void)) r);CHKERRQ(ierr); 1322 PetscFunctionReturn(0); 1323 } 1324 1325 /*@ 1326 PetscWeakFormGetNumFields - Returns the number of fields 1327 1328 Not collective 1329 1330 Input Parameter: 1331 . wf - The PetscWeakForm object 1332 1333 Output Parameter: 1334 . Nf - The nubmer of fields 1335 1336 Level: beginner 1337 1338 .seealso: PetscWeakFormSetNumFields(), PetscWeakFormCreate() 1339 @*/ 1340 PetscErrorCode PetscWeakFormGetNumFields(PetscWeakForm wf, PetscInt *Nf) 1341 { 1342 PetscFunctionBegin; 1343 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 1344 PetscValidPointer(Nf, 2); 1345 *Nf = wf->Nf; 1346 PetscFunctionReturn(0); 1347 } 1348 1349 /*@ 1350 PetscWeakFormSetNumFields - Sets the number of fields 1351 1352 Not collective 1353 1354 Input Parameters: 1355 + wf - The PetscWeakForm object 1356 - Nf - The number of fields 1357 1358 Level: beginner 1359 1360 .seealso: PetscWeakFormGetNumFields(), PetscWeakFormCreate() 1361 @*/ 1362 PetscErrorCode PetscWeakFormSetNumFields(PetscWeakForm wf, PetscInt Nf) 1363 { 1364 PetscFunctionBegin; 1365 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 1366 wf->Nf = Nf; 1367 PetscFunctionReturn(0); 1368 } 1369 1370 /*@ 1371 PetscWeakFormDestroy - Destroys a PetscWeakForm object 1372 1373 Collective on wf 1374 1375 Input Parameter: 1376 . wf - the PetscWeakForm object to destroy 1377 1378 Level: developer 1379 1380 .seealso PetscWeakFormCreate(), PetscWeakFormView() 1381 @*/ 1382 PetscErrorCode PetscWeakFormDestroy(PetscWeakForm *wf) 1383 { 1384 PetscInt f; 1385 PetscErrorCode ierr; 1386 1387 PetscFunctionBegin; 1388 if (!*wf) PetscFunctionReturn(0); 1389 PetscValidHeaderSpecific((*wf), PETSCWEAKFORM_CLASSID, 1); 1390 1391 if (--((PetscObject)(*wf))->refct > 0) {*wf = NULL; PetscFunctionReturn(0);} 1392 ((PetscObject) (*wf))->refct = 0; 1393 ierr = PetscChunkBufferDestroy(&(*wf)->funcs);CHKERRQ(ierr); 1394 for (f = 0; f < PETSC_NUM_WF; ++f) {ierr = PetscHMapFormDestroy(&(*wf)->form[f]);CHKERRQ(ierr);} 1395 ierr = PetscFree((*wf)->form);CHKERRQ(ierr); 1396 ierr = PetscHeaderDestroy(wf);CHKERRQ(ierr); 1397 PetscFunctionReturn(0); 1398 } 1399 1400 static PetscErrorCode PetscWeakFormViewTable_Ascii(PetscWeakForm wf, PetscViewer viewer, PetscBool splitField, const char tableName[], PetscHMapForm map) 1401 { 1402 PetscInt Nf = wf->Nf, Nk, k; 1403 PetscErrorCode ierr; 1404 1405 PetscFunctionBegin; 1406 ierr = PetscHMapFormGetSize(map, &Nk);CHKERRQ(ierr); 1407 if (Nk) { 1408 PetscFormKey *keys; 1409 void (**funcs)(void); 1410 const char *name; 1411 PetscInt off = 0, n, i; 1412 1413 ierr = PetscMalloc1(Nk, &keys);CHKERRQ(ierr); 1414 ierr = PetscHMapFormGetKeys(map, &off, keys);CHKERRQ(ierr); 1415 ierr = PetscViewerASCIIPrintf(viewer, "%s\n", tableName);CHKERRQ(ierr); 1416 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 1417 for (k = 0; k < Nk; ++k) { 1418 if (keys[k].label) { 1419 ierr = PetscObjectGetName((PetscObject) keys[k].label, &name);CHKERRQ(ierr); 1420 ierr = PetscViewerASCIIPrintf(viewer, "(%s, %D) ", name, keys[k].value);CHKERRQ(ierr); 1421 } else {ierr = PetscViewerASCIIPrintf(viewer, "");CHKERRQ(ierr);} 1422 ierr = PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);CHKERRQ(ierr); 1423 if (splitField) {ierr = PetscViewerASCIIPrintf(viewer, "(%D, %D) ", keys[k].field/Nf, keys[k].field%Nf);CHKERRQ(ierr);} 1424 else {ierr = PetscViewerASCIIPrintf(viewer, "(%D) ", keys[k].field);CHKERRQ(ierr);} 1425 ierr = PetscWeakFormGetFunction_Private(wf, map, keys[k].label, keys[k].value, keys[k].field, keys[k].part, &n, &funcs);CHKERRQ(ierr); 1426 for (i = 0; i < n; ++i) { 1427 char *fname; 1428 1429 if (i > 0) {ierr = PetscViewerASCIIPrintf(viewer, ", ");CHKERRQ(ierr);} 1430 ierr = PetscDLAddr(funcs[i], &fname);CHKERRQ(ierr); 1431 if (fname) {ierr = PetscViewerASCIIPrintf(viewer, "%s", fname);CHKERRQ(ierr);} 1432 else {ierr = PetscViewerASCIIPrintf(viewer, "%p", funcs[i]);CHKERRQ(ierr);} 1433 ierr = PetscFree(fname);CHKERRQ(ierr); 1434 } 1435 ierr = PetscViewerASCIIPrintf(viewer, "\n");CHKERRQ(ierr); 1436 ierr = PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);CHKERRQ(ierr); 1437 } 1438 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 1439 ierr = PetscFree(keys);CHKERRQ(ierr); 1440 } 1441 PetscFunctionReturn(0); 1442 } 1443 1444 static PetscErrorCode PetscWeakFormView_Ascii(PetscWeakForm wf, PetscViewer viewer) 1445 { 1446 PetscViewerFormat format; 1447 PetscInt f; 1448 PetscErrorCode ierr; 1449 1450 PetscFunctionBegin; 1451 ierr = PetscViewerGetFormat(viewer, &format);CHKERRQ(ierr); 1452 ierr = PetscViewerASCIIPrintf(viewer, "Weak Form System with %d fields\n", wf->Nf);CHKERRQ(ierr); 1453 ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr); 1454 for (f = 0; f < PETSC_NUM_WF; ++f) { 1455 ierr = PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_FALSE, PetscWeakFormKinds[f], wf->form[f]);CHKERRQ(ierr); 1456 } 1457 ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr); 1458 PetscFunctionReturn(0); 1459 } 1460 1461 /*@C 1462 PetscWeakFormView - Views a PetscWeakForm 1463 1464 Collective on wf 1465 1466 Input Parameter: 1467 + wf - the PetscWeakForm object to view 1468 - v - the viewer 1469 1470 Level: developer 1471 1472 .seealso PetscWeakFormDestroy(), PetscWeakFormCreate() 1473 @*/ 1474 PetscErrorCode PetscWeakFormView(PetscWeakForm wf, PetscViewer v) 1475 { 1476 PetscBool iascii; 1477 PetscErrorCode ierr; 1478 1479 PetscFunctionBegin; 1480 PetscValidHeaderSpecific(wf, PETSCWEAKFORM_CLASSID, 1); 1481 if (!v) {ierr = PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject) wf), &v);CHKERRQ(ierr);} 1482 else {PetscValidHeaderSpecific(v, PETSC_VIEWER_CLASSID, 2);} 1483 ierr = PetscObjectTypeCompare((PetscObject) v, PETSCVIEWERASCII, &iascii);CHKERRQ(ierr); 1484 if (iascii) {ierr = PetscWeakFormView_Ascii(wf, v);CHKERRQ(ierr);} 1485 if (wf->ops->view) {ierr = (*wf->ops->view)(wf, v);CHKERRQ(ierr);} 1486 PetscFunctionReturn(0); 1487 } 1488 1489 /*@ 1490 PetscWeakFormCreate - Creates an empty PetscWeakForm object. 1491 1492 Collective 1493 1494 Input Parameter: 1495 . comm - The communicator for the PetscWeakForm object 1496 1497 Output Parameter: 1498 . wf - The PetscWeakForm object 1499 1500 Level: beginner 1501 1502 .seealso: PetscDS, PetscWeakFormDestroy() 1503 @*/ 1504 PetscErrorCode PetscWeakFormCreate(MPI_Comm comm, PetscWeakForm *wf) 1505 { 1506 PetscWeakForm p; 1507 PetscInt f; 1508 PetscErrorCode ierr; 1509 1510 PetscFunctionBegin; 1511 PetscValidPointer(wf, 2); 1512 *wf = NULL; 1513 ierr = PetscDSInitializePackage();CHKERRQ(ierr); 1514 1515 ierr = PetscHeaderCreate(p, PETSCWEAKFORM_CLASSID, "PetscWeakForm", "Weak Form System", "PetscWeakForm", comm, PetscWeakFormDestroy, PetscWeakFormView);CHKERRQ(ierr); 1516 1517 p->Nf = 0; 1518 ierr = PetscChunkBufferCreate(sizeof(&PetscWeakFormCreate), 2, &p->funcs);CHKERRQ(ierr); 1519 ierr = PetscMalloc1(PETSC_NUM_WF, &p->form);CHKERRQ(ierr); 1520 for (f = 0; f < PETSC_NUM_WF; ++f) {ierr = PetscHMapFormCreate(&p->form[f]);CHKERRQ(ierr);} 1521 *wf = p; 1522 PetscFunctionReturn(0); 1523 } 1524