[CIG-SHORT] convergence troubles: drucker-prager material with fault

[Brad Aagaard]

On 12/05/2013 03:52 PM, Eric Lindsey wrote:
> On Thu, Dec 5, 2013 at 8:55 AM, Matthew Knepley <knepley at mcs.anl.gov> wrote:
>
>> On Wed, Dec 4, 2013 at 11:23 PM, Eric Lindsey <elindsey at ucsd.edu> wrote:
>>
>>> Hi guys,
>>>
>>> I'm having some trouble getting the model with a fault in an
>>> elastoplastic rheology (plane strain) to converge. I'm pretty sure I've
>>> fixed all the boundary conditions to behave correctly now. The LU method
>>> fails entirely once yielding occurs, and the approximate method
>>> (multiplicative fieldsplit using ml and jacobi for the two domains) does
>>> not converge. More details:
>>>
>>> The first thing I tried was using (solver04.cfg) from the tutorial, which
>>> is a schur fieldsplit with LU and the custom Pylith preconditioner. In the
>>> case of a fault in a purely elastic material, this procedure is working
>>> fine. The nonlinear solve does take over 100 iterations, I'm not sure if
>>> this is just a function of my mesh; anyway I'm not too worried about
>>> optimization yet.
>>>
>>
>> Is this a through going fault? The linear solve is exact, so we would not
>> expect this.
>>
>
> The fault does extend all the way through the domain; there are traction
> boundary conditions on those sides of the domain that match the stresses on
> the fault. The fault itself has static friction -- so I think some
> nonlinear iterations are required to solve for the slip that satisfies the
> mohr-coulomb condition on the fault. The linear solves do finish in 1
> iteration, but I think the Jacobi preconditioner for the fault
> (fieldsplit_1) is pretty weak, which explains the number of nonlinear
> iterations. Unless I've misunderstood this.

When the fault is not slipping and the domain is elastic, the number of 
SNES iterations should be 1. Once the fault starts slipping, then for a 
through-going fault, the number of iterations will be significant 
because the slip/traction estimate is setup to converge quickly for a 
fault that is not through-going.

>>> In the plastic case, whenever stresses near the bends in the fault begin
>>> to exceed the yield criterion (I increase them slowly to this value over
>>> several time steps), I get a Zero Pivot error. But maybe this is not
>>> unexpected for a plastic rheology?
>>>
>>
>> Using this option
>>
>> fs_fieldsplit_0_pc_factor_shift_type = nonzero
>>
>> will prevent a 0-pviot, but the preconditioner becomes weaker.
>>
>
>>     Matt
>>
>
> Thanks - I think this helped... It did prevent the zero pivot in the
> plastic case, but now one of the inner linear solves fails to converge
> (after 10,000 iterations), eventually leading to Pylith giving up with the
> message:
>
> WARNING! Fault opening with nonzero traction., v_fault: 175, opening:
> 1.41881e-09, normal traction: -6103.07
>
> The normal traction is still negative and this doesn't happen in the
> elastic case... I think this must be a numerical issue. So there is still
> an issue with my solver once plastic yielding starts, I think.

You shouldn't need to allow so many linear iterations. Notice that the 
opening is small. This suggests that your solver tolerances are not 
setup properly. Please refer to the discussion of solver tolerances for 
simulations with friction in the June tutorial (Session V, slides 9-10). 
You may be able to loosen the absolute tolerances a couple orders of 
magnitude from those mentioned in the slides (start with the elastic case).

Brad