Naïve question, but is your problem really symmetric? For example, if you do FEM, how do you impose Dirichlet boundary conditions? If you can reproduce this on a smaller example and are at liberty to share the Mat, IS, and RHS Vec, feel free to send them at [email protected], I can have a look.
Thanks, Pierre > On 6 Sep 2023, at 8:14 AM, Carl-Johan Thore <[email protected]> wrote: > > > Thanks for this! > I now get convergence with sbaij and fieldsplit. However, as can be seen in > the attached file, it requires around three times as many outer iterations as > with aij to reach the same point (to within rtol). Switching from > -pc_fieldsplit_schur_fact_type LOWER to FULL (which is the "correct" choice?) > helps quite a bit, but still sbaij takes many more iterations. So it seems > that your sbaij-implementation works but that I'm doing something wrong when > setting up the fieldsplit pc, or that some of the subsolvers doesn't work > properly with sbaij. I've tried bjacobi as you had in your logs, but not > hypre yet. But anyway, one doesn't expect the matrix format to have this much > impact on convergence? Any other suggestions? > /Carl-Johan > >> On Mon, Sep 4, 2023 at 9:31 PM Pierre Jolivet <[email protected]> wrote: >> The branch should now be good to go >> (https://gitlab.com/petsc/petsc/-/merge_requests/6841). >> Sorry, I made a mistake before, hence the error on PetscObjectQuery(). >> I’m not sure the code will be covered by the pipeline, but I have tested >> this on a Raviart—Thomas discretization with PCFIELDSPLIT. >> You’ll see in the attached logs that: >> 1) the numerics match >> 2) in the SBAIJ case, PCFIELDSPLIT extract the (non-symmetric) A_{01} block >> from the global (symmetric) A and we get the A_{10} block cheaply by just >> using MatCreateHermitianTranspose(), instead of calling another time >> MatCreateSubMatrix() >> Please let me know if you have some time to test the branch and whether it >> fails or succeeds on your test cases. >> >> Also, I do not agree with what Hong said. >> Sometimes, the assembly of a coefficient can be more expensive than the >> communication of the said coefficient. >> So they are instances where SBAIJ would be more efficient than AIJ even if >> it would require more communication, it is not a black and white picture. >> >> Thanks, >> Pierre >> >> >>>> On 28 Aug 2023, at 12:12 PM, Pierre Jolivet <[email protected]> wrote: >>>> >>>> >>>> On 28 Aug 2023, at 6:50 PM, Carl-Johan Thore <[email protected]> >>>> wrote: >>>> >>>> I've tried the new files, and with them, PCFIELDSPLIT now gets set up >>>> without crashes (but the setup is significantly slower than for MATAIJ) >>> >>> I’ll be back from Japan at the end of this week, my schedule is too packed >>> to get anything done in the meantime. >>> But I’ll let you know when things are working properly (last I checked, I >>> think it was working, but I may have forgotten about a corner case or two). >>> But yes, though one would except things to be faster and less memory >>> intensive with SBAIJ, it’s unfortunately not always the case. >>> >>> Thanks, >>> Pierre >>> >>>> Unfortunately I still get errors later in the process: >>>> >>>> [0]PETSC ERROR: Null argument, when expecting valid pointer >>>> [0]PETSC ERROR: Null Pointer: Parameter # 1 >>>> [0]PETSC ERROR: Petsc Development GIT revision: v3.19.4-1023-ga6d78fcba1d >>>> GIT Date: 2023-08-22 20:32:33 -0400 >>>> [0]PETSC ERROR: Configure options -f --with-fortran-bindings=0 --with-cuda >>>> --with-cusp --download-scalapack --download-hdf5 --download-zlib >>>> --download-mumps --download-parmetis --download-metis --download-ptscotch >>>> --download-hypre --download-spai >>>> [0]PETSC ERROR: #1 PetscObjectQuery() at >>>> /mnt/c/mathware/petsc/petsc-v3-19-4/src/sys/objects/inherit.c:742 >>>> [0]PETSC ERROR: #2 MatCreateSubMatrix_MPISBAIJ() at >>>> /mnt/c/mathware/petsc/petsc-v3-19-4/src/mat/impls/sbaij/mpi/mpisbaij.c:1414 >>>> [0]PETSC ERROR: #3 MatCreateSubMatrix() at >>>> /mnt/c/mathware/petsc/petsc-v3-19-4/src/mat/interface/matrix.c:8476 >>>> [0]PETSC ERROR: #4 PCSetUp_FieldSplit() at >>>> /mnt/c/mathware/petsc/petsc-v3-19-4/src/ksp/pc/impls/fieldsplit/fieldsplit.c:826 >>>> [0]PETSC ERROR: #5 PCSetUp() at >>>> /mnt/c/mathware/petsc/petsc-v3-19-4/src/ksp/pc/interface/precon.c:1069 >>>> [0]PETSC ERROR: #6 KSPSetUp() at >>>> /mnt/c/mathware/petsc/petsc-v3-19-4/src/ksp/ksp/interface/itfunc.c:415 >>>> >>>> The code I'm running here works without any problems for MATAIJ. To run it >>>> with MATSBAIJ I've simply used the command-line option >>>> -dm_mat_type sbaij >>>> >>>> >>>> Kind regards, >>>> Carl-Johan >>>> >>>> >>>>> On Sat, Aug 26, 2023 at 5:21 PM Pierre Jolivet via petsc-users >>>>> <[email protected]> wrote: >>>>> >>>>> >>>>>>> On 27 Aug 2023, at 12:14 AM, Carl-Johan Thore <[email protected]> >>>>>>> wrote: >>>>>>> >>>>>>> “Well, your A00 and A11 will possibly be SBAIJ also, so you’ll end up >>>>>>> with the same issue.” >>>>>>> I’m not sure I follow. Does PCFIELDSPLIT extract further submatrices >>>>>>> from these blocks, or is there >>>>>>> somewhere else in the code that things will go wrong? >>>>>> >>>>>> Ah, no, you are right, in that case it should work. >>>>>> >>>>>> For the MATNEST I was thinking to get some savings from the >>>>>> block-symmetry at least >>>>>> even if symmetry in A00 and A11 cannot be exploited; using SBAIJ for >>>>>> them would just be a >>>>>> (pretty big) bonus. >>>>>> >>>>>> “I’ll rebase on top of main and try to get it integrated if it could be >>>>>> useful to you (but I’m traveling >>>>>> right now so everything gets done more slowly, sorry).” >>>>>> Sound great, Thanks again! >>>>> >>>>> The MR is there https://gitlab.com/petsc/petsc/-/merge_requests/6841. >>>>> I need to add a new code path in MatCreateRedundantMatrix() to make sure >>>>> the resulting Mat is indeed SBAIJ, but that is orthogonal to the >>>>> PCFIELDSPLIT issue. >>>>> The branch should be usable in its current state. >>>>> >>>>> Thanks, >>>>> Pierre >>>>> >>>>>> >>>>>> From: Pierre Jolivet <[email protected]> >>>>>> Sent: Saturday, August 26, 2023 4:36 PM >>>>>> To: Carl-Johan Thore <[email protected]> >>>>>> Cc: Carl-Johan Thore <[email protected]>; [email protected] >>>>>> Subject: Re: [petsc-users] PCFIELDSPLIT with MATSBAIJ >>>>>> >>>>>> >>>>>> >>>>>> >>>>>> On 26 Aug 2023, at 11:16 PM, Carl-Johan Thore <[email protected]> >>>>>> wrote: >>>>>> >>>>>> "(Sadly) MATSBAIJ is extremely broken, in particular, it cannot be used >>>>>> to retrieve rectangular blocks in MatCreateSubMatrices, thus you cannot >>>>>> get the A01 and A10 blocks in PCFIELDSPLIT. >>>>>> I have a branch that fixes this, but I haven’t rebased in a while (and >>>>>> I’m AFK right now), would you want me to rebase and give it a go, or >>>>>> must you stick to a release tarball?" >>>>>> >>>>>> Ok, would be great if you could look at this! I don't need to stick to >>>>>> any particular branch. >>>>>> >>>>>> Do you think MATNEST could be an alternative here? >>>>>> >>>>>> Well, your A00 and A11 will possibly be SBAIJ also, so you’ll end up >>>>>> with the same issue. >>>>>> I’m using both approaches (monolithic SBAIJ or Nest + SBAIJ), it was >>>>>> crashing but I think it was thoroughly fixed in >>>>>> https://gitlab.com/petsc/petsc/-/commits/jolivet/feature-matcreatesubmatrices-rectangular-sbaij/ >>>>>> It is ugly code on top of ugly code, so I didn’t try to get it >>>>>> integrated and just used the branch locally, and then moved to some >>>>>> other stuff. >>>>>> I’ll rebase on top of main and try to get it integrated if it could be >>>>>> useful to you (but I’m traveling right now so everything gets done more >>>>>> slowly, sorry). >>>>>> >>>>>> Thanks, >>>>>> Pierre >>>>>> >>>>>> >>>>>> My matrix is >>>>>> [A00 A01; >>>>>> A01^t A11] >>>>>> so perhaps with MATNEST I can make use of the block-symmetry at least, >>>>>> and then use MATSBAIJ for >>>>>> A00 and A11 if it's possible to combine matrix types which the manual >>>>>> seems to imply. >>>>>> >>>>>> Kind regards >>>>>> Carl-Johan >>>>>> >>>>>> >>>>>> >>>>>> On 26 Aug 2023, at 10:09 PM, Carl-Johan Thore <[email protected]> >>>>>> wrote: >>>>>> >>>>>> Hi, >>>>>> >>>>>> I'm trying to use PCFIELDSPLIT with MATSBAIJ in PETSc 3.19.4. >>>>>> According to the manual "[t]he fieldsplit preconditioner cannot >>>>>> currently be used with the MATBAIJ or MATSBAIJ data formats if the >>>>>> blocksize is larger than 1". Since my blocksize is exactly 1 it would >>>>>> seem that I can use PCFIELDSPLIT. But this fails with "PETSC ERROR: For >>>>>> symmetric format, iscol must equal isrow" >>>>>> from MatCreateSubMatrix_MPISBAIJ. Tracing backwards one ends up in >>>>>> fieldsplit.c at >>>>>> >>>>>> /* extract the A01 and A10 matrices */ ilink = jac->head; >>>>>> PetscCall(ISComplement(ilink->is_col, rstart, rend, &ccis)); if >>>>>> (jac->offdiag_use_amat) { PetscCall(MatCreateSubMatrix(pc->mat, >>>>>> ilink->is, ccis, MAT_INITIAL_MATRIX, &jac->B)); } else { >>>>>> PetscCall(MatCreateSubMatrix(pc->pmat, ilink->is, ccis, >>>>>> MAT_INITIAL_MATRIX, &jac->B)); } >>>>>> >>>>>> This, since my A01 and A10 are not square, seems to explain why iscol is >>>>>> not equal to isrow. >>>>>> From this I gather that it is in general NOT possible to use >>>>>> PCFIELDSPLIT with MATSBAIJ even with block size 1? >>>>>> >>>>>> Kind regards, >>>>>> Carl-Johan >>>>> >>> >> > > <experiment_with_sbaij.txt>
