- density_sensitivityName of the density_sensitivity variable.
C++ Type:VariableName
Unit:(no unit assumed)
Controllable:No
Description:Name of the density_sensitivity variable.
- design_densityDesign density variable name.
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Design density variable name.
- volume_fractionVolume Fraction
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Volume Fraction
Density Update
Compute updated densities based on sensitivities using an optimality criteria method to keep the volume constraint satisified.
Description
The DensityUpdate
class is an ElementUserObject
that calculates and updates the design densities based on the filtered sensitivities. The density update is performed using the Solid Isotropic Material with Penalization (SIMP) method. The class requires the design density variable, the filtered design density variable, the density sensitivity variable, the penalty power, and the volume fraction as inputs.
Example Input File
An example of how to use the DensityUpdate
class in an input file:
listing test/tests/materials/compliance_sensitivity/2d_mbb.i block=UserObjects/update
Input Parameters
- bisection_lower_bound0Lower bound for the bisection algorithm.
Default:0
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Lower bound for the bisection algorithm.
- bisection_upper_bound1e+16Upper bound for the bisection algorithm.
Default:1e+16
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Upper bound for the bisection algorithm.
- blockThe list of blocks (ids or names) that this object will be applied
C++ Type:std::vector<SubdomainName>
Unit:(no unit assumed)
Controllable:No
Description:The list of blocks (ids or names) that this object will be applied
- execute_onTIMESTEP_ENDThe list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html.
Default:TIMESTEP_END
C++ Type:ExecFlagEnum
Unit:(no unit assumed)
Options:XFEM_MARK, FORWARD, ADJOINT, HOMOGENEOUS_FORWARD, ADJOINT_TIMESTEP_BEGIN, ADJOINT_TIMESTEP_END, NONE, INITIAL, LINEAR, NONLINEAR_CONVERGENCE, NONLINEAR, POSTCHECK, TIMESTEP_END, TIMESTEP_BEGIN, MULTIAPP_FIXED_POINT_END, MULTIAPP_FIXED_POINT_BEGIN, FINAL, CUSTOM
Controllable:No
Description:The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html.
- 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
Unit:(no unit assumed)
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.
Optional Parameters
- allow_duplicate_execution_on_initialFalseIn the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).
Default:False
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:In the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).
- control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Unit:(no unit assumed)
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
- enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Unit:(no unit assumed)
Controllable:Yes
Description:Set the enabled status of the MooseObject.
- execution_order_group0Execution order groups are executed in increasing order (e.g., the lowest number is executed first). Note that negative group numbers may be used to execute groups before the default (0) group. Please refer to the user object documentation for ordering of user object execution within a group.
Default:0
C++ Type:int
Unit:(no unit assumed)
Controllable:No
Description:Execution order groups are executed in increasing order (e.g., the lowest number is executed first). Note that negative group numbers may be used to execute groups before the default (0) group. Please refer to the user object documentation for ordering of user object execution within a group.
- force_postauxFalseForces the UserObject to be executed in POSTAUX
Default:False
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:Forces the UserObject to be executed in POSTAUX
- force_preauxFalseForces the UserObject to be executed in PREAUX
Default:False
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:Forces the UserObject to be executed in PREAUX
- force_preicFalseForces the UserObject to be executed in PREIC during initial setup
Default:False
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:Forces the UserObject to be executed in PREIC during initial setup
- implicitTrueDetermines whether this object is calculated using an implicit or explicit form
Default:True
C++ Type:bool
Unit:(no unit assumed)
Controllable:No
Description:Determines whether this object is calculated using an implicit or explicit form
- seed0The seed for the master random number generator
Default:0
C++ Type:unsigned int
Unit:(no unit assumed)
Controllable:No
Description:The seed for the master random number generator
- 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
Unit:(no unit assumed)
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
Input Files
- (modules/combined/test/tests/optimization/compliance_sensitivity/2d_mbb.i)
- (modules/combined/examples/optimization/2d_mbb.i)
- (modules/combined/test/tests/optimization/optimization_density_update/top_opt_2d.i)
- (modules/combined/examples/optimization/3d_mbb.i)
- (modules/combined/test/tests/optimization/optimization_density_update/top_opt_2d_pde_filter.i)
- (modules/combined/test/tests/optimization/thermal_sensitivity/2d_root.i)
- (modules/combined/test/tests/optimization/optimization_density_update/top_opt_3d_pde_filter.i)
- (modules/combined/test/tests/optimization/compliance_sensitivity/2d_mbb_pde.i)
- (modules/combined/examples/optimization/helmholtz_multimat_strip.i)
- (modules/combined/examples/optimization/2d_mbb_pde_amr.i)
- (modules/optimization/test/tests/simp/2d.i)
- (modules/combined/test/tests/optimization/compliance_sensitivity/3d_mbb.i)
- (modules/combined/examples/optimization/thermomechanical/thermomechanical_main.i)
- (modules/combined/examples/optimization/multi-load/single_main.i)
- (modules/combined/test/tests/optimization/compliance_sensitivity/2d_mbb_pde_amr.i)
- (modules/combined/test/tests/optimization/optimization_density_update/top_opt_3d.i)
- (modules/combined/test/tests/optimization/compliance_sensitivity/2d_mmb_2material.i)
- (modules/combined/examples/optimization/2d_mbb_pde.i)
- (modules/combined/examples/optimization/helmholtz_multimat_nostrip.i)
- (modules/combined/examples/optimization/multi-load/square_main.i)
- (modules/combined/test/tests/optimization/compliance_sensitivity/2d_mmb_2material_cost.i)
(modules/combined/test/tests/optimization/compliance_sensitivity/2d_mbb.i)
vol_frac = 0.5
E0 = 1
Emin = 1e-8
power = 2
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 2
nx = 150
ny = 50
xmin = 0
xmax = 30
ymin = 0
ymax = 10
[]
[node]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = hold
nodes = 0
[]
[push]
type = ExtraNodesetGenerator
input = node
new_boundary = push
coord = '30 10 0'
[]
[]
[Variables]
[disp_x]
[]
[disp_y]
[]
[]
[AuxVariables]
[Dc]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = disp_y
boundary = hold
value = 0.0
[]
[no_y]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[]
[NodalKernels]
[push]
type = NodalGravity
variable = disp_y
boundary = push
gravity_value = -1
mass = 1
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'mat_den'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys
incremental = false
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[rad_avg]
type = RadialAverage
radius = 1.2
weights = linear
prop_name = sensitivity
execute_on = TIMESTEP_END
force_preaux = true
[]
[update]
type = DensityUpdate
density_sensitivity = Dc
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
[]
[calc_sense]
type = SensitivityFilter
density_sensitivity = Dc
design_density = mat_den
filter_UO = rad_avg
execute_on = TIMESTEP_END
force_postaux = true
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
nl_abs_tol = 1e-8
dt = 1.0
num_steps = 70
[]
[Outputs]
[out]
type = CSV
execute_on = 'TIMESTEP_END'
[]
print_linear_residuals = false
[]
[Postprocessors]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[sensitivity]
type = ElementIntegralMaterialProperty
mat_prop = sensitivity
[]
[]
(modules/combined/examples/optimization/2d_mbb.i)
vol_frac = 0.5
E0 = 1
Emin = 1e-8
power = 2
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 2
nx = 150
ny = 50
xmin = 0
xmax = 30
ymin = 0
ymax = 10
[]
[node]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = pull
nodes = 0
[]
[push]
type = ExtraNodesetGenerator
input = node
new_boundary = push
coord = '30 10 0'
[]
[]
[AuxVariables]
[Emin]
family = MONOMIAL
order = CONSTANT
initial_condition = ${Emin}
[]
[power]
family = MONOMIAL
order = CONSTANT
initial_condition = ${power}
[]
[E0]
family = MONOMIAL
order = CONSTANT
initial_condition = ${E0}
[]
[Dc]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = disp_y
boundary = pull
value = 0.0
[]
[no_y]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[]
[NodalKernels]
[pull]
type = NodalGravity
variable = disp_y
boundary = push
gravity_value = -1
mass = 1
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'Emin mat_den power E0'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys
incremental = false
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[rad_avg]
type = RadialAverage
radius = 1.2
weights = linear
prop_name = sensitivity
execute_on = TIMESTEP_END
force_preaux = true
[]
[update]
type = DensityUpdate
density_sensitivity = Dc
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
[]
[calc_sense]
type = SensitivityFilter
density_sensitivity = Dc
design_density = mat_den
filter_UO = rad_avg
execute_on = TIMESTEP_END
force_postaux = true
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
nl_abs_tol = 1e-8
dt = 1.0
num_steps = 70
[]
[Outputs]
[out]
type = Exodus
execute_on = 'TIMESTEP_END'
[]
print_linear_residuals = false
[]
[Postprocessors]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[]
(modules/combined/test/tests/optimization/optimization_density_update/top_opt_2d.i)
vol_frac = 0.4
E0 = 1e5
Emin = 1e-4
power = 2
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 2
nx = 40
ny = 20
xmin = 0
xmax = 20
ymin = 0
ymax = 10
[]
[node]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = pull
nodes = 0
[]
[]
[AuxVariables]
[sensitivity]
family = MONOMIAL
order = FIRST
initial_condition = -1.0
[AuxKernel]
type = MaterialRealAux
variable = sensitivity
property = sensitivity
execute_on = LINEAR
[]
[]
[compliance]
family = MONOMIAL
order = CONSTANT
[]
[Dc]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[no_y]
type = DirichletBC
variable = disp_y
boundary = right
value = 0.0
[]
[]
[NodalKernels]
[pull]
type = NodalGravity
variable = disp_y
boundary = pull
gravity_value = -1
mass = 1
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'mat_den'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys
incremental = false
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[rad_avg]
type = RadialAverage
radius = 0.5
weights = constant
prop_name = sensitivity
execute_on = TIMESTEP_END
force_preaux = true
execution_order_group = -1
[]
[update]
type = DensityUpdate
density_sensitivity = Dc
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
[]
[calc_sense]
type = SensitivityFilter
density_sensitivity = Dc
design_density = mat_den
filter_UO = rad_avg
execute_on = TIMESTEP_END
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type '
petsc_options_value = 'lu'
nl_abs_tol = 1e-10
start_time = 0.0
dt = 1.0
num_steps = 50
[]
[Outputs]
[out]
type = Exodus
time_step_interval = 10
[]
[]
(modules/combined/examples/optimization/3d_mbb.i)
vol_frac = 0.5
E0 = 1
Emin = 1e-8
power = 2
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 3
nx = 60
ny = 20
nz = 20
xmin = 0
xmax = 30
ymin = 0
ymax = 10
zmin = 0
zmax = 10
[]
[node]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = hold_y
coord = '0 0 0; 0 0 10'
[]
[push]
type = ExtraNodesetGenerator
input = node
new_boundary = push
coord = '30 10 5'
[]
[]
[Variables]
[disp_z]
[]
[Dc]
initial_condition = -1.0
[]
[]
[AuxVariables]
[Emin]
family = MONOMIAL
order = CONSTANT
initial_condition = ${Emin}
[]
[power]
family = MONOMIAL
order = CONSTANT
initial_condition = ${power}
[]
[E0]
family = MONOMIAL
order = CONSTANT
initial_condition = ${E0}
[]
[sensitivity]
family = MONOMIAL
order = FIRST
initial_condition = -1.0
[AuxKernel]
type = MaterialRealAux
variable = sensitivity
property = sensitivity
execute_on = LINEAR
[]
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[Dc_elem]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[AuxKernel]
type = SelfAux
variable = Dc_elem
v = Dc
execute_on = 'TIMESTEP_END'
[]
[]
[mat_den_nodal]
family = L2_LAGRANGE
order = FIRST
initial_condition = ${vol_frac}
[AuxKernel]
type = SelfAux
execute_on = TIMESTEP_END
variable = mat_den_nodal
v = mat_den
[]
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[Kernels]
[diffusion]
type = FunctionDiffusion
variable = Dc
function = 0.15 # radius coeff
[]
[potential]
type = Reaction
variable = Dc
[]
[source]
type = CoupledForce
variable = Dc
v = sensitivity
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = disp_y
boundary = hold_y
value = 0.0
[]
[no_y]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[boundary_penalty]
type = ADRobinBC
variable = Dc
boundary = 'left top front back'
coefficient = 10
[]
[boundary_penalty_right]
type = ADRobinBC
variable = Dc
boundary = 'right'
coefficient = 10
[]
[]
[NodalKernels]
[push]
type = NodalGravity
variable = disp_y
boundary = push
gravity_value = -1
mass = 1
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'Emin mat_den power E0'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys
incremental = false
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[update]
type = DensityUpdate
density_sensitivity = Dc_elem
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
force_postaux = true
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
line_search = none
nl_abs_tol = 1e-4
l_max_its = 200
start_time = 0.0
dt = 1.0
num_steps = 70
[]
[Outputs]
[out]
type = Exodus
execute_on = 'INITIAL TIMESTEP_END'
[]
print_linear_residuals = false
[]
[Postprocessors]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[]
[Controls]
[first_period]
type = TimePeriod
start_time = 0.0
end_time = 10
enable_objects = 'BCs::boundary_penalty_right'
execute_on = 'initial timestep_begin'
[]
[]
(modules/combined/test/tests/optimization/optimization_density_update/top_opt_2d_pde_filter.i)
vol_frac = 0.4
E0 = 1e5
Emin = 1e-4
power = 2
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 2
nx = 40
ny = 20
xmin = 0
xmax = 20
ymin = 0
ymax = 10
[]
[node]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = pull
nodes = 0
[]
[]
[Variables]
[Dc]
initial_condition = -1.0
[]
[]
[AuxVariables]
[sensitivity]
family = MONOMIAL
order = FIRST
initial_condition = -1.0
[AuxKernel]
type = MaterialRealAux
variable = sensitivity
property = sensitivity
execute_on = LINEAR
[]
[]
[compliance]
family = MONOMIAL
order = CONSTANT
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[Dc_elem]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[AuxKernel]
type = SelfAux
variable = Dc_elem
v = Dc
execute_on = 'TIMESTEP_END'
[]
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[Kernels]
[diffusion]
type = FunctionDiffusion
variable = Dc
function = 0.05
[]
[potential]
type = Reaction
variable = Dc
[]
[source]
type = CoupledForce
variable = Dc
v = sensitivity
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[no_y]
type = DirichletBC
variable = disp_y
boundary = right
value = 0.0
[]
[boundary_penalty]
type = ADRobinBC
variable = Dc
boundary = 'left top'
coefficient = 10
[]
[]
[NodalKernels]
[pull]
type = NodalGravity
variable = disp_y
boundary = pull
gravity_value = -1
mass = 1
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'mat_den'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys
incremental = false
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[update]
type = DensityUpdate
density_sensitivity = Dc_elem
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
force_postaux = true
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type '
petsc_options_value = 'lu'
nl_abs_tol = 1e-10
line_search = none
dt = 1.0
num_steps = 30
[]
[Outputs]
[out]
type = Exodus
time_step_interval = 10
[]
[]
(modules/combined/test/tests/optimization/thermal_sensitivity/2d_root.i)
vol_frac = 0.5
E0 = 1
Emin = 1e-4
power = 1
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 2
nx = 20
ny = 20
xmin = 0
xmax = 40
ymin = 0
ymax = 40
[]
[]
[Variables]
[T]
initial_condition = 100
[]
[]
[AuxVariables]
[Dc]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[]
[Kernels]
[heat]
type = HeatConduction
diffusion_coefficient = k
variable = T
[]
[heat_source]
type = HeatSource
function = 1e-2
variable = T
[]
[]
[DiracKernels]
[src]
type = ConstantPointSource
variable = T
point = '0 5 0'
value = 10
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = T
boundary = 'right top bottom'
value = 0.0
[]
[]
[Materials]
[k]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = k
[]
[dc]
type = ThermalSensitivity
temperature = T
design_density = mat_den
thermal_conductivity = k
[]
#only needed for objective function output in postprocessor
[thermal_compliance]
type = ThermalCompliance
temperature = T
thermal_conductivity = k
[]
[]
[UserObjects]
[rad_avg]
type = RadialAverage
radius = 3
weights = linear
prop_name = thermal_sensitivity
execute_on = TIMESTEP_END
force_preaux = true
[]
[update]
type = DensityUpdate
density_sensitivity = Dc
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
[]
[calc_sense]
type = SensitivityFilter
density_sensitivity = Dc
design_density = mat_den
filter_UO = rad_avg
execute_on = TIMESTEP_END
force_postaux = true
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
nl_abs_tol = 1e-8
dt = 1.0
dtmin = 1.0
num_steps = 20
[]
[Outputs]
[out]
type = CSV
execute_on = 'FINAL'
[]
print_linear_residuals = false
[]
[Postprocessors]
[mesh_volume]
type = VolumePostprocessor
execute_on = 'initial timestep_end'
[]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[vol_frac]
type = ParsedPostprocessor
expression = 'total_vol / mesh_volume'
pp_names = 'total_vol mesh_volume'
[]
[sensitivity]
type = ElementIntegralMaterialProperty
mat_prop = thermal_sensitivity
[]
[objective_thermal]
type = ElementIntegralMaterialProperty
mat_prop = thermal_compliance
execute_on = 'INITIAL TIMESTEP_END'
[]
[]
(modules/combined/test/tests/optimization/optimization_density_update/top_opt_3d_pde_filter.i)
vol_frac = 0.4
E0 = 1e5
Emin = 1e-4
power = 2
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 3
nx = 24
ny = 12
nz = 12
xmin = 0
xmax = 20
ymin = 0
ymax = 10
zmin = 0
zmax = 10
[]
[middle_bottom_left_edge]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = pull
coord = '0 0 5'
[]
[]
[Variables]
[Dc]
initial_condition = -1.0
[]
[]
[AuxVariables]
[sensitivity]
family = MONOMIAL
order = FIRST
initial_condition = -1.0
[AuxKernel]
type = MaterialRealAux
variable = sensitivity
property = sensitivity
execute_on = LINEAR
[]
[]
[compliance]
family = MONOMIAL
order = CONSTANT
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[Dc_elem]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[AuxKernel]
type = SelfAux
variable = Dc_elem
v = Dc
execute_on = 'TIMESTEP_END'
[]
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[Kernels]
[diffusion]
type = FunctionDiffusion
variable = Dc
function = 0.05
[]
[potential]
type = Reaction
variable = Dc
[]
[source]
type = CoupledForce
variable = Dc
v = sensitivity
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[no_y]
type = DirichletBC
variable = disp_y
boundary = right
value = 0.0
[]
[no_z]
type = DirichletBC
variable = disp_z
boundary = right
value = 0.0
[]
[boundary_penalty]
type = ADRobinBC
variable = Dc
boundary = 'left top front back'
coefficient = 10
[]
[]
[NodalKernels]
[pull]
type = NodalGravity
variable = disp_y
boundary = pull
gravity_value = -1
mass = 1
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'mat_den'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys
incremental = false
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[update]
type = DensityUpdate
density_sensitivity = Dc_elem
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type '
petsc_options_value = 'lu'
nl_abs_tol = 1e-10
line_search = none
dt = 1.0
num_steps = 10
[]
[Outputs]
[out]
type = Exodus
time_step_interval = 10
[]
[]
(modules/combined/test/tests/optimization/compliance_sensitivity/2d_mbb_pde.i)
vol_frac = 0.5
E0 = 1
Emin = 1e-8
power = 3
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 2
nx = 150
ny = 50
xmin = 0
xmax = 30
ymin = 0
ymax = 10
[]
[node]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = hold_y
nodes = 0
[]
[push]
type = ExtraNodesetGenerator
input = node
new_boundary = push
coord = '30 10 0'
[]
[]
[Variables]
[disp_x]
[]
[disp_y]
[]
[Dc]
initial_condition = -1.0
[]
[]
[AuxVariables]
[sensitivity]
family = MONOMIAL
order = FIRST
initial_condition = -1.0
[AuxKernel]
type = MaterialRealAux
variable = sensitivity
property = sensitivity
execute_on = LINEAR
[]
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[Dc_elem]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[AuxKernel]
type = SelfAux
variable = Dc_elem
v = Dc
execute_on = 'TIMESTEP_END'
[]
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[Kernels]
[diffusion]
type = FunctionDiffusion
variable = Dc
function = 0.15 # radius coeff
[]
[potential]
type = Reaction
variable = Dc
[]
[source]
type = CoupledForce
variable = Dc
v = sensitivity
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = disp_y
boundary = hold_y
value = 0.0
[]
[no_y]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[boundary_penalty]
type = ADRobinBC
variable = Dc
boundary = 'left top'
coefficient = 10
[]
[boundary_penalty_right]
type = ADRobinBC
variable = Dc
boundary = 'right'
coefficient = 10
[]
[]
[NodalKernels]
[push]
type = NodalGravity
variable = disp_y
boundary = push
gravity_value = -1
mass = 1
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'mat_den'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys
incremental = false
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[update]
type = DensityUpdate
density_sensitivity = Dc_elem
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
force_postaux = true
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
line_search = none
nl_abs_tol = 1e-4
l_max_its = 200
start_time = 0.0
dt = 1.0
num_steps = 70
[]
[Outputs]
[out]
type = CSV
execute_on = 'INITIAL TIMESTEP_END'
[]
print_linear_residuals = false
[]
[Postprocessors]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[sensitivity]
type = ElementIntegralMaterialProperty
mat_prop = sensitivity
[]
[]
[Controls]
[first_period]
type = TimePeriod
start_time = 0.0
end_time = 10
enable_objects = 'BCs::boundary_penalty_right'
execute_on = 'initial timestep_begin'
[]
[]
(modules/combined/examples/optimization/helmholtz_multimat_strip.i)
vol_frac = 0.35
power = 1.1
Emin = 1.0e-6
Ess = 0.475 # ss
Et = 1.0 # w
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
# final_generator = 'MoveRight'
[Bottom]
type = GeneratedMeshGenerator
dim = 2
nx = 320
ny = 30
xmin = 0
xmax = 150
ymin = 0
ymax = 15
[]
[RenameBottom]
type = RenameBoundaryGenerator
input = Bottom
old_boundary = 'top bottom right left'
new_boundary = 'top_bottom bottom_bottom right_bottom left_bottom'
[]
[Middle]
type = GeneratedMeshGenerator
dim = 2
nx = 320
ny = 6
xmin = 0
xmax = 150
ymin = 0
ymax = 3
[]
[MoveMiddle]
type = TransformGenerator
input = Middle
transform = TRANSLATE
vector_value = '0 15 0'
[]
[RenameMiddle]
type = RenameBoundaryGenerator
input = MoveMiddle
old_boundary = 'top bottom right left'
new_boundary = 'top_middle bottom_middle right_middle left_middle'
[]
[Top]
type = GeneratedMeshGenerator
dim = 2
nx = 320
ny = 30
xmin = 0
xmax = 150
ymin = 0
ymax = 15
[]
[MoveTop]
type = TransformGenerator
input = Top
transform = TRANSLATE
vector_value = '0 18 0'
[]
[RenameTop]
type = RenameBoundaryGenerator
input = MoveTop
old_boundary = 'top bottom right left'
new_boundary = 'top_top bottom_top right_top left_top'
[]
[bottom_gen]
type = ParsedSubdomainMeshGenerator
input = RenameBottom
combinatorial_geometry = 'y <= 15'
block_id = 1
[]
[middle_gen]
type = ParsedSubdomainMeshGenerator
input = RenameMiddle
combinatorial_geometry = 'y <= 18 & y > 15'
block_id = 2
[]
[top_gen]
type = ParsedSubdomainMeshGenerator
input = RenameTop
combinatorial_geometry = 'y > 18'
block_id = 3
[]
[stitch]
type = StitchedMeshGenerator
inputs = 'bottom_gen middle_gen top_gen'
stitch_boundaries_pairs = 'top_bottom bottom_middle; top_middle bottom_top'
[]
[left_load]
type = ExtraNodesetGenerator
input = stitch
new_boundary = left_load
coord = '37.5 33 0'
[]
[right_load]
type = ExtraNodesetGenerator
input = left_load
new_boundary = right_load
coord = '112.5 33 0'
[]
[left_support]
type = ExtraNodesetGenerator
input = right_load
new_boundary = left_support
coord = '0 0 0'
[]
[right_support]
type = ExtraNodesetGenerator
input = left_support
new_boundary = right_support
coord = '150 0 0'
[]
[]
[Variables]
[disp_x]
[]
[disp_y]
[]
[Dc]
initial_condition = -1.0
[]
[]
[AuxVariables]
[Cc]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[sensitivity]
family = MONOMIAL
order = FIRST
initial_condition = -1.0
[AuxKernel]
type = MaterialRealAux
variable = sensitivity
property = sensitivity
execute_on = LINEAR
[]
block = '1 2 3'
[]
[mat_den_nodal]
family = L2_LAGRANGE
order = FIRST
initial_condition = ${vol_frac}
[AuxKernel]
type = SelfAux
execute_on = TIMESTEP_END
variable = mat_den_nodal
v = mat_den
[]
[]
[Dc_elem]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[AuxKernel]
type = SelfAux
variable = Dc_elem
v = Dc
execute_on = 'TIMESTEP_END'
[]
[]
[]
[Kernels]
[diffusion]
type = FunctionDiffusion
variable = Dc
function = 4.0
[]
[potential]
type = Reaction
variable = Dc
[]
[source]
type = CoupledForce
variable = Dc
v = sensitivity
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[BCs]
[no_y]
type = DirichletBC
variable = disp_y
boundary = left_support
value = 0.0
[]
[no_x]
type = DirichletBC
variable = disp_x
boundary = left_support
value = 0.0
[]
[no_y_right]
type = DirichletBC
variable = disp_y
boundary = right_support
value = 0.0
[]
[boundary_penalty]
type = ADRobinBC
variable = Dc
boundary = 'bottom_bottom right_bottom left_bottom top_top right_top left_top left_middle '
'right_middle'
coefficient = 10
[]
[]
[NodalKernels]
[left_down]
type = NodalGravity
variable = disp_y
boundary = left_load
gravity_value = -1e-3
mass = 1
[]
[right_down]
type = NodalGravity
variable = disp_y
boundary = right_load
gravity_value = -1e-3
mass = 1
[]
[]
[Materials]
[sensitivity]
type = ParsedMaterial
property_name = 'sensitivity'
block = '2'
expression = '0'
[]
[elasticity_tensor_one]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys_one
poissons_ratio = poissons_ratio
args = 'mat_den'
block = '1'
[]
[elasticity_tensor_three]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys_three
poissons_ratio = poissons_ratio
args = 'mat_den'
block = '3'
[]
[elasticity_tensor_two]
type = ComputeIsotropicElasticityTensor
youngs_modulus = 1.0
poissons_ratio = 0.3
block = '2'
[]
# One: Tungsten
[E_phys_one]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${Et}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys_one
block = '1'
outputs = 'exodus'
[]
# Three: SS316
[E_phys_three]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${Ess}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys_three
block = '3'
outputs = 'exodus'
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc_one]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys_one
block = '1'
[]
[dc_three]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys_three
block = '3'
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[update_one]
type = DensityUpdate
density_sensitivity = Dc_elem
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
block = '1'
[]
[update_three]
type = DensityUpdate
density_sensitivity = Dc_elem
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
block = '3'
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
nl_abs_tol = 1e-10
dt = 1.0
num_steps = 90
[]
[Outputs]
exodus = true
[out]
type = CSV
execute_on = 'TIMESTEP_END'
[]
print_linear_residuals = false
[]
[Postprocessors]
[mesh_volume]
type = VolumePostprocessor
execute_on = 'initial timestep_end'
[]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[vol_frac]
type = ParsedPostprocessor
expression = 'total_vol / mesh_volume'
pp_names = 'total_vol mesh_volume'
[]
[sensitivity]
type = ElementIntegralMaterialProperty
mat_prop = sensitivity
block = '1 3'
[]
[objective_one]
type = ElementIntegralMaterialProperty
mat_prop = strain_energy_density
execute_on = 'INITIAL TIMESTEP_END'
block = '1'
[]
[objective_three]
type = ElementIntegralMaterialProperty
mat_prop = strain_energy_density
execute_on = 'INITIAL TIMESTEP_END'
block = '3'
[]
[]
(modules/combined/examples/optimization/2d_mbb_pde_amr.i)
vol_frac = 0.5
E0 = 1
Emin = 1e-8
power = 2
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 2
nx = 30
ny = 10
xmin = 0
xmax = 30
ymin = 0
ymax = 10
[]
[node]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = pull
nodes = 0
[]
[push]
type = ExtraNodesetGenerator
input = node
new_boundary = push
coord = '30 10 0'
[]
[]
[Variables]
[Dc]
initial_condition = -1.0
[]
[]
[AuxVariables]
[Emin]
family = MONOMIAL
order = CONSTANT
initial_condition = ${Emin}
[]
[power]
family = MONOMIAL
order = CONSTANT
initial_condition = ${power}
[]
[E0]
family = MONOMIAL
order = CONSTANT
initial_condition = ${E0}
[]
[sensitivity]
family = MONOMIAL
order = FIRST
initial_condition = -1.0
[AuxKernel]
type = MaterialRealAux
variable = sensitivity
property = sensitivity
execute_on = LINEAR
[]
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[Dc_elem]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[AuxKernel]
type = SelfAux
variable = Dc_elem
v = Dc
execute_on = 'TIMESTEP_END'
[]
[]
[mat_den_nodal]
family = L2_LAGRANGE
order = FIRST
initial_condition = ${vol_frac}
[AuxKernel]
type = SelfAux
execute_on = TIMESTEP_END
variable = mat_den_nodal
v = mat_den
[]
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[Kernels]
[diffusion]
type = FunctionDiffusion
variable = Dc
function = 0.15 # radius coeff
[]
[potential]
type = Reaction
variable = Dc
[]
[source]
type = CoupledForce
variable = Dc
v = sensitivity
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = disp_y
boundary = pull
value = 0.0
[]
[no_y]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[boundary_penalty]
type = ADRobinBC
variable = Dc
boundary = 'left top'
coefficient = 10
[]
[boundary_penalty_right]
type = ADRobinBC
variable = Dc
boundary = 'right'
coefficient = 10
[]
[]
[NodalKernels]
[pull]
type = NodalGravity
variable = disp_y
boundary = push
gravity_value = -1
mass = 1
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'Emin mat_den power E0'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys
incremental = false
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[update]
type = DensityUpdate
density_sensitivity = Dc_elem
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
force_postaux = true
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
line_search = none
nl_abs_tol = 1e-4
l_max_its = 200
start_time = 0.0
dt = 1.0
num_steps = 70
[]
[Outputs]
[out]
type = Exodus
execute_on = 'INITIAL TIMESTEP_END'
[]
print_linear_residuals = false
[]
[Postprocessors]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[]
[Controls]
[first_period]
type = TimePeriod
start_time = 0.0
end_time = 40
enable_objects = 'BCs::boundary_penalty_right'
execute_on = 'initial timestep_begin'
[]
[]
[Adaptivity]
max_h_level = 2
recompute_markers_during_cycles = true
interval = 1
cycles_per_step = 1
marker = density_marker
[Indicators]
[density_jump]
type = ValueJumpIndicator
variable = mat_den_nodal
[]
[]
[Markers]
[density_marker]
type = ErrorToleranceMarker
indicator = density_jump
coarsen = 0.1
refine = 0.1
[]
[]
[]
(modules/optimization/test/tests/simp/2d.i)
vol_frac = 0.2
[Mesh]
[planet]
type = ConcentricCircleMeshGenerator
has_outer_square = false
radii = 1
num_sectors = 10
rings = 2
preserve_volumes = false
[]
[moon]
type = ConcentricCircleMeshGenerator
has_outer_square = false
radii = 0.5
num_sectors = 8
rings = 2
preserve_volumes = false
[]
[combine]
type = CombinerGenerator
inputs = 'planet moon'
positions = '0 0 0 -1.5 -0.5 0'
[]
[]
[AuxVariables]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = 0.1
[]
[Dc]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[]
[]
[Variables]
[u]
[]
[v]
[]
[]
[Kernels]
[diff_u]
type = Diffusion
variable = u
[]
[dt_u]
type = TimeDerivative
variable = u
[]
[diff_v]
type = Diffusion
variable = v
[]
[dt_v]
type = TimeDerivative
variable = v
[]
[]
[Materials]
[thermal_cond]
type = GenericFunctionMaterial
prop_values = '-1.4*abs(y)-2.7*abs(x)'
prop_names = thermal_cond
outputs = 'exodus'
[]
[thermal_compliance_sensitivity]
type = GenericFunctionMaterial
prop_values = '-3*abs(y)-1.5*abs(x)'
prop_names = thermal_sensitivity
outputs = 'exodus'
[]
[]
[BCs]
[flux_u]
type = DirichletBC
variable = u
boundary = outer
value = 3.0
[]
[flux_v]
type = DirichletBC
variable = v
boundary = outer
value = 7.0
[]
[]
[UserObjects]
[rad_avg]
type = RadialAverage
radius = 0.1
weights = linear
prop_name = thermal_sensitivity
execute_on = TIMESTEP_END
force_preaux = true
[]
[update]
type = DensityUpdate
density_sensitivity = Dc
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
[]
[calc_sense]
type = SensitivityFilter
density_sensitivity = Dc
design_density = mat_den
filter_UO = rad_avg
execute_on = TIMESTEP_END
force_postaux = true
[]
[]
[Executioner]
type = Transient
dt = 0.1
num_steps = 15
nl_rel_tol = 1e-04
[]
[Outputs]
exodus = true
[]
(modules/combined/test/tests/optimization/compliance_sensitivity/3d_mbb.i)
vol_frac = 0.5
E0 = 1
Emin = 1e-8
power = 3
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 3
nx = 30
ny = 10
nz = 10
xmin = 0
xmax = 30
ymin = 0
ymax = 10
zmin = 0
zmax = 10
[]
[node]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = hold_y
coord = '0 0 0; 0 0 10'
[]
[push]
type = ExtraNodesetGenerator
input = node
new_boundary = push
coord = '30 10 5'
[]
[]
[Variables]
[disp_x]
[]
[disp_y]
[]
[disp_z]
[]
[Dc]
initial_condition = -1.0
[]
[]
[AuxVariables]
[sensitivity]
family = MONOMIAL
order = FIRST
initial_condition = -1.0
[AuxKernel]
type = MaterialRealAux
variable = sensitivity
property = sensitivity
execute_on = LINEAR
[]
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[Dc_elem]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[AuxKernel]
type = SelfAux
variable = Dc_elem
v = Dc
execute_on = 'TIMESTEP_END'
[]
[]
[mat_den_nodal]
family = L2_LAGRANGE
order = FIRST
initial_condition = ${vol_frac}
[AuxKernel]
type = SelfAux
execute_on = TIMESTEP_END
variable = mat_den_nodal
v = mat_den
[]
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[Kernels]
[diffusion]
type = FunctionDiffusion
variable = Dc
function = 0.15 # radius coeff
[]
[potential]
type = Reaction
variable = Dc
[]
[source]
type = CoupledForce
variable = Dc
v = sensitivity
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = disp_y
boundary = hold_y
value = 0.0
[]
[no_y]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[boundary_penalty]
type = ADRobinBC
variable = Dc
boundary = 'left top front back'
coefficient = 10
[]
[boundary_penalty_right]
type = ADRobinBC
variable = Dc
boundary = 'right'
coefficient = 10
[]
[]
[NodalKernels]
[push]
type = NodalGravity
variable = disp_y
boundary = push
gravity_value = -1
mass = 1
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'mat_den'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys
incremental = false
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[update]
type = DensityUpdate
density_sensitivity = Dc_elem
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
force_postaux = true
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
line_search = none
nl_abs_tol = 1e-4
l_max_its = 200
start_time = 0.0
dt = 1.0
num_steps = 2
[]
[Outputs]
[out]
type = CSV
execute_on = 'INITIAL TIMESTEP_END'
[]
print_linear_residuals = false
[]
[Postprocessors]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[sensitivity]
type = ElementIntegralMaterialProperty
mat_prop = sensitivity
[]
[]
[Controls]
[first_period]
type = TimePeriod
start_time = 0.0
end_time = 10
enable_objects = 'BCs::boundary_penalty_right'
execute_on = 'initial timestep_begin'
[]
[]
(modules/combined/examples/optimization/thermomechanical/thermomechanical_main.i)
vol_frac = 0.4
power = 2.0
E0 = 1.0e-6
E1 = 1.0
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 2
nx = 40
ny = 40
xmin = 0
xmax = 40
ymin = 0
ymax = 40
[]
[node]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = hold
nodes = 0
[]
[push_left]
type = ExtraNodesetGenerator
input = node
new_boundary = push_left
coord = '16 0 0'
[]
[push_center]
type = ExtraNodesetGenerator
input = push_left
new_boundary = push_center
coord = '24 0 0'
[]
[extra]
type = SideSetsFromBoundingBoxGenerator
input = push_center
bottom_left = '-0.01 17.999 0'
top_right = '5 22.001 0'
boundary_new = n1
boundaries_old = left
[]
[dirichlet_bc]
type = SideSetsFromNodeSetsGenerator
input = extra
[]
[]
[Variables]
[disp_x]
[]
[disp_y]
[]
[]
[AuxVariables]
[mat_den]
family = MONOMIAL
order = FIRST
initial_condition = 0.02
[]
[sensitivity_one]
family = MONOMIAL
order = FIRST
initial_condition = -1.0
[]
[sensitivity_two]
family = MONOMIAL
order = FIRST
initial_condition = -1.0
[]
[total_sensitivity]
family = MONOMIAL
order = FIRST
initial_condition = -1.0
[]
[]
[AuxKernels]
[total_sensitivity]
type = ParsedAux
variable = total_sensitivity
expression = '(1-1.0e-7)*sensitivity_one + 1.0e-7*sensitivity_two'
coupled_variables = 'sensitivity_one sensitivity_two'
execute_on = 'LINEAR TIMESTEP_END'
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[BCs]
[no_y]
type = DirichletBC
variable = disp_y
boundary = hold
value = 0.0
[]
[no_x_symm]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'mat_den'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${E1} + (mat_den ^ ${power}) * (${E1}-${E0})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
# We do filtering in the subapps
[update]
type = DensityUpdate
density_sensitivity = total_sensitivity
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = MULTIAPP_FIXED_POINT_BEGIN
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
nl_abs_tol = 1e-8
dt = 1.0
num_steps = 2
[]
[Outputs]
[out]
type = CSV
execute_on = 'TIMESTEP_END'
[]
print_linear_residuals = false
exodus = true
[]
[Postprocessors]
[mesh_volume]
type = VolumePostprocessor
execute_on = 'initial timestep_end'
[]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[vol_frac]
type = ParsedPostprocessor
expression = 'total_vol / mesh_volume'
pp_names = 'total_vol mesh_volume'
[]
[sensitivity]
type = ElementIntegralVariablePostprocessor
variable = total_sensitivity
[]
[]
[MultiApps]
[sub_app_one]
type = TransientMultiApp
input_files = structural_sub.i
[]
[sub_app_two]
type = TransientMultiApp
input_files = thermal_sub.i
[]
[]
[Transfers]
# First SUB-APP: STRUCTURAL
# To subapp densities
[subapp_one_density]
type = MultiAppCopyTransfer
to_multi_app = sub_app_one
source_variable = mat_den # Here
variable = mat_den
[]
# From subapp sensitivity
[subapp_one_sensitivity]
type = MultiAppCopyTransfer
from_multi_app = sub_app_one
source_variable = Dc # sensitivity_var
variable = sensitivity_one # Here
[]
# Second SUB-APP: HEAT CONDUCTIVITY
# To subapp densities
[subapp_two_density]
type = MultiAppCopyTransfer
to_multi_app = sub_app_two
source_variable = mat_den # Here
variable = mat_den
[]
# From subapp sensitivity
[subapp_two_sensitivity]
type = MultiAppCopyTransfer
from_multi_app = sub_app_two
source_variable = Tc # sensitivity_var
variable = sensitivity_two # Here
[]
[]
(modules/combined/examples/optimization/multi-load/single_main.i)
vol_frac = 0.3
power = 1.1
E0 = 1.0
Emin = 1.0e-6
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
# final_generator = 'MoveRight'
[Bottom]
type = GeneratedMeshGenerator
dim = 2
nx = 80
ny = 40
xmin = 0
xmax = 150
ymin = 0
ymax = 75
[]
[left_load]
type = ExtraNodesetGenerator
input = Bottom
new_boundary = left_load
coord = '37.5 75 0'
[]
[right_load]
type = ExtraNodesetGenerator
input = left_load
new_boundary = right_load
coord = '112.5 75 0'
[]
[left_support]
type = ExtraNodesetGenerator
input = right_load
new_boundary = left_support
coord = '0 0 0'
[]
[right_support]
type = ExtraNodesetGenerator
input = left_support
new_boundary = right_support
coord = '150 0 0'
[]
[]
[Variables]
[disp_x]
[]
[disp_y]
[]
[]
[AuxVariables]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = 0.02
[]
[sensitivity_one]
family = MONOMIAL
order = SECOND
initial_condition = -1.0
[]
[sensitivity_two]
family = MONOMIAL
order = SECOND
initial_condition = -1.0
[]
[total_sensitivity]
family = MONOMIAL
order = SECOND
initial_condition = -1.0
[]
[]
[AuxKernels]
[total_sensitivity]
type = ParsedAux
variable = total_sensitivity
expression = '0.5*sensitivity_one + 0.5*sensitivity_two'
coupled_variables = 'sensitivity_one sensitivity_two'
execute_on = 'LINEAR TIMESTEP_END'
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[BCs]
[no_y]
type = DirichletBC
variable = disp_y
boundary = left_support
value = 0.0
[]
[no_x]
type = DirichletBC
variable = disp_x
boundary = left_support
value = 0.0
[]
[no_y_right]
type = DirichletBC
variable = disp_y
boundary = right_support
value = 0.0
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'mat_den'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.0
[]
[stress]
type = ComputeLinearElasticStress
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
# We do filtering in the subapps
[update]
type = DensityUpdate
density_sensitivity = total_sensitivity
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = MULTIAPP_FIXED_POINT_BEGIN
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
nl_abs_tol = 1e-10
dt = 1.0
num_steps = 25
[]
[Outputs]
exodus = true
[out]
type = CSV
execute_on = 'TIMESTEP_END'
[]
print_linear_residuals = false
[]
[Postprocessors]
[mesh_volume]
type = VolumePostprocessor
execute_on = 'initial timestep_end'
[]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[vol_frac]
type = ParsedPostprocessor
expression = 'total_vol / mesh_volume'
pp_names = 'total_vol mesh_volume'
[]
[sensitivity]
type = ElementIntegralVariablePostprocessor
variable = total_sensitivity
[]
[]
[MultiApps]
[sub_app_one]
type = TransientMultiApp
input_files = single_subapp_one.i
[]
[sub_app_two]
type = TransientMultiApp
input_files = single_subapp_two.i
[]
[]
[Transfers]
# First SUB-APP
# To subapp densities
[subapp_one_density]
type = MultiAppCopyTransfer
to_multi_app = sub_app_one
source_variable = mat_den # Here
variable = mat_den
[]
# From subapp sensitivity
[subapp_one_sensitivity]
type = MultiAppCopyTransfer
from_multi_app = sub_app_one
source_variable = Dc # sensitivity_var
variable = sensitivity_one # Here
[]
# Second SUB-APP
# To subapp densities
[subapp_two_density]
type = MultiAppCopyTransfer
to_multi_app = sub_app_two
source_variable = mat_den # Here
variable = mat_den
[]
# From subapp sensitivity
[subapp_two_sensitivity]
type = MultiAppCopyTransfer
from_multi_app = sub_app_two
source_variable = Dc # sensitivity_var
variable = sensitivity_two # Here
[]
[]
(modules/combined/test/tests/optimization/compliance_sensitivity/2d_mbb_pde_amr.i)
vol_frac = 0.5
E0 = 1
Emin = 1e-8
power = 3
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 2
nx = 30
ny = 10
xmin = 0
xmax = 30
ymin = 0
ymax = 10
[]
[node]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = hold
nodes = 0
[]
[push]
type = ExtraNodesetGenerator
input = node
new_boundary = push
coord = '30 10 0'
[]
[]
[Variables]
[disp_x]
[]
[disp_y]
[]
[Dc]
initial_condition = -1.0
[]
[]
[AuxVariables]
[sensitivity]
family = MONOMIAL
order = FIRST
initial_condition = -1.0
[AuxKernel]
type = MaterialRealAux
variable = sensitivity
property = sensitivity
execute_on = LINEAR
[]
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[Dc_elem]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[AuxKernel]
type = SelfAux
variable = Dc_elem
v = Dc
execute_on = 'TIMESTEP_END'
[]
[]
[mat_den_nodal]
family = L2_LAGRANGE
order = FIRST
initial_condition = ${vol_frac}
[AuxKernel]
type = SelfAux
execute_on = TIMESTEP_END
variable = mat_den_nodal
v = mat_den
[]
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[Kernels]
[diffusion]
type = FunctionDiffusion
variable = Dc
function = 0.15 # radius coeff
[]
[potential]
type = Reaction
variable = Dc
[]
[source]
type = CoupledForce
variable = Dc
v = sensitivity
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = disp_y
boundary = hold
value = 0.0
[]
[no_y]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[boundary_penalty]
type = ADRobinBC
variable = Dc
boundary = 'left top'
coefficient = 10
[]
[boundary_penalty_right]
type = ADRobinBC
variable = Dc
boundary = 'right'
coefficient = 10
[]
[]
[NodalKernels]
[push]
type = NodalGravity
variable = disp_y
boundary = push
gravity_value = -1
mass = 1
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'mat_den'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys
incremental = false
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[update]
type = DensityUpdate
density_sensitivity = Dc_elem
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
force_postaux = true
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
line_search = none
nl_abs_tol = 1e-4
l_max_its = 200
start_time = 0.0
dt = 1.0
num_steps = 70
[]
[Outputs]
[out]
type = CSV
execute_on = 'INITIAL TIMESTEP_END'
[]
print_linear_residuals = false
[]
[Postprocessors]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[sensitivity]
type = ElementIntegralMaterialProperty
mat_prop = sensitivity
[]
[]
[Controls]
[first_period]
type = TimePeriod
start_time = 0.0
end_time = 40
enable_objects = 'BCs::boundary_penalty_right'
execute_on = 'initial timestep_begin'
[]
[]
[Adaptivity]
max_h_level = 2
recompute_markers_during_cycles = true
interval = 1
cycles_per_step = 1
marker = density_marker
[Indicators]
[density_jump]
type = ValueJumpIndicator
variable = mat_den_nodal
[]
[]
[Markers]
[density_marker]
type = ErrorToleranceMarker
indicator = density_jump
coarsen = 0.1
refine = 0.1
[]
[]
[]
(modules/combined/test/tests/optimization/optimization_density_update/top_opt_3d.i)
vol_frac = 0.5
E0 = 1e5
Emin = 1e-2
power = 2
[GlobalParams]
displacements = 'disp_x disp_y disp_z'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 3
nx = 24
ny = 12
nz = 12
xmin = 0
xmax = 20
ymin = 0
ymax = 10
zmin = 0
zmax = 10
[]
[middle_bottom_left_edge]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = pull
coord = '0 0 5'
[]
[]
[AuxVariables]
[compliance]
family = MONOMIAL
order = CONSTANT
[]
[Dc]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[no_y]
type = DirichletBC
variable = disp_y
boundary = right
value = 0.0
[]
[no_z]
type = DirichletBC
variable = disp_z
boundary = right
value = 0.0
[]
[]
[NodalKernels]
[pull]
type = NodalGravity
variable = disp_y
boundary = pull
gravity_value = -1
mass = 1
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'mat_den'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys
incremental = false
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[rad_avg]
type = RadialAverage
radius = 0.5
weights = constant
prop_name = sensitivity
execute_on = TIMESTEP_END
execution_order_group = -1
[]
[update]
type = DensityUpdate
density_sensitivity = Dc
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
[]
[calc_sense]
type = SensitivityFilter
density_sensitivity = Dc
design_density = mat_den
filter_UO = rad_avg
execute_on = TIMESTEP_END
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu '
nl_abs_tol = 1e-10
l_max_its = 200
start_time = 0.0
dt = 1.0
num_steps = 10
[]
[Outputs]
[out]
type = Exodus
time_step_interval = 10
[]
[]
(modules/combined/test/tests/optimization/compliance_sensitivity/2d_mmb_2material.i)
vol_frac = 0.5
power = 1
E0 = 1e-5
E1 = 0.6
E2 = 1.0
rho0 = 0.0
rho1 = 0.4
rho2 = 1.0
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 2
nx = 150
ny = 50
xmin = 0
xmax = 30
ymin = 0
ymax = 10
[]
[node]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = hold
nodes = 0
[]
[push]
type = ExtraNodesetGenerator
input = node
new_boundary = push
coord = '30 10 0'
[]
[]
[Variables]
[disp_x]
[]
[disp_y]
[]
[]
[AuxVariables]
[Dc]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
# initial_condition = ${vol_frac}
[]
[]
[ICs]
[mat_den]
type = RandomIC
seed = 5
variable = mat_den
max = '${fparse vol_frac+0.15}'
min = '${fparse vol_frac-0.15}'
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = disp_y
boundary = hold
value = 0.0
[]
[no_y]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[]
[NodalKernels]
[push]
type = NodalGravity
variable = disp_y
boundary = push
gravity_value = -1
mass = 1
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'mat_den'
[]
[E_phys]
type = DerivativeParsedMaterial
# ordered multimaterial simp
expression = "A1:=(${E0}-${E1})/(${rho0}^${power}-${rho1}^${power}); "
"B1:=${E0}-A1*${rho0}^${power}; E1:=A1*mat_den^${power}+B1; "
"A2:=(${E1}-${E2})/(${rho1}^${power}-${rho2}^${power}); "
"B2:=${E1}-A2*${rho1}^${power}; E2:=A2*mat_den^${power}+B2; "
"if(mat_den<${rho1},E1,E2)"
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc]
type = ComplianceSensitivity2
design_density = mat_den
youngs_modulus = E_phys
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[rad_avg]
type = RadialAverage
radius = 1.2
weights = linear
prop_name = sensitivity
execute_on = TIMESTEP_END
force_preaux = true
[]
[update]
type = DensityUpdate
density_sensitivity = Dc
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
[]
[calc_sense]
type = SensitivityFilter
density_sensitivity = Dc
design_density = mat_den
filter_UO = rad_avg
execute_on = TIMESTEP_END
force_postaux = true
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
nl_abs_tol = 1e-8
dt = 1.0
num_steps = 70
[]
[Outputs]
exodus = true
[out]
type = CSV
execute_on = 'TIMESTEP_END'
[]
print_linear_residuals = false
[]
[Postprocessors]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[sensitivity]
type = ElementIntegralMaterialProperty
mat_prop = sensitivity
[]
[]
(modules/combined/examples/optimization/2d_mbb_pde.i)
vol_frac = 0.5
E0 = 1
Emin = 1e-8
power = 2
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 2
nx = 150
ny = 50
xmin = 0
xmax = 30
ymin = 0
ymax = 10
[]
[node]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = hold_y
nodes = 0
[]
[push]
type = ExtraNodesetGenerator
input = node
new_boundary = push
coord = '30 10 0'
[]
[]
[Variables]
[Dc]
initial_condition = -1.0
[]
[]
[AuxVariables]
[Emin]
family = MONOMIAL
order = CONSTANT
initial_condition = ${Emin}
[]
[power]
family = MONOMIAL
order = CONSTANT
initial_condition = ${power}
[]
[E0]
family = MONOMIAL
order = CONSTANT
initial_condition = ${E0}
[]
[sensitivity]
family = MONOMIAL
order = FIRST
initial_condition = -1.0
[AuxKernel]
type = MaterialRealAux
variable = sensitivity
property = sensitivity
execute_on = LINEAR
[]
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[Dc_elem]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[AuxKernel]
type = SelfAux
variable = Dc_elem
v = Dc
execute_on = 'TIMESTEP_END'
[]
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[Kernels]
[diffusion]
type = FunctionDiffusion
variable = Dc
function = 0.15 # radius coeff
[]
[potential]
type = Reaction
variable = Dc
[]
[source]
type = CoupledForce
variable = Dc
v = sensitivity
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = disp_y
boundary = hold_y
value = 0.0
[]
[no_y]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[boundary_penalty]
type = ADRobinBC
variable = Dc
boundary = 'left top'
coefficient = 10
[]
[boundary_penalty_right]
type = ADRobinBC
variable = Dc
boundary = 'right'
coefficient = 10
[]
[]
[NodalKernels]
[push]
type = NodalGravity
variable = disp_y
boundary = push
gravity_value = -1
mass = 1
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'Emin mat_den power E0'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys
incremental = false
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[update]
type = DensityUpdate
density_sensitivity = Dc_elem
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
force_postaux = true
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
line_search = none
nl_abs_tol = 1e-4
l_max_its = 200
start_time = 0.0
dt = 1.0
num_steps = 70
[]
[Outputs]
[out]
type = Exodus
execute_on = 'INITIAL TIMESTEP_END'
[]
print_linear_residuals = false
[]
[Postprocessors]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[]
[Controls]
[first_period]
type = TimePeriod
start_time = 0.0
end_time = 10
enable_objects = 'BCs::boundary_penalty_right'
execute_on = 'initial timestep_begin'
[]
[]
(modules/combined/examples/optimization/helmholtz_multimat_nostrip.i)
vol_frac = 0.35
power = 1.1
Emin = 1.0e-6
Ess = 0.475 # ss
Et = 1.0 # w
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
# final_generator = 'MoveRight'
[Bottom]
type = GeneratedMeshGenerator
dim = 2
nx = 320
ny = 30
xmin = 0
xmax = 150
ymin = 0
ymax = 15
[]
[RenameBottom]
type = RenameBoundaryGenerator
input = Bottom
old_boundary = 'top bottom right left'
new_boundary = 'top_bottom bottom_bottom right_bottom left_bottom'
[]
[Top]
type = GeneratedMeshGenerator
dim = 2
nx = 320
ny = 30
xmin = 0
xmax = 150
ymin = 0
ymax = 15
[]
[MoveTop]
type = TransformGenerator
input = Top
transform = TRANSLATE
vector_value = '0 15 0'
[]
[RenameTop]
type = RenameBoundaryGenerator
input = MoveTop
old_boundary = 'top bottom right left'
new_boundary = 'top_top bottom_top right_top left_top'
[]
[bottom_gen]
type = ParsedSubdomainMeshGenerator
input = RenameBottom
combinatorial_geometry = 'y <= 15'
block_id = 1
[]
[top_gen]
type = ParsedSubdomainMeshGenerator
input = RenameTop
combinatorial_geometry = 'y > 15'
block_id = 3
[]
[stitch]
type = StitchedMeshGenerator
inputs = 'bottom_gen top_gen'
stitch_boundaries_pairs = 'top_bottom bottom_top'
[]
[left_load]
type = ExtraNodesetGenerator
input = stitch
new_boundary = left_load
coord = '37.5 30 0'
[]
[right_load]
type = ExtraNodesetGenerator
input = left_load
new_boundary = right_load
coord = '112.5 30 0'
[]
[left_support]
type = ExtraNodesetGenerator
input = right_load
new_boundary = left_support
coord = '0 0 0'
[]
[right_support]
type = ExtraNodesetGenerator
input = left_support
new_boundary = right_support
coord = '150 0 0'
[]
[]
[Variables]
[disp_x]
[]
[disp_y]
[]
[Dc]
initial_condition = -1.0
[]
[]
[AuxVariables]
[Cc]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[sensitivity]
family = MONOMIAL
order = FIRST
initial_condition = -1.0
[AuxKernel]
type = MaterialRealAux
variable = sensitivity
property = sensitivity
execute_on = LINEAR
[]
[]
[mat_den_nodal]
family = L2_LAGRANGE
order = FIRST
initial_condition = ${vol_frac}
[AuxKernel]
type = SelfAux
execute_on = TIMESTEP_END
variable = mat_den_nodal
v = mat_den
[]
[]
[Dc_elem]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[AuxKernel]
type = SelfAux
variable = Dc_elem
v = Dc
execute_on = 'TIMESTEP_END'
[]
[]
[]
[Kernels]
[diffusion]
type = FunctionDiffusion
variable = Dc
function = 4.0
[]
[potential]
type = Reaction
variable = Dc
[]
[source]
type = CoupledForce
variable = Dc
v = sensitivity
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[BCs]
[no_y]
type = DirichletBC
variable = disp_y
boundary = left_support
value = 0.0
[]
[no_x]
type = DirichletBC
variable = disp_x
boundary = left_support
value = 0.0
[]
[no_y_right]
type = DirichletBC
variable = disp_y
boundary = right_support
value = 0.0
[]
[boundary_penalty]
type = ADRobinBC
variable = Dc
boundary = 'bottom_bottom right_bottom left_bottom top_top right_top left_top'
coefficient = 10
[]
[]
[NodalKernels]
[left_down]
type = NodalGravity
variable = disp_y
boundary = left_load
gravity_value = -1e-3
mass = 1
[]
[right_down]
type = NodalGravity
variable = disp_y
boundary = right_load
gravity_value = -1e-3
mass = 1
[]
[]
[Materials]
[elasticity_tensor_one]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys_one
poissons_ratio = poissons_ratio
args = 'mat_den'
block = '1'
[]
[elasticity_tensor_three]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys_three
poissons_ratio = poissons_ratio
args = 'mat_den'
block = '3'
[]
# One: Tungsten
[E_phys_one]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${Et}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys_one
block = '1'
outputs = 'exodus'
[]
# Three: SS316
[E_phys_three]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${Ess}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys_three
block = '3'
outputs = 'exodus'
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc_one]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys_one
block = '1'
[]
[dc_three]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys_three
block = '3'
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[update_one]
type = DensityUpdate
density_sensitivity = Dc_elem
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
block = '1'
[]
[update_three]
type = DensityUpdate
density_sensitivity = Dc_elem
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
block = '3'
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
nl_abs_tol = 1e-10
dt = 1.0
num_steps = 90
[]
[Outputs]
exodus = true
[out]
type = CSV
execute_on = 'TIMESTEP_END'
[]
print_linear_residuals = false
[]
[Postprocessors]
[mesh_volume]
type = VolumePostprocessor
execute_on = 'initial timestep_end'
[]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[vol_frac]
type = ParsedPostprocessor
expression = 'total_vol / mesh_volume'
pp_names = 'total_vol mesh_volume'
[]
[sensitivity]
type = ElementIntegralMaterialProperty
mat_prop = sensitivity
block = '1 3'
[]
[objective_one]
type = ElementIntegralMaterialProperty
mat_prop = strain_energy_density
execute_on = 'INITIAL TIMESTEP_END'
block = '1'
[]
[objective_three]
type = ElementIntegralMaterialProperty
mat_prop = strain_energy_density
execute_on = 'INITIAL TIMESTEP_END'
block = '3'
[]
[]
(modules/combined/examples/optimization/multi-load/square_main.i)
# This example is intended to reproduce a 2D example with opposing horizontal
# loads (see [1]). This test has an undefined solution if reguar SIMP is applied.
# Using multi-loads SIMP, on the other hand, generates a structure that optimizes
# the response to both loads individually,
# [1]. Lat. Am. j. solids struct. 12 (5), May 2015
# Topological derivative-based topology optimization of structures subject to multiple load-cases
vol_frac = 0.5
power = 1.0
E0 = 1.0
Emin = 1.0e-6
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
[Bottom]
type = GeneratedMeshGenerator
dim = 2
nx = 100
ny = 100
xmin = 0
xmax = 150
ymin = 0
ymax = 150
[]
[left_load]
type = ExtraNodesetGenerator
input = Bottom
new_boundary = left_load
coord = '0 150 0'
[]
[right_load]
type = ExtraNodesetGenerator
input = left_load
new_boundary = right_load
coord = '150 150 0'
[]
[left_support]
type = ExtraNodesetGenerator
input = right_load
new_boundary = left_support
coord = '0 0 0'
[]
[right_support]
type = ExtraNodesetGenerator
input = left_support
new_boundary = right_support
coord = '150 0 0'
[]
[]
[Variables]
[disp_x]
[]
[disp_y]
[]
[]
[AuxVariables]
[mat_den]
family = MONOMIAL
order = CONSTANT
[]
[sensitivity_one]
family = MONOMIAL
order = SECOND
initial_condition = -1.0
[]
[sensitivity_two]
family = MONOMIAL
order = SECOND
initial_condition = -1.0
[]
[total_sensitivity]
family = MONOMIAL
order = SECOND
initial_condition = -1.0
[]
[]
[ICs]
[mat_den]
type = RandomIC
seed = 7
variable = mat_den
max = '${fparse vol_frac+0.35}'
min = '${fparse vol_frac-0.35}'
[]
[]
[AuxKernels]
[total_sensitivity]
type = ParsedAux
variable = total_sensitivity
expression = '0.5*sensitivity_one + 0.5*sensitivity_two'
coupled_variables = 'sensitivity_one sensitivity_two'
execute_on = 'LINEAR TIMESTEP_END'
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[BCs]
[no_y]
type = DirichletBC
variable = disp_y
boundary = left_support
value = 0.0
[]
[no_x]
type = DirichletBC
variable = disp_x
boundary = left_support
value = 0.0
[]
[no_y_right]
type = DirichletBC
variable = disp_y
boundary = right_support
value = 0.0
[]
[no_x_right]
type = DirichletBC
variable = disp_x
boundary = right_support
value = 0.0
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'mat_den'
[]
[E_phys]
type = DerivativeParsedMaterial
# Emin + (density^penal) * (E0 - Emin)
expression = '${Emin} + (mat_den ^ ${power}) * (${E0}-${Emin})'
coupled_variables = 'mat_den'
property_name = E_phys
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[update]
type = DensityUpdate
density_sensitivity = total_sensitivity
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = MULTIAPP_FIXED_POINT_BEGIN
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
nl_abs_tol = 1e-10
dt = 1.0
num_steps = 10
[]
[Outputs]
exodus = true
[out]
type = CSV
execute_on = 'TIMESTEP_END'
[]
print_linear_residuals = false
[]
[Postprocessors]
[mesh_volume]
type = VolumePostprocessor
execute_on = 'initial timestep_end'
[]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[vol_frac]
type = ParsedPostprocessor
expression = 'total_vol / mesh_volume'
pp_names = 'total_vol mesh_volume'
[]
[sensitivity]
type = ElementIntegralVariablePostprocessor
variable = total_sensitivity
[]
[]
[MultiApps]
[sub_app_one]
type = TransientMultiApp
input_files = square_subapp_one.i
[]
[sub_app_two]
type = TransientMultiApp
input_files = square_subapp_two.i
[]
[]
[Transfers]
# First SUB-APP
# To subapp densities
[subapp_one_density]
type = MultiAppCopyTransfer
to_multi_app = sub_app_one
source_variable = mat_den # Here
variable = mat_den
[]
# From subapp sensitivity
[subapp_one_sensitivity]
type = MultiAppCopyTransfer
from_multi_app = sub_app_one
source_variable = Dc # sensitivity_var
variable = sensitivity_one # Here
[]
# Second SUB-APP
# To subapp densities
[subapp_two_density]
type = MultiAppCopyTransfer
to_multi_app = sub_app_two
source_variable = mat_den # Here
variable = mat_den
[]
# From subapp sensitivity
[subapp_two_sensitivity]
type = MultiAppCopyTransfer
from_multi_app = sub_app_two
source_variable = Dc # sensitivity_var
variable = sensitivity_two # Here
[]
[]
(modules/combined/test/tests/optimization/compliance_sensitivity/2d_mmb_2material_cost.i)
vol_frac = 0.5
power = 3
E0 = 1.0e-6
E1 = 0.3
E2 = 1.0
rho0 = 1.0e-6
rho1 = 0.3
rho2 = 1.0
C0 = 1.0e-6
C1 = 0.5
C2 = 1.0
[GlobalParams]
displacements = 'disp_x disp_y'
[]
[Mesh]
[MeshGenerator]
type = GeneratedMeshGenerator
dim = 2
nx = 150
ny = 50
xmin = 0
xmax = 30
ymin = 0
ymax = 10
[]
[node]
type = ExtraNodesetGenerator
input = MeshGenerator
new_boundary = hold
nodes = 0
[]
[push]
type = ExtraNodesetGenerator
input = node
new_boundary = push
coord = '30 10 0'
[]
[]
[Variables]
[disp_x]
[]
[disp_y]
[]
[]
[AuxVariables]
[Dc]
family = MONOMIAL
order = CONSTANT
initial_condition = -1.0
[]
[mat_den]
family = MONOMIAL
order = CONSTANT
initial_condition = ${vol_frac}
[]
[]
[Physics/SolidMechanics/QuasiStatic]
[all]
strain = SMALL
add_variables = true
incremental = false
[]
[]
[BCs]
[no_x]
type = DirichletBC
variable = disp_y
boundary = hold
value = 0.0
[]
[no_y]
type = DirichletBC
variable = disp_x
boundary = right
value = 0.0
[]
[]
[NodalKernels]
[push]
type = NodalGravity
variable = disp_y
boundary = push
gravity_value = -1
mass = 1
[]
[]
[Materials]
[elasticity_tensor]
type = ComputeVariableIsotropicElasticityTensor
youngs_modulus = E_phys
poissons_ratio = poissons_ratio
args = 'mat_den'
[]
[E_phys]
type = DerivativeParsedMaterial
# ordered multimaterial simp
expression = "A1:=(${E0}-${E1})/(${rho0}^${power}-${rho1}^${power}); "
"B1:=${E0}-A1*${rho0}^${power}; E1:=A1*mat_den^${power}+B1; "
"A2:=(${E1}-${E2})/(${rho1}^${power}-${rho2}^${power}); "
"B2:=${E1}-A2*${rho1}^${power}; E2:=A2*mat_den^${power}+B2; "
"if(mat_den<${rho1},E1,E2)"
coupled_variables = 'mat_den'
property_name = E_phys
[]
[Cost]
type = DerivativeParsedMaterial
# ordered multimaterial simp
expression = "A1:=(${C0}-${C1})/(${rho0}^${power}-${rho1}^${power}); "
"B1:=${C0}-A1*${rho0}^${power}; C1:=A1*mat_den^${power}+B1; "
"A2:=(${C1}-${C2})/(${rho1}^${power}-${rho2}^${power}); "
"B2:=${C1}-A2*${rho1}^${power}; C2:=A2*mat_den^${power}+B2; "
"if(mat_den<${rho1},C1,C2)"
coupled_variables = 'mat_den'
property_name = Cost
[]
[poissons_ratio]
type = GenericConstantMaterial
prop_names = poissons_ratio
prop_values = 0.3
[]
[stress]
type = ComputeLinearElasticStress
[]
[dc]
type = ComplianceSensitivity
design_density = mat_den
youngs_modulus = E_phys
[]
[cc]
type = CostSensitivity
design_density = mat_den
cost = Cost
outputs = 'exodus'
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[UserObjects]
[rad_avg]
type = RadialAverage
radius = 1.2
weights = linear
prop_name = sensitivity
execute_on = TIMESTEP_END
force_preaux = true
[]
[update]
type = DensityUpdate
density_sensitivity = Dc
design_density = mat_den
volume_fraction = ${vol_frac}
execute_on = TIMESTEP_BEGIN
[]
[calc_sense]
type = SensitivityFilter
density_sensitivity = Dc
design_density = mat_den
filter_UO = rad_avg
execute_on = TIMESTEP_END
force_postaux = true
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
petsc_options_iname = '-pc_type -pc_factor_mat_solver_package'
petsc_options_value = 'lu superlu_dist'
nl_abs_tol = 1e-8
dt = 1.0
num_steps = 70
[]
[Outputs]
exodus = true
[out]
type = CSV
execute_on = 'TIMESTEP_END'
[]
print_linear_residuals = false
[]
[Postprocessors]
[total_vol]
type = ElementIntegralVariablePostprocessor
variable = mat_den
execute_on = 'INITIAL TIMESTEP_END'
[]
[sensitivity]
type = ElementIntegralMaterialProperty
mat_prop = sensitivity
[]
[]