Actual source code: dtweakform.c

petsc-main 2021-04-20
Report Typos and Errors
  1: #include <petsc/private/petscdsimpl.h>

  3: PetscClassId PETSCWEAKFORM_CLASSID = 0;

  5: static PetscErrorCode PetscChunkBufferCreate(size_t unitbytes, size_t expected, PetscChunkBuffer **buffer)
  6: {

 10:   PetscNew(buffer);
 11:   PetscCalloc1(expected*unitbytes, &(*buffer)->array);
 12:   (*buffer)->size      = expected;
 13:   (*buffer)->unitbytes = unitbytes;
 14:   (*buffer)->alloc     = expected*unitbytes;
 15:   return(0);
 16: }

 18: static PetscErrorCode PetscChunkBufferDuplicate(PetscChunkBuffer *buffer, PetscChunkBuffer **bufferNew)
 19: {

 23:   PetscNew(bufferNew);
 24:   PetscCalloc1(buffer->size*buffer->unitbytes, &(*bufferNew)->array);
 25:   PetscMemcpy((*bufferNew)->array, buffer->array, buffer->size*buffer->unitbytes);
 26:   (*bufferNew)->size      = buffer->size;
 27:   (*bufferNew)->unitbytes = buffer->unitbytes;
 28:   (*bufferNew)->alloc     = buffer->size*buffer->unitbytes;
 29:   return(0);
 30: }

 32: static PetscErrorCode PetscChunkBufferDestroy(PetscChunkBuffer **buffer)
 33: {

 37:   PetscFree((*buffer)->array);
 38:   PetscFree(*buffer);
 39:   return(0);
 40: }

 42: static PetscErrorCode PetscChunkBufferCreateChunk(PetscChunkBuffer *buffer, PetscInt size, PetscChunk *chunk)
 43: {

 47:   if ((buffer->size + size)*buffer->unitbytes > buffer->alloc) {
 48:     char *tmp;

 50:     if (!buffer->alloc) buffer->alloc = (buffer->size + size)*buffer->unitbytes;
 51:     while ((buffer->size + size)*buffer->unitbytes > buffer->alloc) buffer->alloc *= 2;
 52:     PetscMalloc(buffer->alloc, &tmp);
 53:     PetscMemcpy(tmp, buffer->array, buffer->size*buffer->unitbytes);
 54:     PetscFree(buffer->array);
 55:     buffer->array = tmp;
 56:   }
 57:   chunk->start    = buffer->size*buffer->unitbytes;
 58:   chunk->size     = size;
 59:   chunk->reserved = size;
 60:   buffer->size   += size;
 61:   return(0);
 62: }

 64: static PetscErrorCode PetscChunkBufferEnlargeChunk(PetscChunkBuffer *buffer, PetscInt size, PetscChunk *chunk)
 65: {
 66:   size_t         siz = size;

 70:   if (chunk->size + size > chunk->reserved) {
 71:     PetscChunk newchunk;
 72:     PetscInt   reserved = chunk->size;

 74:     /* TODO Here if we had a chunk list, we could update them all to reclaim unused space */
 75:     while (reserved < chunk->size+size) reserved *= 2;
 76:     PetscChunkBufferCreateChunk(buffer, (size_t) reserved, &newchunk);
 77:     newchunk.size = chunk->size+size;
 78:     PetscMemcpy(&buffer->array[newchunk.start], &buffer->array[chunk->start], chunk->size * buffer->unitbytes);
 79:     *chunk = newchunk;
 80:   } else {
 81:     chunk->size += siz;
 82:   }
 83:   return(0);
 84: }

 86: /*@C
 87:   PetscHashFormKeySort - Sorts an array of PetscHashFormKey in place in increasing order.

 89:   Not Collective

 91:   Input Parameters:
 92: + n - number of values
 93: - X - array of PetscHashFormKey

 95:   Level: intermediate

 97: .seealso: PetscIntSortSemiOrdered(), PetscSortInt()
 98: @*/
 99: PetscErrorCode PetscHashFormKeySort(PetscInt n, PetscHashFormKey arr[])
100: {

104:   if (n <= 1) return(0);
106:   PetscTimSort(n, arr, sizeof(PetscHashFormKey), Compare_PetscHashFormKey_Private, NULL);
107:   return(0);
108: }

110: PetscErrorCode PetscWeakFormGetFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt *n, void (***func)())
111: {
112:   PetscHashFormKey key;
113:   PetscChunk       chunk;
114:   PetscErrorCode   ierr;

117:   key.label = label; key.value = value; key.field = f;
118:   PetscHMapFormGet(ht, key, &chunk);
119:   if (chunk.size < 0) {*n = 0;          *func = NULL;}
120:   else                {*n = chunk.size; *func = (void (**)()) &wf->funcs->array[chunk.start];}
121:   return(0);
122: }

124: /* A NULL argument for func causes this to clear the key */
125: PetscErrorCode PetscWeakFormSetFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt n, void (**func)())
126: {
127:   PetscHashFormKey key;
128:   PetscChunk       chunk;
129:   PetscInt         i;
130:   PetscErrorCode   ierr;

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

151: PetscErrorCode PetscWeakFormAddFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, void (*func)())
152: {
153:   PetscHashFormKey key;
154:   PetscChunk       chunk;
155:   PetscErrorCode   ierr;

158:   if (!func) return(0);
159:   key.label = label; key.value = value; key.field = f;
160:   PetscHMapFormGet(ht, key, &chunk);
161:   if (chunk.size < 0) {
162:     PetscChunkBufferCreateChunk(wf->funcs, 1, &chunk);
163:     PetscHMapFormSet(ht, key, chunk);
164:     ((void (**)()) &wf->funcs->array[chunk.start])[0] = func;
165:   } else {
166:     PetscChunkBufferEnlargeChunk(wf->funcs, 1, &chunk);
167:     PetscHMapFormSet(ht, key, chunk);
168:     ((void (**)()) &wf->funcs->array[chunk.start])[chunk.size-1] = func;
169:   }
170:   return(0);
171: }

173: PetscErrorCode PetscWeakFormGetIndexFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt ind, void (**func)())
174: {
175:   PetscHashFormKey key;
176:   PetscChunk       chunk;
177:   PetscErrorCode   ierr;

180:   key.label = label; key.value = value; key.field = f;
181:   PetscHMapFormGet(ht, key, &chunk);
182:   if (chunk.size < 0) {*func = NULL;}
183:   else {
184:     if (ind >= chunk.size) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Index %D not in [0, %D)", ind, chunk.size);
185:     *func = ((void (**)()) &wf->funcs->array[chunk.start])[ind];
186:   }
187:   return(0);
188: }

190: /* A NULL argument for func causes this to clear the slot, and if there is nothing else, clear the key */
191: PetscErrorCode PetscWeakFormSetIndexFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt ind, void (*func)())
192: {
193:   PetscHashFormKey key;
194:   PetscChunk       chunk;
195:   PetscErrorCode   ierr;

198:   key.label = label; key.value = value; key.field = f;
199:   PetscHMapFormGet(ht, key, &chunk);
200:   if (chunk.size < 0) {
201:     if (!func) return(0);
202:     PetscChunkBufferCreateChunk(wf->funcs, ind+1, &chunk);
203:     PetscHMapFormSet(ht, key, chunk);
204:   } else if (!func && !ind && chunk.size == 1) {
205:     PetscHMapFormDel(ht, key);
206:     return(0);
207:   } else if (chunk.size <= ind) {
208:     PetscChunkBufferEnlargeChunk(wf->funcs, ind - chunk.size + 1, &chunk);
209:     PetscHMapFormSet(ht, key, chunk);
210:   }
211:   ((void (**)()) &wf->funcs->array[chunk.start])[ind] = func;
212:   return(0);
213: }

215: /*@
216:   PetscWeakFormCopy - Copy the pointwise functions to another PetscWeakForm

218:   Not Collective

220:   Input Parameter:
221: . wf - The original PetscWeakForm

223:   Output Parameter:
224: . wfNew - The copy PetscWeakForm

226:   Level: intermediate

228: .seealso: PetscWeakFormCreate(), PetscWeakFormDestroy()
229: @*/
230: PetscErrorCode PetscWeakFormCopy(PetscWeakForm wf, PetscWeakForm wfNew)
231: {

235:   wfNew->Nf = wf->Nf;
236:   PetscChunkBufferDestroy(&wfNew->funcs);
237:   PetscChunkBufferDuplicate(wf->funcs, &wfNew->funcs);
238:   PetscHMapFormDestroy(&wfNew->obj);
239:   PetscHMapFormDuplicate(wf->obj, &wfNew->obj);
240:   PetscHMapFormDestroy(&wfNew->f0);
241:   PetscHMapFormDuplicate(wf->f0, &wfNew->f0);
242:   PetscHMapFormDestroy(&wfNew->f1);
243:   PetscHMapFormDuplicate(wf->f1, &wfNew->f1);
244:   PetscHMapFormDestroy(&wfNew->g0);
245:   PetscHMapFormDuplicate(wf->g0, &wfNew->g0);
246:   PetscHMapFormDestroy(&wfNew->g1);
247:   PetscHMapFormDuplicate(wf->g1, &wfNew->g1);
248:   PetscHMapFormDestroy(&wfNew->g2);
249:   PetscHMapFormDuplicate(wf->g2, &wfNew->g2);
250:   PetscHMapFormDestroy(&wfNew->g3);
251:   PetscHMapFormDuplicate(wf->g3, &wfNew->g3);
252:   PetscHMapFormDestroy(&wfNew->gp0);
253:   PetscHMapFormDuplicate(wf->gp0, &wfNew->gp0);
254:   PetscHMapFormDestroy(&wfNew->gp1);
255:   PetscHMapFormDuplicate(wf->gp1, &wfNew->gp1);
256:   PetscHMapFormDestroy(&wfNew->gp2);
257:   PetscHMapFormDuplicate(wf->gp2, &wfNew->gp2);
258:   PetscHMapFormDestroy(&wfNew->gp3);
259:   PetscHMapFormDuplicate(wf->gp3, &wfNew->gp3);
260:   PetscHMapFormDestroy(&wfNew->gt0);
261:   PetscHMapFormDuplicate(wf->gt0, &wfNew->gt0);
262:   PetscHMapFormDestroy(&wfNew->gt1);
263:   PetscHMapFormDuplicate(wf->gt1, &wfNew->gt1);
264:   PetscHMapFormDestroy(&wfNew->gt2);
265:   PetscHMapFormDuplicate(wf->gt2, &wfNew->gt2);
266:   PetscHMapFormDestroy(&wfNew->gt3);
267:   PetscHMapFormDuplicate(wf->gt3, &wfNew->gt3);
268:   PetscHMapFormDestroy(&wfNew->bdf0);
269:   PetscHMapFormDuplicate(wf->bdf0, &wfNew->bdf0);
270:   PetscHMapFormDestroy(&wfNew->bdf1);
271:   PetscHMapFormDuplicate(wf->bdf1, &wfNew->bdf1);
272:   PetscHMapFormDestroy(&wfNew->bdg0);
273:   PetscHMapFormDuplicate(wf->bdg0, &wfNew->bdg0);
274:   PetscHMapFormDestroy(&wfNew->bdg1);
275:   PetscHMapFormDuplicate(wf->bdg1, &wfNew->bdg1);
276:   PetscHMapFormDestroy(&wfNew->bdg2);
277:   PetscHMapFormDuplicate(wf->bdg2, &wfNew->bdg2);
278:   PetscHMapFormDestroy(&wfNew->bdg3);
279:   PetscHMapFormDuplicate(wf->bdg3, &wfNew->bdg3);
280:   PetscHMapFormDestroy(&wfNew->bdgp0);
281:   PetscHMapFormDuplicate(wf->bdgp0, &wfNew->bdgp0);
282:   PetscHMapFormDestroy(&wfNew->bdgp1);
283:   PetscHMapFormDuplicate(wf->bdgp1, &wfNew->bdgp1);
284:   PetscHMapFormDestroy(&wfNew->bdgp2);
285:   PetscHMapFormDuplicate(wf->bdgp2, &wfNew->bdgp2);
286:   PetscHMapFormDestroy(&wfNew->bdgp3);
287:   PetscHMapFormDuplicate(wf->bdgp3, &wfNew->bdgp3);
288:   PetscHMapFormDestroy(&wfNew->r);
289:   PetscHMapFormDuplicate(wf->r, &wfNew->r);
290:   return(0);
291: }

293: static PetscErrorCode PetscWeakFormRewriteKeys_Internal(PetscWeakForm wf, PetscHMapForm hmap, DMLabel label, PetscInt Nv, const PetscInt values[])
294: {
295:   PetscHashFormKey *keys;
296:   PetscInt          n, i, v, off = 0;
297:   PetscErrorCode    ierr;

300:   PetscHMapFormGetSize(hmap, &n);
301:   PetscMalloc1(n, &keys);
302:   PetscHMapFormGetKeys(hmap, &off, keys);
303:   for (i = 0; i < n; ++i) {
304:     if (keys[i].label == label) {
305:       PetscBool clear = PETSC_TRUE;
306:       void   (**funcs)();
307:       PetscInt  Nf;

309:       PetscWeakFormGetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, &Nf, &funcs);
310:       for (v = 0; v < Nv; ++v) {
311:         PetscWeakFormSetFunction_Private(wf, hmap, keys[i].label, values[v], keys[i].field, Nf, funcs);
312:         if (values[v] == keys[i].value) clear = PETSC_FALSE;
313:       }
314:       if (clear) {PetscWeakFormSetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, 0, NULL);}
315:     }
316:   }
317:   PetscFree(keys);
318:   return(0);
319: }

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

324:   Not Collective

326:   Input Parameters:
327: + wf     - The original PetscWeakForm
328: . label  - The label to change keys for
329: . Nv     - The number of new label values
330: - values - The set of new values to relabel keys with

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

334:   Level: intermediate

336: .seealso: PetscWeakFormCreate(), PetscWeakFormDestroy()
337: @*/
338: PetscErrorCode PetscWeakFormRewriteKeys(PetscWeakForm wf, DMLabel label, PetscInt Nv, const PetscInt values[])
339: {

343:   PetscWeakFormRewriteKeys_Internal(wf, wf->obj,   label, Nv, values);
344:   PetscWeakFormRewriteKeys_Internal(wf, wf->f0,    label, Nv, values);
345:   PetscWeakFormRewriteKeys_Internal(wf, wf->f1,    label, Nv, values);
346:   PetscWeakFormRewriteKeys_Internal(wf, wf->g0,    label, Nv, values);
347:   PetscWeakFormRewriteKeys_Internal(wf, wf->g1,    label, Nv, values);
348:   PetscWeakFormRewriteKeys_Internal(wf, wf->g2,    label, Nv, values);
349:   PetscWeakFormRewriteKeys_Internal(wf, wf->g3,    label, Nv, values);
350:   PetscWeakFormRewriteKeys_Internal(wf, wf->gp0,   label, Nv, values);
351:   PetscWeakFormRewriteKeys_Internal(wf, wf->gp1,   label, Nv, values);
352:   PetscWeakFormRewriteKeys_Internal(wf, wf->gp2,   label, Nv, values);
353:   PetscWeakFormRewriteKeys_Internal(wf, wf->gp3,   label, Nv, values);
354:   PetscWeakFormRewriteKeys_Internal(wf, wf->gt0,   label, Nv, values);
355:   PetscWeakFormRewriteKeys_Internal(wf, wf->gt1,   label, Nv, values);
356:   PetscWeakFormRewriteKeys_Internal(wf, wf->gt2,   label, Nv, values);
357:   PetscWeakFormRewriteKeys_Internal(wf, wf->gt3,   label, Nv, values);
358:   PetscWeakFormRewriteKeys_Internal(wf, wf->bdf0,  label, Nv, values);
359:   PetscWeakFormRewriteKeys_Internal(wf, wf->bdf1,  label, Nv, values);
360:   PetscWeakFormRewriteKeys_Internal(wf, wf->bdg0,  label, Nv, values);
361:   PetscWeakFormRewriteKeys_Internal(wf, wf->bdg1,  label, Nv, values);
362:   PetscWeakFormRewriteKeys_Internal(wf, wf->bdg2,  label, Nv, values);
363:   PetscWeakFormRewriteKeys_Internal(wf, wf->bdg3,  label, Nv, values);
364:   PetscWeakFormRewriteKeys_Internal(wf, wf->bdgp0, label, Nv, values);
365:   PetscWeakFormRewriteKeys_Internal(wf, wf->bdgp1, label, Nv, values);
366:   PetscWeakFormRewriteKeys_Internal(wf, wf->bdgp2, label, Nv, values);
367:   PetscWeakFormRewriteKeys_Internal(wf, wf->bdgp3, label, Nv, values);
368:   PetscWeakFormRewriteKeys_Internal(wf, wf->r,     label, Nv, values);
369:   return(0);
370: }

372: PetscErrorCode PetscWeakFormGetObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt *n,
373:                                          void (***obj)(PetscInt, PetscInt, PetscInt,
374:                                                        const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
375:                                                        const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
376:                                                        PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
377: {

381:   PetscWeakFormGetFunction_Private(wf, wf->obj, label, val, f, n, (void (***)(void)) obj);
382:   return(0);
383: }

385: PetscErrorCode PetscWeakFormSetObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, 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: {

394:   PetscWeakFormSetFunction_Private(wf, wf->obj, label, val, f, n, (void (**)(void)) obj);
395:   return(0);
396: }

398: PetscErrorCode PetscWeakFormAddObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f,
399:                                          void (*obj)(PetscInt, PetscInt, PetscInt,
400:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
401:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
402:                                                      PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
403: {

407:   PetscWeakFormAddFunction_Private(wf, wf->obj, label, val, f, (void (*)(void)) obj);
408:   return(0);
409: }

411: PetscErrorCode PetscWeakFormGetIndexObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt ind,
412:                                               void (**obj)(PetscInt, PetscInt, PetscInt,
413:                                                            const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
414:                                                            const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
415:                                                            PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
416: {

420:   PetscWeakFormGetIndexFunction_Private(wf, wf->obj, label, val, f, ind, (void (**)(void)) obj);
421:   return(0);
422: }

424: PetscErrorCode PetscWeakFormSetIndexObjective(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt ind,
425:                                               void (*obj)(PetscInt, PetscInt, PetscInt,
426:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
427:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
428:                                                           PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
429: {

433:   PetscWeakFormSetIndexFunction_Private(wf, wf->obj, label, val, f, ind, (void (*)(void)) obj);
434:   return(0);
435: }

437: PetscErrorCode PetscWeakFormGetResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f,
438:                                         PetscInt *n0,
439:                                         void (***f0)(PetscInt, PetscInt, PetscInt,
440:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
441:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
442:                                                      PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
443:                                         PetscInt *n1,
444:                                         void (***f1)(PetscInt, PetscInt, PetscInt,
445:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
446:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
447:                                                      PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
448: {

452:   PetscWeakFormGetFunction_Private(wf, wf->f0, label, val, f, n0, (void (***)(void)) f0);
453:   PetscWeakFormGetFunction_Private(wf, wf->f1, label, val, f, n1, (void (***)(void)) f1);
454:   return(0);
455: }

457: PetscErrorCode PetscWeakFormAddResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f,
458:                                         void (*f0)(PetscInt, PetscInt, PetscInt,
459:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
460:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
461:                                                    PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
462:                                         void (*f1)(PetscInt, PetscInt, PetscInt,
463:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
464:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
465:                                                    PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
466: {

470:   PetscWeakFormAddFunction_Private(wf, wf->f0, label, val, f, (void (*)(void)) f0);
471:   PetscWeakFormAddFunction_Private(wf, wf->f1, label, val, f, (void (*)(void)) f1);
472:   return(0);
473: }

475: PetscErrorCode PetscWeakFormSetResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f,
476:                                         PetscInt n0,
477:                                         void (**f0)(PetscInt, PetscInt, PetscInt,
478:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
479:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
480:                                                    PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
481:                                         PetscInt n1,
482:                                         void (**f1)(PetscInt, PetscInt, PetscInt,
483:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
484:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
485:                                                    PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
486: {

490:   PetscWeakFormSetFunction_Private(wf, wf->f0, label, val, f, n0, (void (**)(void)) f0);
491:   PetscWeakFormSetFunction_Private(wf, wf->f1, label, val, f, n1, (void (**)(void)) f1);
492:   return(0);
493: }

495: PetscErrorCode PetscWeakFormSetIndexResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f,
496:                                         PetscInt i0,
497:                                         void (*f0)(PetscInt, PetscInt, PetscInt,
498:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
499:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
500:                                                    PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
501:                                         PetscInt i1,
502:                                         void (*f1)(PetscInt, PetscInt, PetscInt,
503:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
504:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
505:                                                    PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
506: {

510:   PetscWeakFormSetIndexFunction_Private(wf, wf->f0, label, val, f, i0, (void (*)(void)) f0);
511:   PetscWeakFormSetIndexFunction_Private(wf, wf->f1, label, val, f, i1, (void (*)(void)) f1);
512:   return(0);
513: }

515: PetscErrorCode PetscWeakFormGetBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f,
516:                                           PetscInt *n0,
517:                                         void (***f0)(PetscInt, PetscInt, PetscInt,
518:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
519:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
520:                                                      PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
521:                                         PetscInt *n1,
522:                                         void (***f1)(PetscInt, PetscInt, PetscInt,
523:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
524:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
525:                                                      PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
526: {

530:   PetscWeakFormGetFunction_Private(wf, wf->bdf0, label, val, f, n0, (void (***)(void)) f0);
531:   PetscWeakFormGetFunction_Private(wf, wf->bdf1, label, val, f, n1, (void (***)(void)) f1);
532:   return(0);
533: }

535: PetscErrorCode PetscWeakFormAddBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f,
536:                                           void (*f0)(PetscInt, PetscInt, PetscInt,
537:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
538:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
539:                                                      PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
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: {

548:   PetscWeakFormAddFunction_Private(wf, wf->bdf0, label, val, f, (void (*)(void)) f0);
549:   PetscWeakFormAddFunction_Private(wf, wf->bdf1, label, val, f, (void (*)(void)) f1);
550:   return(0);
551: }

553: PetscErrorCode PetscWeakFormSetBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f,
554:                                           PetscInt n0,
555:                                           void (**f0)(PetscInt, PetscInt, PetscInt,
556:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
557:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
558:                                                      PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
559:                                           PetscInt n1,
560:                                           void (**f1)(PetscInt, PetscInt, PetscInt,
561:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
562:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
563:                                                      PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
564: {

568:   PetscWeakFormSetFunction_Private(wf, wf->bdf0, label, val, f, n0, (void (**)(void)) f0);
569:   PetscWeakFormSetFunction_Private(wf, wf->bdf1, label, val, f, n1, (void (**)(void)) f1);
570:   return(0);
571: }

573: PetscErrorCode PetscWeakFormSetIndexBdResidual(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f,
574:                                           PetscInt i0,
575:                                           void (*f0)(PetscInt, PetscInt, PetscInt,
576:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
577:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
578:                                                      PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
579:                                           PetscInt i1,
580:                                           void (*f1)(PetscInt, PetscInt, PetscInt,
581:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
582:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
583:                                                      PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
584: {

588:   PetscWeakFormSetIndexFunction_Private(wf, wf->bdf0, label, val, f, i0, (void (*)(void)) f0);
589:   PetscWeakFormSetIndexFunction_Private(wf, wf->bdf1, label, val, f, i1, (void (*)(void)) f1);
590:   return(0);
591: }

593: PetscErrorCode PetscWeakFormHasJacobian(PetscWeakForm wf, PetscBool *hasJac)
594: {
595:   PetscInt       n0, n1, n2, n3;

601:   PetscHMapFormGetSize(wf->g0, &n0);
602:   PetscHMapFormGetSize(wf->g1, &n1);
603:   PetscHMapFormGetSize(wf->g2, &n2);
604:   PetscHMapFormGetSize(wf->g3, &n3);
605:   *hasJac = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE;
606:   return(0);
607: }

609: PetscErrorCode PetscWeakFormGetJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
610:                                         PetscInt *n0,
611:                                         void (***g0)(PetscInt, PetscInt, PetscInt,
612:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
613:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
614:                                                      PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
615:                                         PetscInt *n1,
616:                                         void (***g1)(PetscInt, PetscInt, PetscInt,
617:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
618:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
619:                                                      PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
620:                                         PetscInt *n2,
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:                                         PetscInt *n3,
626:                                         void (***g3)(PetscInt, PetscInt, PetscInt,
627:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
628:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
629:                                                      PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
630: {
631:   PetscInt       find = f*wf->Nf + g;

635:   PetscWeakFormGetFunction_Private(wf, wf->g0, label, val, find, n0, (void (***)(void)) g0);
636:   PetscWeakFormGetFunction_Private(wf, wf->g1, label, val, find, n1, (void (***)(void)) g1);
637:   PetscWeakFormGetFunction_Private(wf, wf->g2, label, val, find, n2, (void (***)(void)) g2);
638:   PetscWeakFormGetFunction_Private(wf, wf->g3, label, val, find, n3, (void (***)(void)) g3);
639:   return(0);
640: }

642: PetscErrorCode PetscWeakFormAddJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
643:                                         void (*g0)(PetscInt, PetscInt, PetscInt,
644:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
645:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
646:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
647:                                         void (*g1)(PetscInt, PetscInt, PetscInt,
648:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
649:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
650:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
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:                                         void (*g3)(PetscInt, PetscInt, PetscInt,
656:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
657:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
658:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
659: {
660:   PetscInt       find = f*wf->Nf + g;

664:   PetscWeakFormAddFunction_Private(wf, wf->g0, label, val, find, (void (*)(void)) g0);
665:   PetscWeakFormAddFunction_Private(wf, wf->g1, label, val, find, (void (*)(void)) g1);
666:   PetscWeakFormAddFunction_Private(wf, wf->g2, label, val, find, (void (*)(void)) g2);
667:   PetscWeakFormAddFunction_Private(wf, wf->g3, label, val, find, (void (*)(void)) g3);
668:   return(0);
669: }

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

697:   PetscWeakFormSetFunction_Private(wf, wf->g0, label, val, find, n0, (void (**)(void)) g0);
698:   PetscWeakFormSetFunction_Private(wf, wf->g1, label, val, find, n1, (void (**)(void)) g1);
699:   PetscWeakFormSetFunction_Private(wf, wf->g2, label, val, find, n2, (void (**)(void)) g2);
700:   PetscWeakFormSetFunction_Private(wf, wf->g3, label, val, find, n3, (void (**)(void)) g3);
701:   return(0);
702: }

704: PetscErrorCode PetscWeakFormSetIndexJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
705:                                         PetscInt i0,
706:                                         void (*g0)(PetscInt, PetscInt, PetscInt,
707:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
708:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
709:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
710:                                         PetscInt i1,
711:                                         void (*g1)(PetscInt, PetscInt, PetscInt,
712:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
713:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
714:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
715:                                         PetscInt i2,
716:                                         void (*g2)(PetscInt, PetscInt, PetscInt,
717:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
718:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
719:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
720:                                         PetscInt i3,
721:                                         void (*g3)(PetscInt, PetscInt, PetscInt,
722:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
723:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
724:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
725: {
726:   PetscInt       find = f*wf->Nf + g;

730:   PetscWeakFormSetIndexFunction_Private(wf, wf->g0, label, val, find, i0, (void (*)(void)) g0);
731:   PetscWeakFormSetIndexFunction_Private(wf, wf->g1, label, val, find, i1, (void (*)(void)) g1);
732:   PetscWeakFormSetIndexFunction_Private(wf, wf->g2, label, val, find, i2, (void (*)(void)) g2);
733:   PetscWeakFormSetIndexFunction_Private(wf, wf->g3, label, val, find, i3, (void (*)(void)) g3);
734:   return(0);
735: }

737: PetscErrorCode PetscWeakFormHasJacobianPreconditioner(PetscWeakForm wf, PetscBool *hasJacPre)
738: {
739:   PetscInt       n0, n1, n2, n3;

745:   PetscHMapFormGetSize(wf->gp0, &n0);
746:   PetscHMapFormGetSize(wf->gp1, &n1);
747:   PetscHMapFormGetSize(wf->gp2, &n2);
748:   PetscHMapFormGetSize(wf->gp3, &n3);
749:   *hasJacPre = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE;
750:   return(0);
751: }

753: PetscErrorCode PetscWeakFormGetJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
754:                                                       PetscInt *n0,
755:                                                       void (***g0)(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:                                                       PetscInt *n1,
760:                                                       void (***g1)(PetscInt, PetscInt, PetscInt,
761:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
762:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
763:                                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
764:                                                       PetscInt *n2,
765:                                                       void (***g2)(PetscInt, PetscInt, PetscInt,
766:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
767:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
768:                                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
769:                                                       PetscInt *n3,
770:                                                       void (***g3)(PetscInt, PetscInt, PetscInt,
771:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
772:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
773:                                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
774: {
775:   PetscInt       find = f*wf->Nf + g;

779:   PetscWeakFormGetFunction_Private(wf, wf->gp0, label, val, find, n0, (void (***)(void)) g0);
780:   PetscWeakFormGetFunction_Private(wf, wf->gp1, label, val, find, n1, (void (***)(void)) g1);
781:   PetscWeakFormGetFunction_Private(wf, wf->gp2, label, val, find, n2, (void (***)(void)) g2);
782:   PetscWeakFormGetFunction_Private(wf, wf->gp3, label, val, find, n3, (void (***)(void)) g3);
783:   return(0);
784: }

786: PetscErrorCode PetscWeakFormAddJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
787:                                         void (*g0)(PetscInt, PetscInt, PetscInt,
788:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
789:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
790:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
791:                                         void (*g1)(PetscInt, PetscInt, PetscInt,
792:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
793:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
794:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
795:                                         void (*g2)(PetscInt, PetscInt, PetscInt,
796:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
797:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
798:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
799:                                         void (*g3)(PetscInt, PetscInt, PetscInt,
800:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
801:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
802:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
803: {
804:   PetscInt       find = f*wf->Nf + g;

808:   PetscWeakFormAddFunction_Private(wf, wf->gp0, label, val, find, (void (*)(void)) g0);
809:   PetscWeakFormAddFunction_Private(wf, wf->gp1, label, val, find, (void (*)(void)) g1);
810:   PetscWeakFormAddFunction_Private(wf, wf->gp2, label, val, find, (void (*)(void)) g2);
811:   PetscWeakFormAddFunction_Private(wf, wf->gp3, label, val, find, (void (*)(void)) g3);
812:   return(0);
813: }

815: PetscErrorCode PetscWeakFormSetJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
816:                                                       PetscInt n0,
817:                                                       void (**g0)(PetscInt, PetscInt, PetscInt,
818:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
819:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
820:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
821:                                                       PetscInt n1,
822:                                                       void (**g1)(PetscInt, PetscInt, PetscInt,
823:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
824:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
825:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
826:                                                       PetscInt n2,
827:                                                       void (**g2)(PetscInt, PetscInt, PetscInt,
828:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
829:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
830:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
831:                                                       PetscInt n3,
832:                                                       void (**g3)(PetscInt, PetscInt, PetscInt,
833:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
834:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
835:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
836: {
837:   PetscInt       find = f*wf->Nf + g;

841:   PetscWeakFormSetFunction_Private(wf, wf->gp0, label, val, find, n0, (void (**)(void)) g0);
842:   PetscWeakFormSetFunction_Private(wf, wf->gp1, label, val, find, n1, (void (**)(void)) g1);
843:   PetscWeakFormSetFunction_Private(wf, wf->gp2, label, val, find, n2, (void (**)(void)) g2);
844:   PetscWeakFormSetFunction_Private(wf, wf->gp3, label, val, find, n3, (void (**)(void)) g3);
845:   return(0);
846: }

848: PetscErrorCode PetscWeakFormSetIndexJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
849:                                                       PetscInt i0,
850:                                                       void (*g0)(PetscInt, PetscInt, PetscInt,
851:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
852:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
853:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
854:                                                       PetscInt i1,
855:                                                       void (*g1)(PetscInt, PetscInt, PetscInt,
856:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
857:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
858:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
859:                                                       PetscInt i2,
860:                                                       void (*g2)(PetscInt, PetscInt, PetscInt,
861:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
862:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
863:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
864:                                                       PetscInt i3,
865:                                                       void (*g3)(PetscInt, PetscInt, PetscInt,
866:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
867:                                                                  const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
868:                                                                  PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
869: {
870:   PetscInt       find = f*wf->Nf + g;

874:   PetscWeakFormSetIndexFunction_Private(wf, wf->gp0, label, val, find, i0, (void (*)(void)) g0);
875:   PetscWeakFormSetIndexFunction_Private(wf, wf->gp1, label, val, find, i1, (void (*)(void)) g1);
876:   PetscWeakFormSetIndexFunction_Private(wf, wf->gp2, label, val, find, i2, (void (*)(void)) g2);
877:   PetscWeakFormSetIndexFunction_Private(wf, wf->gp3, label, val, find, i3, (void (*)(void)) g3);
878:   return(0);
879: }

881: PetscErrorCode PetscWeakFormHasBdJacobian(PetscWeakForm wf, PetscBool *hasJac)
882: {
883:   PetscInt       n0, n1, n2, n3;

889:   PetscHMapFormGetSize(wf->bdg0, &n0);
890:   PetscHMapFormGetSize(wf->bdg1, &n1);
891:   PetscHMapFormGetSize(wf->bdg2, &n2);
892:   PetscHMapFormGetSize(wf->bdg3, &n3);
893:   *hasJac = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE;
894:   return(0);
895: }

897: PetscErrorCode PetscWeakFormGetBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
898:                                           PetscInt *n0,
899:                                           void (***g0)(PetscInt, PetscInt, PetscInt,
900:                                                        const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
901:                                                        const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
902:                                                        PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
903:                                           PetscInt *n1,
904:                                           void (***g1)(PetscInt, PetscInt, PetscInt,
905:                                                        const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
906:                                                        const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
907:                                                        PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
908:                                           PetscInt *n2,
909:                                           void (***g2)(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 *n3,
914:                                           void (***g3)(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: {
919:   PetscInt       find = f*wf->Nf + g;

923:   PetscWeakFormGetFunction_Private(wf, wf->bdg0, label, val, find, n0, (void (***)(void)) g0);
924:   PetscWeakFormGetFunction_Private(wf, wf->bdg1, label, val, find, n1, (void (***)(void)) g1);
925:   PetscWeakFormGetFunction_Private(wf, wf->bdg2, label, val, find, n2, (void (***)(void)) g2);
926:   PetscWeakFormGetFunction_Private(wf, wf->bdg3, label, val, find, n3, (void (***)(void)) g3);
927:   return(0);
928: }

930: PetscErrorCode PetscWeakFormAddBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
931:                                           void (*g0)(PetscInt, PetscInt, PetscInt,
932:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
933:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
934:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
935:                                           void (*g1)(PetscInt, PetscInt, PetscInt,
936:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
937:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
938:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
939:                                           void (*g2)(PetscInt, PetscInt, PetscInt,
940:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
941:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
942:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
943:                                           void (*g3)(PetscInt, PetscInt, PetscInt,
944:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
945:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
946:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
947: {
948:   PetscInt       find = f*wf->Nf + g;

952:   PetscWeakFormAddFunction_Private(wf, wf->bdg0, label, val, find, (void (*)(void)) g0);
953:   PetscWeakFormAddFunction_Private(wf, wf->bdg1, label, val, find, (void (*)(void)) g1);
954:   PetscWeakFormAddFunction_Private(wf, wf->bdg2, label, val, find, (void (*)(void)) g2);
955:   PetscWeakFormAddFunction_Private(wf, wf->bdg3, label, val, find, (void (*)(void)) g3);
956:   return(0);
957: }

959: PetscErrorCode PetscWeakFormSetBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
960:                                           PetscInt n0,
961:                                           void (**g0)(PetscInt, PetscInt, PetscInt,
962:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
963:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
964:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
965:                                           PetscInt n1,
966:                                           void (**g1)(PetscInt, PetscInt, PetscInt,
967:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
968:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
969:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
970:                                           PetscInt n2,
971:                                           void (**g2)(PetscInt, PetscInt, PetscInt,
972:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
973:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
974:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
975:                                           PetscInt n3,
976:                                           void (**g3)(PetscInt, PetscInt, PetscInt,
977:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
978:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
979:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
980: {
981:   PetscInt       find = f*wf->Nf + g;

985:   PetscWeakFormSetFunction_Private(wf, wf->bdg0, label, val, find, n0, (void (**)(void)) g0);
986:   PetscWeakFormSetFunction_Private(wf, wf->bdg1, label, val, find, n1, (void (**)(void)) g1);
987:   PetscWeakFormSetFunction_Private(wf, wf->bdg2, label, val, find, n2, (void (**)(void)) g2);
988:   PetscWeakFormSetFunction_Private(wf, wf->bdg3, label, val, find, n3, (void (**)(void)) g3);
989:   return(0);
990: }

992: PetscErrorCode PetscWeakFormSetIndexBdJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
993:                                           PetscInt i0,
994:                                           void (*g0)(PetscInt, PetscInt, PetscInt,
995:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
996:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
997:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
998:                                           PetscInt i1,
999:                                           void (*g1)(PetscInt, PetscInt, PetscInt,
1000:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1001:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1002:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1003:                                           PetscInt i2,
1004:                                           void (*g2)(PetscInt, PetscInt, PetscInt,
1005:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1006:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1007:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1008:                                           PetscInt i3,
1009:                                           void (*g3)(PetscInt, PetscInt, PetscInt,
1010:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1011:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1012:                                                      PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1013: {
1014:   PetscInt       find = f*wf->Nf + g;

1018:   PetscWeakFormSetIndexFunction_Private(wf, wf->bdg0, label, val, find, i0, (void (*)(void)) g0);
1019:   PetscWeakFormSetIndexFunction_Private(wf, wf->bdg1, label, val, find, i1, (void (*)(void)) g1);
1020:   PetscWeakFormSetIndexFunction_Private(wf, wf->bdg2, label, val, find, i2, (void (*)(void)) g2);
1021:   PetscWeakFormSetIndexFunction_Private(wf, wf->bdg3, label, val, find, i3, (void (*)(void)) g3);
1022:   return(0);
1023: }

1025: PetscErrorCode PetscWeakFormHasBdJacobianPreconditioner(PetscWeakForm wf, PetscBool *hasJacPre)
1026: {
1027:   PetscInt       n0, n1, n2, n3;

1033:   PetscHMapFormGetSize(wf->bdgp0, &n0);
1034:   PetscHMapFormGetSize(wf->bdgp1, &n1);
1035:   PetscHMapFormGetSize(wf->bdgp2, &n2);
1036:   PetscHMapFormGetSize(wf->bdgp3, &n3);
1037:   *hasJacPre = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE;
1038:   return(0);
1039: }

1041: PetscErrorCode PetscWeakFormGetBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
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;

1067:   PetscWeakFormGetFunction_Private(wf, wf->bdgp0, label, val, find, n0, (void (***)(void)) g0);
1068:   PetscWeakFormGetFunction_Private(wf, wf->bdgp1, label, val, find, n1, (void (***)(void)) g1);
1069:   PetscWeakFormGetFunction_Private(wf, wf->bdgp2, label, val, find, n2, (void (***)(void)) g2);
1070:   PetscWeakFormGetFunction_Private(wf, wf->bdgp3, label, val, find, n3, (void (***)(void)) g3);
1071:   return(0);
1072: }

1074: PetscErrorCode PetscWeakFormAddBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
1075:                                                         void (*g0)(PetscInt, PetscInt, PetscInt,
1076:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1077:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1078:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
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:                                                         void (*g2)(PetscInt, PetscInt, PetscInt,
1084:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1085:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1086:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1087:                                                         void (*g3)(PetscInt, PetscInt, PetscInt,
1088:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1089:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1090:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1091: {
1092:   PetscInt       find = f*wf->Nf + g;

1096:   PetscWeakFormAddFunction_Private(wf, wf->bdgp0, label, val, find, (void (*)(void)) g0);
1097:   PetscWeakFormAddFunction_Private(wf, wf->bdgp1, label, val, find, (void (*)(void)) g1);
1098:   PetscWeakFormAddFunction_Private(wf, wf->bdgp2, label, val, find, (void (*)(void)) g2);
1099:   PetscWeakFormAddFunction_Private(wf, wf->bdgp3, label, val, find, (void (*)(void)) g3);
1100:   return(0);
1101: }

1103: PetscErrorCode PetscWeakFormSetBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
1104:                                                         PetscInt n0,
1105:                                                         void (**g0)(PetscInt, PetscInt, PetscInt,
1106:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1107:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1108:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1109:                                                         PetscInt n1,
1110:                                                         void (**g1)(PetscInt, PetscInt, PetscInt,
1111:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1112:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1113:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1114:                                                         PetscInt n2,
1115:                                                         void (**g2)(PetscInt, PetscInt, PetscInt,
1116:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1117:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1118:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1119:                                                         PetscInt n3,
1120:                                                         void (**g3)(PetscInt, PetscInt, PetscInt,
1121:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1122:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1123:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1124: {
1125:   PetscInt       find = f*wf->Nf + g;

1129:   PetscWeakFormSetFunction_Private(wf, wf->bdgp0, label, val, find, n0, (void (**)(void)) g0);
1130:   PetscWeakFormSetFunction_Private(wf, wf->bdgp1, label, val, find, n1, (void (**)(void)) g1);
1131:   PetscWeakFormSetFunction_Private(wf, wf->bdgp2, label, val, find, n2, (void (**)(void)) g2);
1132:   PetscWeakFormSetFunction_Private(wf, wf->bdgp3, label, val, find, n3, (void (**)(void)) g3);
1133:   return(0);
1134: }

1136: PetscErrorCode PetscWeakFormSetIndexBdJacobianPreconditioner(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
1137:                                                         PetscInt i0,
1138:                                                         void (*g0)(PetscInt, PetscInt, PetscInt,
1139:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1140:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1141:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1142:                                                         PetscInt i1,
1143:                                                         void (*g1)(PetscInt, PetscInt, PetscInt,
1144:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1145:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1146:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1147:                                                         PetscInt i2,
1148:                                                         void (*g2)(PetscInt, PetscInt, PetscInt,
1149:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1150:                                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1151:                                                                    PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1152:                                                         PetscInt i3,
1153:                                                         void (*g3)(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[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1157: {
1158:   PetscInt       find = f*wf->Nf + g;

1162:   PetscWeakFormSetIndexFunction_Private(wf, wf->bdgp0, label, val, find, i0, (void (*)(void)) g0);
1163:   PetscWeakFormSetIndexFunction_Private(wf, wf->bdgp1, label, val, find, i1, (void (*)(void)) g1);
1164:   PetscWeakFormSetIndexFunction_Private(wf, wf->bdgp2, label, val, find, i2, (void (*)(void)) g2);
1165:   PetscWeakFormSetIndexFunction_Private(wf, wf->bdgp3, label, val, find, i3, (void (*)(void)) g3);
1166:   return(0);
1167: }

1169: PetscErrorCode PetscWeakFormHasDynamicJacobian(PetscWeakForm wf, PetscBool *hasDynJac)
1170: {
1171:   PetscInt       n0, n1, n2, n3;

1177:   PetscHMapFormGetSize(wf->gt0, &n0);
1178:   PetscHMapFormGetSize(wf->gt1, &n1);
1179:   PetscHMapFormGetSize(wf->gt2, &n2);
1180:   PetscHMapFormGetSize(wf->gt3, &n3);
1181:   *hasDynJac = n0+n1+n2+n3 ? PETSC_TRUE : PETSC_FALSE;
1182:   return(0);
1183: }

1185: PetscErrorCode PetscWeakFormGetDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
1186:                                         PetscInt *n0,
1187:                                         void (***g0)(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 *n1,
1192:                                         void (***g1)(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:                                         PetscInt *n2,
1197:                                         void (***g2)(PetscInt, PetscInt, PetscInt,
1198:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1199:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1200:                                                      PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1201:                                         PetscInt *n3,
1202:                                         void (***g3)(PetscInt, PetscInt, PetscInt,
1203:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1204:                                                      const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1205:                                                      PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1206: {
1207:   PetscInt       find = f*wf->Nf + g;

1211:   PetscWeakFormGetFunction_Private(wf, wf->gt0, label, val, find, n0, (void (***)(void)) g0);
1212:   PetscWeakFormGetFunction_Private(wf, wf->gt1, label, val, find, n1, (void (***)(void)) g1);
1213:   PetscWeakFormGetFunction_Private(wf, wf->gt2, label, val, find, n2, (void (***)(void)) g2);
1214:   PetscWeakFormGetFunction_Private(wf, wf->gt3, label, val, find, n3, (void (***)(void)) g3);
1215:   return(0);
1216: }

1218: PetscErrorCode PetscWeakFormAddDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
1219:                                         void (*g0)(PetscInt, PetscInt, PetscInt,
1220:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1221:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1222:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1223:                                         void (*g1)(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:                                         void (*g2)(PetscInt, PetscInt, PetscInt,
1228:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1229:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1230:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1231:                                         void (*g3)(PetscInt, PetscInt, PetscInt,
1232:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1233:                                                    const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1234:                                                    PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1235: {
1236:   PetscInt       find = f*wf->Nf + g;

1240:   PetscWeakFormAddFunction_Private(wf, wf->gt0, label, val, find, (void (*)(void)) g0);
1241:   PetscWeakFormAddFunction_Private(wf, wf->gt1, label, val, find, (void (*)(void)) g1);
1242:   PetscWeakFormAddFunction_Private(wf, wf->gt2, label, val, find, (void (*)(void)) g2);
1243:   PetscWeakFormAddFunction_Private(wf, wf->gt3, label, val, find, (void (*)(void)) g3);
1244:   return(0);
1245: }

1247: PetscErrorCode PetscWeakFormSetDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
1248:                                                PetscInt n0,
1249:                                                void (**g0)(PetscInt, PetscInt, PetscInt,
1250:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1251:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1252:                                                           PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1253:                                                PetscInt n1,
1254:                                                void (**g1)(PetscInt, PetscInt, PetscInt,
1255:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1256:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1257:                                                           PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1258:                                                PetscInt n2,
1259:                                                void (**g2)(PetscInt, PetscInt, PetscInt,
1260:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1261:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1262:                                                           PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1263:                                                PetscInt n3,
1264:                                                void (**g3)(PetscInt, PetscInt, PetscInt,
1265:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1266:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1267:                                                           PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1268: {
1269:   PetscInt       find = f*wf->Nf + g;

1273:   PetscWeakFormSetFunction_Private(wf, wf->gt0, label, val, find, n0, (void (**)(void)) g0);
1274:   PetscWeakFormSetFunction_Private(wf, wf->gt1, label, val, find, n1, (void (**)(void)) g1);
1275:   PetscWeakFormSetFunction_Private(wf, wf->gt2, label, val, find, n2, (void (**)(void)) g2);
1276:   PetscWeakFormSetFunction_Private(wf, wf->gt3, label, val, find, n3, (void (**)(void)) g3);
1277:   return(0);
1278: }

1280: PetscErrorCode PetscWeakFormSetIndexDynamicJacobian(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt g,
1281:                                                PetscInt i0,
1282:                                                void (*g0)(PetscInt, PetscInt, PetscInt,
1283:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1284:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1285:                                                           PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1286:                                                PetscInt i1,
1287:                                                void (*g1)(PetscInt, PetscInt, PetscInt,
1288:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1289:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1290:                                                           PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1291:                                                PetscInt i2,
1292:                                                void (*g2)(PetscInt, PetscInt, PetscInt,
1293:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1294:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1295:                                                           PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
1296:                                                PetscInt i3,
1297:                                                void (*g3)(PetscInt, PetscInt, PetscInt,
1298:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1299:                                                           const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
1300:                                                           PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1301: {
1302:   PetscInt       find = f*wf->Nf + g;

1306:   PetscWeakFormSetIndexFunction_Private(wf, wf->gt0, label, val, find, i0, (void (*)(void)) g0);
1307:   PetscWeakFormSetIndexFunction_Private(wf, wf->gt1, label, val, find, i1, (void (*)(void)) g1);
1308:   PetscWeakFormSetIndexFunction_Private(wf, wf->gt2, label, val, find, i2, (void (*)(void)) g2);
1309:   PetscWeakFormSetIndexFunction_Private(wf, wf->gt3, label, val, find, i3, (void (*)(void)) g3);
1310:   return(0);
1311: }

1313: PetscErrorCode PetscWeakFormGetRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt *n,
1314:                                              void (***r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *))
1315: {

1319:   PetscWeakFormGetFunction_Private(wf, wf->r, label, val, f, n, (void (***)(void)) r);
1320:   return(0);
1321: }

1323: PetscErrorCode PetscWeakFormSetRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f,
1324:                                              PetscInt n,
1325:                                              void (**r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *))
1326: {

1330:   PetscWeakFormSetFunction_Private(wf, wf->r, label, val, f, n, (void (**)(void)) r);
1331:   return(0);
1332: }

1334: PetscErrorCode PetscWeakFormSetIndexRiemannSolver(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f,
1335:                                                   PetscInt i,
1336:                                                   void (*r)(PetscInt, PetscInt, const PetscReal[], const PetscReal[], const PetscScalar[], const PetscScalar[], PetscInt, const PetscScalar[], PetscScalar[], void *))
1337: {

1341:   PetscWeakFormSetIndexFunction_Private(wf, wf->r, label, val, f, i, (void (*)(void)) r);
1342:   return(0);
1343: }

1345: /*@
1346:   PetscWeakFormGetNumFields - Returns the number of fields

1348:   Not collective

1350:   Input Parameter:
1351: . wf - The PetscWeakForm object

1353:   Output Parameter:
1354: . Nf - The nubmer of fields

1356:   Level: beginner

1358: .seealso: PetscWeakFormSetNumFields(), PetscWeakFormCreate()
1359: @*/
1360: PetscErrorCode PetscWeakFormGetNumFields(PetscWeakForm wf, PetscInt *Nf)
1361: {
1365:   *Nf = wf->Nf;
1366:   return(0);
1367: }

1369: /*@
1370:   PetscWeakFormSetNumFields - Sets the number of fields

1372:   Not collective

1374:   Input Parameters:
1375: + wf - The PetscWeakForm object
1376: - Nf - The number of fields

1378:   Level: beginner

1380: .seealso: PetscWeakFormGetNumFields(), PetscWeakFormCreate()
1381: @*/
1382: PetscErrorCode PetscWeakFormSetNumFields(PetscWeakForm wf, PetscInt Nf)
1383: {
1386:   wf->Nf = Nf;
1387:   return(0);
1388: }

1390: /*@
1391:   PetscWeakFormDestroy - Destroys a PetscWeakForm object

1393:   Collective on wf

1395:   Input Parameter:
1396: . wf - the PetscWeakForm object to destroy

1398:   Level: developer

1400: .seealso PetscWeakFormCreate(), PetscWeakFormView()
1401: @*/
1402: PetscErrorCode PetscWeakFormDestroy(PetscWeakForm *wf)
1403: {

1407:   if (!*wf) return(0);

1410:   if (--((PetscObject)(*wf))->refct > 0) {*wf = NULL; return(0);}
1411:   ((PetscObject) (*wf))->refct = 0;
1412:   PetscChunkBufferDestroy(&(*wf)->funcs);
1413:   PetscHMapFormDestroy(&(*wf)->obj);
1414:   PetscHMapFormDestroy(&(*wf)->f0);
1415:   PetscHMapFormDestroy(&(*wf)->f1);
1416:   PetscHMapFormDestroy(&(*wf)->g0);
1417:   PetscHMapFormDestroy(&(*wf)->g1);
1418:   PetscHMapFormDestroy(&(*wf)->g2);
1419:   PetscHMapFormDestroy(&(*wf)->g3);
1420:   PetscHMapFormDestroy(&(*wf)->gp0);
1421:   PetscHMapFormDestroy(&(*wf)->gp1);
1422:   PetscHMapFormDestroy(&(*wf)->gp2);
1423:   PetscHMapFormDestroy(&(*wf)->gp3);
1424:   PetscHMapFormDestroy(&(*wf)->gt0);
1425:   PetscHMapFormDestroy(&(*wf)->gt1);
1426:   PetscHMapFormDestroy(&(*wf)->gt2);
1427:   PetscHMapFormDestroy(&(*wf)->gt3);
1428:   PetscHMapFormDestroy(&(*wf)->bdf0);
1429:   PetscHMapFormDestroy(&(*wf)->bdf1);
1430:   PetscHMapFormDestroy(&(*wf)->bdg0);
1431:   PetscHMapFormDestroy(&(*wf)->bdg1);
1432:   PetscHMapFormDestroy(&(*wf)->bdg2);
1433:   PetscHMapFormDestroy(&(*wf)->bdg3);
1434:   PetscHMapFormDestroy(&(*wf)->bdgp0);
1435:   PetscHMapFormDestroy(&(*wf)->bdgp1);
1436:   PetscHMapFormDestroy(&(*wf)->bdgp2);
1437:   PetscHMapFormDestroy(&(*wf)->bdgp3);
1438:   PetscHMapFormDestroy(&(*wf)->r);
1439:   PetscHeaderDestroy(wf);
1440:   return(0);
1441: }

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

1449:   PetscHMapFormGetSize(map, &Nk);
1450:   if (Nk) {
1451:     PetscHashFormKey *keys;
1452:     void           (**funcs)(void);
1453:     const char       *name;
1454:     PetscInt          off = 0, n, i;

1456:     PetscMalloc1(Nk, &keys);
1457:     PetscHMapFormGetKeys(map, &off, keys);
1458:     PetscViewerASCIIPrintf(viewer, "%s\n", tableName);
1459:     PetscViewerASCIIPushTab(viewer);
1460:     for (k = 0; k < Nk; ++k) {
1461:       if (keys[k].label) {
1462:         PetscObjectGetName((PetscObject) keys[k].label, &name);
1463:         PetscViewerASCIIPrintf(viewer, "(%s, %D) ", name, keys[k].value);
1464:       } else {PetscViewerASCIIPrintf(viewer, "");}
1465:       PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
1466:       if (splitField) {PetscViewerASCIIPrintf(viewer, "(%D, %D) ", keys[k].field/Nf, keys[k].field%Nf);}
1467:       else            {PetscViewerASCIIPrintf(viewer, "(%D) ", keys[k].field);}
1468:       PetscWeakFormGetFunction_Private(wf, map, keys[k].label, keys[k].value, keys[k].field, &n, &funcs);
1469:       for (i = 0; i < n; ++i) {
1470:         if (i > 0) {PetscViewerASCIIPrintf(viewer, ", ");}
1471:         PetscDLAddr(funcs[i], &name);
1472:         if (name) {PetscViewerASCIIPrintf(viewer, "%s", name);}
1473:         else      {PetscViewerASCIIPrintf(viewer, "%p", funcs[i]);}
1474:       }
1475:       PetscViewerASCIIPrintf(viewer, "\n");
1476:       PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
1477:     }
1478:     PetscViewerASCIIPopTab(viewer);
1479:     PetscFree(keys);
1480:   }
1481:   return(0);
1482: }

1484: static PetscErrorCode PetscWeakFormView_Ascii(PetscWeakForm wf, PetscViewer viewer)
1485: {
1486:   PetscViewerFormat format;
1487:   PetscErrorCode    ierr;

1490:   PetscViewerGetFormat(viewer, &format);
1491:   PetscViewerASCIIPrintf(viewer, "Weak Form System with %d fields\n", wf->Nf);
1492:   PetscViewerASCIIPushTab(viewer);
1493:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_FALSE, "Objective", wf->obj);
1494:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_FALSE, "Residual f0", wf->f0);
1495:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_FALSE, "Residual f1", wf->f1);
1496:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_FALSE, "Boundary Residual f0", wf->bdf0);
1497:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_FALSE, "Boundary Residual f1", wf->bdf1);
1498:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Jacobian g0", wf->g0);
1499:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Jacobian g1", wf->g1);
1500:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Jacobian g2", wf->g2);
1501:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Jacobian g3", wf->g3);
1502:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Jacobian Preconditioner g0", wf->gp0);
1503:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Jacobian Preconditioner g1", wf->gp1);
1504:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Jacobian Preconditioner g2", wf->gp2);
1505:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Jacobian Preconditioner g3", wf->gp3);
1506:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Dynamic Jacobian g0", wf->gt0);
1507:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Dynamic Jacobian g1", wf->gt1);
1508:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Dynamic Jacobian g2", wf->gt2);
1509:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Dynamic Jacobian g3", wf->gt3);
1510:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Boundary Jacobian g0", wf->bdg0);
1511:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Boundary Jacobian g1", wf->bdg1);
1512:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Boundary Jacobian g2", wf->bdg2);
1513:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Boundary Jacobian g3", wf->bdg3);
1514:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Boundary Jacobian Preconditioner g0", wf->bdgp0);
1515:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Boundary Jacobian Preconditioner g1", wf->bdgp1);
1516:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Boundary Jacobian Preconditioner g2", wf->bdgp2);
1517:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE,  "Boundary Jacobian Preconditioner g3", wf->bdgp3);
1518:   PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_FALSE, "Riemann Solver", wf->r);
1519:   PetscViewerASCIIPopTab(viewer);
1520:   return(0);
1521: }

1523: /*@C
1524:   PetscWeakFormView - Views a PetscWeakForm

1526:   Collective on wf

1528:   Input Parameter:
1529: + wf - the PetscWeakForm object to view
1530: - v  - the viewer

1532:   Level: developer

1534: .seealso PetscWeakFormDestroy(), PetscWeakFormCreate()
1535: @*/
1536: PetscErrorCode PetscWeakFormView(PetscWeakForm wf, PetscViewer v)
1537: {
1538:   PetscBool      iascii;

1543:   if (!v) {PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject) wf), &v);}
1545:   PetscObjectTypeCompare((PetscObject) v, PETSCVIEWERASCII, &iascii);
1546:   if (iascii) {PetscWeakFormView_Ascii(wf, v);}
1547:   if (wf->ops->view) {(*wf->ops->view)(wf, v);}
1548:   return(0);
1549: }

1551: /*@
1552:   PetscWeakFormCreate - Creates an empty PetscWeakForm object.

1554:   Collective

1556:   Input Parameter:
1557: . comm - The communicator for the PetscWeakForm object

1559:   Output Parameter:
1560: . wf - The PetscWeakForm object

1562:   Level: beginner

1564: .seealso: PetscDS, PetscWeakFormDestroy()
1565: @*/
1566: PetscErrorCode PetscWeakFormCreate(MPI_Comm comm, PetscWeakForm *wf)
1567: {
1568:   PetscWeakForm  p;

1573:   *wf  = NULL;
1574:   PetscDSInitializePackage();

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

1578:   p->Nf = 0;
1579:   PetscChunkBufferCreate(sizeof(&PetscWeakFormCreate), 2, &p->funcs);
1580:   PetscHMapFormCreate(&p->obj);
1581:   PetscHMapFormCreate(&p->f0);
1582:   PetscHMapFormCreate(&p->f1);
1583:   PetscHMapFormCreate(&p->g0);
1584:   PetscHMapFormCreate(&p->g1);
1585:   PetscHMapFormCreate(&p->g2);
1586:   PetscHMapFormCreate(&p->g3);
1587:   PetscHMapFormCreate(&p->gp0);
1588:   PetscHMapFormCreate(&p->gp1);
1589:   PetscHMapFormCreate(&p->gp2);
1590:   PetscHMapFormCreate(&p->gp3);
1591:   PetscHMapFormCreate(&p->gt0);
1592:   PetscHMapFormCreate(&p->gt1);
1593:   PetscHMapFormCreate(&p->gt2);
1594:   PetscHMapFormCreate(&p->gt3);
1595:   PetscHMapFormCreate(&p->bdf0);
1596:   PetscHMapFormCreate(&p->bdf1);
1597:   PetscHMapFormCreate(&p->bdg0);
1598:   PetscHMapFormCreate(&p->bdg1);
1599:   PetscHMapFormCreate(&p->bdg2);
1600:   PetscHMapFormCreate(&p->bdg3);
1601:   PetscHMapFormCreate(&p->bdgp0);
1602:   PetscHMapFormCreate(&p->bdgp1);
1603:   PetscHMapFormCreate(&p->bdgp2);
1604:   PetscHMapFormCreate(&p->bdgp3);
1605:   PetscHMapFormCreate(&p->r);
1606:   *wf = p;
1607:   return(0);
1608: }