# Convergence issues with rate and state

**URL:** <https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143>\
**Category:** PyLith\
**Created:** [December 22, 2021, 10:51pm UTC](https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143 "2021-12-22T22:51:05Z")\
**Posts on this page:** 13\
**Page:** 1

<div class="post-metadata">

**Author:** ![LauraBM](https://avatars.discourse-cdn.com/v4/letter/l/3bc359/32.png) [@LauraBM](https://community.geodynamics.org/u/LauraBM)\
**Post date:** [December 22, 2021, 10:51pm UTC](https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143/1 "2021-12-22T22:51:05Z")

</div>

Hello,

I am trying to conduct a simulation with the rate and state model but I have convergence problems. I found in the v2.2.2 manual (section 6.4.5.1) that “The solution of these two linear systems gives the increment in slip assuming all the degrees of freedom except those imme-  
diately adjacent to the fault remain fixed. In real applications where the deformation associated with fault slip is localized  
around the fault, this provides good enough approximations so that the nonlinear solver converges quickly. In problems where  
deformation associated with slip on the fault is not localized (as in the case in some of the example problems), the increment in  
slip computed by solving these two linear systems is not a good approximation and the nonlinear solve requires a large number  
of iterations.”  
I do not know whether this could be the issue (deformation is not localized on the fault surface in my case). I have checked the examples and YouTube tutorials and applied the suggested settings:  
pc\_type = asm  
sub\_pc\_factor\_shift\_type = nonzero  
ksp\_rtol = 1.0e-30  
ksp\_atol = 1.0e-12  
snes\_rtol = 1.0e-30  
snes\_atol = 1.0e-10  
zero\_tolerance = 1.0e-11  
zero\_tolerance\_normal = 1.0e-11  
friction.linear\_slip\_rate = 1.0e-9  
friction\_pc\_type = asm

friction\_sub\_pc\_factor\_shift\_type = nonzero

friction\_ksp\_max\_it = 100

friction\_ksp\_gmres\_restart = 30

friction\_ksp\_error\_if\_not\_converged = true

friction\_ksp\_monitor = true

#friction\_ksp\_view = true

friction\_ksp\_converged\_reason = true

The linear solver takes thousands of iterations every time step, and eventually the non linear solver fails to converge. The model has about 30k elements.  
Are there any specific settings for cases where deformation is not localized in the fault plane?

Thanks,

Laura

---

<div class="post-metadata">

**Author:** ![baagaard](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/baagaard/32/615_2.png) [@baagaard](https://community.geodynamics.org/u/baagaard)\
**Post date:** [December 24, 2021, 12:59am UTC](https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143/3 "2021-12-24T00:59:36Z")

</div>

It looks like you are using the `asm` preconditioner. It is known not to work well for problems with faults. As discussed in PyLith v2.2.2 manual section 4.1 and the 2013 Aagard, Knepley, and Williams JGR paper, the custom fault preconditioner with algebriac multigrid preconditioning on the displacement block performs much better. See Table 4.4 with the corresponding settings provided in `share/settings/solver_fault_fieldsplit.cfg`.

---

<div class="post-metadata">

**Author:** ![LauraBM](https://avatars.discourse-cdn.com/v4/letter/l/3bc359/32.png) [@LauraBM](https://community.geodynamics.org/u/LauraBM)\
**Post date:** [December 26, 2021, 1:57pm UTC](https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143/4 "2021-12-26T13:57:39Z")

</div>

Thanks for your reply. I’ve tried those settings before and I also get convergence problems (at an early time step). Could you please explain what does fs\_\* stand for?. Also, is it OK to use friction\_\* settings combined with fs\_\* settings ? When I set the \*\_error\_if\_not\_converged & \* \_monitor to true, the friction solves converge within a few iterations, but SNES solver does not converge.  
Is there a way to display which node has the highest residual?

Thanks,

Laura

---

<div class="post-metadata">

**Author:** ![baagaard](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/baagaard/32/615_2.png) [@baagaard](https://community.geodynamics.org/u/baagaard)\
**Post date:** [December 27, 2021, 2:47am UTC](https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143/5 "2021-12-27T02:47:39Z")

</div>

`fs_` stands for field split.

PyLith v2.2.2 with friction uses two solvers. The “outer” solve is for the entire domain, while the “inner” solve is for the fault friction. The `friction_` solve settings apply only to the inner solve.

When trying to improve the performance, the first thing to do is to get the outer linear solve to converge reasonably fast. This should happen with the field split solver and custom fault preconditioner that I mentioned. The inner (friction) solve should always converge quickly with the ASM preconditioner.

Diagnosing failure of the nonlinear solve to converge can be tricky. You should turn on the SNES monitor to `snes_monitor = true` to see what the residual is doing. Is it decreasing but doing so very slowly, does it approach a constant value, or does it increase to a large value?

---

<div class="post-metadata">

**Author:** ![LauraBM](https://avatars.discourse-cdn.com/v4/letter/l/3bc359/32.png) [@LauraBM](https://community.geodynamics.org/u/LauraBM)\
**Post date:** [December 29, 2021, 10:41pm UTC](https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143/6 "2021-12-29T22:41:01Z")

</div>

Thanks for your answer. I’ve turned on the SNES monitor and the residual is very large the first time (0 SNES Function norm 4.877127797937e+05), then it gets stuck at about 5 e-7 until it eventually fails. I wonder why it fails at a number of iterations smaller than snes\_max\_it (message is Nonlinear solve did not converge due to DIVERGED\_FUNCTION\_COUNT iterations 258).  
I’ve set the cohesion to 100 MPa to make sure there is no slip on the fault, and the problem persists, so it seems the problem is outside the fault, but I don’t know where.

Laura

---

<div class="post-metadata">

**Author:** ![baagaard](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/baagaard/32/615_2.png) [@baagaard](https://community.geodynamics.org/u/baagaard)\
**Post date:** [December 30, 2021, 3:34am UTC](https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143/7 "2021-12-30T03:34:51Z")

</div>

If there is no slip on the fault and the bulk rheology is linear elastic, then the solution should converge in 1 SNES iteration. I would start by checking the mesh quality. Make sure there are no cells with condition numbers greater than about 2-5. I would then omit the fault from the simulation settings (you should change your preconditioner settings to simply use `pc_type = ml` and not use the field split and custom fault preconditioner settings). If the problem works without the fault, then put it back in and run a problem with a prescribed slip of zero (use the “fault” preconditioning settings). If prescribed slip works, then use a very large compressive normal traction like 1.0e+12 Pa and make sure the fault doesn’t slip and that the solution converges in one SNES iteration. If that works, slowly adjust the parameters to get to the desired ones. Hopefully, this will lead you to identify what is going wrong.

---

<div class="post-metadata">

**Author:** ![LauraBM](https://avatars.discourse-cdn.com/v4/letter/l/3bc359/32.png) [@LauraBM](https://community.geodynamics.org/u/LauraBM)\
**Post date:** [December 30, 2021, 11:58pm UTC](https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143/8 "2021-12-30T23:58:21Z")

</div>

Ok thanks, I will try that.

Laura

---

<div class="post-metadata">

**Author:** ![LauraBM](https://avatars.discourse-cdn.com/v4/letter/l/3bc359/32.png) [@LauraBM](https://community.geodynamics.org/u/LauraBM)\
**Post date:** [January 8, 2022, 12:00am UTC](https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143/9 "2022-01-08T00:00:46Z")

</div>

Hello,

I was wondering how much tight should the absolute tolerances of ksp and snes be when using rate and state. I know there are constraints between zero\_tolerance, snes\_atol and ksp\_atol (and also the linear slip rate), but can we just shift all the values?  
Also, the value of the reference slip rate V0 can significantly affect the value of the computed friction coefficient. How do you estimate it?

Thanks,

Laura

---

<div class="post-metadata">

**Author:** ![baagaard](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/baagaard/32/615_2.png) [@baagaard](https://community.geodynamics.org/u/baagaard)\
**Post date:** [January 10, 2022, 4:30pm UTC](https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143/10 "2022-01-10T16:30:37Z")

</div>

Yes, you can just shift/scale the KSP absolute tolerance, SNES absolute tolerance, and zero\_tolerance values. Note that these are all nondimensional values. That is, they apply directly to the internal values solver values that have been nondimensionalized.

The reference slip rate V0 in the rate-state friction model should be set according to the desired friction model values. The value should be chosen in accordance with the coefficient of friction and the characteristic slip distance.

---

<div class="post-metadata">

**Author:** ![LauraBM](https://avatars.discourse-cdn.com/v4/letter/l/3bc359/32.png) [@LauraBM](https://community.geodynamics.org/u/LauraBM)\
**Post date:** [January 11, 2022, 11:30pm UTC](https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143/11 "2022-01-11T23:30:58Z")

</div>

Thanks for your answer. In choosing the reference slip rate, should one also consider the expected values of the slip rate V (and thus the time step size), given the log dependence?  
As suggested, I checked the condition number of the mesh, it is between 1 and about 2.5. When I try the field split settings that you suggested, the simulation is extremely slow and it does not converge in a reasonable time (the KSP within SNES errates, the norms decrease and then increase and so on, globally decreasing but after thousands of iterations). The simulation runs with ASM. The conceptual model is basically a box with rollers on x=0, y=0 and zbottom, and Neuman boundary conditions on the other three boundaries. Gravity is on, and there is a quasi-vertical fault in the middle.

Thanks again,

Laura

---

<div class="post-metadata">

**Author:** ![baagaard](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/baagaard/32/615_2.png) [@baagaard](https://community.geodynamics.org/u/baagaard)\
**Post date:** [January 19, 2022, 6:44pm UTC](https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143/12 "2022-01-19T18:44:07Z")

</div>

If the KSP (linear) solver does not converge in a reasonable number of iterations (usually less than 100) with the field split settings, then something is wrong in the model setup. Is this true for prescribed slip?

---

<div class="post-metadata">

**Author:** ![LauraBM](https://avatars.discourse-cdn.com/v4/letter/l/3bc359/32.png) [@LauraBM](https://community.geodynamics.org/u/LauraBM)\
**Post date:** [February 10, 2022, 10:58pm UTC](https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143/13 "2022-02-10T22:58:17Z")

</div>

Yes, it also happens with prescribed slip

---

<div class="post-metadata">

**Author:** ![baagaard](https://yyz2.discourse-cdn.com/flex036/user_avatar/community.geodynamics.org/baagaard/32/615_2.png) [@baagaard](https://community.geodynamics.org/u/baagaard)\
**Post date:** [February 12, 2022, 2:39am UTC](https://community.geodynamics.org/t/convergence-issues-with-rate-and-state/2143/14 "2022-02-12T02:39:28Z")

</div>

This suggests there is something wrong or unusual in the model setup. If you provide all of the input files for prescribed slip, we can take a look as see if we can find the problem. Along with the input files, please include a diagram / sketch showing the problem you are trying to solve.
