I have modified step13 in hex8 (slip-weakening fault friction in a quasistatic problem) to build a spontaneous rupture model for my own fault. I have encountered an error regarding the tractions on the fault. The fault is planar with a fine mesh. Here’s what pylith is returning:
“terminating with uncaught exception of type std::runtime_error: WARNING! Fault opening with nonzero traction., v_fault: 1735, opening: 9.60479e-07, normal traction: -1.05722e-05”
I tried attaching the json file but I get an error saying “new users cannot upload attachments”. I can provide it via email if that would help solve the problem.
You need to adjust the solver tolerances. See the discussion of solver tolerances in the slides of the 2017 CDM tutorial 2-D Subduction Zone: Fault Friction.
I have tried using the tolerances recommended in the slides but continue to get the same error. I have also tried increasing the tolerances while respecting the suggestions in the tutorial (keeping ksp_atol and snes atol larger than zero_tolerance) with no luck. Increasing the tolerance of the zero_tolerance_normal to 10^-3 made it crash much later, yet it still ultimately brought up the same error.
Solver tolerances of 10^-3 are almost surely too large for a reasonable solution. Can you post the .cfg files you are using when the solver tolerances are around 10^-8?
Note: You should keep the ksp_rtol and snes_rtol values at 10^-20
I have changed the tolerances to your recommended values (10^-20) and keep getting the same error. I get it with at least three different scripts. I attach the .cfg file for one of them (modified after step13 in hex8). I also attach the pylith.cfg file.
You are using an absolute tolerance of 1.0e-20 for the KSP (linear) solve. This should be 1-2 orders of magnitude smaller than the SNES absolute tolerance (1.0e-9). So try something like ksp_atol = 1.0e-11.
The normal traction and opening are still quite small, so I would try increasing the zero_tolerance_normal to try to keep the fault closed. Sometimes small asymmetries result in very small amounts of fault opening that are difficult to avoid.
See examples/2d/subduction/step06.cfg for an example of how to set zero_tolerance_normal.
I would use something like
[pylithapp.timedependent.interfaces.fault]
# Zero tolerance for locking fault (no slip). Must be larger than ksp tolerance.
zero_tolerance = 1.0e-8
zero_tolerance_normal = 1.0e-6
I followed the script on step06 and tried the values you suggested and others (while keeping values larger than ksp tolerance) but saw no difference. The updated log error is attached.
Are the asymmetries a problem from the mesh? I have run other problems (both coseismic only and including post-seismic relaxation) on this mesh with no issues.
Brad took a look at it and was able to figure out the issue. The buried edges of the fault were not specified and that is what was causing the ‘opening tractions’ error.
I have similar error (WARNING! Fault opening with nonzero traction), in modeling non-planar fault. I am interested in the discussion and in the solutions. But I don’t understand clearly with the last solution (related to the buried edges of the fault).
I would like to ask, what do you mean by “the buried edges of the fault were not specified” in this case (example Step 13)?
When you have a nonplanar fault, slip will induce fault opening because the geometry of the fault surface is incompatible with slip along that surface. The workaround is to increase the zero tolerance for fault opening so that PyLith treats a small amount of opening as zero. If you need to keep increasing the zero_tolerance_normal to larger and larger values, it means your fault is too jagged and you need to smooth the fault geometry, refine your mesh, or both.
[pylithapp.timedependent.interfaces.fault]
# Zero tolerance for locking fault (no slip). Must be larger than ksp tolerance.
zero_tolerance = 1.0e-8
zero_tolerance_normal = 1.0e-6
When I opened the results,
(1) I saw the fault folded at the tips (blue circle top fig), which I did not modeled. Is it ok to have such folded fault?
(2) I exaggerated the block’s displacement to 10000, an opening is still observed on some part of the fault (green circle bottom fig). How do I check the opening value at this portion of fault?
(3) I wonder, can I check/ confirm whether the fault is already slip or not yet, from the traction values on the fault? How do I check it? Kindly advice.
Figures are attached below. The fault specifications are: cohesion 2 MPa, friction 0.5, running 0 to 200 year, axial and shear displacement boundary condition).
The faults getting “folded” at the tips looks like it is associated with not specifying the buries edges of the fault. This is necessary to prevent the buried edges of the fault extending arbitrarily one face into the domain. Refer to Figures 6.3 and 6.4 in the PyLith v2.2.2 manual for diagrams related to this issue. A discussion of specifying buried edges in a .cfg is included in Section 7.18.7 (excerpt from step03.cfg).
A nonplanar fault will generate opening. The amount of opening depends on the curvature and amount of slip. The zero tolerance for opening controls how much opening is considered to be sufficient to create a free surface, that is the two sides are no longer in contact. If the opening is less than the zero tolerance, the code treats the two sides as being in contact and applies the associated physics. If the opening is greater than the zero tolerance, the code treats the two sides as not being in contact and they become free surfaces.
The best way to monitor the amount of slip is to look at the components of slip vector. You can look at the traction components at the same time steps to associate the slip with the effect on the tractions.
I have tried with buried edges, and the fault looks good (no folded edge). Thank you.
I found some interesting anomaly feature (shown in circle in the top figure) that only exist in certain friction and cohesion parameter. However, If I ran with a new mesh (finer mesh) the feature diverge into two patches of anomaly which geologically hard to explain (shown in circle in the bottom figure).
I would like to ask, how I differentiate whether the anomaly is due to deformation result or due to an artefact because of parameter modeling/ mesh? Kindly advice. Thank you.
