Subduction 2D Kinematic Prescription and errors with new mesh

Hi all, I have two questions regarding a modified subduction-2D model.

I want to impose a slip ramp function to create a slow slip event where creep starts at an initiation time and stops after a given time. Can the KinSrcRamp kinematic prescription accommodate spatially varying slip along the interface? So far, I have only been able to make this kinematic prescription work by applying a uniform ramp along my slipping 2D surface, but not a spatially varying slip (for example, the slab interface is locked until 410km and from 410-460km depth, I am applying the slip ramp prescription).

I am also working on modifying the original slab in the PyLith example, taking inspiration from the reverse-2d example, where a fault and splay fault cut through a portion of the model domain, but the line extends to the model domain boundary. I am trying to do something similar, where I have a line (fault) cutting the entire slab but only a portion of the fault slip, and I receive a zero pivot in LU factorization error (see error file in zip folder). I did look at this previous discussion post, and it relates to the subduction-3d example. I don’t believe my error pertains to fault orientation. I have attached the necessary files to run this intraslab fault example with the fault that cuts the entire slab as one line (mesh_tri_lithos_intraslab.msh) and the new fault geometry where there is an extra point midway through the slab representing the slipping surface for comparison (mesh_tri_lithos_intraslab60.msh). Thanks for the guidance!
PyLith Discussion forum files.zip (1015.0 KB)

Spatially varying slip

I didn’t see any spatial database in the zip file with your slip distribution. For a simple depth dependence to the slip, create a SimpleDB with points spaced over the depth. This corresponds to data-dim=1 (variation in the field along a line). For KinSrcRamp in 2D, the spatial database would have values for initiation_time, rise_time, impulse_duration, final_slip_opening (zero), and final_slip_left_lateral (this one would vary with depth). Make sure the header information agrees with the data in terms of the number of points, spatial dimension, etc.

Zero pivot error

I am assuming the zero pivot error is from when you use mesh_tri_lithos_intraslab60.msh. The fault_edge physical group for mesh_tri_lithos_intraslab60.msh only has a point at the bottom of the fault. Because the fault is completely buried, you need a point from each end in the fault_edge. You do have this for mesh_tri_lithos_intraslab.msh.

Spatially varying slip
I didn’t include the files as I was first curious whether spatially varying slip was possible. I have included the files to run this example. I did include the correct header information and the other recommended components for the spatial database files and I receive this SNESSolve error (see text file). I was able to make the spatial database files work when using a uniform slip value across all depths. Thanks in advance for the help!
KinSrcRamp_files.zip (15.3 KB)

Zero pivot error
I included both points in my fault_edge and that solved the issue, thanks! I didn’t realize this fault is buried but it makes sense now.

Tip for all users seeking help: When sending parameter and spatial database files, remove settings that are commented out and check to make sure all comments are consistent with the current simulation. Comments that describe what you are trying to do in the parameter files help make your intention clear and make it easier to verify that the parameters match your intention.

I don’t see any problems in your spatial database files. I suspect that the issue is the nondimensionalization. The default time scale is 1 year and the default length scale is 1 km. As a result, imposing 0.04 m of total slip over 155 days is going to give very small displacements starting out. The length scale is appropriate for the size of the cells and domain. I would reduce the time scale for nondimensionalization to 1 day. You will also likely need to use smaller KSP and SNES tolerances because the deformation is very small.

Some recommendations:

  1. Start with a KinSrcStep time function and impose the variable slip with values on the order of 1 m for the maximum slip. This will make it easier to debug the spatial variation in slip.
  2. Update your KinSrcRamp settings so that the slip is about 0.1 m per time step. This should be okay with the default time scale and current solver tolerances.
  3. Update the default time scale and solver tolerances.

Settings for the time scale and solver tolerances:

[pylithapp.problem]
normalizer = spatialdata.units.NondimElasticQuasistatic
normalizer.relaxation_time = 1*day

[pylithapp.petsc]
# The absolute tolerances (atol) may need adjustments depending on the magnitude of the displacement.
# These atol values are 3 orders of magnitude smaller than the defaults.
ksp_rtol = 1.0e-20
ksp_atol = 1.0e-15
snes_rtol = 1.0e-20
snes_atol = 1.0e-12

Thanks for the recommendations Brad!

  1. I successfully ran this model with KinSrcStep, so spatial variation in slip does work.
  2. With this slip constrain in mind, I increased the amount of slip to about 30m total and I still receive the same SNESSolve error: SNESSolve has not converged due to Nan or Inf norm... I even increased the slip to hundreds of meters with the same error.
  3. I updated the default time and solver tolerances to the above for a final slip of 30m but the error still persists. I increased the tolerances again (below) just in case but no change. I’ve included the updated files below as well.
ksp_rtol = 1.0e-30
ksp_atol = 1.0e-20
snes_rtol = 1.0e-30
snes_atol = 1.0e-20

KinSrcRamp_files.zip (14.7 KB)

Please provide the mesh_tri_lithos_pull.msh file and any other files required to run the simulation.

Sorry about that, I’ve included the mesh now and other required files, thanks!
KinSrcRamp_files.zip (514.8 KB)

It looks like there is an error in the KinSrcRamp implementation. It has not been used as much as KinSrcStep or KinSrcTimeHistory. I will dig deeper sometime in the next few days. You should be able to switch to KinSrcTimeHistory and create the behavior you want.

The KinSrcRamp implementation for zero final slip results in a divide by zero, which leads to NaN in the residual. I have fixed the issue in the PyLith development branch (main). The bugfix will be included in the next release. In the meantime, the workaround is to use the KinSrcTimeHistory implementation.

Thanks for the tip, I tried implementing the KinSrcTimeHistory prescription, and I am receiving a querying error after 155 days of prescribed slip. Is there an explicit way to tell PyLith that there shouldn’t be any more slipping after the time history file period?

KinSrcTimeHistory.zip (16.4 KB)

Just add a point at the end of the time history with the same amount of slip. PyLith will interpolate between the points.

#TIME HISTORY ascii
TimeHistory {
  num-points = 3
  time-units = day
}
  0.0  0.0
155.0  1.0
999.0  1.0