Actual source code: zasmf.c

  1: #include <petsc/private/ftnimpl.h>
  2: #include <petscksp.h>

  4: #if defined(PETSC_HAVE_FORTRAN_CAPS)
  5:   #define pcasmgetsubksp_           PCASMGETSUBKSP
  6:   #define pcasmrestoresubksp_       PCASMRESTORESUBKSP
  7:   #define pcasmgetlocalsubmatrices_ PCASMGETLOCALSUBMATRICES
  8:   #define pcasmgetlocalsubdomains_  PCASMGETLOCALSUBDOMAINS
  9:   #define pcasmcreatesubdomains_    PCASMCREATESUBDOMAINS
 10:   #define pcasmdestroysubdomains_   PCASMDESTROYSUBDOMAINS
 11:   #define pcasmcreatesubdomains2d_  PCASMCREATESUBDOMAINS2D
 12: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 13:   #define pcasmgetsubksp_           pcasmgetsubksp
 14:   #define pcasmrestoresubksp_       pcasmrestoresubksp
 15:   #define pcasmgetlocalsubmatrices_ pcasmgetlocalsubmatrices
 16:   #define pcasmgetlocalsubdomains_  pcasmgetlocalsubdomains
 17:   #define pcasmcreatesubdomains_    pcasmcreatesubdomains
 18:   #define pcasmdestroysubdomains_   pcasmdestroysubdomains
 19:   #define pcasmcreatesubdomains2d_  pcasmcreatesubdomains2d
 20: #endif

 22: PETSC_EXTERN void pcasmcreatesubdomains_(Mat *mat, PetscInt n, F90Array1d *is, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd1))
 23: {
 24:   IS *insubs;

 26:   CHKFORTRANNULLOBJECT(is);
 27:   *ierr = PCASMCreateSubdomains(*mat, n, &insubs);
 28:   if (*ierr) return;
 29:   if (insubs) *ierr = F90Array1dCreate(insubs, MPIU_FORTRANADDR, 1, n, is PETSC_F90_2PTR_PARAM(ptrd1));
 30: }

 32: PETSC_EXTERN void pcasmgetlocalsubmatrices_(PC *pc, PetscInt *n, F90Array1d *mat, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
 33: {
 34:   PetscInt nloc;
 35:   Mat     *tmat;

 37:   CHKFORTRANNULLOBJECT(mat);
 38:   CHKFORTRANNULLINTEGER(n);
 39:   *ierr = PCASMGetLocalSubmatrices(*pc, &nloc, &tmat);
 40:   if (n) *n = nloc;
 41:   if (mat) *ierr = F90Array1dCreate(tmat, MPIU_FORTRANADDR, 1, nloc, mat PETSC_F90_2PTR_PARAM(ptrd));
 42: }

 44: PETSC_EXTERN void pcasmgetlocalsubdomains_(PC *pc, PetscInt *n, F90Array1d *is, F90Array1d *is_local, int *ierr PETSC_F90_2PTR_PROTO(ptrd1) PETSC_F90_2PTR_PROTO(ptrd2))
 45: {
 46:   PetscInt nloc;
 47:   IS      *tis, *tis_local;

 49:   CHKFORTRANNULLOBJECT(is);
 50:   CHKFORTRANNULLOBJECT(is_local);
 51:   CHKFORTRANNULLINTEGER(n);
 52:   *ierr = PCASMGetLocalSubdomains(*pc, &nloc, &tis, &tis_local);
 53:   if (*ierr) return;
 54:   if (n) *n = nloc;
 55:   if (is) *ierr = F90Array1dCreate(tis, MPIU_FORTRANADDR, 1, nloc, is PETSC_F90_2PTR_PARAM(ptrd1));
 56:   if (*ierr) return;
 57:   if (is_local) *ierr = F90Array1dCreate(tis_local, MPIU_FORTRANADDR, 1, nloc, is_local PETSC_F90_2PTR_PARAM(ptrd2));
 58:   if (*ierr) return;
 59: }

 61: PETSC_EXTERN void pcasmdestroysubdomains_(PetscInt *n, F90Array1d *is1, F90Array1d *is2, int *ierr PETSC_F90_2PTR_PROTO(ptrd1) PETSC_F90_2PTR_PROTO(ptrd2))
 62: {
 63:   IS *isa, *isb;

 65:   *ierr = F90Array1dAccess(is1, MPIU_FORTRANADDR, (void **)&isa PETSC_F90_2PTR_PARAM(ptrd1));
 66:   if (*ierr) return;
 67:   *ierr = F90Array1dAccess(is2, MPIU_FORTRANADDR, (void **)&isb PETSC_F90_2PTR_PARAM(ptrd2));
 68:   if (*ierr) return;
 69:   *ierr = PCASMDestroySubdomains(*n, &isa, &isb);
 70:   if (*ierr) return;
 71:   *ierr = F90Array1dDestroy(is1, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd1));
 72:   if (*ierr) return;
 73:   *ierr = F90Array1dDestroy(is2, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd2));
 74:   if (*ierr) return;
 75: }

 77: PETSC_EXTERN void pcasmcreatesubdomains2d_(PetscInt *m, PetscInt *n, PetscInt *M, PetscInt *N, PetscInt *dof, PetscInt *overlap, PetscInt *Nsub, F90Array1d *is1, F90Array1d *is2, int *ierr PETSC_F90_2PTR_PROTO(ptrd1) PETSC_F90_2PTR_PROTO(ptrd2))
 78: {
 79:   IS *iis, *iisl;

 81:   *ierr = PCASMCreateSubdomains2D(*m, *n, *M, *N, *dof, *overlap, Nsub, &iis, &iisl);
 82:   if (*ierr) return;
 83:   *ierr = F90Array1dCreate(iis, MPIU_FORTRANADDR, 1, *Nsub, is1 PETSC_F90_2PTR_PARAM(ptrd1));
 84:   if (*ierr) return;
 85:   *ierr = F90Array1dCreate(iisl, MPIU_FORTRANADDR, 1, *Nsub, is2 PETSC_F90_2PTR_PARAM(ptrd2));
 86:   if (*ierr) return;
 87: }

 89: PETSC_EXTERN void pcasmgetsubksp_(PC *pc, PetscInt *n_local, PetscInt *first_local, F90Array1d *ksp, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
 90: {
 91:   KSP     *tksp;
 92:   PetscInt nloc, flocal;

 94:   CHKFORTRANNULLINTEGER(n_local);
 95:   CHKFORTRANNULLINTEGER(first_local);
 96:   *ierr = PCASMGetSubKSP(*pc, &nloc, &flocal, &tksp);
 97:   if (n_local) *n_local = nloc;
 98:   if (first_local) *first_local = flocal;
 99:   *ierr = F90Array1dCreate(tksp, MPIU_FORTRANADDR, 1, nloc, ksp PETSC_F90_2PTR_PARAM(ptrd));
100: }

102: PETSC_EXTERN void pcasmrestoresubksp_(PC *pc, PetscInt *n_local, PetscInt *first_local, F90Array1d *ksp, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd))
103: {
104:   *ierr = F90Array1dDestroy(ksp, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd));
105: }