Actual source code: zdsf.c
1: #include <petsc/private/fortranimpl.h>
2: #include <petscds.h>
3: #include <petscviewer.h>
5: #if defined(PETSC_HAVE_FORTRAN_CAPS)
6: #define petscdsviewfromoptions_ PETSCDSVIEWFROMOPTIONS
7: #define petscdsview_ PETSCDSVIEW
8: #define petscdssetcontext_ PETSCDSSETCONTEXT
9: #define petscdssetriemannsolver_ PETSCDSSETRIEMANNSOLVER
10: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
11: #define petscdsviewfromoptions_ petscdsviewfromoptions
12: #define petscdsview_ petscdsview
13: #define petscdssetcontext_ petscdssetcontext
14: #define petscdssetriemannsolver_ petscdssetriemannsolver
15: #endif
17: static PetscFortranCallbackId riemannsolver;
19: // We can't use PetscObjectUseFortranCallback() because this function returns void
20: static void ourriemannsolver(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscInt numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx)
21: {
22: void (*func)(PetscInt *dim, PetscInt *Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], const PetscInt *numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx);
23: void *_ctx;
24: PetscCallAbort(PETSC_COMM_SELF, PetscObjectGetFortranCallback((PetscObject)ctx, PETSC_FORTRAN_CALLBACK_CLASS, riemannsolver, (PetscVoidFn **)&func, &_ctx));
25: if (func) { (*func)(&dim, &Nf, x, n, uL, uR, &numConstants, constants, flux, _ctx); }
26: }
28: PETSC_EXTERN void petscdsviewfromoptions_(PetscDS *ao, PetscObject obj, char *type, PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T len)
29: {
30: char *t;
32: FIXCHAR(type, len, t);
33: CHKFORTRANNULLOBJECT(obj);
34: *ierr = PetscDSViewFromOptions(*ao, obj, t);
35: if (*ierr) return;
36: FREECHAR(type, t);
37: }
39: PETSC_EXTERN void petscdsview_(PetscDS *prob, PetscViewer *vin, PetscErrorCode *ierr)
40: {
41: PetscViewer v;
42: PetscPatchDefaultViewers_Fortran(vin, v);
43: *ierr = PetscDSView(*prob, v);
44: if (*ierr) return;
45: }
47: PETSC_EXTERN void petscdssetcontext_(PetscDS *prob, PetscInt *f, void *ctx, PetscErrorCode *ierr)
48: {
49: *ierr = PetscDSSetContext(*prob, *f, *prob);
50: if (*ierr) return;
51: }
53: PETSC_EXTERN void petscdssetriemannsolver_(PetscDS *prob, PetscInt *f, void (*rs)(PetscInt *, PetscInt *, PetscReal *, PetscReal *, PetscScalar *, PetscScalar *, PetscInt *, PetscScalar *, PetscScalar *, void *, PetscErrorCode *), PetscErrorCode *ierr)
54: {
55: *ierr = PetscObjectSetFortranCallback((PetscObject)*prob, PETSC_FORTRAN_CALLBACK_CLASS, &riemannsolver, (PetscVoidFn *)rs, NULL);
56: if (*ierr) return;
57: *ierr = PetscDSSetRiemannSolver(*prob, *f, ourriemannsolver);
58: if (*ierr) return;
59: }