# SNES for bordered block matrix **URL:** \ **Category:** General Discussion\ **Created:** [December 16, 2022, 4:24pm UTC](https://community.freefem.org/t/snes-for-bordered-block-matrix/2193 "2022-12-16T16:24:48Z")\ **Posts on this page:** 6\ **Page:** 2
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [December 19, 2022, 8:54pm UTC](https://community.freefem.org/t/snes-for-bordered-block-matrix/2193/23 "2022-12-19T20:54:47Z")
Even with one process, I don’t get a deterministic behavior. ```auto finding initial point initial point found. moving to continuation part starting loop: 1/5. Current nu = 0.1. Predicted nu = 0.09973 0 SNES Function norm 2.805428653365e-07 nu = 0.09973 1 SNES Function norm 1.695141030744e-07 nu = 0.0997301 2 SNES Function norm 6.278092520575e-13 nu = 0.0997308 3 SNES Function norm 1.354552924947e-14 Nonlinear solve converged due to CONVERGED_SNORM_RELATIVE iterations 3 starting loop: 2/5. Current nu = 0.0997308. Predicted nu = 0.0994608 0 SNES Function norm 1.891297254561e-07 nu = 0.0994608 1 SNES Function norm 7.815465447968e-13 nu = 0.0994616 2 SNES Function norm 1.375679948496e-14 Nonlinear solve converged due to CONVERGED_SNORM_RELATIVE iterations 2 starting loop: 3/5. Current nu = 0.0994616. Predicted nu = 0.0991916 0 SNES Function norm 2.822903166950e-07 nu = 0.0991916 1 SNES Function norm 1.709241490889e-07 nu = 0.0991917 2 SNES Function norm 6.395819877626e-13 nu = 0.0991924 3 SNES Function norm 1.349931255097e-14 Nonlinear solve converged due to CONVERGED_SNORM_RELATIVE iterations 3 ``` vs. ```auto finding initial point initial point found. moving to continuation part starting loop: 1/5. Current nu = 0.1. Predicted nu = 0.09973 0 SNES Function norm 2.805428653365e-07 nu = 0.09973 1 SNES Function norm 1.695141030744e-07 nu = 0.0997301 2 SNES Function norm 6.278092520575e-13 nu = 0.0997308 3 SNES Function norm 1.354552924947e-14 Nonlinear solve converged due to CONVERGED_SNORM_RELATIVE iterations 3 starting loop: 2/5. Current nu = 0.0997308. Predicted nu = 0.0994608 0 SNES Function norm 2.081688271474e-03 nu = 0.0994608 1 SNES Function norm 5.759118650479e-09 nu = 0.0994141 2 SNES Function norm 1.409543210835e-14 Nonlinear solve converged due to CONVERGED_FNORM_RELATIVE iterations 2 starting loop: 3/5. Current nu = 0.0994141. Predicted nu = 0.0991441 0 SNES Function norm 1.900550339798e-07 nu = 0.0991441 Nonlinear solve did not converge due to DIVERGED_LINE_SEARCH iterations 0 starting loop: 4/5. Current nu = 0.0991441. Predicted nu = 0.0988741 0 SNES Function norm 3.809037608138e-07 ``` If you find where it starts to behave differently, I could try to make sense of this, but I cannot debug the full code right now, sorry. ---
**Author:** ![cmd](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/cmd/32/67_2.png) [@cmd](https://community.freefem.org/u/cmd)\ **Post date:** [December 19, 2022, 9:18pm UTC](https://community.freefem.org/t/snes-for-bordered-block-matrix/2193/24 "2022-12-19T21:18:39Z")
Interesting. Thank you. Appreciate your help. Will follow up if I find any clues. ---
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [December 20, 2022, 8:12am UTC](https://community.freefem.org/t/snes-for-bordered-block-matrix/2193/25 "2022-12-20T08:12:29Z")
Ah, OK, I see at least one simple error. You are calling `J(yq,yq)`, but `yq` is the result of a `KSPSolve()`, i.e., it’s in PETSc numbering. For `J(.,.)` you need to feed vectors in FreeFEM numbering. So, either: - replace `KSPSolve()` by `^-1` or add a call to `ChangeNumbering(..., inverse = true, exchange = true)`, or, - compute the scalar product yourself with a vector in PETSc numbering (`tmp = yq' * yq; mpiAllReduce(...);`). ---
**Author:** ![cmd](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/cmd/32/67_2.png) [@cmd](https://community.freefem.org/u/cmd)\ **Post date:** [December 20, 2022, 9:11am UTC](https://community.freefem.org/t/snes-for-bordered-block-matrix/2193/26 "2022-12-20T09:11:52Z")
Thank you for pointing out that issue. This does not solve the non-determinism I’m seeing, but I think it must be related to some similar bug I have involving a mix between PETSc and FreeFEM numbering. I noticed a couple other mix-ups as well, and I think I’m getting close to a working code. I should be able to figure it out from here, but thank you so much for all your help! 🙂 ---
**Author:** ![cmd](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/cmd/32/67_2.png) [@cmd](https://community.freefem.org/u/cmd)\ **Post date:** [February 7, 2023, 7:08pm UTC](https://community.freefem.org/t/snes-for-bordered-block-matrix/2193/27 "2023-02-07T19:08:56Z")
Hi there Pierre, I finally found the time to get back to this, and I got both the “simple” 2x2 block system (N+1)x(N+1) and the more complex 3x3 block system (N+N+1)x(N+N+1) to work. For once, I don’t have any problems or questions for you, and instead just wanted to report back to say thank you, and all is well! 🙂 Here’s the route I took: In **the 2x2 case** , it was relatively straighforward using `-pc_fieldsplit_type schur` as you showed me. Here we have a block structure like ```auto Mat JlPM(JlP.n, mpirank == 0 ? 1 : 0), yqPM(yqP.n, mpirank == 0 ? 1 : 0); Ja = [[J,JlPM],[yqPM',-1.0]]; ``` which can be solved using ```auto set(Ja, sparams = "-ksp_type preonly -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_precondition self ", setup = 1); set(J, sparams = " -fieldsplit_0_ksp_type preonly -fieldsplit_0_pc_type lu ", prefix = "fieldsplit_0_", parent = Ja); ``` However, **the 3x3 case** was a bit trickier. Eventually, I gave up on dealing with it directly as a 3x3 and did it like this: ```auto Mat JJ = [[J, 0],[H, J]]; Mat JlPM(2*JlP.n, mpirank == 0 ? 1 : 0), qmPM(2*qmP.n, mpirank == 0 ? 1 : 0); Mat JJa = [[JJ,JlPM],[qqmPM',0]]; ``` with ```auto func real[int] blockLU(real[int]& qqm) { PetscScalar[int] out1, out2, in1 = qqm(0:J.n-1), in2 = qqm(J.n:2*J.n-1); KSPSolve(J, in1, out1); MatMult(H, out1, out2); in2 -= out2; KSPSolve(J, in2, out2); PetscScalar[int] outPETSc = [out1, out2]; return outPETSc; } set(JJa, sparams = "-ksp_type preonly -pc_type fieldsplit -pc_fieldsplit_type schur -pc_fieldsplit_schur_precondition self ", setup = 1); set(JJ, sparams = "-fieldsplit_0_ksp_type preonly -fieldsplit_0_pc_type shell ", precon = blockLU, prefix = "fieldsplit_0_", parent = JJa); set(J, sparams = "-fieldsplit_ksp_type preonly -fieldsplit_pc_type lu", parent = JJ); ``` Thanks again! ---
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [February 7, 2023, 9:13pm UTC](https://community.freefem.org/t/snes-for-bordered-block-matrix/2193/28 "2023-02-07T21:13:59Z")
Eheh, a nest of nest is something completely valid which can be useful to breakdown complexity, as seen right here. Nicely done! [Previous page](https://community.freefem.org/t/snes-for-bordered-block-matrix/2193.md?page=1)