1 #include <petsc/private/ftnimpl.h> 2 #include <petscmat.h> 3 4 #if defined(PETSC_HAVE_FORTRAN_CAPS) 5 #define matcreatenest_ MATCREATENEST 6 #define matnestsetsubmats_ MATNESTSETSUBMATS 7 #define matnestgetsubmats_ MATNESTGETSUBMATS 8 #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 9 #define matcreatenest_ matcreatenest 10 #define matnestsetsubmats_ matnestsetsubmats 11 #define matnestgetsubmats_ matnestgetsubmats 12 #endif 13 14 PETSC_EXTERN void matcreatenest_(MPI_Fint *comm, PetscInt *nr, IS is_row[], PetscInt *nc, IS is_col[], Mat a[], Mat *B, PetscErrorCode *ierr) 15 { 16 Mat *m, *tmp; 17 PetscInt i; 18 19 CHKFORTRANNULLOBJECT(is_row); 20 CHKFORTRANNULLOBJECT(is_col); 21 22 *ierr = PetscMalloc1((*nr) * (*nc), &m); 23 if (*ierr) return; 24 for (i = 0; i < (*nr) * (*nc); i++) { 25 tmp = &a[i]; 26 CHKFORTRANNULLOBJECT(tmp); 27 if (a[i] == (Mat)-2 || a[i] == (Mat)-3) { 28 (void)PetscError(MPI_Comm_f2c(*comm), __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_WRONG, PETSC_ERROR_INITIAL, "Use PETSC_NULL_MAT for missing blocks"); 29 *ierr = PETSC_ERR_ARG_WRONG; 30 return; 31 } 32 m[i] = (tmp == NULL ? NULL : a[i]); 33 } 34 *ierr = MatCreateNest(MPI_Comm_f2c(*comm), *nr, is_row, *nc, is_col, m, B); 35 if (*ierr) return; 36 *ierr = PetscFree(m); 37 } 38 39 PETSC_EXTERN void matnestsetsubmats_(Mat *B, PetscInt *nr, IS is_row[], PetscInt *nc, IS is_col[], Mat a[], PetscErrorCode *ierr) 40 { 41 Mat *m, *tmp; 42 PetscInt i; 43 MPI_Comm comm; 44 45 CHKFORTRANNULLOBJECT(is_row); 46 CHKFORTRANNULLOBJECT(is_col); 47 48 *ierr = PetscMalloc1((*nr) * (*nc), &m); 49 if (*ierr) return; 50 for (i = 0; i < (*nr) * (*nc); i++) { 51 tmp = &a[i]; 52 CHKFORTRANNULLOBJECT(tmp); 53 if (a[i] == (Mat)-2 || a[i] == (Mat)-3) { 54 *ierr = PetscObjectGetComm((PetscObject)*B, &comm); 55 if (*ierr) return; 56 (void)PetscError(comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_WRONG, PETSC_ERROR_INITIAL, "Use PETSC_NULL_MAT for missing blocks"); 57 *ierr = PETSC_ERR_ARG_WRONG; 58 return; 59 } 60 m[i] = (tmp == NULL ? NULL : a[i]); 61 } 62 *ierr = MatNestSetSubMats(*B, *nr, is_row, *nc, is_col, m); 63 if (*ierr) return; 64 *ierr = PetscFree(m); 65 } 66 67 PETSC_EXTERN void matnestgetsubmats_(Mat *A, PetscInt *M, PetscInt *N, Mat *sub, PetscErrorCode *ierr) 68 { 69 PetscInt i, j, m, n; 70 Mat **mat; 71 72 CHKFORTRANNULLINTEGER(M); 73 CHKFORTRANNULLINTEGER(N); 74 CHKFORTRANNULLOBJECT(sub); 75 76 *ierr = MatNestGetSubMats(*A, &m, &n, &mat); 77 78 if (M) *M = m; 79 if (N) *N = n; 80 if (sub) { 81 for (i = 0; i < m; i++) { 82 for (j = 0; j < n; j++) { 83 if (mat[i][j]) { 84 sub[j + n * i] = mat[i][j]; 85 } else { 86 sub[j + n * i] = (Mat)-1; 87 } 88 } 89 } 90 } 91 } 92