1 #include <petscmat.h> 2 #include <petsc/private/ftnimpl.h> 3 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define matmpiaijgetseqaij_ MATMPIAIJGETSEQAIJ 6 #define matmpiaijrestoreseqaij_ MATMPIAIJRESTORESEQAIJ 7 #define matdensegetarray1d_ MATDENSEGETARRAY1D 8 #define matdenserestorearray1d_ MATDENSERESTOREARRAY1D 9 #define matdensegetarrayread1d_ MATDENSEGETARRAYREAD1D 10 #define matdenserestorearrayread1d_ MATDENSERESTOREARRAYREAD1D 11 #define matdensegetarraywrite1d_ MATDENSEGETARRAYWRITE1D 12 #define matdenserestorearraywrite1d_ MATDENSERESTOREARRAYWRITE1D 13 #define matdensegetarray2d_ MATDENSEGETARRAY2D 14 #define matdenserestorearray2d_ MATDENSERESTOREARRAY2D 15 #define matdensegetarrayread2d_ MATDENSEGETARRAYREAD2D 16 #define matdenserestorearrayread2d_ MATDENSERESTOREARRAYREAD2D 17 #define matdensegetarraywrite2d_ MATDENSEGETARRAYWRITE2D 18 #define matdenserestorearraywrite2d_ MATDENSERESTOREARRAYWRITE2D 19 #define matdensegetcolumn_ MATDENSEGETCOLUMN 20 #define matdenserestorecolumn_ MATDENSERESTORECOLUMN 21 #define matseqaijgetarray_ MATSEQAIJGETARRAY 22 #define matseqaijrestorearray_ MATSEQAIJRESTOREARRAY 23 #define matgetghosts_ MATGETGHOSTS 24 #define matgetrowij_ MATGETROWIJ 25 #define matrestorerowij_ MATRESTOREROWIJ 26 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 27 #define matmpiaijgetseqaij_ matmpiaijgetseqaij 28 #define matmpiaijrestoreseqaij_ matmpiaijrestoreseqaij 29 #define matdensegetarray1d_ matdensegetarray1d 30 #define matdenserestorearray1d_ matdenserestorearray1d 31 #define matdensegetarrayread1d_ matdensegetarrayread1d 32 #define matdenserestorearrayread1d_ matdenserestorearrayread1d 33 #define matdensegetarraywrite1d_ matdensegetarraywrite1d 34 #define matdenserestorearraywrite1d_ matdenserestorearraywrite1d 35 #define matdensegetarray2d_ matdensegetarray2d 36 #define matdenserestorearray2d_ matdenserestorearray2d 37 #define matdensegetarrayread2d_ matdensegetarrayread2d 38 #define matdenserestorearrayread2d_ matdenserestorearrayread2d 39 #define matdensegetarraywrite2d_ matdensegetarraywrite2d 40 #define matdenserestorearraywrite2d_ matdenserestorearraywrite2d 41 #define matdensegetcolumn_ matdensegetcolumn 42 #define matdenserestorecolumn_ matdenserestorecolumn 43 #define matseqaijgetarray_ matseqaijgetarray 44 #define matseqaijrestorearray_ matseqaijrestorearray 45 #define matgetghosts_ matgetghosts 46 #define matgetrowij_ matgetrowij 47 #define matrestorerowij_ matrestorerowij 48 #endif 49 50 PETSC_EXTERN void matgetghosts_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 51 { 52 const PetscInt *ghosts; 53 PetscInt N; 54 55 *ierr = MatGetGhosts(*mat, &N, &ghosts); 56 if (*ierr) return; 57 *ierr = F90Array1dCreate((PetscInt *)ghosts, MPIU_INT, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd)); 58 } 59 PETSC_EXTERN void matdensegetarray2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 60 { 61 PetscScalar *fa; 62 PetscInt m, N, lda; 63 *ierr = MatDenseGetArray(*mat, &fa); 64 if (*ierr) return; 65 *ierr = MatGetLocalSize(*mat, &m, NULL); 66 if (*ierr) return; 67 *ierr = MatGetSize(*mat, NULL, &N); 68 if (*ierr) return; 69 *ierr = MatDenseGetLDA(*mat, &lda); 70 if (*ierr) return; 71 if (m != lda) { // TODO: add F90Array2dLDACreate() 72 *ierr = PetscError(((PetscObject)*mat)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "Array lda %" PetscInt_FMT " must match number of local rows %" PetscInt_FMT, lda, m); 73 return; 74 } 75 *ierr = F90Array2dCreate(fa, MPIU_SCALAR, 1, m, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd)); 76 } 77 PETSC_EXTERN void matdenserestorearray2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 78 { 79 PetscScalar *fa; 80 *ierr = F90Array2dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 81 if (*ierr) return; 82 *ierr = F90Array2dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 83 if (*ierr) return; 84 *ierr = MatDenseRestoreArray(*mat, &fa); 85 } 86 PETSC_EXTERN void matdensegetarrayread2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 87 { 88 const PetscScalar *fa; 89 PetscInt m, N, lda; 90 *ierr = MatDenseGetArrayRead(*mat, &fa); 91 if (*ierr) return; 92 *ierr = MatGetLocalSize(*mat, &m, NULL); 93 if (*ierr) return; 94 *ierr = MatGetSize(*mat, NULL, &N); 95 if (*ierr) return; 96 *ierr = MatDenseGetLDA(*mat, &lda); 97 if (*ierr) return; 98 if (m != lda) { // TODO: add F90Array2dLDACreate() 99 *ierr = PetscError(((PetscObject)*mat)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "Array lda %" PetscInt_FMT " must match number of local rows %" PetscInt_FMT, lda, m); 100 return; 101 } 102 *ierr = F90Array2dCreate((void **)fa, MPIU_SCALAR, 1, m, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd)); 103 } 104 PETSC_EXTERN void matdenserestorearrayread2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 105 { 106 const PetscScalar *fa; 107 *ierr = F90Array2dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 108 if (*ierr) return; 109 *ierr = F90Array2dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 110 if (*ierr) return; 111 *ierr = MatDenseRestoreArrayRead(*mat, &fa); 112 } 113 PETSC_EXTERN void matdensegetarraywrite2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 114 { 115 PetscScalar *fa; 116 PetscInt m, N, lda; 117 *ierr = MatDenseGetArrayWrite(*mat, &fa); 118 if (*ierr) return; 119 *ierr = MatGetLocalSize(*mat, &m, NULL); 120 if (*ierr) return; 121 *ierr = MatGetSize(*mat, NULL, &N); 122 if (*ierr) return; 123 *ierr = MatDenseGetLDA(*mat, &lda); 124 if (*ierr) return; 125 if (m != lda) { // TODO: add F90Array2dLDACreate() 126 *ierr = PetscError(((PetscObject)*mat)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "Array lda %" PetscInt_FMT " must match number of local rows %" PetscInt_FMT, lda, m); 127 return; 128 } 129 *ierr = F90Array2dCreate(fa, MPIU_SCALAR, 1, m, 1, N, ptr PETSC_F90_2PTR_PARAM(ptrd)); 130 } 131 PETSC_EXTERN void matdenserestorearraywrite2d_(Mat *mat, F90Array2d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 132 { 133 PetscScalar *fa; 134 *ierr = F90Array2dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 135 if (*ierr) return; 136 *ierr = F90Array2dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 137 if (*ierr) return; 138 *ierr = MatDenseRestoreArrayWrite(*mat, &fa); 139 } 140 PETSC_EXTERN void matdensegetarray1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 141 { 142 PetscScalar *fa; 143 PetscInt m, N, lda; 144 *ierr = MatDenseGetArray(*mat, &fa); 145 if (*ierr) return; 146 *ierr = MatGetLocalSize(*mat, &m, NULL); 147 if (*ierr) return; 148 *ierr = MatGetSize(*mat, NULL, &N); 149 if (*ierr) return; 150 *ierr = MatDenseGetLDA(*mat, &lda); 151 if (*ierr) return; 152 if (m != lda) { // TODO: add F90Array1dLDACreate() 153 *ierr = PetscError(((PetscObject)*mat)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "Array lda %" PetscInt_FMT " must match number of local rows %" PetscInt_FMT, lda, m); 154 return; 155 } 156 *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, m * N, ptr PETSC_F90_2PTR_PARAM(ptrd)); 157 } 158 PETSC_EXTERN void matdenserestorearray1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 159 { 160 PetscScalar *fa; 161 *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 162 if (*ierr) return; 163 *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 164 if (*ierr) return; 165 *ierr = MatDenseRestoreArray(*mat, &fa); 166 } 167 PETSC_EXTERN void matdensegetarrayread1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 168 { 169 const PetscScalar *fa; 170 PetscInt m, N, lda; 171 *ierr = MatDenseGetArrayRead(*mat, &fa); 172 if (*ierr) return; 173 *ierr = MatGetLocalSize(*mat, &m, NULL); 174 if (*ierr) return; 175 *ierr = MatGetSize(*mat, NULL, &N); 176 if (*ierr) return; 177 *ierr = MatDenseGetLDA(*mat, &lda); 178 if (*ierr) return; 179 if (m != lda) { // TODO: add F90Array1dLDACreate() 180 *ierr = PetscError(((PetscObject)*mat)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "Array lda %" PetscInt_FMT " must match number of local rows %" PetscInt_FMT, lda, m); 181 return; 182 } 183 *ierr = F90Array1dCreate((void **)fa, MPIU_SCALAR, 1, m * N, ptr PETSC_F90_2PTR_PARAM(ptrd)); 184 } 185 PETSC_EXTERN void matdenserestorearrayread1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 186 { 187 const PetscScalar *fa; 188 *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 189 if (*ierr) return; 190 *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 191 if (*ierr) return; 192 *ierr = MatDenseRestoreArrayRead(*mat, &fa); 193 } 194 PETSC_EXTERN void matdensegetarraywrite1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 195 { 196 PetscScalar *fa; 197 PetscInt m, N, lda; 198 *ierr = MatDenseGetArrayWrite(*mat, &fa); 199 if (*ierr) return; 200 *ierr = MatGetLocalSize(*mat, &m, NULL); 201 if (*ierr) return; 202 *ierr = MatGetSize(*mat, NULL, &N); 203 if (*ierr) return; 204 *ierr = MatDenseGetLDA(*mat, &lda); 205 if (*ierr) return; 206 if (m != lda) { // TODO: add F90Array1dLDACreate() 207 *ierr = PetscError(((PetscObject)*mat)->comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_BADPTR, PETSC_ERROR_INITIAL, "Array lda %" PetscInt_FMT " must match number of local rows %" PetscInt_FMT, lda, m); 208 return; 209 } 210 *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, m * N, ptr PETSC_F90_2PTR_PARAM(ptrd)); 211 } 212 PETSC_EXTERN void matdenserestorearraywrite1d_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 213 { 214 PetscScalar *fa; 215 *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 216 if (*ierr) return; 217 *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 218 if (*ierr) return; 219 *ierr = MatDenseRestoreArrayWrite(*mat, &fa); 220 } 221 PETSC_EXTERN void matdensegetcolumn_(Mat *mat, PetscInt *col, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 222 { 223 PetscScalar *fa; 224 PetscInt m; 225 *ierr = MatDenseGetColumn(*mat, *col, &fa); 226 if (*ierr) return; 227 *ierr = MatGetLocalSize(*mat, &m, NULL); 228 if (*ierr) return; 229 *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, m, ptr PETSC_F90_2PTR_PARAM(ptrd)); 230 } 231 PETSC_EXTERN void matdenserestorecolumn_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 232 { 233 PetscScalar *fa; 234 *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 235 if (*ierr) return; 236 *ierr = MatDenseRestoreColumn(*mat, &fa); 237 } 238 PETSC_EXTERN void matseqaijgetarray_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 239 { 240 PetscScalar *fa; 241 PetscInt m, n; 242 *ierr = MatSeqAIJGetArray(*mat, &fa); 243 if (*ierr) return; 244 *ierr = MatGetLocalSize(*mat, &m, &n); 245 if (*ierr) return; 246 *ierr = F90Array1dCreate(fa, MPIU_SCALAR, 1, m * n, ptr PETSC_F90_2PTR_PARAM(ptrd)); 247 } 248 PETSC_EXTERN void matseqaijrestorearray_(Mat *mat, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 249 { 250 PetscScalar *fa; 251 *ierr = F90Array1dAccess(ptr, MPIU_SCALAR, (void **)&fa PETSC_F90_2PTR_PARAM(ptrd)); 252 if (*ierr) return; 253 *ierr = F90Array1dDestroy(ptr, MPIU_SCALAR PETSC_F90_2PTR_PARAM(ptrd)); 254 if (*ierr) return; 255 *ierr = MatSeqAIJRestoreArray(*mat, &fa); 256 } 257 PETSC_EXTERN void matgetrowij_(Mat *B, PetscInt *shift, PetscBool *sym, PetscBool *blockcompressed, PetscInt *n, F90Array1d *ia, F90Array1d *ja, PetscBool *done, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(iad) PETSC_F90_2PTR_PROTO(jad)) 258 { 259 const PetscInt *IA, *JA; 260 *ierr = MatGetRowIJ(*B, *shift, *sym, *blockcompressed, n, &IA, &JA, done); 261 if (*ierr) return; 262 if (!*done) return; 263 *ierr = F90Array1dCreate((PetscInt *)IA, MPIU_INT, 1, *n + 1, ia PETSC_F90_2PTR_PARAM(iad)); 264 *ierr = F90Array1dCreate((PetscInt *)JA, MPIU_INT, 1, IA[*n], ja PETSC_F90_2PTR_PARAM(jad)); 265 } 266 PETSC_EXTERN void matrestorerowij_(Mat *B, PetscInt *shift, PetscBool *sym, PetscBool *blockcompressed, PetscInt *n, F90Array1d *ia, F90Array1d *ja, PetscBool *done, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(iad) PETSC_F90_2PTR_PROTO(jad)) 267 { 268 const PetscInt *IA, *JA; 269 *ierr = F90Array1dAccess(ia, MPIU_INT, (void **)&IA PETSC_F90_2PTR_PARAM(iad)); 270 if (*ierr) return; 271 *ierr = F90Array1dDestroy(ia, MPIU_INT PETSC_F90_2PTR_PARAM(iad)); 272 if (*ierr) return; 273 *ierr = F90Array1dAccess(ja, MPIU_INT, (void **)&JA PETSC_F90_2PTR_PARAM(jad)); 274 if (*ierr) return; 275 *ierr = F90Array1dDestroy(ja, MPIU_INT PETSC_F90_2PTR_PARAM(jad)); 276 if (*ierr) return; 277 *ierr = MatRestoreRowIJ(*B, *shift, *sym, *blockcompressed, n, &IA, &JA, done); 278 } 279 PETSC_EXTERN void matmpiaijgetseqaij_(Mat *mat, Mat *A, Mat *B, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 280 { 281 const PetscInt *fa; 282 PetscInt n; 283 *ierr = MatMPIAIJGetSeqAIJ(*mat, A, B, &fa); 284 if (*ierr) return; 285 *ierr = MatGetLocalSize(*B, NULL, &n); 286 if (*ierr) return; 287 *ierr = F90Array1dCreate((void *)fa, MPIU_INT, 1, n, ptr PETSC_F90_2PTR_PARAM(ptrd)); 288 } 289 PETSC_EXTERN void matmpiaijrestoreseqaij_(Mat *mat, Mat *A, Mat *B, F90Array1d *ptr, int *ierr PETSC_F90_2PTR_PROTO(ptrd)) 290 { 291 *ierr = F90Array1dDestroy(ptr, MPIU_INT PETSC_F90_2PTR_PARAM(ptrd)); 292 if (*ierr) return; 293 } 294