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;
127: PetscInt i;
129: PetscFunctionBegin;
130: key.label = label;
131: key.value = value;
132: key.field = f;
133: key.part = part;
134: if (!func) {
135: PetscCall(PetscHMapFormDel(ht, key));
136: PetscFunctionReturn(PETSC_SUCCESS);
137: } else PetscCall(PetscHMapFormGet(ht, key, &chunk));
138: if (chunk.size < 0) {
139: PetscCall(PetscChunkBufferCreateChunk(wf->funcs, n, &chunk));
140: PetscCall(PetscHMapFormSet(ht, key, chunk));
141: } else if (chunk.size <= n) {
142: PetscCall(PetscChunkBufferEnlargeChunk(wf->funcs, n - chunk.size, &chunk));
143: PetscCall(PetscHMapFormSet(ht, key, chunk));
144: }
145: for (i = 0; i < n; ++i) ((PetscVoidFn **)&wf->funcs->array[chunk.start])[i] = func[i];
146: PetscFunctionReturn(PETSC_SUCCESS);
147: }
149: static PetscErrorCode PetscWeakFormAddFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscVoidFn *func)
150: {
151: PetscFormKey key;
152: PetscChunk chunk;
154: PetscFunctionBegin;
155: if (!func) PetscFunctionReturn(PETSC_SUCCESS);
156: key.label = label;
157: key.value = value;
158: key.field = f;
159: key.part = part;
160: PetscCall(PetscHMapFormGet(ht, key, &chunk));
161: if (chunk.size < 0) {
162: PetscCall(PetscChunkBufferCreateChunk(wf->funcs, 1, &chunk));
163: PetscCall(PetscHMapFormSet(ht, key, chunk));
164: ((PetscVoidFn **)&wf->funcs->array[chunk.start])[0] = func;
165: } else {
166: PetscCall(PetscChunkBufferEnlargeChunk(wf->funcs, 1, &chunk));
167: PetscCall(PetscHMapFormSet(ht, key, chunk));
168: ((PetscVoidFn **)&wf->funcs->array[chunk.start])[chunk.size - 1] = func;
169: }
170: PetscFunctionReturn(PETSC_SUCCESS);
171: }
173: static PetscErrorCode PetscWeakFormGetIndexFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt ind, PetscVoidFn **func)
174: {
175: PetscFormKey key;
176: PetscChunk chunk;
178: PetscFunctionBegin;
179: key.label = label;
180: key.value = value;
181: key.field = f;
182: key.part = part;
183: PetscCall(PetscHMapFormGet(ht, key, &chunk));
184: if (chunk.size < 0) {
185: *func = NULL;
186: } else {
187: PetscCheck(ind < chunk.size, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Index %" PetscInt_FMT " not in [0, %" PetscCount_FMT ")", ind, chunk.size);
188: *func = ((PetscVoidFn **)&wf->funcs->array[chunk.start])[ind];
189: }
190: PetscFunctionReturn(PETSC_SUCCESS);
191: }
193: /* Ignore a NULL func */
194: static PetscErrorCode PetscWeakFormSetIndexFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt ind, PetscVoidFn *func)
195: {
196: PetscFormKey key;
197: PetscChunk chunk;
199: PetscFunctionBegin;
200: if (!func) PetscFunctionReturn(PETSC_SUCCESS);
201: key.label = label;
202: key.value = value;
203: key.field = f;
204: key.part = part;
205: PetscCall(PetscHMapFormGet(ht, key, &chunk));
206: if (chunk.size < 0) {
207: PetscCall(PetscChunkBufferCreateChunk(wf->funcs, ind + 1, &chunk));
208: PetscCall(PetscHMapFormSet(ht, key, chunk));
209: } else if (chunk.size <= ind) {
210: PetscCall(PetscChunkBufferEnlargeChunk(wf->funcs, ind - chunk.size + 1, &chunk));
211: PetscCall(PetscHMapFormSet(ht, key, chunk));
212: }
213: ((PetscVoidFn **)&wf->funcs->array[chunk.start])[ind] = func;
214: PetscFunctionReturn(PETSC_SUCCESS);
215: }
217: static PetscErrorCode PetscWeakFormClearIndexFunction_Private(PetscWeakForm wf, PetscHMapForm ht, DMLabel label, PetscInt value, PetscInt f, PetscInt part, PetscInt ind)
218: {
219: PetscFormKey key;
220: PetscChunk chunk;
222: PetscFunctionBegin;
223: key.label = label;
224: key.value = value;
225: key.field = f;
226: key.part = part;
227: PetscCall(PetscHMapFormGet(ht, key, &chunk));
228: if (chunk.size < 0) PetscFunctionReturn(PETSC_SUCCESS);
229: else if (!ind && chunk.size == 1) {
230: PetscCall(PetscHMapFormDel(ht, key));
231: PetscFunctionReturn(PETSC_SUCCESS);
232: } else if (chunk.size <= ind) PetscFunctionReturn(PETSC_SUCCESS);
233: ((PetscVoidFn **)&wf->funcs->array[chunk.start])[ind] = NULL;
234: PetscFunctionReturn(PETSC_SUCCESS);
235: }
237: /*@
238: PetscWeakFormCopy - Copy the pointwise functions to another `PetscWeakForm`
240: Not Collective
242: Input Parameter:
243: . wf - The original `PetscWeakForm`
245: Output Parameter:
246: . wfNew - The copy of the `PetscWeakForm`
248: Level: intermediate
250: .seealso: `PetscWeakForm`, `PetscWeakFormCreate()`, `PetscWeakFormDestroy()`
251: @*/
252: PetscErrorCode PetscWeakFormCopy(PetscWeakForm wf, PetscWeakForm wfNew)
253: {
254: PetscInt f;
256: PetscFunctionBegin;
257: wfNew->Nf = wf->Nf;
258: PetscCall(PetscChunkBufferDestroy(&wfNew->funcs));
259: PetscCall(PetscChunkBufferDuplicate(wf->funcs, &wfNew->funcs));
260: for (f = 0; f < PETSC_NUM_WF; ++f) {
261: PetscCall(PetscHMapFormDestroy(&wfNew->form[f]));
262: PetscCall(PetscHMapFormDuplicate(wf->form[f], &wfNew->form[f]));
263: }
264: PetscFunctionReturn(PETSC_SUCCESS);
265: }
267: /*@
268: PetscWeakFormClear - Clear all functions from the `PetscWeakForm`
270: Not Collective
272: Input Parameter:
273: . wf - The original `PetscWeakForm`
275: Level: intermediate
277: .seealso: `PetscWeakForm`, `PetscWeakFormCopy()`, `PetscWeakFormCreate()`, `PetscWeakFormDestroy()`
278: @*/
279: PetscErrorCode PetscWeakFormClear(PetscWeakForm wf)
280: {
281: PetscInt f;
283: PetscFunctionBegin;
284: for (f = 0; f < PETSC_NUM_WF; ++f) PetscCall(PetscHMapFormClear(wf->form[f]));
285: PetscFunctionReturn(PETSC_SUCCESS);
286: }
288: static PetscErrorCode PetscWeakFormRewriteKeys_Internal(PetscWeakForm wf, PetscHMapForm hmap, DMLabel label, PetscInt Nv, const PetscInt values[])
289: {
290: PetscFormKey *keys;
291: PetscVoidFn **tmpfuncs;
292: PetscInt n, off = 0, maxNf = 0;
294: PetscFunctionBegin;
295: PetscCall(PetscHMapFormGetSize(hmap, &n));
296: PetscCall(PetscMalloc1(n, &keys));
297: PetscCall(PetscHMapFormGetKeys(hmap, &off, keys));
298: // Need to make a copy since SetFunction() can invalidate the storage
299: for (PetscInt i = 0; i < n; ++i) {
300: if (keys[i].label == label) {
301: PetscVoidFn **funcs;
302: PetscInt Nf;
304: PetscCall(PetscWeakFormGetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &Nf, &funcs));
305: maxNf = PetscMax(maxNf, Nf);
306: }
307: }
308: PetscCall(PetscMalloc1(maxNf, &tmpfuncs));
309: for (PetscInt i = 0; i < n; ++i) {
310: if (keys[i].label == label) {
311: PetscBool clear = PETSC_TRUE;
312: PetscVoidFn **funcs;
313: PetscInt Nf;
315: PetscCall(PetscWeakFormGetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &Nf, &funcs));
316: for (PetscInt f = 0; f < Nf; ++f) tmpfuncs[f] = funcs[f];
317: for (PetscInt v = 0; v < Nv; ++v) {
318: PetscCall(PetscWeakFormSetFunction_Private(wf, hmap, keys[i].label, values[v], keys[i].field, keys[i].part, Nf, tmpfuncs));
319: if (values[v] == keys[i].value) clear = PETSC_FALSE;
320: }
321: if (clear) PetscCall(PetscWeakFormSetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, 0, NULL));
322: }
323: }
324: PetscCall(PetscFree(tmpfuncs));
325: PetscCall(PetscFree(keys));
326: PetscFunctionReturn(PETSC_SUCCESS);
327: }
329: /*@
330: PetscWeakFormRewriteKeys - Change any key on the given label to use the new set of label values
332: Not Collective
334: Input Parameters:
335: + wf - The original `PetscWeakForm`
336: . label - The label to change keys for
337: . Nv - The number of new label values
338: - values - The set of new values to relabel keys with
340: Level: intermediate
342: Note:
343: This is used internally when boundary label values are specified from the command line.
345: .seealso: `PetscWeakForm`, `DMLabel`, `PetscWeakFormReplaceLabel()`, `PetscWeakFormCreate()`, `PetscWeakFormDestroy()`
346: @*/
347: PetscErrorCode PetscWeakFormRewriteKeys(PetscWeakForm wf, DMLabel label, PetscInt Nv, const PetscInt values[])
348: {
349: PetscInt f;
351: PetscFunctionBegin;
352: for (f = 0; f < PETSC_NUM_WF; ++f) PetscCall(PetscWeakFormRewriteKeys_Internal(wf, wf->form[f], label, Nv, values));
353: PetscFunctionReturn(PETSC_SUCCESS);
354: }
356: static PetscErrorCode PetscWeakFormReplaceLabel_Internal(PetscWeakForm wf, PetscHMapForm hmap, DMLabel label)
357: {
358: PetscFormKey *keys;
359: PetscInt n, i, off = 0, maxFuncs = 0;
360: PetscVoidFn **tmpf;
361: const char *name = NULL;
363: PetscFunctionBegin;
364: if (label) PetscCall(PetscObjectGetName((PetscObject)label, &name));
365: PetscCall(PetscHMapFormGetSize(hmap, &n));
366: PetscCall(PetscMalloc1(n, &keys));
367: PetscCall(PetscHMapFormGetKeys(hmap, &off, keys));
368: for (i = 0; i < n; ++i) {
369: PetscBool match = PETSC_FALSE;
370: const char *lname = NULL;
372: if (label == keys[i].label) continue;
373: if (keys[i].label) PetscCall(PetscObjectGetName((PetscObject)keys[i].label, &lname));
374: PetscCall(PetscStrcmp(name, lname, &match));
375: if ((!name && !lname) || match) {
376: PetscVoidFn **funcs;
377: PetscInt Nf;
379: PetscCall(PetscWeakFormGetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &Nf, &funcs));
380: maxFuncs = PetscMax(maxFuncs, Nf);
381: }
382: }
383: /* Need temp space because chunk buffer can be reallocated in SetFunction() call */
384: PetscCall(PetscMalloc1(maxFuncs, &tmpf));
385: for (i = 0; i < n; ++i) {
386: PetscBool match = PETSC_FALSE;
387: const char *lname = NULL;
389: if (label == keys[i].label) continue;
390: if (keys[i].label) PetscCall(PetscObjectGetName((PetscObject)keys[i].label, &lname));
391: PetscCall(PetscStrcmp(name, lname, &match));
392: if ((!name && !lname) || match) {
393: PetscVoidFn **funcs;
394: PetscInt Nf, j;
396: PetscCall(PetscWeakFormGetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &Nf, &funcs));
397: for (j = 0; j < Nf; ++j) tmpf[j] = funcs[j];
398: PetscCall(PetscWeakFormSetFunction_Private(wf, hmap, label, keys[i].value, keys[i].field, keys[i].part, Nf, tmpf));
399: PetscCall(PetscWeakFormSetFunction_Private(wf, hmap, keys[i].label, keys[i].value, keys[i].field, keys[i].part, 0, NULL));
400: }
401: }
402: PetscCall(PetscFree(tmpf));
403: PetscCall(PetscFree(keys));
404: PetscFunctionReturn(PETSC_SUCCESS);
405: }
407: /*@
408: PetscWeakFormReplaceLabel - Change any key on a label of the same name to use the new label
410: Not Collective
412: Input Parameters:
413: + wf - The original `PetscWeakForm`
414: - label - The label to change keys for
416: Level: intermediate
418: Note:
419: This is used internally when meshes are modified
421: .seealso: `PetscWeakForm`, `DMLabel`, `PetscWeakFormRewriteKeys()`, `PetscWeakFormCreate()`, `PetscWeakFormDestroy()`
422: @*/
423: PetscErrorCode PetscWeakFormReplaceLabel(PetscWeakForm wf, DMLabel label)
424: {
425: PetscInt f;
427: PetscFunctionBegin;
428: for (f = 0; f < PETSC_NUM_WF; ++f) PetscCall(PetscWeakFormReplaceLabel_Internal(wf, wf->form[f], label));
429: PetscFunctionReturn(PETSC_SUCCESS);
430: }
432: PetscErrorCode PetscWeakFormClearIndex(PetscWeakForm wf, DMLabel label, PetscInt val, PetscInt f, PetscInt part, PetscWeakFormKind kind, PetscInt ind)
433: {
434: PetscFunctionBegin;
435: PetscCall(PetscWeakFormClearIndexFunction_Private(wf, wf->form[kind], label, val, f, part, ind));
436: PetscFunctionReturn(PETSC_SUCCESS);
437: }
439: 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[]))
440: {
441: PetscFunctionBegin;
442: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, n, (void (***)(void))obj));
443: PetscFunctionReturn(PETSC_SUCCESS);
444: }
446: 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[]))
447: {
448: PetscFunctionBegin;
449: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, n, (PetscVoidFn **)obj));
450: PetscFunctionReturn(PETSC_SUCCESS);
451: }
453: 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[]))
454: {
455: PetscFunctionBegin;
456: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, (PetscVoidFn *)obj));
457: PetscFunctionReturn(PETSC_SUCCESS);
458: }
460: 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[]))
461: {
462: PetscFunctionBegin;
463: PetscCall(PetscWeakFormGetIndexFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, ind, (PetscVoidFn **)obj));
464: PetscFunctionReturn(PETSC_SUCCESS);
465: }
467: 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[]))
468: {
469: PetscFunctionBegin;
470: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_OBJECTIVE], label, val, f, part, ind, (PetscVoidFn *)obj));
471: PetscFunctionReturn(PETSC_SUCCESS);
472: }
474: 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[]))
475: {
476: PetscFunctionBegin;
477: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, n0, (void (***)(void))f0));
478: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, n1, (void (***)(void))f1));
479: PetscFunctionReturn(PETSC_SUCCESS);
480: }
482: 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[]))
483: {
484: PetscFunctionBegin;
485: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, (PetscVoidFn *)f0));
486: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, (PetscVoidFn *)f1));
487: PetscFunctionReturn(PETSC_SUCCESS);
488: }
490: 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[]))
491: {
492: PetscFunctionBegin;
493: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, n0, (PetscVoidFn **)f0));
494: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, n1, (PetscVoidFn **)f1));
495: PetscFunctionReturn(PETSC_SUCCESS);
496: }
498: 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[]))
499: {
500: PetscFunctionBegin;
501: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_F0], label, val, f, part, i0, (PetscVoidFn *)f0));
502: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_F1], label, val, f, part, i1, (PetscVoidFn *)f1));
503: PetscFunctionReturn(PETSC_SUCCESS);
504: }
506: 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[]))
507: {
508: PetscFunctionBegin;
509: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, n0, (void (***)(void))f0));
510: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, n1, (void (***)(void))f1));
511: PetscFunctionReturn(PETSC_SUCCESS);
512: }
514: 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[]))
515: {
516: PetscFunctionBegin;
517: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, (PetscVoidFn *)f0));
518: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, (PetscVoidFn *)f1));
519: PetscFunctionReturn(PETSC_SUCCESS);
520: }
522: 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[]))
523: {
524: PetscFunctionBegin;
525: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, n0, (PetscVoidFn **)f0));
526: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, n1, (PetscVoidFn **)f1));
527: PetscFunctionReturn(PETSC_SUCCESS);
528: }
530: 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[]))
531: {
532: PetscFunctionBegin;
533: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDF0], label, val, f, part, i0, (PetscVoidFn *)f0));
534: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDF1], label, val, f, part, i1, (PetscVoidFn *)f1));
535: PetscFunctionReturn(PETSC_SUCCESS);
536: }
538: PetscErrorCode PetscWeakFormHasJacobian(PetscWeakForm wf, PetscBool *hasJac)
539: {
540: PetscInt n0, n1, n2, n3;
542: PetscFunctionBegin;
544: PetscAssertPointer(hasJac, 2);
545: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_G0], &n0));
546: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_G1], &n1));
547: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_G2], &n2));
548: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_G3], &n3));
549: *hasJac = n0 + n1 + n2 + n3 ? PETSC_TRUE : PETSC_FALSE;
550: PetscFunctionReturn(PETSC_SUCCESS);
551: }
553: 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[]))
554: {
555: PetscInt find = f * wf->Nf + g;
557: PetscFunctionBegin;
558: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, n0, (void (***)(void))g0));
559: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, n1, (void (***)(void))g1));
560: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, n2, (void (***)(void))g2));
561: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, n3, (void (***)(void))g3));
562: PetscFunctionReturn(PETSC_SUCCESS);
563: }
565: 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[]))
566: {
567: PetscInt find = f * wf->Nf + g;
569: PetscFunctionBegin;
570: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, (PetscVoidFn *)g0));
571: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, (PetscVoidFn *)g1));
572: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, (PetscVoidFn *)g2));
573: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, (PetscVoidFn *)g3));
574: PetscFunctionReturn(PETSC_SUCCESS);
575: }
577: 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[]))
578: {
579: PetscInt find = f * wf->Nf + g;
581: PetscFunctionBegin;
582: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, n0, (PetscVoidFn **)g0));
583: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, n1, (PetscVoidFn **)g1));
584: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, n2, (PetscVoidFn **)g2));
585: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, n3, (PetscVoidFn **)g3));
586: PetscFunctionReturn(PETSC_SUCCESS);
587: }
589: 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[]))
590: {
591: PetscInt find = f * wf->Nf + g;
593: PetscFunctionBegin;
594: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G0], label, val, find, part, i0, (PetscVoidFn *)g0));
595: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G1], label, val, find, part, i1, (PetscVoidFn *)g1));
596: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G2], label, val, find, part, i2, (PetscVoidFn *)g2));
597: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_G3], label, val, find, part, i3, (PetscVoidFn *)g3));
598: PetscFunctionReturn(PETSC_SUCCESS);
599: }
601: PetscErrorCode PetscWeakFormHasJacobianPreconditioner(PetscWeakForm wf, PetscBool *hasJacPre)
602: {
603: PetscInt n0, n1, n2, n3;
605: PetscFunctionBegin;
607: PetscAssertPointer(hasJacPre, 2);
608: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GP0], &n0));
609: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GP1], &n1));
610: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GP2], &n2));
611: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GP3], &n3));
612: *hasJacPre = n0 + n1 + n2 + n3 ? PETSC_TRUE : PETSC_FALSE;
613: PetscFunctionReturn(PETSC_SUCCESS);
614: }
616: 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[]))
617: {
618: PetscInt find = f * wf->Nf + g;
620: PetscFunctionBegin;
621: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, n0, (void (***)(void))g0));
622: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, n1, (void (***)(void))g1));
623: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, n2, (void (***)(void))g2));
624: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, n3, (void (***)(void))g3));
625: PetscFunctionReturn(PETSC_SUCCESS);
626: }
628: 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[]))
629: {
630: PetscInt find = f * wf->Nf + g;
632: PetscFunctionBegin;
633: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, (PetscVoidFn *)g0));
634: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, (PetscVoidFn *)g1));
635: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, (PetscVoidFn *)g2));
636: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, (PetscVoidFn *)g3));
637: PetscFunctionReturn(PETSC_SUCCESS);
638: }
640: 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[]))
641: {
642: PetscInt find = f * wf->Nf + g;
644: PetscFunctionBegin;
645: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, n0, (PetscVoidFn **)g0));
646: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, n1, (PetscVoidFn **)g1));
647: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, n2, (PetscVoidFn **)g2));
648: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, n3, (PetscVoidFn **)g3));
649: PetscFunctionReturn(PETSC_SUCCESS);
650: }
652: 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[]))
653: {
654: PetscInt find = f * wf->Nf + g;
656: PetscFunctionBegin;
657: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP0], label, val, find, part, i0, (PetscVoidFn *)g0));
658: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP1], label, val, find, part, i1, (PetscVoidFn *)g1));
659: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP2], label, val, find, part, i2, (PetscVoidFn *)g2));
660: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GP3], label, val, find, part, i3, (PetscVoidFn *)g3));
661: PetscFunctionReturn(PETSC_SUCCESS);
662: }
664: PetscErrorCode PetscWeakFormHasBdJacobian(PetscWeakForm wf, PetscBool *hasJac)
665: {
666: PetscInt n0, n1, n2, n3;
668: PetscFunctionBegin;
670: PetscAssertPointer(hasJac, 2);
671: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDG0], &n0));
672: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDG1], &n1));
673: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDG2], &n2));
674: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDG3], &n3));
675: *hasJac = n0 + n1 + n2 + n3 ? PETSC_TRUE : PETSC_FALSE;
676: PetscFunctionReturn(PETSC_SUCCESS);
677: }
679: 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[]))
680: {
681: PetscInt find = f * wf->Nf + g;
683: PetscFunctionBegin;
684: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, n0, (void (***)(void))g0));
685: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, n1, (void (***)(void))g1));
686: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, n2, (void (***)(void))g2));
687: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, n3, (void (***)(void))g3));
688: PetscFunctionReturn(PETSC_SUCCESS);
689: }
691: 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[]))
692: {
693: PetscInt find = f * wf->Nf + g;
695: PetscFunctionBegin;
696: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, (PetscVoidFn *)g0));
697: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, (PetscVoidFn *)g1));
698: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, (PetscVoidFn *)g2));
699: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, (PetscVoidFn *)g3));
700: PetscFunctionReturn(PETSC_SUCCESS);
701: }
703: 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[]))
704: {
705: PetscInt find = f * wf->Nf + g;
707: PetscFunctionBegin;
708: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, n0, (PetscVoidFn **)g0));
709: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, n1, (PetscVoidFn **)g1));
710: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, n2, (PetscVoidFn **)g2));
711: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, n3, (PetscVoidFn **)g3));
712: PetscFunctionReturn(PETSC_SUCCESS);
713: }
715: 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[]))
716: {
717: PetscInt find = f * wf->Nf + g;
719: PetscFunctionBegin;
720: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG0], label, val, find, part, i0, (PetscVoidFn *)g0));
721: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG1], label, val, find, part, i1, (PetscVoidFn *)g1));
722: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG2], label, val, find, part, i2, (PetscVoidFn *)g2));
723: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDG3], label, val, find, part, i3, (PetscVoidFn *)g3));
724: PetscFunctionReturn(PETSC_SUCCESS);
725: }
727: PetscErrorCode PetscWeakFormHasBdJacobianPreconditioner(PetscWeakForm wf, PetscBool *hasJacPre)
728: {
729: PetscInt n0, n1, n2, n3;
731: PetscFunctionBegin;
733: PetscAssertPointer(hasJacPre, 2);
734: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP0], &n0));
735: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP1], &n1));
736: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP2], &n2));
737: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_BDGP3], &n3));
738: *hasJacPre = n0 + n1 + n2 + n3 ? PETSC_TRUE : PETSC_FALSE;
739: PetscFunctionReturn(PETSC_SUCCESS);
740: }
742: 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[]))
743: {
744: PetscInt find = f * wf->Nf + g;
746: PetscFunctionBegin;
747: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, n0, (void (***)(void))g0));
748: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, n1, (void (***)(void))g1));
749: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, n2, (void (***)(void))g2));
750: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, n3, (void (***)(void))g3));
751: PetscFunctionReturn(PETSC_SUCCESS);
752: }
754: 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[]))
755: {
756: PetscInt find = f * wf->Nf + g;
758: PetscFunctionBegin;
759: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, (PetscVoidFn *)g0));
760: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, (PetscVoidFn *)g1));
761: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, (PetscVoidFn *)g2));
762: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, (PetscVoidFn *)g3));
763: PetscFunctionReturn(PETSC_SUCCESS);
764: }
766: 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[]))
767: {
768: PetscInt find = f * wf->Nf + g;
770: PetscFunctionBegin;
771: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, n0, (PetscVoidFn **)g0));
772: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, n1, (PetscVoidFn **)g1));
773: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, n2, (PetscVoidFn **)g2));
774: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, n3, (PetscVoidFn **)g3));
775: PetscFunctionReturn(PETSC_SUCCESS);
776: }
778: 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[]))
779: {
780: PetscInt find = f * wf->Nf + g;
782: PetscFunctionBegin;
783: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP0], label, val, find, part, i0, (PetscVoidFn *)g0));
784: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP1], label, val, find, part, i1, (PetscVoidFn *)g1));
785: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP2], label, val, find, part, i2, (PetscVoidFn *)g2));
786: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_BDGP3], label, val, find, part, i3, (PetscVoidFn *)g3));
787: PetscFunctionReturn(PETSC_SUCCESS);
788: }
790: PetscErrorCode PetscWeakFormHasDynamicJacobian(PetscWeakForm wf, PetscBool *hasDynJac)
791: {
792: PetscInt n0, n1, n2, n3;
794: PetscFunctionBegin;
796: PetscAssertPointer(hasDynJac, 2);
797: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GT0], &n0));
798: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GT1], &n1));
799: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GT2], &n2));
800: PetscCall(PetscHMapFormGetSize(wf->form[PETSC_WF_GT3], &n3));
801: *hasDynJac = n0 + n1 + n2 + n3 ? PETSC_TRUE : PETSC_FALSE;
802: PetscFunctionReturn(PETSC_SUCCESS);
803: }
805: 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[]))
806: {
807: PetscInt find = f * wf->Nf + g;
809: PetscFunctionBegin;
810: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, n0, (void (***)(void))g0));
811: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, n1, (void (***)(void))g1));
812: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, n2, (void (***)(void))g2));
813: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, n3, (void (***)(void))g3));
814: PetscFunctionReturn(PETSC_SUCCESS);
815: }
817: 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[]))
818: {
819: PetscInt find = f * wf->Nf + g;
821: PetscFunctionBegin;
822: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, (PetscVoidFn *)g0));
823: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, (PetscVoidFn *)g1));
824: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, (PetscVoidFn *)g2));
825: PetscCall(PetscWeakFormAddFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, (PetscVoidFn *)g3));
826: PetscFunctionReturn(PETSC_SUCCESS);
827: }
829: 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[]))
830: {
831: PetscInt find = f * wf->Nf + g;
833: PetscFunctionBegin;
834: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, n0, (PetscVoidFn **)g0));
835: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, n1, (PetscVoidFn **)g1));
836: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, n2, (PetscVoidFn **)g2));
837: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, n3, (PetscVoidFn **)g3));
838: PetscFunctionReturn(PETSC_SUCCESS);
839: }
841: 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[]))
842: {
843: PetscInt find = f * wf->Nf + g;
845: PetscFunctionBegin;
846: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT0], label, val, find, part, i0, (PetscVoidFn *)g0));
847: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT1], label, val, find, part, i1, (PetscVoidFn *)g1));
848: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT2], label, val, find, part, i2, (PetscVoidFn *)g2));
849: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_GT3], label, val, find, part, i3, (PetscVoidFn *)g3));
850: PetscFunctionReturn(PETSC_SUCCESS);
851: }
853: 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 *))
854: {
855: PetscFunctionBegin;
856: PetscCall(PetscWeakFormGetFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, n, (void (***)(void))r));
857: PetscFunctionReturn(PETSC_SUCCESS);
858: }
860: 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 *))
861: {
862: PetscFunctionBegin;
863: PetscCall(PetscWeakFormSetFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, n, (PetscVoidFn **)r));
864: PetscFunctionReturn(PETSC_SUCCESS);
865: }
867: 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 *))
868: {
869: PetscFunctionBegin;
870: PetscCall(PetscWeakFormSetIndexFunction_Private(wf, wf->form[PETSC_WF_R], label, val, f, part, i, (PetscVoidFn *)r));
871: PetscFunctionReturn(PETSC_SUCCESS);
872: }
874: /*@
875: PetscWeakFormGetNumFields - Returns the number of fields in a `PetscWeakForm`
877: Not Collective
879: Input Parameter:
880: . wf - The `PetscWeakForm` object
882: Output Parameter:
883: . Nf - The number of fields
885: Level: beginner
887: .seealso: `PetscWeakForm`, `PetscWeakFormSetNumFields()`, `PetscWeakFormCreate()`
888: @*/
889: PetscErrorCode PetscWeakFormGetNumFields(PetscWeakForm wf, PetscInt *Nf)
890: {
891: PetscFunctionBegin;
893: PetscAssertPointer(Nf, 2);
894: *Nf = wf->Nf;
895: PetscFunctionReturn(PETSC_SUCCESS);
896: }
898: /*@
899: PetscWeakFormSetNumFields - Sets the number of fields
901: Not Collective
903: Input Parameters:
904: + wf - The `PetscWeakForm` object
905: - Nf - The number of fields
907: Level: beginner
909: .seealso: `PetscWeakForm`, `PetscWeakFormGetNumFields()`, `PetscWeakFormCreate()`
910: @*/
911: PetscErrorCode PetscWeakFormSetNumFields(PetscWeakForm wf, PetscInt Nf)
912: {
913: PetscFunctionBegin;
915: wf->Nf = Nf;
916: PetscFunctionReturn(PETSC_SUCCESS);
917: }
919: /*@
920: PetscWeakFormDestroy - Destroys a `PetscWeakForm` object
922: Collective
924: Input Parameter:
925: . wf - the `PetscWeakForm` object to destroy
927: Level: developer
929: .seealso: `PetscWeakForm`, `PetscWeakFormCreate()`, `PetscWeakFormView()`
930: @*/
931: PetscErrorCode PetscWeakFormDestroy(PetscWeakForm *wf)
932: {
933: PetscInt f;
935: PetscFunctionBegin;
936: if (!*wf) PetscFunctionReturn(PETSC_SUCCESS);
939: if (--((PetscObject)*wf)->refct > 0) {
940: *wf = NULL;
941: PetscFunctionReturn(PETSC_SUCCESS);
942: }
943: ((PetscObject)*wf)->refct = 0;
944: PetscCall(PetscChunkBufferDestroy(&(*wf)->funcs));
945: for (f = 0; f < PETSC_NUM_WF; ++f) PetscCall(PetscHMapFormDestroy(&(*wf)->form[f]));
946: PetscCall(PetscFree((*wf)->form));
947: PetscCall(PetscHeaderDestroy(wf));
948: PetscFunctionReturn(PETSC_SUCCESS);
949: }
951: static PetscErrorCode PetscWeakFormViewTable_Ascii(PetscWeakForm wf, PetscViewer viewer, PetscBool splitField, const char tableName[], PetscHMapForm map)
952: {
953: PetscInt Nf = wf->Nf, Nk, k;
955: PetscFunctionBegin;
956: PetscCall(PetscHMapFormGetSize(map, &Nk));
957: if (Nk) {
958: PetscFormKey *keys;
959: PetscVoidFn **funcs = NULL;
960: const char **names;
961: PetscInt *values, *idx1, *idx2, *idx;
962: PetscBool showPart = PETSC_FALSE, showPointer = PETSC_FALSE;
963: PetscInt off = 0;
965: PetscCall(PetscMalloc6(Nk, &keys, Nk, &names, Nk, &values, Nk, &idx1, Nk, &idx2, Nk, &idx));
966: PetscCall(PetscHMapFormGetKeys(map, &off, keys));
967: /* Sort keys by label name and value */
968: {
969: /* First sort values */
970: for (k = 0; k < Nk; ++k) {
971: values[k] = keys[k].value;
972: idx1[k] = k;
973: }
974: PetscCall(PetscSortIntWithPermutation(Nk, values, idx1));
975: /* If the string sort is stable, it will be sorted correctly overall */
976: for (k = 0; k < Nk; ++k) {
977: if (keys[idx1[k]].label) PetscCall(PetscObjectGetName((PetscObject)keys[idx1[k]].label, &names[k]));
978: else names[k] = "";
979: idx2[k] = k;
980: }
981: PetscCall(PetscSortStrWithPermutation(Nk, names, idx2));
982: for (k = 0; k < Nk; ++k) {
983: if (keys[k].label) PetscCall(PetscObjectGetName((PetscObject)keys[k].label, &names[k]));
984: else names[k] = "";
985: idx[k] = idx1[idx2[k]];
986: }
987: }
988: PetscCall(PetscViewerASCIIPrintf(viewer, "%s\n", tableName));
989: PetscCall(PetscViewerASCIIPushTab(viewer));
990: for (k = 0; k < Nk; ++k) {
991: if (keys[k].part != 0) showPart = PETSC_TRUE;
992: }
993: for (k = 0; k < Nk; ++k) {
994: const PetscInt i = idx[k];
995: PetscInt n, f;
997: if (keys[i].label) {
998: if (showPointer) PetscCall(PetscViewerASCIIPrintf(viewer, "(%s:%p, %" PetscInt_FMT ") ", names[i], (void *)keys[i].label, keys[i].value));
999: else PetscCall(PetscViewerASCIIPrintf(viewer, "(%s, %" PetscInt_FMT ") ", names[i], keys[i].value));
1000: }
1001: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
1002: if (splitField) PetscCall(PetscViewerASCIIPrintf(viewer, "(%" PetscInt_FMT ", %" PetscInt_FMT ") ", keys[i].field / Nf, keys[i].field % Nf));
1003: else PetscCall(PetscViewerASCIIPrintf(viewer, "(%" PetscInt_FMT ") ", keys[i].field));
1004: if (showPart) PetscCall(PetscViewerASCIIPrintf(viewer, "(%" PetscInt_FMT ") ", keys[i].part));
1005: PetscCall(PetscWeakFormGetFunction_Private(wf, map, keys[i].label, keys[i].value, keys[i].field, keys[i].part, &n, &funcs));
1006: for (f = 0; f < n; ++f) {
1007: char *fname;
1008: size_t len, l;
1010: if (f > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", "));
1011: PetscCall(PetscDLAddr(funcs[f], &fname));
1012: if (fname) {
1013: /* Eliminate argument types */
1014: PetscCall(PetscStrlen(fname, &len));
1015: for (l = 0; l < len; ++l)
1016: if (fname[l] == '(') {
1017: fname[l] = '\0';
1018: break;
1019: }
1020: PetscCall(PetscViewerASCIIPrintf(viewer, "%s", fname));
1021: } else if (showPointer) {
1022: #if defined(__clang__)
1023: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat-pedantic")
1024: #elif defined(__GNUC__) || defined(__GNUG__)
1025: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat")
1026: #endif
1027: PetscCall(PetscViewerASCIIPrintf(viewer, "%p", funcs[f]));
1028: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_END()
1029: }
1030: PetscCall(PetscFree(fname));
1031: }
1032: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
1033: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
1034: }
1035: PetscCall(PetscViewerASCIIPopTab(viewer));
1036: PetscCall(PetscFree6(keys, names, values, idx1, idx2, idx));
1037: }
1038: PetscFunctionReturn(PETSC_SUCCESS);
1039: }
1041: static PetscErrorCode PetscWeakFormView_Ascii(PetscWeakForm wf, PetscViewer viewer)
1042: {
1043: PetscViewerFormat format;
1044: PetscInt f;
1046: PetscFunctionBegin;
1047: PetscCall(PetscViewerGetFormat(viewer, &format));
1048: PetscCall(PetscViewerASCIIPrintf(viewer, "Weak Form System with %" PetscInt_FMT " fields\n", wf->Nf));
1049: PetscCall(PetscViewerASCIIPushTab(viewer));
1050: for (f = 0; f < PETSC_NUM_WF; ++f) PetscCall(PetscWeakFormViewTable_Ascii(wf, viewer, PETSC_TRUE, PetscWeakFormKinds[f], wf->form[f]));
1051: PetscCall(PetscViewerASCIIPopTab(viewer));
1052: PetscFunctionReturn(PETSC_SUCCESS);
1053: }
1055: /*@
1056: PetscWeakFormView - Views a `PetscWeakForm`
1058: Collective
1060: Input Parameters:
1061: + wf - the `PetscWeakForm` object to view
1062: - v - the viewer
1064: Level: developer
1066: .seealso: `PetscViewer`, `PetscWeakForm`, `PetscWeakFormDestroy()`, `PetscWeakFormCreate()`
1067: @*/
1068: PetscErrorCode PetscWeakFormView(PetscWeakForm wf, PetscViewer v)
1069: {
1070: PetscBool isascii;
1072: PetscFunctionBegin;
1074: if (!v) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)wf), &v));
1076: PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &isascii));
1077: if (isascii) PetscCall(PetscWeakFormView_Ascii(wf, v));
1078: PetscTryTypeMethod(wf, view, v);
1079: PetscFunctionReturn(PETSC_SUCCESS);
1080: }
1082: /*@
1083: PetscWeakFormCreate - Creates an empty `PetscWeakForm` object.
1085: Collective
1087: Input Parameter:
1088: . comm - The communicator for the `PetscWeakForm` object
1090: Output Parameter:
1091: . wf - The `PetscWeakForm` object
1093: Level: beginner
1095: .seealso: `PetscWeakForm`, `PetscDS`, `PetscWeakFormDestroy()`
1096: @*/
1097: PetscErrorCode PetscWeakFormCreate(MPI_Comm comm, PetscWeakForm *wf)
1098: {
1099: PetscWeakForm p;
1100: PetscInt f;
1102: PetscFunctionBegin;
1103: PetscAssertPointer(wf, 2);
1104: PetscCall(PetscDSInitializePackage());
1106: PetscCall(PetscHeaderCreate(p, PETSCWEAKFORM_CLASSID, "PetscWeakForm", "Weak Form System", "PetscWeakForm", comm, PetscWeakFormDestroy, PetscWeakFormView));
1107: p->Nf = 0;
1108: PetscCall(PetscChunkBufferCreate(sizeof(&PetscWeakFormCreate), 2, &p->funcs));
1109: PetscCall(PetscMalloc1(PETSC_NUM_WF, &p->form));
1110: for (f = 0; f < PETSC_NUM_WF; ++f) PetscCall(PetscHMapFormCreate(&p->form[f]));
1111: *wf = p;
1112: PetscFunctionReturn(PETSC_SUCCESS);
1113: }