[petsc-users] Understanding index sets for PCGASM

[Barry Smith]

   Thank you for the test code. I have a fix in the branch  barry/2023-04-29/fix-pcasmcreatesubdomains2d <https://gitlab.com/petsc/petsc/-/commits/barry/2023-04-29/fix-pcasmcreatesubdomains2d> with merge request https://gitlab.com/petsc/petsc/-/merge_requests/6394

   The functions did not have proper Fortran stubs and interfaces so I had to provide them manually in the new branch.

   Use

   git fetch
   git checkout barry/2023-04-29/fix-pcasmcreatesubdomains2d <https://gitlab.com/petsc/petsc/-/commits/barry/2023-04-29/fix-pcasmcreatesubdomains2d>
   ./configure etc

   Your now working test code is in src/ksp/ksp/tests/ex71f.F90  I had to change things slightly and I updated the error handling for the latest version.

   Please let us know if you have any later questions.

  Barry




> On Apr 28, 2023, at 12:07 PM, LEONARDO MUTTI <leonardo.mutti01 at universitadipavia.it> wrote:
> 
> Hello. I am having a hard time understanding the index sets to feed PCGASMSetSubdomains, and I am working in Fortran (as a PETSc novice). To get more intuition on how the IS objects behave I tried the following minimal (non) working example, which should tile a 16x16 matrix into 16 square, non-overlapping submatrices:
> 
> #include <petsc/finclude/petscmat.h>
> #include <petsc/finclude/petscksp.h>
> #include <petsc/finclude/petscpc.h>
>       USE petscmat
>       USE petscksp
>       USE petscpc
>       
>       Mat   :: A
>       PetscInt :: M, NSubx, dof, overlap, NSub
>       INTEGER :: I,J
>       PetscErrorCode      :: ierr
>       PetscScalar :: v
>       KSP            :: ksp
>       PC             :: pc
>       IS :: subdomains_IS, inflated_IS
>       
>       call PetscInitialize(PETSC_NULL_CHARACTER , ierr)
>       
> !-----Create a dummy matrix
>       M = 16
>       call MatCreateAIJ(MPI_COMM_WORLD, PETSC_DECIDE, PETSC_DECIDE,
>      & M, M,
>      & PETSC_DEFAULT_INTEGER, PETSC_NULL_INTEGER,
>      & PETSC_DEFAULT_INTEGER, PETSC_NULL_INTEGER,
>      & A, ierr)
>       
>       DO I=1,M
>          DO J=1,M
>             v = I*J
>             CALL MatSetValue (A,I-1,J-1,v,
>      &                     INSERT_VALUES , ierr)
>          END DO
>       END DO
>       
>       call MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY , ierr)
>       call MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY , ierr)
>       
> !-----Create KSP and PC
>       call KSPCreate(PETSC_COMM_WORLD,ksp, ierr)
>       call KSPSetOperators(ksp,A,A, ierr)
>       call KSPSetType(ksp,"bcgs",ierr)
>       call KSPGetPC(ksp,pc,ierr)
>       call KSPSetUp(ksp, ierr)
>       call PCSetType(pc,PCGASM, ierr)
>       call PCSetUp(pc , ierr)
>       
> !-----GASM setup
>       NSubx = 4
>       dof = 1
>       overlap = 0
>       
>       call PCGASMCreateSubdomains2D(pc, 
>      &      M, M,
>      &      NSubx, NSubx,
>      &      dof, overlap,
>      &      NSub, subdomains_IS, inflated_IS, ierr)
> 
>       call ISView(subdomains_IS, PETSC_VIEWER_STDOUT_WORLD, ierr)
>       
>       call KSPDestroy(ksp, ierr)  
>       call PetscFinalize(ierr)
> 
> Running this on one processor, I get NSub = 4.
> If PCASM and PCASMCreateSubdomains2D are used instead, I get NSub = 16 as expected.
> Moreover, I get in the end "forrtl: severe (157): Program Exception - access violation". So:
> 1) why do I get two different results with ASM, and GASM?
> 2) why do I get access violation and how can I solve this?
> In fact, in C, subdomains_IS, inflated_IS should pointers to IS objects. As I see on the Fortran interface, the arguments to PCGASMCreateSubdomains2D are IS objects:
> 
>        subroutine PCGASMCreateSubdomains2D(a,b,c,d,e,f,g,h,i,j,z)
>        import tPC,tIS
>        PC a ! PC
>        PetscInt b ! PetscInt
>        PetscInt c ! PetscInt
>        PetscInt d ! PetscInt
>        PetscInt e ! PetscInt
>        PetscInt f ! PetscInt
>        PetscInt g ! PetscInt
>        PetscInt h ! PetscInt
>        IS i ! IS
>        IS j ! IS
>        PetscErrorCode z
>        end subroutine PCGASMCreateSubdomains2D
> Thus:
> 3) what should be inside e.g., subdomains_IS? I expect it to contain, for every created subdomain, the list of rows and columns defining the subblock in the matrix, am I right?
> 
> Context: I have a block-tridiagonal system arising from space-time finite elements, and I want to solve it with GMRES+PCGASM preconditioner, where each overlapping submatrix is on the diagonal and of size 3x3 blocks (and spanning multiple processes). This is PETSc 3.17.1 on Windows.
> 
> Thanks in advance,
> Leonardo

-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20230429/86f224ef/attachment.html>