- absolute_tolerance1e-11Absolute convergence tolerance for Newton iteration
Default:1e-11
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Absolute convergence tolerance for Newton iteration
- acceptable_multiplier10Factor applied to relative and absolute tolerance for acceptable convergence if iterations are no longer making progress
Default:10
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Factor applied to relative and absolute tolerance for acceptable convergence if iterations are no longer making progress
- adaptive_substeppingFalseUse adaptive substepping, where the number of substeps is successively doubled until the return mapping model successfully converges or the maximum number of substeps is reached.
Default:False
C++ Type:bool
Controllable:No
Description:Use adaptive substepping, where the number of substeps is successively doubled until the return mapping model successfully converges or the maximum number of substeps is reached.
- automatic_differentiation_return_mappingFalseWhether to use automatic differentiation to compute the derivative.
Default:False
C++ Type:bool
Controllable:No
Description:Whether to use automatic differentiation to compute the derivative.
- base_nameOptional parameter that defines a prefix for all material properties related to this stress update model. This allows for multiple models of the same type to be used without naming conflicts.
C++ Type:std::string
Controllable:No
Description:Optional parameter that defines a prefix for all material properties related to this stress update model. This allows for multiple models of the same type to be used without naming conflicts.
- blockThe list of blocks (ids or names) that this object will be applied
C++ Type:std::vector<SubdomainName>
Controllable:No
Description:The list of blocks (ids or names) that this object will be applied
- boundaryThe list of boundaries (ids or names) from the mesh where this object applies
C++ Type:std::vector<BoundaryName>
Controllable:No
Description:The list of boundaries (ids or names) from the mesh where this object applies
- constant_onNONEWhen ELEMENT, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps.When SUBDOMAIN, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps. Evaluations on element qps will be skipped
Default:NONE
C++ Type:MooseEnum
Controllable:No
Description:When ELEMENT, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps.When SUBDOMAIN, MOOSE will only call computeQpProperties() for the 0th quadrature point, and then copy that value to the other qps. Evaluations on element qps will be skipped
- declare_suffixAn optional suffix parameter that can be appended to any declared 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 declared properties. The suffix will be prepended with a '_' character.
- hardening_constantHardening slope
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Hardening slope
- hardening_functionTrue stress as a function of plastic strain
C++ Type:FunctionName
Unit:(no unit assumed)
Controllable:No
Description:True stress as a function of plastic strain
- max_inelastic_increment0.0001The maximum inelastic strain increment allowed in a time step
Default:0.0001
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:The maximum inelastic strain increment allowed in a time step
- maximum_number_substeps25The maximum number of substeps allowed before cutting the time step.
Default:25
C++ Type:unsigned int
Controllable:No
Description:The maximum number of substeps allowed before cutting the time step.
- relative_tolerance1e-08Relative convergence tolerance for Newton iteration
Default:1e-08
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Relative convergence tolerance for Newton iteration
- temperatureCoupled Temperature
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Coupled Temperature
- use_substep_integration_errorFalseIf true, it establishes a substep size that will yield, at most,the creep numerical integration error given by substep_strain_tolerance.
Default:False
C++ Type:bool
Controllable:No
Description:If true, it establishes a substep size that will yield, at most,the creep numerical integration error given by substep_strain_tolerance.
- use_substeppingNONEWhether and how to use substepping
Default:NONE
C++ Type:MooseEnum
Controllable:No
Description:Whether and how to use substepping
- yield_stressThe point at which plastic strain begins accumulating
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:The point at which plastic strain begins accumulating
- yield_stress_functionYield stress as a function of temperature
C++ Type:FunctionName
Unit:(no unit assumed)
Controllable:No
Description:Yield stress as a function of temperature
Isotropic Plasticity Stress Update
This class uses the discrete material in a radial return isotropic plasticity model. This class is one of the basic radial return constitutive models, yet it can be used in conjunction with other creep and plasticity materials for more complex simulations.
Description
In this numerical approach, a trial stress is calculated at the start of each simulation time increment; the trial stress calculation assumed all of the new strain increment is elastic strain:
The algorithms checks to see if the trial stress state is outside of the yield surface, as shown in the figure to the right. If the stress state is outside of the yield surface, the algorithm recomputes the scalar effective inelastic strain required to return the stress state to the yield surface. This approach is given the name Radial Return because the yield surface used is the von Mises yield surface: in the deviatoric stress space, this yield surface has the shape of a circle, and the scalar inelastic strain is assumed to always be directed at the circle center.
Recompute Iterations on the Effective Plastic Strain Increment
The recompute radial return materials each individually calculate, using the Newton Method, the amount of effective inelastic strain required to return the stress state to the yield surface.
where the change in the iterative effective inelastic strain is defined as the yield surface over the derivative of the yield surface with respect to the inelastic strain increment.
In isotropic linear hardening plasticity, with the hardening function , the effective plastic strain increment has the form: where is the isotropic shear modulus, and is the scalar von Mises trial stress.
This class calculates an effective trial stress, an effective scalar plastic strain increment, and the derivative of the scalar effective plastic strain increment; these values are passed to the RadialReturnStressUpdate to compute the radial return stress increment. This isotropic plasticity class also computes the plastic strain as a stateful material property.
This class is based on the implicit integration algorithm in Dunne and Petrinic (2005) pg. 146–149.
The ADIsotropicPlasticityStressUpdate version of this class uses forward mode automatic differentiation to provide all necessary material property derivatives to assemble a perfect Jacobian (this replaces the approximated tangent operator).
Example Input File Syntax
[./isotropic_plasticity]
type = IsotropicPlasticityStressUpdate
yield_stress = 50.0
hardening_function = hf
[../]IsotropicPlasticityStressUpdate must be run in conjunction with the inelastic strain return mapping stress calculator as shown below:
[./radial_return_stress]
type = ComputeMultipleInelasticStress
tangent_operator = elastic
inelastic_models = 'isotropic_plasticity'
[../]Input Parameters
- apply_strainTrueFlag to apply strain. Used for testing.
Default:True
C++ Type:bool
Controllable:No
Description:Flag to apply strain. Used for testing.
- 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.
- effective_inelastic_strain_nameeffective_plastic_strainName of the material property that stores the effective inelastic strain
Default:effective_plastic_strain
C++ Type:std::string
Controllable:No
Description:Name of the material property that stores the effective inelastic strain
- enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:Yes
Description:Set the enabled status of the MooseObject.
- implicitTrueDetermines whether this object is calculated using an implicit or explicit form
Default:True
C++ Type:bool
Controllable:No
Description:Determines whether this object is calculated using an implicit or explicit form
- search_methodnearest_node_connected_sidesChoice of search algorithm. All options begin by finding the nearest node in the primary boundary to a query point in the secondary boundary. In the default nearest_node_connected_sides algorithm, primary boundary elements are searched iff that nearest node is one of their nodes. This is fast to determine via a pregenerated node-to-elem map and is robust on conforming meshes. In the optional all_proximate_sides algorithm, primary boundary elements are searched iff they touch that nearest node, even if they are not topologically connected to it. This is more CPU-intensive but is necessary for robustness on any boundary surfaces which has disconnections (such as Flex IGA meshes) or non-conformity (such as hanging nodes in adaptively h-refined meshes).
Default:nearest_node_connected_sides
C++ Type:MooseEnum
Options:nearest_node_connected_sides, all_proximate_sides
Controllable:No
Description:Choice of search algorithm. All options begin by finding the nearest node in the primary boundary to a query point in the secondary boundary. In the default nearest_node_connected_sides algorithm, primary boundary elements are searched iff that nearest node is one of their nodes. This is fast to determine via a pregenerated node-to-elem map and is robust on conforming meshes. In the optional all_proximate_sides algorithm, primary boundary elements are searched iff they touch that nearest node, even if they are not topologically connected to it. This is more CPU-intensive but is necessary for robustness on any boundary surfaces which has disconnections (such as Flex IGA meshes) or non-conformity (such as hanging nodes in adaptively h-refined meshes).
- seed0The seed for the master random number generator
Default:0
C++ Type:unsigned int
Controllable:No
Description:The seed for the master random number generator
- substep_strain_tolerance0.1Maximum ratio of the initial elastic strain increment at start of the return mapping solve to the maximum inelastic strain allowable in a single substep. Reduce this value to increase the number of substeps
Default:0.1
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Maximum ratio of the initial elastic strain increment at start of the return mapping solve to the maximum inelastic strain allowable in a single substep. Reduce this value to increase the number of substeps
- 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
- internal_solve_full_iteration_historyFalseSet true to output full internal Newton iteration history at times determined by `internal_solve_output_on`. If false, only a summary is output.
Default:False
C++ Type:bool
Controllable:No
Description:Set true to output full internal Newton iteration history at times determined by `internal_solve_output_on`. If false, only a summary is output.
- internal_solve_output_onon_errorWhen to output internal Newton solve information
Default:on_error
C++ Type:MooseEnum
Options:never, on_error, always
Controllable:No
Description:When to output internal Newton solve information
Debug Parameters
- output_propertiesList of material properties, from this material, to output (outputs must also be defined to an output type)
C++ Type:std::vector<std::string>
Controllable:No
Description:List of material properties, from this material, to output (outputs must also be defined to an output type)
- outputsnone Vector of output names where you would like to restrict the output of variables(s) associated with this object
Default:none
C++ Type:std::vector<OutputName>
Controllable:No
Description:Vector of output names where you would like to restrict the output of variables(s) associated with this object
Outputs 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/combined/test/tests/combined_plasticity_temperature/plasticity_temperature_dep_yield.i)
- (modules/combined/test/tests/inelastic_strain/elas_plas/elas_plas_nl1_cycle.i)
- (modules/solid_mechanics/test/tests/combined_creep_plasticity/combined_stress_prescribed.i)
- (modules/combined/test/tests/combined_plasticity_temperature/ad_plasticity_temperature_dep_yield.i)
- (modules/solid_mechanics/test/tests/recompute_radial_return/isotropic_plasticity_incremental_strain.i)
- (modules/solid_mechanics/test/tests/strain_energy_density/rate_incr_model_elas_plas.i)
- (modules/solid_mechanics/test/tests/combined_creep_plasticity/combined_creep_plasticity.i)
- (modules/combined/test/tests/inelastic_strain/elas_plas/elas_plas_nl1.i)
- (modules/solid_mechanics/test/tests/combined_creep_plasticity/combined_creep_plasticity_start_time.i)
- (modules/solid_mechanics/test/tests/power_law_creep/composite_power_law_creep_plasticity.i)
- (modules/solid_mechanics/test/tests/power_law_creep/composite_power_law_creep_onePhaseMulti.i)
- (modules/solid_mechanics/test/tests/strain_energy_density/incr_model_elas_plas.i)
- (modules/solid_mechanics/test/tests/recompute_radial_return/cp_affine_plasticity.i)
- (modules/solid_mechanics/test/tests/combined_creep_plasticity/creepWithPlasticity.i)
- (modules/combined/test/tests/phase_field_fracture/crack2d_computeCrackedStress_finitestrain_plastic.i)
- (modules/solid_mechanics/test/tests/recompute_radial_return/affine_plasticity.i)
- (modules/solid_mechanics/test/tests/recompute_radial_return/isotropic_plasticity_finite_strain.i)
- (modules/solid_mechanics/test/tests/recompute_radial_return/isotropic_plasticity_errors.i)
- (modules/solid_mechanics/test/tests/material_limit_time_step/elas_plas/nafems_nl1_lim.i)
- (modules/solid_mechanics/test/tests/power_law_creep/cp_power_law_creep.i)
References
- Fionn Dunne and Nik Petrinic.
Introduction to Computational Plasticity.
Oxford University Press on Demand, 2005.[BibTeX]
@book{dunne2005introduction, author = "Dunne, Fionn and Petrinic, Nik", title = "Introduction to Computational Plasticity", year = "2005", publisher = "Oxford University Press on Demand" }
(modules/solid_mechanics/test/tests/recompute_radial_return/isotropic_plasticity_finite_strain.i)
# This simulation uses the piece-wise linear strain hardening model
# with the incremental small strain formulation; incremental small strain
# is required to produce the strain_increment for the DiscreteRadialReturnStressIncrement
# class, which handles the calculation of the stress increment to return
# to the yield surface in a J2 (isotropic) plasticity problem.
#
# This test assumes a Poissons ratio of 0.3 and applies a displacement loading
# condition on the top in the y direction.
#
# An identical problem was run in Abaqus on a similar 1 element mesh and was used
# to verify the SolidMechanics solution; this SolidMechanics code matches the
# SolidMechanics solution.
#
# Mechanical strain is the sum of the elastic and plastic strains but is different
# from total strain in cases with eigen strains, e.g. thermal strain.
[Mesh]
file = 1x1x1cube.e
[]
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Functions]
[./top_pull]
type = ParsedFunction
expression = t*(0.0625)
[../]
[./hf]
type = PiecewiseLinear
x = '0 0.001 0.003 0.023'
y = '50 52 54 56'
[../]
[]
[Physics/SolidMechanics/QuasiStatic]
[./all]
strain = FINITE
add_variables = true
generate_output = 'stress_yy plastic_strain_xx plastic_strain_yy plastic_strain_zz'
[../]
[]
[BCs]
[./y_pull_function]
type = FunctionDirichletBC
variable = disp_y
boundary = 5
function = top_pull
[../]
[./x_bot]
type = DirichletBC
variable = disp_x
boundary = 4
value = 0.0
[../]
[./y_bot]
type = DirichletBC
variable = disp_y
boundary = 3
value = 0.0
[../]
[./z_bot]
type = DirichletBC
variable = disp_z
boundary = 2
value = 0.0
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
youngs_modulus = 2.1e5
poissons_ratio = 0.3
[../]
[./isotropic_plasticity]
type = IsotropicPlasticityStressUpdate
yield_stress = 50.0
hardening_function = hf
[../]
[./radial_return_stress]
type = ComputeMultipleInelasticStress
tangent_operator = elastic
inelastic_models = 'isotropic_plasticity'
[../]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options = '-snes_ksp_ew'
petsc_options_iname = '-ksp_gmres_restart'
petsc_options_value = '101'
line_search = 'none'
l_max_its = 50
nl_max_its = 50
nl_rel_tol = 1e-12
nl_abs_tol = 1e-10
l_tol = 1e-9
start_time = 0.0
end_time = 0.075
dt = 0.00125
dtmin = 0.0001
[]
[Outputs]
[./out]
type = Exodus
elemental_as_nodal = true
[../]
[]
(modules/solid_mechanics/test/tests/recompute_radial_return/isotropic_plasticity_finite_strain.i)
# This simulation uses the piece-wise linear strain hardening model
# with the incremental small strain formulation; incremental small strain
# is required to produce the strain_increment for the DiscreteRadialReturnStressIncrement
# class, which handles the calculation of the stress increment to return
# to the yield surface in a J2 (isotropic) plasticity problem.
#
# This test assumes a Poissons ratio of 0.3 and applies a displacement loading
# condition on the top in the y direction.
#
# An identical problem was run in Abaqus on a similar 1 element mesh and was used
# to verify the SolidMechanics solution; this SolidMechanics code matches the
# SolidMechanics solution.
#
# Mechanical strain is the sum of the elastic and plastic strains but is different
# from total strain in cases with eigen strains, e.g. thermal strain.
[Mesh]
file = 1x1x1cube.e
[]
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Functions]
[./top_pull]
type = ParsedFunction
expression = t*(0.0625)
[../]
[./hf]
type = PiecewiseLinear
x = '0 0.001 0.003 0.023'
y = '50 52 54 56'
[../]
[]
[Physics/SolidMechanics/QuasiStatic]
[./all]
strain = FINITE
add_variables = true
generate_output = 'stress_yy plastic_strain_xx plastic_strain_yy plastic_strain_zz'
[../]
[]
[BCs]
[./y_pull_function]
type = FunctionDirichletBC
variable = disp_y
boundary = 5
function = top_pull
[../]
[./x_bot]
type = DirichletBC
variable = disp_x
boundary = 4
value = 0.0
[../]
[./y_bot]
type = DirichletBC
variable = disp_y
boundary = 3
value = 0.0
[../]
[./z_bot]
type = DirichletBC
variable = disp_z
boundary = 2
value = 0.0
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
youngs_modulus = 2.1e5
poissons_ratio = 0.3
[../]
[./isotropic_plasticity]
type = IsotropicPlasticityStressUpdate
yield_stress = 50.0
hardening_function = hf
[../]
[./radial_return_stress]
type = ComputeMultipleInelasticStress
tangent_operator = elastic
inelastic_models = 'isotropic_plasticity'
[../]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options = '-snes_ksp_ew'
petsc_options_iname = '-ksp_gmres_restart'
petsc_options_value = '101'
line_search = 'none'
l_max_its = 50
nl_max_its = 50
nl_rel_tol = 1e-12
nl_abs_tol = 1e-10
l_tol = 1e-9
start_time = 0.0
end_time = 0.075
dt = 0.00125
dtmin = 0.0001
[]
[Outputs]
[./out]
type = Exodus
elemental_as_nodal = true
[../]
[]
(modules/combined/test/tests/combined_plasticity_temperature/plasticity_temperature_dep_yield.i)
#
# This is a test of the piece-wise linear strain hardening model using the
# small strain formulation. This test exercises the temperature-dependent
# yield stress.
#
# Test procedure:
# 1. The element is pulled to and then beyond the yield stress for a given
# temperature.
# 2. The displacement is then constant while the temperature increases and
# the yield stress decreases. This results in a lower stress with more
# plastic strain.
# 3. The temperature decreases beyond its original value giving a higher
# yield stress. The displacement increases, causing increases stress to
# the new yield stress.
# 4. The temperature and yield stress are constant with increasing
# displacement giving a constant stress and more plastic strain.
#
# Plotting total_strain_yy on the x axis and stress_yy on the y axis shows
# the stress history in a clear way.
#
# s |
# t | *****
# r | *
# e | ***** *
# s | * * *
# s | * *
# |*
# +------------------
# total strain
#
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
[]
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Physics/SolidMechanics/QuasiStatic]
[./all]
strain = FINITE
incremental = true
add_variables = true
generate_output = 'stress_yy plastic_strain_xx plastic_strain_yy plastic_strain_zz'
[../]
[]
[Variables]
[./temp]
order = FIRST
family = LAGRANGE
[../]
[]
[Functions]
[./top_pull]
type = PiecewiseLinear
x = '0 1 2 4 5 6'
y = '0 0.025 0.05 0.05 0.06 0.085'
[../]
[./yield]
type = PiecewiseLinear
x = '400 500 600'
y = '6e3 5e3 4e3'
[../]
[./temp]
type = PiecewiseLinear
x = '0 1 2 3 4'
y = '500 500 500 600 400'
[../]
[]
[Kernels]
[./heat]
type = HeatConduction
variable = temp
[../]
[]
[BCs]
[./y_pull_function]
type = FunctionDirichletBC
variable = disp_y
boundary = top
function = top_pull
[../]
[./x_bot]
type = DirichletBC
variable = disp_x
boundary = left
value = 0.0
[../]
[./y_bot]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[../]
[./z_bot]
type = DirichletBC
variable = disp_z
boundary = back
value = 0.0
[../]
[./temp]
type = FunctionDirichletBC
variable = temp
function = temp
boundary = left
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
block = 0
youngs_modulus = 2.0e5
poissons_ratio = 0.3
[../]
[./creep_plas]
type = ComputeMultipleInelasticStress
tangent_operator = elastic
block = 0
inelastic_models = 'plasticity'
max_iterations = 50
absolute_tolerance = 1e-05
[../]
[./plasticity]
type = IsotropicPlasticityStressUpdate
block = 0
hardening_constant = 0
yield_stress_function = yield
temperature = temp
[../]
[./heat_conduction]
type = HeatConductionMaterial
block = 0
specific_heat = 1
thermal_conductivity = 1
[../]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options = '-snes_ksp_ew'
petsc_options_iname = '-ksp_gmres_restart -pc_type -pc_hypre_type -pc_hypre_boomeramg_max_iter'
petsc_options_value = '201 hypre boomeramg 4'
line_search = 'none'
l_max_its = 100
nl_max_its = 100
nl_rel_tol = 1e-12
nl_abs_tol = 1e-10
l_tol = 1e-9
start_time = 0.0
end_time = 6
dt = 0.1
[]
[Outputs]
exodus = true
[]
(modules/combined/test/tests/inelastic_strain/elas_plas/elas_plas_nl1_cycle.i)
#
# Test for effective strain calculation.
# Boundary conditions from NAFEMS test NL1
#
#
# This is not a verification test. The boundary conditions are applied such
# that the first step generates only elastic stresses. The rest of the load
# steps generate cycles of tension and compression in the axial (i.e., y-axis)
# direction. The axial stresses and strains also cycle, however the effective
# plastic strain increases in value throughout the analysis.
#
[GlobalParams]
order = FIRST
family = LAGRANGE
volumetric_locking_correction = true
displacements = 'disp_x disp_y'
[]
[Mesh]
file = one_elem2.e
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./vonmises]
order = CONSTANT
family = MONOMIAL
[../]
[./pressure]
order = CONSTANT
family = MONOMIAL
[../]
[./elastic_strain_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./elastic_strain_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./elastic_strain_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./plastic_strain_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./plastic_strain_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./plastic_strain_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./tot_strain_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./tot_strain_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./tot_strain_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./eff_plastic_strain]
order = CONSTANT
family = MONOMIAL
[../]
[]
[Kernels]
[./TensorMechanics]
use_displaced_mesh = true
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
execute_on = timestep_end
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
execute_on = timestep_end
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
execute_on = timestep_end
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
execute_on = timestep_end
[../]
[./vonmises]
type = RankTwoScalarAux
rank_two_tensor = stress
variable = vonmises
scalar_type = VonMisesStress
execute_on = timestep_end
[../]
[./pressure]
type = RankTwoScalarAux
rank_two_tensor = stress
variable = pressure
scalar_type = Hydrostatic
execute_on = timestep_end
[../]
[./elastic_strain_xx]
type = RankTwoAux
rank_two_tensor = elastic_strain
variable = elastic_strain_xx
index_i = 0
index_j = 0
execute_on = timestep_end
[../]
[./elastic_strain_yy]
type = RankTwoAux
rank_two_tensor = elastic_strain
variable = elastic_strain_yy
index_i = 1
index_j = 1
execute_on = timestep_end
[../]
[./elastic_strain_zz]
type = RankTwoAux
rank_two_tensor = elastic_strain
variable = elastic_strain_zz
index_i = 2
index_j = 2
execute_on = timestep_end
[../]
[./plastic_strain_xx]
type = RankTwoAux
rank_two_tensor = plastic_strain
variable = plastic_strain_xx
index_i = 0
index_j = 0
execute_on = timestep_end
[../]
[./plastic_strain_yy]
type = RankTwoAux
rank_two_tensor = plastic_strain
variable = plastic_strain_yy
index_i = 1
index_j = 1
execute_on = timestep_end
[../]
[./plastic_strain_zz]
type = RankTwoAux
rank_two_tensor = plastic_strain
variable = plastic_strain_zz
index_i = 2
index_j = 2
execute_on = timestep_end
[../]
[./tot_strain_xx]
type = RankTwoAux
rank_two_tensor = total_strain
variable = tot_strain_xx
index_i = 0
index_j = 0
execute_on = timestep_end
[../]
[./tot_strain_yy]
type = RankTwoAux
rank_two_tensor = total_strain
variable = tot_strain_yy
index_i = 1
index_j = 1
execute_on = timestep_end
[../]
[./tot_strain_zz]
type = RankTwoAux
rank_two_tensor = total_strain
variable = tot_strain_zz
index_i = 2
index_j = 2
execute_on = timestep_end
[../]
[./eff_plastic_strain]
type = MaterialRealAux
property = effective_plastic_strain
variable = eff_plastic_strain
[../]
[]
[Functions]
[./appl_dispy]
type = PiecewiseLinear
x = '0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0 9.0 10.0 11.0 12.0'
y = '0.0 0.208e-4 0.50e-4 1.00e-4 0.784e-4 0.50e-4 0.0 0.216e-4 0.5e-4 1.0e-4 0.785e-4 0.50e-4 0.0'
[../]
[]
[BCs]
[./side_x]
type = DirichletBC
variable = disp_x
boundary = 101
value = 0.0
[../]
[./origin_x]
type = DirichletBC
variable = disp_x
boundary = 103
value = 0.0
[../]
[./bot_y]
type = DirichletBC
variable = disp_y
boundary = 102
value = 0.0
[../]
[./origin_y]
type = DirichletBC
variable = disp_y
boundary = 103
value = 0.0
[../]
[./top_y]
type = FunctionDirichletBC
variable = disp_y
boundary = 1
function = appl_dispy
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
block = 1
youngs_modulus = 250e9
poissons_ratio = 0.25
[../]
[./strain]
type = ComputePlaneFiniteStrain
block = 1
[../]
[./stress]
type = ComputeMultipleInelasticStress
inelastic_models = 'isoplas'
block = 1
[../]
[./isoplas]
type = IsotropicPlasticityStressUpdate
yield_stress = 5e6
hardening_constant = 0.0
relative_tolerance = 1e-20
absolute_tolerance = 1e-8
[../]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
nl_rel_tol = 1e-10
nl_abs_tol = 1e-12
l_tol = 1e-4
l_max_its = 100
nl_max_its = 20
dt = 1.0
start_time = 0.0
num_steps = 100
end_time = 12.0
[]
[Postprocessors]
[./stress_xx]
type = ElementAverageValue
variable = stress_xx
[../]
[./stress_yy]
type = ElementAverageValue
variable = stress_yy
[../]
[./stress_zz]
type = ElementAverageValue
variable = stress_zz
[../]
[./stress_xy]
type = ElementAverageValue
variable = stress_xy
[../]
[./vonmises]
type = ElementAverageValue
variable = vonmises
[../]
[./pressure]
type = ElementAverageValue
variable = pressure
[../]
[./el_strain_xx]
type = ElementAverageValue
variable = elastic_strain_xx
[../]
[./el_strain_yy]
type = ElementAverageValue
variable = elastic_strain_yy
[../]
[./el_strain_zz]
type = ElementAverageValue
variable = elastic_strain_zz
[../]
[./pl_strain_xx]
type = ElementAverageValue
variable = plastic_strain_xx
[../]
[./pl_strain_yy]
type = ElementAverageValue
variable = plastic_strain_yy
[../]
[./pl_strain_zz]
type = ElementAverageValue
variable = plastic_strain_zz
[../]
[./eff_plastic_strain]
type = ElementAverageValue
variable = eff_plastic_strain
[../]
[./tot_strain_xx]
type = ElementAverageValue
variable = tot_strain_xx
[../]
[./tot_strain_yy]
type = ElementAverageValue
variable = tot_strain_yy
[../]
[./tot_strain_zz]
type = ElementAverageValue
variable = tot_strain_zz
[../]
[./disp_x1]
type = NodalVariableValue
nodeid = 0
variable = disp_x
[../]
[./disp_x4]
type = NodalVariableValue
nodeid = 3
variable = disp_x
[../]
[./disp_y1]
type = NodalVariableValue
nodeid = 0
variable = disp_y
[../]
[./disp_y4]
type = NodalVariableValue
nodeid = 3
variable = disp_y
[../]
[./_dt]
type = TimestepSize
[../]
[]
[Outputs]
exodus = true
[./console]
type = Console
output_linear = true
[../]
[]
(modules/solid_mechanics/test/tests/combined_creep_plasticity/combined_stress_prescribed.i)
#
# 1x1x1 unit cube with time-varying pressure on top face
#
# The problem is a one-dimensional creep analysis. The top face has a
# pressure load that is a function of time. The creep strain can be
# calculated analytically. There is no practical active linear
# isotropic plasticity because the yield stress for the plasticity
# model is set to 1e30 MPa, which will not be reached in this
# regression test.
#
# The analytic solution to this problem is:
#
# d ec
# ---- = a*S^b with S = c*t^d
# dt
#
# d ec = a*c^b*t^(b*d) dt
#
# a*c^b
# ec = ----- t^(b*d+1)
# b*d+1
#
# where S = stress
# ec = creep strain
# t = time
# a = constant
# b = constant
# c = constant
# d = constant
#
# With a = 3e-24,
# b = 4,
# c = 1,
# d = 1/2, and
# t = 32400
# we have
#
# S = t^(1/2) = 180
#
# ec = 1e-24*t^3 = 3.4012224e-11
#
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
[]
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
volumetric_locking_correction = true
[]
[Physics/SolidMechanics/QuasiStatic]
[./all]
strain = FINITE
incremental = true
add_variables = true
generate_output = 'stress_yy creep_strain_yy'
[../]
[]
[Functions]
[./pressure]
type = ParsedFunction
expression = 'sqrt(t)'
[../]
[./dts]
type = PiecewiseLinear
y = '1e-2 1e-1 1e0 1e1 1e2'
x = '0 7e-1 7e0 7e1 1e2'
[../]
[]
[BCs]
[./top_pressure]
type = Pressure
variable = disp_y
boundary = top
function = pressure
[../]
[./u_bottom_fix]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[../]
[./u_yz_fix]
type = DirichletBC
variable = disp_x
boundary = left
value = 0.0
[../]
[./u_xy_fix]
type = DirichletBC
variable = disp_z
boundary = back
value = 0.0
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
youngs_modulus = 2.8e7
poissons_ratio = 0.3
[../]
[./creep_plas]
type = ComputeMultipleInelasticStress
inelastic_models = 'creep plas'
tangent_operator = elastic
[../]
[./creep]
type = PowerLawCreepStressUpdate
coefficient = 3.0e-24
n_exponent = 4
m_exponent = 0
activation_energy = 0
[../]
[./plas]
type = IsotropicPlasticityStressUpdate
hardening_constant = 1
yield_stress = 1e30
[../]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = ' lu superlu_dist'
line_search = 'none'
l_max_its = 100
nl_max_its = 100
nl_rel_tol = 1e-10
nl_abs_tol = 1e-7
l_tol = 1e-6
start_time = 0.0
end_time = 32400
dt = 1e-2
[./TimeStepper]
type = FunctionDT
function = dts
[../]
[]
[Postprocessors]
[./timestep]
type = TimestepSize
[../]
[]
[Outputs]
exodus = true
[]
(modules/combined/test/tests/combined_plasticity_temperature/ad_plasticity_temperature_dep_yield.i)
#
# This is a test of the piece-wise linear strain hardening model using the
# small strain formulation. This test exercises the temperature-dependent
# yield stress.
#
# Test procedure:
# 1. The element is pulled to and then beyond the yield stress for a given
# temperature.
# 2. The displacement is then constant while the temperature increases and
# the yield stress decreases. This results in a lower stress with more
# plastic strain.
# 3. The temperature decreases beyond its original value giving a higher
# yield stress. The displacement increases, causing increases stress to
# the new yield stress.
# 4. The temperature and yield stress are constant with increasing
# displacement giving a constant stress and more plastic strain.
#
# Plotting total_strain_yy on the x axis and stress_yy on the y axis shows
# the stress history in a clear way.
#
# s |
# t | *****
# r | *
# e | ***** *
# s | * * *
# s | * *
# |*
# +------------------
# total strain
#
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
[]
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Physics/SolidMechanics/QuasiStatic]
[./all]
strain = FINITE
incremental = true
add_variables = true
generate_output = 'stress_yy plastic_strain_xx plastic_strain_yy plastic_strain_zz'
use_automatic_differentiation = true
[../]
[]
[Variables]
[./temp]
order = FIRST
family = LAGRANGE
[../]
[]
[Functions]
[./top_pull]
type = PiecewiseLinear
x = '0 1 2 4 5 6'
y = '0 0.025 0.05 0.05 0.06 0.085'
[../]
[./yield]
type = PiecewiseLinear
x = '400 500 600'
y = '6e3 5e3 4e3'
[../]
[./temp]
type = PiecewiseLinear
x = '0 1 2 3 4'
y = '500 500 500 600 400'
[../]
[]
[Kernels]
[./heat]
type = ADHeatConduction
variable = temp
[../]
[]
[BCs]
[./y_pull_function]
type = FunctionDirichletBC
variable = disp_y
boundary = top
function = top_pull
[../]
[./x_bot]
type = DirichletBC
variable = disp_x
boundary = left
value = 0.0
[../]
[./y_bot]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[../]
[./z_bot]
type = DirichletBC
variable = disp_z
boundary = back
value = 0.0
[../]
[./temp]
type = FunctionDirichletBC
variable = temp
function = temp
boundary = left
[../]
[]
[Materials]
[./elasticity_tensor]
type = ADComputeIsotropicElasticityTensor
block = 0
youngs_modulus = 2.0e5
poissons_ratio = 0.3
[../]
[./creep_plas]
type = ADComputeMultipleInelasticStress
block = 0
inelastic_models = 'plasticity'
max_iterations = 50
absolute_tolerance = 1e-05
[../]
[./plasticity]
type = ADIsotropicPlasticityStressUpdate
block = 0
hardening_constant = 0
yield_stress_function = yield
temperature = temp
[../]
[./heat_conduction]
type = ADHeatConductionMaterial
block = 0
specific_heat = 1
thermal_conductivity = 1
[../]
[]
[Executioner]
type = Transient
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
line_search = 'none'
l_max_its = 100
nl_max_its = 100
nl_rel_tol = 1e-12
nl_abs_tol = 1e-10
l_tol = 1e-9
start_time = 0.0
end_time = 6
dt = 0.1
[]
[Outputs]
exodus = true
[]
(modules/solid_mechanics/test/tests/recompute_radial_return/isotropic_plasticity_incremental_strain.i)
# This simulation uses the piece-wise linear strain hardening model
# with the incremental small strain formulation; incremental small strain
# is required to produce the strain_increment for the DiscreteRadialReturnStressIncrement
# class, which handles the calculation of the stress increment to return
# to the yield surface in a J2 (isotropic) plasticity problem.
#
# This test assumes a Poissons ratio of zero and applies a displacement loading
# condition on the top in the y direction while fixing the displacement in the x
# and z directions; thus, only the normal stress and the normal strains in the
# y direction are compared in this problem.
#
# A similar problem was run in Abaqus on a similar 1 element mesh and was used
# to verify the SolidMechanics solution; this SolidMechanics code matches the
# SolidMechanics solution.
#
# Mechanical strain is the sum of the elastic and plastic strains but is different
# from total strain in cases with eigen strains, e.g. thermal strain.
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
[]
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Functions]
[./top_pull]
type = ParsedFunction
expression = t*(0.01)
[../]
[./hf]
type = PiecewiseLinear
x = '0 0.00004 0.0001 0.1'
y = '50 54 56 60'
[../]
[]
[Physics/SolidMechanics/QuasiStatic]
[./all]
strain = SMALL
incremental = true
add_variables = true
generate_output = 'stress_yy plastic_strain_xx plastic_strain_yy plastic_strain_zz'
[../]
[]
[BCs]
[./y_pull_function]
type = FunctionDirichletBC
variable = disp_y
boundary = top
function = top_pull
[../]
[./x_sides]
type = DirichletBC
variable = disp_x
boundary = 'left right'
value = 0.0
[../]
[./y_bot]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[../]
[./z_sides]
type = DirichletBC
variable = disp_z
boundary = 'back front'
value = 0.0
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
youngs_modulus = 2.5e5
poissons_ratio = 0.0
[../]
[./isotropic_plasticity]
type = IsotropicPlasticityStressUpdate
yield_stress = 25.
hardening_constant = 1000.0
[../]
[./radial_return_stress]
type = ComputeMultipleInelasticStress
tangent_operator = elastic
inelastic_models = 'isotropic_plasticity'
[../]
[]
[Executioner]
type = Transient
#Preconditioned JFNK (default)
solve_type = 'PJFNK'
petsc_options = '-snes_ksp_ew'
petsc_options_iname = '-ksp_gmres_restart'
petsc_options_value = '101'
line_search = 'none'
l_max_its = 20
nl_max_its = 20
nl_rel_tol = 1e-10
nl_abs_tol = 1e-12
l_tol = 1e-9
start_time = 0.0
end_time = 0.01875
dt = 0.00125
dtmin = 0.0001
[]
[Outputs]
exodus = true
print_linear_residuals = true
perf_graph = true
[]
(modules/solid_mechanics/test/tests/strain_energy_density/rate_incr_model_elas_plas.i)
# Single element test to check the strain energy density calculation
[GlobalParams]
displacements = 'disp_x disp_y'
volumetric_locking_correction = true
[]
[Mesh]
type = GeneratedMesh
dim = 2
nx = 1
ny = 1
xmin = 0
xmax = 1
ymin = 0
ymax = 2
[]
[AuxVariables]
[./SED]
order = CONSTANT
family = MONOMIAL
[../]
[]
[Functions]
[./rampConstantUp]
type = PiecewiseLinear
x = '0. 1.'
y = '0. 1.'
scale_factor = -100
[../]
[./ramp_disp_y]
type = PiecewiseLinear
x = '0. 1. 2.'
y = '0. 6.8e-6 1.36e-5'
[../]
[]
[Physics/SolidMechanics/QuasiStatic]
[./master]
strain = SMALL
add_variables = true
incremental = true
generate_output = 'stress_xx stress_yy stress_zz vonmises_stress elastic_strain_xx elastic_strain_yy elastic_strain_zz plastic_strain_xx plastic_strain_yy plastic_strain_zz strain_xx strain_yy strain_zz'
planar_formulation = PLANE_STRAIN
[../]
[]
[AuxKernels]
[./SED]
type = MaterialRealAux
variable = SED
property = strain_energy_density
execute_on = timestep_end
[../]
[]
[BCs]
[./no_x]
type = DirichletBC
variable = disp_x
preset = false
boundary = 'left'
value = 0.0
[../]
[./no_y]
type = DirichletBC
variable = disp_y
preset = false
boundary = 'bottom'
value = 0.0
[../]
[./top_disp]
type = FunctionDirichletBC
variable = disp_y
preset = false
boundary = 'top'
function = ramp_disp_y
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
youngs_modulus = 30e+6
poissons_ratio = 0.3
[../]
[./elastic_stress]
type = ComputeMultipleInelasticStress
inelastic_models = 'isoplas'
[../]
[./isoplas]
type = IsotropicPlasticityStressUpdate
yield_stress = 1e2
hardening_constant = 0.0
[../]
[./strain_energy_density]
type = StrainEnergyDensity
incremental = true
[../]
[]
[Executioner]
type = Transient
petsc_options_iname = '-ksp_gmres_restart -pc_type -pc_hypre_type -pc_hypre_boomeramg_max_iter'
petsc_options_value = '201 hypre boomeramg 4'
line_search = 'none'
l_max_its = 50
nl_max_its = 20
nl_abs_tol = 3e-7
nl_rel_tol = 1e-12
l_tol = 1e-2
start_time = 0.0
dt = 1
end_time = 2
num_steps = 2
[]
[Postprocessors]
[./epxx]
type = ElementalVariableValue
variable = elastic_strain_xx
elementid = 0
[../]
[./epyy]
type = ElementalVariableValue
variable = elastic_strain_yy
elementid = 0
[../]
[./epzz]
type = ElementalVariableValue
variable = elastic_strain_zz
elementid = 0
[../]
[./eplxx]
type = ElementalVariableValue
variable = plastic_strain_xx
elementid = 0
[../]
[./eplyy]
type = ElementalVariableValue
variable = plastic_strain_yy
elementid = 0
[../]
[./eplzz]
type = ElementalVariableValue
variable = plastic_strain_zz
elementid = 0
[../]
[./etxx]
type = ElementalVariableValue
variable = strain_xx
elementid = 0
[../]
[./etyy]
type = ElementalVariableValue
variable = strain_yy
elementid = 0
[../]
[./etzz]
type = ElementalVariableValue
variable = strain_zz
elementid = 0
[../]
[./sigxx]
type = ElementAverageValue
variable = stress_xx
[../]
[./sigyy]
type = ElementAverageValue
variable = stress_yy
[../]
[./sigzz]
type = ElementAverageValue
variable = stress_zz
[../]
[./SED]
type = ElementAverageValue
variable = SED
[../]
[]
[Outputs]
csv = true
[]
(modules/solid_mechanics/test/tests/combined_creep_plasticity/combined_creep_plasticity.i)
#
# This test is Example 2 from "A Consistent Formulation for the Integration
# of Combined Plasticity and Creep" by P. Duxbury, et al., Int J Numerical
# Methods in Engineering, Vol. 37, pp. 1277-1295, 1994.
#
# The problem is a one-dimensional bar which is loaded from yield to a value of twice
# the initial yield stress and then unloaded to return to the original stress. The
# bar must harden to the required yield stress during the load ramp, with no
# further yielding during unloading. The initial yield stress (sigma_0) is prescribed
# as 20 with a plastic strain hardening of 100. The mesh is a 1x1x1 cube with symmetry
# boundary conditions on three planes to provide a uniaxial stress field.
#
# In the PowerLawCreep model, the creep strain rate is defined by:
#
# edot = A(sigma)**n * exp(-Q/(RT)) * t**m
#
# The creep law specified in the paper, however, defines the creep strain rate as:
#
# edot = Ao * mo * (sigma)**n * t**(mo-1)
# with the creep parameters given by
# Ao = 1e-7
# mo = 0.5
# n = 5
#
# thus, input parameters for the test were specified as:
# A = Ao * mo = 1e-7 * 0.5 = 0.5e-7
# m = mo-1 = -0.5
# n = 5
# Q = 0
#
# The variation of load P with time is:
# P = 20 + 20t 0 < t < 1
# P = 40 - 40(t-1) 1 < t 1.5
#
# The analytic solution for total strain during the loading period 0 < t < 1 is:
#
# e_tot = (sigma_0 + 20*t)/E + 0.2*t + A * t**0.5 * sigma_0**n * [ 1 + (5/3)*t +
# + 2*t**2 + (10/7)*t**3 + (5/9)**t**4 + (1/11)*t**5 ]
#
# and during the unloading period 1 < t < 1.5:
#
# e_tot = (sigma_1 - 40*(t-1))/E + 0.2 + (4672/693) * A * sigma_0**n +
# A * sigma_0**n * [ t**0.5 * ( 32 - (80/3)*t + 16*t**2 - (40/7)*t**3
# + (10/9)*t**4 - (1/11)*t**5 ) - (11531/693) ]
#
# where sigma_1 is the stress at time t = 1.
#
# Assuming a Young's modulus (E) of 1000 and using the parameters defined above:
#
# e_tot(1) = 2.39734
# e_tot(1.5) = 3.16813
#
#
# The numerically computed solution is:
#
# e_tot(1) = 2.39718 (~0.006% error)
# e_tot(1.5) = 3.15555 (~0.40% error)
#
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
[]
[Physics/SolidMechanics/QuasiStatic]
[./all]
strain = FINITE
incremental = true
add_variables = true
generate_output = 'stress_yy elastic_strain_yy creep_strain_yy plastic_strain_yy'
[../]
[]
[Functions]
[./top_pull]
type = PiecewiseLinear
x = ' 0 1 1.5'
y = '-20 -40 -20'
[../]
[./dts]
type = PiecewiseLinear
x = '0 0.5 1.0 1.5'
y = '0.015 0.015 0.005 0.005'
[../]
[]
[BCs]
[./u_top_pull]
type = Pressure
variable = disp_y
boundary = top
factor = 1
function = top_pull
[../]
[./u_bottom_fix]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[../]
[./u_yz_fix]
type = DirichletBC
variable = disp_x
boundary = left
value = 0.0
[../]
[./u_xy_fix]
type = DirichletBC
variable = disp_z
boundary = back
value = 0.0
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
block = 0
youngs_modulus = 1e3
poissons_ratio = 0.3
[../]
[./creep_plas]
type = ComputeMultipleInelasticStress
block = 0
tangent_operator = elastic
inelastic_models = 'creep plas'
max_iterations = 50
absolute_tolerance = 1e-05
combined_inelastic_strain_weights = '0.0 1.0'
[../]
[./creep]
type = PowerLawCreepStressUpdate
block = 0
coefficient = 0.5e-7
n_exponent = 5
m_exponent = -0.5
activation_energy = 0
[../]
[./plas]
type = IsotropicPlasticityStressUpdate
block = 0
hardening_constant = 100
yield_stress = 20
[../]
[]
[Executioner]
type = Transient
#Preconditioned JFNK (default)
solve_type = 'PJFNK'
petsc_options = '-snes_ksp'
petsc_options_iname = '-ksp_gmres_restart'
petsc_options_value = '101'
line_search = 'none'
l_max_its = 20
nl_max_its = 6
nl_rel_tol = 1e-6
nl_abs_tol = 1e-10
l_tol = 1e-5
start_time = 0.0
end_time = 1.5
[./TimeStepper]
type = FunctionDT
function = dts
[../]
[]
[Outputs]
exodus = true
[]
(modules/combined/test/tests/inelastic_strain/elas_plas/elas_plas_nl1.i)
#
# Test for effective strain calculation.
# Boundary conditions from NAFEMS test NL1
#
# This is not a verification test. The boundary conditions are applied such
# that the first step generates only elastic stresses. The second and third
# steps generate plastic deformation and the effective strain should be
# increasing throughout the run.
#
[GlobalParams]
order = FIRST
family = LAGRANGE
volumetric_locking_correction = true
displacements = 'disp_x disp_y'
[]
[Mesh]
file = one_elem2.e
[]
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./vonmises]
order = CONSTANT
family = MONOMIAL
[../]
[./pressure]
order = CONSTANT
family = MONOMIAL
[../]
[./elastic_strain_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./elastic_strain_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./elastic_strain_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./plastic_strain_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./plastic_strain_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./plastic_strain_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./tot_strain_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./tot_strain_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./tot_strain_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./eff_plastic_strain]
order = CONSTANT
family = MONOMIAL
[../]
[]
[Kernels]
[./TensorMechanics]
use_displaced_mesh = true
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
execute_on = timestep_end
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
execute_on = timestep_end
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
execute_on = timestep_end
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
execute_on = timestep_end
[../]
[./vonmises]
type = RankTwoScalarAux
rank_two_tensor = stress
variable = vonmises
scalar_type = VonMisesStress
execute_on = timestep_end
[../]
[./pressure]
type = RankTwoScalarAux
rank_two_tensor = stress
variable = pressure
scalar_type = Hydrostatic
execute_on = timestep_end
[../]
[./elastic_strain_xx]
type = RankTwoAux
rank_two_tensor = elastic_strain
variable = elastic_strain_xx
index_i = 0
index_j = 0
execute_on = timestep_end
[../]
[./elastic_strain_yy]
type = RankTwoAux
rank_two_tensor = elastic_strain
variable = elastic_strain_yy
index_i = 1
index_j = 1
execute_on = timestep_end
[../]
[./elastic_strain_zz]
type = RankTwoAux
rank_two_tensor = elastic_strain
variable = elastic_strain_zz
index_i = 2
index_j = 2
execute_on = timestep_end
[../]
[./plastic_strain_xx]
type = RankTwoAux
rank_two_tensor = plastic_strain
variable = plastic_strain_xx
index_i = 0
index_j = 0
execute_on = timestep_end
[../]
[./plastic_strain_yy]
type = RankTwoAux
rank_two_tensor = plastic_strain
variable = plastic_strain_yy
index_i = 1
index_j = 1
execute_on = timestep_end
[../]
[./plastic_strain_zz]
type = RankTwoAux
rank_two_tensor = plastic_strain
variable = plastic_strain_zz
index_i = 2
index_j = 2
execute_on = timestep_end
[../]
[./tot_strain_xx]
type = RankTwoAux
rank_two_tensor = total_strain
variable = tot_strain_xx
index_i = 0
index_j = 0
execute_on = timestep_end
[../]
[./tot_strain_yy]
type = RankTwoAux
rank_two_tensor = total_strain
variable = tot_strain_yy
index_i = 1
index_j = 1
execute_on = timestep_end
[../]
[./tot_strain_zz]
type = RankTwoAux
rank_two_tensor = total_strain
variable = tot_strain_zz
index_i = 2
index_j = 2
execute_on = timestep_end
[../]
[./eff_plastic_strain]
type = MaterialRealAux
property = effective_plastic_strain
variable = eff_plastic_strain
[../]
[]
[Functions]
[./appl_dispy]
type = PiecewiseLinear
x = '0 1.0 2.0 3.0'
y = '0.0 0.208e-4 0.50e-4 1.00e-4'
[../]
[]
[BCs]
[./side_x]
type = DirichletBC
variable = disp_x
boundary = 101
value = 0.0
[../]
[./origin_x]
type = DirichletBC
variable = disp_x
boundary = 103
value = 0.0
[../]
[./bot_y]
type = DirichletBC
variable = disp_y
boundary = 102
value = 0.0
[../]
[./origin_y]
type = DirichletBC
variable = disp_y
boundary = 103
value = 0.0
[../]
[./top_y]
type = FunctionDirichletBC
variable = disp_y
boundary = 1
function = appl_dispy
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
block = 1
youngs_modulus = 250e9
poissons_ratio = 0.25
[../]
[./strain]
type = ComputePlaneFiniteStrain
block = 1
[../]
[./stress]
type = ComputeMultipleInelasticStress
inelastic_models = 'isoplas'
block = 1
[../]
[./isoplas]
type = IsotropicPlasticityStressUpdate
yield_stress = 5e6
hardening_constant = 0.0
[../]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
nl_rel_tol = 1e-10
nl_abs_tol = 1e-12
l_tol = 1e-4
l_max_its = 100
nl_max_its = 20
dt = 1.0
start_time = 0.0
num_steps = 100
end_time = 3.0
[] # Executioner
[Postprocessors]
[./stress_xx]
type = ElementAverageValue
variable = stress_xx
[../]
[./stress_yy]
type = ElementAverageValue
variable = stress_yy
[../]
[./stress_zz]
type = ElementAverageValue
variable = stress_zz
[../]
[./stress_xy]
type = ElementAverageValue
variable = stress_xy
[../]
[./vonmises]
type = ElementAverageValue
variable = vonmises
[../]
[./pressure]
type = ElementAverageValue
variable = pressure
[../]
[./el_strain_xx]
type = ElementAverageValue
variable = elastic_strain_xx
[../]
[./el_strain_yy]
type = ElementAverageValue
variable = elastic_strain_yy
[../]
[./el_strain_zz]
type = ElementAverageValue
variable = elastic_strain_zz
[../]
[./pl_strain_xx]
type = ElementAverageValue
variable = plastic_strain_xx
[../]
[./pl_strain_yy]
type = ElementAverageValue
variable = plastic_strain_yy
[../]
[./pl_strain_zz]
type = ElementAverageValue
variable = plastic_strain_zz
[../]
[./eff_plastic_strain]
type = ElementAverageValue
variable = eff_plastic_strain
[../]
[./tot_strain_xx]
type = ElementAverageValue
variable = tot_strain_xx
[../]
[./tot_strain_yy]
type = ElementAverageValue
variable = tot_strain_yy
[../]
[./tot_strain_zz]
type = ElementAverageValue
variable = tot_strain_zz
[../]
[./disp_x1]
type = NodalVariableValue
nodeid = 0
variable = disp_x
[../]
[./disp_x4]
type = NodalVariableValue
nodeid = 3
variable = disp_x
[../]
[./disp_y1]
type = NodalVariableValue
nodeid = 0
variable = disp_y
[../]
[./disp_y4]
type = NodalVariableValue
nodeid = 3
variable = disp_y
[../]
[./_dt]
type = TimestepSize
[../]
[]
[Outputs]
exodus = true
[./console]
type = Console
output_linear = true
[../]
[] # Outputs
(modules/solid_mechanics/test/tests/combined_creep_plasticity/combined_creep_plasticity_start_time.i)
#
# This test is Example 2 from "A Consistent Formulation for the Integration
# of Combined Plasticity and Creep" by P. Duxbury, et al., Int J Numerical
# Methods in Engineering, Vol. 37, pp. 1277-1295, 1994.
#
# The problem is a one-dimensional bar which is loaded from yield to a value of twice
# the initial yield stress and then unloaded to return to the original stress. The
# bar must harden to the required yield stress during the load ramp, with no
# further yielding during unloading. The initial yield stress (sigma_0) is prescribed
# as 20 with a plastic strain hardening of 100. The mesh is a 1x1x1 cube with symmetry
# boundary conditions on three planes to provide a uniaxial stress field.
#
# In the PowerLawCreep model, the creep strain rate is defined by:
#
# edot = A(sigma)**n * exp(-Q/(RT)) * t**m
#
# The creep law specified in the paper, however, defines the creep strain rate as:
#
# edot = Ao * mo * (sigma)**n * t**(mo-1)
# with the creep parameters given by
# Ao = 1e-7
# mo = 0.5
# n = 5
#
# thus, input parameters for the test were specified as:
# A = Ao * mo = 1e-7 * 0.5 = 0.5e-7
# m = mo-1 = -0.5
# n = 5
# Q = 0
#
# The variation of load P with time is:
# P = 20 + 20t 0 < t < 1
# P = 40 - 40(t-1) 1 < t 1.5
#
# The analytic solution for total strain during the loading period 0 < t < 1 is:
#
# e_tot = (sigma_0 + 20*t)/E + 0.2*t + A * t**0.5 * sigma_0**n * [ 1 + (5/3)*t +
# + 2*t**2 + (10/7)*t**3 + (5/9)**t**4 + (1/11)*t**5 ]
#
# and during the unloading period 1 < t < 1.5:
#
# e_tot = (sigma_1 - 40*(t-1))/E + 0.2 + (4672/693) * A * sigma_0**n +
# A * sigma_0**n * [ t**0.5 * ( 32 - (80/3)*t + 16*t**2 - (40/7)*t**3
# + (10/9)*t**4 - (1/11)*t**5 ) - (11531/693) ]
#
# where sigma_1 is the stress at time t = 1.
#
# Assuming a Young's modulus (E) of 1000 and using the parameters defined above:
#
# e_tot(1) = 2.39734
# e_tot(1.5) = 3.16813
#
#
# The numerically computed solution is:
#
# e_tot(1) = 2.39718 (~0.006% error)
# e_tot(1.5) = 3.15555 (~0.40% error)
#
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
[]
[Physics/SolidMechanics/QuasiStatic]
[./all]
strain = FINITE
incremental = true
add_variables = true
generate_output = 'stress_yy elastic_strain_yy creep_strain_yy plastic_strain_yy'
[../]
[]
[Functions]
[./top_pull]
type = PiecewiseLinear
x = ' 10 11 11.5'
y = '-20 -40 -20'
[../]
[./dts]
type = PiecewiseLinear
x = '10 10.5 11.0 11.5'
y = '0.015 0.015 0.005 0.005'
[../]
[]
[BCs]
[./u_top_pull]
type = Pressure
variable = disp_y
boundary = top
factor = 1
function = top_pull
[../]
[./u_bottom_fix]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[../]
[./u_yz_fix]
type = DirichletBC
variable = disp_x
boundary = left
value = 0.0
[../]
[./u_xy_fix]
type = DirichletBC
variable = disp_z
boundary = back
value = 0.0
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
block = 0
youngs_modulus = 1e3
poissons_ratio = 0.3
[../]
[./creep_plas]
type = ComputeMultipleInelasticStress
block = 0
tangent_operator = elastic
inelastic_models = 'creep plas'
max_iterations = 50
absolute_tolerance = 1e-05
combined_inelastic_strain_weights = '0.0 1.0'
[../]
[./creep]
type = PowerLawCreepStressUpdate
block = 0
coefficient = 0.5e-7
n_exponent = 5
m_exponent = -0.5
activation_energy = 0
start_time = 10
[../]
[./plas]
type = IsotropicPlasticityStressUpdate
block = 0
hardening_constant = 100
yield_stress = 20
[../]
[]
[Executioner]
type = Transient
#Preconditioned JFNK (default)
solve_type = 'PJFNK'
petsc_options = '-snes_ksp'
petsc_options_iname = '-ksp_gmres_restart'
petsc_options_value = '101'
line_search = 'none'
l_max_its = 20
nl_max_its = 6
nl_rel_tol = 1e-6
nl_abs_tol = 1e-10
l_tol = 1e-5
start_time = 10.0
end_time = 11.5
[./TimeStepper]
type = FunctionDT
function = dts
[../]
[]
[Outputs]
exodus = true
[]
(modules/solid_mechanics/test/tests/power_law_creep/composite_power_law_creep_plasticity.i)
# 1x1x1 unit cube with uniform pressure on top face and 2 phases with different materials
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 6
zmax = 1
xmax = 1
ymax = 1
[]
[Variables]
[temp]
order = FIRST
family = LAGRANGE
initial_condition = 1000.0
[]
[]
[ICs]
[phase1IC]
type = BoundingBoxIC
x1 = -1
x2 = 1.5
y1 = -1
y2 = 1.5
z1 = -1
z2 = 1.0
inside = 1
outside = 0
variable = phase1
int_width=0.01
[]
[phase2IC]
type = BoundingBoxIC
x1 = -1
x2 = 1.5
y1 = -1
y2 = 1.5
z1 = -1
z2 = 1.0
inside = 0
outside = 1
variable = phase2
int_width=0.01
[]
[]
[AuxVariables]
[phase1]
[]
[phase2]
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = FINITE
incremental = true
add_variables = true
generate_output = 'stress_yy creep_strain_xx creep_strain_yy creep_strain_zz elastic_strain_yy'
[]
[]
[Functions]
[top_pull]
type = PiecewiseLinear
x = '0 1'
y = '1 1'
[]
[]
[Kernels]
[heat]
type = Diffusion
variable = temp
[]
[heat_ie]
type = TimeDerivative
variable = temp
[]
[]
[BCs]
[u_top_pull]
type = Pressure
variable = disp_y
boundary = top
factor = -10.0e6
function = top_pull
[]
[u_bottom_fix]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[]
[u_yz_fix]
type = DirichletBC
variable = disp_x
boundary = left
value = 0.0
[]
[u_xy_fix]
type = DirichletBC
variable = disp_z
boundary = back
value = 0.0
[]
[temp_fix]
type = DirichletBC
variable = temp
boundary = 'bottom top'
value = 1000.0
[]
[]
[Materials]
[elasticity_tensor1]
type = ComputeIsotropicElasticityTensor
base_name = C1
youngs_modulus = 2e11
poissons_ratio = 0.3
[]
[elasticity_tensor2]
type = ComputeIsotropicElasticityTensor
base_name = C2
youngs_modulus = 2e11
poissons_ratio = 0.3
[]
[h1]
type = ParsedMaterial
property_name = h1
coupled_variables = phase1
expression = '0.5*tanh(20*(phase1-0.5))+0.5'
[]
[h2]
type = ParsedMaterial
property_name = h2
coupled_variables = phase2
expression = '0.5*tanh(20*(phase2-0.5))+0.5'
[]
[./C]
type = CompositeElasticityTensor
coupled_variables = 'phase1 phase2'
tensors = 'C1 C2'
weights = 'h1 h2'
[../]
[radial_return_stress]
type = ComputeMultipleInelasticStress
inelastic_models = 'power_law_creep plas'
tangent_operator = elastic
[]
[power_law_creep]
type = CompositePowerLawCreepStressUpdate
coefficient = '1.0e-15 2.0e-18'
n_exponent = '4 5'
activation_energy = '3.0e5 3.5e5'
switching_functions = 'h1 h2'
temperature = temp
[]
[./plas]
type = IsotropicPlasticityStressUpdate
hardening_constant = 1
yield_stress = 1e30
[../]
[]
[VectorPostprocessors]
[./soln]
type = LineValueSampler
warn_discontinuous_face_values = false
sort_by = x
variable = 'disp_x disp_y disp_z creep_strain_xx creep_strain_yy creep_strain_zz'
start_point = '0 0 0.0'
end_point = '1.0 1.0 1.0'
num_points = 5
outputs = tests
[../]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options = '-snes_ksp'
petsc_options_iname = '-ksp_gmres_restart'
petsc_options_value = '101'
line_search = 'none'
l_max_its = 20
nl_max_its = 20
nl_rel_tol = 1.0e-9
nl_abs_tol = 1.0e-9
l_tol = 1e-10
start_time = 0.0
end_time = 1.0
num_steps = 10
dt = 0.1
[]
[Outputs]
exodus = false
[./tests]
type = CSV
execute_on = final
[../]
[]
(modules/solid_mechanics/test/tests/power_law_creep/composite_power_law_creep_onePhaseMulti.i)
# 1x1x1 unit cube with uniform pressure on top face and 2 phases with different materials
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 6
zmax = 1
xmax = 1
ymax = 1
[]
[Variables]
[temp]
order = FIRST
family = LAGRANGE
initial_condition = 1000.0
[]
[]
[ICs]
[phase1IC]
type = BoundingBoxIC
x1 = -1
x2 = 1.5
y1 = -1
y2 = 1.5
z1 = -1
z2 = 1.0
inside = 1
outside = 0
variable = phase1
int_width=0.01
[]
[phase2IC]
type = BoundingBoxIC
x1 = -1
x2 = 1.5
y1 = -1
y2 = 1.5
z1 = -1
z2 = 1.0
inside = 0
outside = 1
variable = phase2
int_width=0.01
[]
[]
[AuxVariables]
[phase1]
[]
[phase2]
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = FINITE
incremental = true
add_variables = true
generate_output = 'stress_yy creep_strain_xx creep_strain_yy creep_strain_zz elastic_strain_yy'
[]
[]
[Functions]
[top_pull]
type = PiecewiseLinear
x = '0 1'
y = '1 1'
[]
[]
[Kernels]
[heat]
type = Diffusion
variable = temp
[]
[heat_ie]
type = TimeDerivative
variable = temp
[]
[]
[BCs]
[u_top_pull]
type = Pressure
variable = disp_y
boundary = top
factor = -10.0e6
function = top_pull
[]
[u_bottom_fix]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[]
[u_yz_fix]
type = DirichletBC
variable = disp_x
boundary = left
value = 0.0
[]
[u_xy_fix]
type = DirichletBC
variable = disp_z
boundary = back
value = 0.0
[]
[temp_fix]
type = DirichletBC
variable = temp
boundary = 'bottom top'
value = 1000.0
[]
[]
[Materials]
[elasticity_tensor1]
type = ComputeIsotropicElasticityTensor
base_name = C1
youngs_modulus = 2e11
poissons_ratio = 0.3
[]
[elasticity_tensor2]
type = ComputeIsotropicElasticityTensor
base_name = C2
youngs_modulus = 2e11
poissons_ratio = 0.3
[]
[h1]
type = ParsedMaterial
property_name = h1
coupled_variables = phase1
expression = '0.5*tanh(20*(phase1-0.5))+0.5'
[]
[h2]
type = ParsedMaterial
property_name = h2
coupled_variables = phase2
expression = '0.5*tanh(20*(phase2-0.5))+0.5'
[]
[./C]
type = CompositeElasticityTensor
coupled_variables = 'phase1 phase2'
tensors = 'C1 C2'
weights = 'h1 h2'
[../]
[radial_return_stress]
type = ComputeMultipleInelasticStress
inelastic_models = 'power_law_creep plas'
tangent_operator = elastic
[]
[power_law_creep]
type = CompositePowerLawCreepStressUpdate
coefficient = '1e-15 5e-20'
n_exponent = '4 5'
activation_energy = '3.0e5 3.5e5'
switching_functions = 'h1 h1'
temperature = temp
[]
[./plas]
type = IsotropicPlasticityStressUpdate
hardening_constant = 1
yield_stress = 1e30
[../]
[]
[VectorPostprocessors]
[./soln]
type = LineValueSampler
warn_discontinuous_face_values = false
sort_by = x
variable = 'disp_x disp_y disp_z creep_strain_xx creep_strain_yy creep_strain_zz'
start_point = '0 0 0.0'
end_point = '1.0 1.0 1.0'
num_points = 5
outputs = tests
[../]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options = '-snes_ksp'
petsc_options_iname = '-ksp_gmres_restart'
petsc_options_value = '101'
line_search = 'none'
l_max_its = 20
nl_max_its = 20
nl_rel_tol = 1.0e-10
nl_abs_tol = 1.0e-10
l_tol = 1e-10
start_time = 0.0
end_time = 1.0
num_steps = 10
dt = 0.1
[]
[Outputs]
exodus = false
[./tests]
type = CSV
execute_on = final
[../]
[]
(modules/solid_mechanics/test/tests/strain_energy_density/incr_model_elas_plas.i)
# Single element test to check the strain energy density calculation
[GlobalParams]
displacements = 'disp_x disp_y'
volumetric_locking_correction = true
[]
[Mesh]
type = GeneratedMesh
dim = 2
nx = 1
ny = 1
xmin = 0
xmax = 1
ymin = 0
ymax = 2
[]
[AuxVariables]
[./SED]
order = CONSTANT
family = MONOMIAL
[../]
[]
[Functions]
[./rampConstantUp]
type = PiecewiseLinear
x = '0. 1.'
y = '0. 1.'
scale_factor = -100
[../]
[./ramp_disp_y]
type = PiecewiseLinear
x = '0. 1. 2.'
y = '0. 6.8e-6 1.36e-5'
[../]
[]
[Physics/SolidMechanics/QuasiStatic]
[./master]
strain = SMALL
add_variables = true
incremental = true
generate_output = 'stress_xx stress_yy stress_zz vonmises_stress elastic_strain_xx elastic_strain_yy elastic_strain_zz plastic_strain_xx plastic_strain_yy plastic_strain_zz strain_xx strain_yy strain_zz'
planar_formulation = PLANE_STRAIN
[../]
[]
[AuxKernels]
[./SED]
type = MaterialRealAux
variable = SED
property = strain_energy_density
execute_on = timestep_end
[../]
[]
[BCs]
[./no_x]
type = DirichletBC
variable = disp_x
preset = false
boundary = 'left'
value = 0.0
[../]
[./no_y]
type = DirichletBC
variable = disp_y
preset = false
boundary = 'bottom'
value = 0.0
[../]
[./top_disp]
type = FunctionDirichletBC
variable = disp_y
preset = false
boundary = 'top'
function = ramp_disp_y
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
youngs_modulus = 30e+6
poissons_ratio = 0.3
[../]
[./elastic_stress]
type = ComputeMultipleInelasticStress
inelastic_models = 'isoplas'
[../]
[./isoplas]
type = IsotropicPlasticityStressUpdate
yield_stress = 1e2
hardening_constant = 0.0
[../]
[./strain_energy_density]
type = StrainEnergyDensity
incremental = true
[../]
[]
[Executioner]
type = Transient
petsc_options_iname = '-ksp_gmres_restart -pc_type -pc_hypre_type -pc_hypre_boomeramg_max_iter'
petsc_options_value = '201 hypre boomeramg 4'
line_search = 'none'
l_max_its = 50
nl_max_its = 20
nl_abs_tol = 3e-7
nl_rel_tol = 1e-12
l_tol = 1e-2
start_time = 0.0
dt = 1
end_time = 2
num_steps = 2
[]
[Postprocessors]
[./epxx]
type = ElementalVariableValue
variable = elastic_strain_xx
elementid = 0
[../]
[./epyy]
type = ElementalVariableValue
variable = elastic_strain_yy
elementid = 0
[../]
[./epzz]
type = ElementalVariableValue
variable = elastic_strain_zz
elementid = 0
[../]
[./eplxx]
type = ElementalVariableValue
variable = plastic_strain_xx
elementid = 0
[../]
[./eplyy]
type = ElementalVariableValue
variable = plastic_strain_yy
elementid = 0
[../]
[./eplzz]
type = ElementalVariableValue
variable = plastic_strain_zz
elementid = 0
[../]
[./etxx]
type = ElementalVariableValue
variable = strain_xx
elementid = 0
[../]
[./etyy]
type = ElementalVariableValue
variable = strain_yy
elementid = 0
[../]
[./etzz]
type = ElementalVariableValue
variable = strain_zz
elementid = 0
[../]
[./sigxx]
type = ElementAverageValue
variable = stress_xx
[../]
[./sigyy]
type = ElementAverageValue
variable = stress_yy
[../]
[./sigzz]
type = ElementAverageValue
variable = stress_zz
[../]
[./SED]
type = ElementAverageValue
variable = SED
[../]
[]
[Outputs]
csv = true
[]
(modules/solid_mechanics/test/tests/recompute_radial_return/cp_affine_plasticity.i)
# Affine Plasticity Test for Transient Stress Eigenvalues with Stationary Eigenvectors
# This test is taken from K. Jamojjala, R. Brannon, A. Sadeghirad, J. Guilkey,
# "Verification tests in solid mechanics," Engineering with Computers, Vol 31.,
# p. 193-213.
# The test involves applying particular strains and expecting particular stresses.
# The material properties are:
# Yield in shear 165 MPa
# Shear modulus 79 GPa
# Poisson's ratio 1/3
# The strains are:
# Time e11 e22 e33
# 0 0 0 0
# 1 -0.003 -0.003 0.006
# 2 -0.0103923 0 0.0103923
# The expected stresses are:
# sigma11:
# -474*t 0 < t <= 0.201
# -95.26 0.201 < t <= 1
# (189.4+0.1704*sqrt(a)-0.003242*a)
# --------------------------------- 1 < t <= 2
# 1+0.00001712*a
# -189.4 t > 2 (paper erroneously gives a positive value)
#
# sigma22:
# -474*t 0 < t <= 0.201
# -95.26 0.201 < t <= 1
# -(76.87+1.443*sqrt(a)-0.001316*a)
# --------------------------------- 1 < t <= 2 (paper gives opposite sign)
# 1+0.00001712*a
# 76.87 t > 2
#
# sigma33:
# 948*t 0 < t <= 0.201
# 190.5 0.201 < t <= 1
# -(112.5-1.272*sqrt(a)-0.001926*a)
# --------------------------------- 1 < t <= 2 (paper has two sign errors here)
# 1+0.00001712*a
# 112.5 t > 2
#
# where a = exp(12.33*t).
#
# Note: If planning to run this case with strain type ComputeFiniteStrain, the
# displacement function must be adjusted. Instead of
# strain = (l - l0)/l0 = (u+l0 - l0)/l0 = u/l0
# with l0=1.0, we would have
# strain = log(l/l0) = log((u+l0)/l0)
# with l0=1.0. So, for strain = -0.003,
# -0.003 = log((u+l0)/l0) ->
# u = exp(-0.003)*l0 - l0 = -0.0029955044966269995.
#
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
block = '0'
[]
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
# This test uses ElementalVariableValue postprocessors on specific
# elements, so element numbering needs to stay unchanged
allow_renumbering = false
[]
[Functions]
[disp_x]
type = PiecewiseLinear
x = '0. 1. 2.'
y = '0. -0.003 -0.0103923'
[]
[disp_y]
type = PiecewiseLinear
x = '0. 1. 2.'
y = '0. -0.003 0.'
[]
[disp_z]
type = PiecewiseLinear
x = '0. 1. 2.'
y = '0. 0.006 0.0103923'
[]
[stress_xx]
type = ParsedFunction
# The paper gives 0.201 as the time at initial yield, but 0.20097635952803425 is the exact value.
# The paper gives -95.26 MPa as the stress at yield, but -95.26279441628823 is the exact value.
# The paper gives 12.33 as the factor in the exponential, but 12.332921390339125 is the exact value.
# 189.409039923814000, 0.170423791206825, -0.003242011311945, 1.711645501845780E-05 - exact values
symbol_names = 'timeAtYield stressAtYield expFac a b c d'
symbol_values = '0.20097635952803425 -95.26279441628823 12.332921390339125 189.409039923814000 0.170423791206825 -0.003242011311945 1.711645501845780E-05'
value = '1e6*
if(t<=timeAtYield, -474*t,
if(t<=1, stressAtYield,
(a+b*sqrt(exp(expFac*t))+c*exp(expFac*t))/(1.0+d*exp(expFac*t))))' # tends to -a
[]
[stress_yy]
type = ParsedFunction
# The paper gives 0.201 as the time at initial yield, but 0.20097635952803425 is the exact value.
# the paper gives -95.26 MPa as the stress at yield, but -95.26279441628823 is the exact value.
# The paper gives 12.33 as the factor in the exponential, but 12.332921390339125 is the exact value.
# -76.867432297315000, -1.442488120272900, 0.001315697947301, 1.711645501845780E-05 - exact values
symbol_names = 'timeAtYield stressAtYield expFac a b c d'
symbol_values = '0.20097635952803425 -95.26279441628823 12.332921390339125 -76.867432297315000 -1.442488120272900 0.001315697947301 1.711645501845780E-05'
value = '1e6*
if(t<=timeAtYield, -474*t,
if(t<=1, stressAtYield,
(a+b*sqrt(exp(expFac*t))+c*exp(expFac*t))/(1.0+d*exp(expFac*t))))' # tends to -a
[]
[stress_zz]
type = ParsedFunction
# The paper gives 0.201 as the time at initial yield, but 0.20097635952803425 is the exact value.
# the paper gives 190.5 MPa as the stress at yield, but 190.52558883257645 is the exact value.
# The paper gives 12.33 as the factor in the exponential, but 12.332921390339125 is the exact value.
# -112.541607626499000, 1.272064329066080, 0.001926313364644, 1.711645501845780E-05 - exact values
symbol_names = 'timeAtYield stressAtYield expFac a b c d'
symbol_values = '0.20097635952803425 190.52558883257645 12.332921390339125 -112.541607626499000 1.272064329066080 0.001926313364644 1.711645501845780E-05'
value = '1e6*
if(t<=timeAtYield, 948*t,
if(t<=1, stressAtYield,
(a+b*sqrt(exp(expFac*t))+c*exp(expFac*t))/(1.0+d*exp(expFac*t))))' # tends to -a
[]
[]
[Variables]
[disp_x]
order = FIRST
family = LAGRANGE
[]
[disp_y]
order = FIRST
family = LAGRANGE
[]
[disp_z]
order = FIRST
family = LAGRANGE
[]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./vonmises]
order = CONSTANT
family = MONOMIAL
[../]
[./plastic_strain_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./plastic_strain_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./plastic_strain_zz]
order = CONSTANT
family = MONOMIAL
[../]
[]
[Kernels]
[SolidMechanics]
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
execute_on = 'timestep_end'
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
execute_on = 'timestep_end'
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
execute_on = 'timestep_end'
[../]
[./vonmises]
type = RankTwoScalarAux
rank_two_tensor = stress
variable = vonmises
scalar_type = vonmisesStress
execute_on = 'timestep_end'
[../]
[./plastic_strain_xx]
type = RankTwoAux
rank_two_tensor = plastic_strain
variable = plastic_strain_xx
index_i = 0
index_j = 0
execute_on = 'timestep_end'
[../]
[./plastic_strain_yy]
type = RankTwoAux
rank_two_tensor = plastic_strain
variable = plastic_strain_yy
index_i = 1
index_j = 1
execute_on = 'timestep_end'
[../]
[./plastic_strain_zz]
type = RankTwoAux
rank_two_tensor = plastic_strain
variable = plastic_strain_zz
index_i = 2
index_j = 2
execute_on = 'timestep_end'
[../]
[]
[BCs]
[fixed_x]
type = DirichletBC
variable = disp_x
boundary = left
value = 0.0
[]
[fixed_y]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[]
[fixed_z]
type = DirichletBC
variable = disp_z
boundary = back
value = 0.0
[]
[disp_x]
type = FunctionDirichletBC
variable = disp_x
boundary = right
function = disp_x
[]
[disp_y]
type = FunctionDirichletBC
variable = disp_y
boundary = top
function = disp_y
[]
[disp_z]
type = FunctionDirichletBC
variable = disp_z
boundary = front
function = disp_z
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeIsotropicElasticityTensor
youngs_modulus = 210666666666.666667
poissons_ratio = 0.3333333333333333
[../]
[./strain]
type = ComputeIncrementalStrain
[../]
[creep]
type = PowerLawCreepStressUpdate
coefficient = 0
n_exponent = 1
m_exponent = 1
activation_energy = 0
temperature = 1
[]
[isotropic_plasticity]
type = IsotropicPlasticityStressUpdate
yield_stress = 285788383.2488647 # = sqrt(3)*165e6 = sqrt(3) * yield in shear
hardening_constant = 0.0
[]
[radial_return_stress]
type = ComputeCreepPlasticityStress
tangent_operator = elastic
creep_model = creep
plasticity_model = isotropic_plasticity
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options = '-snes_ksp_ew'
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
nl_abs_tol = 1e-10
l_max_its = 20
start_time = 0.0
dt = 0.01 # use 0.0001 for a nearly exact match
end_time = 2.0
[]
[Postprocessors]
[analytic_xx]
type = FunctionValuePostprocessor
function = stress_xx
[]
[analytic_yy]
type = FunctionValuePostprocessor
function = stress_yy
[]
[analytic_zz]
type = FunctionValuePostprocessor
function = stress_zz
[]
[stress_xx]
type = ElementalVariableValue
variable = stress_xx
elementid = 0
[]
[stress_yy]
type = ElementalVariableValue
variable = stress_yy
elementid = 0
[]
[stress_zz]
type = ElementalVariableValue
variable = stress_zz
elementid = 0
[]
[stress_xx_l2_error]
type = ElementL2Error
variable = stress_xx
function = stress_xx
[]
[stress_yy_l2_error]
type = ElementL2Error
variable = stress_yy
function = stress_yy
[]
[stress_zz_l2_error]
type = ElementL2Error
variable = stress_zz
function = stress_zz
[]
[]
[Outputs]
exodus = true
[]
(modules/solid_mechanics/test/tests/combined_creep_plasticity/creepWithPlasticity.i)
#
# This test is Example 2 from "A Consistent Formulation for the Integration
# of Combined Plasticity and Creep" by P. Duxbury, et al., Int J Numerical
# Methods in Engineering, Vol. 37, pp. 1277-1295, 1994.
#
# The problem is a one-dimensional bar which is loaded from yield to a value of twice
# the initial yield stress and then unloaded to return to the original stress. The
# bar must harden to the required yield stress during the load ramp, with no
# further yielding during unloading. The initial yield stress (sigma_0) is prescribed
# as 20 with a plastic strain hardening of 100. The mesh is a 1x1x1 cube with symmetry
# boundary conditions on three planes to provide a uniaxial stress field.
#
# In the PowerLawCreep model, the creep strain rate is defined by:
#
# edot = A(sigma)**n * exp(-Q/(RT)) * t**m
#
# The creep law specified in the paper, however, defines the creep strain rate as:
#
# edot = Ao * mo * (sigma)**n * t**(mo-1)
# with the creep parameters given by
# Ao = 1e-7
# mo = 0.5
# n = 5
#
# thus, input parameters for the test were specified as:
# A = Ao * mo = 1e-7 * 0.5 = 0.5e-7
# m = mo-1 = -0.5
# n = 5
# Q = 0
#
# The variation of load P with time is:
# P = 20 + 20t 0 < t < 1
# P = 40 - 40(t-1) 1 < t 1.5
#
# The analytic solution for total strain during the loading period 0 < t < 1 is:
#
# e_tot = (sigma_0 + 20*t)/E + 0.2*t + A * t**0.5 * sigma_0**n * [ 1 + (5/3)*t +
# + 2*t**2 + (10/7)*t**3 + (5/9)**t**4 + (1/11)*t**5 ]
#
# and during the unloading period 1 < t < 1.5:
#
# e_tot = (sigma_1 - 40*(t-1))/E + 0.2 + (4672/693) * A * sigma_0**n +
# A * sigma_0**n * [ t**0.5 * ( 32 - (80/3)*t + 16*t**2 - (40/7)*t**3
# + (10/9)*t**4 - (1/11)*t**5 ) - (11531/693) ]
#
# where sigma_1 is the stress at time t = 1.
#
# Assuming a Young's modulus (E) of 1000 and using the parameters defined above:
#
# e_tot(1) = 2.39734
# e_tot(1.5) = 3.16813
#
#
# The numerically computed solution is:
#
# e_tot(1) = 2.39718 (~0.006% error)
# e_tot(1.5) = 3.15555 (~0.40% error)
#
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = FINITE
incremental = true
add_variables = true
generate_output = 'stress_yy elastic_strain_yy creep_strain_yy plastic_strain_yy'
[]
[]
[Functions]
[top_pull]
type = PiecewiseLinear
x = ' 0 1 1.5'
y = '-20 -40 -20'
[]
[dts]
type = PiecewiseLinear
x = '0 0.5 1.0 1.5'
y = '0.015 0.015 0.005 0.005'
[]
[]
[BCs]
[u_top_pull]
type = Pressure
variable = disp_y
boundary = top
factor = 1
function = top_pull
[]
[u_bottom_fix]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[]
[u_yz_fix]
type = DirichletBC
variable = disp_x
boundary = left
value = 0.0
[]
[u_xy_fix]
type = DirichletBC
variable = disp_z
boundary = back
value = 0.0
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeIsotropicElasticityTensor
block = 0
youngs_modulus = 1e3
poissons_ratio = 0.3
[]
[creep_plas]
type = ComputeCreepPlasticityStress
block = 0
tangent_operator = elastic
creep_model = creep
plasticity_model = plasticity
max_iterations = 50
relative_tolerance = 1e-8
absolute_tolerance = 1e-8
[]
[creep]
type = PowerLawCreepStressUpdate
block = 0
coefficient = 0.5e-7
n_exponent = 5
m_exponent = -0.5
activation_energy = 0
temperature = 1
[]
[plasticity]
type = IsotropicPlasticityStressUpdate
block = 0
yield_stress = 20
hardening_constant = 100
[]
[]
[Executioner]
type = Transient
#Preconditioned JFNK (default)
solve_type = 'PJFNK'
petsc_options = '-snes_ksp'
petsc_options_iname = '-ksp_gmres_restart'
petsc_options_value = '101'
line_search = 'none'
l_max_its = 20
nl_max_its = 6
nl_rel_tol = 1e-6
nl_abs_tol = 1e-10
l_tol = 1e-5
start_time = 0.0
end_time = 1.5
[TimeStepper]
type = FunctionDT
function = dts
[]
[]
[Outputs]
exodus = true
[]
(modules/combined/test/tests/phase_field_fracture/crack2d_computeCrackedStress_finitestrain_plastic.i)
#This input uses PhaseField-Nonconserved Action to add phase field fracture bulk rate kernels
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
nx = 40
ny = 20
ymax = 0.5
[]
[./noncrack]
type = BoundingBoxNodeSetGenerator
new_boundary = noncrack
bottom_left = '0.5 0 0'
top_right = '1 0 0'
input = gen
[../]
[]
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[AuxVariables]
[./strain_yy]
family = MONOMIAL
order = CONSTANT
[../]
[./elastic_strain_yy]
family = MONOMIAL
order = CONSTANT
[../]
[./plastic_strain_yy]
family = MONOMIAL
order = CONSTANT
[../]
[./uncracked_stress_yy]
family = MONOMIAL
order = CONSTANT
[../]
[]
[Physics]
[./SolidMechanics]
[./QuasiStatic]
[./All]
add_variables = true
strain = FINITE
planar_formulation = PLANE_STRAIN
additional_generate_output = 'stress_yy vonmises_stress'
strain_base_name = uncracked
[../]
[../]
[../]
[]
[Modules]
[./PhaseField]
[./Nonconserved]
[./c]
free_energy = E_el
kappa = kappa_op
mobility = L
[../]
[../]
[../]
[]
[Kernels]
[./solid_x]
type = PhaseFieldFractureMechanicsOffDiag
variable = disp_x
component = 0
c = c
[../]
[./solid_y]
type = PhaseFieldFractureMechanicsOffDiag
variable = disp_y
component = 1
c = c
[../]
[./off_disp]
type = AllenCahnElasticEnergyOffDiag
variable = c
displacements = 'disp_x disp_y'
mob_name = L
[../]
[]
[AuxKernels]
[./strain_yy]
type = RankTwoAux
variable = strain_yy
rank_two_tensor = uncracked_mechanical_strain
index_i = 1
index_j = 1
execute_on = TIMESTEP_END
[../]
[./elastic_strain_yy]
type = RankTwoAux
variable = elastic_strain_yy
rank_two_tensor = uncracked_elastic_strain
index_i = 1
index_j = 1
execute_on = TIMESTEP_END
[../]
[./plastic_strain_yy]
type = RankTwoAux
variable = plastic_strain_yy
rank_two_tensor = uncracked_plastic_strain
index_i = 1
index_j = 1
execute_on = TIMESTEP_END
[../]
[./uncracked_stress_yy]
type = RankTwoAux
variable = uncracked_stress_yy
rank_two_tensor = uncracked_stress
index_i = 1
index_j = 1
execute_on = TIMESTEP_END
[../]
[]
[BCs]
[./ydisp]
type = FunctionDirichletBC
variable = disp_y
boundary = top
function = 't'
[../]
[./yfix]
type = DirichletBC
variable = disp_y
boundary = noncrack
value = 0
[../]
[./xfix]
type = DirichletBC
variable = disp_x
boundary = right
value = 0
[../]
[]
[Functions]
[./hf]
type = PiecewiseLinear
x = '0 0.001 0.003 0.023'
y = '0.85 1.0 1.25 1.5'
[../]
[]
[Materials]
[./pfbulkmat]
type = GenericConstantMaterial
prop_names = 'gc_prop l visco'
prop_values = '1e-3 0.05 5e-3'
[../]
[./elasticity_tensor]
type = ComputeElasticityTensor
C_ijkl = '120.0 80.0'
fill_method = symmetric_isotropic
base_name = uncracked
[../]
[./isotropic_plasticity]
type = IsotropicPlasticityStressUpdate
yield_stress = 0.85
hardening_function = hf
base_name = uncracked
[../]
[./radial_return_stress]
type = ComputeMultipleInelasticStress
tangent_operator = elastic
inelastic_models = 'isotropic_plasticity'
base_name = uncracked
[../]
[./cracked_stress]
type = ComputeCrackedStress
c = c
F_name = E_el
use_current_history_variable = true
uncracked_base_name = uncracked
finite_strain_model = true
[../]
[]
[Postprocessors]
[./av_stress_yy]
type = ElementAverageValue
variable = stress_yy
[../]
[./av_strain_yy]
type = SideAverageValue
variable = disp_y
boundary = top
[../]
[./av_uncracked_stress_yy]
type = ElementAverageValue
variable = uncracked_stress_yy
[../]
[./max_c]
type = ElementExtremeValue
variable = c
[../]
[]
[Preconditioning]
[./smp]
type = SMP
full = true
[../]
[]
[Executioner]
type = Transient
solve_type = PJFNK
petsc_options_iname = '-pc_type -pc_factor_mat_solving_package'
petsc_options_value = 'lu superlu_dist'
nl_rel_tol = 1e-8
l_tol = 1e-4
l_max_its = 100
nl_max_its = 10
dt = 2.0e-5
num_steps = 2
[]
[Outputs]
exodus = true
[]
(modules/solid_mechanics/test/tests/recompute_radial_return/affine_plasticity.i)
# Affine Plasticity Test for Transient Stress Eigenvalues with Stationary Eigenvectors
# This test is taken from K. Jamojjala, R. Brannon, A. Sadeghirad, J. Guilkey,
# "Verification tests in solid mechanics," Engineering with Computers, Vol 31.,
# p. 193-213.
# The test involves applying particular strains and expecting particular stresses.
# The material properties are:
# Yield in shear 165 MPa
# Shear modulus 79 GPa
# Poisson's ratio 1/3
# The strains are:
# Time e11 e22 e33
# 0 0 0 0
# 1 -0.003 -0.003 0.006
# 2 -0.0103923 0 0.0103923
# The expected stresses are:
# sigma11:
# -474*t 0 < t <= 0.201
# -95.26 0.201 < t <= 1
# (189.4+0.1704*sqrt(a)-0.003242*a)
# --------------------------------- 1 < t <= 2
# 1+0.00001712*a
# -189.4 t > 2 (paper erroneously gives a positive value)
#
# sigma22:
# -474*t 0 < t <= 0.201
# -95.26 0.201 < t <= 1
# -(76.87+1.443*sqrt(a)-0.001316*a)
# --------------------------------- 1 < t <= 2 (paper gives opposite sign)
# 1+0.00001712*a
# 76.87 t > 2
#
# sigma33:
# 948*t 0 < t <= 0.201
# 190.5 0.201 < t <= 1
# -(112.5-1.272*sqrt(a)-0.001926*a)
# --------------------------------- 1 < t <= 2 (paper has two sign errors here)
# 1+0.00001712*a
# 112.5 t > 2
#
# where a = exp(12.33*t).
#
# Note: If planning to run this case with strain type ComputeFiniteStrain, the
# displacement function must be adjusted. Instead of
# strain = (l - l0)/l0 = (u+l0 - l0)/l0 = u/l0
# with l0=1.0, we would have
# strain = log(l/l0) = log((u+l0)/l0)
# with l0=1.0. So, for strain = -0.003,
# -0.003 = log((u+l0)/l0) ->
# u = exp(-0.003)*l0 - l0 = -0.0029955044966269995.
#
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
block = '0'
[]
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
# This test uses ElementalVariableValue postprocessors on specific
# elements, so element numbering needs to stay unchanged
allow_renumbering = false
[]
[Functions]
[./disp_x]
type = PiecewiseLinear
x = '0. 1. 2.'
y = '0. -0.003 -0.0103923'
[../]
[./disp_y]
type = PiecewiseLinear
x = '0. 1. 2.'
y = '0. -0.003 0.'
[../]
[./disp_z]
type = PiecewiseLinear
x = '0. 1. 2.'
y = '0. 0.006 0.0103923'
[../]
[./stress_xx]
type = ParsedFunction
# The paper gives 0.201 as the time at initial yield, but 0.20097635952803425 is the exact value.
# The paper gives -95.26 MPa as the stress at yield, but -95.26279441628823 is the exact value.
# The paper gives 12.33 as the factor in the exponential, but 12.332921390339125 is the exact value.
# 189.409039923814000, 0.170423791206825, -0.003242011311945, 1.711645501845780E-05 - exact values
symbol_names = 'timeAtYield stressAtYield expFac a b c d'
symbol_values = '0.20097635952803425 -95.26279441628823 12.332921390339125 189.409039923814000 0.170423791206825 -0.003242011311945 1.711645501845780E-05'
value = '1e6*
if(t<=timeAtYield, -474*t,
if(t<=1, stressAtYield,
(a+b*sqrt(exp(expFac*t))+c*exp(expFac*t))/(1.0+d*exp(expFac*t))))' # tends to -a
[../]
[./stress_yy]
type = ParsedFunction
# The paper gives 0.201 as the time at initial yield, but 0.20097635952803425 is the exact value.
# the paper gives -95.26 MPa as the stress at yield, but -95.26279441628823 is the exact value.
# The paper gives 12.33 as the factor in the exponential, but 12.332921390339125 is the exact value.
# -76.867432297315000, -1.442488120272900, 0.001315697947301, 1.711645501845780E-05 - exact values
symbol_names = 'timeAtYield stressAtYield expFac a b c d'
symbol_values = '0.20097635952803425 -95.26279441628823 12.332921390339125 -76.867432297315000 -1.442488120272900 0.001315697947301 1.711645501845780E-05'
value = '1e6*
if(t<=timeAtYield, -474*t,
if(t<=1, stressAtYield,
(a+b*sqrt(exp(expFac*t))+c*exp(expFac*t))/(1.0+d*exp(expFac*t))))' # tends to -a
[../]
[./stress_zz]
type = ParsedFunction
# The paper gives 0.201 as the time at initial yield, but 0.20097635952803425 is the exact value.
# the paper gives 190.5 MPa as the stress at yield, but 190.52558883257645 is the exact value.
# The paper gives 12.33 as the factor in the exponential, but 12.332921390339125 is the exact value.
# -112.541607626499000, 1.272064329066080, 0.001926313364644, 1.711645501845780E-05 - exact values
symbol_names = 'timeAtYield stressAtYield expFac a b c d'
symbol_values = '0.20097635952803425 190.52558883257645 12.332921390339125 -112.541607626499000 1.272064329066080 0.001926313364644 1.711645501845780E-05'
value = '1e6*
if(t<=timeAtYield, 948*t,
if(t<=1, stressAtYield,
(a+b*sqrt(exp(expFac*t))+c*exp(expFac*t))/(1.0+d*exp(expFac*t))))' # tends to -a
[../]
[]
[Variables]
[./disp_x]
order = FIRST
family = LAGRANGE
[../]
[./disp_y]
order = FIRST
family = LAGRANGE
[../]
[./disp_z]
order = FIRST
family = LAGRANGE
[../]
[]
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./vonmises]
order = CONSTANT
family = MONOMIAL
[../]
[./plastic_strain_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./plastic_strain_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./plastic_strain_zz]
order = CONSTANT
family = MONOMIAL
[../]
[]
[Kernels]
[SolidMechanics]
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
execute_on = 'timestep_end'
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
execute_on = 'timestep_end'
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
execute_on = 'timestep_end'
[../]
[./vonmises]
type = RankTwoScalarAux
rank_two_tensor = stress
variable = vonmises
scalar_type = vonmisesStress
execute_on = 'timestep_end'
[../]
[./plastic_strain_xx]
type = RankTwoAux
rank_two_tensor = plastic_strain
variable = plastic_strain_xx
index_i = 0
index_j = 0
execute_on = 'timestep_end'
[../]
[./plastic_strain_yy]
type = RankTwoAux
rank_two_tensor = plastic_strain
variable = plastic_strain_yy
index_i = 1
index_j = 1
execute_on = 'timestep_end'
[../]
[./plastic_strain_zz]
type = RankTwoAux
rank_two_tensor = plastic_strain
variable = plastic_strain_zz
index_i = 2
index_j = 2
execute_on = 'timestep_end'
[../]
[]
[BCs]
[./fixed_x]
type = DirichletBC
variable = disp_x
boundary = left
value = 0.0
[../]
[./fixed_y]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[../]
[./fixed_z]
type = DirichletBC
variable = disp_z
boundary = back
value = 0.0
[../]
[./disp_x]
type = FunctionDirichletBC
variable = disp_x
boundary = right
function = disp_x
[../]
[./disp_y]
type = FunctionDirichletBC
variable = disp_y
boundary = top
function = disp_y
[../]
[./disp_z]
type = FunctionDirichletBC
variable = disp_z
boundary = front
function = disp_z
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
youngs_modulus = 210666666666.666667
poissons_ratio = 0.3333333333333333
[../]
[./strain]
type = ComputeIncrementalStrain
[../]
[./isotropic_plasticity]
type = IsotropicPlasticityStressUpdate
yield_stress = 285788383.2488647 # = sqrt(3)*165e6 = sqrt(3) * yield in shear
hardening_constant = 0.0
[../]
[./radial_return_stress]
type = ComputeMultipleInelasticStress
tangent_operator = elastic
inelastic_models = 'isotropic_plasticity'
[../]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options = '-snes_ksp_ew'
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
nl_abs_tol = 1e-10
l_max_its = 20
start_time = 0.0
dt = 0.01 # use 0.0001 for a nearly exact match
end_time = 2.0
[]
[Postprocessors]
[./analytic_xx]
type = FunctionValuePostprocessor
function = stress_xx
[../]
[./analytic_yy]
type = FunctionValuePostprocessor
function = stress_yy
[../]
[./analytic_zz]
type = FunctionValuePostprocessor
function = stress_zz
[../]
[./stress_xx]
type = ElementalVariableValue
variable = stress_xx
elementid = 0
[../]
[./stress_yy]
type = ElementalVariableValue
variable = stress_yy
elementid = 0
[../]
[./stress_zz]
type = ElementalVariableValue
variable = stress_zz
elementid = 0
[../]
[./stress_xx_l2_error]
type = ElementL2Error
variable = stress_xx
function = stress_xx
[../]
[./stress_yy_l2_error]
type = ElementL2Error
variable = stress_yy
function = stress_yy
[../]
[./stress_zz_l2_error]
type = ElementL2Error
variable = stress_zz
function = stress_zz
[../]
[]
[Outputs]
exodus = true
[]
(modules/solid_mechanics/test/tests/recompute_radial_return/isotropic_plasticity_finite_strain.i)
# This simulation uses the piece-wise linear strain hardening model
# with the incremental small strain formulation; incremental small strain
# is required to produce the strain_increment for the DiscreteRadialReturnStressIncrement
# class, which handles the calculation of the stress increment to return
# to the yield surface in a J2 (isotropic) plasticity problem.
#
# This test assumes a Poissons ratio of 0.3 and applies a displacement loading
# condition on the top in the y direction.
#
# An identical problem was run in Abaqus on a similar 1 element mesh and was used
# to verify the SolidMechanics solution; this SolidMechanics code matches the
# SolidMechanics solution.
#
# Mechanical strain is the sum of the elastic and plastic strains but is different
# from total strain in cases with eigen strains, e.g. thermal strain.
[Mesh]
file = 1x1x1cube.e
[]
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Functions]
[./top_pull]
type = ParsedFunction
expression = t*(0.0625)
[../]
[./hf]
type = PiecewiseLinear
x = '0 0.001 0.003 0.023'
y = '50 52 54 56'
[../]
[]
[Physics/SolidMechanics/QuasiStatic]
[./all]
strain = FINITE
add_variables = true
generate_output = 'stress_yy plastic_strain_xx plastic_strain_yy plastic_strain_zz'
[../]
[]
[BCs]
[./y_pull_function]
type = FunctionDirichletBC
variable = disp_y
boundary = 5
function = top_pull
[../]
[./x_bot]
type = DirichletBC
variable = disp_x
boundary = 4
value = 0.0
[../]
[./y_bot]
type = DirichletBC
variable = disp_y
boundary = 3
value = 0.0
[../]
[./z_bot]
type = DirichletBC
variable = disp_z
boundary = 2
value = 0.0
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
youngs_modulus = 2.1e5
poissons_ratio = 0.3
[../]
[./isotropic_plasticity]
type = IsotropicPlasticityStressUpdate
yield_stress = 50.0
hardening_function = hf
[../]
[./radial_return_stress]
type = ComputeMultipleInelasticStress
tangent_operator = elastic
inelastic_models = 'isotropic_plasticity'
[../]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options = '-snes_ksp_ew'
petsc_options_iname = '-ksp_gmres_restart'
petsc_options_value = '101'
line_search = 'none'
l_max_its = 50
nl_max_its = 50
nl_rel_tol = 1e-12
nl_abs_tol = 1e-10
l_tol = 1e-9
start_time = 0.0
end_time = 0.075
dt = 0.00125
dtmin = 0.0001
[]
[Outputs]
[./out]
type = Exodus
elemental_as_nodal = true
[../]
[]
(modules/solid_mechanics/test/tests/recompute_radial_return/isotropic_plasticity_errors.i)
# This simulation uses the piece-wise linear strain hardening model
# with the incremental small strain formulation; incremental small strain
# is required to produce the strain_increment for the DiscreteRadialReturnStressIncrement
# class, which handles the calculation of the stress increment to return
# to the yield surface in a J2 (isotropic) plasticity problem.
#
# This test is used to check the error messages in the discrete radial return
# model DiscreteRRIsotropicPlasticity; cli_args are used to check all of the
# error messages in the DiscreteRRIsotropicPlasticity model.
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
[]
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Functions]
[./top_pull]
type = ParsedFunction
expression = t*(0.0625)
[../]
[./harden_func]
type = PiecewiseLinear
x = '0 0.0003 0.0007 0.0009'
y = '50 52 54 56'
[../]
[]
[Physics/SolidMechanics/QuasiStatic]
[./all]
strain = SMALL
incremental = true
add_variables = true
generate_output = 'stress_yy plastic_strain_xx plastic_strain_yy plastic_strain_zz'
[../]
[]
[BCs]
[./y_pull_function]
type = FunctionDirichletBC
variable = disp_y
boundary = top
function = top_pull
[../]
[./x_bot]
type = DirichletBC
variable = disp_x
boundary = left
value = 0.0
[../]
[./y_bot]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[../]
[./z_bot]
type = DirichletBC
variable = disp_z
boundary = back
value = 0.0
[../]
[]
[Materials]
[./elasticity_tensor]
[../]
[./isotropic_plasticity]
type = IsotropicPlasticityStressUpdate
relative_tolerance = 1e-25
absolute_tolerance = 1e-5
[../]
[./radial_return_stress]
type = ComputeMultipleInelasticStress
tangent_operator = elastic
inelastic_models = 'isotropic_plasticity'
[../]
[]
[Executioner]
type = Transient
#Preconditioned JFNK (default)
solve_type = 'PJFNK'
petsc_options = '-snes_ksp_ew'
petsc_options_iname = '-ksp_gmres_restart'
petsc_options_value = '101'
line_search = 'none'
l_max_its = 100
nl_max_its = 100
nl_rel_tol = 1e-18
nl_abs_tol = 1e-10
l_tol = 1e-12
start_time = 0.0
end_time = 0.025
dt = 0.00125
dtmin = 0.0001
[]
[Outputs]
[./out]
type = Exodus
elemental_as_nodal = true
[../]
[]
(modules/solid_mechanics/test/tests/material_limit_time_step/elas_plas/nafems_nl1_lim.i)
#
# Tests material model IsotropicPlasticity with material based time stepper
# Boundary conditions from NAFEMS test NL1
#
[GlobalParams]
order = FIRST
family = LAGRANGE
volumetric_locking_correction = true
displacements = 'disp_x disp_y'
[]
[Mesh]#Comment
file = one_elem2.e
[] # Mesh
[Variables]
[./disp_x]
[../]
[./disp_y]
[../]
[] # Variables
[AuxVariables]
[./stress_xx]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_zz]
order = CONSTANT
family = MONOMIAL
[../]
[./stress_xy]
order = CONSTANT
family = MONOMIAL
[../]
[./vonmises]
order = CONSTANT
family = MONOMIAL
[../]
[./elastic_strain_yy]
order = CONSTANT
family = MONOMIAL
[../]
[./plastic_strain_eff]
order = CONSTANT
family = MONOMIAL
[../]
[./tot_strain_yy]
order = CONSTANT
family = MONOMIAL
[../]
[] # AuxVariables
[Kernels]
[SolidMechanics]
use_displaced_mesh = true
[../]
[]
[AuxKernels]
[./stress_xx]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xx
index_i = 0
index_j = 0
[../]
[./stress_yy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_yy
index_i = 1
index_j = 1
[../]
[./stress_zz]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_zz
index_i = 2
index_j = 2
[../]
[./stress_xy]
type = RankTwoAux
rank_two_tensor = stress
variable = stress_xy
index_i = 0
index_j = 1
[../]
[./vonmises]
type = RankTwoScalarAux
rank_two_tensor = stress
variable = vonmises
scalar_type = VonMisesStress
execute_on = timestep_end
[../]
[./elastic_strain_yy]
type = RankTwoAux
rank_two_tensor = elastic_strain
variable = elastic_strain_yy
index_i = 1
index_j = 1
[../]
[./plastic_strain_eff]
type = MaterialRealAux
property = effective_plastic_strain
variable = plastic_strain_eff
[../]
[./tot_strain_yy]
type = RankTwoAux
rank_two_tensor = total_strain
variable = tot_strain_yy
index_i = 1
index_j = 1
[../]
[] # AuxKernels
[Functions]
[./appl_dispx]
type = PiecewiseLinear
x = '0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0'
y = '0.0 0.25e-4 0.50e-4 0.50e-4 0.50e-4 0.25e-4 0.0 0.0 0.0'
[../]
[./appl_dispy]
type = PiecewiseLinear
x = '0 1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0'
y = '0.0 0.0 0.0 0.25e-4 0.50e-4 0.50e-4 0.50e-4 0.25e-4 0.0 '
[../]
[]
[BCs]
[./side_x]
type = DirichletBC
variable = disp_x
boundary = 101
value = 0.0
[../]
[./origin_x]
type = DirichletBC
variable = disp_x
boundary = 103
value = 0.0
[../]
[./bot_y]
type = DirichletBC
variable = disp_y
boundary = 102
value = 0.0
[../]
[./origin_y]
type = DirichletBC
variable = disp_y
boundary = 103
value = 0.0
[../]
[./top_y]
type = FunctionDirichletBC
variable = disp_y
boundary = 1
function = appl_dispy
[../]
[./right_x]
type = FunctionDirichletBC
variable = disp_x
boundary = 2
function = appl_dispx
[../]
[]
[Materials]
[./elasticity_tensor]
type = ComputeIsotropicElasticityTensor
block = 1
youngs_modulus = 250e9
poissons_ratio = 0.25
[../]
[./strain]
type = ComputePlaneFiniteStrain
block = 1
[../]
[./stress]
type = ComputeMultipleInelasticStress
inelastic_models = 'isoplas'
block = 1
[../]
[./isoplas]
type = IsotropicPlasticityStressUpdate
yield_stress = 5e6
hardening_constant = 0.0
relative_tolerance = 1e-20
absolute_tolerance = 1e-8
max_inelastic_increment = 0.000001
[../]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
nl_rel_tol = 1e-10
nl_abs_tol = 1e-12
l_tol = 1e-4
l_max_its = 100
nl_max_its = 20
[./TimeStepper]
type = IterationAdaptiveDT
dt = 0.1
time_t = '1.0 2.0 3.0 4.0 5.0 6.0 7.0 8.0'
time_dt = '0.1 0.1 0.1 0.1 0.1 0.1 0.1 0.1'
optimal_iterations = 30
iteration_window = 9
growth_factor = 2.0
cutback_factor = 0.5
timestep_limiting_postprocessor = matl_ts_min
[../]
start_time = 0.0
num_steps = 1000
end_time = 8.0
[] # Executioner
[Postprocessors]
[./matl_ts_min]
type = MaterialTimeStepPostprocessor
[../]
[./stress_xx]
type = ElementAverageValue
variable = stress_xx
[../]
[./stress_yy]
type = ElementAverageValue
variable = stress_yy
[../]
[./stress_zz]
type = ElementAverageValue
variable = stress_zz
[../]
[./vonmises]
type = ElementAverageValue
variable = vonmises
[../]
[./el_strain_yy]
type = ElementAverageValue
variable = elastic_strain_yy
[../]
[./plas_strain_eff]
type = ElementAverageValue
variable = plastic_strain_eff
[../]
[./tot_strain_yy]
type = ElementAverageValue
variable = tot_strain_yy
[../]
[./disp_x1]
type = NodalVariableValue
nodeid = 0
variable = disp_x
[../]
[./disp_x4]
type = NodalVariableValue
nodeid = 3
variable = disp_x
[../]
[./disp_y1]
type = NodalVariableValue
nodeid = 0
variable = disp_y
[../]
[./disp_y4]
type = NodalVariableValue
nodeid = 3
variable = disp_y
[../]
[./_dt]
type = TimestepSize
[../]
[]
[Outputs]
exodus = true
csv = true
[./console]
type = Console
output_linear = true
[../]
[] # Outputs
(modules/solid_mechanics/test/tests/power_law_creep/cp_power_law_creep.i)
# 1x1x1 unit cube with uniform pressure on top face
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Mesh]
type = GeneratedMesh
dim = 3
nx = 1
ny = 1
nz = 1
[]
[Variables]
[temp]
order = FIRST
family = LAGRANGE
initial_condition = 1000.0
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = FINITE
incremental = true
add_variables = true
generate_output = 'stress_yy creep_strain_xx creep_strain_yy creep_strain_zz elastic_strain_yy'
[]
[]
[Functions]
[top_pull]
type = PiecewiseLinear
x = '0 1'
y = '1 1'
[]
[]
[Kernels]
[heat]
type = Diffusion
variable = temp
[]
[heat_ie]
type = TimeDerivative
variable = temp
[]
[]
[BCs]
[u_top_pull]
type = Pressure
variable = disp_y
boundary = top
factor = -10.0e6
function = top_pull
[]
[u_bottom_fix]
type = DirichletBC
variable = disp_y
boundary = bottom
value = 0.0
[]
[u_yz_fix]
type = DirichletBC
variable = disp_x
boundary = left
value = 0.0
[]
[u_xy_fix]
type = DirichletBC
variable = disp_z
boundary = back
value = 0.0
[]
[temp_fix]
type = DirichletBC
variable = temp
boundary = 'bottom top'
value = 1000.0
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeIsotropicElasticityTensor
youngs_modulus = 2e11
poissons_ratio = 0.3
[]
[radial_return_stress]
type = ComputeCreepPlasticityStress
creep_model = power_law_creep
plasticity_model = isotropic_plasticity
tangent_operator = elastic
[]
[power_law_creep]
type = PowerLawCreepStressUpdate
coefficient = 1.0e-15
n_exponent = 4
activation_energy = 3.0e5
temperature = temp
[]
[isotropic_plasticity]
type = IsotropicPlasticityStressUpdate
yield_stress = 1e30
hardening_constant = 0.0
[]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options = '-snes_ksp'
petsc_options_iname = '-ksp_gmres_restart'
petsc_options_value = '101'
line_search = 'none'
l_max_its = 20
nl_max_its = 20
nl_rel_tol = 1e-6
nl_abs_tol = 1e-6
l_tol = 1e-5
start_time = 0.0
end_time = 1.0
num_steps = 10
dt = 0.1
[]
[Outputs]
exodus = true
[]