16dd63270SBarry Smith #include <petsc/private/ftnimpl.h> 245c38901SJed Brown #include <petscmat.h> 345c38901SJed Brown 445c38901SJed Brown #if defined(PETSC_HAVE_FORTRAN_CAPS) 545c38901SJed Brown #define matcreatenest_ MATCREATENEST 6*58ad77e8SBarry Smith #define matnestsetsubmats_ MATNESTSETSUBMATS 7ffa9b3b1SVincent Le Chenadec #define matnestgetsubmats_ MATNESTGETSUBMATS 845c38901SJed Brown #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 945c38901SJed Brown #define matcreatenest_ matcreatenest 10*58ad77e8SBarry Smith #define matnestsetsubmats_ matnestsetsubmats 11ffa9b3b1SVincent Le Chenadec #define matnestgetsubmats_ matnestgetsubmats 1245c38901SJed Brown #endif 1345c38901SJed Brown 14*58ad77e8SBarry Smith PETSC_EXTERN void matcreatenest_(MPI_Fint *comm, PetscInt *nr, IS is_row[], PetscInt *nc, IS is_col[], Mat a[], Mat *B, PetscErrorCode *ierr) 1545c38901SJed Brown { 162f6eced2SAlex Fikl Mat *m, *tmp; 172f6eced2SAlex Fikl PetscInt i; 182f6eced2SAlex Fikl 1945c38901SJed Brown CHKFORTRANNULLOBJECT(is_row); 2045c38901SJed Brown CHKFORTRANNULLOBJECT(is_col); 212f6eced2SAlex Fikl 225975b3b6SBarry Smith *ierr = PetscMalloc1((*nr) * (*nc), &m); 235975b3b6SBarry Smith if (*ierr) return; 242f6eced2SAlex Fikl for (i = 0; i < (*nr) * (*nc); i++) { 25f4f49eeaSPierre Jolivet tmp = &a[i]; 262f6eced2SAlex Fikl CHKFORTRANNULLOBJECT(tmp); 27*58ad77e8SBarry Smith if (a[i] == (Mat)-2 || a[i] == (Mat)-3) { 28*58ad77e8SBarry Smith (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*58ad77e8SBarry Smith *ierr = PETSC_ERR_ARG_WRONG; 30*58ad77e8SBarry Smith return; 31*58ad77e8SBarry Smith } 322f6eced2SAlex Fikl m[i] = (tmp == NULL ? NULL : a[i]); 332f6eced2SAlex Fikl } 345975b3b6SBarry Smith *ierr = MatCreateNest(MPI_Comm_f2c(*comm), *nr, is_row, *nc, is_col, m, B); 355975b3b6SBarry Smith if (*ierr) return; 362f6eced2SAlex Fikl *ierr = PetscFree(m); 3745c38901SJed Brown } 383a4d7b9aSSatish Balay 39*58ad77e8SBarry Smith PETSC_EXTERN void matnestsetsubmats_(Mat *B, PetscInt *nr, IS is_row[], PetscInt *nc, IS is_col[], Mat a[], PetscErrorCode *ierr) 40*58ad77e8SBarry Smith { 41*58ad77e8SBarry Smith Mat *m, *tmp; 42*58ad77e8SBarry Smith PetscInt i; 43*58ad77e8SBarry Smith MPI_Comm comm; 44*58ad77e8SBarry Smith 45*58ad77e8SBarry Smith CHKFORTRANNULLOBJECT(is_row); 46*58ad77e8SBarry Smith CHKFORTRANNULLOBJECT(is_col); 47*58ad77e8SBarry Smith 48*58ad77e8SBarry Smith *ierr = PetscMalloc1((*nr) * (*nc), &m); 49*58ad77e8SBarry Smith if (*ierr) return; 50*58ad77e8SBarry Smith for (i = 0; i < (*nr) * (*nc); i++) { 51*58ad77e8SBarry Smith tmp = &a[i]; 52*58ad77e8SBarry Smith CHKFORTRANNULLOBJECT(tmp); 53*58ad77e8SBarry Smith if (a[i] == (Mat)-2 || a[i] == (Mat)-3) { 54*58ad77e8SBarry Smith *ierr = PetscObjectGetComm((PetscObject)*B, &comm); 55*58ad77e8SBarry Smith if (*ierr) return; 56*58ad77e8SBarry Smith (void)PetscError(comm, __LINE__, PETSC_FUNCTION_NAME, __FILE__, PETSC_ERR_ARG_WRONG, PETSC_ERROR_INITIAL, "Use PETSC_NULL_MAT for missing blocks"); 57*58ad77e8SBarry Smith *ierr = PETSC_ERR_ARG_WRONG; 58*58ad77e8SBarry Smith return; 59*58ad77e8SBarry Smith } 60*58ad77e8SBarry Smith m[i] = (tmp == NULL ? NULL : a[i]); 61*58ad77e8SBarry Smith } 62*58ad77e8SBarry Smith *ierr = MatNestSetSubMats(*B, *nr, is_row, *nc, is_col, m); 63*58ad77e8SBarry Smith if (*ierr) return; 64*58ad77e8SBarry Smith *ierr = PetscFree(m); 65*58ad77e8SBarry Smith } 66*58ad77e8SBarry Smith 67*58ad77e8SBarry Smith PETSC_EXTERN void matnestgetsubmats_(Mat *A, PetscInt *M, PetscInt *N, Mat *sub, PetscErrorCode *ierr) 68ffa9b3b1SVincent Le Chenadec { 69351962e3SVincent Le Chenadec PetscInt i, j, m, n; 70ffa9b3b1SVincent Le Chenadec Mat **mat; 71351962e3SVincent Le Chenadec 72351962e3SVincent Le Chenadec CHKFORTRANNULLINTEGER(M); 73351962e3SVincent Le Chenadec CHKFORTRANNULLINTEGER(N); 74351962e3SVincent Le Chenadec CHKFORTRANNULLOBJECT(sub); 75351962e3SVincent Le Chenadec 76351962e3SVincent Le Chenadec *ierr = MatNestGetSubMats(*A, &m, &n, &mat); 77351962e3SVincent Le Chenadec 785975b3b6SBarry Smith if (M) { *M = m; } 795975b3b6SBarry Smith if (N) { *N = n; } 80351962e3SVincent Le Chenadec if (sub) { 81351962e3SVincent Le Chenadec for (i = 0; i < m; i++) { 82351962e3SVincent Le Chenadec for (j = 0; j < n; j++) { 832f6eced2SAlex Fikl if (mat[i][j]) { 84351962e3SVincent Le Chenadec sub[j + n * i] = mat[i][j]; 852f6eced2SAlex Fikl } else { 861c8b34f3SBarry Smith sub[j + n * i] = (Mat)-1; 872f6eced2SAlex Fikl } 88351962e3SVincent Le Chenadec } 89ffa9b3b1SVincent Le Chenadec } 90ffa9b3b1SVincent Le Chenadec } 91ffa9b3b1SVincent Le Chenadec } 92