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