1 #include "private/fortranimpl.h" 2 #include "petscmat.h" 3 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define matgetfactor_ MATGETFACTOR 6 #define matgetrowij_ MATGETROWIJ 7 #define matrestorerowij_ MATRESTOREROWIJ 8 #define matgetrow_ MATGETROW 9 #define matrestorerow_ MATRESTOREROW 10 #define matview_ MATVIEW 11 #define matgetarray_ MATGETARRAY 12 #define matrestorearray_ MATRESTOREARRAY 13 #define matconvert_ MATCONVERT 14 #define matgetsubmatrices_ MATGETSUBMATRICES 15 #define matzerorows_ MATZEROROWS 16 #define matzerorowsis_ MATZEROROWSIS 17 #define matzerorowslocal_ MATZEROROWSLOCAL 18 #define matzerorowslocalis_ MATZEROROWSLOCALIS 19 #define matsetoptionsprefix_ MATSETOPTIONSPREFIX 20 #define matgetvecs_ MATGETVECS 21 #define matnullspaceremove_ MATNULLSPACEREMOVE 22 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 23 #define matgetfactor_ matgetfactor 24 #define matgetvecs_ matgetvecs 25 #define matgetrowij_ matgetrowij 26 #define matrestorerowij_ matrestorerowij 27 #define matgetrow_ matgetrow 28 #define matrestorerow_ matrestorerow 29 #define matview_ matview 30 #define matgetarray_ matgetarray 31 #define matrestorearray_ matrestorearray 32 #define matconvert_ matconvert 33 #define matgetsubmatrices_ matgetsubmatrices 34 #define matzerorows_ matzerorows 35 #define matzerorowsis_ matzerorowsis 36 #define matzerorowslocal_ matzerorowslocal 37 #define matzerorowslocalis_ matzerorowslocalis 38 #define matsetoptionsprefix_ matsetoptionsprefix 39 #define matnullspaceremove_ matnullspaceremove 40 #endif 41 42 EXTERN_C_BEGIN 43 44 void PETSC_STDCALL matgetvecs_(Mat *mat,Vec *right,Vec *left, int *ierr ) 45 { 46 CHKFORTRANNULLOBJECT(right); 47 CHKFORTRANNULLOBJECT(left); 48 *ierr = MatGetVecs(*mat,right,left); 49 } 50 51 void PETSC_STDCALL matgetrowij_(Mat *B,PetscInt *shift,PetscTruth *sym,PetscTruth *blockcompressed,PetscInt *n,PetscInt *ia,size_t *iia, 52 PetscInt *ja,size_t *jja,PetscTruth *done,PetscErrorCode *ierr) 53 { 54 PetscInt *IA,*JA; 55 *ierr = MatGetRowIJ(*B,*shift,*sym,*blockcompressed,n,&IA,&JA,done);if (*ierr) return; 56 *iia = PetscIntAddressToFortran(ia,IA); 57 *jja = PetscIntAddressToFortran(ja,JA); 58 } 59 60 void PETSC_STDCALL matrestorerowij_(Mat *B,PetscInt *shift,PetscTruth *sym,PetscTruth *blockcompressed, PetscInt *n,PetscInt *ia,size_t *iia, 61 PetscInt *ja,size_t *jja,PetscTruth *done,PetscErrorCode *ierr) 62 { 63 PetscInt *IA = PetscIntAddressFromFortran(ia,*iia),*JA = PetscIntAddressFromFortran(ja,*jja); 64 *ierr = MatRestoreRowIJ(*B,*shift,*sym,*blockcompressed,n,&IA,&JA,done); 65 } 66 67 /* 68 This is a poor way of storing the column and value pointers 69 generated by MatGetRow() to be returned with MatRestoreRow() 70 but there is not natural,good place else to store them. Hence 71 Fortran programmers can only have one outstanding MatGetRows() 72 at a time. 73 */ 74 static PetscErrorCode matgetrowactive = 0; 75 static const PetscInt *my_ocols = 0; 76 static const PetscScalar *my_ovals = 0; 77 78 void PETSC_STDCALL matgetrow_(Mat *mat,PetscInt *row,PetscInt *ncols,PetscInt *cols,PetscScalar *vals,PetscErrorCode *ierr) 79 { 80 const PetscInt **oocols = &my_ocols; 81 const PetscScalar **oovals = &my_ovals; 82 83 if (matgetrowactive) { 84 PetscError(__LINE__,"MatGetRow_Fortran",__FILE__,__SDIR__,1,0, 85 "Cannot have two MatGetRow() active simultaneously\n\ 86 call MatRestoreRow() before calling MatGetRow() a second time"); 87 *ierr = 1; 88 return; 89 } 90 91 CHKFORTRANNULLINTEGER(cols); if (!cols) oocols = PETSC_NULL; 92 CHKFORTRANNULLSCALAR(vals); if (!vals) oovals = PETSC_NULL; 93 94 *ierr = MatGetRow(*mat,*row,ncols,oocols,oovals); 95 if (*ierr) return; 96 97 if (oocols) { *ierr = PetscMemcpy(cols,my_ocols,(*ncols)*sizeof(PetscInt)); if (*ierr) return;} 98 if (oovals) { *ierr = PetscMemcpy(vals,my_ovals,(*ncols)*sizeof(PetscScalar)); if (*ierr) return; } 99 matgetrowactive = 1; 100 } 101 102 void PETSC_STDCALL matrestorerow_(Mat *mat,PetscInt *row,PetscInt *ncols,PetscInt *cols,PetscScalar *vals,PetscErrorCode *ierr) 103 { 104 const PetscInt **oocols = &my_ocols; 105 const PetscScalar **oovals = &my_ovals; 106 if (!matgetrowactive) { 107 PetscError(__LINE__,"MatRestoreRow_Fortran",__FILE__,__SDIR__,1,0, 108 "Must call MatGetRow() first"); 109 *ierr = 1; 110 return; 111 } 112 CHKFORTRANNULLINTEGER(cols); if (!cols) oocols = PETSC_NULL; 113 CHKFORTRANNULLSCALAR(vals); if (!vals) oovals = PETSC_NULL; 114 115 *ierr = MatRestoreRow(*mat,*row,ncols,oocols,oovals); 116 matgetrowactive = 0; 117 } 118 119 void PETSC_STDCALL matview_(Mat *mat,PetscViewer *vin,PetscErrorCode *ierr) 120 { 121 PetscViewer v; 122 PetscPatchDefaultViewers_Fortran(vin,v); 123 *ierr = MatView(*mat,v); 124 } 125 126 void PETSC_STDCALL matgetarray_(Mat *mat,PetscScalar *fa,size_t *ia,PetscErrorCode *ierr) 127 { 128 PetscScalar *mm; 129 PetscInt m,n; 130 131 *ierr = MatGetArray(*mat,&mm); if (*ierr) return; 132 *ierr = MatGetSize(*mat,&m,&n); if (*ierr) return; 133 *ierr = PetscScalarAddressToFortran((PetscObject)*mat,1,fa,mm,m*n,ia); if (*ierr) return; 134 } 135 136 void PETSC_STDCALL matrestorearray_(Mat *mat,PetscScalar *fa,size_t *ia,PetscErrorCode *ierr) 137 { 138 PetscScalar *lx; 139 PetscInt m,n; 140 141 *ierr = MatGetSize(*mat,&m,&n); if (*ierr) return; 142 *ierr = PetscScalarAddressFromFortran((PetscObject)*mat,fa,*ia,m*n,&lx);if (*ierr) return; 143 *ierr = MatRestoreArray(*mat,&lx);if (*ierr) return; 144 } 145 146 void PETSC_STDCALL matgetfactor_(Mat *mat,CHAR outtype PETSC_MIXED_LEN(len),MatFactorType ftype,Mat *M,PetscErrorCode *ierr PETSC_END_LEN(len)) 147 { 148 char *t; 149 FIXCHAR(outtype,len,t); 150 *ierr = MatGetFactor(*mat,t,ftype,M); 151 FREECHAR(outtype,t); 152 } 153 154 void PETSC_STDCALL matconvert_(Mat *mat,CHAR outtype PETSC_MIXED_LEN(len),MatReuse *reuse,Mat *M,PetscErrorCode *ierr PETSC_END_LEN(len)) 155 { 156 char *t; 157 FIXCHAR(outtype,len,t); 158 *ierr = MatConvert(*mat,t,*reuse,M); 159 FREECHAR(outtype,t); 160 } 161 162 /* 163 MatGetSubmatrices() is slightly different from C since the 164 Fortran provides the array to hold the submatrix objects,while in C that 165 array is allocated by the MatGetSubmatrices() 166 */ 167 void PETSC_STDCALL matgetsubmatrices_(Mat *mat,PetscInt *n,IS *isrow,IS *iscol,MatReuse *scall,Mat *smat,PetscErrorCode *ierr) 168 { 169 Mat *lsmat; 170 PetscInt i; 171 172 if (*scall == MAT_INITIAL_MATRIX) { 173 *ierr = MatGetSubMatrices(*mat,*n,isrow,iscol,*scall,&lsmat); 174 for (i=0; i<*n; i++) { 175 smat[i] = lsmat[i]; 176 } 177 *ierr = PetscFree(lsmat); 178 } else { 179 *ierr = MatGetSubMatrices(*mat,*n,isrow,iscol,*scall,&smat); 180 } 181 } 182 183 void PETSC_STDCALL matzerorows_(Mat *mat,PetscInt *numRows,PetscInt *rows,PetscScalar *diag,PetscErrorCode *ierr) 184 { 185 *ierr = MatZeroRows(*mat,*numRows,rows,*diag); 186 } 187 188 void PETSC_STDCALL matzerorowsis_(Mat *mat,IS *is,PetscScalar *diag,PetscErrorCode *ierr) 189 { 190 *ierr = MatZeroRowsIS(*mat,*is,*diag); 191 } 192 193 void PETSC_STDCALL matzerorowslocal_(Mat *mat,PetscInt *numRows,PetscInt *rows,PetscScalar *diag,PetscErrorCode *ierr) 194 { 195 *ierr = MatZeroRowsLocal(*mat,*numRows,rows,*diag); 196 } 197 198 void PETSC_STDCALL matzerorowslocalis_(Mat *mat,IS *is,PetscScalar *diag,PetscErrorCode *ierr) 199 { 200 *ierr = MatZeroRowsLocalIS(*mat,*is,*diag); 201 } 202 203 204 void PETSC_STDCALL matsetoptionsprefix_(Mat *mat,CHAR prefix PETSC_MIXED_LEN(len), 205 PetscErrorCode *ierr PETSC_END_LEN(len)) 206 { 207 char *t; 208 209 FIXCHAR(prefix,len,t); 210 *ierr = MatSetOptionsPrefix(*mat,t); 211 FREECHAR(prefix,t); 212 } 213 214 void PETSC_STDCALL matnullspaceremove_(MatNullSpace *sp,Vec *vec,Vec *out,PetscErrorCode *ierr) 215 { 216 CHKFORTRANNULLOBJECT(out); 217 *ierr = MatNullSpaceRemove(*sp,*vec,out); 218 } 219 220 EXTERN_C_END 221