[cig-commits] r20230 - short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction
brad at geodynamics.org
brad at geodynamics.org
Wed May 30 12:45:46 PDT 2012
Author: brad
Date: 2012-05-30 12:45:45 -0700 (Wed, 30 May 2012)
New Revision: 20230
Added:
short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/afterslip_tractions.py
short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/afterslip_tractions.spatialdb
short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/step04.cfg
Modified:
short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/README
short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/step01.cfg
Log:
Added afterslip example.
Modified: short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/README
===================================================================
--- short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/README 2012-05-30 19:14:57 UTC (rev 20229)
+++ short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/README 2012-05-30 19:45:45 UTC (rev 20230)
@@ -86,6 +86,18 @@
pylith step03.cfg
+Step 4. Friction controlled afterslip
+
+ This simulation uses the stress changes associated with the the
+ coseismic deformation in Step 1 in a simulation of afterslip
+ governed by static friction. The afterslip occurs at the down-dip
+ extent of rupture where the coseismic slip increases the shear
+ tractions.
+
+ Run the simulation via the following command:
+ pylith step04.cfg
+
+
Suggestions variations
The list below includes some suggested modifications to the problem
Added: short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/afterslip_tractions.py
===================================================================
--- short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/afterslip_tractions.py (rev 0)
+++ short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/afterslip_tractions.py 2012-05-30 19:45:45 UTC (rev 20230)
@@ -0,0 +1,50 @@
+import numpy
+import h5py
+
+# Load in change in tractions from coseismic simulation
+filename = "output/step01-fault.h5"
+
+h5 = h5py.File(filename, 'r', driver="sec2")
+vertices = h5['geometry/vertices'][:]
+tractions_change = h5['vertex_fields/traction_change'][0,:,:]
+h5.close()
+
+# Parameters for tractions associated with background stress field
+# Nominal density
+density = 2900.0
+gacc = 9.80665
+coef_friction = 0.6
+
+# Background normal tractions are compressive (less than zero because
+# y is negative)
+tractions_bg_normal = density*gacc*(vertices[:,1])
+
+# Background shear tractions are reverse ("right-lateral" is negative)
+# because the normal tractions are negative.
+tractions_bg_shear = coef_friction*tractions_bg_normal
+
+tractions_shear = tractions_bg_shear + tractions_change[:,0]
+tractions_normal = tractions_bg_normal + tractions_change[:,1]
+
+filename = "afterslip_tractions.spatialdb"
+
+# Create coordinate system for output
+from spatialdata.geocoords.CSCart import CSCart
+cs = CSCart()
+cs._configure()
+cs.setSpaceDim(2)
+
+# Create writer for spatialdata base file
+from spatialdata.spatialdb.SimpleIOAscii import SimpleIOAscii
+writer = SimpleIOAscii()
+writer.inventory.filename = filename
+writer._configure()
+writer.write({'points': vertices,
+ 'coordsys': cs,
+ 'data_dim': 1,
+ 'values': [{'name': "traction-shear",
+ 'units': "Pa",
+ 'data': tractions_shear},
+ {'name': "traction-normal",
+ 'units': "Pa",
+ 'data': tractions_normal}]})
Property changes on: short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/afterslip_tractions.py
___________________________________________________________________
Name: svn:executable
+ *
Added: short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/afterslip_tractions.spatialdb
===================================================================
--- short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/afterslip_tractions.spatialdb (rev 0)
+++ short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/afterslip_tractions.spatialdb 2012-05-30 19:45:45 UTC (rev 20230)
@@ -0,0 +1,61 @@
+#SPATIAL.ascii 1
+SimpleDB {
+ num-values = 2
+ value-names = traction-shear traction-normal
+ value-units = Pa Pa
+ num-locs = 48
+ data-dim = 1
+ space-dim = 2
+ cs-data = cartesian {
+ to-meters = 1
+ space-dim = 2
+}
+}
+ -7.900000e+03 -4.000000e+04 -6.806167e+08 -1.138178e+09
+ 6.217879e+04 -2.151659e+04 -3.827731e+08 -6.118344e+08
+ 7.537575e+03 -3.538001e+04 -6.190239e+08 -1.007624e+09
+ 2.301448e+04 -3.089404e+04 -5.233966e+08 -8.808484e+08
+ 3.078374e+04 -2.876014e+04 -4.775233e+08 -8.191033e+08
+ 3.858097e+04 -2.673089e+04 -4.527056e+08 -7.602052e+08
+ 7.809440e+04 -1.900824e+04 -3.214286e+08 -5.416606e+08
+ 7.011929e+04 -2.015270e+04 -3.414270e+08 -5.751686e+08
+ 8.609358e+04 -1.804554e+04 -3.066258e+08 -5.143901e+08
+ 1.500204e+05 -1.009572e+04 -1.617417e+08 -2.825305e+08
+ 1.578557e+05 -8.221434e+03 -1.649012e+08 -2.353765e+08
+ 1.021271e+05 -1.643855e+04 -2.770489e+08 -4.681209e+08
+ 1.101485e+05 -1.568141e+04 -2.654716e+08 -4.457931e+08
+ 1.181656e+05 -1.488009e+04 -2.429373e+08 -4.230659e+08
+ 1.261718e+05 -1.397798e+04 -2.318145e+08 -3.965166e+08
+ 1.341587e+05 -1.291802e+04 -2.164727e+08 -3.665086e+08
+ 1.421138e+05 -1.164290e+04 -1.745752e+08 -3.305236e+08
+ 1.656000e+05 -6.000000e+03 -1.210554e+08 -1.759786e+08
+ -1.827295e+02 -3.768493e+04 -6.498480e+08 -1.072043e+09
+ 1.526778e+04 -3.310856e+04 -5.656255e+08 -9.436834e+08
+ 4.641061e+04 -2.483068e+04 -4.189697e+08 -7.046004e+08
+ 5.427594e+04 -2.308426e+04 -3.959394e+08 -6.559993e+08
+ 9.410690e+04 -1.720783e+04 -2.902652e+08 -4.900463e+08
+ -1.253533e+05 -7.672410e+04 -1.309863e+09 -2.181960e+09
+ -2.391216e+05 -1.279221e+05 -2.182994e+09 -3.638046e+09
+ -6.921401e+04 -5.834591e+04 -9.976098e+08 -1.660803e+09
+ -5.071555e+05 -2.902538e+05 -4.952758e+09 -8.254536e+09
+ -6.000000e+05 -3.400000e+05 -5.801607e+09 -9.669356e+09
+ -3.707411e+05 -2.057980e+05 -3.511744e+09 -5.852682e+09
+ -3.197051e+05 -1.730797e+05 -2.953471e+09 -4.922257e+09
+ -2.760049e+05 -1.479104e+05 -2.524034e+09 -4.206485e+09
+ -2.080485e+05 -1.120422e+05 -1.912077e+09 -3.186437e+09
+ -1.410140e+05 -8.254640e+04 -1.409049e+09 -2.347527e+09
+ -9.212276e+04 -6.549934e+04 -1.118838e+09 -1.862407e+09
+ -1.012707e+05 -6.846624e+04 -1.169246e+09 -1.947230e+09
+ -5.390095e+04 -5.370808e+04 -9.240258e+08 -1.527955e+09
+ -6.156041e+04 -5.601732e+04 -9.693735e+08 -1.593212e+09
+ -3.090161e+04 -4.685024e+04 -8.042750e+08 -1.332769e+09
+ -2.323282e+04 -4.457215e+04 -7.585454e+08 -1.267450e+09
+ -1.818044e+05 -9.963948e+04 -1.700513e+09 -2.833707e+09
+ -8.449740e+04 -6.308054e+04 -1.078065e+09 -1.794688e+09
+ -7.686028e+04 -6.069846e+04 -1.037555e+09 -1.725958e+09
+ -4.623705e+04 -5.141357e+04 -8.685688e+08 -1.461076e+09
+ -3.857014e+04 -4.912917e+04 -8.259365e+08 -1.397715e+09
+ -1.556516e+04 -4.229027e+04 -7.172898e+08 -1.202252e+09
+ -4.316281e+05 -2.458566e+05 -4.195244e+09 -6.991855e+09
+ -1.596584e+05 -9.001740e+04 -1.536412e+09 -2.560041e+09
+ -1.122342e+05 -7.213547e+04 -1.231812e+09 -2.051539e+09
Modified: short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/step01.cfg
===================================================================
--- short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/step01.cfg 2012-05-30 19:14:57 UTC (rev 20229)
+++ short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/step01.cfg 2012-05-30 19:45:45 UTC (rev 20230)
@@ -35,13 +35,12 @@
# boundary conditions
# ----------------------------------------------------------------------
[pylithapp.timedependent]
-# Set bc to an array of 5 boundary conditions:
+# Set bc to an array of 4 boundary conditions:
# 'boundary_east_crust'
# 'boundary_east_mantle'
# 'boundary_west'
-# 'boundary_bottom_slab'
# 'boundary_bottom_mantle'
-bc = [boundary_east_crust,boundary_east_mantle,boundary_west,boundary_bottom_slab,boundary_bottom_mantle]
+bc = [boundary_east_crust,boundary_east_mantle,boundary_west,boundary_bottom_mantle]
# For all boundaries, we fix the displacement normal to the boundary
# (roller boundary condition) by retaining the default ZeroDispDB,
@@ -49,12 +48,13 @@
#
# The label corresponds to the name of the nodeset in CUBIT.
-# East boundaries
+# East boundary (crust)
[pylithapp.timedependent.bc.boundary_east_crust]
bc_dof = [0]
label = bndry_east_crust
db_initial.label = Dirichlet BC on east boundary (crust)
+# East boundary (mantle)
[pylithapp.timedependent.bc.boundary_east_mantle]
bc_dof = [0]
label = bndry_east_mantle
@@ -66,12 +66,7 @@
label = bndry_west
db_initial.label = Dirichlet BC on west boundary
-# Bottom boundaries
-[pylithapp.timedependent.bc.boundary_bottom_slab]
-bc_dof = [1]
-label = bndry_bot_slab
-db_initial.label = Dirichlet BC on bottom boundary (slab)
-
+# Bottom boundary (mantle)
[pylithapp.timedependent.bc.boundary_bottom_mantle]
bc_dof = [1]
label = bndry_bot_mantle
@@ -90,7 +85,7 @@
# Set the parameters for the fault interface condition.
[pylithapp.timedependent.interfaces.fault]
# The label corresponds to the name of the nodeset in CUBIT.
-label = fault_coseismic
+label = fault_slabtop
# We must define the quadrature information for fault cells.
# The fault cells are 1D (line).
Added: short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/step04.cfg
===================================================================
--- short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/step04.cfg (rev 0)
+++ short/3D/PyLith/branches/v1.7-trunk/examples/2d/subduction/step04.cfg 2012-05-30 19:45:45 UTC (rev 20230)
@@ -0,0 +1,176 @@
+[pylithapp]
+
+# ----------------------------------------------------------------------
+# PROBLEM DESCRIPTION
+# ----------------------------------------------------------------------
+#
+# This simulation computes afterslip following the coseismic rupture
+# (imposed in step01) with slip occurring where the fault tractions
+# increase. We assume initial tractions on the fault surface that are
+# 2.0 MPa below the friction threshold before adding in changes in
+# fault tractions from coseismic slip.
+#
+# ----------------------------------------------------------------------
+# RUNNING THE SIMULATON
+# ----------------------------------------------------------------------
+#
+# This is not a self-contained simulation configuration file. This
+# file specifies only the boundary conditions and earthquake
+# parameters for the simulation. The general quasi-static and mesh
+# parameters are specificed in the pylithapp.cfg file which PyLith
+# reads by default.
+#
+# To run the simulation:
+# pylith step04.cfg
+#
+# Output will be directed to directory output.
+
+# ----------------------------------------------------------------------
+# problem
+# ----------------------------------------------------------------------
+[pylithapp.timedependent.formulation.time_step]
+total_time = 0.0*year
+dt = 5.0*year
+
+[pylithapp.timedependent]
+# For this problem that uses friction we must switch to a nonlinear solver.
+implicit.solver = pylith.problems.SolverNonlinear
+
+# ----------------------------------------------------------------------
+# boundary conditions
+# ----------------------------------------------------------------------
+[pylithapp.timedependent]
+# Set bc to an array of 4 boundary conditions:
+# 'boundary_east_crust'
+# 'boundary_east_mantle'
+# 'boundary_west'
+# 'boundary_bottom_mantle'
+bc = [boundary_east_crust,boundary_east_mantle,boundary_west,boundary_bottom_mantle]
+
+# For all boundaries, we fix the displacement normal to the boundary
+# (roller boundary condition) by retaining the default ZeroDispDB,
+# which specifies a zero value.
+#
+# The label corresponds to the name of the nodeset in CUBIT.
+
+# East boundary (crust)
+[pylithapp.timedependent.bc.boundary_east_crust]
+bc_dof = [0]
+label = bndry_east_crust
+db_initial.label = Dirichlet BC on east boundary (crust)
+
+# East boundary (mantle)
+[pylithapp.timedependent.bc.boundary_east_mantle]
+bc_dof = [0]
+label = bndry_east_mantle
+db_initial.label = Dirichlet BC on east boundary (mantle)
+
+# West boundary
+[pylithapp.timedependent.bc.boundary_west]
+bc_dof = [0]
+label = bndry_west
+db_initial.label = Dirichlet BC on west boundary
+
+# Bottom boundary (mantle)
+[pylithapp.timedependent.bc.boundary_bottom_mantle]
+bc_dof = [1]
+label = bndry_bot_mantle
+db_initial.label = Dirichlet BC on bottom boundary (mantle)
+
+
+# ----------------------------------------------------------------------
+# faults
+# ----------------------------------------------------------------------
+[pylithapp.timedependent]
+interfaces = [fault]
+
+# Set the type of fault interface condition.
+[pylithapp.timedependent.interfaces]
+fault = pylith.faults.FaultCohesiveDyn
+
+# Set the parameters for the fault interface condition.
+[pylithapp.timedependent.interfaces.fault]
+# The label corresponds to the name of the nodeset in CUBIT.
+label = fault_slabtop
+
+# We must define the quadrature information for fault cells.
+# The fault cells are 1D (line).
+quadrature.cell = pylith.feassemble.FIATSimplex
+quadrature.cell.dimension = 1
+
+# Specify zero tolerance for detecting slip.
+zero_tolerance = 1.0e-9
+
+# Friction model
+friction = pylith.friction.StaticFriction
+friction.label = Static friction
+
+friction.db_properties = spatialdata.spatialdb.UniformDB
+friction.db_properties.label = Static friction
+friction.db_properties.values = [friction-coefficient, cohesion]
+friction.db_properties.data = [0.6, 2.0*MPa]
+
+# Initial tractions
+traction_perturbation = pylith.faults.TractPerturbation
+
+[pylithapp.timedependent.interfaces.fault.traction_perturbation]
+db_initial = spatialdata.spatialdb.SimpleDB
+db_initial.label = Initial fault tractions
+db_initial.iohandler.filename = afterslip_tractions.spatialdb
+db_initial.query_type = nearest
+
+# ----------------------------------------------------------------------
+# output
+# ----------------------------------------------------------------------
+# Domain
+[pylithapp.problem.formulation.output.domain]
+writer.filename = output/step04.h5
+
+# Ground surface
+[pylithapp.problem.formulation.output.subdomain]
+writer.filename = output/step04-groundsurf.h5
+
+# Fault
+[pylithapp.problem.interfaces.fault.output]
+writer = pylith.meshio.DataWriterHDF5SubSubMesh
+writer.filename = output/step04-fault.h5
+
+# Materials
+[pylithapp.timedependent.materials.continent_crust.output]
+writer.filename = output/step04-concrust.h5
+
+[pylithapp.timedependent.materials.continent_mantle.output]
+writer.filename = output/step04-conmantle.h5
+
+[pylithapp.timedependent.materials.ocean_crust.output]
+writer.filename = output/step04-oceancrust.h5
+
+[pylithapp.timedependent.materials.ocean_mantle.output]
+writer.filename = output/step04-oceanmantle.h5
+
+
+# ----------------------------------------------------------------------
+# PETSc
+# ----------------------------------------------------------------------
+# NOTE: There are additional settings specific to fault friction.
+[pylithapp.petsc]
+
+# 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 = 25
+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
+
+# Convergence parameters. Tighten tolerances to be less than zero
+# tolerance for slip. Increase number of iterations.
+ksp_rtol = 1.0e-12
+ksp_atol = 1.0e-13
+ksp_max_it = 400
+
+# End of file
More information about the CIG-COMMITS
mailing list