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