- cohesionA SolidMechanicsHardening UserObject that defines hardening of the cohesion. Physically the cohesion should not be negative.
C++ Type:UserObjectName
Controllable:No
Description:A SolidMechanicsHardening UserObject that defines hardening of the cohesion. Physically the cohesion should not be negative.
- dilation_angleA SolidMechanicsHardening UserObject that defines hardening of the dilation angle (in radians). Usually the dilation angle is not greater than the friction angle, and it is between 0 and 90deg.
C++ Type:UserObjectName
Controllable:No
Description:A SolidMechanicsHardening UserObject that defines hardening of the dilation angle (in radians). Usually the dilation angle is not greater than the friction angle, and it is between 0 and 90deg.
- friction_angleA SolidMechanicsHardening UserObject that defines hardening of the friction angle (in radians). Physically the friction angle should be between 0 and 90deg.
C++ Type:UserObjectName
Controllable:No
Description:A SolidMechanicsHardening UserObject that defines hardening of the friction angle (in radians). Physically the friction angle should be between 0 and 90deg.
- internal_constraint_toleranceThe Newton-Raphson process is only deemed converged if the internal constraint is less than this.
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:The Newton-Raphson process is only deemed converged if the internal constraint is less than this.
- mc_tip_smootherSmoothing parameter: the cone vertex at mean = cohesion*cot(friction_angle), will be smoothed by the given amount. Typical value is 0.1*cohesion
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Smoothing parameter: the cone vertex at mean = cohesion*cot(friction_angle), will be smoothed by the given amount. Typical value is 0.1*cohesion
- yield_function_toleranceIf the yield function is less than this amount, the (stress, internal parameter) are deemed admissible.
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:If the yield function is less than this amount, the (stress, internal parameter) are deemed admissible.
SolidMechanicsPlasticMohrCoulomb
buildconstruction:Undocumented Class
The SolidMechanicsPlasticMohrCoulomb has not been documented. The content listed below should be used as a starting point for documenting the class, which includes the typical automatic documentation associated with a MooseObject; however, what is contained is ultimately determined by what is necessary to make the documentation clear for users.
Non-associative Mohr-Coulomb plasticity with hardening/softening
Overview
Example Input File Syntax
Input Parameters
- cap_rate0For the 'cap' tip_scheme, this controls how quickly the cap degenerates to a hemisphere: small values mean a slow degeneration to a hemisphere (and zero means the 'cap' will be totally inactive). Typical value is 1/tensile_strength
Default:0
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:For the 'cap' tip_scheme, this controls how quickly the cap degenerates to a hemisphere: small values mean a slow degeneration to a hemisphere (and zero means the 'cap' will be totally inactive). Typical value is 1/tensile_strength
- cap_start0For the 'cap' tip_scheme, smoothing is performed in the stress_mean > cap_start region
Default:0
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:For the 'cap' tip_scheme, smoothing is performed in the stress_mean > cap_start region
- mc_edge_smoother25Smoothing parameter: the edges of the cone are smoothed by the given amount.
Default:25
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Smoothing parameter: the edges of the cone are smoothed by the given amount.
- mc_lode_cutoffIf the second invariant of stress is less than this amount, the Lode angle is assumed to be zero. This is to gaurd against precision-loss problems, and this parameter should be set small. Default = 0.00001*((yield_Function_tolerance)^2)
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:If the second invariant of stress is less than this amount, the Lode angle is assumed to be zero. This is to gaurd against precision-loss problems, and this parameter should be set small. Default = 0.00001*((yield_Function_tolerance)^2)
- tip_schemehyperbolicScheme by which the pyramid's tip will be smoothed.
Default:hyperbolic
C++ Type:MooseEnum
Options:hyperbolic, cap
Controllable:No
Description:Scheme by which the pyramid's tip will be smoothed.
Optional Parameters
- allow_duplicate_execution_on_initialFalseIn the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).
Default:False
C++ Type:bool
Controllable:No
Description:In the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).
- execute_onTIMESTEP_ENDThe list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html.
Default:TIMESTEP_END
C++ Type:ExecFlagEnum
Options:XFEM_MARK, FORWARD, ADJOINT, HOMOGENEOUS_FORWARD, ADJOINT_TIMESTEP_BEGIN, ADJOINT_TIMESTEP_END, NONE, INITIAL, LINEAR, LINEAR_CONVERGENCE, NONLINEAR, NONLINEAR_CONVERGENCE, POSTCHECK, TIMESTEP_END, TIMESTEP_BEGIN, MULTIAPP_FIXED_POINT_END, MULTIAPP_FIXED_POINT_BEGIN, MULTIAPP_FIXED_POINT_CONVERGENCE, FINAL, CUSTOM, TRANSFER
Controllable:No
Description:The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html.
- execution_order_group0Execution order groups are executed in increasing order (e.g., the lowest number is executed first). Note that negative group numbers may be used to execute groups before the default (0) group. Please refer to the user object documentation for ordering of user object execution within a group.
Default:0
C++ Type:int
Controllable:No
Description:Execution order groups are executed in increasing order (e.g., the lowest number is executed first). Note that negative group numbers may be used to execute groups before the default (0) group. Please refer to the user object documentation for ordering of user object execution within a group.
- force_postauxFalseForces the UserObject to be executed in POSTAUX
Default:False
C++ Type:bool
Controllable:No
Description:Forces the UserObject to be executed in POSTAUX
- force_preauxFalseForces the UserObject to be executed in PREAUX
Default:False
C++ Type:bool
Controllable:No
Description:Forces the UserObject to be executed in PREAUX
- force_preicFalseForces the UserObject to be executed in PREIC during initial setup
Default:False
C++ Type:bool
Controllable:No
Description:Forces the UserObject to be executed in PREIC during initial setup
Execution Scheduling Parameters
- control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
- enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:Yes
Description:Set the enabled status of the MooseObject.
- use_displaced_meshFalseWhether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Default:False
C++ Type:bool
Controllable:No
Description:Whether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Advanced Parameters
- prop_getter_suffixAn optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:An optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
- use_interpolated_stateFalseFor the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Default:False
C++ Type:bool
Controllable:No
Description:For the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Material Property Retrieval Parameters
Input Files
- (modules/solid_mechanics/test/tests/mohr_coulomb/except4.i)
- (modules/solid_mechanics/test/tests/multi/rock1.i)
- (modules/solid_mechanics/test/tests/initial_stress/mc_tensile.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/small_deform4.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/many_deforms_cap.i)
- (modules/solid_mechanics/test/tests/multi/paper3.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/small_deform1_uo.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/uni_axial3.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/uni_axial1.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/small_deform_hard_cubic.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/small_deform2.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/uni_axial2.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/except3.i)
- (modules/solid_mechanics/test/tests/jacobian/cto12.i)
- (modules/solid_mechanics/test/tests/jacobian/cto13.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/small_deform2_small_strain.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/small_deform3.i)
- (modules/solid_mechanics/test/tests/notched_plastic_block/biaxial_abbo.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/small_deform5.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/small_deform_hard3.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/random.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/except2.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/except1.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/small_deform_hard2.i)
- (modules/solid_mechanics/test/tests/multi/mc_wpt_1.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/small_deform_hard1.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/uni_axial1_small_strain.i)
- (modules/solid_mechanics/test/tests/multi/special_rock1.i)
- (modules/solid_mechanics/test/tests/mohr_coulomb/small_deform6.i)
(modules/solid_mechanics/test/tests/mohr_coulomb/except4.i)
# checking for exception error messages on the edge smoothing
# here edge_smoother=5deg, which means the friction_angle must be <= 35.747
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '1E-6*x'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '1E-6*y'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '1E-6*z'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10
[../]
[./mc_phi]
type = SolidMechanicsHardeningExponential
value_0 = 0.52359878 # 30deg
value_residual = 0.62831853 # 36deg
rate = 3000.0
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 1
mc_edge_smoother = 5
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
debug_jac_at_stress = '10 0 0 0 10 0 0 0 10'
debug_jac_at_pm = 1
debug_jac_at_intnl = 1
debug_stress_change = 1E-5
debug_pm_change = 1E-6
debug_intnl_change = 1E-6
[../]
[]
[Executioner]
end_time = 1
dt = 1
type = Transient
[]
[Outputs]
file_base = except4
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/multi/rock1.i)
# Plasticity models:
# Mohr-Coulomb with cohesion = 40MPa, friction angle = 35deg, dilation angle = 10deg
# Tensile with strength = 1MPa
# WeakPlaneShear with cohesion = 1MPa, friction angle = 25deg, dilation angle = 25deg
# WeakPlaneTensile with strength = 0.01MPa
#
# Lame lambda = 1GPa. Lame mu = 1.3GPa
#
# A line of elements is perturbed randomly, and return to the yield surface at each quadpoint is checked
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1000
ny = 1234
nz = 1
xmin = 0
xmax = 1000
ymin = 0
ymax = 1234
zmin = 0
zmax = 1
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[ICs]
[./x]
type = RandomIC
min = -0.1
max = 0.1
variable = disp_x
[../]
[./y]
type = RandomIC
min = -0.1
max = 0.1
variable = disp_y
[../]
[./z]
type = RandomIC
min = -0.1
max = 0.1
variable = disp_z
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '0'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '0'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '0'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./f0]
order = CONSTANT
family = MONOMIAL
[../]
[./f1]
order = CONSTANT
family = MONOMIAL
[../]
[./f2]
order = CONSTANT
family = MONOMIAL
[../]
[./f3]
order = CONSTANT
family = MONOMIAL
[../]
[./int0]
order = CONSTANT
family = MONOMIAL
[../]
[./int1]
order = CONSTANT
family = MONOMIAL
[../]
[./int2]
order = CONSTANT
family = MONOMIAL
[../]
[./int3]
order = CONSTANT
family = MONOMIAL
[../]
[./iter]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./f0]
type = MaterialStdVectorAux
property = plastic_yield_function
index = 0
variable = f0
[../]
[./f1]
type = MaterialStdVectorAux
property = plastic_yield_function
index = 1
variable = f1
[../]
[./f2]
type = MaterialStdVectorAux
property = plastic_yield_function
index = 2
variable = f2
[../]
[./f3]
type = MaterialStdVectorAux
property = plastic_yield_function
index = 3
variable = f3
[../]
[./int0]
type = MaterialStdVectorAux
property = plastic_internal_parameter
factor = 1E6
index = 0
variable = int0
[../]
[./int1]
type = MaterialStdVectorAux
property = plastic_internal_parameter
factor = 1E6
index = 1
variable = int1
[../]
[./int2]
type = MaterialStdVectorAux
property = plastic_internal_parameter
factor = 1E6
index = 2
variable = int2
[../]
[./int3]
type = MaterialStdVectorAux
property = plastic_internal_parameter
factor = 1E6
index = 3
variable = int3
[../]
[./iter]
type = MaterialRealAux
property = plastic_NR_iterations
variable = iter
[../]
[]
[Postprocessors]
[./raw_f0]
type = ElementExtremeValue
variable = f0
outputs = console
[../]
[./raw_f1]
type = ElementExtremeValue
variable = f1
outputs = console
[../]
[./raw_f2]
type = ElementExtremeValue
variable = f2
outputs = console
[../]
[./raw_f3]
type = ElementExtremeValue
variable = f3
outputs = console
[../]
[./iter]
type = ElementExtremeValue
variable = iter
outputs = console
[../]
[./f0]
type = FunctionValuePostprocessor
function = should_be_zero0_fcn
[../]
[./f1]
type = FunctionValuePostprocessor
function = should_be_zero1_fcn
[../]
[./f2]
type = FunctionValuePostprocessor
function = should_be_zero2_fcn
[../]
[./f3]
type = FunctionValuePostprocessor
function = should_be_zero3_fcn
[../]
[]
[Functions]
[./should_be_zero0_fcn]
type = ParsedFunction
expression = 'if(a<1E-1,0,a)'
symbol_names = 'a'
symbol_values = 'raw_f0'
[../]
[./should_be_zero1_fcn]
type = ParsedFunction
expression = 'if(a<1E-1,0,a)'
symbol_names = 'a'
symbol_values = 'raw_f1'
[../]
[./should_be_zero2_fcn]
type = ParsedFunction
expression = 'if(a<1E-1,0,a)'
symbol_names = 'a'
symbol_values = 'raw_f2'
[../]
[./should_be_zero3_fcn]
type = ParsedFunction
expression = 'if(a<1E-1,0,a)'
symbol_names = 'a'
symbol_values = 'raw_f3'
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 4E7
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 35
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 10
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4E6
yield_function_tolerance = 1.0E-1
internal_constraint_tolerance = 1.0E-7
[../]
[./ts]
type = SolidMechanicsHardeningConstant
value = 1E6
[../]
[./tensile]
type = SolidMechanicsPlasticTensile
tensile_strength = ts
tensile_tip_smoother = 1E5
yield_function_tolerance = 1.0E-1
internal_constraint_tolerance = 1.0E-7
[../]
[./coh]
type = SolidMechanicsHardeningConstant
value = 1E6
[../]
[./tanphi]
type = SolidMechanicsHardeningConstant
value = 0.46630766
[../]
[./tanpsi]
type = SolidMechanicsHardeningConstant
value = 0.46630766
[../]
[./wps]
type = SolidMechanicsPlasticWeakPlaneShear
cohesion = coh
tan_friction_angle = tanphi
tan_dilation_angle = tanpsi
smoother = 1E5
yield_function_tolerance = 1.0E-1
internal_constraint_tolerance = 1.0E-7
[../]
[./str]
type = SolidMechanicsHardeningConstant
value = 0.01E6
[../]
[./wpt]
type = SolidMechanicsPlasticWeakPlaneTensile
tensile_strength = str
yield_function_tolerance = 1.0E-1
internal_constraint_tolerance = 1.0E-7
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '1E9 1.3E9'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./multi]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-7
plastic_models = 'mc tensile wps wpt'
deactivation_scheme = 'optimized_to_safe_to_dumb'
max_NR_iterations = 20
min_stepsize = 1E-4
max_stepsize_for_dumb = 1E-3
debug_fspb = crash
debug_jac_at_stress = '10 0 0 0 10 0 0 0 10'
debug_jac_at_pm = '1 1 1 1'
debug_jac_at_intnl = '1 1 1 1'
debug_stress_change = 1E1
debug_pm_change = '1E-6 1E-6 1E-6 1E-6'
debug_intnl_change = '1E-6 1E-6 1E-6 1E-6'
[../]
[]
[Executioner]
end_time = 1
dt = 1
type = Transient
[]
[Outputs]
file_base = rock1
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/initial_stress/mc_tensile.i)
# In this example, an initial stress is applied that
# is inadmissible, and the return-map algorithm must be
# used to return to the yield surface before any other
# computations can be carried out.
# In this case, the return-map algorithm must subdivide
# the initial stress, otherwise it does not converge.
# This test is testing that subdivision process.
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./bottomx]
type = DirichletBC
variable = disp_x
boundary = 'back'
value = 0.0
[../]
[./bottomy]
type = DirichletBC
variable = disp_y
boundary = 'back'
value = 0.0
[../]
[./bottomz]
type = DirichletBC
variable = disp_z
boundary = 'back'
value = 0.0
[../]
[./topx]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front'
function = '2*t-1'
[../]
[./topy]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front'
function = 't-1'
[../]
[./topz]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front'
function = 't-1'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./iter]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./iter]
type = MaterialRealAux
property = plastic_NR_iterations
variable = iter
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[./iter]
type = PointValue
point = '0 0 0'
variable = iter
outputs = console
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 1E5
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 60
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4.0
mc_edge_smoother = 25
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-9
[../]
[./str]
type = SolidMechanicsHardeningConstant
value = 1
[../]
[./pt]
type = SolidMechanicsPlasticTensile
tensile_strength = str
yield_function_tolerance = 1E-3
tensile_tip_smoother = 0.05
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
displacements = 'disp_x disp_y disp_z'
eigenstrain_names = ini_stress
[../]
[./ini_stress]
type = ComputeEigenstrainFromInitialStress
initial_stress = '8E6 4E6 -18E6 4E6 -40E6 -2E6 -18E6 -2E6 -34E6'
eigenstrain_name = ini_stress
[../]
[./mc]
type = ComputeMultiPlasticityStress
ep_plastic_tolerance = 1E-9
plastic_models = 'pt mc'
deactivation_scheme = safe
max_NR_iterations = 100
min_stepsize = 0.1
[../]
[]
[Executioner]
end_time = 2
dt = 1
type = Transient
[]
[Outputs]
file_base = mc_tensile
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/small_deform4.i)
# apply repeated stretches in z direction, and smaller stretches in the x and y directions
# so that sigma_II = sigma_III,
# which means that lode angle = -30deg.
# The allows yield surface in meridional plane to be mapped out
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '0.25E-6*x*sin(t)'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '0.25E-6*y*sin(t)'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '1E-6*z*t'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./mc_int]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./mc_int_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_internal_parameter
variable = mc_int
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[./internal]
type = PointValue
point = '0 0 0'
variable = mc_int
[../]
[./f]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 50
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningExponential
value_0 = 0
value_residual = 0.8726646 # 50deg
rate = 3000.0
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4
mc_edge_smoother = 20
yield_function_tolerance = 1E-8
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
[../]
[]
[Executioner]
end_time = 30
dt = 1
type = Transient
[]
[Outputs]
file_base = small_deform4
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/many_deforms_cap.i)
# apply many large deformations, checking that the algorithm returns correctly to
# the yield surface each time
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./bottomx]
type = DirichletBC
variable = disp_x
boundary = back
value = 0.0
[../]
[./bottomy]
type = DirichletBC
variable = disp_y
boundary = back
value = 0.0
[../]
[./bottomz]
type = DirichletBC
variable = disp_z
boundary = back
value = 0.0
[../]
[./topx]
type = FunctionDirichletBC
variable = disp_x
boundary = front
function = '(sin(0.05*t)+x)/1E0'
[../]
[./topy]
type = FunctionDirichletBC
variable = disp_y
boundary = front
function = '(cos(0.04*t)+x*y)/1E0'
[../]
[./topz]
type = FunctionDirichletBC
variable = disp_z
boundary = front
function = 't/1E2'
[../]
[]
[AuxVariables]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[]
[Postprocessors]
[./yield_fcn_at_zero]
type = PointValue
point = '0 0 0'
variable = yield_fcn
outputs = 'console'
[../]
[./should_be_zero]
type = FunctionValuePostprocessor
function = should_be_zero_fcn
[../]
[]
[Functions]
[./should_be_zero_fcn]
type = ParsedFunction
expression = 'if(a<1E-3,0,a)'
symbol_names = 'a'
symbol_values = 'yield_fcn_at_zero'
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 1E3
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 30
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
tip_scheme = cap
mc_tip_smoother = 0.0
cap_start = 1000
cap_rate = 1E-3
mc_edge_smoother = 10
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-6
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
max_NR_iterations = 1000
ep_plastic_tolerance = 1E-6
plastic_models = mc
debug_fspb = crash
deactivation_scheme = safe
[../]
[]
[Executioner]
end_time = 1000
dt = 1
type = Transient
[]
[Outputs]
file_base = many_deforms_cap
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/multi/paper3.i)
# This runs the third example models described in the 'MultiSurface' plasticity paper
# Just change the deactivation_scheme
#
# Plasticity models:
# Mohr-Coulomb with cohesion = 40MPa, friction angle = 35deg, dilation angle = 5deg
# Tensile with strength = 1MPa
# WeakPlaneTensile with strength = 1000Pa
# WeakPlaneShear with cohesion = 0.1MPa and friction angle = 25, dilation angle = 5deg
#
# Lame lambda = 1.2GPa. Lame mu = 1.2GPa (Young = 3GPa, poisson = 0.5)
#
# A line of elements is perturbed randomly, and return to the yield surface at each quadpoint is checked
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1000
ny = 125
nz = 1
xmin = 0
xmax = 1000
ymin = 0
ymax = 125
zmin = 0
zmax = 1
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[GlobalParams]
volumetric_locking_correction=true
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[ICs]
[./x]
type = RandomIC
min = -0.1
max = 0.1
variable = disp_x
[../]
[./y]
type = RandomIC
min = -0.1
max = 0.1
variable = disp_y
[../]
[./z]
type = RandomIC
min = -0.1
max = 0.1
variable = disp_z
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '0'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '0'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '0'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./linesearch]
order = CONSTANT
family = MONOMIAL
[../]
[./ld]
order = CONSTANT
family = MONOMIAL
[../]
[./constr_added]
order = CONSTANT
family = MONOMIAL
[../]
[./iter]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./linesearch]
type = MaterialRealAux
property = plastic_linesearch_needed
variable = linesearch
[../]
[./ld]
type = MaterialRealAux
property = plastic_linear_dependence_encountered
variable = ld
[../]
[./constr_added]
type = MaterialRealAux
property = plastic_constraints_added
variable = constr_added
[../]
[./iter]
type = MaterialRealAux
property = plastic_NR_iterations
variable = iter
[../]
[]
[Postprocessors]
[./max_iter]
type = ElementExtremeValue
variable = iter
outputs = console
[../]
[./av_iter]
type = ElementAverageValue
variable = iter
outputs = 'console csv'
[../]
[./av_linesearch]
type = ElementAverageValue
variable = linesearch
outputs = 'console csv'
[../]
[./av_ld]
type = ElementAverageValue
variable = ld
outputs = 'console csv'
[../]
[./av_constr_added]
type = ElementAverageValue
variable = constr_added
outputs = 'console csv'
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 4E7
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 35
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulombMulti
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
yield_function_tolerance = 1.0
shift = 1.0
internal_constraint_tolerance = 1.0E-7
[../]
[./mc_smooth]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4E6
yield_function_tolerance = 1.0
internal_constraint_tolerance = 1.0E-7
[../]
[./ts]
type = SolidMechanicsHardeningConstant
value = 1E6
[../]
[./tensile]
type = SolidMechanicsPlasticTensileMulti
tensile_strength = ts
yield_function_tolerance = 1.0
shift = 1.0
internal_constraint_tolerance = 1.0E-7
use_custom_returnMap = false
use_custom_cto = false
[../]
[./tensile_smooth]
type = SolidMechanicsPlasticTensile
tensile_strength = ts
tensile_tip_smoother = 1E5
yield_function_tolerance = 1.0
internal_constraint_tolerance = 1.0E-7
[../]
[./wpt_str]
type = SolidMechanicsHardeningConstant
value = 1.0E3
[../]
[./wpt]
type = SolidMechanicsPlasticWeakPlaneTensile
tensile_strength = wpt_str
yield_function_tolerance = 1.0
internal_constraint_tolerance = 1.0E-7
[../]
[./wps_c]
type = SolidMechanicsHardeningConstant
value = 1.0E5
[../]
[./wps_tan_phi]
type = SolidMechanicsHardeningConstant
value = 0.466
[../]
[./wps_tan_psi]
type = SolidMechanicsHardeningConstant
value = 0.087
[../]
[./wps]
type = SolidMechanicsPlasticWeakPlaneShear
cohesion = wps_c
tan_friction_angle = wps_tan_phi
tan_dilation_angle = wps_tan_psi
smoother = 1.0E4
yield_function_tolerance = 1.0
internal_constraint_tolerance = 1.0E-7
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '1.2E9 1.2E9'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./multi]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-7
plastic_models = 'tensile_smooth mc_smooth wpt wps'
max_NR_iterations = 30
specialIC = 'none'
deactivation_scheme = 'optimized'
min_stepsize = 1E-6
max_stepsize_for_dumb = 1E-2
debug_fspb = crash
debug_jac_at_stress = '10 0 0 0 10 0 0 0 10'
debug_jac_at_pm = '1 1 1 1'
debug_jac_at_intnl = '1 1 1 1'
debug_stress_change = 1E1
debug_pm_change = '1E-6 1E-6 1E-6 1E-6'
debug_intnl_change = '1E-6 1E-6 1E-6 1E-6'
[../]
[]
[Executioner]
end_time = 1
dt = 1
type = Transient
[]
[Outputs]
file_base = paper3
exodus = false
csv = true
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/small_deform1_uo.i)
# apply uniform stretch in x, y and z directions.
# With cohesion = 10, friction_angle = 60deg, tip_smoother = 4, the
# algorithm should return to
# sigma_m = (10*Cos(60) - 4)/Sin(60) = 1.1547
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '1E-6*x'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '1E-6*y'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '1E-6*z'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./f]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./f]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = f
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[./f]
type = PointValue
point = '0 0 0'
variable = f
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 60
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4
mc_edge_smoother = 25
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
debug_jac_at_stress = '10 0 0 0 10 0 0 0 10'
debug_jac_at_pm = 1
debug_jac_at_intnl = 1
debug_stress_change = 1E-5
debug_pm_change = 1E-6
debug_intnl_change = 1E-6
[../]
[]
[Executioner]
end_time = 1
dt = 1
type = Transient
[]
[Outputs]
file_base = small_deform1_uo
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/uni_axial3.i)
[Mesh]
type = FileMesh
file = quarter_hole.e
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./zmin_zzero]
type = DirichletBC
variable = disp_z
boundary = 'zmin'
value = '0'
[../]
[./xmin_xzero]
type = DirichletBC
variable = disp_x
boundary = 'xmin'
value = '0'
[../]
[./ymin_yzero]
type = DirichletBC
variable = disp_y
boundary = 'ymin'
value = '0'
[../]
[./ymax_disp]
type = FunctionDirichletBC
variable = disp_y
boundary = 'ymax'
function = '-1E-4*t'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./mc_int]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./mc_int_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_internal_parameter
variable = mc_int
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0.005 0.02 0.002'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0.005 0.02 0.002'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0.005 0.02 0.002'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0.005 0.02 0.002'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0.005 0.02 0.002'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0.005 0.02 0.002'
variable = stress_zz
[../]
[./f]
type = PointValue
point = '0.005 0.02 0.002'
variable = yield_fcn
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10E6
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 40
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 40
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 0.01E6
mc_edge_smoother = 29
yield_function_tolerance = 1E-5
internal_constraint_tolerance = 1E-11
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 1
fill_method = symmetric_isotropic
C_ijkl = '0 5E9' # young = 10Gpa, poisson = 0.0
[../]
[./strain]
type = ComputeIncrementalStrain
block = 1
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 1
ep_plastic_tolerance = 1E-11
plastic_models = mc
max_NR_iterations = 1000
debug_fspb = crash
[../]
[]
# Preconditioning and Executioner options kindly provided by Andrea
[Preconditioning]
[./andy]
type = SMP
full = true
[../]
[]
[Executioner]
end_time = 1.05
dt = 0.1
solve_type = NEWTON
type = Transient
nl_abs_tol = 1E-10
nl_rel_tol = 1E-12
l_tol = 1E-2
l_max_its = 50
nl_max_its = 400
petsc_options_iname = '-pc_type -pc_asm_overlap -sub_pc_type -ksp_type -ksp_gmres_restart'
petsc_options_value = ' asm 2 lu gmres 200'
[]
[Outputs]
file_base = uni_axial3
exodus = true
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/uni_axial1.i)
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
# back = zmin
# front = zmax
# bottom = ymin
# top = ymax
# left = xmin
# right = xmax
[./xmin_xzero]
type = DirichletBC
variable = disp_x
boundary = 'left'
value = '0'
[../]
[./ymin_yzero]
type = DirichletBC
variable = disp_y
boundary = 'bottom'
value = '0'
[../]
[./zmin_zzero]
type = DirichletBC
variable = disp_z
boundary = 'back'
value = '0'
[../]
[./zmax_disp]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front'
function = '-1E-3*t'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./mc_int]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./mc_int_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_internal_parameter
variable = mc_int
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[./mc_int]
type = PointValue
point = '0 0 0'
variable = mc_int
[../]
[./f]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10E6
[../]
[./mc_phi]
type = SolidMechanicsHardeningExponential
value_0 = 0
value_residual = 0.6981317 # 40deg
rate = 10000
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 0
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 0
mc_edge_smoother = 25
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-10
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '5.77E10 3.85E10' # young = 100Gpa, poisson = 0.3
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-10
plastic_models = mc
max_NR_iterations = 1000
debug_fspb = crash
[../]
[]
[Preconditioning]
[./andy]
type = SMP
full = true
[../]
[]
[Executioner]
end_time = 0.5
dt = 0.05
solve_type = PJFNK # cannot use NEWTON because we are using ComputeFiniteStrain, and hence the Jacobian contributions will not be correct, even though ComputeMultiPlasticityStress will compute the correct consistent tangent operator for small strains
type = Transient
line_search = 'none'
nl_rel_tol = 1E-10
l_tol = 1E-3
l_max_its = 200
nl_max_its = 10
petsc_options_iname = '-pc_type -pc_asm_overlap -sub_pc_type -ksp_type -ksp_gmres_restart'
petsc_options_value = ' asm 2 lu gmres 200'
[]
[Outputs]
file_base = uni_axial1
exodus = true
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/small_deform_hard_cubic.i)
# apply uniform stretches in x, y and z directions.
# let cohesion = 10, cohesion_residual = 2, cohesion_limit = 0.0003
# With cohesion = C, friction_angle = 60deg, tip_smoother = 4, the
# algorithm should return to
# sigma_m = (C*Cos(60) - 4)/Sin(60)
# This allows checking of the relationship for C
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '1E-6*x*t'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '1E-6*y*t'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '1E-6*z*t'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./mc_int]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./mc_int_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_internal_parameter
variable = mc_int
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[./internal]
type = PointValue
point = '0 0 0'
variable = mc_int
[../]
[./f]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningCubic
value_0 = 10
value_residual = 2
internal_limit = 0.0003
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 60
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4
mc_edge_smoother = 25
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
debug_jac_at_stress = '10 1 2 1 10 3 2 3 10'
debug_jac_at_pm = 1
debug_jac_at_intnl = 1E-4
debug_stress_change = 1E-5
debug_pm_change = 1E-6
debug_intnl_change = 1E-8
[../]
[]
[Executioner]
end_time = 10
dt = 0.25
type = Transient
[]
[Outputs]
file_base = small_deform_hard_cubic
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/small_deform2.i)
# apply repeated stretches in x, y and z directions, so that mean_stress = 0
# This maps out the yield surface in the octahedral plane for zero mean stress
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '2E-6*x*sin(t)'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '2E-6*y*sin(2*t)'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '-2E-6*z*(sin(t)+sin(2*t))'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[./f]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[./f0]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[./f1]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 20
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 0
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4
mc_edge_smoother = 20
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
[../]
[]
[Executioner]
end_time = 100
dt = 1
type = Transient
[]
[Outputs]
file_base = small_deform2
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/uni_axial2.i)
[Mesh]
type = FileMesh
file = quarter_hole.e
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./zmin_zzero]
type = DirichletBC
variable = disp_z
boundary = 'zmin'
value = '0'
[../]
[./xmin_xzero]
type = DirichletBC
variable = disp_x
boundary = 'xmin'
value = '0'
[../]
[./ymin_yzero]
type = DirichletBC
variable = disp_y
boundary = 'ymin'
value = '0'
[../]
[./ymax_disp]
type = FunctionDirichletBC
variable = disp_y
boundary = 'ymax'
function = '-1E-4*t'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./mc_int]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./mc_int_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_internal_parameter
variable = mc_int
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0.005 0.02 0.002'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0.005 0.02 0.002'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0.005 0.02 0.002'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0.005 0.02 0.002'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0.005 0.02 0.002'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0.005 0.02 0.002'
variable = stress_zz
[../]
[./f]
type = PointValue
point = '0.005 0.02 0.002'
variable = yield_fcn
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10E6
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 2
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 2
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 0.01E6
mc_edge_smoother = 29
yield_function_tolerance = 1E-5
internal_constraint_tolerance = 1E-11
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 1
fill_method = symmetric_isotropic
C_ijkl = '0 5E9' # young = 10Gpa, poisson = 0.0
[../]
[./strain]
type = ComputeIncrementalStrain
block = 1
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 1
ep_plastic_tolerance = 1E-11
plastic_models = mc
max_NR_iterations = 1000
debug_fspb = crash
[../]
[]
# Preconditioning and Executioner options kindly provided by Andrea
[Preconditioning]
[./andy]
type = SMP
full = true
[../]
[]
[Executioner]
end_time = 0.5
dt = 0.1
solve_type = NEWTON
type = Transient
l_tol = 1E-2
nl_abs_tol = 1E-9
nl_rel_tol = 1E-11
l_max_its = 200
nl_max_its = 400
petsc_options_iname = '-pc_type -pc_asm_overlap -sub_pc_type -ksp_type -ksp_gmres_restart'
petsc_options_value = ' asm 2 lu gmres 200'
[]
[Outputs]
file_base = uni_axial2
exodus = true
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/except3.i)
# checking for exception error messages on the edge smoothing
# here edge_smoother=5deg, which means the friction_angle must be <= 35.747
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '1E-6*x'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '1E-6*y'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '1E-6*z'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 36
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 1
mc_edge_smoother = 5
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
debug_jac_at_stress = '10 0 0 0 10 0 0 0 10'
debug_jac_at_pm = 1
debug_jac_at_intnl = 1
debug_stress_change = 1E-5
debug_pm_change = 1E-6
debug_intnl_change = 1E-6
[../]
[]
[Executioner]
end_time = 1
dt = 1
type = Transient
[]
[Outputs]
file_base = except3
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/jacobian/cto12.i)
# checking jacobian for nonlinear plasticity (single surface, smoothed MohrCoulomb)
# note: must have min_stepsize=1 otherwise the nonlinearities compound and make the jacobian more inaccurate
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[GlobalParams]
block = 0
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 60
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4
mc_edge_smoother = 25
yield_function_tolerance = 1E-11
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1'
[../]
[./strain]
type = ComputeIncrementalStrain
displacements = 'disp_x disp_y disp_z'
eigenstrain_names = ini_stress
[../]
[./ini_stress]
type = ComputeEigenstrainFromInitialStress
initial_stress = '3 0 0 0 3 0 0 0 1.5'
eigenstrain_name = ini_stress
[../]
[./mc]
type = ComputeMultiPlasticityStress
ep_plastic_tolerance = 1E-11
plastic_models = mc
tangent_operator = nonlinear
min_stepsize = 1
[../]
[]
[Preconditioning]
[./andy]
type = SMP
full = true
petsc_options_iname = '-ksp_type -pc_type -snes_atol -snes_rtol -snes_max_it -snes_type'
petsc_options_value = 'bcgs bjacobi 1E-15 1E-10 10000 test'
[../]
[]
[Executioner]
type = Transient
solve_type = Newton
[]
(modules/solid_mechanics/test/tests/jacobian/cto13.i)
# checking jacobian for nonlinear plasticity (single surface, smoothed MohrCoulomb)
# note: must have min_stepsize=1 otherwise the nonlinearities compound and make the jacobian more inaccurate
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[GlobalParams]
block = 0
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 60
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4
mc_edge_smoother = 25
yield_function_tolerance = 1E-11
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1'
[../]
[./strain]
type = ComputeIncrementalStrain
displacements = 'disp_x disp_y disp_z'
eigenstrain_names = ini_stress
[../]
[./ini_stress]
type = ComputeEigenstrainFromInitialStress
initial_stress = '6 5 4 5 7 2 4 2 2'
eigenstrain_name = ini_stress
[../]
[./mc]
type = ComputeMultiPlasticityStress
ep_plastic_tolerance = 1E-11
plastic_models = mc
tangent_operator = nonlinear
min_stepsize = 1
[../]
[]
[Preconditioning]
[./andy]
type = SMP
full = true
petsc_options_iname = '-ksp_type -pc_type -snes_atol -snes_rtol -snes_max_it -snes_type'
petsc_options_value = 'bcgs bjacobi 1E-15 1E-10 10000 test'
[../]
[]
[Executioner]
type = Transient
solve_type = Newton
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/small_deform2_small_strain.i)
# apply repeated stretches in x, y and z directions, so that mean_stress = 0
# This maps out the yield surface in the octahedral plane for zero mean stress
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '2E-6*x*sin(t)'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '2E-6*y*sin(2*t)'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '-2E-6*z*(sin(t)+sin(2*t))'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[./f]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[./f0]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[./f1]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 20
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 0
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4
mc_edge_smoother = 20
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeIncrementalStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
[../]
[]
[Executioner]
end_time = 100
dt = 1
type = Transient
[]
[Outputs]
file_base = small_deform2_small_strain
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/small_deform3.i)
# apply repeated stretches in x z directions, and smaller stretches along the y direction,
# so that sigma_I = sigma_II,
# which means that lode angle = 30deg.
# The allows yield surface in meridional plane to be mapped out
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '1E-6*x*t'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '0.25E-6*y*sin(t)'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '1E-6*z*t'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./mc_int]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./mc_int_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_internal_parameter
variable = mc_int
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[./internal]
type = PointValue
point = '0 0 0'
variable = mc_int
[../]
[./f]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 50
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningExponential
value_0 = 0
value_residual = 0.8726646 # 50deg
rate = 3000.0
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4
mc_edge_smoother = 20
yield_function_tolerance = 1E-8
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
[../]
[]
[Executioner]
end_time = 30
dt = 1
type = Transient
[]
[Outputs]
file_base = small_deform3
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/notched_plastic_block/biaxial_abbo.i)
# Uses an Abbo et al smoothed version of Mohr-Coulomb (via SolidMechanicsPlasticMohrCoulomb and ComputeMultiPlasticityStress) to simulate the following problem.
# A cubical block is notched around its equator.
# All of its outer surfaces have roller BCs, but the notched region is free to move as needed
# The block is initialised with a high hydrostatic tensile stress
# Without the notch, the BCs do not allow contraction of the block, and this stress configuration is admissible
# With the notch, however, the interior parts of the block are free to move in order to relieve stress, and this causes plastic failure
# The top surface is then pulled upwards (the bottom is fixed because of the roller BCs)
# This causes more failure
[Mesh]
[generated_mesh]
type = GeneratedMeshGenerator
dim = 3
nx = 9
ny = 9
nz = 9
xmin = 0
xmax = 0.1
ymin = 0
ymax = 0.1
zmin = 0
zmax = 0.1
[]
[block_to_remove_xmin]
type = SubdomainBoundingBoxGenerator
bottom_left = '-0.01 -0.01 0.045'
top_right = '0.01 0.11 0.055'
location = INSIDE
block_id = 1
input = generated_mesh
[]
[block_to_remove_xmax]
type = SubdomainBoundingBoxGenerator
bottom_left = '0.09 -0.01 0.045'
top_right = '0.11 0.11 0.055'
location = INSIDE
block_id = 1
input = block_to_remove_xmin
[]
[block_to_remove_ymin]
type = SubdomainBoundingBoxGenerator
bottom_left = '-0.01 -0.01 0.045'
top_right = '0.11 0.01 0.055'
location = INSIDE
block_id = 1
input = block_to_remove_xmax
[]
[block_to_remove_ymax]
type = SubdomainBoundingBoxGenerator
bottom_left = '-0.01 0.09 0.045'
top_right = '0.11 0.11 0.055'
location = INSIDE
block_id = 1
input = block_to_remove_ymin
[]
[remove_block]
type = BlockDeletionGenerator
block = 1
input = block_to_remove_ymax
[]
[]
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Physics/SolidMechanics/QuasiStatic]
[./all]
add_variables = true
incremental = true
generate_output = 'max_principal_stress mid_principal_stress min_principal_stress stress_zz'
eigenstrain_names = ini_stress
[../]
[]
[Postprocessors]
[./uz]
type = PointValue
point = '0 0 0.1'
use_displaced_mesh = false
variable = disp_z
[../]
[./s_zz]
type = ElementAverageValue
use_displaced_mesh = false
variable = stress_zz
[../]
[./num_res]
type = NumResidualEvaluations
[../]
[./nr_its] # num_iters is the average number of NR iterations encountered per element in this timestep
type = ElementAverageValue
variable = num_iters
[../]
[./max_nr_its] # num_iters is the average number of NR iterations encountered in the element in this timestep, so we must get max(max_nr_its) to obtain the max number of iterations
type = ElementExtremeValue
variable = num_iters
[../]
[./runtime]
type = PerfGraphData
data_type = TOTAL
section_name = 'Root'
[../]
[]
[BCs]
# back=zmin, front=zmax, bottom=ymin, top=ymax, left=xmin, right=xmax
[./xmin_xzero]
type = DirichletBC
variable = disp_x
boundary = left
value = 0.0
[../]
[./xmax_xzero]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[../]
[./ymin_yzero]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[../]
[./ymax_yzero]
type = DirichletBC
variable = disp_y
boundary = top
value = 0.0
[../]
[./zmin_zzero]
type = DirichletBC
variable = disp_z
boundary = back
value = '0'
[../]
[./zmax_disp]
type = FunctionDirichletBC
variable = disp_z
boundary = front
function = '1E-6*max(t,0)'
[../]
[]
[AuxVariables]
[./mc_int]
order = CONSTANT
family = MONOMIAL
[../]
[./num_iters]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./mc_int_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_internal_parameter
variable = mc_int
[../]
[./num_iters_auxk]
type = MaterialRealAux
property = plastic_NR_iterations
variable = num_iters
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 5E6
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 35
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 10
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 0.02E6
mc_edge_smoother = 29
yield_function_tolerance = 1E-5
internal_constraint_tolerance = 1E-11
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
youngs_modulus = 16E9
poissons_ratio = 0.25
[../]
[./mc]
type = ComputeMultiPlasticityStress
ep_plastic_tolerance = 1E-11
plastic_models = mc
max_NR_iterations = 1000
debug_fspb = crash
[../]
[./strain_from_initial_stress]
type = ComputeEigenstrainFromInitialStress
initial_stress = '6E6 0 0 0 6E6 0 0 0 6E6'
eigenstrain_name = ini_stress
[../]
[]
[Preconditioning]
[./andy]
type = SMP
full = true
[../]
[]
[Executioner]
start_time = -1
end_time = 10
dt = 1
solve_type = NEWTON
type = Transient
l_tol = 1E-2
nl_abs_tol = 1E-5
nl_rel_tol = 1E-7
l_max_its = 200
nl_max_its = 400
petsc_options_iname = '-pc_type -pc_asm_overlap -sub_pc_type -ksp_type -ksp_gmres_restart'
petsc_options_value = ' asm 2 lu gmres 200'
[]
[Outputs]
file_base = biaxial_abbo
perf_graph = true
exodus = false
csv = true
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/small_deform5.i)
# apply repeated stretches in x z directions, and smaller stretches along the y direction,
# so that sigma_I = sigma_II,
# which means that lode angle = 30deg.
# The allows yield surface in meridional plane to be mapped out
# Use 'cap' smoothing
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '1E-6*x*t'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '0.9E-6*y*sin(t)'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '1E-6*z*t'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./mc_int]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[./iter]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./mc_int_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_internal_parameter
variable = mc_int
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[./iter_auxk]
type = MaterialRealAux
property = plastic_NR_iterations
variable = iter
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[./internal]
type = PointValue
point = '0 0 0'
variable = mc_int
[../]
[./f]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[./iter]
type = PointValue
point = '0 0 0'
variable = iter
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 50
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningExponential
value_0 = 0
value_residual = 0.8726646 # 50deg
rate = 3000.0
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
tip_scheme = cap
mc_tip_smoother = 0
cap_start = 3
cap_rate = 0.8
mc_edge_smoother = 20
yield_function_tolerance = 1E-8
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
[../]
[]
[Executioner]
end_time = 150
dt = 5
type = Transient
[]
[Outputs]
file_base = small_deform5
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/small_deform_hard3.i)
# apply repeated stretches in x z directions, and smaller stretches along the y direction,
# so that sigma_I = sigma_II,
# which means that lode angle = 30deg.
# The allows yield surface in meridional plane to be mapped out
#
# friction_angle = 50deg, friction_angle_residual=51deg, friction_angle_rate = 1E7 (huge)
# cohesion = 10, cohesion_residual = 9.9, cohesion_rate = 1E7 (huge)
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '1E-6*x*t'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '0.25E-6*y*sin(t)'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '1E-6*z*t'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./mc_int]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./mc_int_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_internal_parameter
variable = mc_int
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[./internal]
type = PointValue
point = '0 0 0'
variable = mc_int
[../]
[./f]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningExponential
value_0 = 10
value_residual = 9.9
rate = 1E7
[../]
[./mc_phi]
type = SolidMechanicsHardeningExponential
value_0 = 0.8726646 # 50deg
value_residual = 0.8901179 # 51deg
rate = 1E7
[../]
[./mc_psi]
type = SolidMechanicsHardeningExponential
value_0 = 0
value_residual = 0.8726646 # 50deg
rate = 3000
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4
mc_edge_smoother = 20
yield_function_tolerance = 1E-6
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
debug_jac_at_stress = '10 1 2 1 11 -3 2 -3 8'
debug_jac_at_pm = 1
debug_jac_at_intnl = 1
debug_stress_change = 1E-5
debug_pm_change = 1E-6
debug_intnl_change = 1E-6
[../]
[]
[Executioner]
end_time = 30
dt = 1
type = Transient
[]
[Outputs]
file_base = small_deform_hard3
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/random.i)
# apply many random large deformations, checking that the algorithm returns correctly to
# the yield surface each time.
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1000
ny = 125
nz = 1
xmin = 0
xmax = 1000
ymin = 0
ymax = 125
zmin = 0
zmax = 1
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[ICs]
[./x]
type = RandomIC
min = -0.1
max = 0.1
variable = disp_x
[../]
[./y]
type = RandomIC
min = -0.1
max = 0.1
variable = disp_y
[../]
[./z]
type = RandomIC
min = -0.1
max = 0.1
variable = disp_z
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '0'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '0'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '0'
[../]
[]
[AuxVariables]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[./iter]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[./iter]
type = MaterialRealAux
property = plastic_NR_iterations
variable = iter
[../]
[]
[Postprocessors]
[./yield_fcn_at_zero]
type = PointValue
point = '0 0 0'
variable = yield_fcn
outputs = 'console'
[../]
[./should_be_zero]
type = FunctionValuePostprocessor
function = should_be_zero_fcn
[../]
[./av_iter]
type = ElementAverageValue
variable = iter
outputs = 'console'
[../]
[]
[Functions]
[./should_be_zero_fcn]
type = ParsedFunction
expression = 'if(a<1E-3,0,a)'
symbol_names = 'a'
symbol_values = 'yield_fcn_at_zero'
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 1E3
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 30
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 0.1E3
mc_edge_smoother = 25
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-6
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0.7E7 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
max_NR_iterations = 1000
ep_plastic_tolerance = 1E-6
min_stepsize = 1E-3
plastic_models = mc
debug_fspb = crash
deactivation_scheme = safe
[../]
[]
[Executioner]
end_time = 1
dt = 1
type = Transient
[]
[Outputs]
file_base = random
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/except2.i)
# checking for exception error messages
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '1E-6*x'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '1E-6*y'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '1E-6*z'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 30
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 1
mc_edge_smoother = 25
mc_lode_cutoff = -1.0E-6
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
debug_jac_at_stress = '10 0 0 0 10 0 0 0 10'
debug_jac_at_pm = 1
debug_jac_at_intnl = 1
debug_stress_change = 1E-5
debug_pm_change = 1E-6
debug_intnl_change = 1E-6
[../]
[]
[Executioner]
end_time = 1
dt = 1
type = Transient
[]
[Outputs]
file_base = except2
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/except1.i)
# checking for exception error messages
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '1E-6*x'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '1E-6*y'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '1E-6*z'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 30
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 45
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 1
mc_edge_smoother = 25
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
debug_jac_at_stress = '10 0 0 0 10 0 0 0 10'
debug_jac_at_pm = 1
debug_jac_at_intnl = 1
debug_stress_change = 1E-5
debug_pm_change = 1E-6
debug_intnl_change = 1E-6
[../]
[]
[Executioner]
end_time = 1
dt = 1
type = Transient
[]
[Outputs]
file_base = except1
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/small_deform_hard2.i)
# apply uniform stretches in x, y and z directions.
# let friction_angle = 60deg, friction_angle_residual=10deg, friction_angle_rate = 0.5E4
# With cohesion = C, friction_angle = phi, tip_smoother = T, the
# algorithm should return to
# sigma_m = (C*Cos(phi) - T)/Sin(phi)
# Or, when T=C,
# phi = 2*pi*n - 2*arctan(sigma_m/C)
# This allows checking of the relationship for phi
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '1E-6*x*t'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '1E-6*y*t'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '1E-6*z*t'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./mc_int]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./mc_int_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_internal_parameter
variable = mc_int
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[./internal]
type = PointValue
point = '0 0 0'
variable = mc_int
[../]
[./f]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10
[../]
[./mc_phi]
type = SolidMechanicsHardeningExponential
value_0 = 1.04719755 # 60deg
value_residual = 0.17453293 # 10deg
rate = 0.5E4
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 10
mc_edge_smoother = 25
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
debug_jac_at_stress = '10 1 2 1 10 3 2 3 10'
debug_jac_at_pm = 1
debug_jac_at_intnl = 1E-3
debug_stress_change = 1E-5
debug_pm_change = 1E-6
debug_intnl_change = 1E-6
[../]
[]
[Executioner]
end_time = 10
dt = 1
type = Transient
[]
[Outputs]
file_base = small_deform_hard2
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/multi/mc_wpt_1.i)
# checking for small deformation
# A single element is stretched by 1E-6m in x,y and z directions.
# stress_zz = Youngs Modulus*Strain = 2E6*1E-6 = 2 Pa
# wpt_tensile_strength is set to 1Pa
# Then the final stress should return to the yeild surface and its value should be 1pa.
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./bottomx]
type = DirichletBC
variable = disp_x
boundary = back
value = 0.0
[../]
[./bottomy]
type = DirichletBC
variable = disp_y
boundary = back
value = 0.0
[../]
[./bottomz]
type = DirichletBC
variable = disp_z
boundary = back
value = 0.0
[../]
[./topx]
type = DirichletBC
variable = disp_x
boundary = front
value = 0E-6
[../]
[./topy]
type = DirichletBC
variable = disp_y
boundary = front
value = 0E-6
[../]
[./topz]
type = DirichletBC
variable = disp_z
boundary = front
value = 1E-6
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 100
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 60
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4
mc_edge_smoother = 25
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-9
[../]
[./str]
type = SolidMechanicsHardeningConstant
value = 1
[../]
[./wpt]
type = SolidMechanicsPlasticWeakPlaneTensile
tensile_strength = str
yield_function_tolerance = 1E-6
internal_constraint_tolerance = 1E-5
[../]
[]
[Materials]
[./mc]
type = FiniteStrainMultiPlasticity
block = 0
disp_x = disp_x
disp_y = disp_y
disp_z = disp_z
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
ep_plastic_tolerance = 1E-9
plastic_models = 'mc wpt'
debug_fspb = crash
debug_jac_at_stress = '10 0 0 0 10 0 0 0 10'
debug_jac_at_pm = 1
debug_jac_at_intnl = 1
debug_stress_change = 1E-5
debug_pm_change = 1E-6
debug_intnl_change = 1E-6
[../]
[]
[Executioner]
end_time = 1
dt = 1
type = Transient
[]
[Outputs]
file_base = mc_wpt_1
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/small_deform_hard1.i)
# apply uniform stretches in x, y and z directions.
# let mc_cohesion = 10, mc_cohesion_residual = 2, mc_cohesion_rate =
# With cohesion = C, friction_angle = 60deg, tip_smoother = 4, the
# algorithm should return to
# sigma_m = (C*Cos(60) - 4)/Sin(60)
# This allows checking of the relationship for C
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '1E-6*x*t'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '1E-6*y*t'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '1E-6*z*t'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./mc_int]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./mc_int_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_internal_parameter
variable = mc_int
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[./internal]
type = PointValue
point = '0 0 0'
variable = mc_int
[../]
[./f]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningExponential
value_0 = 10
value_residual = 2
rate = 1E4
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 60
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4
mc_edge_smoother = 25
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
debug_jac_at_stress = '10 1 2 1 10 3 2 3 10'
debug_jac_at_pm = 1
debug_jac_at_intnl = 1E-4
debug_stress_change = 1E-5
debug_pm_change = 1E-6
debug_intnl_change = 1E-8
[../]
[]
[Executioner]
end_time = 10
dt = 1
type = Transient
[]
[Outputs]
file_base = small_deform_hard1
exodus = false
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/uni_axial1_small_strain.i)
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
# back = zmin
# front = zmax
# bottom = ymin
# top = ymax
# left = xmin
# right = xmax
[./xmin_xzero]
type = DirichletBC
variable = disp_x
boundary = 'left'
value = '0'
[../]
[./ymin_yzero]
type = DirichletBC
variable = disp_y
boundary = 'bottom'
value = '0'
[../]
[./zmin_zzero]
type = DirichletBC
variable = disp_z
boundary = 'back'
value = '0'
[../]
[./zmax_disp]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front'
function = '-1E-3*t'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./mc_int]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./mc_int_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_internal_parameter
variable = mc_int
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[./mc_int]
type = PointValue
point = '0 0 0'
variable = mc_int
[../]
[./f]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10E6
[../]
[./mc_phi]
type = SolidMechanicsHardeningExponential
value_0 = 0
value_residual = 0.6981317 # 40deg
rate = 10000
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 0
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 0
mc_edge_smoother = 25
yield_function_tolerance = 1E-3
internal_constraint_tolerance = 1E-10
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '5.77E10 3.85E10' # young = 100Gpa, poisson = 0.3
[../]
[./strain]
type = ComputeIncrementalStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-10
plastic_models = mc
max_NR_iterations = 1000
debug_fspb = crash
[../]
[]
[Preconditioning]
[./andy]
type = SMP
full = true
[../]
[]
[Executioner]
end_time = 0.5
dt = 0.05
solve_type = NEWTON
type = Transient
line_search = 'none'
nl_rel_tol = 1E-10
l_tol = 1E-3
l_max_its = 200
nl_max_its = 10
petsc_options_iname = '-pc_type -pc_asm_overlap -sub_pc_type -ksp_type -ksp_gmres_restart'
petsc_options_value = ' asm 2 lu gmres 200'
[]
[Outputs]
file_base = uni_axial1_small_strain
exodus = true
[./csv]
type = CSV
[../]
[]
(modules/solid_mechanics/test/tests/multi/special_rock1.i)
# Plasticity models:
# Mohr-Coulomb with cohesion = 40MPa, friction angle = 35deg, dilation angle = 5deg
# Tensile with strength = 1MPa
#
# Lame lambda = 1GPa. Lame mu = 1.3GPa
#
# A line of elements is perturbed randomly, and return to the yield surface at each quadpoint is checked
#
# NOTE: The yield function tolerances here are set at 100-times what i would usually use
# This is because otherwise the test fails on the 'pearcey' architecture.
# This is because identical stress tensors yield slightly different eigenvalues
# (and hence return-map residuals) on 'pearcey' than elsewhere, which results in
# a different number of NR iterations are needed to return to the yield surface.
# This is presumably because of compiler internals, or the BLAS routines being
# optimised differently or something similar.
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1000
ny = 125
nz = 1
xmin = 0
xmax = 1000
ymin = 0
ymax = 125
zmin = 0
zmax = 1
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[GlobalParams]
volumetric_locking_correction=true
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[ICs]
[./x]
type = RandomIC
min = -0.1
max = 0.1
variable = disp_x
[../]
[./y]
type = RandomIC
min = -0.1
max = 0.1
variable = disp_y
[../]
[./z]
type = RandomIC
min = -0.1
max = 0.1
variable = disp_z
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '0'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '0'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '0'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./linesearch]
order = CONSTANT
family = MONOMIAL
[../]
[./ld]
order = CONSTANT
family = MONOMIAL
[../]
[./constr_added]
order = CONSTANT
family = MONOMIAL
[../]
[./iter]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./linesearch]
type = MaterialRealAux
property = plastic_linesearch_needed
variable = linesearch
[../]
[./ld]
type = MaterialRealAux
property = plastic_linear_dependence_encountered
variable = ld
[../]
[./constr_added]
type = MaterialRealAux
property = plastic_constraints_added
variable = constr_added
[../]
[./iter]
type = MaterialRealAux
property = plastic_NR_iterations
variable = iter
[../]
[]
[Postprocessors]
[./max_iter]
type = ElementExtremeValue
variable = iter
outputs = console
[../]
[./av_linesearch]
type = ElementAverageValue
variable = linesearch
outputs = 'console csv'
[../]
[./av_ld]
type = ElementAverageValue
variable = ld
outputs = 'console csv'
[../]
[./av_constr_added]
type = ElementAverageValue
variable = constr_added
outputs = 'console csv'
[../]
[./av_iter]
type = ElementAverageValue
variable = iter
outputs = 'console csv'
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 4E7
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 35
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningConstant
value = 5
convert_to_radians = true
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulombMulti
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
use_custom_returnMap = false
yield_function_tolerance = 1.0E+2 # Note larger value
shift = 1.0E+2 # Note larger value
internal_constraint_tolerance = 1.0E-7
[../]
[./mc_smooth]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
mc_tip_smoother = 4E6
yield_function_tolerance = 1.0
internal_constraint_tolerance = 1.0E-7
[../]
[./ts]
type = SolidMechanicsHardeningConstant
value = 1E6
[../]
[./tensile]
type = SolidMechanicsPlasticTensileMulti
tensile_strength = ts
yield_function_tolerance = 1.0E+2 # Note larger value
shift = 1.0E+2 # Note larger value
internal_constraint_tolerance = 1.0E-7
use_custom_returnMap = false
use_custom_cto = false
[../]
[./tensile_smooth]
type = SolidMechanicsPlasticTensile
tensile_strength = ts
tensile_tip_smoother = 1E5
yield_function_tolerance = 1.0
internal_constraint_tolerance = 1.0E-7
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '1.0E9 1.3E9'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./multi]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-5 # Note larger value, to match the larger yield_function_tolerances
plastic_models = 'tensile mc'
max_NR_iterations = 5
specialIC = 'rock'
deactivation_scheme = 'safe'
min_stepsize = 1
max_stepsize_for_dumb = 1
debug_fspb = crash
debug_jac_at_stress = '10 0 0 0 10 0 0 0 10'
debug_jac_at_pm = '1 1 1 1'
debug_jac_at_intnl = '1 1 1 1'
debug_stress_change = 1E1
debug_pm_change = '1E-6 1E-6 1E-6 1E-6'
debug_intnl_change = '1E-6 1E-6 1E-6 1E-6'
[../]
[]
[Executioner]
end_time = 1
dt = 1
type = Transient
[]
[Outputs]
file_base = special_rock1
exodus = false
csv = true
[]
(modules/solid_mechanics/test/tests/mohr_coulomb/small_deform6.i)
# apply repeated stretches in z direction, and smaller stretches in the x and y directions
# so that sigma_II = sigma_III,
# which means that lode angle = -30deg.
# The allows yield surface in meridional plane to be mapped out
# Using cap smoothing
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
xmin = -0.5
xmax = 0.5
ymin = -0.5
ymax = 0.5
zmin = -0.5
zmax = 0.5
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[./disp_z]
[../]
[]
[Kernels]
[SolidMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
[]
[BCs]
[./x]
type = FunctionDirichletBC
variable = disp_x
boundary = 'front back'
function = '0.9E-6*x*sin(t)'
[../]
[./y]
type = FunctionDirichletBC
variable = disp_y
boundary = 'front back'
function = '0.9E-6*y*sin(t)'
[../]
[./z]
type = FunctionDirichletBC
variable = disp_z
boundary = 'front back'
function = '1E-6*z*t'
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./mc_int]
order = CONSTANT
family = MONOMIAL
[../]
[./yield_fcn]
order = CONSTANT
family = MONOMIAL
[../]
[./iter]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./stress_xz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xz
index_i = 0
index_j = 2
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_yz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yz
index_i = 1
index_j = 2
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./mc_int_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_internal_parameter
variable = mc_int
[../]
[./yield_fcn_auxk]
type = MaterialStdVectorAux
index = 0
property = plastic_yield_function
variable = yield_fcn
[../]
[./iter_auxk]
type = MaterialRealAux
property = plastic_NR_iterations
variable = iter
[../]
[]
[Postprocessors]
[./s_xx]
type = PointValue
point = '0 0 0'
variable = stress_xx
[../]
[./s_xy]
type = PointValue
point = '0 0 0'
variable = stress_xy
[../]
[./s_xz]
type = PointValue
point = '0 0 0'
variable = stress_xz
[../]
[./s_yy]
type = PointValue
point = '0 0 0'
variable = stress_yy
[../]
[./s_yz]
type = PointValue
point = '0 0 0'
variable = stress_yz
[../]
[./s_zz]
type = PointValue
point = '0 0 0'
variable = stress_zz
[../]
[./internal]
type = PointValue
point = '0 0 0'
variable = mc_int
[../]
[./f]
type = PointValue
point = '0 0 0'
variable = yield_fcn
[../]
[./iter]
type = PointValue
point = '0 0 0'
variable = iter
[../]
[]
[UserObjects]
[./mc_coh]
type = SolidMechanicsHardeningConstant
value = 10
[../]
[./mc_phi]
type = SolidMechanicsHardeningConstant
value = 50
convert_to_radians = true
[../]
[./mc_psi]
type = SolidMechanicsHardeningExponential
value_0 = 0
value_residual = 0.8726646 # 50deg
rate = 3000.0
[../]
[./mc]
type = SolidMechanicsPlasticMohrCoulomb
cohesion = mc_coh
friction_angle = mc_phi
dilation_angle = mc_psi
tip_scheme = cap
mc_tip_smoother = 0
cap_start = 3
cap_rate = 0.8
mc_edge_smoother = 20
yield_function_tolerance = 1E-8
internal_constraint_tolerance = 1E-9
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeElasticityTensor
block = 0
fill_method = symmetric_isotropic
C_ijkl = '0 1E7'
[../]
[./strain]
type = ComputeFiniteStrain
block = 0
displacements = 'disp_x disp_y disp_z'
[../]
[./mc]
type = ComputeMultiPlasticityStress
block = 0
ep_plastic_tolerance = 1E-9
plastic_models = mc
debug_fspb = crash
[../]
[]
[Executioner]
end_time = 30
dt = 1
type = Transient
[]
[Outputs]
file_base = small_deform6
exodus = false
[./csv]
type = CSV
[../]
[]