1 #include <petsc/private/ftnimpl.h>
2 #include <petsc/private/matimpl.h>
3
4 /* Declare these pointer types instead of void* for clarity, but do not include petscts.h so that this code does have an actual reverse dependency. */
5 typedef struct _p_TS *TS;
6 typedef struct _p_SNES *SNES;
7
8 #if defined(PETSC_HAVE_FORTRAN_CAPS)
9 #define matfdcoloringsetfunctionts_ MATFDCOLORINGSETFUNCTIONTS
10 #define matfdcoloringsetfunction_ MATFDCOLORINGSETFUNCTION
11 #define matfdcoloringgetperturbedcolumns_ MATFDCOLORINGGETPERTURBEDCOLUMNS
12 #define matfdcoloringrestoreperturbedcolumns_ MATFDCOLORINGRESTOREPERTURBEDCOLUMNS
13 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
14 #define matfdcoloringsetfunctionts_ matfdcoloringsetfunctionts
15 #define matfdcoloringsetfunction_ matfdcoloringsetfunction
16 #define matfdcoloringgetperturbedcolumns_ matfdcoloringgetperturbedcolumns
17 #define matfdcoloringrestoreperturbedcolumns_ matfdcoloringrestoreperturbedcolumns
18 #endif
19
matfdcoloringgetperturbedcolumns_(MatFDColoring * x,PetscInt * len,F90Array1d * ptr,int * __ierr PETSC_F90_2PTR_PROTO (ptrd))20 PETSC_EXTERN void matfdcoloringgetperturbedcolumns_(MatFDColoring *x, PetscInt *len, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
21 {
22 const PetscInt *fa;
23
24 *__ierr = MatFDColoringGetPerturbedColumns(*x, len, &fa);
25 if (*__ierr) return;
26 *__ierr = F90Array1dCreate((void *)fa, MPIU_INT, 1, *len, ptr PETSC_F90_2PTR_PARAM(ptrd));
27 }
matfdcoloringrestoreperturbedcolumns_(MatFDColoring * x,PetscInt * len,F90Array1d * ptr,int * __ierr PETSC_F90_2PTR_PROTO (ptrd))28 PETSC_EXTERN void matfdcoloringrestoreperturbedcolumns_(MatFDColoring *x, PetscInt *len, F90Array1d *ptr, int *__ierr PETSC_F90_2PTR_PROTO(ptrd))
29 {
30 *__ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd));
31 }
32
33 /* These are not extern C because they are passed into non-extern C user level functions */
ourmatfdcoloringfunctionts(TS ts,PetscReal t,Vec x,Vec y,MatFDColoring fd)34 static PetscErrorCode ourmatfdcoloringfunctionts(TS ts, PetscReal t, Vec x, Vec y, MatFDColoring fd)
35 {
36 PetscErrorCode ierr = PETSC_SUCCESS;
37 (*(void (*)(TS *, PetscReal *, Vec *, Vec *, void *, PetscErrorCode *))fd->ftn_func_pointer)(&ts, &t, &x, &y, fd->ftn_func_cntx, &ierr);
38 return ierr;
39 }
40
ourmatfdcoloringfunctionsnes(SNES snes,Vec x,Vec y,MatFDColoring fd)41 static PetscErrorCode ourmatfdcoloringfunctionsnes(SNES snes, Vec x, Vec y, MatFDColoring fd)
42 {
43 PetscErrorCode ierr = PETSC_SUCCESS;
44 (*(void (*)(SNES *, Vec *, Vec *, void *, PetscErrorCode *))fd->ftn_func_pointer)(&snes, &x, &y, fd->ftn_func_cntx, &ierr);
45 return ierr;
46 }
47
48 /*
49 MatFDColoringSetFunction sticks the Fortran function and its context into the MatFDColoring structure and passes the MatFDColoring object
50 in as the function context. ourmafdcoloringfunctionsnes() and ourmatfdcoloringfunctionts() then access the function and its context from the
51 MatFDColoring that is passed in. This is the same way that fortran_func_pointers is used in PETSc objects.
52
53 NOTE: FORTRAN USER CANNOT PUT IN A NEW J OR B currently.
54 */
55
matfdcoloringsetfunctionts_(MatFDColoring * fd,void (* f)(TS *,double *,Vec *,Vec *,void *,PetscErrorCode *),PetscCtx ctx,PetscErrorCode * ierr)56 PETSC_EXTERN void matfdcoloringsetfunctionts_(MatFDColoring *fd, void (*f)(TS *, double *, Vec *, Vec *, void *, PetscErrorCode *), PetscCtx ctx, PetscErrorCode *ierr)
57 {
58 (*fd)->ftn_func_pointer = (PetscFortranCallbackFn *)f;
59 (*fd)->ftn_func_cntx = ctx;
60
61 *ierr = MatFDColoringSetFunction(*fd, (MatFDColoringFn *)(PetscVoidFn *)ourmatfdcoloringfunctionts, *fd);
62 }
63
matfdcoloringsetfunction_(MatFDColoring * fd,void (* f)(SNES *,Vec *,Vec *,void *,PetscErrorCode *),PetscCtx ctx,PetscErrorCode * ierr)64 PETSC_EXTERN void matfdcoloringsetfunction_(MatFDColoring *fd, void (*f)(SNES *, Vec *, Vec *, void *, PetscErrorCode *), PetscCtx ctx, PetscErrorCode *ierr)
65 {
66 (*fd)->ftn_func_pointer = (PetscFortranCallbackFn *)f;
67 (*fd)->ftn_func_cntx = ctx;
68
69 *ierr = MatFDColoringSetFunction(*fd, (MatFDColoringFn *)ourmatfdcoloringfunctionsnes, *fd);
70 }
71