Actual source code: zasmf.c

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

  4: #if defined(PETSC_HAVE_FORTRAN_CAPS)
  5:   #define pcasmgetsubksp1_          PCASMGETSUBKSP1
  6:   #define pcasmgetsubksp2_          PCASMGETSUBKSP2
  7:   #define pcasmgetsubksp3_          PCASMGETSUBKSP3
  8:   #define pcasmgetsubksp4_          PCASMGETSUBKSP4
  9:   #define pcasmgetsubksp5_          PCASMGETSUBKSP5
 10:   #define pcasmgetsubksp6_          PCASMGETSUBKSP6
 11:   #define pcasmgetsubksp7_          PCASMGETSUBKSP7
 12:   #define pcasmgetsubksp8_          PCASMGETSUBKSP8
 13:   #define pcasmgetlocalsubmatrices_ PCASMGETLOCALSUBMATRICES
 14:   #define pcasmgetlocalsubdomains_  PCASMGETLOCALSUBDOMAINS
 15:   #define pcasmcreatesubdomains_    PCASMCREATESUBDOMAINS
 16:   #define pcasmdestroysubdomains_   PCASMDESTROYSUBDOMAINS
 17:   #define pcasmcreatesubdomains2d_  PCASMCREATESUBDOMAINS2D
 18: #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE)
 19:   #define pcasmgetsubksp1_          pcasmgetsubksp1
 20:   #define pcasmgetsubksp2_          pcasmgetsubksp2
 21:   #define pcasmgetsubksp3_          pcasmgetsubksp3
 22:   #define pcasmgetsubksp4_          pcasmgetsubksp4
 23:   #define pcasmgetsubksp5_          pcasmgetsubksp5
 24:   #define pcasmgetsubksp6_          pcasmgetsubksp6
 25:   #define pcasmgetsubksp7_          pcasmgetsubksp7
 26:   #define pcasmgetsubksp8_          pcasmgetsubksp8
 27:   #define pcasmgetlocalsubmatrices_ pcasmgetlocalsubmatrices
 28:   #define pcasmgetlocalsubdomains_  pcasmgetlocalsubdomains
 29:   #define pcasmcreatesubdomains_    pcasmcreatesubdomains
 30:   #define pcasmdestroysubdomains_   pcasmdestroysubdomains
 31:   #define pcasmcreatesubdomains2d_  pcasmcreatesubdomains2d
 32: #endif

 34: PETSC_EXTERN void pcasmcreatesubdomains2d_(PetscInt *m, PetscInt *n, PetscInt *M, PetscInt *N, PetscInt *dof, PetscInt *overlap, PetscInt *Nsub, IS *is, IS *isl, int *ierr)
 35: {
 36:   IS *iis, *iisl;
 37:   *ierr = PCASMCreateSubdomains2D(*m, *n, *M, *N, *dof, *overlap, Nsub, &iis, &iisl);
 38:   if (*ierr) return;
 39:   *ierr = PetscMemcpy(is, iis, *Nsub * sizeof(IS));
 40:   if (*ierr) return;
 41:   *ierr = PetscMemcpy(isl, iisl, *Nsub * sizeof(IS));
 42:   if (*ierr) return;
 43:   *ierr = PetscFree(iis);
 44:   if (*ierr) return;
 45:   *ierr = PetscFree(iisl);
 46: }

 48: PETSC_EXTERN void pcasmcreatesubdomains_(Mat *mat, PetscInt *n, IS *subs, PetscErrorCode *ierr)
 49: {
 50:   PetscInt i;
 51:   IS      *insubs;

 53:   *ierr = PCASMCreateSubdomains(*mat, *n, &insubs);
 54:   if (*ierr) return;
 55:   for (i = 0; i < *n; i++) subs[i] = insubs[i];
 56:   *ierr = PetscFree(insubs);
 57: }

 59: PETSC_EXTERN void pcasmdestroysubdomains_(PetscInt *n, IS *subs, IS *isubs, PetscErrorCode *ierr)
 60: {
 61:   PetscInt i;

 63:   for (i = 0; i < *n; i++) {
 64:     *ierr = ISDestroy(&subs[i]);
 65:     if (*ierr) return;
 66:     *ierr = ISDestroy(&isubs[i]);
 67:     if (*ierr) return;
 68:   }
 69: }

 71: PETSC_EXTERN void pcasmgetsubksp1_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr)
 72: {
 73:   KSP     *tksp;
 74:   PetscInt i, nloc;
 75:   CHKFORTRANNULLINTEGER(n_local);
 76:   CHKFORTRANNULLINTEGER(first_local);
 77:   CHKFORTRANNULLOBJECT(ksp);
 78:   *ierr = PCASMGetSubKSP(*pc, &nloc, first_local, &tksp);
 79:   if (n_local) *n_local = nloc;
 80:   if (ksp) {
 81:     for (i = 0; i < nloc; i++) ksp[i] = tksp[i];
 82:   }
 83: }

 85: PETSC_EXTERN void pcasmgetsubksp2_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr)
 86: {
 87:   pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr);
 88: }

 90: PETSC_EXTERN void pcasmgetsubksp3_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr)
 91: {
 92:   pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr);
 93: }

 95: PETSC_EXTERN void pcasmgetsubksp4_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr)
 96: {
 97:   pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr);
 98: }

100: PETSC_EXTERN void pcasmgetsubksp5_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr)
101: {
102:   pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr);
103: }

105: PETSC_EXTERN void pcasmgetsubksp6_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr)
106: {
107:   pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr);
108: }

110: PETSC_EXTERN void pcasmgetsubksp7_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr)
111: {
112:   pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr);
113: }

115: PETSC_EXTERN void pcasmgetsubksp8_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr)
116: {
117:   pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr);
118: }

120: PETSC_EXTERN void pcasmgetlocalsubmatrices_(PC *pc, PetscInt *n, Mat *mat, PetscErrorCode *ierr)
121: {
122:   PetscInt nloc, i;
123:   Mat     *tmat;
124:   CHKFORTRANNULLOBJECT(mat);
125:   CHKFORTRANNULLINTEGER(n);
126:   *ierr = PCASMGetLocalSubmatrices(*pc, &nloc, &tmat);
127:   if (n) *n = nloc;
128:   if (mat) {
129:     for (i = 0; i < nloc; i++) mat[i] = tmat[i];
130:   }
131: }
132: PETSC_EXTERN void pcasmgetlocalsubdomains_(PC *pc, PetscInt *n, IS *is, IS *is_local, PetscErrorCode *ierr)
133: {
134:   PetscInt nloc, i;
135:   IS      *tis, *tis_local;
136:   CHKFORTRANNULLOBJECT(is);
137:   CHKFORTRANNULLOBJECT(is_local);
138:   CHKFORTRANNULLINTEGER(n);
139:   *ierr = PCASMGetLocalSubdomains(*pc, &nloc, &tis, &tis_local);
140:   if (n) *n = nloc;
141:   if (is) {
142:     for (i = 0; i < nloc; i++) is[i] = tis[i];
143:   }
144:   if (is_local && tis_local) {
145:     for (i = 0; i < nloc; i++) is_local[i] = tis_local[i];
146:   }
147: }