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 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 20 #define pcasmgetsubksp1_ pcasmgetsubksp1 21 #define pcasmgetsubksp2_ pcasmgetsubksp2 22 #define pcasmgetsubksp3_ pcasmgetsubksp3 23 #define pcasmgetsubksp4_ pcasmgetsubksp4 24 #define pcasmgetsubksp5_ pcasmgetsubksp5 25 #define pcasmgetsubksp6_ pcasmgetsubksp6 26 #define pcasmgetsubksp7_ pcasmgetsubksp7 27 #define pcasmgetsubksp8_ pcasmgetsubksp8 28 #define pcasmsetlocalsubdomains_ pcasmsetlocalsubdomains 29 #define pcasmsetglobalsubdomains_ pcasmsetglobalsubdomains 30 #define pcasmgetlocalsubmatrices_ pcasmgetlocalsubmatrices 31 #define pcasmgetlocalsubdomains_ pcasmgetlocalsubdomains 32 #define pcasmcreatesubdomains_ pcasmcreatesubdomains 33 #define pcasmdestroysubdomains_ pcasmdestroysubdomains 34 #endif 35 36 PETSC_EXTERN void pcasmcreatesubdomains_(Mat *mat, PetscInt *n, IS *subs, PetscErrorCode *ierr) 37 { 38 PetscInt i; 39 IS *insubs; 40 41 *ierr = PCASMCreateSubdomains(*mat, *n, &insubs); 42 if (*ierr) return; 43 for (i = 0; i < *n; i++) subs[i] = insubs[i]; 44 *ierr = PetscFree(insubs); 45 } 46 47 PETSC_EXTERN void pcasmdestroysubdomains_(Mat *mat, PetscInt *n, IS *subs, PetscErrorCode *ierr) 48 { 49 PetscInt i; 50 51 for (i = 0; i < *n; i++) { 52 *ierr = ISDestroy(&subs[i]); 53 if (*ierr) return; 54 } 55 } 56 57 PETSC_EXTERN void pcasmgetsubksp1_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 58 { 59 KSP *tksp; 60 PetscInt i, nloc; 61 CHKFORTRANNULLINTEGER(n_local); 62 CHKFORTRANNULLINTEGER(first_local); 63 CHKFORTRANNULLOBJECT(ksp); 64 *ierr = PCASMGetSubKSP(*pc, &nloc, first_local, &tksp); 65 if (n_local) *n_local = nloc; 66 if (ksp) { 67 for (i = 0; i < nloc; i++) ksp[i] = tksp[i]; 68 } 69 } 70 71 PETSC_EXTERN void pcasmgetsubksp2_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 72 { 73 pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr); 74 } 75 76 PETSC_EXTERN void pcasmgetsubksp3_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 77 { 78 pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr); 79 } 80 81 PETSC_EXTERN void pcasmgetsubksp4_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 82 { 83 pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr); 84 } 85 86 PETSC_EXTERN void pcasmgetsubksp5_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 87 { 88 pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr); 89 } 90 91 PETSC_EXTERN void pcasmgetsubksp6_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 92 { 93 pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr); 94 } 95 96 PETSC_EXTERN void pcasmgetsubksp7_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 97 { 98 pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr); 99 } 100 101 PETSC_EXTERN void pcasmgetsubksp8_(PC *pc, PetscInt *n_local, PetscInt *first_local, KSP *ksp, PetscErrorCode *ierr) 102 { 103 pcasmgetsubksp1_(pc, n_local, first_local, ksp, ierr); 104 } 105 106 PETSC_EXTERN void pcasmsetlocalsubdomains_(PC *pc, PetscInt *n, IS *is, IS *is_local, PetscErrorCode *ierr) 107 { 108 CHKFORTRANNULLOBJECT(is); 109 CHKFORTRANNULLOBJECT(is_local); 110 *ierr = PCASMSetLocalSubdomains(*pc, *n, is, is_local); 111 } 112 113 PETSC_EXTERN void pcasmsettotalsubdomains_(PC *pc, PetscInt *N, IS *is, IS *is_local, PetscErrorCode *ierr) 114 { 115 CHKFORTRANNULLOBJECT(is); 116 CHKFORTRANNULLOBJECT(is_local); 117 *ierr = PCASMSetTotalSubdomains(*pc, *N, is, is_local); 118 } 119 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 } 148