Python ElasticMaterialApp.ElasticMaterialApp类代码示例

  def __init__(self, name="powerlaw3dtimedep"):
    """
    Constructor.
    """
    ElasticMaterialApp.__init__(self, name)

    # import pdb
    # pdb.set_trace()

    numLocs = 2

    self.dimension = dimension
    self.numLocs = numLocs

    self.dbPropertyValues = ["density", "vs", "vp",
                             "reference-strain-rate", "reference-stress",
                             "power-law-exponent"]
    self.propertyValues = ["density", "mu", "lambda",
                           "reference-strain-rate", "reference-stress",
                           "power-law-exponent"]
    self.numPropertyValues = numpy.array([1, 1, 1, 1, 1, 1], dtype=numpy.int32)

    self.dbStateVarValues = ["viscous-strain-xx",
                             "viscous-strain-yy",
                             "viscous-strain-zz",
                             "viscous-strain-xy",
                             "viscous-strain-yz",
                             "viscous-strain-xz",
                             "stress-xx",
                             "stress-yy",
                             "stress-zz",
                             "stress-xy",
                             "stress-yz",
                             "stress-xz",
                             ]
    self.stateVarValues = ["viscous-strain", "stress"]
    self.numStateVarValues = numpy.array([6, 6], dtype=numpy.int32)

    self.alpha = 0.5
    self.dt = 2.0e5

    densityA = 2500.0
    vsA = 3000.0
    vpA = vsA*3**0.5
    powerLawCoeffA = 1.0/3.0e18
    refStrainRateA = 1.0e-6
    powerLawExponentA = 1.0
    strainA = [1.1e-4, 1.2e-4, 1.3e-4, 1.4e-4, 1.5e-4, 1.6e-4]
    initialStressA = [2.1e4, 2.2e4, 2.3e4, 2.4e4, 2.5e4, 2.6e4]
    initialStrainA = [3.6e-5, 3.5e-5, 3.4e-5, 3.3e-5, 3.2e-5, 3.1e-5]
    muA = vsA*vsA*densityA
    lambdaA = vpA*vpA*densityA - 2.0*muA
    viscosityCoeffA = (1.0/((3.0**0.5)**(powerLawExponentA + 1.0) \
                            * powerLawCoeffA))**(1.0/powerLawExponentA)
    refStressA = viscosityCoeffA * \
                 (2.0 * refStrainRateA) ** (1.0/powerLawExponentA)

    densityB = 2000.0
    vsB = 1200.0
    vpB = vsB*3**0.5
    powerLawCoeffB = 1.0/9.0e36
    refStrainRateB = 1.0e-6
    powerLawExponentB = 3.0
    strainB = [4.1e-4, 4.2e-4, 4.3e-4, 4.4e-4, 4.5e-4, 4.6e-4]
    initialStressB = [5.1e4, 5.2e4, 5.3e4, 5.4e4, 5.5e4, 5.6e4]
    initialStrainB = [6.1e-5, 6.2e-5, 6.3e-5, 6.6e-5, 6.5e-5, 6.4e-5]
    muB = vsB*vsB*densityB
    lambdaB = vpB*vpB*densityB - 2.0*muB
    viscosityCoeffB = (1.0/((3.0**0.5)**(powerLawExponentB + 1.0) \
                            * powerLawCoeffB))**(1.0/powerLawExponentB)
    refStressB = viscosityCoeffB * \
                 (2.0 * refStrainRateB) ** (1.0/powerLawExponentB)

    self.lengthScale = 1.0e+3
    self.pressureScale = muA
    self.timeScale = 1.0
    self.densityScale = muA / (self.lengthScale / self.timeScale)**2
    self.strainRateScale = 1.0/self.timeScale

    self.dbProperties = numpy.array([ [densityA, vsA, vpA, \
                                       refStrainRateA, refStressA, \
                                       powerLawExponentA],
                                      [densityB, vsB, vpB, \
                                       refStrainRateB, refStressB, \
                                       powerLawExponentB] ], 
                                    dtype=numpy.float64)
    self.properties = numpy.array([ [densityA, muA, lambdaA, \
                                     refStrainRateA, refStressA, \
                                     powerLawExponentA],
                                    [densityB, muB, lambdaB, \
                                     refStrainRateB, refStressB, \
                                     powerLawExponentB] ],
                                  dtype=numpy.float64)

    # TEMPORARY, need to determine how to use initial state variables
    self.dbStateVars = numpy.zeros( (numLocs, tensorSize),
                                    dtype=numpy.float64)

    mu0 = self.pressureScale
    density0 = self.densityScale
#.........这里部分代码省略.........

  def __init__(self, name="genmaxwellplanestrainelastic"):
    """
    Constructor.
    """
    ElasticMaterialApp.__init__(self, name)

    # import pdb
    # pdb.set_trace()
    
    numLocs = 2

    self.dimension = dimension
    self.numLocs = numLocs

    self.dbPropertyValues = ["density", "vs", "vp",
                             "shear-ratio-1", "shear-ratio-2", "shear-ratio-3",
                             "viscosity-1", "viscosity-2", "viscosity-3"]
    self.numPropertyValues = numpy.array([1, 1, 1, 1, 1, 1, 1, 1, 1],
                                         dtype=numpy.int32)

    self.dbStateVarValues = ["stress-zz-initial",
                             "total-strain-xx",
                             "total-strain-yy",
                             "total-strain-xy",
                             "viscous-strain-1-xx",
                             "viscous-strain-1-yy",
                             "viscous-strain-1-zz",
                             "viscous-strain-1-xy",
                             "viscous-strain-2-xx",
                             "viscous-strain-2-yy",
                             "viscous-strain-2-zz",
                             "viscous-strain-2-xy",
                             "viscous-strain-3-xx",
                             "viscous-strain-3-yy",
                             "viscous-strain-3-zz",
                             "viscous-strain-3-xy"
                             ]
    self.numStateVarValues = numpy.array([1, tensorSize, 4, 4, 4],
                                         dtype=numpy.int32)

    densityA = 2500.0
    vsA = 3000.0
    vpA = vsA*3**0.5
    shearRatioA = [0.5, 0.1, 0.2]
    viscosityA = [1.0e18, 1.0e17, 1.0e19]
    strainA = [1.1e-4, 1.2e-4, 1.4e-4]
    initialStressA = [2.1e4, 2.2e4, 2.4e4]
    initialStrainA = [3.1e-5, 3.2e-5, 3.4e-5]
    muA = vsA*vsA*densityA
    lambdaA = vpA*vpA*densityA - 2.0*muA
    stressInitialZZA = numpy.array([1.5e4], dtype=numpy.float64)


    densityB = 2000.0
    vsB = 1200.0
    vpB = vsB*3**0.5
    shearRatioB = [0.2, 0.2, 0.2]
    viscosityB = [1.0e18, 1.0e19, 1.0e20]
    strainB = [4.1e-4, 4.2e-4, 4.4e-4]
    initialStressB = [5.1e4, 5.2e4, 5.4e4]
    initialStrainB = [6.1e-5, 6.2e-5, 6.4e-5]
    muB = vsB*vsB*densityB
    lambdaB = vpB*vpB*densityB - 2.0*muB
    stressInitialZZB = numpy.array([4.5e4], dtype=numpy.float64)

    maxwellTimeA = [1.0e30, 1.0e30, 1.0e30]
    maxwellTimeB = [1.0e30, 1.0e30, 1.0e30]
    for i in xrange(numMaxwellModels):
      if shearRatioA[i] != 0.0:
        maxwellTimeA[i] = viscosityA[i]/(muA*shearRatioA[i])
      if shearRatioB[i] != 0.0:
        maxwellTimeB[i] = viscosityB[i]/(muB*shearRatioB[i])

    self.lengthScale = 1.0e+3
    self.pressureScale = muA
    self.timeScale = 1.0
    self.densityScale = muA / (self.lengthScale / self.timeScale)**2

    propA = [densityA, vsA, vpA] + shearRatioA + viscosityA
    propB = [densityB, vsB, vpB] + shearRatioB + viscosityB
    self.dbProperties = numpy.array([propA, propB], dtype=numpy.float64)
    propA = [densityA, muA, lambdaA] + shearRatioA + maxwellTimeA
    propB = [densityB, muB, lambdaB] + shearRatioB + maxwellTimeB
    self.properties = numpy.array([propA, propB], dtype=numpy.float64)

    # TEMPORARY, need to determine how to use initial state variables
    # At present, only the first (stressInitialZZ) is being used.
    self.dbStateVars = numpy.zeros( (numLocs,
                                     1 + tensorSize + 4 * numMaxwellModels),
                                    dtype=numpy.float64)
    self.dbStateVars[0, 0] = stressInitialZZA
    self.dbStateVars[1, 0] = stressInitialZZB
    
    self.stateVars = numpy.zeros( (numLocs,
                                   1 + tensorSize + 4 * numMaxwellModels),
                                  dtype=numpy.float64)
    self.stateVars[0, 0] = stressInitialZZA
    self.stateVars[1, 0] = stressInitialZZB

    mu0 = self.pressureScale
#.........这里部分代码省略.........

  def __init__(self, name="druckerpragerplanestrainelastic"):
    """
    Constructor.
    """
    ElasticMaterialApp.__init__(self, name)

    # import pdb
    # pdb.set_trace()
    numLocs = 2

    self.dimension = dimension
    self.numLocs = numLocs

    self.dbPropertyValues = ["density", "vs", "vp",
                             "friction-angle", "cohesion",
                             "dilatation-angle"]
    self.numPropertyValues = numpy.array([1, 1, 1, 1, 1, 1], dtype=numpy.int32)

    self.dbStateVarValues = ["stress-zz-initial",
                             "plastic-strain-xx",
                             "plastic-strain-yy",
                             "plastic-strain-zz",
                             "plastic-strain-xy"
                             ]
    self.numStateVarValues = numpy.array([1, 4], dtype=numpy.int32)

    densityA = 2500.0
    vsA = 3000.0
    vpA = vsA*3**0.5
    # First case has different values for friction angle and dilatation angle.
    frictionAngleA = math.radians(30.0)
    dilatationAngleA = math.radians(20.0)
    cohesionA = 3.0e5
    strainA = [1.1e-4, 1.2e-4, 1.4e-4]
    initialStressA = [2.1e4, 2.2e4, 2.4e4]
    initialStrainA = [3.1e-4, 3.2e-4, 3.4e-4]
    muA = vsA*vsA*densityA
    lambdaA = vpA*vpA*densityA - 2.0*muA
    stressZZInitialA = 1.075e+4

    denomFrictionA = math.sqrt(3.0) * (3.0 - math.sin(frictionAngleA))
    denomDilatationA = math.sqrt(3.0) * (3.0 - math.sin(dilatationAngleA))
    alphaYieldA = 2.0 * math.sin(frictionAngleA)/denomFrictionA
    betaA = 6.0 * cohesionA * math.cos(frictionAngleA)/denomFrictionA
    alphaFlowA = 2.0 * math.sin(dilatationAngleA)/denomDilatationA
    
    densityB = 2000.0
    vsB = 1200.0
    vpB = vsB*3**0.5
    # Second case has same values for friction angle and dilatation angle.
    frictionAngleB = math.radians(25.0)
    dilatationAngleB = math.radians(25.0)
    cohesionB = 1.0e5
    strainB = [4.1e-4, 4.2e-4, 4.4e-4]
    initialStressB = [5.1e4, 5.2e4, 5.4e4]
    initialStrainB = [6.1e-4, 6.2e-4, 6.4e-4]
    muB = vsB*vsB*densityB
    lambdaB = vpB*vpB*densityB - 2.0*muB
    denomFrictionB = math.sqrt(3.0) * (3.0 - math.sin(frictionAngleB))
    denomDilatationB = math.sqrt(3.0) * (3.0 - math.sin(dilatationAngleB))
    alphaYieldB = 2.0 * math.sin(frictionAngleB)/denomFrictionB
    betaB = 6.0 * cohesionB * math.cos(frictionAngleB)/denomFrictionB
    alphaFlowB = 2.0 * math.sin(dilatationAngleB)/denomDilatationB
    stressZZInitialB = 2.575e+4

    self.lengthScale = 1.0e+3
    self.pressureScale = muA
    self.timeScale = 1.0
    self.densityScale = muA / (self.lengthScale / self.timeScale)**2

    self.dbProperties = numpy.array([ [densityA, vsA, vpA, \
                                       frictionAngleA, cohesionA, \
                                       dilatationAngleA],
                                      [densityB, vsB, vpB, \
                                       frictionAngleB, cohesionB, \
                                       dilatationAngleB] ], 
                                    dtype=numpy.float64)
    self.properties = numpy.array([ [densityA, muA, lambdaA, \
                                     alphaYieldA, betaA, \
                                     alphaFlowA],
                                    [densityB, muB, lambdaB, \
                                     alphaYieldB, betaB, \
                                     alphaFlowB] ],
                                     dtype=numpy.float64)

    # TEMPORARY, need to determine how to use initial state variables
    self.dbStateVars = numpy.zeros( (numLocs, 1 + 4), dtype=numpy.float64)
    self.dbStateVars[0, 0] = stressZZInitialA
    self.dbStateVars[1, 0] = stressZZInitialB

    self.stateVars = numpy.zeros( (numLocs, 1 + 4), dtype=numpy.float64)
    self.stateVars[0, 0] = stressZZInitialA
    self.stateVars[1, 0] = stressZZInitialB

    mu0 = self.pressureScale
    density0 = self.densityScale
    self.propertiesNondim = \
        numpy.array([ [densityA/density0, muA/mu0, lambdaA/mu0, \
                       alphaYieldA, betaA/mu0, \
                       alphaFlowA],
#.........这里部分代码省略.........

  def __init__(self, name="maxwellisotropic3dtimedep"):
    """
    Constructor.
    """
    ElasticMaterialApp.__init__(self, name)

    # import pdb
    # pdb.set_trace()

    numLocs = 2

    self.dimension = dimension
    self.numLocs = numLocs

    self.dbPropertyValues = ["density", "vs", "vp", "viscosity"]
    self.propertyValues = ["density", "mu", "lambda", "maxwellTime"]
    self.numPropertyValues = numpy.array([1, 1, 1, 1], dtype=numpy.int32)

    self.dbStateVarValues = ["total-strain-xx",
                             "total-strain-yy",
                             "total-strain-zz",
                             "total-strain-xy",
                             "total-strain-yz",
                             "total-strain-xz",
                             "viscous-strain-xx",
                             "viscous-strain-yy",
                             "viscous-strain-zz",
                             "viscous-strain-xy",
                             "viscous-strain-yz",
                             "viscous-strain-xz",
                             ]
    self.stateVarValues = ["total-strain", "viscous-strain"]
    self.numStateVarValues = numpy.array([6, 6], dtype=numpy.int32)

    self.dt = 2.0e5

    densityA = 2500.0
    vsA = 3000.0
    vpA = vsA*3**0.5
    viscosityA = 1.0e18
    strainA = [1.1e-4, 1.2e-4, 1.3e-4, 1.4e-4, 1.5e-4, 1.6e-4]
    initialStressA = [2.1e4, 2.2e4, 2.3e4, 2.4e4, 2.5e4, 2.6e4]
    initialStrainA = [3.6e-5, 3.5e-5, 3.4e-5, 3.3e-5, 3.2e-5, 3.1e-5]
    muA = vsA*vsA*densityA
    lambdaA = vpA*vpA*densityA - 2.0*muA
    maxwellTimeA = viscosityA / muA
    meanStrainA = (strainA[0] + strainA[1] + strainA[2])/3.0

    densityB = 2000.0
    vsB = 1200.0
    vpB = vsB*3**0.5
    viscosityB = 1.0e19
    strainB = [4.1e-4, 4.2e-4, 4.3e-4, 4.4e-4, 4.5e-4, 4.6e-4]
    initialStressB = [5.1e4, 5.2e4, 5.3e4, 5.4e4, 5.5e4, 5.6e4]
    initialStrainB = [6.1e-5, 6.2e-5, 6.3e-5, 6.6e-5, 6.5e-5, 6.4e-5]
    muB = vsB*vsB*densityB
    lambdaB = vpB*vpB*densityB - 2.0*muB
    maxwellTimeB = viscosityB / muB
    meanStrainB = (strainB[0] + strainB[1] + strainB[2])/3.0

    diag = numpy.array([1.0, 1.0, 1.0, 0.0, 0.0, 0.0],
                       dtype=numpy.float64)

    self.lengthScale = 1.0e+3
    self.pressureScale = muA
    self.timeScale = 1.0
    self.densityScale = muA / (self.lengthScale / self.timeScale)**2

    self.dbProperties = numpy.array([ [densityA, vsA, vpA, viscosityA],
                                      [densityB, vsB, vpB, viscosityB] ], 
                                    dtype=numpy.float64)
    self.properties = numpy.array([ [densityA, muA, lambdaA, maxwellTimeA],
                                    [densityB, muB, lambdaB, maxwellTimeB] ],
                                     dtype=numpy.float64)

    # TEMPORARY, need to determine how to use initial state variables
    self.dbStateVars = numpy.zeros( (numLocs, tensorSize),
                                    dtype=numpy.float64)

    mu0 = self.pressureScale
    density0 = self.densityScale
    time0 = self.timeScale
    self.propertiesNondim = \
        numpy.array([ [densityA/density0, muA/mu0, lambdaA/mu0, \
                       maxwellTimeA/time0],
                      [densityB/density0, muB/mu0, lambdaB/mu0, \
                       maxwellTimeB/time0] ],
                    dtype=numpy.float64)

    self.initialStress = numpy.array([initialStressA,
                                      initialStressB],
                                    dtype=numpy.float64)
    self.initialStrain = numpy.array([initialStrainA,
                                      initialStrainB],
                                    dtype=numpy.float64)
    
    self.density = numpy.array([densityA,
                                densityB],
                               dtype=numpy.float64)

#.........这里部分代码省略.........

  def __init__(self, name="druckerprager3dtimedep"):
    """
    Constructor.
    """
    ElasticMaterialApp.__init__(self, name)

    numLocs = 2

    self.dimension = dimension
    self.numLocs = numLocs

    self.dbPropertyValues = ["density", "vs", "vp",
                             "friction-angle", "cohesion",
                             "dilatation-angle"]
    self.propertyValues = ["density", "mu", "lambda",
                           "alpha_yield", "beta", "alpha_flow"]
    self.numPropertyValues = numpy.array([1, 1, 1, 1, 1, 1], dtype=numpy.int32)

    self.dbStateVarValues = ["plastic-strain-xx",
                             "plastic-strain-yy",
                             "plastic-strain-zz",
                             "plastic-strain-xy",
                             "plastic-strain-yz",
                             "plastic-strain-xz"
                             ]
    self.stateVarValues = ["plastic-strain"]
    self.numStateVarValues = numpy.array([6], dtype=numpy.int32)

    self.dt = 2.0e5

    densityA = 2500.0
    vsA = 3000.0
    vpA = vsA*3**0.5
    # First case has same values for friction angle and dilatation angle.
    frictionAngleA = math.radians(30.0)
    dilatationAngleA = math.radians(20.0)
    cohesionA = 3.0e5
    strainA = [-2.1e-4, 1.2e-4, 1.3e-4, 1.1e-5, 1.1e-5, 1.1e-5]
    initialStressA = [2.1e4, 2.2e4, 2.3e4, 2.4e4, 2.5e4, 2.6e4]
    initialStrainA = [3.6e-5, 3.5e-5, 3.4e-5, 3.3e-5, 3.2e-5, 3.1e-5]
    # initialStressA = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0]
    # initialStrainA = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0]
    muA = vsA*vsA*densityA
    lambdaA = vpA*vpA*densityA - 2.0*muA

    denomFrictionA = math.sqrt(3.0) * (3.0 - math.sin(frictionAngleA))
    denomDilatationA = math.sqrt(3.0) * (3.0 - math.sin(dilatationAngleA))
    alphaYieldA = 2.0 * math.sin(frictionAngleA)/denomFrictionA
    betaA = 6.0 * cohesionA * math.cos(frictionAngleA)/denomFrictionA
    alphaFlowA = 2.0 * math.sin(dilatationAngleA)/denomDilatationA
    
    densityB = 2000.0
    vsB = 1200.0
    vpB = vsB*3**0.5
    # Second case has different values for friction angle and dilatation angle.
    frictionAngleB = math.radians(25.0)
    dilatationAngleB = math.radians(25.0)
    # frictionAngleB = 0.0
    # dilatationAngleB = 0.0
    cohesionB = 1.0e4
    strainB = [4.1e-4, 4.2e-4, 4.3e-4, 1.4e-4, 1.5e-4, 1.6e-4]
    initialStressB = [5.6e4, 5.5e4, 5.4e4, 5.3e4, 5.2e4, 5.1e4]
    initialStrainB = [6.6e-5, 6.5e-5, 6.4e-5, 6.3e-5, 6.2e-5, 6.1e-5]
    # initialStressB = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0]
    # initialStrainB = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0]
    muB = vsB*vsB*densityB
    lambdaB = vpB*vpB*densityB - 2.0*muB
    denomFrictionB = math.sqrt(3.0) * (3.0 - math.sin(frictionAngleB))
    denomDilatationB = math.sqrt(3.0) * (3.0 - math.sin(dilatationAngleB))
    alphaYieldB = 2.0 * math.sin(frictionAngleB)/denomFrictionB
    betaB = 6.0 * cohesionB * math.cos(frictionAngleB)/denomFrictionB
    alphaFlowB = 2.0 * math.sin(dilatationAngleB)/denomDilatationB

    self.lengthScale = 1.0e+3
    self.pressureScale = muA
    self.timeScale = 1.0
    self.densityScale = muA / (self.lengthScale / self.timeScale)**2

    self.dbProperties = numpy.array([ [densityA, vsA, vpA, \
                                       frictionAngleA, cohesionA, \
                                       dilatationAngleA],
                                      [densityB, vsB, vpB, \
                                       frictionAngleB, cohesionB, \
                                       dilatationAngleB] ], 
                                    dtype=numpy.float64)
    self.properties = numpy.array([ [densityA, muA, lambdaA, \
                                     alphaYieldA, betaA, \
                                     alphaFlowA],
                                    [densityB, muB, lambdaB, \
                                     alphaYieldB, betaB, \
                                     alphaFlowB] ],
                                  dtype=numpy.float64)

    # TEMPORARY, need to determine how to use initial state variables
    self.dbStateVars = numpy.zeros( (numLocs, tensorSize), dtype=numpy.float64)

    mu0 = self.pressureScale
    density0 = self.densityScale
    self.propertiesNondim = \
        numpy.array([ [densityA/density0, muA/mu0, lambdaA/mu0, \
#.........这里部分代码省略.........

  def __init__(self, name="elasticisotropic3d"):
    """
    Constructor.
    """
    ElasticMaterialApp.__init__(self, name)

    numLocs = 2

    self.dimension = dimension
    self.numLocs = numLocs
    
    self.dbPropertyValues = ["density", "vs", "vp"]    
    self.propertyValues = ["density", "mu", "lambda"]
    self.numPropertyValues = numpy.array([1, 1, 1], dtype=numpy.int32)

    densityA = 2500.0
    vsA = 3000.0
    vpA = vsA*3**0.5
    strainA = [1.1e-4, 1.2e-4, 1.3e-4, 1.4e-4, 1.5e-4, 1.6e-4]
    initialStressA = [2.1e4, 2.2e4, 2.3e4, 2.4e4, 2.5e4, 2.6e4]
    initialStrainA = [3.1e-4, 3.2e-4, 3.3e-4, 3.4e-4, 3.5e-6, 3.6e-4]
    muA = vsA*vsA*densityA
    lambdaA = vpA*vpA*densityA - 2.0*muA
    
    densityB = 2000.0
    vsB = 1200.0
    vpB = vsB*3**0.5
    strainB = [4.1e-4, 4.2e-4, 4.3e-4, 4.4e-4, 4.5e-4, 4.6e-4]
    initialStressB = [5.1e4, 5.2e4, 5.3e4, 5.4e4, 5.5e4, 5.6e4]
    initialStrainB = [6.1e-4, 6.2e-4, 6.3e-4, 6.4e-4, 6.5e-6, 6.6e-4]
    muB = vsB*vsB*densityB
    lambdaB = vpB*vpB*densityB - 2.0*muB

    self.lengthScale = 1.0e+3
    self.pressureScale = muA
    self.timeScale = 1.0
    self.densityScale = muA / (self.lengthScale / self.timeScale)**2
    
    self.dbProperties = numpy.array([ [densityA, vsA, vpA],
                                      [densityB, vsB, vpB] ], 
                                    dtype=numpy.float64)
    self.properties = numpy.array([ [densityA, muA, lambdaA],
                                    [densityB, muB, lambdaB] ],
                                     dtype=numpy.float64)

    mu0 = self.pressureScale
    density0 = self.densityScale
    self.propertiesNondim = \
        numpy.array([ [densityA/density0, muA/mu0, lambdaA/mu0],
                      [densityB/density0, muB/mu0, lambdaB/mu0] ],
                    dtype=numpy.float64)

    self.initialStress = numpy.array([initialStressA,
                                      initialStressB],
                                    dtype=numpy.float64)
    self.initialStrain = numpy.array([initialStrainA,
                                      initialStrainB],
                                    dtype=numpy.float64)
    
    self.density = numpy.array([densityA,
                                densityB],
                               dtype=numpy.float64)

    self.strain = numpy.array([strainA,
                               strainB],
                               dtype=numpy.float64)
    
    stress = numpy.zeros( (numLocs, tensorSize), dtype=numpy.float64)
    elasticConsts = numpy.zeros( (numLocs, numElasticConsts),
                                 dtype=numpy.float64)

    (elasticConsts[0,:], stress[0,:]) = \
        self._calcStress(strainA, densityA, muA, lambdaA,
                         initialStressA, initialStrainA)
    (elasticConsts[1,:], stress[1,:]) = \
        self._calcStress(strainB, densityB, muB, lambdaB,
                         initialStressB, initialStrainB)

    self.stress = stress
    self.elasticConsts = elasticConsts

    self.dtStableExplicit = 1000.0 / vpA

    return

    def __init__(self, name="powerlaw3delastic"):
        """
    Constructor.
    """
        ElasticMaterialApp.__init__(self, name)

        # import pdb
        # pdb.set_trace()
        numLocs = 2

        self.dimension = dimension
        self.numLocs = numLocs

        self.dbPropertyValues = [
            "density",
            "vs",
            "vp",
            "reference-strain-rate",
            "reference-stress",
            "power-law-exponent",
        ]
        self.numPropertyValues = numpy.array([1, 1, 1, 1, 1, 1], dtype=numpy.int32)

        self.dbStateVarValues = [
            "viscous-strain-xx",
            "viscous-strain-yy",
            "viscous-strain-zz",
            "viscous-strain-xy",
            "viscous-strain-yz",
            "viscous-strain-xz",
            "stress-xx",
            "stress-yy",
            "stress-zz",
            "stress-xy",
            "stress-yz",
            "stress-xz",
        ]
        self.numStateVarValues = numpy.array([6, 6], dtype=numpy.int32)

        densityA = 2500.0
        vsA = 3000.0
        vpA = vsA * 3 ** 0.5
        # Derive new values in based on previous value for power-law coefficient
        # and viscosity coefficient.
        powerLawCoeffA = 1.0 / 3.0e18
        refStrainRateA = 1.0e-6
        powerLawExponentA = 1.0
        strainA = [1.1e-4, 1.2e-4, 1.3e-4, 1.4e-4, 1.5e-4, 1.6e-4]
        initialStressA = [2.1e4, 2.2e4, 2.3e4, 2.4e4, 2.5e4, 2.6e4]
        initialStrainA = [3.1e-4, 3.2e-4, 3.3e-4, 3.4e-4, 3.5e-4, 3.6e-4]
        muA = vsA * vsA * densityA
        lambdaA = vpA * vpA * densityA - 2.0 * muA

        viscosityCoeffA = (1.0 / ((3.0 ** 0.5) ** (powerLawExponentA + 1.0) * powerLawCoeffA)) ** (
            1.0 / powerLawExponentA
        )
        refStressA = viscosityCoeffA * (2.0 * refStrainRateA) ** (1.0 / powerLawExponentA)
        # refStressA = (refStrainRateA/powerLawCoeffA)**(1.0/powerLawExponentA)

        densityB = 2000.0
        vsB = 1200.0
        vpB = vsB * 3 ** 0.5
        powerLawCoeffB = 1.0 / 9.0e30
        refStrainRateB = 1.0e-6
        powerLawExponentB = 3.0
        strainB = [4.1e-4, 4.2e-4, 4.3e-4, 4.4e-4, 4.5e-4, 4.6e-4]
        initialStressB = [5.1e4, 5.2e4, 5.3e4, 5.4e4, 5.5e4, 5.6e4]
        initialStrainB = [6.1e-4, 6.2e-4, 6.3e-4, 6.4e-4, 6.5e-4, 6.6e-4]
        muB = vsB * vsB * densityB
        lambdaB = vpB * vpB * densityB - 2.0 * muB
        viscosityCoeffB = (1.0 / ((3.0 ** 0.5) ** (powerLawExponentB + 1.0) * powerLawCoeffB)) ** (
            1.0 / powerLawExponentB
        )
        refStressB = viscosityCoeffB * (2.0 * refStrainRateB) ** (1.0 / powerLawExponentB)
        # refStressB = (refStrainRateB/powerLawCoeffB)**(1.0/powerLawExponentB)

        self.lengthScale = 1.0e3
        self.pressureScale = muA
        self.timeScale = 1.0
        self.densityScale = muA / (self.lengthScale / self.timeScale) ** 2
        self.strainRateScale = 1.0 / self.timeScale

        self.dbProperties = numpy.array(
            [
                [densityA, vsA, vpA, refStrainRateA, refStressA, powerLawExponentA],
                [densityB, vsB, vpB, refStrainRateB, refStressB, powerLawExponentB],
            ],
            dtype=numpy.float64,
        )
        self.properties = numpy.array(
            [
                [densityA, muA, lambdaA, refStrainRateA, refStressA, powerLawExponentA],
                [densityB, muB, lambdaB, refStrainRateB, refStressB, powerLawExponentB],
            ],
            dtype=numpy.float64,
        )

        # TEMPORARY, need to determine how to use initial state variables
        self.dbStateVars = numpy.zeros((numLocs, tensorSize + tensorSize), dtype=numpy.float64)
        self.stateVars = numpy.zeros((numLocs, tensorSize + tensorSize), dtype=numpy.float64)
#.........这里部分代码省略.........

    def __init__(self, name="genmaxwellqpqsisotropic3dtimedep"):
        """
    Constructor.
    """
        ElasticMaterialApp.__init__(self, name)

        numLocs = 2
        self.dt = 2.0e5
        self.dimension = dimension
        self.numLocs = numLocs

        self.dbPropertyValues = [
            "density",
            "vs",
            "vp",
            "shear-ratio-1",
            "shear-ratio-2",
            "shear-ratio-3",
            "shear-viscosity-1",
            "shear-viscosity-2",
            "shear-viscosity-3",
            "bulk-ratio-1",
            "bulk-ratio-2",
            "bulk-ratio-3",
            "bulk-viscosity-1",
            "bulk-viscosity-2",
            "bulk-viscosity-3",
        ]
        self.numPropertyValues = numpy.array([1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1], dtype=numpy.int32)

        self.dbStateVarValues = [
            "total-strain-xx",
            "total-strain-yy",
            "total-strain-zz",
            "total-strain-xy",
            "total-strain-yz",
            "total-strain-xz",
            "viscousdeviatoric-strain-1-xx",
            "viscousdeviatoric-strain-1-yy",
            "viscousdeviatoric-strain-1-zz",
            "viscousdeviatoric-strain-1-xy",
            "viscousdeviatoric-strain-1-yz",
            "viscousdeviatoric-strain-1-xz",
            "viscousdeviatoric-strain-2-xx",
            "viscousdeviatoric-strain-2-yy",
            "viscousdeviatoric-strain-2-zz",
            "viscousdeviatoric-strain-2-xy",
            "viscousdeviatoric-strain-2-yz",
            "viscousdeviatoric-strain-2-xz",
            "viscousdeviatoric-strain-3-xx",
            "viscousdeviatoric-strain-3-yy",
            "viscousdeviatoric-strain-3-zz",
            "viscousdeviatoric-strain-3-xy",
            "viscousdeviatoric-strain-3-yz",
            "viscousdeviatoric-strain-3-xz",
            "viscous-mean-strain-1",
            "viscous-mean-strain-2",
            "viscous-mean-strain-3",
        ]
        #    self.numStateVarValues = numpy.array([tensorSize]*(1+numMaxwellModels) + numMaxwellModels,
        #                                         dtype=numpy.int32)

        self.numStateVarValues = numpy.array(
            [1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1], dtype=numpy.int32
        )

        densityA = 2500.0
        vsA = 3000.0
        vpA = vsA * 3 ** 0.5
        shearRatioA = [0.5, 0.1, 0.2]
        shearViscosityA = [1.0e18, 1.0e17, 1.0e19]
        bulkRatioA = [0.4, 0.3, 0.1]
        bulkViscosityA = [2.0e18, 2.0e17, 2.0e19]
        strainA = [1.1e-4, 2.2e-4, 3.3e-4, 4.4e-4, 5.5e-4, 6.6e-4]
        initialStressA = [2.1e4, 2.2e4, 2.3e4, 2.4e4, 2.5e4, 2.6e4]
        # initialStrainA = [3.1e-4, 3.2e-4, 3.3e-4, 3.4e-4, 3.5e-4, 3.6e-4]
        muA = vsA * vsA * densityA
        lambdaA = vpA * vpA * densityA - 2.0 * muA
        kA = lambdaA + 2 / 3.0 * muA

        densityB = 2000.0
        vsB = 1200.0
        vpB = vsB * 3 ** 0.5
        shearRatioB = [0.2, 0.2, 0.2]
        shearViscosityB = [1.0e18, 1.0e19, 1.0e20]
        bulkRatioB = [0.2, 0.2, 0.2]
        bulkViscosityB = [1.0e18, 1.0e19, 1.0e20]
        strainB = [1.2e-4, 2.3e-4, 3.4e-4, 4.5e-4, 5.6e-4, 6.7e-4]
        initialStressB = [5.1e4, 5.2e4, 5.3e4, 5.4e4, 5.5e4, 5.6e4]
        # initialStrainB = [6.1e-4, 6.2e-4, 6.3e-4, 6.4e-4, 6.5e-4, 6.6e-4]
        muB = vsB * vsB * densityB
        lambdaB = vpB * vpB * densityB - 2.0 * muB
        kB = lambdaB + 2 / 3.0 * muB

        # TEMPORARY, need to determine how to use initial strain
        initialStrainA = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0]
        initialStrainB = [0.0, 0.0, 0.0, 0.0, 0.0, 0.0]

        # Simplest approach for now is to assume this is the first step after the
        # elastic solution. In that case, both the total strain from the last step
#.........这里部分代码省略.........

  def __init__(self, name="genmaxwellisotropic3dtimedep"):
    """
    Constructor.
    """
    ElasticMaterialApp.__init__(self, name)

    # import pdb
    # pdb.set_trace()

    numLocs = 2
    self.dt = 2.0e5
    self.dimension = dimension
    self.numLocs = numLocs

    self.dbPropertyValues = ["density", "vs", "vp",
                             "shear-ratio-1", "shear-ratio-2", "shear-ratio-3",
                             "viscosity-1", "viscosity-2", "viscosity-3"]
    self.numPropertyValues = numpy.array([1, 1, 1, 1, 1, 1, 1, 1, 1],
                                         dtype=numpy.int32)

    self.dbStateVarValues = ["total-strain-xx",
                             "total-strain-yy",
                             "total-strain-zz",
                             "total-strain-xy",
                             "total-strain-yz",
                             "total-strain-xz",
                             "viscous-strain-1-xx",
                             "viscous-strain-1-yy",
                             "viscous-strain-1-zz",
                             "viscous-strain-1-xy",
                             "viscous-strain-1-yz",
                             "viscous-strain-1-xz",
                             "viscous-strain-2-xx",
                             "viscous-strain-2-yy",
                             "viscous-strain-2-zz",
                             "viscous-strain-2-xy",
                             "viscous-strain-2-yz",
                             "viscous-strain-2-xz",
                             "viscous-strain-3-xx",
                             "viscous-strain-3-yy",
                             "viscous-strain-3-zz",
                             "viscous-strain-3-xy",
                             "viscous-strain-3-yz",
                             "viscous-strain-3-xz",
                             ]
    self.numStateVarValues = numpy.array([tensorSize]*(1+numMaxwellModels),
                                         dtype=numpy.int32)

    densityA = 2500.0
    vsA = 3000.0
    vpA = vsA*3**0.5
    shearRatioA = [0.5, 0.1, 0.2]
    viscosityA = [1.0e18, 1.0e17, 1.0e19]
    strainA = [1.1e-4, 2.2e-4, 3.3e-4, 4.4e-4, 5.5e-4, 6.6e-4]
    initialStressA = [2.1e4, 2.2e4, 2.3e4, 2.4e4, 2.5e4, 2.6e4]
    initialStrainA = [3.1e-5, 3.2e-5, 3.3e-5, 3.4e-5, 3.5e-5, 3.6e-5]
    muA = vsA*vsA*densityA
    lambdaA = vpA*vpA*densityA - 2.0*muA
    meanStrainA = (strainA[0] + strainA[1] + strainA[2])/3.0

    densityB = 2000.0
    vsB = 1200.0
    vpB = vsB*3**0.5
    shearRatioB = [0.2, 0.2, 0.2]
    viscosityB = [1.0e18, 1.0e19, 1.0e20]
    strainB = [1.2e-4, 2.3e-4, 3.4e-4, 4.5e-4, 5.6e-4, 6.7e-4]
    initialStressB = [5.1e4, 5.2e4, 5.3e4, 5.4e4, 5.5e4, 5.6e4]
    initialStrainB = [6.1e-5, 6.2e-5, 6.3e-5, 6.4e-5, 6.5e-5, 6.6e-5]
    muB = vsB*vsB*densityB
    lambdaB = vpB*vpB*densityB - 2.0*muB
    meanStrainB = (strainB[0] + strainB[1] + strainB[2])/3.0

    # Compute Maxwell time and viscous strain.
    diag = numpy.array([1.0, 1.0, 1.0, 0.0, 0.0, 0.0], dtype=numpy.float64)
    strainTA = numpy.array(strainA, dtype=numpy.float64)
    strainTB = numpy.array(strainB, dtype=numpy.float64)

    maxwellTimeA = [1.0e30, 1.0e30, 1.0e30]
    maxwellTimeB = [1.0e30, 1.0e30, 1.0e30]
    visStrainA = numpy.zeros( (numMaxwellModels, tensorSize),
                              dtype=numpy.float64)
    visStrainB = numpy.zeros( (numMaxwellModels, tensorSize),
                              dtype=numpy.float64)
    for imodel in xrange(numMaxwellModels):
      if shearRatioA[imodel] != 0.0:
        maxwellTimeA[imodel] = viscosityA[imodel]/(muA*shearRatioA[imodel])
        visStrainA[imodel,:] = strainTA[:] - diag[:] * meanStrainA
      if shearRatioB[imodel] != 0.0:
        maxwellTimeB[imodel] = viscosityB[imodel]/(muB*shearRatioB[imodel])
        visStrainB[imodel,:] = strainTB[:] - diag[:] * meanStrainB

    dbPropA = [densityA, vsA, vpA] + shearRatioA + viscosityA
    dbPropB = [densityB, vsB, vpB] + shearRatioB + viscosityB
    self.dbProperties = numpy.array([dbPropA, dbPropB], dtype=numpy.float64)
    propA = [densityA, muA, lambdaA] + shearRatioA + maxwellTimeA
    propB = [densityB, muB, lambdaB] + shearRatioB + maxwellTimeB
    self.properties = numpy.array([propA, propB], dtype=numpy.float64)

    # TEMPORARY, need to determine how to use initial state variables
    self.dbStateVars = numpy.zeros(
#.........这里部分代码省略.........

  def __init__(self, name="maxwellisotropic3delastic"):
    """
    Constructor.
    """
    ElasticMaterialApp.__init__(self, name)

    numLocs = 2

    self.dimension = dimension
    self.numLocs = numLocs

    self.dbPropertyValues = ["density", "vs", "vp", "viscosity"]
    self.numPropertyValues = numpy.array([1, 1, 1, 1], dtype=numpy.int32)

    self.dbStateVarValues = ["total-strain-xx",
                             "total-strain-yy",
                             "total-strain-zz",
                             "total-strain-xy",
                             "total-strain-yz",
                             "total-strain-xz",
                             "viscous-strain-xx",
                             "viscous-strain-yy",
                             "viscous-strain-zz",
                             "viscous-strain-xy",
                             "viscous-strain-yz",
                             "viscous-strain-xz",
                             ]
    self.numStateVarValues = numpy.array([6, 6], dtype=numpy.int32)

    densityA = 2500.0
    vsA = 3000.0
    vpA = vsA*3**0.5
    viscosityA = 1.0e18
    strainA = [1.1e-4, 1.2e-4, 1.3e-4, 1.4e-4, 1.5e-4, 1.6e-4]
    initialStressA = [2.1e4, 2.2e4, 2.3e4, 2.4e4, 2.5e4, 2.6e4]
    initialStrainA = [3.1e-5, 3.2e-5, 3.3e-5, 3.4e-5, 3.5e-5, 3.6e-5]
    muA = vsA*vsA*densityA
    lambdaA = vpA*vpA*densityA - 2.0*muA
    maxwellTimeA = viscosityA / muA
    
    densityB = 2000.0
    vsB = 1200.0
    vpB = vsB*3**0.5
    viscosityB = 1.0e18
    strainB = [4.1e-4, 4.2e-4, 4.3e-4, 4.4e-4, 4.5e-4, 4.6e-4]
    initialStressB = [5.1e4, 5.2e4, 5.3e4, 5.4e4, 5.5e4, 5.6e4]
    initialStrainB = [6.1e-5, 6.2e-5, 6.3e-5, 6.4e-5, 6.5e-5, 6.6e-5]
    muB = vsB*vsB*densityB
    lambdaB = vpB*vpB*densityB - 2.0*muB
    maxwellTimeB = viscosityB / muB

    self.lengthScale = 1.0e+3
    self.pressureScale = muA
    self.timeScale = 1.0
    self.densityScale = muA / (self.lengthScale / self.timeScale)**2

    self.dbProperties = numpy.array([ [densityA, vsA, vpA, viscosityA],
                                      [densityB, vsB, vpB, viscosityB] ], 
                                    dtype=numpy.float64)
    self.properties = numpy.array([ [densityA, muA, lambdaA, maxwellTimeA],
                                    [densityB, muB, lambdaB, maxwellTimeB] ],
                                     dtype=numpy.float64)

    # TEMPORARY, need to determine how to use initial state variables
    self.dbStateVars = numpy.zeros( (numLocs, tensorSize+tensorSize),
                                    dtype=numpy.float64)
    self.stateVars = numpy.zeros( (numLocs, tensorSize+tensorSize),
                                  dtype=numpy.float64)

    mu0 = self.pressureScale
    density0 = self.densityScale
    time0 = self.timeScale
    self.propertiesNondim = \
        numpy.array([ [densityA/density0, muA/mu0, lambdaA/mu0, maxwellTimeA/time0],
                      [densityB/density0, muB/mu0, lambdaB/mu0, maxwellTimeB/time0] ],
                    dtype=numpy.float64)

    self.stateVarsNondim = self.stateVars # no scaling

    self.initialStress = numpy.array([initialStressA,
                                      initialStressB],
                                    dtype=numpy.float64)
    self.initialStrain = numpy.array([initialStrainA,
                                      initialStrainB],
                                    dtype=numpy.float64)
    
    self.density = numpy.array([densityA,
                                densityB],
                               dtype=numpy.float64)

    self.strain = numpy.array([strainA,
                               strainB],
                               dtype=numpy.float64)
    
    self.stress = numpy.zeros( (numLocs, tensorSize), dtype=numpy.float64)
    self.elasticConsts = numpy.zeros( (numLocs, numElasticConsts), \
                                      dtype=numpy.float64)
    self.stateVarsUpdated = numpy.zeros( (numLocs, tensorSize + tensorSize), \
                                         dtype=numpy.float64)

#.........这里部分代码省略.........

  def __init__(self, name="genmaxwellbulkisotropic3delastic"):
    """
    Constructor.
    """
    ElasticMaterialApp.__init__(self, name)

    numLocs = 2
    self.dt = 2.0e5
    self.dimension = dimension
    self.numLocs = numLocs

    self.dbPropertyValues = ["density", "vs", "vp",
                             "shear-ratio-1", 
                             "shear-ratio-2", 
                             "shear-ratio-3",
                             "shear-viscosity-1", 
                             "shear-viscosity-2", 
                             "shear-viscosity-3",
                             "bulk-ratio-1", 
                             "bulk-ratio-2", 
                             "bulk-ratio-3",
                             "bulk-viscosity-1", 
                             "bulk-viscosity-2", 
                             "bulk-viscosity-3",
                             ]
    self.numPropertyValues = numpy.array([1, 1, 1,
                                          1, 1, 1,
                                          1, 1, 1,
                                          1, 1, 1,
                                          1, 1, 1], dtype=numpy.int32)

    self.dbStateVarValues = ["total-strain-xx",
                             "total-strain-yy",
                             "total-strain-zz",
                             "total-strain-xy",
                             "total-strain-yz",
                             "total-strain-xz",
                             "viscous-deviatoric-strain-1-xx",
                             "viscous-deviatoric-strain-1-yy",
                             "viscous-deviatoric-strain-1-zz",
                             "viscous-deviatoric-strain-1-xy",
                             "viscous-deviatoric-strain-1-yz",
                             "viscous-deviatoric-strain-1-xz",
                             "viscous-deviatoric-strain-2-xx",
                             "viscous-deviatoric-strain-2-yy",
                             "viscous-deviatoric-strain-2-zz",
                             "viscous-deviatoric-strain-2-xy",
                             "viscous-deviatoric-strain-2-yz",
                             "viscous-deviatoric-strain-2-xz",
                             "viscous-deviatoric-strain-3-xx",
                             "viscous-deviatoric-strain-3-yy",
                             "viscous-deviatoric-strain-3-zz",
                             "viscous-deviatoric-strain-3-xy",
                             "viscous-deviatoric-strain-3-yz",
                             "viscous-deviatoric-strain-3-xz",
                             "viscous-mean-strain-1",
                             "viscous-mean-strain-2",
                             "viscous-mean-strain-3",
                             ]

    self.numStateVarValues = numpy.array([1, 1, 1, 1, 1, 1,
                                          1, 1, 1, 1, 1, 1,
                                          1, 1, 1, 1, 1, 1,
                                          1, 1, 1, 1, 1, 1,
                                          1, 1, 1], dtype=numpy.int32)

    densityA = 2500.0
    vsA = 3000.0
    vpA = vsA*3**0.5
    shearRatioA = [0.5, 0.1, 0.2]
    shearViscosityA = [1.0e+18, 1.0e+17, 1.0e+19]
    bulkRatioA = [0.4, 0.3, 0.1] 
    bulkViscosityA = [2.0e+18, 2.0e+17, 2.0e+19]
    strainA = [1.1e-4, 2.2e-4, 3.3e-4, 4.4e-4, 5.5e-4, 6.6e-4]
    initialStressA = [2.1e4, 2.2e4, 2.3e4, 2.4e4, 2.5e4, 2.6e4]
    initialStrainA = [3.1e-4, 3.2e-4, 3.3e-4, 3.4e-4, 3.5e-4, 3.6e-4]
    muA = vsA*vsA*densityA
    lambdaA = vpA*vpA*densityA - 2.0*muA
    kA = lambdaA + 2/3.0*muA

    densityB = 2000.0
    vsB = 1200.0
    vpB = vsB*3**0.5
    shearRatioB = [0.2, 0.2, 0.2]
    shearViscosityB = [1.0e18, 1.0e19, 1.0e20]
    bulkRatioB = [0.2, 0.2, 0.2] 
    bulkViscosityB = [2.0e18, 2.0e19, 2.0e20]
    strainB = [1.2e-4, 2.3e-4, 3.4e-4, 4.5e-4, 5.6e-4, 6.7e-4]
    initialStressB = [5.1e4, 5.2e4, 5.3e4, 5.4e4, 5.5e4, 5.6e4]
    initialStrainB = [6.1e-4, 6.2e-4, 6.3e-4, 6.4e-4, 6.5e-4, 6.6e-4]
    muB = vsB*vsB*densityB
    lambdaB = vpB*vpB*densityB - 2.0*muB
    kB = lambdaB + 2/3.0*muB

    maxwellTimeA = [0.0, 0.0, 0.0]
    maxwellTimeB = [0.0, 0.0, 0.0]
    for i in xrange(numMaxwellModels):
      if shearRatioA[i] != 0.0:
        maxwellTimeA[i] = shearViscosityA[i]/muA
      if shearRatioB[i] != 0.0:
#.........这里部分代码省略.........

  def __init__(self, name="maxwellplanestrainelastic"):
    """
    Constructor.
    """
    ElasticMaterialApp.__init__(self, name)

    # import pdb
    # pdb.set_trace()

    numLocs = 2

    self.dimension = dimension
    self.numLocs = numLocs

    self.dbPropertyValues = ["density", "vs", "vp", "viscosity"]
    self.numPropertyValues = numpy.array([1, 1, 1, 1], dtype=numpy.int32)

    self.dbStateVarValues = ["stress-zz-initial",
                             "total-strain-xx",
                             "total-strain-yy",
                             "total-strain-xy",
                             "viscous-strain-xx",
                             "viscous-strain-yy",
                             "viscous-strain-zz",
                             "viscous-strain-xy"
                             ]
    self.numStateVarValues = numpy.array([1, 3, 4], dtype=numpy.int32)
    
    densityA = 2500.0
    vsA = 3000.0
    vpA = vsA*3**0.5
    viscosityA = 1.0e18
    strainA = [1.1e-4, 1.2e-4, 1.4e-4]
    initialStressA = [2.1e4, 2.2e4, 2.4e4]
    initialStrainA = [3.1e-5, 3.2e-5, 3.4e-5]
    muA = vsA*vsA*densityA
    lambdaA = vpA*vpA*densityA - 2.0*muA
    maxwellTimeA = viscosityA / muA
    stressInitialZZA = numpy.array([1.5e4], dtype=numpy.float64)
    
    densityB = 2000.0
    vsB = 1200.0
    vpB = vsB*3**0.5
    viscosityB = 1.0e18
    strainB = [4.1e-4, 4.2e-4, 4.4e-4]
    initialStressB = [5.1e4, 5.2e4, 5.4e4]
    initialStrainB = [6.1e-5, 6.2e-5, 6.4e-5]
    muB = vsB*vsB*densityB
    lambdaB = vpB*vpB*densityB - 2.0*muB
    maxwellTimeB = viscosityB / muB
    stressInitialZZB = numpy.array([4.5e4], dtype=numpy.float64)

    self.lengthScale = 1.0e+3
    self.pressureScale = muA
    self.timeScale = 1.0
    self.densityScale = muA / (self.lengthScale / self.timeScale)**2

    self.dbProperties = numpy.array([ [densityA, vsA, vpA, viscosityA],
                                      [densityB, vsB, vpB, viscosityB] ], 
                                    dtype=numpy.float64)
    self.properties = numpy.array([ [densityA, muA, lambdaA, maxwellTimeA],
                                    [densityB, muB, lambdaB, maxwellTimeB] ],
                                     dtype=numpy.float64)

    # TEMPORARY, need to determine how to use initial state variables
    # At present, only the first (stressInitialZZ) is being used.
    self.dbStateVars = numpy.zeros( (numLocs, 1 + tensorSize + 4),
                                    dtype=numpy.float64)
    self.dbStateVars[0, 0] = stressInitialZZA
    self.dbStateVars[1, 0] = stressInitialZZB
    
    self.stateVars = numpy.zeros( (numLocs, 1 + tensorSize + 4),
                                  dtype=numpy.float64)
    self.stateVars[0, 0] = stressInitialZZA
    self.stateVars[1, 0] = stressInitialZZB

    mu0 = self.pressureScale
    density0 = self.densityScale
    time0 = self.timeScale
    self.propertiesNondim = \
        numpy.array([ [densityA/density0, muA/mu0, lambdaA/mu0,
                       maxwellTimeA/time0],
                      [densityB/density0, muB/mu0, lambdaB/mu0,
                       maxwellTimeB/time0] ],
                    dtype=numpy.float64)

    stressInitialZZANondim = stressInitialZZA/mu0
    stressInitialZZBNondim = stressInitialZZB/mu0

    self.stateVarsNondim = numpy.zeros( (numLocs, 1 + tensorSize + 4),
                                        dtype=numpy.float64)

    self.stateVarsNondim[0, 0] = stressInitialZZANondim
    self.stateVarsNondim[1, 0] = stressInitialZZBNondim

    self.initialStress = numpy.array([initialStressA,
                                      initialStressB],
                                    dtype=numpy.float64)
    self.initialStrain = numpy.array([initialStrainA,
                                      initialStrainB],
#.........这里部分代码省略.........

  def __init__(self, name="powerlawplanestrainelastic"):
    """
    Constructor.
    """
    ElasticMaterialApp.__init__(self, name)

    # import pdb
    # pdb.set_trace()
    numLocs = 2

    self.dimension = dimension
    self.numLocs = numLocs

    self.dbPropertyValues = ["density", "vs", "vp",
                             "reference-strain-rate", "reference-stress",
                             "power-law-exponent"]
    self.numPropertyValues = numpy.array([1, 1, 1, 1, 1, 1], dtype=numpy.int32)

    self.dbStateVarValues = ["stress-zz-initial",
                             "viscous-strain-xx",
                             "viscous-strain-yy",
                             "viscous-strain-zz",
                             "viscous-strain-xy",
                             "stress4-xx",
                             "stress4-yy",
                             "stress4-zz",
                             "stress4-xy"
                             ]
    self.numStateVarValues = numpy.array([1, 4, 4], dtype=numpy.int32)

    densityA = 2500.0
    vsA = 3000.0
    vpA = vsA*3**0.5
    # Derive new values in based on previous value for power-law coefficient
    # and viscosity coefficient.
    powerLawCoeffA = 1.0/3.0e18
    refStrainRateA = 1.0e-6
    powerLawExponentA = 1.0
    strainA = [1.1e-4, 1.2e-4, 1.4e-4]
    initialStressA = [2.1e4, 2.2e4, 2.4e4]
    initialStrainA = [3.1e-4, 3.2e-4, 3.4e-4]
    muA = vsA*vsA*densityA
    lambdaA = vpA*vpA*densityA - 2.0*muA
    stressZZInitialA = numpy.array([1.0754e4], dtype=numpy.float64)

    viscosityCoeffA = (1.0/((3.0**0.5)**(powerLawExponentA + 1.0) \
                            * powerLawCoeffA))**(1.0/powerLawExponentA)
    refStressA = viscosityCoeffA * \
                 (2.0 * refStrainRateA) ** (1.0/powerLawExponentA)
    # refStressA = (refStrainRateA/powerLawCoeffA)**(1.0/powerLawExponentA)
    
    densityB = 2000.0
    vsB = 1200.0
    vpB = vsB*3**0.5
    powerLawCoeffB = 1.0/9.0e30
    refStrainRateB = 1.0e-6
    powerLawExponentB = 3.0
    strainB = [4.1e-4, 4.2e-4, 4.4e-4]
    initialStressB = [5.1e4, 5.2e4, 5.4e4]
    initialStrainB = [6.1e-4, 6.2e-4, 6.4e-4]
    muB = vsB*vsB*densityB
    lambdaB = vpB*vpB*densityB - 2.0*muB
    viscosityCoeffB = (1.0/((3.0**0.5)**(powerLawExponentB + 1.0) \
                            * powerLawCoeffB))**(1.0/powerLawExponentB)
    refStressB = viscosityCoeffB * \
                 (2.0 * refStrainRateB) ** (1.0/powerLawExponentB)
    # refStressB = (refStrainRateB/powerLawCoeffB)**(1.0/powerLawExponentB)
    stressZZInitialB = numpy.array([2.575e4], dtype=numpy.float64)

    self.lengthScale = 1.0e+3
    self.pressureScale = muA
    self.timeScale = 1.0
    self.densityScale = muA / (self.lengthScale / self.timeScale)**2
    self.strainRateScale = 1.0/self.timeScale

    self.dbProperties = numpy.array([ [densityA, vsA, vpA, \
                                       refStrainRateA, refStressA, \
                                       powerLawExponentA],
                                      [densityB, vsB, vpB, \
                                       refStrainRateB, refStressB, \
                                       powerLawExponentB] ], 
                                    dtype=numpy.float64)
    self.properties = numpy.array([ [densityA, muA, lambdaA, \
                                     refStrainRateA, refStressA, \
                                     powerLawExponentA],
                                    [densityB, muB, lambdaB, \
                                     refStrainRateB, refStressB, \
                                     powerLawExponentB] ],
                                     dtype=numpy.float64)

    # TEMPORARY, need to determine how to use initial state variables
    self.dbStateVars = numpy.zeros( (numLocs, 1 + 4 + 4),
                                    dtype=numpy.float64)
    self.dbStateVars[0,0] = stressZZInitialA
    self.dbStateVars[1,0] = stressZZInitialB

    self.stateVars = numpy.zeros( (numLocs, 1 + 4 + 4), dtype=numpy.float64)
    self.stateVars[0,0] = stressZZInitialA
    self.stateVars[1,0] = stressZZInitialB

#.........这里部分代码省略.........