1 #include <petsc/private/fortranimpl.h> 2 #include <petscksp.h> 3 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 pcasmsetlocalsubdomains_ PCASMSETLOCALSUBDOMAINS 14 #define pcasmsetglobalsubdomains_ PCASMSETGLOBALSUBDOMAINS 15 #define pcasmgetlocalsubmatrices_ PCASMGETLOCALSUBMATRICES 16 #define pcasmgetlocalsubdomains_ PCASMGETLOCALSUBDOMAINS 17 #define pcasmcreatesubdomains_ PCASMCREATESUBDOMAINS 18 #define pcasmdestroysubdomains_ PCASMDESTROYSUBDOMAINS 19 #define pcasmcreatesubdomains2d_ PCASMCREATESUBDOMAINS2D 20 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 21 #define pcasmgetsubksp1_ pcasmgetsubksp1 22 #define pcasmgetsubksp2_ pcasmgetsubksp2 23 #define pcasmgetsubksp3_ pcasmgetsubksp3 24 #define pcasmgetsubksp4_ pcasmgetsubksp4 25 #define pcasmgetsubksp5_ pcasmgetsubksp5 26 #define pcasmgetsubksp6_ pcasmgetsubksp6 27 #define pcasmgetsubksp7_ pcasmgetsubksp7 28 #define pcasmgetsubksp8_ pcasmgetsubksp8 29 #define pcasmsetlocalsubdomains_ pcasmsetlocalsubdomains 30 #define pcasmsetglobalsubdomains_ pcasmsetglobalsubdomains 31 #define pcasmgetlocalsubmatrices_ pcasmgetlocalsubmatrices 32 #define pcasmgetlocalsubdomains_ pcasmgetlocalsubdomains 33 #define pcasmcreatesubdomains_ pcasmcreatesubdomains 34 #define pcasmdestroysubdomains_ pcasmdestroysubdomains 35 #define pcasmcreatesubdomains2d_ pcasmcreatesubdomains2d 36 #endif 37 38 PETSC_EXTERN void pcasmcreatesubdomains2d_(PetscInt *m, PetscInt *n, PetscInt *M, PetscInt *N, PetscInt *dof, PetscInt *overlap, PetscInt *Nsub, IS *is, IS *isl, int *ierr) 39 { 40 IS *iis, *iisl; 41 *ierr = PCASMCreateSubdomains2D(*m, *n, *M, *N, *dof, *overlap, Nsub, &iis, &iisl); 42 if (*ierr) return; 43 *ierr = PetscMemcpy(is, iis, *Nsub * sizeof(IS)); 44 if (*ierr) return; 45 *ierr = PetscMemcpy(isl, iisl, *Nsub * sizeof(IS)); 46 if (*ierr) return; 47 *ierr = PetscFree(iis); 48 if (*ierr) return; 49 *ierr = PetscFree(iisl); 50 } 51 52 PETSC_EXTERN void pcasmcreatesubdomains_(Mat *mat, PetscInt *n, IS *subs, PetscErrorCode *ierr) 53 { 54 PetscInt i; 55 IS *insubs; 56 57 *ierr = PCASMCreateSubdomains(*mat, *n, &insubs); 58 if (*ierr) return; 59 for (i = 0; i < *n; i++) subs[i] = insubs[i]; 60 *ierr = PetscFree(insubs); 61 } 62 63 PETSC_EXTERN void pcasmdestroysubdomains_(PetscInt *n, IS *subs, IS *isubs, PetscErrorCode *ierr) 64 { 65 PetscInt i; 66 67 for (i = 0; i < *n; i++) { 68 *ierr = ISDestroy(&subs[i]); 69 if (*ierr) return; 70 *ierr = ISDestroy(&isubs[i]); 71 if (*ierr) return; 72 } 73 } 74 75 PETSC_EXTERN void pcasmgetsubksp1_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 76 { 77 KSP *tksp; 78 PetscInt i, nloc; 79 CHKFORTRANNULLINTEGER(n_local); 80 CHKFORTRANNULLINTEGER(first_local); 81 CHKFORTRANNULLOBJECT(ksp); 82 *ierr = PCASMGetSubKSP(*pc, &nloc, first_local, &tksp); 83 if (n_local) *n_local = nloc; 84 if (ksp) { 85 for (i = 0; i < nloc; i++) ksp[i] = tksp[i]; 86 } 87 } 88 89 PETSC_EXTERN void pcasmgetsubksp2_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 90 { 91 pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr); 92 } 93 94 PETSC_EXTERN void pcasmgetsubksp3_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 95 { 96 pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr); 97 } 98 99 PETSC_EXTERN void pcasmgetsubksp4_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 100 { 101 pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr); 102 } 103 104 PETSC_EXTERN void pcasmgetsubksp5_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 105 { 106 pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr); 107 } 108 109 PETSC_EXTERN void pcasmgetsubksp6_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 110 { 111 pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr); 112 } 113 114 PETSC_EXTERN void pcasmgetsubksp7_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 115 { 116 pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr); 117 } 118 119 PETSC_EXTERN void pcasmgetsubksp8_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 120 { 121 pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr); 122 } 123 124 PETSC_EXTERN void pcasmsetlocalsubdomains_(PC *pc, PetscInt *n, IS *is, IS *is_local, PetscErrorCode *ierr) 125 { 126 CHKFORTRANNULLOBJECT(is); 127 CHKFORTRANNULLOBJECT(is_local); 128 *ierr = PCASMSetLocalSubdomains(*pc, *n, is, is_local); 129 } 130 131 PETSC_EXTERN void pcasmsettotalsubdomains_(PC *pc, PetscInt *N, IS *is, IS *is_local, PetscErrorCode *ierr) 132 { 133 CHKFORTRANNULLOBJECT(is); 134 CHKFORTRANNULLOBJECT(is_local); 135 *ierr = PCASMSetTotalSubdomains(*pc, *N, is, is_local); 136 } 137 138 PETSC_EXTERN void pcasmgetlocalsubmatrices_(PC *pc, PetscInt *n, Mat *mat, PetscErrorCode *ierr) 139 { 140 PetscInt nloc, i; 141 Mat *tmat; 142 CHKFORTRANNULLOBJECT(mat); 143 CHKFORTRANNULLINTEGER(n); 144 *ierr = PCASMGetLocalSubmatrices(*pc, &nloc, &tmat); 145 if (n) *n = nloc; 146 if (mat) { 147 for (i = 0; i < nloc; i++) mat[i] = tmat[i]; 148 } 149 } 150 PETSC_EXTERN void pcasmgetlocalsubdomains_(PC *pc, PetscInt *n, IS *is, IS *is_local, PetscErrorCode *ierr) 151 { 152 PetscInt nloc, i; 153 IS *tis, *tis_local; 154 CHKFORTRANNULLOBJECT(is); 155 CHKFORTRANNULLOBJECT(is_local); 156 CHKFORTRANNULLINTEGER(n); 157 *ierr = PCASMGetLocalSubdomains(*pc, &nloc, &tis, &tis_local); 158 if (n) *n = nloc; 159 if (is) { 160 for (i = 0; i < nloc; i++) is[i] = tis[i]; 161 } 162 if (is_local && tis_local) { 163 for (i = 0; i < nloc; i++) is_local[i] = tis_local[i]; 164 } 165 } 166