- default_value-1Value to return if target value is not found on line and error_if_not_found is false.
Default:-1
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Value to return if target value is not found on line and error_if_not_found is false.
- depth36Maximum number of bisections to perform.
Default:36
C++ Type:unsigned int
Controllable:No
Description:Maximum number of bisections to perform.
- end_pointEnd point of the sampling line.
C++ Type:libMesh::Point
Controllable:No
Description:End point of the sampling line.
- error_if_not_foundTrueIf true, stop with error if target value is not found on the line. If false, return default_value.
Default:True
C++ Type:bool
Controllable:No
Description:If true, stop with error if target value is not found on the line. If false, return default_value.
- start_pointStart point of the sampling line.
C++ Type:libMesh::Point
Controllable:No
Description:Start point of the sampling line.
- targetTarget value to locate.
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Target value to locate.
- tol1e-10Stop search if a value is found that is equal to the target with this tolerance applied.
Default:1e-10
C++ Type:double
Unit:(no unit assumed)
Controllable:No
Description:Stop search if a value is found that is equal to the target with this tolerance applied.
- vVariable to inspect
C++ Type:std::vector<VariableName>
Unit:(no unit assumed)
Controllable:No
Description:Variable to inspect
FindValueOnLine
The FindValueOnLine Postprocessor is for sampling a line through a mesh or on a boundary looking for the point at which a specific value occurs. The values of the line must behave monotonically or an indeterminate solution may be found.
Find a specific target value along a sampling line. The variable values along the line should change monotonically. The target value is searched using a bisection algorithm.
Input Parameters
- allow_duplicate_execution_on_initialFalseIn the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).
Default:False
C++ Type:bool
Controllable:No
Description:In the case where this UserObject is depended upon by an initial condition, allow it to be executed twice during the initial setup (once before the IC and again after mesh adaptivity (if applicable).
- execute_onTIMESTEP_ENDThe list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html.
Default:TIMESTEP_END
C++ Type:ExecFlagEnum
Options:XFEM_MARK, FORWARD, ADJOINT, HOMOGENEOUS_FORWARD, ADJOINT_TIMESTEP_BEGIN, ADJOINT_TIMESTEP_END, NONE, INITIAL, LINEAR, LINEAR_CONVERGENCE, NONLINEAR, NONLINEAR_CONVERGENCE, POSTCHECK, TIMESTEP_END, TIMESTEP_BEGIN, MULTIAPP_FIXED_POINT_END, MULTIAPP_FIXED_POINT_BEGIN, MULTIAPP_FIXED_POINT_CONVERGENCE, FINAL, CUSTOM, TRANSFER
Controllable:No
Description:The list of flag(s) indicating when this object should be executed. For a description of each flag, see https://mooseframework.inl.gov/source/interfaces/SetupInterface.html.
- execution_order_group0Execution order groups are executed in increasing order (e.g., the lowest number is executed first). Note that negative group numbers may be used to execute groups before the default (0) group. Please refer to the user object documentation for ordering of user object execution within a group.
Default:0
C++ Type:int
Controllable:No
Description:Execution order groups are executed in increasing order (e.g., the lowest number is executed first). Note that negative group numbers may be used to execute groups before the default (0) group. Please refer to the user object documentation for ordering of user object execution within a group.
- force_postauxFalseForces the UserObject to be executed in POSTAUX
Default:False
C++ Type:bool
Controllable:No
Description:Forces the UserObject to be executed in POSTAUX
- force_preauxFalseForces the UserObject to be executed in PREAUX
Default:False
C++ Type:bool
Controllable:No
Description:Forces the UserObject to be executed in PREAUX
- force_preicFalseForces the UserObject to be executed in PREIC during initial setup
Default:False
C++ Type:bool
Controllable:No
Description:Forces the UserObject to be executed in PREIC during initial setup
Execution Scheduling Parameters
- control_tagsAdds user-defined labels for accessing object parameters via control logic.
C++ Type:std::vector<std::string>
Controllable:No
Description:Adds user-defined labels for accessing object parameters via control logic.
- enableTrueSet the enabled status of the MooseObject.
Default:True
C++ Type:bool
Controllable:Yes
Description:Set the enabled status of the MooseObject.
- outputsVector of output names where you would like to restrict the output of variables(s) associated with this object
C++ Type:std::vector<OutputName>
Controllable:No
Description:Vector of output names where you would like to restrict the output of variables(s) associated with this object
- use_displaced_meshFalseWhether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Default:False
C++ Type:bool
Controllable:No
Description:Whether or not this object should use the displaced mesh for computation. Note that in the case this is true but no displacements are provided in the Mesh block the undisplaced mesh will still be used.
Advanced Parameters
- prop_getter_suffixAn optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
C++ Type:MaterialPropertyName
Unit:(no unit assumed)
Controllable:No
Description:An optional suffix parameter that can be appended to any attempt to retrieve/get material properties. The suffix will be prepended with a '_' character.
- use_interpolated_stateFalseFor the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Default:False
C++ Type:bool
Controllable:No
Description:For the old and older state use projected material properties interpolated at the quadrature points. To set up projection use the ProjectedStatefulMaterialStorageAction.
Material Property Retrieval Parameters
Input Files
- (modules/phase_field/examples/anisotropic_interfaces/echebarria_iso.i)
- (modules/phase_field/test/tests/phase_field_contact_angle/contact_angle_verification.i)
- (modules/combined/examples/phase_field-mechanics/kks_mechanics_KHS.i)
- (modules/combined/examples/publications/rapid_dev/fig7b.i)
- (modules/combined/examples/publications/rapid_dev/fig7a.i)
- (modules/combined/test/tests/phase_field_contact_angle/contact_angle.i)
- (test/tests/postprocessors/find_value_on_line/findvalueonline.i)
References
No citations exist within this document.(modules/phase_field/examples/anisotropic_interfaces/echebarria_iso.i)
#This example implements an isotropic, isothermal version of the
#binary alloy solidification model of Echebarria et al.,
#Physical Review E, 70, 061604 (2004). The governing equations are (132)-(133)
#Temperature gradient, pulling velocity, and interfacial energy anisotropy
#are not included.
#The sinusoidal perturbation at the interface decays appproximately
#exponentially with approximate decay constant 1.55e-4, in agreement
#with linear stability analysis
[Mesh]
type = GeneratedMesh
dim = 2
nx = 120
ny = 300
xmin = 0
xmax = 96
ymin = 0
ymax = 240
[]
[GlobalParams]
enable_jit = true
derivative_order = 2
[]
[Variables]
[phi]
[]
[U]
[]
[]
[AuxVariables]
[c]
[]
[]
[ICs]
[phi_IC]
type = FunctionIC
variable = phi
function = ic_func_phi
[]
[U_IC]
type = FunctionIC
variable = U
function = ic_func_U
[]
[]
[Functions]
[ic_func_phi]
type = ParsedFunction
symbol_names = 'midpoint lambda A'
symbol_values = '40 96 8'
expression = 'tanh(-(y - (midpoint + A * cos(2 * pi * x / lambda))) / sqrt(2))'
[]
[ic_func_U]
type = ParsedFunction
expression = '0'
[]
[]
[Kernels]
# Order parameter phi
[AC_dphi_dt]
type = SusceptibilityTimeDerivative
variable = phi
f_name = dphi_dt_pre
[]
[AC_grad]
type = MatDiffusion
variable = phi
diffusivity = as_sq
[]
[AC_floc]
type = AllenCahn
variable = phi
f_name = f_loc
mob_name = L
coupled_variables = 'U'
[]
#dimensionless supersaturation U
[diff_dU_dt]
type = SusceptibilityTimeDerivative
variable = U
f_name = dU_dt_pre
coupled_variables = 'phi'
[]
[diff_grad]
type = MatDiffusion
variable = U
diffusivity = D_interp
args = 'phi'
[]
[diff_antitrapping]
type = AntitrappingCurrent
variable = U
v = phi
f_name = antitrap_pre
[]
[diff_dphidt]
type = CoupledSusceptibilityTimeDerivative
variable = U
v = phi
f_name = coupled_pre
[]
[]
[AuxKernels]
[c_aux]
type = ParsedAux
variable = c
constant_names = 'c_l k'
constant_expressions = '0.33 0.1712'
coupled_variables = 'phi U'
expression = '(1 + (1-k) * U) / 2 * c_l * (1+k - (1-k)*phi)'
[]
[]
[Materials]
[dphi_dt_pre_material]
type = DerivativeParsedMaterial
property_name = dphi_dt_pre
material_property_names = 'as_sq(phi) k'
expression = '(1-(1-k)*0) * as_sq'
[]
[as_sq_material]
type = DerivativeParsedMaterial
property_name = as_sq
expression = '1'
[]
[f_loc_material]
type = DerivativeParsedMaterial
property_name = f_loc
coupled_variables = 'phi U'
constant_names = 'a1'
constant_expressions = '5*sqrt(2)/8'
material_property_names = 'epsilon'
expression = '-phi^2/2 + phi^4/4 + a1 * epsilon * (phi - 2*phi^3/3 + phi^5/5) * U'
[]
[dU_dt_pre_material]
type = DerivativeParsedMaterial
property_name = dU_dt_pre
coupled_variables = 'phi'
material_property_names = 'k'
expression = '(1+k)/2 - (1-k)/2 * phi'
[]
[D_interp_material]
type = DerivativeParsedMaterial
property_name = D_interp
coupled_variables = 'phi'
material_property_names = 'epsilon'
constant_names = ' a1 a2'
constant_expressions = '5*sqrt(2)/8 0.6267'
expression = 'a1 * a2 * epsilon * (1-phi)/2'
# expression = 'a1 * a2 * epsilon'
output_properties = 'D_interp'
outputs = 'exodus'
[]
[antitrap_pre_material]
type = DerivativeParsedMaterial
property_name = antitrap_pre
coupled_variables = 'U'
material_property_names = 'k'
expression = '1/(2*sqrt(2)) * (1 + (1-k) * U)'
[]
[coupled_pre_material]
type = DerivativeParsedMaterial
property_name = coupled_pre
coupled_variables = 'U'
material_property_names = 'k'
expression = '- (1 + (1-k) * U) / 2'
[]
[const]
type = GenericConstantMaterial
prop_names = ' L k epsilon'
prop_values = '1.0 0.1712 30 '
[]
[]
[Postprocessors]
[int_position]
type = FindValueOnLine
start_point = '0 0 0'
end_point = '0 100 0'
v = phi
target = 0
tol = 1e-8
[]
[]
[Preconditioning]
[SMP]
type = SMP
full = true
[]
[]
[Executioner]
type = Transient
scheme = bdf2
petsc_options_iname = '-pc_type -ksp_gmres_restart -sub_pc_type -pc_asm_overlap'
petsc_options_value = 'asm 31 lu 1'
l_tol = 1.0e-3
nl_rel_tol = 1.0e-8
nl_abs_tol = 1e-9
end_time = 1e9
[TimeStepper]
type = IterationAdaptiveDT
dt = 1
optimal_iterations = 8
iteration_window = 2
[]
[]
[Outputs]
exodus = true
csv = true
[]
(modules/phase_field/test/tests/phase_field_contact_angle/contact_angle_verification.i)
sigma = 25e-3 #10e-3 #25e-3 #surface tension coefficient
epsilon = 1e-6 #width parameter
nu = 1e-4#mobility parameter
contactangle = 2.61799#0.523599#1.0472
lambda = ${fparse 3*sigma*epsilon/(2*sqrt(2))}
prefactor_phi = ${fparse nu*lambda/(epsilon*epsilon)}
prefactor_psi = ${fparse -epsilon*epsilon}
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = 0.2e-3
ymin = 0
ymax = 0.2e-3
nx = 20
ny = 20
elem_type = QUAD9
[]
[]
[ICs]
[pf_ic]
type = BoundingBoxIC
variable = pf
x1 = 0.1e-3
y1 = -0.1e-3
x2 = 0.3e-03
y2 = 0.3e-3
inside = 1
outside = -1
int_width = ${fparse 2*sqrt(2)*epsilon}
[]
[velocity]
type = VectorConstantIC
x_value = 0.0
y_value = 0.0
variable = velocity
[]
[]
[Variables]
[pf]
family = LAGRANGE
order = second
[]
[auxpf]
family = LAGRANGE
order = second
[]
[velocity]
family = LAGRANGE_VEC
[]
[]
[Kernels]
[velocity_timederivative]
type = ADVectorTimeDerivative
variable = velocity
[]
[phasefield_timederivative]
type = ADTimeDerivative
variable = pf
[]
[phasefield_supg]
type = ADPhaseFieldTimeDerivativeSUPG
velocity = velocity
variable = pf
[]
[phasefield_laplacian]
type=ADPrefactorLaplacianSplit
variable = pf
c = auxpf
prefactor = ${prefactor_phi}
[]
[Auxphasefield_firstorder]
type=ADReaction
variable = auxpf
rate = 1.0
[]
[Auxphasefield_laplacian]
type=ADPrefactorLaplacianSplit
variable = auxpf
c = pf
prefactor=${prefactor_psi}
[]
[Auxphasefield_doublewell]
type=ADPhaseFieldCoupledDoubleWellPotential
variable = auxpf
c = pf
prefactor=-1.0
[]
[]
[BCs]
[velocity]
type = VectorFunctionDirichletBC
variable = velocity
boundary = 'left right top bottom'
function_x = 0.0
[]
[ContactangleBC]
type=ADPhaseFieldContactAngleBC
variable = auxpf
pf = pf
epsilon = ${epsilon}
lambda=${lambda}
sigma=${sigma}
contactangle=${contactangle}
boundary = 'top bottom'
[]
[]
[Materials]
[rho]
type = ADPhaseFieldTwoPhaseMaterial
prop_name = rho
prop_value_1 = 1000
prop_value_2 = 840
pf = pf
# outputs = exodus
[]
[mu]
type = ADPhaseFieldTwoPhaseMaterial
prop_name = mu
prop_value_1 = 1e-3
prop_value_2 = 7.6e-3
pf = pf
# outputs = exodus
[]
[]
[Postprocessors]
[contact_angle_top]
type = ObtainAvgContactAngle
boundary = top
pf=pf
execute_on = 'timestep_end'
[]
[x_position]
type = FindValueOnLine
start_point = '0 0.0001 0'
end_point ='0.0002 0.0001 0'
v = pf
target = 0.0
tol = 1e-6
[]
[]
[Preconditioning]
[SMP]
type = SMP
full = true
[]
[]
[Adaptivity]
initial_steps = 2
initial_marker = phase_marker
marker = phase_marker
max_h_level = 4
[Markers]
[phase_marker]
type = ValueRangeMarker
lower_bound = -0.99
upper_bound = 0.99
variable = pf
[]
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
start_time = 0
num_steps = 5
dtmax = 0.25
[TimeStepper]
type = IterationAdaptiveDT
dt = 1e-10
iteration_window = 2
optimal_iterations = 10
growth_factor = 2
cutback_factor = 0.5
[]
# petsc_options_iname = '-pc_type -ksp_gmres_restart -pc_factor_mat_solver_type -pc_factor_shift_type -pc_factor_shift_amount'
# petsc_options_value = 'lu 50 superlu_dist NONZERO 1e-15'
petsc_options_iname = '-pc_type -pc_factor_shift_type -pc_factor_shift_amount'
petsc_options_value = 'lu NONZERO 1e-15 '
line_search = 'none'
nl_rel_tol = 1e-6
nl_abs_tol = 1e-8
nl_max_its = 20
nl_forced_its = 3
l_tol = 1e-6
l_max_its = 20
[]
[Outputs]
[csv]
type = CSV
time_step_interval = 1
[]
[]
(modules/combined/examples/phase_field-mechanics/kks_mechanics_KHS.i)
# KKS phase-field model coupled with elasticity using Khachaturyan's scheme as
# described in L.K. Aagesen et al., Computational Materials Science, 140, 10-21 (2017)
# Original run #170403a
[Mesh]
type = GeneratedMesh
dim = 3
nx = 640
ny = 1
nz = 1
xmin = -10
xmax = 10
ymin = 0
ymax = 0.03125
zmin = 0
zmax = 0.03125
elem_type = HEX8
[]
[Variables]
# order parameter
[./eta]
order = FIRST
family = LAGRANGE
[../]
# solute concentration
[./c]
order = FIRST
family = LAGRANGE
[../]
# chemical potential
[./w]
order = FIRST
family = LAGRANGE
[../]
# solute phase concentration (matrix)
[./cm]
order = FIRST
family = LAGRANGE
[../]
# solute phase concentration (precipitate)
[./cp]
order = FIRST
family = LAGRANGE
[../]
[./disp_x]
order = FIRST
family = LAGRANGE
[../]
[./disp_y]
order = FIRST
family = LAGRANGE
[../]
[./disp_z]
order = FIRST
family = LAGRANGE
[../]
[]
[ICs]
[./eta_ic]
variable = eta
type = FunctionIC
function = ic_func_eta
block = 0
[../]
[./c_ic]
variable = c
type = FunctionIC
function = ic_func_c
block = 0
[../]
[./w_ic]
variable = w
type = ConstantIC
value = 0.00991
block = 0
[../]
[./cm_ic]
variable = cm
type = ConstantIC
value = 0.131
block = 0
[../]
[./cp_ic]
variable = cp
type = ConstantIC
value = 0.236
block = 0
[../]
[]
[Functions]
[./ic_func_eta]
type = ParsedFunction
expression = '0.5*(1.0+tanh((x)/delta_eta/sqrt(2.0)))'
symbol_names = 'delta_eta'
symbol_values = '0.8034'
[../]
[./ic_func_c]
type = ParsedFunction
expression = '0.2389*(0.5*(1.0+tanh(x/delta/sqrt(2.0))))^3*(6*(0.5*(1.0+tanh(x/delta/sqrt(2.0))))^2-15*(0.5*(1.0+tanh(x/delta/sqrt(2.0))))+10)+0.1339*(1-(0.5*(1.0+tanh(x/delta/sqrt(2.0))))^3*(6*(0.5*(1.0+tanh(x/delta/sqrt(2.0))))^2-15*(0.5*(1.0+tanh(x/delta/sqrt(2.0))))+10))'
symbol_names = 'delta'
symbol_values = '0.8034'
[../]
[./psi_eq_int]
type = ParsedFunction
expression = 'volume*psi_alpha'
symbol_names = 'volume psi_alpha'
symbol_values = 'volume psi_alpha'
[../]
[./gamma]
type = ParsedFunction
expression = '(psi_int - psi_eq_int) / dy / dz'
symbol_names = 'psi_int psi_eq_int dy dz'
symbol_values = 'psi_int psi_eq_int 0.03125 0.03125'
[../]
[]
[AuxVariables]
[./sigma11]
order = CONSTANT
family = MONOMIAL
[../]
[./sigma22]
order = CONSTANT
family = MONOMIAL
[../]
[./sigma33]
order = CONSTANT
family = MONOMIAL
[../]
[./e11]
order = CONSTANT
family = MONOMIAL
[../]
[./e12]
order = CONSTANT
family = MONOMIAL
[../]
[./e22]
order = CONSTANT
family = MONOMIAL
[../]
[./e33]
order = CONSTANT
family = MONOMIAL
[../]
[./e_el11]
order = CONSTANT
family = MONOMIAL
[../]
[./e_el12]
order = CONSTANT
family = MONOMIAL
[../]
[./e_el22]
order = CONSTANT
family = MONOMIAL
[../]
[./f_el]
order = CONSTANT
family = MONOMIAL
[../]
[./eigen_strain00]
order = CONSTANT
family = MONOMIAL
[../]
[./Fglobal]
order = CONSTANT
family = MONOMIAL
[../]
[./psi]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./matl_sigma11]
type = RankTwoAux
rank_two_tensor = stress
index_i = 0
index_j = 0
variable = sigma11
[../]
[./matl_sigma22]
type = RankTwoAux
rank_two_tensor = stress
index_i = 1
index_j = 1
variable = sigma22
[../]
[./matl_sigma33]
type = RankTwoAux
rank_two_tensor = stress
index_i = 2
index_j = 2
variable = sigma33
[../]
[./matl_e11]
type = RankTwoAux
rank_two_tensor = total_strain
index_i = 0
index_j = 0
variable = e11
[../]
[./f_el]
type = MaterialRealAux
variable = f_el
property = f_el_mat
execute_on = timestep_end
[../]
[./GlobalFreeEnergy]
variable = Fglobal
type = KKSGlobalFreeEnergy
fa_name = fm
fb_name = fp
w = 0.0264
kappa_names = kappa
interfacial_vars = eta
[../]
[./psi_potential]
variable = psi
type = ParsedAux
coupled_variables = 'Fglobal w c f_el sigma11 e11'
expression = 'Fglobal - w*c + f_el - sigma11*e11'
[../]
[]
[BCs]
[./left_x]
type = DirichletBC
variable = disp_x
boundary = left
value = 0
[../]
[./right_x]
type = DirichletBC
variable = disp_x
boundary = right
value = 0
[../]
[./front_y]
type = DirichletBC
variable = disp_y
boundary = front
value = 0
[../]
[./back_y]
type = DirichletBC
variable = disp_y
boundary = back
value = 0
[../]
[./top_z]
type = DirichletBC
variable = disp_z
boundary = top
value = 0
[../]
[./bottom_z]
type = DirichletBC
variable = disp_z
boundary = bottom
value = 0
[../]
[]
[Materials]
# Chemical free energy of the matrix
[./fm]
type = DerivativeParsedMaterial
property_name = fm
coupled_variables = 'cm'
expression = '6.55*(cm-0.13)^2'
[../]
# Chemical Free energy of the precipitate phase
[./fp]
type = DerivativeParsedMaterial
property_name = fp
coupled_variables = 'cp'
expression = '6.55*(cp-0.235)^2'
[../]
# Elastic energy of the precipitate
[./elastic_free_energy_p]
type = ElasticEnergyMaterial
f_name = f_el_mat
args = 'eta'
outputs = exodus
[../]
# h(eta)
[./h_eta]
type = SwitchingFunctionMaterial
h_order = HIGH
eta = eta
[../]
# 1- h(eta), putting in function explicitly
[./one_minus_h_eta_explicit]
type = DerivativeParsedMaterial
property_name = one_minus_h_explicit
coupled_variables = eta
expression = 1-eta^3*(6*eta^2-15*eta+10)
outputs = exodus
[../]
# g(eta)
[./g_eta]
type = BarrierFunctionMaterial
g_order = SIMPLE
eta = eta
[../]
# constant properties
[./constants]
type = GenericConstantMaterial
prop_names = 'M L kappa misfit'
prop_values = '0.7 0.7 0.01704 0.00377'
[../]
#Mechanical properties
[./Stiffness_matrix]
type = ComputeElasticityTensor
base_name = C_matrix
C_ijkl = '103.3 74.25 74.25 103.3 74.25 103.3 46.75 46.75 46.75'
fill_method = symmetric9
[../]
[./Stiffness_ppt]
type = ComputeElasticityTensor
C_ijkl = '100.7 71.45 71.45 100.7 71.45 100.7 50.10 50.10 50.10'
base_name = C_ppt
fill_method = symmetric9
[../]
[./C]
type = CompositeElasticityTensor
args = eta
tensors = 'C_matrix C_ppt'
weights = 'one_minus_h_explicit h'
[../]
[./stress]
type = ComputeLinearElasticStress
[../]
[./strain]
type = ComputeSmallStrain
displacements = 'disp_x disp_y disp_z'
eigenstrain_names = 'eigenstrain_ppt'
[../]
[./eigen_strain]
type = ComputeVariableEigenstrain
eigen_base = '0.00377 0.00377 0.00377 0 0 0'
prefactor = h
args = eta
eigenstrain_name = 'eigenstrain_ppt'
[../]
[]
[Kernels]
[./TensorMechanics]
displacements = 'disp_x disp_y disp_z'
[../]
# enforce c = (1-h(eta))*cm + h(eta)*cp
[./PhaseConc]
type = KKSPhaseConcentration
ca = cm
variable = cp
c = c
eta = eta
[../]
# enforce pointwise equality of chemical potentials
[./ChemPotVacancies]
type = KKSPhaseChemicalPotential
variable = cm
cb = cp
fa_name = fm
fb_name = fp
[../]
#
# Cahn-Hilliard Equation
#
[./CHBulk]
type = KKSSplitCHCRes
variable = c
ca = cm
fa_name = fm
w = w
[../]
[./dcdt]
type = CoupledTimeDerivative
variable = w
v = c
[../]
[./ckernel]
type = SplitCHWRes
mob_name = M
variable = w
[../]
#
# Allen-Cahn Equation
#
[./ACBulkF]
type = KKSACBulkF
variable = eta
fa_name = fm
fb_name = fp
w = 0.0264
args = 'cp cm'
[../]
[./ACBulkC]
type = KKSACBulkC
variable = eta
ca = cm
cb = cp
fa_name = fm
[../]
[./ACBulk_el] #This adds df_el/deta for strain interpolation
type = AllenCahn
variable = eta
f_name = f_el_mat
[../]
[./ACInterface]
type = ACInterface
variable = eta
kappa_name = kappa
[../]
[./detadt]
type = TimeDerivative
variable = eta
[../]
[]
[Executioner]
type = Transient
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type -sub_pc_type -sub_pc_factor_shift_type'
petsc_options_value = 'asm ilu nonzero'
l_max_its = 30
nl_max_its = 10
l_tol = 1.0e-4
nl_rel_tol = 1.0e-8
nl_abs_tol = 1.0e-11
num_steps = 200
[./TimeStepper]
type = SolutionTimeAdaptiveDT
dt = 0.5
[../]
[]
[Postprocessors]
[./f_el_int]
type = ElementIntegralMaterialProperty
mat_prop = f_el_mat
[../]
[./c_alpha]
type = SideAverageValue
boundary = left
variable = c
[../]
[./c_beta]
type = SideAverageValue
boundary = right
variable = c
[../]
[./e11_alpha]
type = SideAverageValue
boundary = left
variable = e11
[../]
[./e11_beta]
type = SideAverageValue
boundary = right
variable = e11
[../]
[./s11_alpha]
type = SideAverageValue
boundary = left
variable = sigma11
[../]
[./s22_alpha]
type = SideAverageValue
boundary = left
variable = sigma22
[../]
[./s33_alpha]
type = SideAverageValue
boundary = left
variable = sigma33
[../]
[./s11_beta]
type = SideAverageValue
boundary = right
variable = sigma11
[../]
[./s22_beta]
type = SideAverageValue
boundary = right
variable = sigma22
[../]
[./s33_beta]
type = SideAverageValue
boundary = right
variable = sigma33
[../]
[./f_el_alpha]
type = SideAverageValue
boundary = left
variable = f_el
[../]
[./f_el_beta]
type = SideAverageValue
boundary = right
variable = f_el
[../]
[./f_c_alpha]
type = SideAverageValue
boundary = left
variable = Fglobal
[../]
[./f_c_beta]
type = SideAverageValue
boundary = right
variable = Fglobal
[../]
[./chem_pot_alpha]
type = SideAverageValue
boundary = left
variable = w
[../]
[./chem_pot_beta]
type = SideAverageValue
boundary = right
variable = w
[../]
[./psi_alpha]
type = SideAverageValue
boundary = left
variable = psi
[../]
[./psi_beta]
type = SideAverageValue
boundary = right
variable = psi
[../]
[./total_energy]
type = ElementIntegralVariablePostprocessor
variable = Fglobal
[../]
# Get simulation cell size from postprocessor
[./volume]
type = ElementIntegralMaterialProperty
mat_prop = 1
[../]
[./psi_eq_int]
type = FunctionValuePostprocessor
function = psi_eq_int
[../]
[./psi_int]
type = ElementIntegralVariablePostprocessor
variable = psi
[../]
[./gamma]
type = FunctionValuePostprocessor
function = gamma
[../]
[./int_position]
type = FindValueOnLine
start_point = '-10 0 0'
end_point = '10 0 0'
v = eta
target = 0.5
[../]
[]
#
# Precondition using handcoded off-diagonal terms
#
[Preconditioning]
[./full]
type = SMP
full = true
[../]
[]
[Outputs]
[./exodus]
type = Exodus
time_step_interval = 20
[../]
checkpoint = true
[./csv]
type = CSV
execute_on = 'final'
[../]
[]
(modules/combined/examples/publications/rapid_dev/fig7b.i)
#
# Fig. 7 input for 10.1016/j.commatsci.2017.02.017
# D. Schwen et al./Computational Materials Science 132 (2017) 36-45
# Dashed black curve (2)
# Eigenstrain is globally applied. Single global elastic free energies.
# Supply the RADIUS parameter (10-35) on the command line to generate data
# for all curves in the plot.
#
[Mesh]
type = GeneratedMesh
dim = 1
nx = 32
xmin = 0
xmax = 100
second_order = true
coord_type = RSPHERICAL
[]
[GlobalParams]
displacements = 'disp_r'
[]
[Functions]
[./diff]
type = ParsedFunction
expression = '${RADIUS}-pos_c'
symbol_names = pos_c
symbol_values = pos_c
[../]
[]
# AuxVars to compute the free energy density for outputting
[AuxVariables]
[./local_energy]
order = CONSTANT
family = MONOMIAL
[../]
[./cross_energy]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./local_free_energy]
type = TotalFreeEnergy
variable = local_energy
interfacial_vars = 'c'
kappa_names = 'kappa_c'
execute_on = 'INITIAL TIMESTEP_END'
[../]
[]
[Variables]
# Solute concentration variable
[./c]
[./InitialCondition]
type = SmoothCircleIC
invalue = 1
outvalue = 0
x1 = 0
y1 = 0
radius = ${RADIUS}
int_width = 3
[../]
[../]
[./w]
[../]
# Phase order parameter
[./eta]
[./InitialCondition]
type = SmoothCircleIC
invalue = 1
outvalue = 0
x1 = 0
y1 = 0
radius = ${RADIUS}
int_width = 3
[../]
[../]
[./Fe_fit]
order = SECOND
[../]
[]
[Physics/SolidMechanics/QuasiStatic/all]
add_variables = true
eigenstrain_names = eigenstrain
[]
[Kernels]
# Split Cahn-Hilliard kernels
[./c_res]
type = SplitCHParsed
variable = c
f_name = F
args = 'eta'
kappa_name = kappa_c
w = w
[../]
[./wres]
type = SplitCHWRes
variable = w
mob_name = M
[../]
[./time]
type = CoupledTimeDerivative
variable = w
v = c
[../]
# Allen-Cahn and Lagrange-multiplier constraint kernels for order parameter 1
[./detadt]
type = TimeDerivative
variable = eta
[../]
[./ACBulk1]
type = AllenCahn
variable = eta
args = 'c'
mob_name = L
f_name = F
[../]
[./ACInterface]
type = ACInterface
variable = eta
mob_name = L
kappa_name = kappa_eta
[../]
[./Fe]
type = MaterialPropertyValue
prop_name = Fe
variable = Fe_fit
[../]
[./autoadjust]
type = MaskedBodyForce
variable = w
function = diff
mask = mask
[../]
[]
[Materials]
# declare a few constants, such as mobilities (L,M) and interface gradient prefactors (kappa*)
[./consts]
type = GenericConstantMaterial
prop_names = 'M L kappa_c kappa_eta'
prop_values = '1.0 1.0 0.5 1'
[../]
# forcing function mask
[./mask]
type = ParsedMaterial
property_name = mask
expression = grad/dt
material_property_names = 'grad dt'
[../]
[./grad]
type = VariableGradientMaterial
variable = c
prop = grad
[../]
[./time]
type = TimeStepMaterial
[../]
# global mechanical properties
[./elasticity_tensor]
type = ComputeElasticityTensor
C_ijkl = '1 1'
fill_method = symmetric_isotropic
[../]
[./stress]
type = ComputeLinearElasticStress
[../]
# eigenstrain as a function of phase
[./eigenstrain]
type = ComputeVariableEigenstrain
eigen_base = '0.05 0.05 0.05 0 0 0'
prefactor = h
args = eta
eigenstrain_name = eigenstrain
[../]
# switching functions
[./switching]
type = SwitchingFunctionMaterial
function_name = h
eta = eta
h_order = SIMPLE
[../]
[./barrier]
type = BarrierFunctionMaterial
eta = eta
[../]
# chemical free energies
[./chemical_free_energy_1]
type = DerivativeParsedMaterial
property_name = Fc1
expression = 'c^2'
coupled_variables = 'c'
derivative_order = 2
[../]
[./chemical_free_energy_2]
type = DerivativeParsedMaterial
property_name = Fc2
expression = '(1-c)^2'
coupled_variables = 'c'
derivative_order = 2
[../]
# global chemical free energy
[./chemical_free_energy]
type = DerivativeTwoPhaseMaterial
f_name = Fc
fa_name = Fc1
fb_name = Fc2
eta = eta
args = 'c'
W = 4
[../]
# global elastic free energy
[./elastic_free_energy]
type = ElasticEnergyMaterial
f_name = Fe
args = 'eta'
output_properties = Fe
outputs = 'all'
derivative_order = 2
[../]
# free energy
[./free_energy]
type = DerivativeSumMaterial
property_name = F
sum_materials = 'Fc Fe'
coupled_variables = 'c eta'
derivative_order = 2
[../]
[]
[BCs]
[./left_r]
type = DirichletBC
variable = disp_r
boundary = 'left'
value = 0
[../]
[]
[Preconditioning]
[./SMP]
type = SMP
full = true
[../]
[]
# We monitor the total free energy and the total solute concentration (should be constant)
[Postprocessors]
[./total_free_energy]
type = ElementIntegralVariablePostprocessor
variable = local_energy
execute_on = 'INITIAL TIMESTEP_END'
outputs = 'table console'
[../]
[./total_solute]
type = ElementIntegralVariablePostprocessor
variable = c
execute_on = 'INITIAL TIMESTEP_END'
outputs = 'table console'
[../]
[./pos_c]
type = FindValueOnLine
start_point = '0 0 0'
end_point = '100 0 0'
v = c
target = 0.582
tol = 1e-8
execute_on = 'INITIAL TIMESTEP_END'
outputs = 'table console'
[../]
[./pos_eta]
type = FindValueOnLine
start_point = '0 0 0'
end_point = '100 0 0'
v = eta
target = 0.5
tol = 1e-8
execute_on = 'INITIAL TIMESTEP_END'
outputs = 'table console'
[../]
[./c_min]
type = ElementExtremeValue
value_type = min
variable = c
execute_on = 'INITIAL TIMESTEP_END'
outputs = 'table console'
[../]
[]
[VectorPostprocessors]
[./line]
type = LineValueSampler
variable = 'Fe_fit c w'
start_point = '0 0 0'
end_point = '100 0 0'
num_points = 5000
sort_by = x
outputs = vpp
execute_on = 'INITIAL TIMESTEP_END'
[../]
[]
[Executioner]
type = Transient
scheme = bdf2
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type -sub_pc_type'
petsc_options_value = 'asm lu'
l_max_its = 30
nl_max_its = 15
l_tol = 1.0e-4
nl_rel_tol = 1.0e-8
nl_abs_tol = 2.0e-9
start_time = 0.0
end_time = 100000.0
[./TimeStepper]
type = IterationAdaptiveDT
optimal_iterations = 8
iteration_window = 1
dt = 1
[../]
[./Adaptivity]
initial_adaptivity = 5
interval = 10
max_h_level = 5
refine_fraction = 0.9
coarsen_fraction = 0.1
[../]
[]
[Outputs]
print_linear_residuals = false
perf_graph = true
execute_on = 'INITIAL TIMESTEP_END'
[./table]
type = CSV
delimiter = ' '
file_base = radius_${RADIUS}/eigenstrain_pp
[../]
[./vpp]
type = CSV
delimiter = ' '
sync_times = '10 50 100 500 1000 5000 10000 50000 100000'
sync_only = true
time_data = true
file_base = radius_${RADIUS}/eigenstrain_vpp
[../]
[]
(modules/combined/examples/publications/rapid_dev/fig7a.i)
#
# Fig. 7 input for 10.1016/j.commatsci.2017.02.017
# D. Schwen et al./Computational Materials Science 132 (2017) 36-45
# Solid gray curve (1)
# Eigenstrain and elastic energies ar computed per phase and then interpolated.
# Supply the RADIUS parameter (10-35) on the command line to generate data
# for all curves in the plot.
#
[Mesh]
type = GeneratedMesh
dim = 1
nx = 32
xmin = 0
xmax = 100
second_order = true
coord_type = RSPHERICAL
[]
[GlobalParams]
displacements = 'disp_r'
[]
[Functions]
[./diff]
type = ParsedFunction
expression = '${RADIUS}-pos_c'
symbol_names = pos_c
symbol_values = pos_c
[../]
[]
# AuxVars to compute the free energy density for outputting
[AuxVariables]
[./local_energy]
order = CONSTANT
family = MONOMIAL
[../]
[./cross_energy]
order = CONSTANT
family = MONOMIAL
[../]
[]
[AuxKernels]
[./local_free_energy]
type = TotalFreeEnergy
variable = local_energy
interfacial_vars = 'c'
kappa_names = 'kappa_c'
execute_on = 'INITIAL TIMESTEP_END'
[../]
[]
[Variables]
# Solute concentration variable
[./c]
[./InitialCondition]
type = SmoothCircleIC
invalue = 1
outvalue = 0
x1 = 0
y1 = 0
radius = ${RADIUS}
int_width = 3
[../]
[../]
[./w]
[../]
# Phase order parameter
[./eta]
[./InitialCondition]
type = SmoothCircleIC
invalue = 1
outvalue = 0
x1 = 0
y1 = 0
radius = ${RADIUS}
int_width = 3
[../]
[../]
# Mesh displacement
[./disp_r]
order = SECOND
[../]
[./Fe_fit]
order = SECOND
[../]
[]
[Kernels]
# Set up stress divergence kernels
[./TensorMechanics]
[../]
# Split Cahn-Hilliard kernels
[./c_res]
type = SplitCHParsed
variable = c
f_name = F
args = 'eta'
kappa_name = kappa_c
w = w
[../]
[./wres]
type = SplitCHWRes
variable = w
mob_name = M
[../]
[./time]
type = CoupledTimeDerivative
variable = w
v = c
[../]
# Allen-Cahn and Lagrange-multiplier constraint kernels for order parameter 1
[./detadt]
type = TimeDerivative
variable = eta
[../]
[./ACBulk1]
type = AllenCahn
variable = eta
args = 'c'
mob_name = L
f_name = F
[../]
[./ACInterface]
type = ACInterface
variable = eta
mob_name = L
kappa_name = kappa_eta
[../]
[./Fe]
type = MaterialPropertyValue
prop_name = Fe
variable = Fe_fit
[../]
[./autoadjust]
type = MaskedBodyForce
variable = w
function = diff
mask = mask
[../]
[]
[Materials]
# declare a few constants, such as mobilities (L,M) and interface gradient prefactors (kappa*)
[./consts]
type = GenericConstantMaterial
prop_names = 'M L kappa_c kappa_eta'
prop_values = '1.0 1.0 0.5 1'
[../]
# forcing function mask
[./mask]
type = ParsedMaterial
property_name = mask
expression = grad/dt
material_property_names = 'grad dt'
[../]
[./grad]
type = VariableGradientMaterial
variable = c
prop = grad
[../]
[./time]
type = TimeStepMaterial
[../]
# global mechanical properties
[./elasticity_tensor_1]
type = ComputeElasticityTensor
C_ijkl = '1 1'
base_name = phase1
fill_method = symmetric_isotropic
[../]
[./elasticity_tensor_2]
type = ComputeElasticityTensor
C_ijkl = '1 1'
base_name = phase2
fill_method = symmetric_isotropic
[../]
[./strain_1]
type = ComputeRSphericalSmallStrain
base_name = phase1
[../]
[./strain_2]
type = ComputeRSphericalSmallStrain
base_name = phase2
eigenstrain_names = eigenstrain
[../]
[./stress_1]
type = ComputeLinearElasticStress
base_name = phase1
[../]
[./stress_2]
type = ComputeLinearElasticStress
base_name = phase2
[../]
# eigenstrain per phase
[./eigenstrain2]
type = ComputeEigenstrain
eigen_base = '0.05 0.05 0.05 0 0 0'
base_name = phase2
eigenstrain_name = eigenstrain
[../]
# switching functions
[./switching]
type = SwitchingFunctionMaterial
function_name = h
eta = eta
h_order = SIMPLE
[../]
[./barrier]
type = BarrierFunctionMaterial
eta = eta
[../]
# chemical free energies
[./chemical_free_energy_1]
type = DerivativeParsedMaterial
property_name = Fc1
expression = 'c^2'
coupled_variables = 'c'
derivative_order = 2
[../]
[./chemical_free_energy_2]
type = DerivativeParsedMaterial
property_name = Fc2
expression = '(1-c)^2'
coupled_variables = 'c'
derivative_order = 2
[../]
# elastic free energies
[./elastic_free_energy_1]
type = ElasticEnergyMaterial
f_name = Fe1
args = ''
base_name = phase1
derivative_order = 2
[../]
[./elastic_free_energy_2]
type = ElasticEnergyMaterial
f_name = Fe2
args = ''
base_name = phase2
derivative_order = 2
[../]
# per phase free energies
[./free_energy_1]
type = DerivativeSumMaterial
property_name = F1
sum_materials = 'Fc1 Fe1'
coupled_variables = 'c'
derivative_order = 2
[../]
[./free_energy_2]
type = DerivativeSumMaterial
property_name = F2
sum_materials = 'Fc2 Fe2'
coupled_variables = 'c'
derivative_order = 2
[../]
# global chemical free energy
[./global_free_energy]
type = DerivativeTwoPhaseMaterial
f_name = F
fa_name = F1
fb_name = F2
eta = eta
args = 'c'
W = 4
[../]
# global stress
[./global_stress]
type = TwoPhaseStressMaterial
base_A = phase1
base_B = phase2
[../]
[./elastic_free_energy]
type = DerivativeTwoPhaseMaterial
f_name = Fe
fa_name = Fe1
fb_name = Fe2
eta = eta
args = 'c'
W = 0
[../]
[]
[BCs]
[./left_r]
type = DirichletBC
variable = disp_r
boundary = 'left'
value = 0
[../]
[]
[Preconditioning]
[./SMP]
type = SMP
full = true
[../]
[]
# We monitor the total free energy and the total solute concentration (should be constant)
[Postprocessors]
[./total_free_energy]
type = ElementIntegralVariablePostprocessor
variable = local_energy
execute_on = 'INITIAL TIMESTEP_END'
outputs = 'table console'
[../]
[./total_solute]
type = ElementIntegralVariablePostprocessor
variable = c
execute_on = 'INITIAL TIMESTEP_END'
outputs = 'table console'
[../]
[./pos_c]
type = FindValueOnLine
start_point = '0 0 0'
end_point = '100 0 0'
v = c
target = 0.582
tol = 1e-8
execute_on = 'INITIAL TIMESTEP_END'
outputs = 'table console'
[../]
[./pos_eta]
type = FindValueOnLine
start_point = '0 0 0'
end_point = '100 0 0'
v = eta
target = 0.5
tol = 1e-8
execute_on = 'INITIAL TIMESTEP_END'
outputs = 'table console'
[../]
[./c_min]
type = ElementExtremeValue
value_type = min
variable = c
execute_on = 'INITIAL TIMESTEP_END'
outputs = 'table console'
[../]
[]
[VectorPostprocessors]
[./line]
type = LineValueSampler
variable = 'Fe_fit c w'
start_point = '0 0 0'
end_point = '100 0 0'
num_points = 5000
sort_by = x
outputs = vpp
execute_on = 'INITIAL TIMESTEP_END'
[../]
[]
[Executioner]
type = Transient
scheme = bdf2
solve_type = 'PJFNK'
petsc_options_iname = '-pc_type -sub_pc_type'
petsc_options_value = 'asm lu'
l_max_its = 30
nl_max_its = 15
l_tol = 1.0e-4
nl_rel_tol = 1.0e-8
nl_abs_tol = 2.0e-9
start_time = 0.0
end_time = 100000.0
[./TimeStepper]
type = IterationAdaptiveDT
optimal_iterations = 7
iteration_window = 1
dt = 1
[../]
[./Adaptivity]
initial_adaptivity = 5
interval = 10
max_h_level = 5
refine_fraction = 0.9
coarsen_fraction = 0.1
[../]
[]
[Outputs]
print_linear_residuals = false
perf_graph = true
execute_on = 'INITIAL TIMESTEP_END'
[./table]
type = CSV
delimiter = ' '
file_base = radius_${RADIUS}/energy_pp
[../]
[./vpp]
type = CSV
delimiter = ' '
sync_times = '10 50 100 500 1000 5000 10000 50000 100000'
sync_only = true
time_data = true
file_base = radius_${RADIUS}/energy_vpp
[../]
[]
(modules/combined/test/tests/phase_field_contact_angle/contact_angle.i)
sigma = 25e-3 #10e-3 #25e-3 #surface tension coefficient
epsilon = 1e-6 #width parameter
nu = 1e-4#mobility parameter
contactangle = 2.61799#0.523599#1.0472
lambda = ${fparse 3*sigma*epsilon/(2*sqrt(2))}
prefactor_phi = ${fparse nu*lambda/(epsilon*epsilon)}
prefactor_psi = ${fparse -epsilon*epsilon}
coeff = ${fparse lambda/(epsilon*epsilon)}
[Mesh]
[gen]
type = GeneratedMeshGenerator
dim = 2
xmin = 0
xmax = 0.2e-3
ymin = 0
ymax = 0.2e-3
nx = 20
ny = 20
elem_type = QUAD9
[]
[]
[ICs]
[pf_ic]
type = BoundingBoxIC
variable = pf
x1 = 0.1e-3
y1 = -0.1e-3
x2 = 0.3e-03
y2 = 0.3e-3
inside = 1
outside = -1
int_width = ${fparse 2*sqrt(2)*epsilon}
[]
[velocity]
type = VectorConstantIC
x_value = 0.0
y_value = 0.0
variable = velocity
[]
[]
[Variables]
[pf]
family = LAGRANGE
order = second
[]
[auxpf]
family = LAGRANGE
order = second
[]
[velocity]
family = LAGRANGE_VEC
[]
[p]
[]
[]
[Kernels]
[mass]
type = INSADMass
variable = p
[]
[mass_pspg]
type = INSADMassPSPG
variable = p
[]
[momentum_viscous]
type = INSADMomentumViscous
variable = velocity
viscous_form = 'traction'
mu_name = 'mu'
[]
[momentum_pressure]
type = INSADMomentumPressure
variable = velocity
pressure = p
integrate_p_by_parts = false
[]
[surface_tension]
type = ADPhaseFieldTwoPhaseSurfaceTension
variable = velocity
pf = pf
auxpf = auxpf
coeff = ${coeff}
[]
[phasefield_timederivative]
type = ADTimeDerivative
variable = pf
[]
[phasefield_supg]
type = ADPhaseFieldTimeDerivativeSUPG
velocity = velocity
variable = pf
[]
[phasefield_laplacian]
type=ADPrefactorLaplacianSplit
variable = pf
c = auxpf
prefactor = ${prefactor_phi}
[]
[Auxphasefield_firstorder]
type=ADReaction
variable = auxpf
rate = 1.0
[]
[Auxphasefield_laplacian]
type=ADPrefactorLaplacianSplit
variable = auxpf
c = pf
prefactor=${prefactor_psi}
[]
[Auxphasefield_doublewell]
type=ADPhaseFieldCoupledDoubleWellPotential
variable = auxpf
c = pf
prefactor=-1.0
[]
[]
[BCs]
[no_slip]
type = VectorFunctionDirichletBC
variable = velocity
boundary = 'top bottom'
[]
[velocity_L]
type = VectorFunctionDirichletBC
variable = velocity
boundary = 'left'
[]
[velocity_R]
type = VectorFunctionDirichletBC
variable = velocity
boundary = 'right'
[]
[ContactangleBC]
type=ADPhaseFieldContactAngleBC
variable = auxpf
pf = pf
epsilon = ${epsilon}
lambda=${lambda}
sigma=${sigma}
contactangle=${contactangle}
boundary = 'top bottom'
[]
[]
[Materials]
[rho]
type = ADPhaseFieldTwoPhaseMaterial
prop_name = rho
prop_value_1 = 1000
prop_value_2 = 840
pf = pf
# outputs = exodus
[]
[mu]
type = ADPhaseFieldTwoPhaseMaterial
prop_name = mu
prop_value_1 = 1e-3
prop_value_2 = 7.6e-3
pf = pf
# outputs = exodus
[]
[ins_mat]
type = INSADTauMaterial
velocity = velocity
pressure = p
alpha = .1
[]
[]
[Postprocessors]
[contact_angle_top]
type = ObtainAvgContactAngle
boundary = top
pf=pf
execute_on = 'timestep_end'
[]
[x_position]
type = FindValueOnLine
start_point = '0 0.0001 0'
end_point ='0.0002 0.0001 0'
v = pf
target = 0.0
tol = 1e-6
[]
[]
[Preconditioning]
[SMP]
type = SMP
full = true
[]
[]
[Adaptivity]
initial_steps = 2
initial_marker = phase_marker
marker = phase_marker
max_h_level = 4
[Markers]
[phase_marker]
type = ValueRangeMarker
lower_bound = -0.99
upper_bound = 0.99
variable = pf
[]
[]
[]
[Executioner]
type = Transient
solve_type = NEWTON
start_time = 0
num_steps = 5
dtmax = 0.25
[TimeStepper]
type = IterationAdaptiveDT
dt = 1e-10
iteration_window = 2
optimal_iterations = 10
growth_factor = 2
cutback_factor = 0.5
[]
# petsc_options_iname = '-pc_type -ksp_gmres_restart -pc_factor_mat_solver_type -pc_factor_shift_type -pc_factor_shift_amount'
# petsc_options_value = 'lu 50 superlu_dist NONZERO 1e-15'
#petsc_options_iname = '-pc_type'
#petsc_options_value = 'lu '
petsc_options_iname = '-pc_type -pc_factor_shift_type -pc_factor_shift_amount'
petsc_options_value = 'lu NONZERO 1e-15 '
line_search = 'none'
nl_rel_tol = 1e-6
nl_abs_tol = 1e-8
nl_max_its = 20
nl_forced_its = 3
l_tol = 1e-6
l_max_its = 20
[]
[Outputs]
[csv]
type = CSV
time_step_interval = 1
[]
[]
(test/tests/postprocessors/find_value_on_line/findvalueonline.i)
[Mesh]
type = GeneratedMesh
dim = 1
nx = 20
xmin = 0
xmax = 10
[]
[Variables]
[./phi]
[./InitialCondition]
type = FunctionIC
function = if(x<1,1-x,0)
[../]
[../]
[]
[Kernels]
[./diff]
type = Diffusion
variable = phi
[../]
[./dt]
type = TimeDerivative
variable = phi
[../]
[]
[BCs]
[./influx]
type = NeumannBC
boundary = left
variable = phi
value = 1
[../]
[./fix]
type = DirichletBC
boundary = right
variable = phi
value = 0
[../]
[]
[Postprocessors]
[./pos]
type = FindValueOnLine
target = 0.5
v = phi
start_point = '0 0 0'
end_point = '10 0 0'
execute_on = 'initial timestep_end'
[../]
[]
[Executioner]
type = Transient
num_steps = 5
dt = 2.5
[]
[Outputs]
csv = true
[]