Actual source code: zmatrixf.c
1: #include <petsc/private/ftnimpl.h>
2: #include <petscmat.h>
3: #include <petscviewer.h>
5: #if defined(PETSC_HAVE_FORTRAN_CAPS)
6: #define matdestroymatrices_ MATDESTROYMATRICES
7: #define matdestroysubmatrices_ MATDESTROYSUBMATRICES
8: #define matcreatesubmatrices_ MATCREATESUBMATRICES
9: #define matcreatesubmatricesmpi_ MATCREATESUBMATRICESMPI
10: #define matnullspacesetfunction_ MATNULLSPACESETFUNCTION
11: #define matfindnonzerorows_ MATFINDNONZEROROWS
12: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
13: #define matdestroymatrices_ matdestroymatrices
14: #define matdestroysubmatrices_ matdestroysubmatrices
15: #define matcreatesubmatrices_ matcreatesubmatrices
16: #define matcreatesubmatricesmpi_ matcreatesubmatricesmpi
17: #define matnullspacesetfunction_ matnullspacesetfunction
18: #define matfindnonzerorows_ matfindnonzerorows
19: #endif
21: static PetscErrorCode ournullfunction(MatNullSpace sp, Vec x, void *ctx)
22: {
23: PetscCallFortranVoidFunction((*(void (*)(MatNullSpace *, Vec *, void *, PetscErrorCode *))(((PetscObject)sp)->fortran_func_pointers[0]))(&sp, &x, ctx, &ierr));
24: return PETSC_SUCCESS;
25: }
27: PETSC_EXTERN void matnullspacesetfunction_(MatNullSpace *sp, PetscErrorCode (*rem)(MatNullSpace, Vec, void *), void *ctx, PetscErrorCode *ierr)
28: {
29: PetscObjectAllocateFortranPointers(*sp, 1);
30: ((PetscObject)*sp)->fortran_func_pointers[0] = (PetscVoidFn *)rem;
32: *ierr = MatNullSpaceSetFunction(*sp, ournullfunction, ctx);
33: }
35: PETSC_EXTERN void matcreatesubmatrices_(Mat *mat, PetscInt *n, IS *isrow, IS *iscol, MatReuse *scall, F90Array1d *ptr, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
36: {
37: Mat *lsmat;
39: if (*scall == MAT_INITIAL_MATRIX) {
40: *ierr = MatCreateSubMatrices(*mat, *n, isrow, iscol, *scall, &lsmat);
41: *ierr = F90Array1dCreate(lsmat, MPIU_FORTRANADDR, 1, *n + 1, ptr PETSC_F90_2PTR_PARAM(ptrd));
42: } else {
43: *ierr = F90Array1dAccess(ptr, MPIU_FORTRANADDR, (void **)&lsmat PETSC_F90_2PTR_PARAM(ptrd));
44: *ierr = MatCreateSubMatrices(*mat, *n, isrow, iscol, *scall, &lsmat);
45: }
46: }
48: PETSC_EXTERN void matcreatesubmatricesmpi_(Mat *mat, PetscInt *n, IS *isrow, IS *iscol, MatReuse *scall, F90Array1d *ptr, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
49: {
50: Mat *lsmat;
52: if (*scall == MAT_INITIAL_MATRIX) {
53: *ierr = MatCreateSubMatricesMPI(*mat, *n, isrow, iscol, *scall, &lsmat);
54: if (*ierr) return;
55: *ierr = F90Array1dCreate(lsmat, MPIU_FORTRANADDR, 1, *n + 1, ptr PETSC_F90_2PTR_PARAM(ptrd));
56: } else {
57: *ierr = F90Array1dAccess(ptr, MPIU_FORTRANADDR, (void **)&lsmat PETSC_F90_2PTR_PARAM(ptrd));
58: if (*ierr) return;
59: *ierr = MatCreateSubMatricesMPI(*mat, *n, isrow, iscol, *scall, &lsmat);
60: }
61: }
63: PETSC_EXTERN void matdestroymatrices_(PetscInt *n, F90Array1d *ptr, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
64: {
65: PetscInt i;
66: Mat *lsmat;
68: *ierr = F90Array1dAccess(ptr, MPIU_FORTRANADDR, (void **)&lsmat PETSC_F90_2PTR_PARAM(ptrd));
69: if (*ierr) return;
70: for (i = 0; i < *n; i++) {
71: PETSC_FORTRAN_OBJECT_F_DESTROYED_TO_C_NULL(&lsmat[i]);
72: *ierr = MatDestroy(&lsmat[i]);
73: if (*ierr) return;
74: }
75: *ierr = F90Array1dDestroy(ptr, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd));
76: if (*ierr) return;
77: *ierr = PetscFree(lsmat);
78: }
80: PETSC_EXTERN void matdestroysubmatrices_(PetscInt *n, F90Array1d *ptr, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
81: {
82: Mat *lsmat;
84: if (*n == 0) return;
85: *ierr = F90Array1dAccess(ptr, MPIU_FORTRANADDR, (void **)&lsmat PETSC_F90_2PTR_PARAM(ptrd));
86: if (*ierr) return;
87: *ierr = MatDestroySubMatrices(*n, &lsmat);
88: if (*ierr) return;
89: *ierr = F90Array1dDestroy(ptr, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd));
90: if (*ierr) return;
91: *ierr = PetscFree(lsmat);
92: }