[CIG-SHORT] Fault direction in 3D fault
Jiangzhi Chen
jzchenjzarthur at gmail.com
Thu May 1 15:26:36 PDT 2014
Hi everyone,
I have a 3D subduction model simulating the interseismic
deformation, where I have both fault_slabtop and fault_slabbot creeping.
I have set fault_slabtop with reverse_slip of 8 cm/yr and fault_slabbot
with reverse_slip of -8 cm/yr, and used the up_dir parameter [0,0,1],
but still the subduction is in wrong direction. I am not sure what is
wrong. The configuration files are attached. Thanks.
Jiangzhi
-------------- next part --------------
A non-text attachment was scrubbed...
Name: screen.png
Type: image/png
Size: 42397 bytes
Desc: not available
URL: <http://lists.geodynamics.org/pipermail/cig-short/attachments/20140501/d1ff5eca/attachment-0001.png>
-------------- next part --------------
[pylithapp]
# This is not a self-contained simulation configuration file. This
# file only specifies the general parameters common to the simulations
# in this directory.
# ----------------------------------------------------------------------
# journal
# ----------------------------------------------------------------------
# Turn on some journals to show progress.
[pylithapp.journal.info]
timedependent = 1
implicit = 1
petsc = 1
solverlinear = 1
meshiocubit = 1
implicitelasticity = 1
faultcohesivekin = 1
fiatsimplex = 1
pylithapp = 1
materials = 1
# ----------------------------------------------------------------------
# mesh_generator
# ----------------------------------------------------------------------
[pylithapp.mesh_generator]
# Change the default mesh reader to the CUBIT reader.
reader = pylith.meshio.MeshIOCubit
[pylithapp.mesh_generator.reader]
filename = mesh3d.exo
coordsys.space_dim = 3
# ----------------------------------------------------------------------
# problem
# ----------------------------------------------------------------------
[pylithapp.timedependent]
dimension = 3
[pylithapp.timedependent.formulation.time_step]
# Define the total time for the simulation and the time step size.
total_time = 200.0*year
dt = 10.0*year
# [pylithapp.timedependent.formulation]
# split_fields = True
# matrix_type = aij
# use_custom_constraint_pc = True
# ----------------------------------------------------------------------
# materials
# ----------------------------------------------------------------------
[pylithapp.timedependent]
# Set materials to an array of 4 materials:
# 'continent_litho'
# 'continent_mantle'
# 'ocean_litho'
# 'ocean_mantle'
materials = [continent_litho,continent_mantle,ocean_litho,ocean_mantle]
[pylithapp.timedependent.materials]
# Set bulk constitutive model for each material.
continent_litho = pylith.materials.ElasticIsotropic3D
ocean_litho = pylith.materials.ElasticIsotropic3D
continent_mantle = pylith.materials.MaxwellIsotropic3D
ocean_mantle = pylith.materials.MaxwellIsotropic3D
# Oceanic mantle
[pylithapp.timedependent.materials.ocean_mantle]
label = Oceanic mantle
# The id corresponds to the block number from CUBIT.
id = 1
db_properties.label = Oceanic mantle properties
db_properties.iohandler.filename = mat_oceanmantle.spatialdb
# We are doing 3D quadrature for a tetrahedron.
quadrature.cell = pylith.feassemble.FIATSimplex
quadrature.cell.dimension = 3
# Continental litho
[pylithapp.timedependent.materials.continent_litho]
label = Continental litho
# The id corresponds to the block number from CUBIT.
id = 2
db_properties.label = Continental litho properties
db_properties.iohandler.filename = mat_conlitho.spatialdb
# We are doing 3D quadrature for a tetrahedron.
quadrature.cell = pylith.feassemble.FIATSimplex
quadrature.cell.dimension = 3
# Continental mantle
[pylithapp.timedependent.materials.continent_mantle]
label = Continental mantle
# The id corresponds to the block number from CUBIT.
id = 3
db_properties.label = Continental mantle properties
db_properties.iohandler.filename = mat_conmantle.spatialdb
# We are doing 3D quadrature for a tetrahedron.
quadrature.cell = pylith.feassemble.FIATSimplex
quadrature.cell.dimension = 3
# Oceanic litho
[pylithapp.timedependent.materials.ocean_litho]
label = Oceanic litho
# The id corresponds to the block number from CUBIT.
id = 4
db_properties.label = Oceanic litho properties
db_properties.iohandler.filename = mat_oceanlitho.spatialdb
# We are doing 3D quadrature for a tetrahedron.
quadrature.cell = pylith.feassemble.FIATSimplex
quadrature.cell.dimension = 3
# ----------------------------------------------------------------------
# output
# ----------------------------------------------------------------------
# Names of output files are set in stepXX.cfg. We consolidate all of the
# output settings that are common to all of the simulations here.
[pylithapp.timedependent.implicit]
# Set the output to an array of 2 output managers.
# We will output the solution over the domain and the ground surface.
output = [domain,subdomain]
# Set subdomain component to OutputSolnSubset (subset of domain).
output.subdomain = pylith.meshio.OutputSolnSubset
# Domain
[pylithapp.problem.formulation.output.domain]
output_freq = time_step
time_step = 9.99999*year
writer = pylith.meshio.DataWriterHDF5Mesh
# Ground surface
[pylithapp.problem.formulation.output.subdomain]
label = groundsurf ; # Name of CUBIT nodeset for ground surface.
writer = pylith.meshio.DataWriterHDF5SubMesh
# Materials
[pylithapp.timedependent.materials.continent_litho.output]
cell_filter = pylith.meshio.CellFilterAvgMesh
output_freq = time_step
time_step = 9.99999*year
writer = pylith.meshio.DataWriterHDF5Mesh
[pylithapp.timedependent.materials.continent_mantle.output]
cell_filter = pylith.meshio.CellFilterAvgMesh
output_freq = time_step
time_step = 9.99999*year
writer = pylith.meshio.DataWriterHDF5Mesh
[pylithapp.timedependent.materials.ocean_litho.output]
cell_filter = pylith.meshio.CellFilterAvgMesh
output_freq = time_step
time_step = 9.99999*year
writer = pylith.meshio.DataWriterHDF5Mesh
[pylithapp.timedependent.materials.ocean_mantle.output]
cell_filter = pylith.meshio.CellFilterAvgMesh
output_freq = time_step
time_step = 9.99999*year
writer = pylith.meshio.DataWriterHDF5Mesh
# ----------------------------------------------------------------------
# PETSc
# ----------------------------------------------------------------------
# Set the solver options.
[pylithapp.timedependent.formulation]
split_fields = True
matrix_type = aij
use_custom_constraint_pc = True
[pylithapp.petsc]
fs_pc_type = fieldsplit
fs_pc_fieldsplit_type = schur
fs_pc_fieldsplit_schur_factorization_type = full
fs_pc_fieldsplit_real_diagonal = true
fs_fieldsplit_0_pc_type = gamg
fs_fieldsplit_0_ksp_type = preonly
fs_fieldsplit_0_ksp_rtol = 1.0e-10
fs_fieldsplit_1_pc_type = none
fs_fieldsplit_1_ksp_type = minres
fs_fieldsplit_1_ksp_rtol = 1.0e-10
# # Preconditioner settings.
# pc_type = asm
# sub_pc_factor_shift_type = nonzero
#
# Convergence parameters.
ksp_rtol = 1.0e-8
ksp_atol = 1.0e-12
ksp_max_it = 1000
ksp_gmres_restart = 50
# Linear solver monitoring options.
ksp_monitor = true
ksp_view = true
ksp_converged_reason = true
# Nonlinear solver monitoring options.
snes_rtol = 1.0e-8
snes_atol = 1.0e-12
snes_max_it = 1000
snes_monitor = true
snes_view = true
snes_converged_reason = true
# PETSc summary -- useful for performance information.
#log_summary = true
# Uncomment to launch gdb when starting PyLith.
# start_in_debugger = true
# End of file
-------------- next part --------------
[pylithapp]
# ----------------------------------------------------------------------
# PROBLEM DESCRIPTION
# ----------------------------------------------------------------------
#
# This simulation involves aseismic creep along the interfaces between
# the subducting oceanic litho and the mantle. The slip rate is a
# constant 8 cm/yr. The slip time is 40 years for quick test.
#
# ----------------------------------------------------------------------
# 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.
#
# ----------------------------------------------------------------------
# problem
# ----------------------------------------------------------------------
[pylithapp.timedependent.formulation.time_step]
total_time = 40.0*year
dt = 10.0*year
# ----------------------------------------------------------------------
# boundary conditions
# ----------------------------------------------------------------------
[pylithapp.timedependent]
# Set bc to an array of 3 boundary conditions:
# 'boundary_west_mantle'
# 'boundary_east'
# 'boundary_bottom_mantle'
# 'boundary_side'
bc = [boundary_west_mantle,boundary_east,boundary_bottom_mantle,boundary_side]
# 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.
# west boundary (mantle)
[pylithapp.timedependent.bc.boundary_west_mantle]
bc_dof = [0]
label = bndry_west_mantle
db_initial.label = Dirichlet BC on west boundary (mantle)
# east boundary
[pylithapp.timedependent.bc.boundary_east]
bc_dof = [0]
label = bndry_east
db_initial.label = Dirichlet BC on east 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)
# front and back sides
[pylithapp.timedependent.bc.boundary_side]
bc_dof = [2]
label = bndry_side
db_initial.label = Dirichlet BC on side boundary (all)
# ----------------------------------------------------------------------
# faults
# ----------------------------------------------------------------------
[pylithapp.timedependent]
interfaces = [fault_slabtop,fault_slabbot]
# Set the type of fault interface condition.
[pylithapp.timedependent.interfaces]
fault_slabtop = pylith.faults.FaultCohesiveKin
fault_slabbot = pylith.faults.FaultCohesiveKin
# Slab top --------------------
[pylithapp.timedependent.interfaces.fault_slabtop]
# The label corresponds to the name of the nodeset in CUBIT.
label = fault_slabtop
id = 100
up_dir = [0, 0, 1]
# We must define the quadrature information for fault cells.
# The fault cells are 2D (surface).
quadrature.cell = pylith.feassemble.FIATSimplex
quadrature.cell.dimension = 2
# Switch to constant slip rate time function.
[pylithapp.timedependent.interfaces.fault_slabtop.eq_srcs.rupture]
slip_function = pylith.faults.ConstRateSlipFn
# The slip time and final slip are defined in spatial databases.
[pylithapp.timedependent.interfaces.fault_slabtop.eq_srcs.rupture.slip_function]
# slip_rate.iohandler.filename = fault_creep_slabtop.spatialdb
# slip_rate.query_type = linear
# slip_rate.label = Final slip
slip_rate = spatialdata.spatialdb.UniformDB
slip_rate.label = Slip rate
slip_rate.values = [left-lateral-slip,reverse-slip,fault-opening]
slip_rate.data = [0.0*cm/year, 8.0*cm/year, 0.0*cm/year]
# Slip time is uniform, so use UniformDB for convenience
slip_time = spatialdata.spatialdb.UniformDB
slip_time.label = Slip time
slip_time.values = [slip-time]
slip_time.data = [0.0*year]
# Slab bottom --------------------
[pylithapp.timedependent.interfaces.fault_slabbot]
# The label corresponds to the name of the nodeset in CUBIT.
label = fault_slabbot
id = 101
# We must define the quadrature information for fault cells.
# The fault cells are 2D (surface).
quadrature.cell = pylith.feassemble.FIATSimplex
quadrature.cell.dimension = 2
# Switch to constant slip rate time function.
[pylithapp.timedependent.interfaces.fault_slabbot.eq_srcs.rupture]
slip_function = pylith.faults.ConstRateSlipFn
# The creep rate and slip time are uniform, so use UniformDB for simplicity.
[pylithapp.timedependent.interfaces.fault_slabbot.eq_srcs.rupture.slip_function]
slip_rate = spatialdata.spatialdb.UniformDB
slip_rate.label = Slip rate
slip_rate.values = [left-lateral-slip,reverse-slip,fault-opening]
slip_rate.data = [0.0*cm/year, -8.0*cm/year, 0.0*cm/year]
# Slip time is uniform, so use UniformDB for convenience
slip_time = spatialdata.spatialdb.UniformDB
slip_time.label = Slip time
slip_time.values = [slip-time]
slip_time.data = [0.0*year]
# ----------------------------------------------------------------------
# output
# ----------------------------------------------------------------------
# Domain
[pylithapp.problem.formulation.output.domain]
writer.filename = interseismic1/interseismic.h5
# Ground surface
[pylithapp.problem.formulation.output.subdomain]
writer.filename = interseismic1/interseismic-groundsurf.h5
# Faults
[pylithapp.problem.interfaces.fault_slabtop.output]
writer = pylith.meshio.DataWriterHDF5SubSubMesh
writer.filename = interseismic1/interseismic-fault-slabtop.h5
[pylithapp.problem.interfaces.fault_slabbot.output]
writer = pylith.meshio.DataWriterHDF5SubSubMesh
writer.filename = interseismic1/interseismic-fault-slabbot.h5
# Materials
[pylithapp.timedependent.materials.continent_litho.output]
writer.filename = interseismic1/interseismic-conlitho.h5
[pylithapp.timedependent.materials.continent_mantle.output]
writer.filename = interseismic1/interseismic-conmantle.h5
[pylithapp.timedependent.materials.ocean_litho.output]
writer.filename = interseismic1/interseismic-oceanlitho.h5
[pylithapp.timedependent.materials.ocean_mantle.output]
writer.filename = interseismic1/interseismic-oceanmantle.h5
# End of file
-------------- next part --------------
// -*- C++ -*- (syntax highlighting)
//
// This spatial database specifies the parameters for the continental
// crust material, which is a linear Maxwell viscoelastic solid.
//
// The viscosity corresponds to a Maxwell time of 200 years.
// Maxwell time = viscosity / shear_modulus
//
#SPATIAL.ascii 1
SimpleDB {
num-values = 4
value-names = density vs vp viscosity
value-units = kg/m**3 m/s m/s Pa*s // units
num-locs = 1 // number of locations
data-dim = 0
space-dim = 3
cs-data = cartesian {
to-meters = 1.0
space-dim = 3
}
}
// Columns are
// (1) x coordinate (m)
// (2) y coordinate (m)
// (3) z coordinate (m)
// (4) density (kg/m^3)
// (5) vs (m/s)
// (6) vp (m/s)
// (7) viscosity (Pa*s)
0.0 0.0 0.0 4000.0 5600.0 10000.0 7e+20
-------------- next part --------------
// -*- C++ -*- (syntax highlighting)
//
// This spatial database specifies the parameters for the oceanic
// crust material, which is linearly elastic.
//
#SPATIAL.ascii 1
SimpleDB {
// number of physical properties
num-values = 3
// Names and units of physical properties
value-names = density vs vp
value-units = kg/m**3 m/s m/s
// 1 location -> data is uniform -> data dimension is 0
num-locs = 1
data-dim = 0
// Problem is in 3-D
space-dim = 3
// Coordinate system (3-D with coordinates in meters)
cs-data = cartesian {
to-meters = 1.0
space-dim = 3
}
}
// Columns are x, y, z, density, vs, vp
0.0 0.0 0.0 3300.0 4500.0 8100.0
-------------- next part --------------
// -*- C++ -*- (syntax highlighting)
//
// This spatial database specifies the parameters for the continental
// crust material, which is a linear Maxwell viscoelastic solid.
//
// The viscosity corresponds to a Maxwell time of 200 years.
// Maxwell time = viscosity / shear_modulus
//
#SPATIAL.ascii 1
SimpleDB {
num-values = 4
value-names = density vs vp viscosity
value-units = kg/m**3 m/s m/s Pa*s // units
num-locs = 1 // number of locations
data-dim = 0
space-dim = 3
cs-data = cartesian {
to-meters = 1.0
space-dim = 3
}
}
// Columns are
// (1) x coordinate (m)
// (2) y coordinate (m)
// (3) z coordinate (m)
// (4) density (kg/m^3)
// (5) vs (m/s)
// (6) vp (m/s)
// (7) viscosity (Pa-s)
0.0 0.0 0 4000.0 5600.0 10000.0 7e+20
-------------- next part --------------
// -*- C++ -*- (syntax highlighting)
//
// This spatial database specifies the parameters for the continental
// lithosphere material, which is linearly elastic.
//
#SPATIAL.ascii 1
SimpleDB {
// number of physical properties
num-values = 3
// Names and units of physical properties
value-names = density vs vp
value-units = kg/m**3 m/s m/s
// data dimension is 3
num-locs = 10
data-dim = 3
// Problem is in 3-D
space-dim = 3
// Coordinate system (3-D with coordinates in meters)
cs-data = cartesian {
to-meters = 1.0e+3
space-dim = 3
}
}
// Columns are x, y, z, density, vs, vp
676.250875 38.896813 -44.036906 2844.578786 4781.189817 7203.462665
651.958375 35.248172 -38.508008 2789.975056 4492.486075 6708.886601
668.759437 59.451211 -49.604582 2903.650795 4904.665122 7490.714456
659.807000 40.145496 -52.226375 2904.611003 4948.966191 7490.144999
545.609500 -441.873844 -29.972080 2775.120259 3814.499402 6302.034496
558.610937 -429.389313 -34.406617 2794.772671 4103.195185 6525.524078
566.326625 -439.211937 -33.212914 2746.492207 4057.181867 6293.626742
560.179188 -446.880750 -38.889770 2808.076413 4163.244930 6641.081019
652.438562 24.011359 -54.045871 2898.740760 4870.841339 7393.448999
672.746688 16.530051 -43.541328 2830.499942 4754.365784 7117.495957
More information about the CIG-SHORT
mailing list