Problem when setting the initial stress of material

Dear Brad,
I had some problem when setting the initial stress of material. I set the initial stress condition by the following commands in pylithapp.cfg.

db_initial_stress = spatialdata.spatialdb.UniformDB
db_initial_stress.label = Initial stress tensor
db_initial_stress.values = [stress-xx,stress-yy,stress-zz,stress-xy,stress-yz,stress-xz] = [-52.875MPa,-38.6250MPa,-40.0MPa,-12.34MPa,0.0MPa,0.0MPa]

The simulation wouldn’t continue after several time steps, and no error was reported. The interface stopped at the first few time steps, as shown in the picture. After removing the above commands, the program will run properly. Could you please tell me where the problem is?

Please provide more information.

  • Diagram of the boundary value problem you are trying to solve showing the geometry of the domain, the boundary conditions, and materials. If it is a time-dependent problem, what is varying in time?
  • Version information for PyLith (for example, the output of pylith --version)

Hi, Brad. Thanks for helping. I have attached the pylithapp.cfg and material properties files. It is a dynamic simulation, and my Pylith version is v2.2.2. (118.7 KB)

Thanks for including the diagrams. They always help!

I suspect there is an issue with fault opening based on the nonplanar fault geometry. On some systems, stdout and stderr get truncated and a simulation ends before it produces the error message. This is more likely when running in parallel. Sometimes just capturing stdout and stderr to a log will enable you to see the error message.

See Nonzero traction error in spontaneous rupture model - #10 by baagaard for a discussion of setting zero tolerances for nonplanar fault geometry.

Thanks for answering. It is the problem about fault opening. So I just need to set a higher value of zero_tolerance_normal? Would there be any problem if I set the value too big? I found that I have to use a value of 5e-3 to make the simulation continue.

The zero_tolerance_normal value will depend on the smoothness of the surface. If you need a relatively large value, it is an indication that you are not accurately resolving the geometry and stress variations. Your fault surface is relatively rough given your cell size. You could probably use a smaller value if the cell size was smaller.

Hi, Brad. I’m using SlipWeakeningTimeStable friction law in my simulation, the nucleation area is controlled by time-weakening law and other parts of fault are controlled by slip-weakening law. In my understanding, the influence of increasing Dc is to increase the fracture energy and thus decrease the rupture speed for the same stress and frictional conditions. But as I increased the Dc in my simulations, the slip of the fault has increased and the slip rate changed little.

Here is slip result of Dc = 0.3m

Dc = 0.5m

All the other conditions are the same. The result I expected was a smaller amount of slip and slip rate. Do you know why that is?

With a large slip-weakening parameter (Dc) you may be resolving the cohesive zone and resolving the rupture front more accurately. This can lead to larger slip.

Hi, Brad. You mean the results with large Dc are more accurate? The slip of the fault can reach as much as 7m when Dc = 0.5m. I have checked my initial stress and friction conditions, the static coefficient and dynamic coefficient are 0.5 and 0.2 respectively. The relative prestress ratio is around 0.5. I think the the parameters are reasonable. Whether I increase the static coefficient or reduce the stress drop slightly, the rupture wouldn’t propagate. The sudden change from large slip to wouldn’t rupture at all really confused me. Could you please tell me how to reduce the slip?

Dc can affect the accuracy of the results. You need to make sure you are resolving the cohesive zone. Once you have an accurate solution, then you can investigate the sensitivity of the rupture to the fault constitutive parameters. Changing Dc can affect the rupture in a variety of ways as can changing the coefficient of friction. For example, you can keep the fracture energy nearly the same by changing both simultaneously.

Hi, Brad. I read some articles about cohesive zone. I calculated the cohesive zone size of my simulation based on the following equation,
The result is 633m. The cell size of fault surface of my model is 200m. That means the cohesive zone was resolved with about three cell size, is it enough for the finite element method?

A cohesive zone of 3-5 times the cell size is usually okay. See the recent SCEC Dynamic Rupture benchmarking papers and results (SCEC/USGS Spontaneous Rupture Code Verification Project) for additional information.

Thanks a lot for helping.

Hi, Brad. I investigated the influence of the constitutive parameters on rupture according to your suggestion. I changed the Dc and static friction coefficient simultaneously with the fracture energy remaining the same. But I found the results are different. For example, I first set Dc to 0.5, static friction coefficient to 0.35, and dynamic friction coefficient to 0.2. The rupture can propagate on fault. Then I changed the Dc to 0.3 and static friction coefficient to 0.45. According to the equation fracture energy = 1/2 * (static stress drop - dynamic stress drop) * Dc, the fracture energy didn’t change. but the latter simulation wouldn’t rupture at all. Both the cohesive zones of two simulations are well resolved and all the other parameters are the same. I found that the simulations with higher static friction coefficient are harder to rupture even though the fracture energy is the same, do you know why that is?

The relationship between the fault constitutive parameters and rupture behavior is quite complex even for slip-weakening friction. In addition to the fracture energy, the strength excess is a parameter that people often examine as well. These issues are discussed in numerous journal papers over the past 30 years. A good starting point is a paper by Martin Galis and others.