Dear Pylith developers,

I’m modeling a fault in a 2D (plane strain) elastic block using FaultCohesiveDyn. I have a horizontal planar fault with constant friction everywhere.

I prescribe initial shear and normal stresses in the domain and then load the top boundary with a shear stress rate of 0.5 MPa per year. I want the simulation to run for 4.6 years and have the slip profiles for several timesteps in-between.

total_time = 4.6*year
max_dt = 0.2*year

The simulation diverges when the slip on the fault is supposed to start.

Here are my settings. I can’t attach the log file, it is too heavy, but I copy the final error. I’ve tried with and without the fault specific preconditioner, with the preconditioned SNES diverges at first timestep. If I decrease the timestep the simulation runs for longer and I get slip along the fault, but it still takes a lot of iterations to converge and I run out of cluster time.

I though this is supposed to be a very basic problem. Is there something very wrong with my settings? Please advise what’s best to do.

# ----------------------------------------------------------------------

# faults

# ----------------------------------------------------------------------

[pylithapp.timedependent.interfaces]

# Change fault to dynamic fault interface.

fault = pylith.faults.FaultCohesiveDyn

[pylithapp.timedependent.interfaces.fault]

zero_tolerance = 1.0e-8

# The label corresponds to the name of the nodeset in CUBIT.

label = fault

edge = fault_edge

# We must define the quadrature information for fault cells.

# The fault cells are 1D (surface).

quadrature.cell = pylith.feassemble.FIATSimplex

quadrature.cell.dimension = 1

# Use rate- and state-friction with the ageing law for evolution of

# the state variable.

friction = pylith.friction.contrib.DoubleSlipWeakeningFrictionNoHeal

friction.label = DoubleSlipWeakeningFrictionNoHeal

friction.db_properties = spatialdata.spatialdb.SimpleDB

friction.db_properties.label = Rate State Ageing

friction.db_properties.iohandler.filename = spatialdb/1_friction_DSW_single_nucl_1m_3m_fault.spatialdb

# ----------------------------------------------------------------------

# PETSc settings

# ----------------------------------------------------------------------

# NOTE: There are additional settings specific to fault friction.

#[pylithapp.timedependent.formulation]

#split_fields = True

#use_custom_constraint_pc = True

#matrix_type = aij

[pylithapp.petsc]

#fs_pc_type = fieldsplit

#fs_pc_use_amat = true

#fs_pc_fieldsplit_type = multiplicative

#fs_fieldsplit_displacement_pc_type = ml

#fs_fieldsplit_lagrange_multiplier_pc_type = jacobi

#fs_fieldsplit_displacement_ksp_type = preonly

#fs_fieldsplit_lagrange_multiplier_ksp_type = preonly

malloc_dump =

# Preconditioner settings.

pc_type = asm

sub_pc_factor_shift_type = nonzero

#pc_type = lu

#pc_factor_mat_solver_type = umfpack

# Convergence parameters.

ksp_rtol = 1.0e-12

ksp_atol = 5.0e-9

ksp_max_it = 3000

ksp_gmres_restart = 50

# Linear solver monitoring options.

ksp_monitor = true

#ksp_view = true

ksp_converged_reason = true

ksp_error_if_not_converged = true

# Nonlinear convergence parameters.

snes_rtol = 1.0e-12

snes_atol = 1.0e-7

snes_max_it = 3000

# Nonlinear solver monitoring options.

snes_monitor = true

snes_linesearch_monitor = true

#snes_view = true

snes_converged_reason = true

snes_error_if_not_converged = true

# PETSc summary – useful for performance information.

log_summary = true

# Uncomment to launch gdb when starting PyLith.

# start_in_debugger = true

# Friction sensitivity solve used to compute the increment in slip

# associated with changes in the Lagrange multiplier imposed by the

# fault constitutive model.

friction_pc_type = asm

friction_sub_pc_factor_shift_type = nonzero

friction_ksp_max_it = 50

friction_ksp_gmres_restart = 30

# Uncomment to view details of friction sensitivity solve.

friction_ksp_monitor = true

friction_ksp_view = true

friction_ksp_converged_reason = true

##
Error:

PC Object: (friction_) 1 MPI processes

type: asm

total subdomain blocks = 1, amount of overlap = 1

restriction/interpolation type - RESTRICT

Local solve is same for all blocks, in the following KSP and PC objects:

KSP Object: (friction_sub_) 1 MPI processes

type: preonly

maximum iterations=10000, initial guess is zero

tolerances: relative=1e-05, absolute=1e-50, divergence=10000.

left preconditioning

using NONE norm type for convergence test

PC Object: (friction_sub_) 1 MPI processes

type: ilu

out-of-place factorization

0 levels of fill

tolerance for zero pivot 2.22045e-14

using diagonal shift to prevent zero pivot [NONZERO]

matrix ordering: natural

factor fill ratio given 1., needed 1.

Factored matrix follows:

Mat Object: 1 MPI processes

type: seqaij

rows=1198, cols=1198

package used to perform factorization: petsc

total: nonzeros=7180, allocated nonzeros=7180

total number of mallocs used during MatSetValues calls =0

using I-node routines: found 599 nodes, limit used is 5

linear system matrix = precond matrix:

Mat Object: 1 MPI processes

type: seqaij

rows=1198, cols=1198

total: nonzeros=7180, allocated nonzeros=7180

total number of mallocs used during MatSetValues calls =0

using I-node routines: found 599 nodes, limit used is 5

linear system matrix = precond matrix:

Mat Object: 1 MPI processes

type: seqaij

rows=1198, cols=1198, bs=2

total: nonzeros=7180, allocated nonzeros=7180

total number of mallocs used during MatSetValues calls =0

using I-node routines: found 599 nodes, limit used is 5

3000 SNES Function norm 1.831416565858e-07

Nonlinear solve did not converge due to DIVERGED_MAX_IT iterations 3000

[0]PETSC ERROR: --------------------- Error Message --------------------------------------------------------------

[0]PETSC ERROR:

[0]PETSC ERROR: SNESSolve has not converged

[0]PETSC ERROR: See http://www.mcs.anl.gov/petsc/documentation/faq.html for trouble shooting.

[0]PETSC ERROR: Petsc Release Version 3.10.2, Jul, 01, 2019

[0]PETSC ERROR: /home/bolee/PYLITH_FROM_SOURCE/pylith/bin/mpinemesis on a named node127 by bolee Wed Nov 23 12:49:57 2022

[0]PETSC ERROR: Configure options --prefix=/home/bolee/PYLITH_FROM_SOURCE/pylith --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-chaco=1 --download-ml=1 --download-f2cblaslapack=1 --with-hdf5=1 --with-hdf5-dir=/home/bolee/PYLITH_FROM_SOURCE/pylith --with-zlib=1 --LIBS=-lz --with-debugging=0 --with-fc=0 CPPFLAGS=“-I/home/bolee/PYLITH_FROM_SOURCE/pylith/include -I/home/bolee/PYLITH_FROM_SOURCE/pylith/include " LDFLAGS=”-L/home/bolee/PYLITH_FROM_SOURCE/pylith/lib -L/home/bolee/PYLITH_FROM_SOURCE/pylith/lib64 -L/home/bolee/PYLITH_FROM_SOURCE/pylith/lib -L/home/bolee/PYLITH_FROM_SOURCE/pylith/lib64 " CFLAGS=“-g -O2” CXXFLAGS=“-g -O2 -DMPICH_IGNORE_CXX_SEEK” FCFLAGS= PETSC_DIR=/home/bolee/PYLITH_FROM_SOURCE/build/petsc-pylith PETSC_ARCH=arch-pylith

[0]PETSC ERROR: #1 SNESSolve() line 4408 in /home/bolee/PYLITH_FROM_SOURCE/build/petsc-pylith/src/snes/interface/snes.c

[0]PETSC ERROR: #2 void pylith::problems::SolverNonlinear::solve(pylith::topology::Field*, pylith::topology::Jacobian*, const pylith::topology::Field&)() line 152 in …/…/…/pylith-2.2.2/libsrc/pylith/problems/SolverNonlinear.cc

Fatal error. Calling MPI_Abort() to abort PyLith application.

Traceback (most recent call last):

File “/home/bolee/PYLITH_FROM_SOURCE/pylith/lib/python2.7/site-packages/pylith/apps/PetscApplication.py”, line 74, in onComputeNodes

self.main(*args, **kwds)

File “/home/bolee/PYLITH_FROM_SOURCE/pylith/lib/python2.7/site-packages/pylith/apps/PyLithApp.py”, line 138, in main

self.problem.run(self)

File “/home/bolee/PYLITH_FROM_SOURCE/pylith/lib/python2.7/site-packages/pylith/problems/TimeDependent.py”, line 203, in run

self.formulation.step(t, dt)

File “/home/bolee/PYLITH_FROM_SOURCE/pylith/lib/python2.7/site-packages/pylith/problems/Implicit.py”, line 212, in step

self.solver.solve(dispIncr, self.jacobian, residual)

File “/home/bolee/PYLITH_FROM_SOURCE/pylith/lib/python2.7/site-packages/pylith/problems/problems.py”, line 186, in solve

def solve(self, *args): return _problems.SolverNonlinear_solve(self, *args)

RuntimeError: Error detected while in PETSc function.

MPI_ABORT was invoked on rank 0 in communicator MPI_COMM_WORLD

with errorcode -1.

##
NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes.

You may or may not see output from other processes, depending on

exactly when Open MPI kills them.

[warn] Epoll MOD(1) on fd 28 failed. Old events were 6; read change was 0 (none); write change was 2 (del): Bad file descriptor

[warn] Epoll MOD(4) on fd 28 failed. Old events were 6; read change was 2 (del); write change was 0 (none): Bad file descriptor

/home/bolee/PYLITH_FROM_SOURCE/pylith/bin/nemesis: mpirun: exit 255

/home/bolee/PYLITH_FROM_SOURCE/pylith/bin/pylith: /home/bolee/PYLITH_FROM_SOURCE/pylith/bin/nemesis: exit 1