I will add the two interfaces you requested today. (Likely you may need more also).
Barry > On May 4, 2023, at 6:01 PM, Matthew Knepley <knep...@gmail.com> wrote: > > On Thu, May 4, 2023 at 1:43 PM LEONARDO MUTTI > <leonardo.mutt...@universitadipavia.it > <mailto:leonardo.mutt...@universitadipavia.it>> wrote: >> Of course, I'll try to explain. >> >> I am solving a parabolic equation with space-time FEM and I want an >> efficient solver/preconditioner for the resulting system. >> The corresponding matrix, call it X, has an e.g. block bi-diagonal >> structure, if the cG(1)-dG(0) method is used (i.e. implicit Euler solved in >> batch). >> Every block-row of X corresponds to a time instant. >> >> I want to introduce parallelism in time by subdividing X into overlapping >> submatrices of e.g 2x2 or 3x3 blocks, along the block diagonal. >> For instance, call X_i the individual blocks. The submatrices would be, for >> various i, (X_{i-1,i-1},X_{i-1,i};X_{i,i-1},X_{i,i}). >> I'd like each submatrix to be solved in parallel, to combine the various >> results together in an ASM like fashion. >> Every submatrix has thus a predecessor and a successor, and it overlaps with >> both, so that as far as I could understand, GASM has to be used in place of >> ASM. > > Yes, ordered that way you need GASM. I wonder if inverting the ordering would > be useful, namely putting the time index on the inside. > Then the blocks would be over all time, but limited space, which is more the > spirit of ASM I think. > > Have you considered waveform relaxation for this problem? > > Thanks, > > Matt > >> Hope this helps. >> Best, >> Leonardo >> >> Il giorno gio 4 mag 2023 alle ore 18:05 Matthew Knepley <knep...@gmail.com >> <mailto:knep...@gmail.com>> ha scritto: >>> On Thu, May 4, 2023 at 11:24 AM LEONARDO MUTTI >>> <leonardo.mutt...@universitadipavia.it >>> <mailto:leonardo.mutt...@universitadipavia.it>> wrote: >>>> Thank you for the help. >>>> Adding to my example: >>>> call PCGASMSetSubdomains(pc,NSub, subdomains_IS, inflated_IS,ierr) >>>> call PCGASMDestroySubdomains(NSub,subdomains_IS,inflated_IS,ierr) >>>> results in: >>>> Error LNK2019 unresolved external symbol PCGASMDESTROYSUBDOMAINS >>>> referenced in function ... >>>> Error LNK2019 unresolved external symbol PCGASMSETSUBDOMAINS >>>> referenced in function ... >>>> I'm not sure if the interfaces are missing or if I have a compilation >>>> problem. >>> >>> I just want to make sure you really want GASM. It sounded like you might >>> able to do what you want just with ASM. >>> Can you tell me again what you want to do overall? >>> >>> Thanks, >>> >>> Matt >>> >>>> Thank you again. >>>> Best, >>>> Leonardo >>>> >>>> Il giorno sab 29 apr 2023 alle ore 20:30 Barry Smith <bsm...@petsc.dev >>>> <mailto:bsm...@petsc.dev>> ha scritto: >>>>> >>>>> 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.mutt...@universitadipavia.it >>>>>> <mailto:leonardo.mutt...@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 >>>>> >>> >>> >>> -- >>> What most experimenters take for granted before they begin their >>> experiments is infinitely more interesting than any results to which their >>> experiments lead. >>> -- Norbert Wiener >>> >>> https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/> > > > -- > What most experimenters take for granted before they begin their experiments > is infinitely more interesting than any results to which their experiments > lead. > -- Norbert Wiener > > https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
