idaholab / moose

Multiphysics Object Oriented Simulation Environment
https://www.mooseframework.org
GNU Lesser General Public License v2.1
1.72k stars 1.04k forks source link

Fix volumetric_locking_correction in tensor mechanics #11220

Closed tophmatthews closed 5 years ago

tophmatthews commented 6 years ago

Rationale

The locking correction in tensor mechanics tends to make unphysical oscillations in displacements, stresses, and strains. The locking correction should be fixed to retain physical behaviors for materials with high poisson's ratio or high creep strains.

Description

The following input deck can be used to visualize stress, strain, and displacements at the top and bottom surface of a 2DRZ rodlet.

Input File ``` [Mesh] type = GeneratedMesh dim = 2 xmax = 1e-3 ymax = 1e-3 nx = 10 ny = 10 second_order = true [] [GlobalParams] displacements = 'disp_x disp_y' volumetric_locking_correction = true # false order = SECOND [] [Problem] type = ReferenceResidualProblem reference_residual_variables = 'saved_x saved_y saved_temp' solution_variables = 'disp_x disp_y temp' coord_type = RZ [] [Variables] [./temp] initial_condition = 300 [../] [] [Kernels] [./temp_dt] type = HeatConductionTimeDerivative specific_heat = 1 density_name = 1 variable = temp save_in = saved_temp [../] [./temp_diff] type = HeatConduction diffusion_coefficient = 1 variable = temp save_in = saved_temp [../] [./source] type = BodyForce variable = temp function = ramp save_in = saved_temp [../] [] [Functions] [./ramp] type = PiecewiseLinear x = '0 10' y = '0 1e9' [../] [] [Modules/TensorMechanics/Master] [./all] strain = FINITE add_variables = true temperature = temp eigenstrain_names = 'reduced' # 'thermal' save_in = 'saved_x saved_y' use_displaced_mesh = true [../] [] [AuxVariables] [./saved_x] [../] [./saved_y] [../] [./saved_temp] [../] [./stress_theta] order = CONSTANT family = MONOMIAL [../] [./strain_radial_constant] order = CONSTANT family = MONOMIAL [../] [./strain_radial_first] order = FIRST family = MONOMIAL [../] [./strain_hoop_constant] order = CONSTANT family = MONOMIAL [../] [./strain_hoop_first] order = FIRST family = MONOMIAL [../] [./stress_hoop_constant] order = CONSTANT family = MONOMIAL [../] [./stress_hoop_first] order = FIRST family = MONOMIAL [../] [./stress_hoop_second] order = SECOND family = MONOMIAL [../] [./stress_radial_constant] order = CONSTANT family = MONOMIAL [../] [./stress_radial_first] order = FIRST family = MONOMIAL [../] [./stress_radial_second] order = SECOND family = MONOMIAL [../] [./inelastic_constant] order = CONSTANT family = MONOMIAL [../] [./inelastic_first] order = FIRST family = MONOMIAL [../] [./strain_constant] order = CONSTANT family = MONOMIAL [../] [./strain_first] order = FIRST family = MONOMIAL [../] [] [AuxKernels] [./stress_theta] type = RankTwoAux rank_two_tensor = stress index_i = 2 index_j = 2 variable = stress_theta execute_on = timestep_end [../] [./stress_hoop_constant] type = RankTwoAux rank_two_tensor = stress index_i = 2 index_j = 2 variable = stress_hoop_constant execute_on = timestep_end [../] [./stress_hoop_first] type = RankTwoAux rank_two_tensor = stress index_i = 2 index_j = 2 variable = stress_hoop_first execute_on = timestep_end [../] [./stress_hoop_second] type = RankTwoAux rank_two_tensor = stress index_i = 2 index_j = 2 variable = stress_hoop_second execute_on = timestep_end [../] [./stress_radial_constant] type = RankTwoAux rank_two_tensor = stress index_i = 0 index_j = 0 variable = stress_radial_constant execute_on = timestep_end [../] [./stress_radial_first] type = RankTwoAux rank_two_tensor = stress index_i = 0 index_j = 0 variable = stress_radial_first execute_on = timestep_end [../] [./stress_radial_second] type = RankTwoAux rank_two_tensor = stress index_i = 0 index_j = 0 variable = stress_radial_second execute_on = timestep_end [../] [./strain_radial_constant] type = RankTwoAux rank_two_tensor = total_strain index_i = 0 index_j = 0 variable = strain_radial_constant execute_on = timestep_end [../] [./strain_radial_first] type = RankTwoAux rank_two_tensor = total_strain index_i = 0 index_j = 0 variable = strain_radial_first execute_on = timestep_end [../] [./strain_hoop_constant] type = RankTwoAux rank_two_tensor = total_strain index_i = 2 index_j = 2 variable = strain_hoop_constant execute_on = timestep_end [../] [./strain_hoop_first] type = RankTwoAux rank_two_tensor = total_strain index_i = 2 index_j = 2 variable = strain_hoop_first execute_on = timestep_end [../] [./inelastic_constant] type = RankTwoScalarAux rank_two_tensor = combined_inelastic_strain variable = inelastic_constant scalar_type = EffectiveStrain execute_on = timestep_end [../] [./inelastic_first] type = RankTwoScalarAux rank_two_tensor = combined_inelastic_strain variable = inelastic_first scalar_type = EffectiveStrain execute_on = timestep_end [../] [./strain_constant] type = RankTwoScalarAux rank_two_tensor = total_strain variable = strain_constant scalar_type = EffectiveStrain execute_on = timestep_end [../] [./strain_first] type = RankTwoScalarAux rank_two_tensor = total_strain variable = strain_first scalar_type = EffectiveStrain execute_on = timestep_end [../] [] [Materials] [./elasticity_tensor] type = ComputeIsotropicElasticityTensor youngs_modulus = 1e10 poissons_ratio = 0.3 [../] [./elastic_strain] type = ComputeMultipleInelasticStress inelastic_models = creep [../] [./creep] type = PowerLawCreepStressUpdate coefficient = 0 # 1e-25 n_exponent = 4 activation_energy = 0 outputs = all [../] [./thexp] type = ComputeThermalExpansionEigenstrain temperature = temp eigenstrain_name = thermal stress_free_temperature = 300 thermal_expansion_coeff = 1e-6 [../] [./eigen] type = ComputeReducedOrderEigenstrain eigenstrain_name = reduced input_eigenstrain_names = thermal [../] [] [BCs] [./temp_surf] type = DirichletBC variable = temp boundary = right value = 300 [../] [./no_disp_x] type = PresetBC variable = disp_x boundary = left value = 0.0 [../] [./no_disp_y] type = PresetBC variable = disp_y boundary = bottom value = 0.0 [../] [] [Preconditioning] [./smp] type = SMP full = true [../] [] [Executioner] type = Transient petsc_options_iname = '-pc_type -pc_factor_mat_solver_package -ksp_gmres_restart' petsc_options_value = 'lu superlu_dist 51' line_search = 'none' nl_rel_tol = 1e-5 nl_max_its = 20 l_max_its = 100 end_time = 10 dt = 1 dtmin = 1e-2 [] [Postprocessors] [./temp_avg] type = ElementAverageValue variable = temp [../] [./strainTheta] type = ElementAverageValue variable = strain_hoop_constant [../] [./stressTheta] type = ElementAverageValue variable = stress_theta [../] [./temp_max] type = ElementExtremeValue variable = temp [../] [] [VectorPostprocessors] [./radial_top] type = LineMaterialRankTwoSampler index_i = 0 index_j = 0 property = stress sort_by = x start = '0 1e-3 0' end = '1e-3 1e-3 0' [../] [./hoop_top] type = LineMaterialRankTwoSampler index_i = 2 index_j = 2 property = stress sort_by = x start = '0 1e-3 0' end = '1e-3 1e-3 0' [../] [./radial_bottom] type = LineMaterialRankTwoSampler index_i = 0 index_j = 0 property = stress sort_by = x start = '0 0 0' end = '1e-3 0 0' [../] [./hoop_bottom] type = LineMaterialRankTwoSampler index_i = 2 index_j = 2 property = stress sort_by = x start = '0 0 0' end = '1e-3 0 0' [../] [./line_top] type = LineValueSampler num_points = 101 sort_by = x start_point = '0 1e-3 0' end_point = '1e-3 1e-3 0' variable = 'stress_hoop_constant stress_hoop_first stress_hoop_second stress_radial_constant stress_radial_first stress_radial_second temp disp_x disp_y effective_creep_strain inelastic_constant inelastic_first strain_radial_constant strain_radial_first strain_hoop_constant strain_hoop_first' [../] [./line_bottom] type = LineValueSampler num_points = 101 sort_by = x start_point = '0 0 0' end_point = '1e-3 0 0' variable = 'stress_hoop_constant stress_hoop_first stress_hoop_second stress_radial_constant stress_radial_first stress_radial_second temp disp_x disp_y effective_creep_strain inelastic_constant inelastic_first strain_radial_constant strain_radial_first strain_hoop_constant strain_hoop_first' [../] [] [Outputs] exodus = true [./csv] type = CSV execute_on = final [../] [] [Debug] show_var_residual_norms = true [] ```

Impact

This will provide a correct locking correction for tensor mechanics stress/strain calculations.

tophmatthews commented 6 years ago

Ping @jasondhales and @bwspenc, following our conversation last Friday.

tophmatthews commented 6 years ago

Also, a segfault occurs when using 2DRZ with volumetric_locking_correction and use_finite_deform_jacobian #11261

jasondhales commented 5 years ago

@sveerara @bwspenc @novasr Here is a ticket where @tophmatthews supplies an input file that shows bad behavior.

tophmatthews commented 5 years ago

Interesting way to "solve" the problem 😆

*** Warning ***
Volumteric locking correction is not required for second order elements. Using volumetric locking with second order elements could cause zigzag patterns in stresses and strains.

In reality, it's second order elements without locking correction looks good. First order still has stress zigzags going the wrong direction with the locking correction: image

Good to see the problems have been addressed. I would almost make this an error, not a warning, but it is what it is. FYI @novasr.