Actual source code: dtweakform.c

  1: #include <petsc/private/petscdsimpl.h>

  3: PetscClassId PETSCWEAKFORM_CLASSID = 0;

  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};

  7: static PetscErrorCode PetscChunkBufferCreate(size_t unitbytes, size_t expected, PetscChunkBuffer **buffer)
  8: {
  9:   PetscNew(buffer);
 10:   PetscCalloc1(expected*unitbytes, &(*buffer)->array);
 11:   (*buffer)->size      = expected;
 12:   (*buffer)->unitbytes = unitbytes;
 13:   (*buffer)->alloc     = expected*unitbytes;
 14:   return 0;
 15: }

 17: static PetscErrorCode PetscChunkBufferDuplicate(PetscChunkBuffer *buffer, PetscChunkBuffer **bufferNew)
 18: {
 19:   PetscNew(bufferNew);
 20:   PetscCalloc1(buffer->size*buffer->unitbytes, &(*bufferNew)->array);
 21:   PetscMemcpy((*bufferNew)->array, buffer->array, buffer->size*buffer->unitbytes);
 22:   (*bufferNew)->size      = buffer->size;
 23:   (*bufferNew)->unitbytes = buffer->unitbytes;
 24:   (*bufferNew)->alloc     = buffer->size*buffer->unitbytes;
 25:   return 0;
 26: }

 28: static PetscErrorCode PetscChunkBufferDestroy(PetscChunkBuffer **buffer)
 29: {
 30:   PetscFree((*buffer)->array);
 31:   PetscFree(*buffer);
 32:   return 0;
 33: }

 35: static PetscErrorCode PetscChunkBufferCreateChunk(PetscChunkBuffer *buffer, PetscInt size, PetscChunk *chunk)
 36: {
 37:   if ((buffer->size + size)*buffer->unitbytes > buffer->alloc) {
 38:     char *tmp;

 40:     if (!buffer->alloc) buffer->alloc = (buffer->size + size)*buffer->unitbytes;
 41:     while ((buffer->size + size)*buffer->unitbytes > buffer->alloc) buffer->alloc *= 2;
 42:     PetscMalloc(buffer->alloc, &tmp);
 43:     PetscMemcpy(tmp, buffer->array, buffer->size*buffer->unitbytes);
 44:     PetscFree(buffer->array);
 45:     buffer->array = tmp;
 46:   }
 47:   chunk->start    = buffer->size*buffer->unitbytes;
 48:   chunk->size     = size;
 49:   chunk->reserved = size;
 50:   buffer->size   += size;
 51:   return 0;
 52: }

 54: static PetscErrorCode PetscChunkBufferEnlargeChunk(PetscChunkBuffer *buffer, PetscInt size, PetscChunk *chunk)
 55: {
 56:   size_t         siz = size;

 58:   if (chunk->size + size > chunk->reserved) {
 59:     PetscChunk newchunk;
 60:     PetscInt   reserved = chunk->size;

 62:     /* TODO Here if we had a chunk list, we could update them all to reclaim unused space */
 63:     while (reserved < chunk->size+size) reserved *= 2;
 64:     PetscChunkBufferCreateChunk(buffer, (size_t) reserved, &newchunk);
 65:     newchunk.size = chunk->size+size;
 66:     PetscMemcpy(&buffer->array[newchunk.start], &buffer->array[chunk->start], chunk->size * buffer->unitbytes);
 67:     *chunk = newchunk;
 68:   } else {
 69:     chunk->size += siz;
 70:   }
 71:   return 0;
 72: }

 74: /*@C
 75:   PetscFormKeySort - Sorts an array of PetscFormKey in place in increasing order.

 77:   Not Collective

 79:   Input Parameters:
 80: + n - number of values
 81: - X - array of PetscFormKey

 83:   Level: intermediate

 85: .seealso: PetscIntSortSemiOrdered(), PetscSortInt()
 86: @*/
 87: PetscErrorCode PetscFormKeySort(PetscInt n, PetscFormKey arr[])
 88: {
 89:   if (n <= 1) return 0;
 91:   PetscTimSort(n, arr, sizeof(PetscFormKey), Compare_PetscFormKey_Private, NULL);
 92:   return 0;
 93: }

 95: PetscErrorCode PetscWeakFormGetFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt *n, void (***func)())
 96: {
 97:   PetscFormKey   key;
 98:   PetscChunk     chunk;

100:   key.label = label; key.value = value; key.field = f; key.part = part;
101:   PetscHMapFormGet(ht, key, &chunk);
102:   if (chunk.size < 0) {*n = 0;          *func = NULL;}
103:   else                {*n = chunk.size; *func = (void (**)()) &wf->funcs->array[chunk.start];}
104:   return 0;
105: }

107: /* A NULL argument for func causes this to clear the key */
108: PetscErrorCode PetscWeakFormSetFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt n, void (**func)())
109: {
110:   PetscFormKey   key;
111:   PetscChunk     chunk;
112:   PetscInt       i;

114:   key.label = label; key.value = value; key.field = f; key.part = part;
115:   if (!func) {
116:     PetscHMapFormDel(ht, key);
117:     return 0;
118:   } else {
119:     PetscHMapFormGet(ht, key, &chunk);
120:   }
121:   if (chunk.size < 0) {
122:     PetscChunkBufferCreateChunk(wf->funcs, n, &chunk);
123:     PetscHMapFormSet(ht, key, chunk);
124:   } else if (chunk.size <= n) {
125:     PetscChunkBufferEnlargeChunk(wf->funcs, n - chunk.size, &chunk);
126:     PetscHMapFormSet(ht, key, chunk);
127:   }
128:   for (i = 0; i < n; ++i) ((void (**)()) &wf->funcs->array[chunk.start])[i] = func[i];
129:   return 0;
130: }

132: PetscErrorCode PetscWeakFormAddFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, void (*func)())
133: {
134:   PetscFormKey   key;
135:   PetscChunk     chunk;

137:   if (!func) return 0;
138:   key.label = label; key.value = value; key.field = f; key.part = part;
139:   PetscHMapFormGet(ht, key, &chunk);
140:   if (chunk.size < 0) {
141:     PetscChunkBufferCreateChunk(wf->funcs, 1, &chunk);
142:     PetscHMapFormSet(ht, key, chunk);
143:     ((void (**)()) &wf->funcs->array[chunk.start])[0] = func;
144:   } else {
145:     PetscChunkBufferEnlargeChunk(wf->funcs, 1, &chunk);
146:     PetscHMapFormSet(ht, key, chunk);
147:     ((void (**)()) &wf->funcs->array[chunk.start])[chunk.size-1] = func;
148:   }
149:   return 0;
150: }

152: PetscErrorCode PetscWeakFormGetIndexFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt ind, void (**func)())
153: {
154:   PetscFormKey key;
155:   PetscChunk       chunk;

157:   key.label = label; key.value = value; key.field = f; key.part = part;
158:   PetscHMapFormGet(ht, key, &chunk);
159:   if (chunk.size < 0) {*func = NULL;}
160:   else {
162:     *func = ((void (**)()) &wf->funcs->array[chunk.start])[ind];
163:   }
164:   return 0;
165: }

167: /* Ignore a NULL func */
168: PetscErrorCode PetscWeakFormSetIndexFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt ind, void (*func)())
169: {
170:   PetscFormKey key;
171:   PetscChunk       chunk;

173:   if (!func) return 0;
174:   key.label = label; key.value = value; key.field = f; key.part = part;
175:   PetscHMapFormGet(ht, key, &chunk);
176:   if (chunk.size < 0) {
177:     PetscChunkBufferCreateChunk(wf->funcs, ind+1, &chunk);
178:     PetscHMapFormSet(ht, key, chunk);
179:   } else if (chunk.size <= ind) {
180:     PetscChunkBufferEnlargeChunk(wf->funcs, ind - chunk.size + 1, &chunk);
181:     PetscHMapFormSet(ht, key, chunk);
182:   }
183:   ((void (**)()) &wf->funcs->array[chunk.start])[ind] = func;
184:   return 0;
185: }

187: PetscErrorCode PetscWeakFormClearIndexFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt ind)
188: {
189:   PetscFormKey key;
190:   PetscChunk       chunk;

192:   key.label = label; key.value = value; key.field = f; key.part = part;
193:   PetscHMapFormGet(ht, key, &chunk);
194:   if (chunk.size < 0) {
195:     return 0;
196:   } else if (!ind && chunk.size == 1) {
197:     PetscHMapFormDel(ht, key);
198:     return 0;
199:   } else if (chunk.size <= ind) {
200:     return 0;
201:   }
202:   ((void (**)()) &wf->funcs->array[chunk.start])[ind] = NULL;
203:   return 0;
204: }

206: /*@
207:   PetscWeakFormCopy - Copy the pointwise functions to another PetscWeakForm

209:   Not Collective

211:   Input Parameter:
212: . wf - The original PetscWeakForm

214:   Output Parameter:
215: . wfNew - The copy PetscWeakForm

217:   Level: intermediate

219: .seealso: PetscWeakFormCreate(), PetscWeakFormDestroy()
220: @*/
221: PetscErrorCode PetscWeakFormCopy(PetscWeakForm wf, PetscWeakForm wfNew)
222: {
223:   PetscInt       f;

225:   wfNew->Nf = wf->Nf;
226:   PetscChunkBufferDestroy(&wfNew->funcs);
227:   PetscChunkBufferDuplicate(wf->funcs, &wfNew->funcs);
228:   for (f = 0; f < PETSC_NUM_WF; ++f) {
229:     PetscHMapFormDestroy(&wfNew->form[f]);
230:     PetscHMapFormDuplicate(wf->form[f], &wfNew->form[f]);
231:   }
232:   return 0;
233: }

235: /*@
236:   PetscWeakFormClear - Clear all functions from the PetscWeakForm

238:   Not Collective

240:   Input Parameter:
241: . wf - The original PetscWeakForm

243:   Level: intermediate

245: .seealso: PetscWeakFormCopy(), PetscWeakFormCreate(), PetscWeakFormDestroy()
246: @*/
247: PetscErrorCode PetscWeakFormClear(PetscWeakForm wf)
248: {
249:   PetscInt       f;

251:   for (f = 0; f < PETSC_NUM_WF; ++f) PetscHMapFormClear(wf->form[f]);
252:   return 0;
253: }

255: static PetscErrorCode PetscWeakFormRewriteKeys_Internal(PetscWeakForm wf, PetscHMapForm hmap, DMLabel label, PetscInt Nv, const PetscInt values[])
256: {
257:   PetscFormKey  *keys;
258:   PetscInt       n, i, v, off = 0;

260:   PetscHMapFormGetSize(hmap, &n);
261:   PetscMalloc1(n, &keys);
262:   PetscHMapFormGetKeys(hmap, &off, keys);
263:   for (i = 0; i < n; ++i) {
264:     if (keys[i].label == label) {
265:       PetscBool clear = PETSC_TRUE;
266:       void   (**funcs)();
267:       PetscInt  Nf;

269:       PetscWeakFormGetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &Nf, &funcs);
270:       for (v = 0; v < Nv; ++v) {
271:         PetscWeakFormSetFunction_Private(wf, hmap, keys[i].label, values[v], keys[i].field, keys[i].part, Nf, funcs);
272:         if (values[v] == keys[i].value) clear = PETSC_FALSE;
273:       }
274:       if (clear) PetscWeakFormSetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, 0, NULL);
275:     }
276:   }
277:   PetscFree(keys);
278:   return 0;
279: }

281: /*@C
282:   PetscWeakFormRewriteKeys - Change any key on the given label to use the new set of label values

284:   Not Collective

286:   Input Parameters:
287: + wf     - The original PetscWeakForm
288: . label  - The label to change keys for
289: . Nv     - The number of new label values
290: - values - The set of new values to relabel keys with

292:   Note: This is used internally when boundary label values are specified from the command line.

294:   Level: intermediate

296: .seealso: PetscWeakFormReplaceLabel(), PetscWeakFormCreate(), PetscWeakFormDestroy()
297: @*/
298: PetscErrorCode PetscWeakFormRewriteKeys(PetscWeakForm wf, DMLabel label, PetscInt Nv, const PetscInt values[])
299: {
300:   PetscInt       f;

302:   for (f = 0; f < PETSC_NUM_WF; ++f) PetscWeakFormRewriteKeys_Internal(wf, wf->form[f], label, Nv, values);
303:   return 0;
304: }

306: static PetscErrorCode PetscWeakFormReplaceLabel_Internal(PetscWeakForm wf, PetscHMapForm hmap, DMLabel label)
307: {
308:   PetscFormKey  *keys;
309:   PetscInt       n, i, off = 0, maxFuncs = 0;
310:   void       (**tmpf)();
311:   const char    *name = NULL;

313:   if (label) PetscObjectGetName((PetscObject) label, &name);
314:   PetscHMapFormGetSize(hmap, &n);
315:   PetscMalloc1(n, &keys);
316:   PetscHMapFormGetKeys(hmap, &off, keys);
317:   for (i = 0; i < n; ++i) {
318:     PetscBool   match = PETSC_FALSE;
319:     const char *lname = NULL;

321:     if (label == keys[i].label) continue;
322:     if (keys[i].label) PetscObjectGetName((PetscObject) keys[i].label, &lname);
323:     PetscStrcmp(name, lname, &match);
324:     if ((!name && !lname) || match) {
325:       void  (**funcs)();
326:       PetscInt Nf;

328:       PetscWeakFormGetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &Nf, &funcs);
329:       maxFuncs = PetscMax(maxFuncs, Nf);
330:     }
331:   }
332:   /* Need temp space because chunk buffer can be reallocated in SetFunction() call */
333:   PetscMalloc1(maxFuncs, &tmpf);
334:   for (i = 0; i < n; ++i) {
335:     PetscBool   match = PETSC_FALSE;
336:     const char *lname = NULL;

338:     if (label == keys[i].label) continue;
339:     if (keys[i].label) PetscObjectGetName((PetscObject) keys[i].label, &lname);
340:     PetscStrcmp(name, lname, &match);
341:     if ((!name && !lname) || match) {
342:       void  (**funcs)();
343:       PetscInt Nf, j;

345:       PetscWeakFormGetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &Nf, &funcs);
346:       for (j = 0; j < Nf; ++j) tmpf[j] = funcs[j];
347:       PetscWeakFormSetFunction_Private(wf, hmap, label,         keys[i].value, keys[i].field, keys[i].part,  Nf,  tmpf);
348:       PetscWeakFormSetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part,  0,   NULL);
349:     }
350:   }
351:   PetscFree(tmpf);
352:   PetscFree(keys);
353:   return 0;
354: }

356: /*@C
357:   PetscWeakFormReplaceLabel - Change any key on a label of the same name to use the new label

359:   Not Collective

361:   Input Parameters:
362: + wf    - The original PetscWeakForm
363: - label - The label to change keys for

365:   Note: This is used internally when meshes are modified

367:   Level: intermediate

369: .seealso: PetscWeakFormRewriteKeys(), PetscWeakFormCreate(), PetscWeakFormDestroy()
370: @*/
371: PetscErrorCode PetscWeakFormReplaceLabel(PetscWeakForm wf, DMLabel label)
372: {
373:   PetscInt       f;

375:   for (f = 0; f < PETSC_NUM_WF; ++f) PetscWeakFormReplaceLabel_Internal(wf, wf->form[f], label);
376:   return 0;
377: }

379: PetscErrorCode PetscWeakFormClearIndex(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscWeakFormKind kind, PetscInt ind)
380: {
381:   PetscWeakFormClearIndexFunction_Private(wf, wf->form[kind], label, val, f, part, ind);
382:   return 0;
383: }

385: PetscErrorCode PetscWeakFormGetObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt *n,
386:                                          void (***obj)(PetscInt, PetscInt, PetscInt,
387:                                                        const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
388:                                                        const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
389:                                                        PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
390: {
391:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, n, (void (***)(void)) obj);
392:   return 0;
393: }

395: PetscErrorCode PetscWeakFormSetObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt n,
396:                                          void (**obj)(PetscInt, PetscInt, PetscInt,
397:                                                       const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const  PetscScalar[],
398:                                                       const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
399:                                                       PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
400: {
401:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, n, (void (**)(void)) obj);
402:   return 0;
403: }

405: PetscErrorCode PetscWeakFormAddObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part,
406:                                          void (*obj)(PetscInt, PetscInt, PetscInt,
407:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
408:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
409:                                                      PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
410: {
411:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, (void (*)(void)) obj);
412:   return 0;
413: }

415: PetscErrorCode PetscWeakFormGetIndexObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt ind,
416:                                               void (**obj)(PetscInt, PetscInt, PetscInt,
417:                                                            const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
418:                                                            const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
419:                                                            PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
420: {
421:   PetscWeakFormGetIndexFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, ind, (void (**)(void)) obj);
422:   return 0;
423: }

425: PetscErrorCode PetscWeakFormSetIndexObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt ind,
426:                                               void (*obj)(PetscInt, PetscInt, PetscInt,
427:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
428:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
429:                                                           PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
430: {
431:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, ind, (void (*)(void)) obj);
432:   return 0;
433: }

435: PetscErrorCode PetscWeakFormGetResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part,
436:                                         PetscInt *n0,
437:                                         void (***f0)(PetscInt, PetscInt, PetscInt,
438:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
439:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
440:                                                      PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
441:                                         PetscInt *n1,
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:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, n0, (void (***)(void)) f0);
448:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, n1, (void (***)(void)) f1);
449:   return 0;
450: }

452: PetscErrorCode PetscWeakFormAddResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part,
453:                                         void (*f0)(PetscInt, PetscInt, PetscInt,
454:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
455:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
456:                                                    PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
457:                                         void (*f1)(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: {
462:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, (void (*)(void)) f0);
463:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, (void (*)(void)) f1);
464:   return 0;
465: }

467: PetscErrorCode PetscWeakFormSetResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part,
468:                                         PetscInt n0,
469:                                         void (**f0)(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:                                         PetscInt n1,
474:                                         void (**f1)(PetscInt, PetscInt, PetscInt,
475:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
476:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
477:                                                    PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
478: {
479:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, n0, (void (**)(void)) f0);
480:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, n1, (void (**)(void)) f1);
481:   return 0;
482: }

484: PetscErrorCode PetscWeakFormSetIndexResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part,
485:                                         PetscInt i0,
486:                                         void (*f0)(PetscInt, PetscInt, PetscInt,
487:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
488:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
489:                                                    PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
490:                                         PetscInt i1,
491:                                         void (*f1)(PetscInt, PetscInt, PetscInt,
492:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
493:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
494:                                                    PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
495: {
496:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, i0, (void (*)(void)) f0);
497:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, i1, (void (*)(void)) f1);
498:   return 0;
499: }

501: PetscErrorCode PetscWeakFormGetBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part,
502:                                           PetscInt *n0,
503:                                         void (***f0)(PetscInt, PetscInt, PetscInt,
504:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
505:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
506:                                                      PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
507:                                         PetscInt *n1,
508:                                         void (***f1)(PetscInt, PetscInt, PetscInt,
509:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
510:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
511:                                                      PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
512: {
513:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, n0, (void (***)(void)) f0);
514:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, n1, (void (***)(void)) f1);
515:   return 0;
516: }

518: PetscErrorCode PetscWeakFormAddBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part,
519:                                           void (*f0)(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[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
523:                                           void (*f1)(PetscInt, PetscInt, PetscInt,
524:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
525:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
526:                                                      PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
527: {
528:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, (void (*)(void)) f0);
529:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, (void (*)(void)) f1);
530:   return 0;
531: }

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:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, n0, (void (**)(void)) f0);
546:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, n1, (void (**)(void)) f1);
547:   return 0;
548: }

550: PetscErrorCode PetscWeakFormSetIndexBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part,
551:                                           PetscInt i0,
552:                                           void (*f0)(PetscInt, PetscInt, PetscInt,
553:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
554:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
555:                                                      PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
556:                                           PetscInt i1,
557:                                           void (*f1)(PetscInt, PetscInt, PetscInt,
558:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
559:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
560:                                                      PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
561: {
562:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, i0, (void (*)(void)) f0);
563:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, i1, (void (*)(void)) f1);
564:   return 0;
565: }

567: PetscErrorCode PetscWeakFormHasJacobian(PetscWeakForm wf, PetscBool *hasJac)
568: {
569:   PetscInt       n0, n1, n2, n3;

573:   PetscHMapFormGetSize(wf->form[PETSC_WF_G0], &n0);
574:   PetscHMapFormGetSize(wf->form[PETSC_WF_G1], &n1);
575:   PetscHMapFormGetSize(wf->form[PETSC_WF_G2], &n2);
576:   PetscHMapFormGetSize(wf->form[PETSC_WF_G3], &n3);
577:   *hasJac = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE;
578:   return 0;
579: }

581: PetscErrorCode PetscWeakFormGetJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
582:                                         PetscInt *n0,
583:                                         void (***g0)(PetscInt, PetscInt, PetscInt,
584:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
585:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
586:                                                      PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
587:                                         PetscInt *n1,
588:                                         void (***g1)(PetscInt, PetscInt, PetscInt,
589:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
590:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
591:                                                      PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
592:                                         PetscInt *n2,
593:                                         void (***g2)(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, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
597:                                         PetscInt *n3,
598:                                         void (***g3)(PetscInt, PetscInt, PetscInt,
599:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
600:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
601:                                                      PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
602: {
603:   PetscInt       find = f*wf->Nf + g;

605:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, n0, (void (***)(void)) g0);
606:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, n1, (void (***)(void)) g1);
607:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, n2, (void (***)(void)) g2);
608:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, n3, (void (***)(void)) g3);
609:   return 0;
610: }

612: PetscErrorCode PetscWeakFormAddJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
613:                                         void (*g0)(PetscInt, PetscInt, PetscInt,
614:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
615:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
616:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
617:                                         void (*g1)(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, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
621:                                         void (*g2)(PetscInt, PetscInt, PetscInt,
622:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
623:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
624:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
625:                                         void (*g3)(PetscInt, PetscInt, PetscInt,
626:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
627:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
628:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
629: {
630:   PetscInt       find = f*wf->Nf + g;

632:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, (void (*)(void)) g0);
633:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, (void (*)(void)) g1);
634:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, (void (*)(void)) g2);
635:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, (void (*)(void)) g3);
636:   return 0;
637: }

639: PetscErrorCode PetscWeakFormSetJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
640:                                         PetscInt n0,
641:                                         void (**g0)(PetscInt, PetscInt, PetscInt,
642:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
643:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
644:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
645:                                         PetscInt n1,
646:                                         void (**g1)(PetscInt, PetscInt, PetscInt,
647:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
648:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
649:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
650:                                         PetscInt n2,
651:                                         void (**g2)(PetscInt, PetscInt, PetscInt,
652:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
653:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
654:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
655:                                         PetscInt n3,
656:                                         void (**g3)(PetscInt, PetscInt, PetscInt,
657:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
658:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
659:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
660: {
661:   PetscInt       find = f*wf->Nf + g;

663:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, n0, (void (**)(void)) g0);
664:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, n1, (void (**)(void)) g1);
665:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, n2, (void (**)(void)) g2);
666:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, n3, (void (**)(void)) g3);
667:   return 0;
668: }

670: PetscErrorCode PetscWeakFormSetIndexJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
671:                                         PetscInt i0,
672:                                         void (*g0)(PetscInt, PetscInt, PetscInt,
673:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
674:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
675:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
676:                                         PetscInt i1,
677:                                         void (*g1)(PetscInt, PetscInt, PetscInt,
678:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
679:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
680:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
681:                                         PetscInt i2,
682:                                         void (*g2)(PetscInt, PetscInt, PetscInt,
683:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
684:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
685:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
686:                                         PetscInt i3,
687:                                         void (*g3)(PetscInt, PetscInt, PetscInt,
688:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
689:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
690:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
691: {
692:   PetscInt       find = f*wf->Nf + g;

694:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, i0, (void (*)(void)) g0);
695:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, i1, (void (*)(void)) g1);
696:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, i2, (void (*)(void)) g2);
697:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, i3, (void (*)(void)) g3);
698:   return 0;
699: }

701: PetscErrorCode PetscWeakFormHasJacobianPreconditioner(PetscWeakForm wf, PetscBool *hasJacPre)
702: {
703:   PetscInt       n0, n1, n2, n3;

707:   PetscHMapFormGetSize(wf->form[PETSC_WF_GP0], &n0);
708:   PetscHMapFormGetSize(wf->form[PETSC_WF_GP1], &n1);
709:   PetscHMapFormGetSize(wf->form[PETSC_WF_GP2], &n2);
710:   PetscHMapFormGetSize(wf->form[PETSC_WF_GP3], &n3);
711:   *hasJacPre = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE;
712:   return 0;
713: }

715: PetscErrorCode PetscWeakFormGetJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
716:                                                       PetscInt *n0,
717:                                                       void (***g0)(PetscInt, PetscInt, PetscInt,
718:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
719:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
720:                                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
721:                                                       PetscInt *n1,
722:                                                       void (***g1)(PetscInt, PetscInt, PetscInt,
723:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
724:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
725:                                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
726:                                                       PetscInt *n2,
727:                                                       void (***g2)(PetscInt, PetscInt, PetscInt,
728:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
729:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
730:                                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
731:                                                       PetscInt *n3,
732:                                                       void (***g3)(PetscInt, PetscInt, PetscInt,
733:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
734:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
735:                                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
736: {
737:   PetscInt       find = f*wf->Nf + g;

739:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, n0, (void (***)(void)) g0);
740:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, n1, (void (***)(void)) g1);
741:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, n2, (void (***)(void)) g2);
742:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, n3, (void (***)(void)) g3);
743:   return 0;
744: }

746: PetscErrorCode PetscWeakFormAddJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
747:                                         void (*g0)(PetscInt, PetscInt, PetscInt,
748:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
749:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
750:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
751:                                         void (*g1)(PetscInt, PetscInt, PetscInt,
752:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
753:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
754:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
755:                                         void (*g2)(PetscInt, PetscInt, PetscInt,
756:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
757:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
758:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
759:                                         void (*g3)(PetscInt, PetscInt, PetscInt,
760:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
761:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
762:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
763: {
764:   PetscInt       find = f*wf->Nf + g;

766:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, (void (*)(void)) g0);
767:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, (void (*)(void)) g1);
768:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, (void (*)(void)) g2);
769:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, (void (*)(void)) g3);
770:   return 0;
771: }

773: PetscErrorCode PetscWeakFormSetJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
774:                                                       PetscInt n0,
775:                                                       void (**g0)(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:                                                       PetscInt n1,
780:                                                       void (**g1)(PetscInt, PetscInt, PetscInt,
781:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
782:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
783:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
784:                                                       PetscInt n2,
785:                                                       void (**g2)(PetscInt, PetscInt, PetscInt,
786:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
787:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
788:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
789:                                                       PetscInt n3,
790:                                                       void (**g3)(PetscInt, PetscInt, PetscInt,
791:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
792:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
793:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
794: {
795:   PetscInt       find = f*wf->Nf + g;

797:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, n0, (void (**)(void)) g0);
798:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, n1, (void (**)(void)) g1);
799:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, n2, (void (**)(void)) g2);
800:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, n3, (void (**)(void)) g3);
801:   return 0;
802: }

804: PetscErrorCode PetscWeakFormSetIndexJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
805:                                                       PetscInt i0,
806:                                                       void (*g0)(PetscInt, PetscInt, PetscInt,
807:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
808:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
809:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
810:                                                       PetscInt i1,
811:                                                       void (*g1)(PetscInt, PetscInt, PetscInt,
812:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
813:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
814:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
815:                                                       PetscInt i2,
816:                                                       void (*g2)(PetscInt, PetscInt, PetscInt,
817:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
818:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
819:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
820:                                                       PetscInt i3,
821:                                                       void (*g3)(PetscInt, PetscInt, PetscInt,
822:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
823:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
824:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
825: {
826:   PetscInt       find = f*wf->Nf + g;

828:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, i0, (void (*)(void)) g0);
829:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, i1, (void (*)(void)) g1);
830:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, i2, (void (*)(void)) g2);
831:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, i3, (void (*)(void)) g3);
832:   return 0;
833: }

835: PetscErrorCode PetscWeakFormHasBdJacobian(PetscWeakForm wf, PetscBool *hasJac)
836: {
837:   PetscInt       n0, n1, n2, n3;

841:   PetscHMapFormGetSize(wf->form[PETSC_WF_BDG0], &n0);
842:   PetscHMapFormGetSize(wf->form[PETSC_WF_BDG1], &n1);
843:   PetscHMapFormGetSize(wf->form[PETSC_WF_BDG2], &n2);
844:   PetscHMapFormGetSize(wf->form[PETSC_WF_BDG3], &n3);
845:   *hasJac = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE;
846:   return 0;
847: }

849: PetscErrorCode PetscWeakFormGetBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
850:                                           PetscInt *n0,
851:                                           void (***g0)(PetscInt, PetscInt, PetscInt,
852:                                                        const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
853:                                                        const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
854:                                                        PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
855:                                           PetscInt *n1,
856:                                           void (***g1)(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[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
860:                                           PetscInt *n2,
861:                                           void (***g2)(PetscInt, PetscInt, PetscInt,
862:                                                        const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
863:                                                        const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
864:                                                        PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
865:                                           PetscInt *n3,
866:                                           void (***g3)(PetscInt, PetscInt, PetscInt,
867:                                                        const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
868:                                                        const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
869:                                                        PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
870: {
871:   PetscInt       find = f*wf->Nf + g;

873:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, n0, (void (***)(void)) g0);
874:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, n1, (void (***)(void)) g1);
875:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, n2, (void (***)(void)) g2);
876:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, n3, (void (***)(void)) g3);
877:   return 0;
878: }

880: PetscErrorCode PetscWeakFormAddBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
881:                                           void (*g0)(PetscInt, PetscInt, PetscInt,
882:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
883:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
884:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
885:                                           void (*g1)(PetscInt, PetscInt, PetscInt,
886:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
887:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
888:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
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:                                           void (*g3)(PetscInt, PetscInt, PetscInt,
894:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
895:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
896:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
897: {
898:   PetscInt       find = f*wf->Nf + g;

900:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, (void (*)(void)) g0);
901:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, (void (*)(void)) g1);
902:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, (void (*)(void)) g2);
903:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, (void (*)(void)) g3);
904:   return 0;
905: }

907: PetscErrorCode PetscWeakFormSetBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
908:                                           PetscInt n0,
909:                                           void (**g0)(PetscInt, PetscInt, PetscInt,
910:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
911:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
912:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
913:                                           PetscInt n1,
914:                                           void (**g1)(PetscInt, PetscInt, PetscInt,
915:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
916:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
917:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
918:                                           PetscInt n2,
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:                                           PetscInt n3,
924:                                           void (**g3)(PetscInt, PetscInt, PetscInt,
925:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
926:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
927:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
928: {
929:   PetscInt       find = f*wf->Nf + g;

931:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, n0, (void (**)(void)) g0);
932:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, n1, (void (**)(void)) g1);
933:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, n2, (void (**)(void)) g2);
934:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, n3, (void (**)(void)) g3);
935:   return 0;
936: }

938: PetscErrorCode PetscWeakFormSetIndexBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
939:                                           PetscInt i0,
940:                                           void (*g0)(PetscInt, PetscInt, PetscInt,
941:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
942:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
943:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
944:                                           PetscInt i1,
945:                                           void (*g1)(PetscInt, PetscInt, PetscInt,
946:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
947:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
948:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
949:                                           PetscInt i2,
950:                                           void (*g2)(PetscInt, PetscInt, PetscInt,
951:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
952:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
953:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
954:                                           PetscInt i3,
955:                                           void (*g3)(PetscInt, PetscInt, PetscInt,
956:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
957:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
958:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
959: {
960:   PetscInt       find = f*wf->Nf + g;

962:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, i0, (void (*)(void)) g0);
963:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, i1, (void (*)(void)) g1);
964:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, i2, (void (*)(void)) g2);
965:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, i3, (void (*)(void)) g3);
966:   return 0;
967: }

969: PetscErrorCode PetscWeakFormHasBdJacobianPreconditioner(PetscWeakForm wf, PetscBool *hasJacPre)
970: {
971:   PetscInt       n0, n1, n2, n3;

975:   PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP0], &n0);
976:   PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP1], &n1);
977:   PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP2], &n2);
978:   PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP3], &n3);
979:   *hasJacPre = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE;
980:   return 0;
981: }

983: PetscErrorCode PetscWeakFormGetBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
984:                                                         PetscInt *n0,
985:                                                         void (***g0)(PetscInt, PetscInt, PetscInt,
986:                                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
987:                                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
988:                                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
989:                                                         PetscInt *n1,
990:                                                         void (***g1)(PetscInt, PetscInt, PetscInt,
991:                                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
992:                                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
993:                                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
994:                                                         PetscInt *n2,
995:                                                         void (***g2)(PetscInt, PetscInt, PetscInt,
996:                                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
997:                                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
998:                                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
999:                                                         PetscInt *n3,
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;

1007:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, n0, (void (***)(void)) g0);
1008:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, n1, (void (***)(void)) g1);
1009:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, n2, (void (***)(void)) g2);
1010:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, n3, (void (***)(void)) g3);
1011:   return 0;
1012: }

1014: PetscErrorCode PetscWeakFormAddBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
1015:                                                         void (*g0)(PetscInt, PetscInt, PetscInt,
1016:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1017:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1018:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1019:                                                         void (*g1)(PetscInt, PetscInt, PetscInt,
1020:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1021:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1022:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1023:                                                         void (*g2)(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:                                                         void (*g3)(PetscInt, PetscInt, PetscInt,
1028:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1029:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1030:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1031: {
1032:   PetscInt       find = f*wf->Nf + g;

1034:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, (void (*)(void)) g0);
1035:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, (void (*)(void)) g1);
1036:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, (void (*)(void)) g2);
1037:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, (void (*)(void)) g3);
1038:   return 0;
1039: }

1041: PetscErrorCode PetscWeakFormSetBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
1042:                                                         PetscInt n0,
1043:                                                         void (**g0)(PetscInt, PetscInt, PetscInt,
1044:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1045:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1046:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1047:                                                         PetscInt n1,
1048:                                                         void (**g1)(PetscInt, PetscInt, PetscInt,
1049:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1050:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1051:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1052:                                                         PetscInt n2,
1053:                                                         void (**g2)(PetscInt, PetscInt, PetscInt,
1054:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1055:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1056:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1057:                                                         PetscInt n3,
1058:                                                         void (**g3)(PetscInt, PetscInt, PetscInt,
1059:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1060:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1061:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1062: {
1063:   PetscInt       find = f*wf->Nf + g;

1065:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, n0, (void (**)(void)) g0);
1066:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, n1, (void (**)(void)) g1);
1067:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, n2, (void (**)(void)) g2);
1068:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, n3, (void (**)(void)) g3);
1069:   return 0;
1070: }

1072: PetscErrorCode PetscWeakFormSetIndexBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
1073:                                                         PetscInt i0,
1074:                                                         void (*g0)(PetscInt, PetscInt, PetscInt,
1075:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1076:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1077:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1078:                                                         PetscInt i1,
1079:                                                         void (*g1)(PetscInt, PetscInt, PetscInt,
1080:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1081:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1082:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1083:                                                         PetscInt i2,
1084:                                                         void (*g2)(PetscInt, PetscInt, PetscInt,
1085:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1086:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1087:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1088:                                                         PetscInt i3,
1089:                                                         void (*g3)(PetscInt, PetscInt, PetscInt,
1090:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1091:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1092:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1093: {
1094:   PetscInt       find = f*wf->Nf + g;

1096:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, i0, (void (*)(void)) g0);
1097:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, i1, (void (*)(void)) g1);
1098:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, i2, (void (*)(void)) g2);
1099:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, i3, (void (*)(void)) g3);
1100:   return 0;
1101: }

1103: PetscErrorCode PetscWeakFormHasDynamicJacobian(PetscWeakForm wf, PetscBool *hasDynJac)
1104: {
1105:   PetscInt       n0, n1, n2, n3;

1109:   PetscHMapFormGetSize(wf->form[PETSC_WF_GT0], &n0);
1110:   PetscHMapFormGetSize(wf->form[PETSC_WF_GT1], &n1);
1111:   PetscHMapFormGetSize(wf->form[PETSC_WF_GT2], &n2);
1112:   PetscHMapFormGetSize(wf->form[PETSC_WF_GT3], &n3);
1113:   *hasDynJac = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE;
1114:   return 0;
1115: }

1117: PetscErrorCode PetscWeakFormGetDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
1118:                                         PetscInt *n0,
1119:                                         void (***g0)(PetscInt, PetscInt, PetscInt,
1120:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1121:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1122:                                                      PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1123:                                         PetscInt *n1,
1124:                                         void (***g1)(PetscInt, PetscInt, PetscInt,
1125:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1126:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1127:                                                      PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1128:                                         PetscInt *n2,
1129:                                         void (***g2)(PetscInt, PetscInt, PetscInt,
1130:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1131:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1132:                                                      PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1133:                                         PetscInt *n3,
1134:                                         void (***g3)(PetscInt, PetscInt, PetscInt,
1135:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1136:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1137:                                                      PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1138: {
1139:   PetscInt       find = f*wf->Nf + g;

1141:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, n0, (void (***)(void)) g0);
1142:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, n1, (void (***)(void)) g1);
1143:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, n2, (void (***)(void)) g2);
1144:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, n3, (void (***)(void)) g3);
1145:   return 0;
1146: }

1148: PetscErrorCode PetscWeakFormAddDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
1149:                                         void (*g0)(PetscInt, PetscInt, PetscInt,
1150:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1151:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1152:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1153:                                         void (*g1)(PetscInt, PetscInt, PetscInt,
1154:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1155:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1156:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1157:                                         void (*g2)(PetscInt, PetscInt, PetscInt,
1158:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1159:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1160:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1161:                                         void (*g3)(PetscInt, PetscInt, PetscInt,
1162:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1163:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1164:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1165: {
1166:   PetscInt       find = f*wf->Nf + g;

1168:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, (void (*)(void)) g0);
1169:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, (void (*)(void)) g1);
1170:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, (void (*)(void)) g2);
1171:   PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, (void (*)(void)) g3);
1172:   return 0;
1173: }

1175: PetscErrorCode PetscWeakFormSetDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
1176:                                                PetscInt n0,
1177:                                                void (**g0)(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 n1,
1182:                                                void (**g1)(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:                                                PetscInt n2,
1187:                                                void (**g2)(PetscInt, PetscInt, PetscInt,
1188:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1189:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1190:                                                           PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1191:                                                PetscInt n3,
1192:                                                void (**g3)(PetscInt, PetscInt, PetscInt,
1193:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1194:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1195:                                                           PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1196: {
1197:   PetscInt       find = f*wf->Nf + g;

1199:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, n0, (void (**)(void)) g0);
1200:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, n1, (void (**)(void)) g1);
1201:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, n2, (void (**)(void)) g2);
1202:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, n3, (void (**)(void)) g3);
1203:   return 0;
1204: }

1206: PetscErrorCode PetscWeakFormSetIndexDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g, PetscInt part,
1207:                                                PetscInt i0,
1208:                                                void (*g0)(PetscInt, PetscInt, PetscInt,
1209:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1210:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1211:                                                           PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1212:                                                PetscInt i1,
1213:                                                void (*g1)(PetscInt, PetscInt, PetscInt,
1214:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1215:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1216:                                                           PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1217:                                                PetscInt i2,
1218:                                                void (*g2)(PetscInt, PetscInt, PetscInt,
1219:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1220:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1221:                                                           PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1222:                                                PetscInt i3,
1223:                                                void (*g3)(PetscInt, PetscInt, PetscInt,
1224:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1225:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1226:                                                           PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1227: {
1228:   PetscInt       find = f*wf->Nf + g;

1230:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, i0, (void (*)(void)) g0);
1231:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, i1, (void (*)(void)) g1);
1232:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, i2, (void (*)(void)) g2);
1233:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, i3, (void (*)(void)) g3);
1234:   return 0;
1235: }

1237: PetscErrorCode PetscWeakFormGetRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscInt *n,
1238:                                              void (***r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *))
1239: {
1240:   PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, n, (void (***)(void)) r);
1241:   return 0;
1242: }

1244: PetscErrorCode PetscWeakFormSetRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part,
1245:                                              PetscInt n,
1246:                                              void (**r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *))
1247: {
1248:   PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, n, (void (**)(void)) r);
1249:   return 0;
1250: }

1252: PetscErrorCode PetscWeakFormSetIndexRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part,
1253:                                                   PetscInt i,
1254:                                                   void (*r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *))
1255: {
1256:   PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, i, (void (*)(void)) r);
1257:   return 0;
1258: }

1260: /*@
1261:   PetscWeakFormGetNumFields - Returns the number of fields

1263:   Not collective

1265:   Input Parameter:
1266: . wf - The PetscWeakForm object

1268:   Output Parameter:
1269: . Nf - The number of fields

1271:   Level: beginner

1273: .seealso: PetscWeakFormSetNumFields(), PetscWeakFormCreate()
1274: @*/
1275: PetscErrorCode PetscWeakFormGetNumFields(PetscWeakForm wf, PetscInt *Nf)
1276: {
1279:   *Nf = wf->Nf;
1280:   return 0;
1281: }

1283: /*@
1284:   PetscWeakFormSetNumFields - Sets the number of fields

1286:   Not collective

1288:   Input Parameters:
1289: + wf - The PetscWeakForm object
1290: - Nf - The number of fields

1292:   Level: beginner

1294: .seealso: PetscWeakFormGetNumFields(), PetscWeakFormCreate()
1295: @*/
1296: PetscErrorCode PetscWeakFormSetNumFields(PetscWeakForm wf, PetscInt Nf)
1297: {
1299:   wf->Nf = Nf;
1300:   return 0;
1301: }

1303: /*@
1304:   PetscWeakFormDestroy - Destroys a PetscWeakForm object

1306:   Collective on wf

1308:   Input Parameter:
1309: . wf - the PetscWeakForm object to destroy

1311:   Level: developer

1313: .seealso PetscWeakFormCreate(), PetscWeakFormView()
1314: @*/
1315: PetscErrorCode PetscWeakFormDestroy(PetscWeakForm *wf)
1316: {
1317:   PetscInt       f;

1319:   if (!*wf) return 0;

1322:   if (--((PetscObject)(*wf))->refct > 0) {*wf = NULL; return 0;}
1323:   ((PetscObject) (*wf))->refct = 0;
1324:   PetscChunkBufferDestroy(&(*wf)->funcs);
1325:   for (f = 0; f < PETSC_NUM_WF; ++f) PetscHMapFormDestroy(&(*wf)->form[f]);
1326:   PetscFree((*wf)->form);
1327:   PetscHeaderDestroy(wf);
1328:   return 0;
1329: }

1331: static PetscErrorCode PetscWeakFormViewTable_Ascii(PetscWeakForm wf, PetscViewer viewer, PetscBool splitField, const char tableName[], PetscHMapForm map)
1332: {
1333:   PetscInt       Nf = wf->Nf, Nk, k;

1335:   PetscHMapFormGetSize(map, &Nk);
1336:   if (Nk) {
1337:     PetscFormKey *keys;
1338:     void       (**funcs)(void);
1339:     const char  **names;
1340:     PetscInt     *values, *idx1, *idx2, *idx;
1341:     PetscBool     showPart = PETSC_FALSE, showPointer = PETSC_FALSE;
1342:     PetscInt      off = 0;

1344:     PetscMalloc6(Nk, &keys, Nk, &names, Nk, &values, Nk, &idx1, Nk, &idx2, Nk, &idx);
1345:     PetscHMapFormGetKeys(map, &off, keys);
1346:     /* Sort keys by label name and value */
1347:     {
1348:       /* First sort values */
1349:       for (k = 0; k < Nk; ++k) {values[k] = keys[k].value; idx1[k] = k;}
1350:       PetscSortIntWithPermutation(Nk, values, idx1);
1351:       /* If the string sort is stable, it will be sorted correctly overall */
1352:       for (k = 0; k < Nk; ++k) {
1353:         if (keys[idx1[k]].label) PetscObjectGetName((PetscObject) keys[idx1[k]].label, &names[k]);
1354:         else                     {names[k] = "";}
1355:         idx2[k] = k;
1356:       }
1357:       PetscSortStrWithPermutation(Nk, names, idx2);
1358:       for (k = 0; k < Nk; ++k) {
1359:         if (keys[k].label) PetscObjectGetName((PetscObject) keys[k].label, &names[k]);
1360:         else               {names[k] = "";}
1361:         idx[k] = idx1[idx2[k]];
1362:       }
1363:     }
1364:     PetscViewerASCIIPrintf(viewer, "%s\n", tableName);
1365:     PetscViewerASCIIPushTab(viewer);
1366:     for (k = 0; k < Nk; ++k) {
1367:       if (keys[k].part != 0) showPart = PETSC_TRUE;
1368:     }
1369:     for (k = 0; k < Nk; ++k) {
1370:       const PetscInt i = idx[k];
1371:       PetscInt       n, f;

1373:       if (keys[i].label) {
1374:         if (showPointer) PetscViewerASCIIPrintf(viewer, "(%s:%p, %D) ", names[i], keys[i].label, keys[i].value);
1375:         else             PetscViewerASCIIPrintf(viewer, "(%s, %D) ", names[i], keys[i].value);
1376:       } else PetscViewerASCIIPrintf(viewer, "");
1377:       PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
1378:       if (splitField) PetscViewerASCIIPrintf(viewer, "(%D, %D) ", keys[i].field/Nf, keys[i].field%Nf);
1379:       else            PetscViewerASCIIPrintf(viewer, "(%D) ", keys[i].field);
1380:       if (showPart)   PetscViewerASCIIPrintf(viewer, "(%D) ", keys[i].part);
1381:       PetscWeakFormGetFunction_Private(wf, map, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &n, &funcs);
1382:       for (f = 0; f < n; ++f) {
1383:         char  *fname;
1384:         size_t len, l;

1386:         if (f > 0) PetscViewerASCIIPrintf(viewer, ", ");
1387:         PetscDLAddr(funcs[f], &fname);
1388:         if (fname) {
1389:           /* Eliminate argument types */
1390:           PetscStrlen(fname, &len);
1391:           for (l = 0; l < len; ++l) if (fname[l] == '(') {fname[l] = '\0'; break;}
1392:           PetscViewerASCIIPrintf(viewer, "%s", fname);
1393:         } else if (showPointer) {
1394:           PetscViewerASCIIPrintf(viewer, "%p", funcs[f]);
1395:         }
1396:         PetscFree(fname);
1397:       }
1398:       PetscViewerASCIIPrintf(viewer, "\n");
1399:       PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
1400:     }
1401:     PetscViewerASCIIPopTab(viewer);
1402:     PetscFree6(keys, names, values, idx1, idx2, idx);
1403:   }
1404:   return 0;
1405: }

1407: static PetscErrorCode PetscWeakFormView_Ascii(PetscWeakForm wf, PetscViewer viewer)
1408: {
1409:   PetscViewerFormat format;
1410:   PetscInt          f;

1412:   PetscViewerGetFormat(viewer, &format);
1413:   PetscViewerASCIIPrintf(viewer, "Weak Form System with %d fields\n", wf->Nf);
1414:   PetscViewerASCIIPushTab(viewer);
1415:   for (f = 0; f < PETSC_NUM_WF; ++f) {
1416:     PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE, PetscWeakFormKinds[f], wf->form[f]);
1417:   }
1418:   PetscViewerASCIIPopTab(viewer);
1419:   return 0;
1420: }

1422: /*@C
1423:   PetscWeakFormView - Views a PetscWeakForm

1425:   Collective on wf

1427:   Input Parameters:
1428: + wf - the PetscWeakForm object to view
1429: - v  - the viewer

1431:   Level: developer

1433: .seealso PetscWeakFormDestroy(), PetscWeakFormCreate()
1434: @*/
1435: PetscErrorCode PetscWeakFormView(PetscWeakForm wf, PetscViewer v)
1436: {
1437:   PetscBool      iascii;

1440:   if (!v) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject) wf), &v);
1442:   PetscObjectTypeCompare((PetscObject) v, PETSCVIEWERASCII, &iascii);
1443:   if (iascii) PetscWeakFormView_Ascii(wf, v);
1444:   if (wf->ops->view) (*wf->ops->view)(wf, v);
1445:   return 0;
1446: }

1448: /*@
1449:   PetscWeakFormCreate - Creates an empty PetscWeakForm object.

1451:   Collective

1453:   Input Parameter:
1454: . comm - The communicator for the PetscWeakForm object

1456:   Output Parameter:
1457: . wf - The PetscWeakForm object

1459:   Level: beginner

1461: .seealso: PetscDS, PetscWeakFormDestroy()
1462: @*/
1463: PetscErrorCode PetscWeakFormCreate(MPI_Comm comm, PetscWeakForm *wf)
1464: {
1465:   PetscWeakForm  p;
1466:   PetscInt       f;

1469:   *wf  = NULL;
1470:   PetscDSInitializePackage();

1472:   PetscHeaderCreate(p, PETSCWEAKFORM_CLASSID, "PetscWeakForm", "Weak Form System", "PetscWeakForm", comm, PetscWeakFormDestroy, PetscWeakFormView);

1474:   p->Nf = 0;
1475:   PetscChunkBufferCreate(sizeof(&PetscWeakFormCreate), 2, &p->funcs);
1476:   PetscMalloc1(PETSC_NUM_WF, &p->form);
1477:   for (f = 0; f < PETSC_NUM_WF; ++f) PetscHMapFormCreate(&p->form[f]);
1478:   *wf = p;
1479:   return 0;
1480: }