Dear Sir,
I would like to ask a question regarding my use of poroelastic simulations in Pylith 4.1.3. When I use the no_faultzone_permeability3.spatialdb file for all materials, the simulation runs successfully. However, when I use the no_faultzone_permeability.spatialdb file for all materials, an error occurs, displaying the following message:
[0]: --------------------- Error Message --------------------------------------------------------------
[0]PETSC ERROR: Residual norm computed by GMRES recursion formula 16442.8 is far from the computed residual norm 41861.1 at restart, residual norm at start of cycle 97030.1
Could you please advise on how to resolve this issue? I found a reference in the manual mentioning an “edge” issue, but my simulation only involves blocks.
Additionally, I would like to know how to output stress, pressure, and other values on the fault surface in Pylith 4.1.3. The syntax seems to differ significantly from that in Pylith 2, and I have not yet found any relevant examples. I apologize for the inconvenience.
file.zip (1.6 MB)
Dear Baagaard, I encountered some issues while trying to export the values on the fault surface. Could you provide some guidance or reference?
Try using these solver settings (these will be the default in v4.2.0) for poroelasticity without a fault:
[pylithapp.petsc]
pc_type = fieldsplit
pc_fieldsplit_type = multiplicative
pc_fieldsplit_0_fields = 2
pc_fieldsplit_1_fields = 1
pc_fieldsplit_2_fields = 0
fieldsplit_trace_strain_pc_type = bjacobi
fieldsplit_pressure_pc_type = bjacobi
fieldsplit_displacement_pc_type = ml
fieldsplit_displacement_ksp_type = gmres
Information available on the fault surface is limited to values that are available on the fault. When prescribing slip, this includes the slip and the change in tractions. Note that PyLith v4.1.3 contains a fault orientation bug the gives erroneous results for the tractions. Refer to Fault orientation bug in PyLith v4.1.3 (and v4.0.0-v4.1.1) for more information.
If you want the traction and fluid pressure at the location where a fault is (without having the fault), then use the material output and compute the fault traction from the stress tensor, making use of the fault orientation. This will need to be done as part of the post processing, such as in a Python script.
Thank you very much. However, when I try to set up the fault surface and output the results, I encounter some issues. I’m not quite sure how to configure the fault surface ID in PyLith v4.1.3, and I don’t understand the differences in output on the surface compared to PyLith2. Are there any specific reference examples available?
And I’m still not quite sure why replacing no_faultzone_permeability.spatialdb results in the above PETSC ERROR. The solver I’m currently using is the same as the one you provided earlier.
Fault output
For information on specifying faults in PyLith v3 and later, I recommend the PyLith tutorials on YoutTube (for example, https://www.youtube.com/watch?v=fcHKL8PqFDk&list=PLNTb1ey7FFk7aaP7meQlevDqahYa-7S1h&index=3) and Horizontal Cross-Section of Strike-Slip Fault (2D) — PyLith 4.2.0dev documentation.
In PyLith v3 and later, in most cases you do not need to do as much for output. In examples/strikeslip-2d, we specify the fields we want output (we can select from slip, traction change, and any solution fields) in pylithapp.cfg:
# Output `slip` and the change in tractions on the fault.
observers.observer.data_fields = [slip, traction_change]
Poroelasticity solver error
Changing the material properties may be creating an ill-conditioned system due to the inconsistency between the material properties and the scales used for nondimensionalization. For a well-conditioned system for poroelasticity, the ratio of the nondimensional permeability to viscosity should be about 1. Assuming the stress and length scales are appropriate for your problem, then you would adjust the time scale so that the nondimensional permeability to viscosity ratio is about 1.
Thank you very much for your reply!
Dear Sir,
I have a question. I want to distinguish between the gravitational loading of the water body and the stress changes caused by water diffusion in the fluid diffusion process. Currently, my idea is to set the diffusion model to 1e-30 or even smaller to approximate a model with no diffusion. However, during the simulation, I keep encountering errors like:
50 SNES Function norm 2.712968822457e-05 Nonlinear solve did not converge due to DIVERGED_MAX_IT iterations 50, or 0 TS dt 798. time 0. 0 SNES Function norm 3.433142042464e-04 Linear solve converged due to CONVERGED_ATOL iterations 48 1 SNES Function norm 2.649311676138e-05 Linear solve did not converge due to DIVERGED_DTOL iterations 0
How should I solve this issue?
Or is there a problem with this idea? Do you think this theory is feasible? Since I need to compare it with the previous diffusion model, I want to keep everything else unchanged and only modify the diffusion coefficient
I suggest looking at the governing equations for poroelasticity and the terms in the Jacobian. You may be introducing very large terms in the Jacobian, resulting in an ill-conditioned system. Rather than adjusting material properties to change a poroelastic problem into an elastic problem, it is better to change the material (governing equations) to get the behavior that you want.
I understand that you mean replacing the poroelastic material with a purely elastic material. What should I do to ensure that everything remains the same as before, except for the diffusion coefficient? The spatial file for the elastic material is significantly different from the one for poroelasticity, and I wonder if there will be major changes in the solving process and other settings as well. Besides, is there a simple replacement method?
Even when I change it to elastic, the error still occurs, and I cannot find a solution in the manual. I would like to ask if the issue is due to an incorrect solver, or if there is a problem elsewhere. Below is my file.
file.zip (6.1 MB)
The error message shows:
[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------
[0]PETSC ERROR: KSPSolve has not converged due to Nan or Inf norm
[0]PETSC ERROR: See FAQ — PETSc 3.22.3 documentation for trouble shooting.
[0]PETSC ERROR: Petsc Development GIT revision: v3.22.2-268-g8f3ed387c6c GIT Date: 2025-01-07 15:38:51 +0000
[0]PETSC ERROR: /public1/home/a8s000279/pylith-4.2.0-linux-x86_64/bin/mpinemesis with 50 MPI process(es) and PETSC_ARCH on m4cm2104.para.bscc by a8s000279 Wed Feb 19 17:24:51 2025
[0]PETSC ERROR: Configure options: --prefix=/opt/pylith/dist --with-c2html=0 --with-x=0 --with-clanguage=C --with-mpicompilers=1 --with-shared-libraries=1 --with-64-bit-points=1 --with-large-file-io=1 --download-parmetis=1 --download-metis=1 --download-f2cblaslapack=1 --download-ml --with-fc=0 --with-hwloc=0 --with-ssl=0 --with-x=0 --with-c2html=0 --with-lgrind=0 --with-hdf5=1 --with-zlib=1 --LIBS=-lz --with-debugging=0 --with-fc=0 CFLAGS=“-g -O2” CXXFLAGS=“-g -O2 -DMPICH_IGNORE_CXX_SEEK” FCFLAGS= CPPFLAGS=“-I/opt/pylith/dist/include -I/opt/pylith/dist/include” LDFLAGS=“-L/opt/pylith/dist/lib -L/opt/pylith/dist/lib64 -L/opt/pylith/dist/lib” PETSC_DIR=/opt/pylith/build/cig/petsc-pylith PETSC_ARCH=arch-pylith
[0]PETSC ERROR: #1 KSPGMRESCycle() at /opt/pylith/build/cig/petsc-pylith/src/ksp/ksp/impls/gmres/gmres.c:111
[0]PETSC ERROR: #2 KSPSolve_GMRES() at /opt/pylith/build/cig/petsc-pylith/src/ksp/ksp/impls/gmres/gmres.c:228
[0]PETSC ERROR: #3 KSPSolve_Private() at /opt/pylith/build/cig/petsc-pylith/src/ksp/ksp/interface/itfunc.c:900
[0]PETSC ERROR: #4 KSPSolve() at /opt/pylith/build/cig/petsc-pylith/src/ksp/ksp/interface/itfunc.c:1075
[0]PETSC ERROR: #5 KSPSetUp_Chebyshev() at /opt/pylith/build/cig/petsc-pylith/src/ksp/ksp/impls/cheby/cheby.c:823
[0]PETSC ERROR: #6 KSPSetUp() at /opt/pylith/build/cig/petsc-pylith/src/ksp/ksp/interface/itfunc.c:377
[0]PETSC ERROR: #7 PCSetUp_MG() at /opt/pylith/build/cig/petsc-pylith/src/ksp/pc/impls/mg/mg.c:1151
[0]PETSC ERROR: #8 PCSetUp_GAMG() at /opt/pylith/build/cig/petsc-pylith/src/ksp/pc/impls/gamg/gamg.c:943
[0]PETSC ERROR: #9 PCSetUp() at /opt/pylith/build/cig/petsc-pylith/src/ksp/pc/interface/precon.c:1071
[0]PETSC ERROR: #10 KSPSetUp() at /opt/pylith/build/cig/petsc-pylith/src/ksp/ksp/interface/itfunc.c:415
[0]PETSC ERROR: #11 KSPSolve_Private() at /opt/pylith/build/cig/petsc-pylith/src/ksp/ksp/interface/itfunc.c:826
[0]PETSC ERROR: #12 KSPSolve() at /opt/pylith/build/cig/petsc-pylith/src/ksp/ksp/interface/itfunc.c:1075
[0]PETSC ERROR: #13 SNESSolve_NEWTONLS() at /opt/pylith/build/cig/petsc-pylith/src/snes/impls/ls/ls.c:221
[0]PETSC ERROR: #14 SNESSolve() at /opt/pylith/build/cig/petsc-pylith/src/snes/interface/snes.c:4834
[0]PETSC ERROR: #15 TSTheta_SNESSolve() at /opt/pylith/build/cig/petsc-pylith/src/ts/impls/implicit/theta/theta.c:174
[0]PETSC ERROR: #16 TSStep_Theta() at /opt/pylith/build/cig/petsc-pylith/src/ts/impls/implicit/theta/theta.c:225
[0]PETSC ERROR: #17 TSStep() at /opt/pylith/build/cig/petsc-pylith/src/ts/interface/ts.c:3426
[0]PETSC ERROR: #18 TSSolve() at /opt/pylith/build/cig/petsc-pylith/src/ts/interface/ts.c:4072
[0]PETSC ERROR: #19 void pylith::problems::TimeDependent::solve()() at …/…/…/pylith-4.2.0/libsrc/pylith/problems/TimeDependent.cc:487
Fatal error. Calling MPI_Abort() to abort PyLith application.
Traceback (most recent call last):
File “/public1/home/a8s000279/pylith-4.2.0-linux-x86_64/lib/python3.12/site-packages/pylith/apps/PetscApplication.py”, line 55, in onComputeNodes
self.main(*args, **kwds)
File “/public1/home/a8s000279/pylith-4.2.0-linux-x86_64/lib/python3.12/site-packages/pylith/apps/PyLithApp.py”, line 114, in main
self.problem.run(self)
File “/public1/home/a8s000279/pylith-4.2.0-linux-x86_64/lib/python3.12/site-packages/pylith/problems/TimeDependent.py”, line 134, in run
ModuleTimeDependent.solve(self)
File “/public1/home/a8s000279/pylith-4.2.0-linux-x86_64/lib/python3.12/site-packages/pylith/problems/problems.py”, line 217, in solve
return _problems.TimeDependent_solve(self)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
RuntimeError: Error detected while in PETSc function.
Abort(-1) on node 28 (rank 28 in comm 0): application called MPI_Abort(MPI_COMM_WORLD, -1) - process 28
[server]: PMIU_parse_keyvals: unexpected key delimiter at character 48 in cmd
/public1/home/a8s000279/pylith-4.2.0-linux-x86_64/bin/nemesis: mpirun: exit 255
/public1/home/a8s000279/pylith-4.2.0-linux-x86_64/bin/pylith: /public1/home/a8s000279/pylith-4.2.0-linux-x86_64/bin/nemesis: exit 1
I suggest setting up the elastic and poroelastic simulations so that they use the same spatial database for materials. A spatial database can have values that are not used, so it can include both the values needed for a linear elastic rheology (Vp, Vs, and density) as well as the values needed for a linear poroelastic rheology. This will make it easier for them to be consistent.
You will need to change the solution field (displacement only), materials, and use only the BCs related to the displacement field.
Sorry, sir, I still don’t quite understand what you mean. My understanding is that in the spatial data file, we can freely combine and configure the parameters of linear elastic or linear poroelastic models, right? For example, if I want to distinguish between the effects of hydraulic loading and hydraulic diffusion, I can modify the original poroelastic model by removing the diffusion coefficient and other related parameters, then run the simulation. Is that correct?
Additionally, I would like to know the cause of the error I encountered. Also, is it possible to include a viscosity coefficient in the spatial data file, meaning a viscoelastic poroelastic model?
You can have as many fields as you want in a spatial database file. PyLith knows what fields it needs for each material and will only query for those. Any additional fields are ignored. This means you can put fields for both linear elasticity and poroelasticity in the same spatial database. This can make it easier to keep the values between the two materials consistent (density, shear modulus, and bulk modulus for the solid; the elastic model parameterizes these using Vp, Vs, and density whereas the poroelastic model uses different fields).
If the solver is not converging, then there is likely an error in setting up the problem. The elastic simulation should be easiest to setup as there are fewer parameters. Consequently, I would start with that one. The solution field should just have the displacement field (default), the material should be Elasticity with an IsotropicLinearElastic bulk rheology (both the default), and BCs should be DirichletTimeDependent or NeumannTimeDependent on the displacement field (default).
I will try to find time to look in detail at your files next week.
Thank you very much for your help, and I look forward to your reply again!
Sir, I would like to ask whether it is possible to remove the diffusion coefficient and add parameters such as Vp, Vs, and viscosity in the poroelastic model’s spatial database file. If the diffusion coefficient is removed, will significant modifications to the pylithapp.cfg file also be required?