16dd63270SBarry Smith #include <petsc/private/ftnimpl.h> 2c6db04a5SJed Brown #include <petscdm.h> 3665c2dedSJed Brown #include <petscviewer.h> 49a42bb27SBarry Smith 59a42bb27SBarry Smith #if defined(PETSC_HAVE_FORTRAN_CAPS) 6b14c0cbaSBlaise Bourdin #define dmcreatesuperdm_ DMCREATESUPERDM 78d9ecca5SBarry Smith #define dmcreatefielddecompositiongetname_ DMCREATEFIELDDECOMPOSITIONGETNAME 88d9ecca5SBarry Smith #define dmcreatefielddecompositiongetisdm_ DMCREATEFIELDDECOMPOSITIONGETISDM 98d9ecca5SBarry Smith #define dmcreatefielddecompositionrestoreisdm_ DMCREATEFIELDDECOMPOSITIONRESTOREISDM 109a42bb27SBarry Smith #elif !defined(PETSC_HAVE_FORTRAN_UNDERSCORE) 11b14c0cbaSBlaise Bourdin #define dmcreatesuperdm_ dmreatesuperdm 128d9ecca5SBarry Smith #define dmcreatefielddecompositiongetname_ dmcreatefielddecompositiongetname 138d9ecca5SBarry Smith #define dmcreatefielddecompositiongetisdm_ dmcreatefielddecompositiongetisdm 148d9ecca5SBarry Smith #define dmcreatefielddecompositionrestoreisdm_ dmcreatefielddecompositionrestoreisdm 159a42bb27SBarry Smith #endif 169a42bb27SBarry Smith 178d9ecca5SBarry Smith PETSC_EXTERN void dmcreatesuperdm_(DM dms[], PetscInt *len, IS ***is, DM *superdm, PetscErrorCode *ierr) 186823f3c5SBlaise Bourdin { 196823f3c5SBlaise Bourdin *ierr = DMCreateSuperDM(dms, *len, *is, superdm); 206823f3c5SBlaise Bourdin } 218d9ecca5SBarry Smith 228d9ecca5SBarry Smith PETSC_EXTERN void dmcreatefielddecompositiongetname_(DM *dm, PetscInt *i, char *name, PetscErrorCode *ierr, PETSC_FORTRAN_CHARLEN_T l_b) 238d9ecca5SBarry Smith { 248d9ecca5SBarry Smith PetscInt n; 258d9ecca5SBarry Smith char **names; 268d9ecca5SBarry Smith *ierr = DMCreateFieldDecomposition(*dm, &n, &names, NULL, NULL); 278d9ecca5SBarry Smith if (*ierr) return; 288d9ecca5SBarry Smith *ierr = PetscStrncpy((char *)name, names[*i - 1], l_b); 298d9ecca5SBarry Smith if (*ierr) return; 308d9ecca5SBarry Smith FIXRETURNCHAR(PETSC_TRUE, name, l_b); 31*ac530a7eSPierre Jolivet for (PetscInt j = 0; j < n; j++) *ierr = PetscFree(names[j]); 328d9ecca5SBarry Smith *ierr = PetscFree(names); 338d9ecca5SBarry Smith } 348d9ecca5SBarry Smith 358d9ecca5SBarry Smith PETSC_EXTERN void dmcreatefielddecompositiongetisdm_(DM *dm, F90Array1d *iss, F90Array1d *dms, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd1) PETSC_F90_2PTR_PROTO(ptrd2)) 368d9ecca5SBarry Smith { 378d9ecca5SBarry Smith PetscInt n; 388d9ecca5SBarry Smith IS *tis; 398d9ecca5SBarry Smith DM *tdm; 408d9ecca5SBarry Smith 418d9ecca5SBarry Smith if (iss && dms) { 428d9ecca5SBarry Smith *ierr = DMCreateFieldDecomposition(*dm, &n, NULL, &tis, &tdm); 438d9ecca5SBarry Smith } else if (iss) { 448d9ecca5SBarry Smith *ierr = DMCreateFieldDecomposition(*dm, &n, NULL, &tis, NULL); 458d9ecca5SBarry Smith } else if (dms) { 468d9ecca5SBarry Smith *ierr = DMCreateFieldDecomposition(*dm, &n, NULL, NULL, &tdm); 478d9ecca5SBarry Smith } 488d9ecca5SBarry Smith if (*ierr) return; 498d9ecca5SBarry Smith if (iss) *ierr = F90Array1dCreate(tis, MPIU_FORTRANADDR, 1, n, iss PETSC_F90_2PTR_PARAM(ptrd1)); 508d9ecca5SBarry Smith if (*ierr) return; 518d9ecca5SBarry Smith if (dms) *ierr = F90Array1dCreate(tdm, MPIU_FORTRANADDR, 1, n, dms PETSC_F90_2PTR_PARAM(ptrd2)); 528d9ecca5SBarry Smith } 538d9ecca5SBarry Smith 548d9ecca5SBarry Smith PETSC_EXTERN void dmcreatefielddecompositionrestoreisdm_(DM *dm, F90Array1d *iss, F90Array1d *dms, PetscErrorCode *ierr PETSC_F90_2PTR_PROTO(ptrd1) PETSC_F90_2PTR_PROTO(ptrd2)) 558d9ecca5SBarry Smith { 568d9ecca5SBarry Smith PetscInt n; 578d9ecca5SBarry Smith 588d9ecca5SBarry Smith *ierr = DMGetNumFields(*dm, &n); 598d9ecca5SBarry Smith if (*ierr) return; 608d9ecca5SBarry Smith if (iss) { 618d9ecca5SBarry Smith IS *tis; 628d9ecca5SBarry Smith *ierr = F90Array1dAccess(iss, MPIU_FORTRANADDR, (void **)&tis PETSC_F90_2PTR_PARAM(ptrd1)); 638d9ecca5SBarry Smith if (*ierr) return; 648d9ecca5SBarry Smith *ierr = F90Array1dDestroy(iss, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd1)); 658d9ecca5SBarry Smith if (*ierr) return; 66*ac530a7eSPierre Jolivet for (PetscInt i = 0; i < n; i++) *ierr = ISDestroy(&tis[i]); 678d9ecca5SBarry Smith *ierr = PetscFree(tis); 688d9ecca5SBarry Smith if (*ierr) return; 698d9ecca5SBarry Smith } 708d9ecca5SBarry Smith if (dms) { 718d9ecca5SBarry Smith DM *tdm; 728d9ecca5SBarry Smith *ierr = F90Array1dAccess(dms, MPIU_FORTRANADDR, (void **)&tdm PETSC_F90_2PTR_PARAM(ptrd2)); 738d9ecca5SBarry Smith if (*ierr) return; 748d9ecca5SBarry Smith *ierr = F90Array1dDestroy(dms, MPIU_FORTRANADDR PETSC_F90_2PTR_PARAM(ptrd2)); 758d9ecca5SBarry Smith if (*ierr) return; 76*ac530a7eSPierre Jolivet for (PetscInt i = 0; i < n; i++) *ierr = DMDestroy(&tdm[i]); 778d9ecca5SBarry Smith *ierr = PetscFree(tdm); 788d9ecca5SBarry Smith if (*ierr) return; 798d9ecca5SBarry Smith } 808d9ecca5SBarry Smith } 81