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