Actual source code: zfieldsplitf.c

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

  4: #if defined(PETSC_HAVE_FORTRAN_CAPS)
  5:   #define pcfieldsplitgetsubksp_      PCFIELDSPLITGETSUBKSP
  6:   #define pcfieldsplitschurgetsubksp_ PCFIELDSPLITSCHURGETSUBKSP
  7:   #define pcfieldsplitsetis_          PCFIELDSPLITSETIS
  8:   #define pcfieldsplitgetis_          PCFIELDSPLITGETIS
  9:   #define pcfieldsplitsetfields_      PCFIELDSPLITSETFIELDS
 10: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 11:   #define pcfieldsplitgetsubksp_      pcfieldsplitgetsubksp
 12:   #define pcfieldsplitschurgetsubksp_ pcfieldsplitschurgetsubksp
 13:   #define pcfieldsplitsetis_          pcfieldsplitsetis
 14:   #define pcfieldsplitgetis_          pcfieldsplitgetis
 15:   #define pcfieldsplitsetfields_      pcfieldsplitsetfields
 16: #endif

 18: PETSC_EXTERN void pcfieldsplitschurgetsubksp_(PC *pc, PetscInt *n_local, KSP *ksp, PetscErrorCode *ierr)
 19: {
 20:   KSP     *tksp;
 21:   PetscInt i, nloc;
 22:   CHKFORTRANNULLINTEGER(n_local);
 23:   *ierr = PCFieldSplitSchurGetSubKSP(*pc, &nloc, &tksp);
 24:   if (*ierr) return;
 25:   if (n_local) *n_local = nloc;
 26:   CHKFORTRANNULLOBJECT(ksp);
 27:   if (ksp) {
 28:     for (i = 0; i < nloc; i++) ksp[i] = tksp[i];
 29:   }
 30:   *ierr = PetscFree(tksp);
 31: }

 33: PETSC_EXTERN void pcfieldsplitgetsubksp_(PC *pc, PetscInt *n_local, KSP *ksp, PetscErrorCode *ierr)
 34: {
 35:   KSP     *tksp;
 36:   PetscInt i, nloc;
 37:   CHKFORTRANNULLINTEGER(n_local);
 38:   *ierr = PCFieldSplitGetSubKSP(*pc, &nloc, &tksp);
 39:   if (*ierr) return;
 40:   if (n_local) *n_local = nloc;
 41:   CHKFORTRANNULLOBJECT(ksp);
 42:   if (ksp) {
 43:     for (i = 0; i < nloc; i++) ksp[i] = tksp[i];
 44:   }
 45:   *ierr = PetscFree(tksp);
 46: }

 48: PETSC_EXTERN void pcfieldsplitsetis_(PC *pc, char *splitname, IS *is, PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T len)
 49: {
 50:   char *t;
 51:   FIXCHAR(splitname, len, t);
 52:   *ierr = PCFieldSplitSetIS(*pc, t, *is);
 53:   if (*ierr) return;
 54:   FREECHAR(splitname, t);
 55: }

 57: PETSC_EXTERN void pcfieldsplitgetis_(PC *pc, char *splitname, IS *is, PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T len)
 58: {
 59:   char *t;
 60:   FIXCHAR(splitname, len, t);
 61:   *ierr = PCFieldSplitGetIS(*pc, t, is);
 62:   if (*ierr) return;
 63:   FREECHAR(splitname, t);
 64: }

 66: PETSC_EXTERN void pcfieldsplitsetfields_(PC *pc, char *splitname, PetscInt *n, const PetscInt *fields, const PetscInt *fields_col, PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T len)
 67: {
 68:   char *t;
 69:   FIXCHAR(splitname, len, t);
 70:   *ierr = PCFieldSplitSetFields(*pc, t, *n, fields, fields_col);
 71:   if (*ierr) return;
 72:   FREECHAR(splitname, t);
 73: }