[cig-commits] r16952 - in short/3D/PyLith/trunk/doc/userguide/tutorials/3dhex8: . friction

willic3 at geodynamics.org willic3 at geodynamics.org
Tue Jun 8 23:09:56 PDT 2010 

Author: willic3
Date: 2010-06-08 23:09:56 -0700 (Tue, 08 Jun 2010)
New Revision: 16952

Added:
   short/3D/PyLith/trunk/doc/userguide/tutorials/3dhex8/friction/figs/
   short/3D/PyLith/trunk/doc/userguide/tutorials/3dhex8/friction/friction.lyx
Modified:
   short/3D/PyLith/trunk/doc/userguide/tutorials/3dhex8/3dhex8.lyx
Log:
Started on friction section.



Modified: short/3D/PyLith/trunk/doc/userguide/tutorials/3dhex8/3dhex8.lyx
===================================================================
--- short/3D/PyLith/trunk/doc/userguide/tutorials/3dhex8/3dhex8.lyx	2010-06-09 05:53:21 UTC (rev 16951)
+++ short/3D/PyLith/trunk/doc/userguide/tutorials/3dhex8/3dhex8.lyx	2010-06-09 06:09:56 UTC (rev 16952)
@@ -563,6 +563,13 @@
 \end_inset
 
 
+\begin_inset CommandInset include
+LatexCommand include
+filename "friction/friction.lyx"
+
+\end_inset
+
+
 \end_layout
 
 \begin_layout Subsection

Added: short/3D/PyLith/trunk/doc/userguide/tutorials/3dhex8/friction/friction.lyx
===================================================================
--- short/3D/PyLith/trunk/doc/userguide/tutorials/3dhex8/friction/friction.lyx	                        (rev 0)
+++ short/3D/PyLith/trunk/doc/userguide/tutorials/3dhex8/friction/friction.lyx	2010-06-09 06:09:56 UTC (rev 16952)
@@ -0,0 +1,2187 @@
+#LyX 1.6.6 created this file. For more info see http://www.lyx.org/
+\lyxformat 345
+\begin_document
+\begin_header
+\textclass book
+\begin_preamble
+
+\end_preamble
+\use_default_options false
+\language english
+\inputencoding latin1
+\font_roman default
+\font_sans default
+\font_typewriter default
+\font_default_family default
+\font_sc false
+\font_osf false
+\font_sf_scale 100
+\font_tt_scale 100
+
+\graphics default
+\paperfontsize default
+\spacing single
+\use_hyperref false
+\papersize default
+\use_geometry true
+\use_amsmath 0
+\use_esint 0
+\cite_engine basic
+\use_bibtopic false
+\paperorientation portrait
+\leftmargin 1in
+\topmargin 1in
+\rightmargin 1in
+\bottommargin 2in
+\secnumdepth 3
+\tocdepth 3
+\paragraph_separation indent
+\defskip medskip
+\quotes_language english
+\papercolumns 1
+\papersides 1
+\paperpagestyle default
+\tracking_changes false
+\output_changes false
+\author "" 
+\author "" 
+\end_header
+
+\begin_body
+
+\begin_layout Subsection
+\begin_inset CommandInset label
+LatexCommand label
+name "sec:Tutorial-3d-hex8-friction"
+
+\end_inset
+
+Tutorial Using Hexahedral Mesh Created by CUBIT - Fault Friction Examples
+\end_layout
+
+\begin_layout Standard
+PyLith features discussed in this tutorial:
+\end_layout
+
+\begin_layout Itemize
+Static fault friction
+\end_layout
+
+\begin_layout Itemize
+Slip-weakening fault friction
+\end_layout
+
+\begin_layout Itemize
+Rate-and-state fault friction
+\end_layout
+
+\begin_layout Itemize
+Nonlinear solver
+\end_layout
+
+\begin_layout Subsubsection
+Overview
+\end_layout
+
+\begin_layout Standard
+This set of examples provides an introduction to using fault friction in
+ static and quasi-static problems with PyLith.
+ The boundary conditions are all either static or quasi-static Dirichlet
+ conditions, and only elastic materials are used.
+ Thus, the only aspects that were not covered in previous examples involve
+ the use of fault friction.
+ All of the fault friction models are nonlinear, and they therefore require
+ the usage of the nonlinear solver.
+ All of the examples are contained in the directory 
+\family typewriter
+examples/3d/hex8
+\family default
+, and the corresponding 
+\family typewriter
+.cfg
+\family default
+ files are 
+\family typewriter
+step10.cfg
+\family default
+, 
+\family typewriter
+step11.cfg
+\family default
+, 
+\family typewriter
+step12.cfg
+\family default
+, 
+\family typewriter
+step13.cfg
+\family default
+, and 
+\family typewriter
+step14.cfg
+\family default
+.
+ Each example may be run as follows:
+\end_layout
+
+\begin_layout LyX-Code
+pylith stepxx.cfg
+\end_layout
+
+\begin_layout Standard
+This will cause PyLith to read the default parameters in 
+\family typewriter
+pylithapp.cfg
+\family default
+, and then override or augment them with the additional parameters in the
+ 
+\family typewriter
+stepxx.cfg
+\family default
+ file.
+ Each 
+\family typewriter
+.cfg
+\family default
+ file is extensively documented, to provide detailed information on the
+ various parameters.
+\end_layout
+
+\begin_layout Subsubsection
+Step10 - Static Friction (Stick) with Static Dirichlet Boundary Conditions
+\end_layout
+
+\begin_layout Standard
+The 
+\family typewriter
+step10.cfg
+\family default
+ file describes a problem that is identical to example step01, except for
+ the presence of a vertical fault with static friction.
+ In this case, the applied displacements are insufficient to cause the fault
+ to slip, so the solution should be identical to the results of example
+ step01.
+ As in previous examples involving faults, we must first provide an array
+ defining the fault interfaces:
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent]
+\end_layout
+
+\begin_layout LyX-Code
+# Set interfaces to an array of 1 fault: 'fault'.
+\end_layout
+
+\begin_layout LyX-Code
+interfaces = [fault]
+\end_layout
+
+\begin_layout Standard
+Since all fault friction models are nonlinear we must also use the nonlinear
+ solver:
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.implicit]
+\end_layout
+
+\begin_layout LyX-Code
+# Fault friction is a nonlinear problem so we need to use the nonlinear
+\end_layout
+
+\begin_layout LyX-Code
+# solver.
+\end_layout
+
+\begin_layout LyX-Code
+solver = pylith.problems.SolverNonlinear
+\end_layout
+
+\begin_layout Standard
+We need to change the fault interface from the default (
+\family typewriter
+FaultCohesiveKin
+\family default
+) to 
+\family typewriter
+FaultCohesiveDyn
+\family default
+, and we must also define the friction model to use:
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.interfaces]
+\end_layout
+
+\begin_layout LyX-Code
+# Change fault to dynamic fault interface.
+\end_layout
+
+\begin_layout LyX-Code
+fault = pylith.faults.FaultCohesiveDyn
+\end_layout
+
+\begin_layout LyX-Code
+
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.interfaces.fault]
+\end_layout
+
+\begin_layout LyX-Code
+# Use the static friction model.
+\end_layout
+
+\begin_layout LyX-Code
+friction = pylith.friction.StaticFriction
+\end_layout
+
+\begin_layout LyX-Code
+
+\end_layout
+
+\begin_layout Standard
+We then change the material type of the lower crust, provide a spatial database
+ from which to obtain the material properties (using the default 
+\family typewriter
+SimpleDB
+\family default
+), and request additional output information for the material:
+\end_layout
+
+\begin_layout LyX-Code
+# Change material type of lower crust to Maxwell viscoelastic.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent]
+\end_layout
+
+\begin_layout LyX-Code
+materials.lower_crust = pylith.materials.MaxwellIsotropic3D
+\end_layout
+
+\begin_layout LyX-Code
+
+\end_layout
+
+\begin_layout LyX-Code
+# Provide a spatial database from which to obtain property values.
+\end_layout
+
+\begin_layout LyX-Code
+# Since there are additional properties and state variables for the Maxwell
+\end_layout
+
+\begin_layout LyX-Code
+# model, we explicitly request that they be output.
+ Properties are named in
+\end_layout
+
+\begin_layout LyX-Code
+# cell_info_fields and state variables are named in cell_data_fields.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.materials.lower_crust]
+\end_layout
+
+\begin_layout LyX-Code
+db_properties.iohandler.filename = spatialdb/mat_maxwell.spatialdb
+\end_layout
+
+\begin_layout LyX-Code
+output.cell_info_fields = [density,mu,lambda,maxwell_time]
+\end_layout
+
+\begin_layout LyX-Code
+output.cell_data_fields = [total_strain,stress,viscous_strain]
+\end_layout
+
+\begin_layout Standard
+Note that the default 
+\family typewriter
+output.cell_info_fields
+\family default
+ are those corresponding to an elastic material (
+\family typewriter
+density
+\family default
+, 
+\family typewriter
+mu
+\family default
+, 
+\family typewriter
+lambda
+\family default
+), and the default 
+\family typewriter
+output.cell_data_fields
+\family default
+ are 
+\family typewriter
+total_strain
+\family default
+ and 
+\family typewriter
+stress
+\family default
+.
+ For materials other than elastic, there are generally additional material
+ properties and state variables, and the appropriate additional fields must
+ be specifically requested for each material type.
+\end_layout
+
+\begin_layout Standard
+This example has no displacements in the elastic solution (t = 0), so we
+ retain the default 
+\family typewriter
+ZeroDispDB
+\family default
+ for all instances of 
+\family typewriter
+db_initial
+\family default
+.
+ To apply the velocity boundary conditions, we must specify 
+\family typewriter
+db_rate
+\family default
+, which is zero by default.
+ We use a 
+\family typewriter
+UniformDB
+\family default
+ to assign the velocities:
+\end_layout
+
+\begin_layout LyX-Code
+# +x face
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.bc.x_pos]
+\end_layout
+
+\begin_layout LyX-Code
+bc_dof = [0, 1]
+\end_layout
+
+\begin_layout LyX-Code
+label = face_xpos
+\end_layout
+
+\begin_layout LyX-Code
+db_initial.label = Dirichlet BC on +x
+\end_layout
+
+\begin_layout LyX-Code
+db_rate = spatialdata.spatialdb.UniformDB
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.label = Dirichlet rate BC on +x
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.values = [displacement-rate-x,displacement-rate-y,rate-start-time]
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.data = [0.0*cm/year,1.0*cm/year,0.0*year]
+\end_layout
+
+\begin_layout LyX-Code
+
+\end_layout
+
+\begin_layout LyX-Code
+# -x face
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.bc.x_neg]
+\end_layout
+
+\begin_layout LyX-Code
+bc_dof = [0, 1]
+\end_layout
+
+\begin_layout LyX-Code
+label = face_xneg
+\end_layout
+
+\begin_layout LyX-Code
+db_initial.label = Dirichlet BC on -x
+\end_layout
+
+\begin_layout LyX-Code
+db_rate = spatialdata.spatialdb.UniformDB
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.label = Dirichlet rate BC on +x
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.values = [displacement-rate-x,displacement-rate-y,rate-start-time]
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.data = [0.0*cm/year,-1.0*cm/year,0.0*year]
+\end_layout
+
+\begin_layout Standard
+Note that 
+\family typewriter
+db_rate
+\family default
+ requires a start time, which allows the condition to be applied at any
+ time during the simulation.
+ For this example, we start the velocity boundary conditions at t = 0.
+\end_layout
+
+\begin_layout Standard
+Finally, we must provide information on VTK output.
+ This is slightly more complicated than the static case, because we must
+ decide the frequency with which output occurs for each output type.
+ We also assign a more user-friendly format to the output file time stamp,
+ and we request that the time stamp is in units of 1 year (rather than the
+ default value of seconds):
+\end_layout
+
+\begin_layout LyX-Code
+# Give basename for VTK domain output of solution over domain.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.problem.formulation.output.domain]
+\end_layout
+
+\begin_layout LyX-Code
+# We specify that output occurs in terms of a given time frequency, and
+\end_layout
+
+\begin_layout LyX-Code
+# ask for output every 40 years.
+ The time stamps of the output files are
+\end_layout
+
+\begin_layout LyX-Code
+# in years (rather than the default of seconds), and we give a format for
+\end_layout
+
+\begin_layout LyX-Code
+# the time stamp.
+\end_layout
+
+\begin_layout LyX-Code
+output_freq = time_step
+\end_layout
+
+\begin_layout LyX-Code
+time_step = 40.0*year
+\end_layout
+
+\begin_layout LyX-Code
+writer.filename = output/step04.vtk
+\end_layout
+
+\begin_layout LyX-Code
+writer.time_format = %04.0f
+\end_layout
+
+\begin_layout LyX-Code
+writer.time_constant = 1.0*year
+\end_layout
+
+\begin_layout LyX-Code
+
+\end_layout
+
+\begin_layout LyX-Code
+# Give basename for VTK domain output of solution over ground surface.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.problem.formulation.output.subdomain]
+\end_layout
+
+\begin_layout LyX-Code
+# Name of nodeset for ground surface.
+\end_layout
+
+\begin_layout LyX-Code
+label = face_zpos
+\end_layout
+
+\begin_layout LyX-Code
+# We keep the default output frequency behavior (skip every n steps), and
+\end_layout
+
+\begin_layout LyX-Code
+# ask to skip 0 steps between output, so that we get output every time step.
+\end_layout
+
+\begin_layout LyX-Code
+skip = 0
+\end_layout
+
+\begin_layout LyX-Code
+writer.filename = output/step04-groundsurf.vtk
+\end_layout
+
+\begin_layout LyX-Code
+writer.time_format = %04.0f
+\end_layout
+
+\begin_layout LyX-Code
+writer.time_constant = 1.0*year
+\end_layout
+
+\begin_layout Standard
+We provide similar output information for the two materials (
+\family typewriter
+upper_crust
+\family default
+ and 
+\family typewriter
+lower_crust
+\family default
+).
+ Note that for the domain output, we requested output in terms of a given
+ time frequency, while for the subdomain we requested output in terms of
+ number of time steps.
+ When we have run the simulation, the output VTK files will be contained
+ in 
+\family typewriter
+examples/3d/hex8/output
+\family default
+ (all with a prefix of 
+\family typewriter
+step04
+\family default
+), and the results may be visualized using a tool such as ParaView or mayavi2.
+ Results using ParaView are shown in 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "fig:step04-displ-t200"
+
+\end_inset
+
+.
+\end_layout
+
+\begin_layout Standard
+\begin_inset Float figure
+wide false
+sideways false
+status open
+
+\begin_layout Plain Layout
+\begin_inset Graphics
+	filename figs/step04-displ-t200.jpg
+	width 10cm
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Plain Layout
+\begin_inset Caption
+
+\begin_layout Plain Layout
+Displacement field for example step04 at t = 200 years visualized using
+ ParaView.
+ The mesh has been distorted by the computed displacements (magnified by
+ 500), and the vectors show the computed displacements.
+\begin_inset CommandInset label
+LatexCommand label
+name "fig:step04-displ-t200"
+
+\end_inset
+
+
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Plain Layout
+
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Subsubsection
+Step05 - Time-varying Dirichlet and Neumann Boundary Conditions
+\end_layout
+
+\begin_layout Standard
+The 
+\family typewriter
+step05.cfg
+\family default
+ file describes a problem with time-varying Dirichlet and Neumann boundary
+ conditions.
+ The example is similar to example step04, with a few important differences:
+\end_layout
+
+\begin_layout Itemize
+The Dirichlet boundary conditions on the negative x-face include an initial
+ displacement (applied in the elastic solution), as well as a constant velocity.
+\end_layout
+
+\begin_layout Itemize
+Neumann (traction) boundary conditions are applied in the negative x-direction
+ on the positive x-face, giving a compressive stress.
+ An initial traction is applied in the elastic solution, and then at t =
+ 100 years it begins decreasing linearly until it reaches zero at the end
+ of the simulation (t = 200 years).
+\end_layout
+
+\begin_layout Standard
+We again use a Maxwell viscoelastic material for the lower crust.
+\end_layout
+
+\begin_layout Standard
+For the boundary conditions, we must first change the boundary condition
+ type for the positive x-face from the default Dirichlet to Neumann:
+\end_layout
+
+\begin_layout LyX-Code
+# +x face -- first change bc type to Neumann
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.bc]
+\end_layout
+
+\begin_layout LyX-Code
+x_pos = pylith.bc.Neumann 
+\end_layout
+
+\begin_layout Standard
+We provide quadrature information for this face as we did for example step02.
+ We then use a 
+\family typewriter
+UniformDB
+\family default
+ for both the initial tractions as well as the traction rates.
+ We provide a start time of 100 years for the traction rates, and use a
+ rate of 0.01 MPa/year, so that by the end of 200 years we have completely
+ cancelled the initial traction of -1 MPa:
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.bc.x_pos]
+\end_layout
+
+\begin_layout LyX-Code
+# Fist specify a UniformDB for the initial tractions, along with the values.
+\end_layout
+
+\begin_layout LyX-Code
+db_initial = spatialdata.spatialdb.UniformDB
+\end_layout
+
+\begin_layout LyX-Code
+db_initial.label = Neumann BC on +x
+\end_layout
+
+\begin_layout LyX-Code
+db_initial.values = [traction-shear-horiz,traction-shear-vert,traction-normal]
+\end_layout
+
+\begin_layout LyX-Code
+db_initial.data = [0.0*MPa,0.0*MPa,-1.0*MPa]
+\end_layout
+
+\begin_layout LyX-Code
+# Provide information on traction rates.
+\end_layout
+
+\begin_layout LyX-Code
+db_rate = spatialdata.spatialdb.UniformDB
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.label = Neumann rate BC on +x
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.values = [traction-rate-shear-horiz,traction-rate-shear-vert,traction-rat
+e-normal,rate-start-time]
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.data = [0.0*MPa/year,0.0*MPa/year,0.01*MPa/year,100.0*year]
+\end_layout
+
+\begin_layout Standard
+The boundary conditions on the negative x-face are analogous, but we are
+ instead using Dirichlet boundary conditions, and the initial displacement
+ is in the same direction as the applied velocities:
+\end_layout
+
+\begin_layout LyX-Code
+# -x face
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.bc.x_neg]
+\end_layout
+
+\begin_layout LyX-Code
+bc_dof = [0, 1]
+\end_layout
+
+\begin_layout LyX-Code
+label = face_xneg
+\end_layout
+
+\begin_layout LyX-Code
+# Initial displacements.
+\end_layout
+
+\begin_layout LyX-Code
+db_initial = spatialdata.spatialdb.UniformDB
+\end_layout
+
+\begin_layout LyX-Code
+db_initial.label = Dirichlet BC on -x
+\end_layout
+
+\begin_layout LyX-Code
+db_initial.values = [displacement-x,displacement-y]
+\end_layout
+
+\begin_layout LyX-Code
+db_initial.data = [0.0*cm,-0.5*cm]
+\end_layout
+
+\begin_layout LyX-Code
+# Velocities.
+\end_layout
+
+\begin_layout LyX-Code
+db_rate = spatialdata.spatialdb.UniformDB
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.label = Dirichlet rate BC on -x
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.values = [displacement-rate-x,displacement-rate-y,rate-start-time]
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.data = [0.0*cm/year,-1.0*cm/year,0.0*year]
+\end_layout
+
+\begin_layout Standard
+The boundary conditions on the negative z-face are supplied in the same
+ manner as for example step04.
+ When we have run the simulation, the output VTK files will be contained
+ in 
+\family typewriter
+examples/3d/hex8/output
+\family default
+ (all with a prefix of 
+\family typewriter
+step05
+\family default
+), and the results may be visualized using a tool such as ParaView or mayavi2.
+ Results using ParaView are shown in 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "fig:step05-displ-t40"
+
+\end_inset
+
+.
+\begin_inset Float figure
+wide false
+sideways false
+status open
+
+\begin_layout Plain Layout
+\begin_inset Graphics
+	filename figs/step05-displ-t40.jpg
+	width 10cm
+
+\end_inset
+
+
+\begin_inset Caption
+
+\begin_layout Plain Layout
+Displacement field for example step05 at t = 40 years visualized using ParaView.
+ The mesh has been distorted by the computed displacements (magnified by
+ 500), and the vectors show the computed displacements.
+\begin_inset CommandInset label
+LatexCommand label
+name "fig:step05-displ-t40"
+
+\end_inset
+
+.
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Plain Layout
+
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Subsubsection
+Step06 - Dirichlet Boundary Conditions with Time-Dependent Kinematic Fault
+ Slip
+\end_layout
+
+\begin_layout Standard
+The 
+\family typewriter
+step06.cfg
+\family default
+ file describes a problem with Dirichlet (displacement) boundary conditions
+ corresponding to zero x and y-displacements applied on the negative and
+ positive x-faces and a vertical fault that includes multiple earthquake
+ ruptures as well as steady fault creep.
+ The upper (locked) portion of the fault has 4 m of left-lateral slip every
+ 200 years, while the lower (creeping) portion of the fault slips at a steady
+ rate of 2 cm/year.
+ The problem bears some similarity to the strike-slip fault model of 
+\begin_inset CommandInset citation
+LatexCommand cite
+key "Savage:Prescott:1978"
+
+\end_inset
+
+, except that the fault creep extends through the viscoelastic portion of
+ the domain, and the far-field displacement boundary conditions are held
+ fixed.
+\end_layout
+
+\begin_layout Standard
+In this example and the remainder of the examples in this section, we change
+ the time stepping behavior from the default 
+\family typewriter
+TimeStepUniform
+\family default
+ to 
+\family typewriter
+TimeStepAdapt
+\family default
+.
+ For adaptive time stepping, we provide the maximum permissible time step
+ size, along with a stability factor.
+ The stability factor controls how the time step size is allowed to deviate
+ from the stable time step size provided by the different materials in the
+ model.
+ A 
+\family typewriter
+stability_factor
+\family default
+ of 1.0 means we should use the stable time step size, while a 
+\family typewriter
+stability_factor
+\family default
+ greater than 1.0 means we want to use a smaller time step size.
+ A 
+\family typewriter
+stability_factor
+\family default
+ less than 1.0 allows time step sizes greater than the stable time step size,
+ which may provide inaccurate results.
+ The adaptive time stepping information is provided as:
+\end_layout
+
+\begin_layout LyX-Code
+# Change time stepping algorithm from uniform time step, to adaptive
+\end_layout
+
+\begin_layout LyX-Code
+# time stepping.
+\end_layout
+
+\begin_layout LyX-Code
+time_step = pylith.problems.TimeStepAdapt
+\end_layout
+
+\begin_layout LyX-Code
+
+\end_layout
+
+\begin_layout LyX-Code
+# Change the total simulation time to 700 years, and set the maximum time
+\end_layout
+
+\begin_layout LyX-Code
+# step size to 10 years.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.implicit.time_step]
+\end_layout
+
+\begin_layout LyX-Code
+total_time = 700.0*year
+\end_layout
+
+\begin_layout LyX-Code
+max_dt = 10.0*year
+\end_layout
+
+\begin_layout LyX-Code
+stability_factor = 1.0 ; use time step equal to stable value from materials
+\end_layout
+
+\begin_layout Standard
+Due to the simplicity of the boundary conditions, we are able to use the
+ default 
+\family typewriter
+ZeroDispBC
+\family default
+ for the positive and negative x-faces, as well as the negative z-face.
+ As for example step03, we define a fault interface, we identify the nodeset
+ corresponding to the fault, and we provide quadrature information for the
+ fault.
+ We then define an array of earthquake sources and provide an origin time
+ for each:
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.interfaces.fault]
+\end_layout
+
+\begin_layout LyX-Code
+# Set earthquake sources to an array consisting of creep and 3 ruptures.
+\end_layout
+
+\begin_layout LyX-Code
+eq_srcs = [creep,one,two,three]
+\end_layout
+
+\begin_layout LyX-Code
+eq_srcs.creep.origin_time = 00.0*year
+\end_layout
+
+\begin_layout LyX-Code
+eq_srcs.one.origin_time = 200.0*year
+\end_layout
+
+\begin_layout LyX-Code
+eq_srcs.two.origin_time = 400.0*year
+\end_layout
+
+\begin_layout LyX-Code
+eq_srcs.three.origin_time = 600.0*year
+\end_layout
+
+\begin_layout Standard
+Note that the creep begins at t = 0 years, while the ruptures (
+\family typewriter
+one
+\family default
+, 
+\family typewriter
+two
+\family default
+, 
+\family typewriter
+three
+\family default
+) occur at regular intervals of 200 years.
+ We retain the default 
+\family typewriter
+StepSlipFn
+\family default
+ for the ruptures.
+ Each of the ruptures has the same amount of slip, and slip occurs simultaneousl
+y for the entire rupture region, so we can use the same 
+\family typewriter
+SimpleDB
+\family default
+ files providing slip and slip time for each rupture:
+\end_layout
+
+\begin_layout LyX-Code
+# Define slip and origin time for first rupture.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.interfaces.fault.eq_srcs.one.slip_function]
+\end_layout
+
+\begin_layout LyX-Code
+slip.iohandler.filename = spatialdb/finalslip_rupture.spatialdb
+\end_layout
+
+\begin_layout LyX-Code
+slip_time.iohandler.filename = spatialdb/sliptime.spatialdb
+\end_layout
+
+\begin_layout LyX-Code
+
+\end_layout
+
+\begin_layout LyX-Code
+# Define slip and origin time for second rupture.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.interfaces.fault.eq_srcs.two.slip_function]
+\end_layout
+
+\begin_layout LyX-Code
+slip.iohandler.filename = spatialdb/finalslip_rupture.spatialdb
+\end_layout
+
+\begin_layout LyX-Code
+slip_time.iohandler.filename = spatialdb/sliptime.spatialdb
+\end_layout
+
+\begin_layout LyX-Code
+
+\end_layout
+
+\begin_layout LyX-Code
+# Define slip and origin time for third rupture.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.interfaces.fault.eq_srcs.three.slip_function]
+\end_layout
+
+\begin_layout LyX-Code
+slip.iohandler.filename = spatialdb/finalslip_rupture.spatialdb
+\end_layout
+
+\begin_layout LyX-Code
+slip_time.iohandler.filename = spatialdb/sliptime.spatialdb
+\end_layout
+
+\begin_layout Standard
+For the creep source, we change the slip function to 
+\family typewriter
+ConstRateSlipFn
+\family default
+, and we use a 
+\family typewriter
+SimpleDB
+\family default
+ for both the slip time and the slip rate:
+\end_layout
+
+\begin_layout LyX-Code
+# Define slip rate and origin time for fault creep.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.interfaces.fault.eq_srcs.creep]
+\end_layout
+
+\begin_layout LyX-Code
+slip_function = pylith.faults.ConstRateSlipFn
+\end_layout
+
+\begin_layout LyX-Code
+slip_function.slip_rate.iohandler.filename = spatialdb/sliprate_creep.spatialdb
+\end_layout
+
+\begin_layout LyX-Code
+slip_function.slip_time.iohandler.filename = spatialdb/sliptime.spatialdb
+\end_layout
+
+\begin_layout Standard
+Note that for all earthquake sources we provide both an 
+\family typewriter
+origin_time
+\family default
+ and a 
+\family typewriter
+slip_function.slip_time
+\family default
+.
+ The first provides the starting time for the entire earthquake source,
+ while the second provides any spatial variation in the slip time with respect
+ to the 
+\family typewriter
+origin_time
+\family default
+ (if any).
+ Since there are multiple earthquake sources of different types, there are
+ a number of additional fault information fields available for output.
+ We therefore request that these additional fields be output to the fault
+ information file:
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.interfaces.fault]
+\end_layout
+
+\begin_layout LyX-Code
+output.vertex_info_fields = [normal_dir,strike_dir,dip_dir,final_slip_creep,final
+_slip_one,final_slip_two,final_slip_three,slip_time_creep,slip_time_one,slip_tim
+e_two,slip_time_three]
+\end_layout
+
+\begin_layout Standard
+This additional information will be contained in file 
+\family typewriter
+step06-fault_info.vtk
+\family default
+.
+ It will contain final slip information for each earthquake source along
+ with slip time information.
+ When we have run the simulation, the output VTK files will be contained
+ in 
+\family typewriter
+examples/3d/hex8/output
+\family default
+ (all with a prefix of 
+\family typewriter
+step06
+\family default
+), and the results may be visualized using a tool such as ParaView or mayavi2.
+ Results using ParaView are shown in 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "fig:step06-displ-t300"
+
+\end_inset
+
+.
+\begin_inset Float figure
+wide false
+sideways false
+status open
+
+\begin_layout Plain Layout
+\begin_inset Graphics
+	filename figs/step06-displ-t300.jpg
+	width 10cm
+
+\end_inset
+
+
+\begin_inset Caption
+
+\begin_layout Plain Layout
+Displacement field for example step06 at t = 300 years visualized using
+ ParaView.
+ The mesh has been distorted by the computed displacements (magnified by
+ 500), and the vectors show the computed displacements.
+\begin_inset CommandInset label
+LatexCommand label
+name "fig:step06-displ-t300"
+
+\end_inset
+
+.
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Subsubsection
+Step07 - Dirichlet Velocity Boundary Conditions with Time-Dependent Kinematic
+ Fault Slip
+\end_layout
+
+\begin_layout Standard
+The 
+\family typewriter
+step07.cfg
+\family default
+ file describes a problem nearly identical to example step06.
+ The only difference is that velocity boundary conditions in the positive
+ and negative y-directions on the positive and negative x-faces, so that
+ the external boundaries keep pace with the average fault slip.
+ This problem is nearly identical to the strike-slip fault model of 
+\begin_inset CommandInset citation
+LatexCommand cite
+key "Savage:Prescott:1978"
+
+\end_inset
+
+, except that the fault creep extends through the viscoelastic portion of
+ the domain.
+\end_layout
+
+\begin_layout Standard
+We use the default 
+\family typewriter
+ZeroDispBC
+\family default
+ for the initial displacements on the positive and negative x-faces, as
+ well as the negative z-face.
+ For the velocities on the positive and negative x-faces, we use a 
+\family typewriter
+UniformDB
+\family default
+:
+\end_layout
+
+\begin_layout LyX-Code
+# +x face
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.bc.x_pos]
+\end_layout
+
+\begin_layout LyX-Code
+bc_dof = [0, 1]
+\end_layout
+
+\begin_layout LyX-Code
+label = face_xpos
+\end_layout
+
+\begin_layout LyX-Code
+db_initial.label = Dirichlet BC on +x
+\end_layout
+
+\begin_layout LyX-Code
+db_rate = spatialdata.spatialdb.UniformDB
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.label = Dirichlet rate BC on +x
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.values = [displacement-rate-x,displacement-rate-y,rate-start-time]
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.data = [0.0*cm/year,1.0*cm/year,0.0*year]
+\end_layout
+
+\begin_layout LyX-Code
+
+\end_layout
+
+\begin_layout LyX-Code
+# -x face
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.bc.x_neg]
+\end_layout
+
+\begin_layout LyX-Code
+bc_dof = [0, 1]
+\end_layout
+
+\begin_layout LyX-Code
+label = face_xneg
+\end_layout
+
+\begin_layout LyX-Code
+db_initial.label = Dirichlet BC on -x
+\end_layout
+
+\begin_layout LyX-Code
+db_rate = spatialdata.spatialdb.UniformDB
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.label = Dirichlet rate BC on +x
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.values = [displacement-rate-x,displacement-rate-y,rate-start-time]
+\end_layout
+
+\begin_layout LyX-Code
+db_rate.data = [0.0*cm/year,-1.0*cm/year,0.0*year]
+\end_layout
+
+\begin_layout Standard
+The fault definition information is identical to example step06.
+ In previous examples, we have just used the default output for the domain
+ and subdomain (ground surface), which just provides the displacements.
+ In many cases, it is also useful to provide the velocities.
+ PyLith is able to provide this information, approximating the velocities
+ for the current time step as the difference between the current displacements
+ and the displacements from the previous time step, divided by the time
+ step size.
+ We can obtain this information by explicitly requesting it in 
+\family typewriter
+vertex_data_fields
+\family default
+:
+\end_layout
+
+\begin_layout LyX-Code
+# Give basename for VTK domain output of solution over domain.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.problem.formulation.output.domain]
+\end_layout
+
+\begin_layout LyX-Code
+# We specify that output occurs in terms of a given time frequency, and
+\end_layout
+
+\begin_layout LyX-Code
+# ask for output every 50 years.
+ The time stamps of the output files are
+\end_layout
+
+\begin_layout LyX-Code
+# in years (rather than the default of seconds), and we give a format for
+\end_layout
+
+\begin_layout LyX-Code
+# the time stamp.
+\end_layout
+
+\begin_layout LyX-Code
+# We also request velocity output in addition to displacements.
+\end_layout
+
+\begin_layout LyX-Code
+vertex_data_fields = [displacement,velocity]
+\end_layout
+
+\begin_layout LyX-Code
+output_freq = time_step
+\end_layout
+
+\begin_layout LyX-Code
+time_step = 50.0*year
+\end_layout
+
+\begin_layout LyX-Code
+writer.filename = output/step07.vtk
+\end_layout
+
+\begin_layout LyX-Code
+writer.time_format = %04.0f
+\end_layout
+
+\begin_layout LyX-Code
+writer.time_constant = 1.0*year
+\end_layout
+
+\begin_layout LyX-Code
+
+\end_layout
+
+\begin_layout LyX-Code
+# Give basename for VTK domain output of solution over ground surface.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.problem.formulation.output.subdomain]
+\end_layout
+
+\begin_layout LyX-Code
+# Name of nodeset for ground surface.
+\end_layout
+
+\begin_layout LyX-Code
+label = face_zpos
+\end_layout
+
+\begin_layout LyX-Code
+# We also request velocity output in addition to displacements.
+\end_layout
+
+\begin_layout LyX-Code
+vertex_data_fields = [displacement,velocity]
+\end_layout
+
+\begin_layout LyX-Code
+# We keep the default output frequency behavior (skip every n steps), and
+\end_layout
+
+\begin_layout LyX-Code
+# ask to skip 0 steps between output, so that we get output every time step.
+\end_layout
+
+\begin_layout LyX-Code
+skip = 0
+\end_layout
+
+\begin_layout LyX-Code
+writer.filename = output/step07-groundsurf.vtk
+\end_layout
+
+\begin_layout LyX-Code
+writer.time_format = %04.0f
+\end_layout
+
+\begin_layout LyX-Code
+writer.time_constant = 1.0*year
+\end_layout
+
+\begin_layout Standard
+When we have run the simulation, the output VTK files will be contained
+ in 
+\family typewriter
+examples/3d/hex8/output
+\family default
+ (all with a prefix of 
+\family typewriter
+step07
+\family default
+), and the results may be visualized using a tool such as ParaView or mayavi2.
+ Results using ParaView are shown in 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "fig:step07-displ-vel-t300"
+
+\end_inset
+
+.
+\begin_inset Float figure
+wide false
+sideways false
+status open
+
+\begin_layout Plain Layout
+\begin_inset Graphics
+	filename figs/step07-displ-vel-t300.jpg
+	width 10cm
+
+\end_inset
+
+
+\begin_inset Caption
+
+\begin_layout Plain Layout
+Displacement field (color contours) and velocity field (vectors) for example
+ step07 at t = 300 years visualized using ParaView.
+ The mesh has been distorted by the computed displacements (magnified by
+ 500), and the vectors show the computed velocities.
+\begin_inset CommandInset label
+LatexCommand label
+name "fig:step07-displ-vel-t300"
+
+\end_inset
+
+.
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Subsubsection
+Step08 - Dirichlet Velocity Boundary Conditions with Time-Dependent Kinematic
+ Fault Slip and Power-Law Rheology
+\begin_inset CommandInset label
+LatexCommand label
+name "sub:Tutorial-Step08-Power-law"
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Standard
+The 
+\family typewriter
+step08.cfg
+\family default
+ file describes a problem that is identical to example step07, except the
+ the lower crust is composed of power-law viscoelastic material.
+ Since the material behavior is now nonlinear, we must use the nonlinear
+ solver:
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent]
+\end_layout
+
+\begin_layout LyX-Code
+# For this problem we must switch to a nonlinear solver.
+\end_layout
+
+\begin_layout LyX-Code
+implicit.solver = pylith.problems.SolverNonlinear
+\end_layout
+
+\begin_layout Standard
+Although we have not discussed the PyLith PETSc settings previously, note
+ that the use of the nonlinear solver may require additional options if
+ we wish to override the defaults.
+ These settings are contained in 
+\family typewriter
+pylithapp.cfg
+\family default
+:
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.petsc]
+\end_layout
+
+\begin_layout LyX-Code
+# Nonlinear solver monitoring options.
+\end_layout
+
+\begin_layout LyX-Code
+snes_rtol = 1.0e-8
+\end_layout
+
+\begin_layout LyX-Code
+snes_atol = 1.0e-12
+\end_layout
+
+\begin_layout LyX-Code
+snes_max_it = 100
+\end_layout
+
+\begin_layout LyX-Code
+snes_monitor = true
+\end_layout
+
+\begin_layout LyX-Code
+snes_view = true
+\end_layout
+
+\begin_layout LyX-Code
+snes_converged_reason = true
+\end_layout
+
+\begin_layout Standard
+These settings are unused unless we are using the nonlinear solver.
+\end_layout
+
+\begin_layout Standard
+When using laboratory results for power-law viscoelastic materials, the
+ properties provided do not generally match those expected by PyLith (see
+ 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "sub:Power-Law-Maxwell-Viscoelastic"
+
+\end_inset
+
+).
+ A utility code, 
+\family typewriter
+powerlaw_gendb.py
+\family default
+, has been provided to simplify the process of using laboratory results
+ with PyLith.
+ This utility code should be installed in the same location as PyLith.
+ An example of how to use it is in 
+\family typewriter
+examples/3d/hex8/spatialdb/powerlaw
+\family default
+.
+ The user must provide a spatial database defining the spatial distribution
+ of laboratory-derived parameters (contained in 
+\family typewriter
+powerlaw_params.spatialdb
+\family default
+), another spatial database defining the temperature field in degrees K
+ (contained in 
+\family typewriter
+temperature.spatialdb
+\family default
+), and a set of points for which values are desired (
+\family typewriter
+powerlaw_points.txt
+\family default
+).
+ The parameters for the code are defined in 
+\family typewriter
+powerlaw_gendb.cfg
+\family default
+.
+ The properties expected by PyLith are 
+\family typewriter
+reference_strain_rate
+\family default
+, 
+\family typewriter
+reference_stress
+\family default
+, and 
+\family typewriter
+power_law_exponent
+\family default
+.
+ The user must specify either 
+\family typewriter
+reference_strain_rate
+\family default
+ or 
+\family typewriter
+reference_stress
+\family default
+ so that 
+\family typewriter
+powerlaw_gendb.py
+\family default
+ can compute the other property.
+ Default values of 1.0e-6 1/s and 1 MPa are provided.
+ A spatial database must then be provided for each laboratory-derived parameter
+ and for the temperature.
+ In this example, the same database was used for all parameters, and a separate
+ database was used to define the temperature distribution.
+ In practice, the user can provide any desired thermal model to provide
+ the spatial database for the temperature.
+ In this example, a simple 1D (vertically-varying) distribution was used.
+ The utility code can be used by simply executing it from the 
+\family typewriter
+examples/3d/hex8/spatialdb/powerlaw
+\family default
+ directory:
+\end_layout
+
+\begin_layout LyX-Code
+powerlaw_gendb.py
+\end_layout
+
+\begin_layout Standard
+This will cause the code to automatically read the parameters in powerlaw_gendb.c
+fg, and the file 
+\family typewriter
+examples/3d/hex8/spatialdb/mat_powerlaw.spatialdb
+\family default
+ will be created.
+\end_layout
+
+\begin_layout Standard
+We first change the material type of the lower crust to 
+\family typewriter
+PowerLaw3D
+\family default
+:
+\end_layout
+
+\begin_layout LyX-Code
+# Change material type of lower crust to power-law viscoelastic.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent]
+\end_layout
+
+\begin_layout LyX-Code
+materials.lower_crust = pylith.materials.PowerLaw3D
+\end_layout
+
+\begin_layout Standard
+In many cases, it is useful to obtain the material properties from two different
+ sources.
+ For example, the elastic properties may come from a seismic velocity model
+ while the viscous properties may be derived from a thermal model.
+ In such a case we can use a 
+\family typewriter
+CompositeDB
+\family default
+, which allows a different spatial database to be used for a subset of the
+ properties.
+ We do this as follows:
+\end_layout
+
+\begin_layout LyX-Code
+# Provide a spatial database from which to obtain property values.
+\end_layout
+
+\begin_layout LyX-Code
+# In this case, we prefer to obtain the power-law properties from one
+\end_layout
+
+\begin_layout LyX-Code
+# database and the elastic properties from another database, so we use
+\end_layout
+
+\begin_layout LyX-Code
+# a CompositeDB.
+ Each part of the CompositeDB is a SimpleDB.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.materials.lower_crust]
+\end_layout
+
+\begin_layout LyX-Code
+db_properties = spatialdata.spatialdb.CompositeDB
+\end_layout
+
+\begin_layout LyX-Code
+db_properties.db_A = spatialdata.spatialdb.SimpleDB
+\end_layout
+
+\begin_layout LyX-Code
+db_properties.db_B = spatialdata.spatialdb.SimpleDB
+\end_layout
+
+\begin_layout Standard
+Each subset of properties comes from a 
+\family typewriter
+SimpleDB
+\family default
+, so we must first define the properties that come from each spatial database
+ and then provide the database filename:
+\end_layout
+
+\begin_layout LyX-Code
+# Provide the values to be obtained from each database and the database
+\end_layout
+
+\begin_layout LyX-Code
+# name.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.materials.lower_crust.db_properties]
+\end_layout
+
+\begin_layout LyX-Code
+values_A = [density,vs,vp]   ; Elastic properties.
+\end_layout
+
+\begin_layout LyX-Code
+db_A.label = Elastic properties
+\end_layout
+
+\begin_layout LyX-Code
+db_A.iohandler.filename = spatialdb/mat_elastic.spatialdb
+\end_layout
+
+\begin_layout LyX-Code
+values_B = [reference-stress,reference-strain-rate,power-law-exponent] 
+  ; Power-law properties.
+\end_layout
+
+\begin_layout LyX-Code
+db_B.label = Power-law properties
+\end_layout
+
+\begin_layout LyX-Code
+db_B.iohandler.filename = spatialdb/mat_powerlaw.spatialdb
+\end_layout
+
+\begin_layout Standard
+The 
+\family typewriter
+PowerLaw3D
+\family default
+ material has additional properties and state variables with respect to
+ the default 
+\family typewriter
+ElasticIsotropic3D
+\family default
+ material, so we request that these properties be written to the 
+\family typewriter
+lower_crust
+\family default
+ material files:
+\end_layout
+
+\begin_layout LyX-Code
+# Since there are additional properties and state variables for the
+\end_layout
+
+\begin_layout LyX-Code
+# power-law model, we explicitly request that they be output.
+ Properties are
+\end_layout
+
+\begin_layout LyX-Code
+# named in cell_info_fields and state variables are named in
+\end_layout
+
+\begin_layout LyX-Code
+# cell_data_fields.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.materials.lower_crust]
+\end_layout
+
+\begin_layout LyX-Code
+output.cell_info_fields = [density,mu,lambda,reference_strain_rate,reference_stre
+ss,power_law_exponent]
+\end_layout
+
+\begin_layout LyX-Code
+output.cell_data_fields = [total_strain,stress,viscous_strain]
+\end_layout
+
+\begin_layout Standard
+When we have run the simulation, the output VTK files will be contained
+ in 
+\family typewriter
+examples/3d/hex8/output
+\family default
+ (all with a prefix of 
+\family typewriter
+step08
+\family default
+), and the results may be visualized using a tool such as ParaView or mayavi2.
+ Results using ParaView are shown in 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "fig:step08-strain-displ-t150"
+
+\end_inset
+
+.
+\begin_inset Float figure
+wide false
+sideways false
+status open
+
+\begin_layout Plain Layout
+\begin_inset Graphics
+	filename figs/step08-strain-displ-t150.jpg
+	width 10cm
+
+\end_inset
+
+
+\begin_inset Caption
+
+\begin_layout Plain Layout
+The XY-component of strain (color contours) and displacement field (vectors)
+ for example step08 at t = 150 years visualized using ParaView.
+ For this visualization, we loaded both the 
+\family typewriter
+step08-lower_crust_txxxx.vtk
+\family default
+ and 
+\family typewriter
+step08-upper_crust_txxxx.vtk
+\family default
+ files to contour the strain field, and superimposed on it the displacement
+ field vectors from 
+\family typewriter
+step08_txxxx.vtk
+\family default
+.
+\begin_inset CommandInset label
+LatexCommand label
+name "fig:step08-strain-displ-t150"
+
+\end_inset
+
+.
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\begin_layout Subsubsection
+Step09 - Dirichlet Velocity Boundary Conditions with Time-Dependent Kinematic
+ Fault Slip and Drucker-Prager Elastoplastic Rheology
+\end_layout
+
+\begin_layout Standard
+The 
+\family typewriter
+step09.cfg
+\family default
+ file describes a problem that is nearly identical to example step08, except
+ the the lower crust is composed of Drucker-Prager elastoplastic material.
+ As for example step08, the material behavior is nonlinear so we again use
+ the nonlinear solver.
+ Since the material is elastoplastic, there is no inherent time-dependence,
+ which means that there is no stable time step size for the material.
+ Depending on the loading conditions, this may make it difficult to obtain
+ convergence.
+ To avoid this, we set the maximum time step size to 5 years rather than
+ the value of 10 years used in example step08:
+\end_layout
+
+\begin_layout LyX-Code
+# Change the total simulation time to 700 years, and set the maximum time
+\end_layout
+
+\begin_layout LyX-Code
+# step size to 5 years.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.implicit.time_step]
+\end_layout
+
+\begin_layout LyX-Code
+total_time = 700.0*year
+\end_layout
+
+\begin_layout LyX-Code
+max_dt = 5.0*year
+\end_layout
+
+\begin_layout LyX-Code
+stability_factor = 1.0 ; use time step equal to stable value from materials
+\end_layout
+
+\begin_layout LyX-Code
+# For this problem we set adapt_skip to zero so that the time step size
+ is
+\end_layout
+
+\begin_layout LyX-Code
+# readjusted every time step.
+\end_layout
+
+\begin_layout LyX-Code
+adapt_skip = 0
+\end_layout
+
+\begin_layout Standard
+We change the material type of the lower crust to 
+\family typewriter
+DruckerPrager3D
+\family default
+, and we again use a 
+\family typewriter
+CompositeDB
+\family default
+ to assign the material properties:
+\end_layout
+
+\begin_layout LyX-Code
+# Change material type of lower crust to Drucker-Prager.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent]
+\end_layout
+
+\begin_layout LyX-Code
+materials.lower_crust = pylith.materials.DruckerPrager3D
+\end_layout
+
+\begin_layout LyX-Code
+
+\end_layout
+
+\begin_layout LyX-Code
+# Provide a spatial database from which to obtain property values.
+\end_layout
+
+\begin_layout LyX-Code
+# In this case, we prefer to obtain the Drucker-Prager properties from one
+\end_layout
+
+\begin_layout LyX-Code
+# database and the elastic properties from another database, so we use
+\end_layout
+
+\begin_layout LyX-Code
+# a CompositeDB.
+ Each part of the CompositeDB is a SimpleDB.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.materials.lower_crust]
+\end_layout
+
+\begin_layout LyX-Code
+db_properties = spatialdata.spatialdb.CompositeDB
+\end_layout
+
+\begin_layout LyX-Code
+db_properties.db_A = spatialdata.spatialdb.SimpleDB
+\end_layout
+
+\begin_layout LyX-Code
+db_properties.db_B = spatialdata.spatialdb.SimpleDB
+\end_layout
+
+\begin_layout Standard
+As for the step08 example, we first define the properties that come from
+ each spatial database and then provide the database filename:
+\end_layout
+
+\begin_layout LyX-Code
+# Provide the values to be obtained from each database and the database
+\end_layout
+
+\begin_layout LyX-Code
+# name.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.materials.lower_crust.db_properties]
+\end_layout
+
+\begin_layout LyX-Code
+values_A = [density,vs,vp]   ; Elastic properties.
+\end_layout
+
+\begin_layout LyX-Code
+db_A.label = Elastic properties
+\end_layout
+
+\begin_layout LyX-Code
+db_A.iohandler.filename = spatialdb/mat_elastic.spatialdb
+\end_layout
+
+\begin_layout LyX-Code
+values_B = [friction-angle,cohesion,dilatation-angle]   ; Drucker-Prager
+ properties.
+\end_layout
+
+\begin_layout LyX-Code
+db_B.label = Drucker-Prager properties
+\end_layout
+
+\begin_layout LyX-Code
+db_B.iohandler.filename = spatialdb/mat_druckerprager.spatialdb
+\end_layout
+
+\begin_layout Standard
+We also request output of the properties and state variables that are unique
+ to the 
+\family typewriter
+DruckerPrager3D
+\family default
+ material:
+\end_layout
+
+\begin_layout LyX-Code
+# Since there are additional properties and state variables for the
+\end_layout
+
+\begin_layout LyX-Code
+# Drucker-Prager model, we explicitly request that they be output.
+\end_layout
+
+\begin_layout LyX-Code
+# Properties are named in cell_info_fields and state variables are named
+ in
+\end_layout
+
+\begin_layout LyX-Code
+# cell_data_fields.
+\end_layout
+
+\begin_layout LyX-Code
+[pylithapp.timedependent.materials.lower_crust]
+\end_layout
+
+\begin_layout LyX-Code
+output.cell_info_fields = [density,mu,lambda,alpha_yield,beta,alpha_flow]
+\end_layout
+
+\begin_layout LyX-Code
+output.cell_data_fields = [total_strain,stress,plastic_strain]
+\end_layout
+
+\begin_layout Standard
+When we have run the simulation, the output VTK files will be contained
+ in 
+\family typewriter
+examples/3d/hex8/output
+\family default
+ (all with a prefix of 
+\family typewriter
+step09
+\family default
+), and the results may be visualized using a tool such as ParaView or mayavi2.
+ Results using ParaView are shown in 
+\begin_inset CommandInset ref
+LatexCommand ref
+reference "fig:step09-strain-displ-t150"
+
+\end_inset
+
+.
+\begin_inset Float figure
+wide false
+sideways false
+status open
+
+\begin_layout Plain Layout
+\begin_inset Graphics
+	filename /Users/charlesw/geoframe/cig/short/3D/PyLith/trunk/doc/userguide/tutorials/3dhex8/quasistatic/figs/step08-strain-displ-t150.jpg
+	width 10cm
+
+\end_inset
+
+
+\begin_inset Caption
+
+\begin_layout Plain Layout
+The XY-component of strain (color contours) and displacement field (vectors)
+ for example step09 at t = 150 years visualized using ParaView.
+ For this visualization, we loaded both the 
+\family typewriter
+step09-lower_crust_txxxx.vtk
+\family default
+ and 
+\family typewriter
+step09-upper_crust_txxxx.vtk
+\family default
+ files to contour the strain field, and superimposed on it the displacement
+ field vectors from 
+\family typewriter
+step09_txxxx.vtk
+\family default
+.
+\begin_inset CommandInset label
+LatexCommand label
+name "fig:step09-strain-displ-t150"
+
+\end_inset
+
+.
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\end_inset
+
+
+\end_layout
+
+\end_body
+\end_document

More information about the CIG-COMMITS mailing list