- vVariable to take the value of.
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Variable to take the value of.
- variableThe name of the variable that this object applies to
C++ Type:AuxVariableName
Unit:(no unit assumed)
Controllable:No
Description:The name of the variable that this object applies to
ProjectionAux
Returns the specified variable as an auxiliary variable with a projection of the source variable. If they are the same type, this amounts to a simple copy.
The ProjectionAux
can be used to make a copy of a field variable, for example to lag them in certain numerical schemes. Many AuxKernels
can use variables as arguments, without modifying them, and store the result in a separate variable. The use of a ProjectionAux
can often be avoided for that reason.
The ProjectionAux
can also be used to project between different finite element variable families and order. The "v" parameter is then the source variable.
The default projection method will attempt to enforce that the two variables have the same value on every quadrature point. If the shape of the source cannot be reproduced by the target variable finite element family and order, the modeler is invited to measure the projection error using a ElementL2Difference postprocessor.
Elemental variables that are discontinuous at nodes are projected to nodal variables by computing nodal values as element-volume-weighted averages of the centroid values of neighbor elements.
The block restriction of the auxkernel, specified using the "block" parameter, is used to select the source variable value as well using the "use_block_restriction_for_source" parameter.
Lower dimensional elements are currently not supported. If you require using lower dimensional elements for projections, please reach out on the MOOSE Discussions forum.
Input Parameters
- _allow_nodal_to_elemental_couplingTrue
Default:True
C++ Type:bool
Controllable:No
- 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
- check_boundary_restrictedTrueWhether to check for multiple element sides on the boundary in the case of a boundary restricted, element aux variable. Setting this to false will allow contribution to a single element's elemental value(s) from multiple boundary sides on the same element (example: when the restricted boundary exists on two or more sides of an element, such as at a corner of a mesh
Default:True
C++ Type:bool
Controllable:No
Description:Whether to check for multiple element sides on the boundary in the case of a boundary restricted, element aux variable. Setting this to false will allow contribution to a single element's elemental value(s) from multiple boundary sides on the same element (example: when the restricted boundary exists on two or more sides of an element, such as at a corner of a mesh
- execute_onLINEAR TIMESTEP_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:LINEAR TIMESTEP_END
C++ Type:ExecFlagEnum
Options:XFEM_MARK, FORWARD, ADJOINT, HOMOGENEOUS_FORWARD, ADJOINT_TIMESTEP_BEGIN, ADJOINT_TIMESTEP_END, NONE, INITIAL, LINEAR, LINEAR_CONVERGENCE, NONLINEAR, NONLINEAR_CONVERGENCE, POSTCHECK, TIMESTEP_END, TIMESTEP_BEGIN, MULTIAPP_FIXED_POINT_END, MULTIAPP_FIXED_POINT_BEGIN, MULTIAPP_FIXED_POINT_CONVERGENCE, FINAL, CUSTOM, PRE_DISPLACE
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.
- ghost_layers1The number of layers of elements to ghost.
Default:1
C++ Type:unsigned short
Controllable:No
Description:The number of layers of elements to ghost.
- use_block_restriction_for_sourceFalseWhether to use the auxkernel block restriction to also restrict the locations selected for source variable values
Default:False
C++ Type:bool
Controllable:No
Description:Whether to use the auxkernel block restriction to also restrict the locations selected for source variable values
Optional Parameters
- control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
- enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:Yes
Description:Set the enabled status of the MooseObject.
- 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
- use_displaced_meshFalseWhether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Default:False
C++ Type:bool
Controllable:No
Description:Whether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Advanced Parameters
- prop_getter_suffixAn optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:An optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
- use_interpolated_stateFalseFor the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Default:False
C++ Type:bool
Controllable:No
Description:For the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Material Property Retrieval Parameters
Input Files
- (test/tests/transfers/coord_transform/both-transformed/copy/sub-app.i)
- (test/tests/auxkernels/nodal_aux_var/multi_update_var_deprecated_test.i)
- (test/tests/problems/verbose_setup/sample.i)
- (test/tests/problems/external_problem/update-ghosted-aux-soln.i)
- (modules/navier_stokes/test/tests/finite_element/ins/cg-dg-hybrid/lid-driven/hybrid-cg-dg.i)
- (test/tests/auxkernels/projection_aux/2d_and_lower_d.i)
- (modules/combined/test/tests/optimization/optimization_density_update/top_opt_2d_pde_filter.i)
- (test/tests/misc/boundary_variable_check/test.i)
- (test/tests/transfers/general_field/nearest_node/duplicated_nearest_node_tests/two_way_many_apps_sub.i)
- (test/tests/auxkernels/projection_aux/2d.i)
- (test/tests/ics/from_exodus_solution/elem_part2.i)
- (test/tests/auxkernels/projection_aux/1d.i)
- (modules/phase_field/test/tests/free_energy_material/MathEBFreeEnergy_split_name.i)
- (test/tests/transfers/coord_transform/both-transformed/interpolation/sub-app.i)
- (test/tests/transfers/coord_transform/both-transformed/mesh-function/sub-app.i)
- (test/tests/outputs/debug/show_execution_auxkernels.i)
- (test/tests/transfers/coord_transform/both-transformed/projection/sub-app.i)
- (test/tests/transfers/coord_transform/both-transformed/nearest-node/sub-app.i)
v
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Variable to take the value of.
block
C++ Type:std::vector<SubdomainName>
Controllable:No
Description:The list of blocks (ids or names) that this object will be applied
use_block_restriction_for_source
Default:False
C++ Type:bool
Controllable:No
Description:Whether to use the auxkernel block restriction to also restrict the locations selected for source variable values
(test/tests/transfers/coord_transform/both-transformed/copy/sub-app.i)
[Mesh]
type = GeneratedMesh
dim = 2
xmin = -1
xmax = 0
ymin = 0
ymax = 1
nx = 10
ny = 10
alpha_rotation = -90
[]
[Variables]
[v][]
[]
[AuxVariables]
[v_elem]
order = CONSTANT
family = MONOMIAL
[]
[w][]
[w_elem]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[v_elem]
type = ProjectionAux
v = v
variable = v_elem
[]
[]
[Kernels]
[diff_v]
type = Diffusion
variable = v
[]
[]
[BCs]
[left_v]
type = DirichletBC
variable = v
boundary = bottom
value = 0
[]
[right_v]
type = DirichletBC
variable = v
boundary = top
value = 1
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_hypre_type'
petsc_options_value = 'hypre boomeramg'
[]
[Outputs]
exodus = true
[]
(test/tests/auxkernels/nodal_aux_var/multi_update_var_deprecated_test.i)
[Mesh]
type = GeneratedMesh
dim = 2
xmin = 0
xmax = 1
ymin = 0
ymax = 1
nx = 2
ny = 2
[]
[Variables]
[u]
order = FIRST
family = LAGRANGE
[]
[]
[AuxVariables]
[tt]
order = FIRST
family = LAGRANGE
initial_condition = 0
[]
[ten]
order = FIRST
family = LAGRANGE
initial_condition = 1
[]
[2k]
order = FIRST
family = LAGRANGE
initial_condition = 2
[]
[]
[Kernels]
[diff]
type = Diffusion
variable = u
[]
[]
[AuxKernels]
[do-no-1]
variable = ten
type = ProjectionAux
v = ten
[]
[do-no-2]
variable = 2k
type = ProjectionAux
v = ten
[]
[all]
variable = tt
type = MultipleUpdateAux
use_deprecated_api = true
u = u
var1 = ten
var2 = 2k
[]
[]
[BCs]
active = 'left right'
[left]
type = DirichletBC
variable = u
boundary = 1
value = 0
[]
[right]
type = DirichletBC
variable = u
boundary = 3
value = 1
[]
[]
[Executioner]
type = Steady
solve_type = 'PJFNK'
[]
[Outputs]
file_base = out_multi_var
exodus = true
[]
(test/tests/problems/verbose_setup/sample.i)
[Mesh]
type = GeneratedMesh
dim = 2
[]
[Variables]
[u]
initial_condition = 3
[]
[]
[Kernels]
[diff]
type = Diffusion
variable = u
[]
[]
[AuxVariables]
[c]
[]
[]
[AuxKernels]
[copy]
type = ProjectionAux
v = u
variable = c
[]
[]
[Materials]
[unused]
type = GenericConstantMaterial
prop_names = 'f1'
prop_values = '2'
[]
[]
[Functions]
[f]
type = ConstantFunction
value = 1
[]
[]
[Problem]
type = FEProblem
solve = false
verbose_setup = true
[]
[Executioner]
type = Steady
[]
(test/tests/problems/external_problem/update-ghosted-aux-soln.i)
[Mesh]
type = GeneratedMesh
dim = 1
nx = 1
[]
[Executioner]
type = Transient
num_steps = 3
[]
[Problem]
type = SyncTestExternalProblem
[]
[AuxVariables]
[copy]
family = MONOMIAL
order = CONSTANT
[]
[]
[AuxKernels]
[copy]
type = ProjectionAux
variable = copy
v = heat_source
[]
[]
[Postprocessors]
[original]
type = PointValue
variable = heat_source
point = '0.0 0.0 0.0'
[]
[copy]
type = PointValue
variable = copy
point = '0.0 0.0 0.0'
[]
[]
[Outputs]
csv = true
[]
(modules/navier_stokes/test/tests/finite_element/ins/cg-dg-hybrid/lid-driven/hybrid-cg-dg.i)
mu = 1
rho = 1
l = 200
U = 1
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = ${l}
ymin = 0
ymax = ${l}
nx = 20
ny = 20
[]
[corner_node]
type = ExtraNodesetGenerator
new_boundary = 'pinned_node'
nodes = '0'
input = gen
[]
[]
[Variables]
[u]
family = MONOMIAL
[]
[v]
family = MONOMIAL
[]
[pressure][]
[]
[Kernels]
[momentum_x_convection]
type = ADConservativeAdvection
variable = u
velocity = 'velocity'
advected_quantity = 'rhou'
[]
[momentum_x_diffusion]
type = MatDiffusion
variable = u
diffusivity = 'mu'
[]
[momentum_x_pressure]
type = PressureGradient
integrate_p_by_parts = false
variable = u
pressure = pressure
component = 0
[]
[momentum_y_convection]
type = ADConservativeAdvection
variable = v
velocity = 'velocity'
advected_quantity = 'rhov'
[]
[momentum_y_diffusion]
type = MatDiffusion
variable = v
diffusivity = 'mu'
[]
[momentum_y_pressure]
type = PressureGradient
integrate_p_by_parts = false
variable = v
pressure = pressure
component = 1
[]
[mass]
type = ADConservativeAdvection
variable = pressure
velocity = velocity
advected_quantity = -1
[]
[]
[DGKernels]
[momentum_x_convection]
type = ADDGAdvection
variable = u
velocity = 'velocity'
advected_quantity = 'rhou'
[]
[momentum_x_diffusion]
type = DGDiffusion
variable = u
sigma = 6
epsilon = -1
diff = 'mu'
[]
[momentum_y_convection]
type = ADDGAdvection
variable = v
velocity = 'velocity'
advected_quantity = 'rhov'
[]
[momentum_y_diffusion]
type = DGDiffusion
variable = v
sigma = 6
epsilon = -1
diff = 'mu'
[]
[]
[BCs]
[u_walls]
type = DGFunctionDiffusionDirichletBC
boundary = 'left bottom right'
variable = u
sigma = 6
epsilon = -1
function = '0'
diff = 'mu'
[]
[v_walls]
type = DGFunctionDiffusionDirichletBC
boundary = 'left bottom right top'
variable = v
sigma = 6
epsilon = -1
function = '0'
diff = 'mu'
[]
[u_top]
type = DGFunctionDiffusionDirichletBC
boundary = 'top'
variable = u
sigma = 6
epsilon = -1
function = '${U}'
diff = 'mu'
[]
[pressure_pin]
type = DirichletBC
variable = pressure
boundary = 'pinned_node'
value = 0
[]
[]
[Materials]
[const]
type = ADGenericConstantMaterial
prop_names = 'rho'
prop_values = '${rho}'
[]
[const_reg]
type = GenericConstantMaterial
prop_names = 'mu'
prop_values = '${mu}'
[]
[vel]
type = ADVectorFromComponentVariablesMaterial
vector_prop_name = 'velocity'
u = u
v = v
[]
[rhou]
type = ADParsedMaterial
property_name = 'rhou'
coupled_variables = 'u'
material_property_names = 'rho'
expression = 'rho*u'
[]
[rhov]
type = ADParsedMaterial
property_name = 'rhov'
coupled_variables = 'v'
material_property_names = 'rho'
expression = 'rho*v'
[]
[]
[AuxVariables]
[vel_x]
family = MONOMIAL
order = CONSTANT
[]
[vel_y]
family = MONOMIAL
order = CONSTANT
[]
[p]
family = MONOMIAL
order = CONSTANT
[]
[]
[AuxKernels]
[vel_x]
type = ProjectionAux
variable = vel_x
v = u
[]
[vel_y]
type = ProjectionAux
variable = vel_y
v = v
[]
[p]
type = ProjectionAux
variable = p
v = pressure
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_factor_shift_type'
petsc_options_value = 'lu NONZERO'
nl_rel_tol = 1e-12
[]
[Outputs]
exodus = true
[]
[Postprocessors]
[Re]
type = ParsedPostprocessor
expression = '${rho} * ${U} * ${l} / ${mu}'
[]
[]
(test/tests/auxkernels/projection_aux/2d_and_lower_d.i)
[Mesh]
[cube]
type = FileMeshGenerator
file = 'mesh/cube.e'
[]
[remove_low_low_D]
type = BlockDeletionGenerator
input = cube
block = 2
[]
second_order = true
allow_renumbering = false
[]
[AuxVariables]
[u_elem]
[InitialCondition]
type = FunctionIC
function = 'x + y * x'
[]
family = MONOMIAL
order = CONSTANT
[]
[u_nodal]
[InitialCondition]
type = FunctionIC
function = 'x + y * x'
[]
[]
[v_high]
order = SECOND
[]
[v_low]
block = 0
[]
[]
[AuxKernels]
[node_to_node_higher_order]
type = ProjectionAux
variable = v_high
v = u_nodal
execute_on = 'INITIAL TIMESTEP_END'
# block restrict the lower D blocks 1 & 2 away
block = 0
[]
[elem_to_nodal]
type = ProjectionAux
variable = v_low
v = u_elem
execute_on = 'INITIAL TIMESTEP_END'
# block restrict the lower D blocks 1 & 2 away
block = 0
[]
[]
[Executioner]
type = Steady
[]
[Problem]
solve = false
[]
[Outputs]
csv = true
exodus = true
show = 'v_high v_low'
[]
(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 = ProjectionAux
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
[]
[]
(test/tests/misc/boundary_variable_check/test.i)
[Problem]
boundary_restricted_elem_integrity_check = true
[]
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 1
nx = 10
xmax = 2
[]
[subdomain1]
input = gen
type = SubdomainBoundingBoxGenerator
bottom_left = '1.0 0 0'
block_id = 1
top_right = '2.0 1.0 0'
[]
[interface]
input = subdomain1
type = SideSetsBetweenSubdomainsGenerator
primary_block = '0'
paired_block = '1'
new_boundary = 'primary0_interface'
[]
[]
[AuxVariables]
[dummy][]
[dummy2]
family = MONOMIAL
order = CONSTANT
block = 1
[]
[dummy3]
family = MONOMIAL
order = CONSTANT
block = 0
[]
[]
[AuxKernels]
active = 'bad'
[bad]
type = ProjectionAux
variable = dummy
v = v
boundary = 'left'
[]
[bad_elemental]
type = ProjectionAux
variable = dummy3
v = dummy2
boundary = 'left'
[]
[]
[Variables]
[u]
block = '0'
[]
[v]
block = '1'
[]
[]
[Kernels]
[diff_u]
type = CoeffParamDiffusion
variable = u
D = 4
block = 0
[]
[diff_v]
type = CoeffParamDiffusion
variable = v
D = 2
block = 1
[]
[]
[InterfaceKernels]
active = 'interface'
[interface]
type = InterfaceDiffusion
variable = u
neighbor_var = v
boundary = primary0_interface
D = 'D'
D_neighbor = 'D'
[]
[penalty_interface]
type = PenaltyInterfaceDiffusion
variable = u
neighbor_var = v
boundary = primary0_interface
penalty = 1e6
[]
[]
[BCs]
active = 'left right middle'
[left]
type = DirichletBC
variable = u
boundary = 'left'
value = 1
[]
[bad]
type = MatchedValueBC
variable = u
boundary = 'left'
v = v
[]
[bad_integrated]
type = CoupledVarNeumannBC
variable = u
boundary = 'left'
v = v
[]
[right]
type = DirichletBC
variable = v
boundary = 'right'
value = 0
[]
[middle]
type = MatchedValueBC
variable = v
boundary = 'primary0_interface'
v = u
[]
[]
[Materials]
[stateful]
type = StatefulMaterial
initial_diffusivity = 1
boundary = primary0_interface
[]
[block0]
type = GenericConstantMaterial
block = '0'
prop_names = 'D'
prop_values = '4'
[]
[block1]
type = GenericConstantMaterial
block = '1'
prop_names = 'D'
prop_values = '2'
[]
[]
[Preconditioning]
[smp]
type = SMP
full = true
[]
[]
[Executioner]
type = Steady
solve_type = NEWTON
[]
[Outputs]
exodus = true
[]
[Postprocessors]
active = ''
[bad]
type = NodalExtremeValue
boundary = 'left'
variable = v
[]
[bad_side]
type = SideDiffusiveFluxIntegral
variable = v
diffusivity = 1
boundary = 'left'
[]
[]
(test/tests/transfers/general_field/nearest_node/duplicated_nearest_node_tests/two_way_many_apps_sub.i)
[Mesh]
type = GeneratedMesh
dim = 2
nx = 10
ny = 10
xmax = 0.2
ymax = 0.2
[]
[Variables]
[u]
[]
[]
[AuxVariables]
[from_parent]
[]
[elemental_from_parent]
order = CONSTANT
family = MONOMIAL
[]
[u_elem]
[]
[]
[Kernels]
[diff]
type = Diffusion
variable = u
[]
[]
[AuxKernels]
# this is done to avoid floating point precision on sending u, with two equidistant points
[copy_over]
type = ProjectionAux
v = u
variable = u_elem
[]
[]
[BCs]
[left]
type = DirichletBC
variable = u
boundary = left
value = 0
[]
[right]
type = DirichletBC
variable = u
boundary = right
value = 1
[]
[]
[Executioner]
type = Transient
num_steps = 1
dt = 1
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type -pc_hypre_type'
petsc_options_value = 'hypre boomeramg'
[]
[Outputs]
exodus = true
[]
(test/tests/auxkernels/projection_aux/2d.i)
[Mesh]
type = GeneratedMesh
dim = 2
xmin = 0
xmax = 1
ymin = 0
ymax = 1
nx = 10
ny = 20
elem_type = QUAD9
[]
[ICs]
active = 'constant_elem constant_nodal'
[constant_elem]
type = ConstantIC
variable = base_elem
value = 4
[]
[constant_nodal]
type = ConstantIC
variable = base_nodal
value = 3.5
[]
[linear_elem]
type = FunctionIC
variable = base_elem
function = 2+2*x-y
[]
[linear_nodal]
type = FunctionIC
variable = base_nodal
function = 3+3*x-y
[]
[cubic_elem]
type = FunctionIC
variable = base_elem
function = 2+2*x*x-3*y*y*y
[]
[cubic_nodal]
type = FunctionIC
variable = base_nodal
function = 3+3*x*x-4*y*y*y
[]
[]
[AuxVariables]
# Families:
# LAGRANGE, MONOMIAL, HERMITE, SCALAR, HIERARCHIC, CLOUGH, XYZ, SZABAB, BERNSTEIN,
# L2_LAGRANGE, L2_HIERARCHIC, RATIONAL_BERNSTEIN, SIDE_HIERARCHIC
# Notes:
# - 'elemental': MONOMIAL, XYZ, L2_LAGRANGE, L2_HIERARCHIC
# - 'nodal': LAGRANGE, HERMITE, HIERARCHIC, CLOUGH, SZABAB, BERNSTEIN, RATIONAL_BERNSTEIN
# - Clough, rational Berstein cannot be created in 2D QUAD9
# - Hermite cannot be created on 2D Tri6
# - Clough, Szabab, Hermite, hierarchic, L2_lagrange, L2_hierarchic, Bernstein cannot be created as constant
[base_elem]
family = MONOMIAL
order = CONSTANT
[]
[base_nodal]
[]
[test_elem_lagrange]
[]
[test_elem_lagrange_high]
order = SECOND
[]
[test_elem_mono]
family = MONOMIAL
order = CONSTANT
[]
[test_elem_mono_high]
family = MONOMIAL
order = SECOND
[]
[test_elem_fv]
type = MooseVariableFVReal
[]
[test_elem_hierarchic]
family = HIERARCHIC
order = FIRST
[]
[test_elem_xyz]
family = XYZ
order = CONSTANT
[]
[test_elem_xyz_high]
family = XYZ
order = SECOND
[]
[test_elem_szabab]
family = SZABAB
order = FIRST
[]
[test_elem_bernstein]
family = BERNSTEIN
order = FIRST
[]
[test_elem_l2_lagrange]
family = L2_LAGRANGE
order = FIRST
[]
[test_elem_l2_lagrange_high]
family = L2_LAGRANGE
order = SECOND
[]
[test_elem_l2_hierarchic]
family = L2_HIERARCHIC
order = FIRST
[]
[test_elem_l2_hierarchic_high]
family = L2_HIERARCHIC
order = SECOND
[]
[test_nodal_lagrange]
[]
[test_nodal_lagrange_high]
order = SECOND
[]
[test_nodal_mono]
family = MONOMIAL
order = CONSTANT
[]
[test_nodal_mono_high]
family = MONOMIAL
order = SECOND
[]
[test_nodal_fv]
type = MooseVariableFVReal
[]
[test_nodal_hierarchic]
family = HIERARCHIC
order = FIRST
[]
[test_nodal_xyz]
family = XYZ
order = CONSTANT
[]
[test_nodal_xyz_high]
family = XYZ
order = SECOND
[]
[test_nodal_szabab]
family = SZABAB
order = FIRST
[]
[test_nodal_bernstein]
family = BERNSTEIN
order = FIRST
[]
[test_nodal_l2_lagrange]
family = L2_LAGRANGE
order = FIRST
[]
[test_nodal_l2_lagrange_high]
family = L2_LAGRANGE
order = SECOND
[]
[test_nodal_l2_hierarchic]
family = L2_HIERARCHIC
order = FIRST
[]
[test_nodal_l2_hierarchic_high]
family = L2_HIERARCHIC
order = SECOND
[]
[]
[AuxKernels]
# Project from constant monomial
[base_elem_proj_lagrange]
type = ProjectionAux
variable = test_elem_lagrange
v = base_elem
[]
[base_elem_proj_lagrange_high]
type = ProjectionAux
variable = test_elem_lagrange_high
v = base_elem
[]
[base_elem_proj_mono]
type = ProjectionAux
variable = test_elem_mono
v = base_elem
[]
[base_elem_proj_mono_high]
type = ProjectionAux
variable = test_elem_mono_high
v = base_elem
[]
[base_elem_proj_fv]
type = ProjectionAux
variable = test_elem_fv
v = base_elem
[]
[base_elem_proj_hierarchic]
type = ProjectionAux
variable = test_elem_hierarchic
v = base_elem
[]
[base_elem_proj_xyz]
type = ProjectionAux
variable = test_elem_xyz
v = base_elem
[]
[base_elem_proj_xyz_high]
type = ProjectionAux
variable = test_elem_xyz_high
v = base_elem
[]
[base_elem_proj_szabab]
type = ProjectionAux
variable = test_elem_szabab
v = base_elem
[]
[base_elem_proj_bernstein]
type = ProjectionAux
variable = test_elem_bernstein
v = base_elem
[]
[base_elem_proj_l2_lagrange]
type = ProjectionAux
variable = test_elem_l2_lagrange
v = base_elem
[]
[base_elem_proj_l2_lagrange_high]
type = ProjectionAux
variable = test_elem_l2_lagrange_high
v = base_elem
[]
[base_elem_proj_l2_hierarchic]
type = ProjectionAux
variable = test_elem_l2_hierarchic
v = base_elem
[]
[base_elem_proj_l2_hierarchic_high]
type = ProjectionAux
variable = test_elem_l2_hierarchic_high
v = base_elem
[]
# Project from constant nodal
[base_nodal_proj_lagrange]
type = ProjectionAux
variable = test_nodal_lagrange
v = base_nodal
[]
[base_nodal_proj_lagrange_high]
type = ProjectionAux
variable = test_nodal_lagrange_high
v = base_nodal
[]
[base_nodal_proj_mono]
type = ProjectionAux
variable = test_nodal_mono
v = base_nodal
[]
[base_nodal_proj_mono_high]
type = ProjectionAux
variable = test_nodal_mono_high
v = base_nodal
[]
[base_nodal_proj_fv]
type = ProjectionAux
variable = test_nodal_fv
v = base_nodal
[]
[base_nodal_proj_hierarchic]
type = ProjectionAux
variable = test_nodal_hierarchic
v = base_nodal
[]
[base_nodal_proj_xyz]
type = ProjectionAux
variable = test_nodal_xyz
v = base_nodal
[]
[base_nodal_proj_xyz_high]
type = ProjectionAux
variable = test_nodal_xyz_high
v = base_nodal
[]
[base_nodal_proj_szabab]
type = ProjectionAux
variable = test_nodal_szabab
v = base_nodal
[]
[base_nodal_proj_bernstein]
type = ProjectionAux
variable = test_nodal_bernstein
v = base_nodal
[]
[base_nodal_proj_l2_lagrange]
type = ProjectionAux
variable = test_nodal_l2_lagrange
v = base_nodal
[]
[base_nodal_proj_l2_lagrange_high]
type = ProjectionAux
variable = test_nodal_l2_lagrange_high
v = base_nodal
[]
[base_nodal_proj_l2_hierarchic]
type = ProjectionAux
variable = test_nodal_l2_hierarchic
v = base_nodal
[]
[base_nodal_proj_l2_hierarchic_high]
type = ProjectionAux
variable = test_nodal_l2_hierarchic_high
v = base_nodal
[]
[]
[Postprocessors]
[base_elem_proj_lagrange]
type = ElementL2Difference
variable = test_elem_lagrange
other_variable = base_elem
[]
[base_elem_proj_lagrange_high]
type = ElementL2Difference
variable = test_elem_lagrange_high
other_variable = base_elem
[]
[base_elem_proj_mono]
type = ElementL2Difference
variable = test_elem_mono
other_variable = base_elem
[]
[base_elem_proj_mono_high]
type = ElementL2Difference
variable = test_elem_mono_high
other_variable = base_elem
[]
[base_elem_proj_fv]
type = ElementL2Difference
variable = test_elem_fv
other_variable = base_elem
[]
[base_elem_proj_hierarchic]
type = ElementL2Difference
variable = test_elem_hierarchic
other_variable = base_elem
[]
[base_elem_proj_xyz]
type = ElementL2Difference
variable = test_elem_xyz
other_variable = base_elem
[]
[base_elem_proj_xyz_high]
type = ElementL2Difference
variable = test_elem_xyz_high
other_variable = base_elem
[]
[base_elem_proj_szabab]
type = ElementL2Difference
variable = test_elem_szabab
other_variable = base_elem
[]
[base_elem_proj_bernstein]
type = ElementL2Difference
variable = test_elem_bernstein
other_variable = base_elem
[]
[base_elem_proj_l2_lagrange]
type = ElementL2Difference
variable = test_elem_l2_lagrange
other_variable = base_elem
[]
[base_elem_proj_l2_lagrange_high]
type = ElementL2Difference
variable = test_elem_l2_lagrange_high
other_variable = base_elem
[]
[base_elem_proj_l2_hierarchic]
type = ElementL2Difference
variable = test_elem_l2_hierarchic
other_variable = base_elem
[]
[base_elem_proj_l2_hierarchic_high]
type = ElementL2Difference
variable = test_elem_l2_hierarchic_high
other_variable = base_elem
[]
[base_nodal_proj_lagrange]
type = ElementL2Difference
variable = test_nodal_lagrange
other_variable = base_nodal
[]
[base_nodal_proj_lagrange_high]
type = ElementL2Difference
variable = test_nodal_lagrange_high
other_variable = base_nodal
[]
[base_nodal_proj_mono]
type = ElementL2Difference
variable = test_nodal_mono
other_variable = base_nodal
[]
[base_nodal_proj_mono_high]
type = ElementL2Difference
variable = test_nodal_mono_high
other_variable = base_nodal
[]
[base_nodal_proj_fv]
type = ElementL2Difference
variable = test_nodal_fv
other_variable = base_nodal
[]
[base_nodal_proj_hierarchic]
type = ElementL2Difference
variable = test_nodal_hierarchic
other_variable = base_nodal
[]
[base_nodal_proj_xyz]
type = ElementL2Difference
variable = test_nodal_xyz
other_variable = base_nodal
[]
[base_nodal_proj_xyz_high]
type = ElementL2Difference
variable = test_nodal_xyz_high
other_variable = base_nodal
[]
[base_nodal_proj_szabab]
type = ElementL2Difference
variable = test_nodal_szabab
other_variable = base_nodal
[]
[base_nodal_proj_bernstein]
type = ElementL2Difference
variable = test_nodal_bernstein
other_variable = base_nodal
[]
[base_nodal_proj_l2_lagrange]
type = ElementL2Difference
variable = test_nodal_l2_lagrange
other_variable = base_nodal
[]
[base_nodal_proj_l2_lagrange_high]
type = ElementL2Difference
variable = test_nodal_l2_lagrange_high
other_variable = base_nodal
[]
[base_nodal_proj_l2_hierarchic]
type = ElementL2Difference
variable = test_nodal_l2_hierarchic
other_variable = base_nodal
[]
[base_nodal_proj_l2_hierarchic_high]
type = ElementL2Difference
variable = test_nodal_l2_hierarchic_high
other_variable = base_nodal
[]
[]
[Outputs]
execute_on = 'timestep_end'
exodus = true
csv = true
[]
[Executioner]
type = Steady
[]
[Problem]
solve = false
[]
(test/tests/ics/from_exodus_solution/elem_part2.i)
# Use the exodus file for restarting the problem:
# - restart elemental aux variable
[Mesh]
[fmg]
type = FileMeshGenerator
file = elem_part1_out.e
use_for_exodus_restart = true
[]
# This problem uses ExodusII_IO::copy_elemental_solution(), which only
# works with ReplicatedMesh
parallel_type = replicated
[]
[Functions]
[exact_fn]
type = ParsedFunction
expression = ((x*x)+(y*y))
[]
[forcing_fn]
type = ParsedFunction
expression = -4
[]
[]
[AuxVariables]
[e]
order = CONSTANT
family = MONOMIAL
initial_from_file_var = e
initial_from_file_timestep = 6
[]
[]
[AuxKernels]
[ak]
type = ProjectionAux
variable = e
v = e
[]
[]
[Variables]
[u]
order = FIRST
family = LAGRANGE
[]
[]
[Kernels]
[diff]
type = Diffusion
variable = u
[]
[ffn]
type = BodyForce
variable = u
function = forcing_fn
[]
[]
[BCs]
[all]
type = FunctionDirichletBC
variable = u
boundary = '0 1 2 3'
function = exact_fn
[]
[]
[Executioner]
type = Steady
solve_type = 'PJFNK'
[]
[Outputs]
exodus = true
[]
(test/tests/auxkernels/projection_aux/1d.i)
[Mesh]
type = GeneratedMesh
dim = 1
xmin = 0
xmax = 1
nx = 10
elem_type = EDGE3
[]
[ICs]
active = 'constant_elem constant_nodal'
[constant_elem]
type = ConstantIC
variable = base_elem
value = 4
[]
[constant_nodal]
type = ConstantIC
variable = base_nodal
value = 3.5
[]
[linear_elem]
type = FunctionIC
variable = base_elem
function = 2+2*x
[]
[linear_nodal]
type = FunctionIC
variable = base_nodal
function = 3+3*x
[]
[quadratic_elem]
type = FunctionIC
variable = base_elem
function = 2+2*x*x
[]
[quadratic_nodal]
type = FunctionIC
variable = base_nodal
function = 3+3*x*x
[]
[]
[AuxVariables]
# Families:
# LAGRANGE, MONOMIAL, HERMITE, SCALAR, HIERARCHIC, CLOUGH, XYZ, SZABAB, BERNSTEIN,
# L2_LAGRANGE, L2_HIERARCHIC, RATIONAL_BERNSTEIN, SIDE_HIERARCHIC
# Notes:
# - 'elemental': MONOMIAL, XYZ, L2_LAGRANGE, L2_HIERARCHIC
# - 'nodal': LAGRANGE, HERMITE, HIERARCHIC, CLOUGH, SZABAB, BERNSTEIN, RATIONAL_BERNSTEIN
# - Hermite, Hierachic, Clough, Bernstein, rational Bernstein cannot be created on 1D EDGE3 elements
# - Hermite, Hierarchic, Clough, Bernstein, rational Bernstein, Szabab, L2_lagrange, L2_hierarchic
# cannot be created with a constant order
[base_elem]
family = MONOMIAL
order = CONSTANT
[]
[base_nodal]
[]
[test_elem_lagrange]
[]
[test_elem_lagrange_high]
order = SECOND
[]
[test_elem_mono]
family = MONOMIAL
order = CONSTANT
[]
[test_elem_mono_high]
family = MONOMIAL
order = SECOND
[]
[test_elem_fv]
type = MooseVariableFVReal
[]
[test_elem_xyz]
family = XYZ
order = CONSTANT
[]
[test_elem_xyz_high]
family = XYZ
order = SECOND
[]
[test_elem_szabab]
family = SZABAB
order = FIRST
[]
[test_elem_l2_lagrange]
family = L2_LAGRANGE
order = FIRST
[]
[test_elem_l2_lagrange_high]
family = L2_LAGRANGE
order = SECOND
[]
[test_elem_l2_hierarchic]
family = L2_HIERARCHIC
order = FIRST
[]
[test_elem_l2_hierarchic_high]
family = L2_HIERARCHIC
order = SECOND
[]
[test_nodal_lagrange]
[]
[test_nodal_lagrange_high]
order = SECOND
[]
[test_nodal_mono]
family = MONOMIAL
order = CONSTANT
[]
[test_nodal_mono_high]
family = MONOMIAL
order = SECOND
[]
[test_nodal_fv]
type = MooseVariableFVReal
[]
[test_nodal_xyz]
family = XYZ
order = CONSTANT
[]
[test_nodal_xyz_high]
family = XYZ
order = SECOND
[]
[test_nodal_szabab]
family = SZABAB
order = FIRST
[]
[test_nodal_l2_lagrange]
family = L2_LAGRANGE
order = FIRST
[]
[test_nodal_l2_lagrange_high]
family = L2_LAGRANGE
order = SECOND
[]
[test_nodal_l2_hierarchic]
family = L2_HIERARCHIC
order = FIRST
[]
[test_nodal_l2_hierarchic_high]
family = L2_HIERARCHIC
order = SECOND
[]
[]
[AuxKernels]
# Project from constant monomial
[base_elem_proj_lagrange]
type = ProjectionAux
variable = test_elem_lagrange
v = base_elem
[]
[base_elem_proj_lagrange_high]
type = ProjectionAux
variable = test_elem_lagrange_high
v = base_elem
[]
[base_elem_proj_mono]
type = ProjectionAux
variable = test_elem_mono
v = base_elem
[]
[base_elem_proj_mono_high]
type = ProjectionAux
variable = test_elem_mono_high
v = base_elem
[]
[base_elem_proj_fv]
type = ProjectionAux
variable = test_elem_fv
v = base_elem
[]
[base_elem_proj_xyz]
type = ProjectionAux
variable = test_elem_xyz
v = base_elem
[]
[base_elem_proj_xyz_high]
type = ProjectionAux
variable = test_elem_xyz_high
v = base_elem
[]
[base_elem_proj_szabab]
type = ProjectionAux
variable = test_elem_szabab
v = base_elem
[]
[base_elem_proj_l2_lagrange]
type = ProjectionAux
variable = test_elem_l2_lagrange
v = base_elem
[]
[base_elem_proj_l2_lagrange_high]
type = ProjectionAux
variable = test_elem_l2_lagrange_high
v = base_elem
[]
[base_elem_proj_l2_hierarchic]
type = ProjectionAux
variable = test_elem_l2_hierarchic
v = base_elem
[]
[base_elem_proj_l2_hierarchic_high]
type = ProjectionAux
variable = test_elem_l2_hierarchic_high
v = base_elem
[]
# Project from constant nodal
[base_nodal_proj_lagrange]
type = ProjectionAux
variable = test_nodal_lagrange
v = base_nodal
[]
[base_nodal_proj_lagrange_high]
type = ProjectionAux
variable = test_nodal_lagrange_high
v = base_nodal
[]
[base_nodal_proj_mono]
type = ProjectionAux
variable = test_nodal_mono
v = base_nodal
[]
[base_nodal_proj_mono_high]
type = ProjectionAux
variable = test_nodal_mono_high
v = base_nodal
[]
[base_nodal_proj_fv]
type = ProjectionAux
variable = test_nodal_fv
v = base_nodal
[]
[base_nodal_proj_xyz]
type = ProjectionAux
variable = test_nodal_xyz
v = base_nodal
[]
[base_nodal_proj_xyz_high]
type = ProjectionAux
variable = test_nodal_xyz_high
v = base_nodal
[]
[base_nodal_proj_szabab]
type = ProjectionAux
variable = test_nodal_szabab
v = base_nodal
[]
[base_nodal_proj_l2_lagrange]
type = ProjectionAux
variable = test_nodal_l2_lagrange
v = base_nodal
[]
[base_nodal_proj_l2_lagrange_high]
type = ProjectionAux
variable = test_nodal_l2_lagrange_high
v = base_nodal
[]
[base_nodal_proj_l2_hierarchic]
type = ProjectionAux
variable = test_nodal_l2_hierarchic
v = base_nodal
[]
[base_nodal_proj_l2_hierarchic_high]
type = ProjectionAux
variable = test_nodal_l2_hierarchic_high
v = base_nodal
[]
[]
[Postprocessors]
[base_elem_proj_lagrange]
type = ElementL2Difference
variable = test_elem_lagrange
other_variable = base_elem
[]
[base_elem_proj_lagrange_high]
type = ElementL2Difference
variable = test_elem_lagrange_high
other_variable = base_elem
[]
[base_elem_proj_mono]
type = ElementL2Difference
variable = test_elem_mono
other_variable = base_elem
[]
[base_elem_proj_mono_high]
type = ElementL2Difference
variable = test_elem_mono_high
other_variable = base_elem
[]
[base_elem_proj_fv]
type = ElementL2Difference
variable = test_elem_fv
other_variable = base_elem
[]
[base_elem_proj_xyz]
type = ElementL2Difference
variable = test_elem_xyz
other_variable = base_elem
[]
[base_elem_proj_xyz_high]
type = ElementL2Difference
variable = test_elem_xyz_high
other_variable = base_elem
[]
[base_elem_proj_szabab]
type = ElementL2Difference
variable = test_elem_szabab
other_variable = base_elem
[]
[base_elem_proj_l2_lagrange]
type = ElementL2Difference
variable = test_elem_l2_lagrange
other_variable = base_elem
[]
[base_elem_proj_l2_lagrange_high]
type = ElementL2Difference
variable = test_elem_l2_lagrange_high
other_variable = base_elem
[]
[base_elem_proj_l2_hierarchic]
type = ElementL2Difference
variable = test_elem_l2_hierarchic
other_variable = base_elem
[]
[base_elem_proj_l2_hierarchic_high]
type = ElementL2Difference
variable = test_elem_l2_hierarchic_high
other_variable = base_elem
[]
# Project from constant nodal
[base_nodal_proj_lagrange]
type = ElementL2Difference
variable = test_nodal_lagrange
other_variable = base_nodal
[]
[base_nodal_proj_lagrange_high]
type = ElementL2Difference
variable = test_nodal_lagrange_high
other_variable = base_nodal
[]
[base_nodal_proj_mono]
type = ElementL2Difference
variable = test_nodal_mono
other_variable = base_nodal
[]
[base_nodal_proj_mono_high]
type = ElementL2Difference
variable = test_nodal_mono_high
other_variable = base_nodal
[]
[base_nodal_proj_fv]
type = ElementL2Difference
variable = test_nodal_fv
other_variable = base_nodal
[]
[base_nodal_proj_xyz]
type = ElementL2Difference
variable = test_nodal_xyz
other_variable = base_nodal
[]
[base_nodal_proj_xyz_high]
type = ElementL2Difference
variable = test_nodal_xyz_high
other_variable = base_nodal
[]
[base_nodal_proj_szabab]
type = ElementL2Difference
variable = test_nodal_szabab
other_variable = base_nodal
[]
[base_nodal_proj_l2_lagrange]
type = ElementL2Difference
variable = test_nodal_l2_lagrange
other_variable = base_nodal
[]
[base_nodal_proj_l2_lagrange_high]
type = ElementL2Difference
variable = test_nodal_l2_lagrange_high
other_variable = base_nodal
[]
[base_nodal_proj_l2_hierarchic]
type = ElementL2Difference
variable = test_nodal_l2_hierarchic
other_variable = base_nodal
[]
[base_nodal_proj_l2_hierarchic_high]
type = ElementL2Difference
variable = test_nodal_l2_hierarchic_high
other_variable = base_nodal
[]
[]
[Outputs]
execute_on = 'timestep_end'
exodus = true
csv = true
[]
[Executioner]
type = Steady
[]
[Problem]
solve = false
[]
(modules/phase_field/test/tests/free_energy_material/MathEBFreeEnergy_split_name.i)
[Mesh]
type = GeneratedMesh
dim = 2
nx = 30
ny = 30
xmin = 0.0
xmax = 30.0
ymin = 0.0
ymax = 30.0
elem_type = QUAD4
[]
[Variables]
[d]
[InitialCondition]
type = CrossIC
x1 = 0.0
x2 = 30.0
y1 = 0.0
y2 = 30.0
[]
[]
[w]
[]
[]
[AuxVariables]
[c]
[]
[]
[AuxKernels]
[c]
type = ProjectionAux
variable = c
v = d
execute_on = 'INITIAL TIMESTEP_END FINAL'
[]
[]
[Preconditioning]
active = 'SMP'
[PBP]
type = PBP
solve_order = 'w d'
preconditioner = 'AMG ASM'
off_diag_row = 'd '
off_diag_column = 'w '
[]
[SMP]
type = SMP
off_diag_row = 'w d'
off_diag_column = 'd w'
[]
[]
[Kernels]
[cres]
type = SplitCHParsed
variable = d
kappa_name = kappa_c
w = w
f_name = F
[]
[wres]
type = SplitCHWRes
variable = w
mob_name = M
[]
[time]
type = CoupledTimeDerivative
variable = w
v = d
[]
[]
[BCs]
[Periodic]
[top_bottom]
primary = 0
secondary = 2
translation = '0 30.0 0'
[]
[left_right]
primary = 1
secondary = 3
translation = '-30.0 0 0'
[]
[]
[]
[Materials]
[constant]
type = GenericConstantMaterial
prop_names = 'M kappa_c'
prop_values = '1.0 2.0'
[]
[free_energy]
type = MathEBFreeEnergy
property_name = F
c = d
[]
[]
[Executioner]
type = Transient
scheme = 'BDF2'
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type'
petsc_options_value = 'lu'
l_max_its = 30
l_tol = 1.0e-3
nl_max_its = 50
nl_rel_tol = 1.0e-10
dt = 10.0
num_steps = 2
[]
[Outputs]
exodus = true
hide = d
[]
(test/tests/transfers/coord_transform/both-transformed/interpolation/sub-app.i)
[Mesh]
type = GeneratedMesh
dim = 2
xmin = -1
xmax = 0
ymin = 0
ymax = 1
nx = 10
ny = 10
alpha_rotation = -90
[]
[Variables]
[v][]
[]
[AuxVariables]
[v_elem]
order = CONSTANT
family = MONOMIAL
[]
[w][]
[w_elem]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[v_elem]
type = ProjectionAux
v = v
variable = v_elem
[]
[]
[Kernels]
[diff_v]
type = Diffusion
variable = v
[]
[]
[BCs]
[left_v]
type = DirichletBC
variable = v
boundary = bottom
value = 0
[]
[right_v]
type = DirichletBC
variable = v
boundary = top
value = 1
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_hypre_type'
petsc_options_value = 'hypre boomeramg'
[]
[Outputs]
exodus = true
[]
(test/tests/transfers/coord_transform/both-transformed/mesh-function/sub-app.i)
[Mesh]
type = GeneratedMesh
dim = 2
xmin = -1
xmax = 0
ymin = 0
ymax = 1
nx = 10
ny = 10
alpha_rotation = -90
[]
[Variables]
[v][]
[]
[AuxVariables]
[v_elem]
order = CONSTANT
family = MONOMIAL
[]
[w][]
[w_elem]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[v_elem]
type = ProjectionAux
v = v
variable = v_elem
[]
[]
[Kernels]
[diff_v]
type = Diffusion
variable = v
[]
[]
[BCs]
[left_v]
type = DirichletBC
variable = v
boundary = bottom
value = 0
[]
[right_v]
type = DirichletBC
variable = v
boundary = top
value = 1
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_hypre_type'
petsc_options_value = 'hypre boomeramg'
[]
[Outputs]
exodus = true
[]
(test/tests/outputs/debug/show_execution_auxkernels.i)
[Mesh]
type = GeneratedMesh
dim = 2
xmin = 0
xmax = 1
ymin = 0
ymax = 1
nx = 10
ny = 10
elem_type = QUAD4
[]
[Debug]
show_execution_order = 'ALWAYS'
[]
[AuxVariables]
[a]
initial_condition = 1
[]
[b]
initial_condition = 2
[]
[c]
initial_condition = 3
[]
[a_elem]
order = CONSTANT
family = MONOMIAL
initial_condition = 1
[]
[b_elem]
order = CONSTANT
family = MONOMIAL
initial_condition = 2
[]
[c_elem]
order = CONSTANT
family = MONOMIAL
initial_condition = 3
[]
[d_elem]
order = CONSTANT
family = MONOMIAL
initial_condition = 3
[]
[]
[Variables]
[u]
order = FIRST
family = LAGRANGE
[]
[]
[Functions]
[exact_fn]
type = ParsedFunction
expression = t
[]
[a_fn]
type = ParsedFunction
expression = t
[]
[b_fn]
type = ParsedFunction
expression = (4-t)/2
[]
[]
[AuxKernels]
# Nodal
# this one needs a and b set, should run last
[c_saux]
type = QuotientAux
variable = c
numerator = a
denominator = b
execute_on = 'initial timestep_end'
[]
# setting b requires a
[b_saux]
type = ProjectionAux
variable = b
v = a
execute_on = 'linear timestep_end'
[]
# Elements
# this one needs a and b set, should run last
[c_saux_elem]
type = QuotientAux
variable = c_elem
numerator = a_elem
denominator = b_elem
execute_on = 'initial timestep_end'
[]
# setting b requires a
[b_saux_elem]
type = ProjectionAux
variable = b_elem
v = a_elem
execute_on = 'linear timestep_end'
[]
# boundary auxkernel
[real_property]
type = MaterialRealAux
variable = d_elem
property = 3
boundary = 'top bottom'
[]
[]
[Kernels]
[ie]
type = TimeDerivative
variable = u
[]
[diff]
type = Diffusion
variable = u
[]
[]
[BCs]
[all]
type = FunctionDirichletBC
variable = u
boundary = '0 1 2 3'
function = exact_fn
[]
[]
[Executioner]
type = Transient
scheme = 'implicit-euler'
solve_type = 'PJFNK'
start_time = 0.0
num_steps = 1
dt = 1
[]
(test/tests/transfers/coord_transform/both-transformed/projection/sub-app.i)
[Mesh]
type = GeneratedMesh
dim = 2
xmin = -1
xmax = 0
ymin = 0
ymax = 1
nx = 10
ny = 10
alpha_rotation = -90
[]
[Variables]
[v][]
[]
[AuxVariables]
[v_elem]
order = CONSTANT
family = MONOMIAL
[]
[w][]
[w_elem]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[v_elem]
type = ProjectionAux
v = v
variable = v_elem
[]
[]
[Kernels]
[diff_v]
type = Diffusion
variable = v
[]
[]
[BCs]
[left_v]
type = DirichletBC
variable = v
boundary = bottom
value = 0
[]
[right_v]
type = DirichletBC
variable = v
boundary = top
value = 1
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_hypre_type'
petsc_options_value = 'hypre boomeramg'
[]
[Outputs]
exodus = true
[]
(test/tests/transfers/coord_transform/both-transformed/nearest-node/sub-app.i)
[Mesh]
type = GeneratedMesh
dim = 2
xmin = -1
xmax = 0
ymin = 0
ymax = 1
nx = 10
ny = 10
alpha_rotation = -90
[]
[Variables]
[v][]
[]
[AuxVariables]
[v_elem]
order = CONSTANT
family = MONOMIAL
[]
[w][]
[w_elem]
order = CONSTANT
family = MONOMIAL
[]
[]
[AuxKernels]
[v_elem]
type = ProjectionAux
v = v
variable = v_elem
[]
[]
[Kernels]
[diff_v]
type = Diffusion
variable = v
[]
[]
[BCs]
[left_v]
type = DirichletBC
variable = v
boundary = bottom
value = 0
[]
[right_v]
type = DirichletBC
variable = v
boundary = top
value = 1
[]
[]
[Executioner]
type = Steady
solve_type = 'NEWTON'
petsc_options_iname = '-pc_type -pc_hypre_type'
petsc_options_value = 'hypre boomeramg'
[]
[Outputs]
exodus = true
[]