Sublattice KKS model
The sublattice Kim-Kim-Suzuki (SLKKS) Schwen et al. (2021) model is an extension of the original KKS model incorporating an additional equal chemical potential constraint among the sublattice concentration in each phase.
As such each component in the system will have corresponding phase concentrations, which are split up into per-sublattice concentrations.
Nomenclature
indexes a component, a phase, and a sublattice in the given phase. is the fraction of sublattice sites in phase .
Sublattice equilibrium
All sublattice pairs and in a given phase are assumed to be always in equilibrium, thus
This condition is enforced by the SLKKSChemicalPotential
kernel.
[chempot1a1b]
type = SLKKSChemicalPotential
variable = Cijk
k = ajk
cs = Cijk'
ks = ajk'
F = Fj
[]
With sublattices in a phase such kernels are required. That leaves one sublattice concentration variable in a given phase to be covered by a kernel.
Phase equilibrium
At the same time the original KKS model chemical potential equilibria still hold
meaning that sublattice chemical potentials (corrected for the sublattice site fractions) must be in equilibrium across phases.
That remaining sublattice concentration will couple to the next phase using a KKSPhaseChemicalPotential
kernel.
[chempot1c2a]
type = KKSPhaseChemicalPotential
variable = Cijk
ka = ajk
fa_name = Fj
cb = Cij'k
kb = aj'k
fb_name = Fj'
[]
Note that in this kernel the ajk
and aj'k
must be true site fractions ranging from 0 to 1.
Mass transport
The global concentration variables govern the mass transport along chemical potential gradients.
(1)
This equation can be implemented using multiple MatDiffusion
kernels in MOOSE. To illustrate this we can re-order the summation to yield
(2)
This means we need to add a MatDiffusion
kernel for sublattice in all phases, each operating on the variable . We use the "v" parameter to specify as the variable to take the gradient of (rather than ). The effective diffusivity passed into the kernel is , which can be provided by a DerivativeParsedMaterial
[D'i]
type = DerivativeParsedMaterial
f_name = D'i
function = Di*hi*ajk
material_property_names = 'Di hi'
constant_names = ajk
constant_expressions = 1/3
[]
Here we assume the site fraction of the current sublattice to be .
The derivation should be check to see if a phase dependent concentration or even sublattice dependent concentration is feasible.
Example
An open TDB file for the Fe-Cr system was graciously provided by Aurélie Jacob at TU Wien, based on the work in Jacob et al. (2018).
The TDB file can be converted into MOOSE DerivativeParsedMaterial
syntax using the free_energy.py
script.
mamba activate pycalphad
cd $MOOSE_DIR/modules/phase_field/examples/slkks
../../../../python/calphad/free_energy.py CrFe_Jacob.tdb BCC_A2 SIGMA
The example directory also contains an Jupyter notebook to visualize the phase diagram for the supplied Fe-Cr thermodynamic database file and the simulation results,
#
# SLKKS two phase example for the BCC and SIGMA phases. The sigma phase contains
# multiple sublattices. Free energy from
# Jacob, Aurelie, Erwin Povoden-Karadeniz, and Ernst Kozeschnik. "Revised thermodynamic
# description of the Fe-Cr system based on an improved sublattice model of the sigma phase."
# Calphad 60 (2018): 16-28.
#
# In this simulation we consider diffusion (Cahn-Hilliard) and phase transformation.
#
# This example requires CrFe_sigma_out_var_0001.csv file, which generated by first
# running the CrFe_sigma.i input file.
[Mesh<<<{"href": "../../../syntax/Mesh/index.html"}>>>]
[gen]
type = GeneratedMeshGenerator<<<{"description": "Create a line, square, or cube mesh with uniformly spaced or biased elements.", "href": "../../../source/meshgenerators/GeneratedMeshGenerator.html"}>>>
dim<<<{"description": "The dimension of the mesh to be generated"}>>> = 2
nx<<<{"description": "Number of elements in the X direction"}>>> = 160
ny<<<{"description": "Number of elements in the Y direction"}>>> = 1
nz<<<{"description": "Number of elements in the Z direction"}>>> = 0
xmin<<<{"description": "Lower X Coordinate of the generated mesh"}>>> = -25
xmax<<<{"description": "Upper X Coordinate of the generated mesh"}>>> = 25
ymin<<<{"description": "Lower Y Coordinate of the generated mesh"}>>> = -2.5
ymax<<<{"description": "Upper Y Coordinate of the generated mesh"}>>> = 2.5
elem_type<<<{"description": "The type of element from libMesh to generate (default: linear element for requested dimension)"}>>> = QUAD4
[]
[]
[AuxVariables<<<{"href": "../../../syntax/AuxVariables/index.html"}>>>]
[Fglobal]
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = CONSTANT
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = MONOMIAL
[]
[]
[Functions<<<{"href": "../../../syntax/Functions/index.html"}>>>]
[sigma_cr0]
type = PiecewiseLinear<<<{"description": "Linearly interpolates between pairs of x-y data", "href": "../../../source/functions/PiecewiseLinear.html"}>>>
data_file<<<{"description": "File holding CSV data"}>>> = CrFe_sigma_out_var_0001.csv
format<<<{"description": "Format of csv data file that is in either in columns or rows"}>>> = columns
x_index_in_file<<<{"description": "The abscissa index in the data file"}>>> = 5
y_index_in_file<<<{"description": "The ordinate index in the data file"}>>> = 2
xy_in_file_only<<<{"description": "If the data file only contains abscissa and ordinate data"}>>> = false
[]
[sigma_cr1]
type = PiecewiseLinear<<<{"description": "Linearly interpolates between pairs of x-y data", "href": "../../../source/functions/PiecewiseLinear.html"}>>>
data_file<<<{"description": "File holding CSV data"}>>> = CrFe_sigma_out_var_0001.csv
format<<<{"description": "Format of csv data file that is in either in columns or rows"}>>> = columns
x_index_in_file<<<{"description": "The abscissa index in the data file"}>>> = 5
y_index_in_file<<<{"description": "The ordinate index in the data file"}>>> = 3
xy_in_file_only<<<{"description": "If the data file only contains abscissa and ordinate data"}>>> = false
[]
[sigma_cr2]
type = PiecewiseLinear<<<{"description": "Linearly interpolates between pairs of x-y data", "href": "../../../source/functions/PiecewiseLinear.html"}>>>
data_file<<<{"description": "File holding CSV data"}>>> = CrFe_sigma_out_var_0001.csv
format<<<{"description": "Format of csv data file that is in either in columns or rows"}>>> = columns
x_index_in_file<<<{"description": "The abscissa index in the data file"}>>> = 5
y_index_in_file<<<{"description": "The ordinate index in the data file"}>>> = 4
xy_in_file_only<<<{"description": "If the data file only contains abscissa and ordinate data"}>>> = false
[]
[]
[Variables<<<{"href": "../../../syntax/Variables/index.html"}>>>]
# order parameters
[eta1]
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = FIRST
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = LAGRANGE
initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 0.5
[]
[eta2]
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = FIRST
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = LAGRANGE
initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 0.5
[]
# solute concentration
[cCr]
order<<<{"description": "Specifies the order of the FE shape function to use for this variable (additional orders not listed are allowed)"}>>> = FIRST
family<<<{"description": "Specifies the family of FE shape functions to use for this variable"}>>> = LAGRANGE
[InitialCondition]
type = FunctionIC
function = '(x+25)/50*0.5+0.1'
[]
[]
# sublattice concentrations
[BCC_CR]
initial_condition<<<{"description": "Specifies a constant initial condition for this variable"}>>> = 0.45
[]
[SIGMA_0CR]
[InitialCondition]
type = CoupledValueFunctionIC
function = sigma_cr0
v = cCr
variable = SIGMA_0CR
[]
[]
[SIGMA_1CR]
[InitialCondition]
type = CoupledValueFunctionIC
function = sigma_cr1
v = cCr
variable = SIGMA_1CR
[]
[]
[SIGMA_2CR]
[InitialCondition]
type = CoupledValueFunctionIC
function = sigma_cr2
v = cCr
variable = SIGMA_2CR
[]
[]
# Lagrange multiplier
[lambda]
[]
[]
[Materials<<<{"href": "../../../syntax/Materials/index.html"}>>>]
# CALPHAD free energies
[F_BCC_A2]
type = DerivativeParsedMaterial<<<{"description": "Parsed Function Material with automatic derivatives.", "href": "../../../source/materials/DerivativeParsedMaterial.html"}>>>
property_name<<<{"description": "Name of the parsed material property"}>>> = F_BCC_A2
outputs<<<{"description": "Vector of output names where you would like to restrict the output of variables(s) associated with this object"}>>> = exodus
output_properties<<<{"description": "List of material properties, from this material, to output (outputs must also be defined to an output type)"}>>> = F_BCC_A2
expression<<<{"description": "Parsed function (see FParser) expression for the parsed material"}>>> = 'BCC_FE:=1-BCC_CR; G := 8.3145*T*(1.0*if(BCC_CR > 1.0e-15,BCC_CR*log(BCC_CR),0) + '
'1.0*if(BCC_FE > 1.0e-15,BCC_FE*plog(BCC_FE,eps),0) + 3.0*if(BCC_VA > '
'1.0e-15,BCC_VA*log(BCC_VA),0))/(BCC_CR + BCC_FE) + 8.3145*T*if(T < '
'548.2*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) - 932.5*BCC_CR*BCC_FE*BCC_VA + '
'311.5*BCC_CR*BCC_VA - '
'1043.0*BCC_FE*BCC_VA,-8.13674105561218e-49*T^15/(0.525599232981783*BCC_CR*BCC_FE*BCC_'
'VA*(BCC_CR - BCC_FE) - 0.894055608820709*BCC_CR*BCC_FE*BCC_VA + '
'0.298657718120805*BCC_CR*BCC_VA - BCC_FE*BCC_VA + 9.58772770853308e-13)^15 - '
'4.65558036243985e-30*T^9/(0.525599232981783*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) - '
'0.894055608820709*BCC_CR*BCC_FE*BCC_VA + 0.298657718120805*BCC_CR*BCC_VA - '
'BCC_FE*BCC_VA + 9.58772770853308e-13)^9 - '
'1.3485349181899e-10*T^3/(0.525599232981783*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) - '
'0.894055608820709*BCC_CR*BCC_FE*BCC_VA + 0.298657718120805*BCC_CR*BCC_VA - '
'BCC_FE*BCC_VA + 9.58772770853308e-13)^3 + 1 - '
'0.905299382744392*(548.2*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) - '
'932.5*BCC_CR*BCC_FE*BCC_VA + 311.5*BCC_CR*BCC_VA - 1043.0*BCC_FE*BCC_VA + '
'1.0e-9)/T,if(T < -548.2*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) + '
'932.5*BCC_CR*BCC_FE*BCC_VA - 311.5*BCC_CR*BCC_VA + '
'1043.0*BCC_FE*BCC_VA,-8.13674105561218e-49*T^15/(-0.525599232981783*BCC_CR*BCC_FE*BCC'
'_VA*(BCC_CR - BCC_FE) + 0.894055608820709*BCC_CR*BCC_FE*BCC_VA - '
'0.298657718120805*BCC_CR*BCC_VA + BCC_FE*BCC_VA + 9.58772770853308e-13)^15 - '
'4.65558036243985e-30*T^9/(-0.525599232981783*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) '
'+ 0.894055608820709*BCC_CR*BCC_FE*BCC_VA - 0.298657718120805*BCC_CR*BCC_VA + '
'BCC_FE*BCC_VA + 9.58772770853308e-13)^9 - '
'1.3485349181899e-10*T^3/(-0.525599232981783*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) + '
'0.894055608820709*BCC_CR*BCC_FE*BCC_VA - 0.298657718120805*BCC_CR*BCC_VA + '
'BCC_FE*BCC_VA + 9.58772770853308e-13)^3 + 1 - '
'0.905299382744392*(-548.2*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) + '
'932.5*BCC_CR*BCC_FE*BCC_VA - 311.5*BCC_CR*BCC_VA + 1043.0*BCC_FE*BCC_VA + '
'1.0e-9)/T,if(T > -548.2*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) + '
'932.5*BCC_CR*BCC_FE*BCC_VA - 311.5*BCC_CR*BCC_VA + 1043.0*BCC_FE*BCC_VA & '
'548.2*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) - 932.5*BCC_CR*BCC_FE*BCC_VA + '
'311.5*BCC_CR*BCC_VA - 1043.0*BCC_FE*BCC_VA < '
'0,-79209031311018.7*(-0.525599232981783*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) + '
'0.894055608820709*BCC_CR*BCC_FE*BCC_VA - 0.298657718120805*BCC_CR*BCC_VA + '
'BCC_FE*BCC_VA + 9.58772770853308e-13)^5/T^5 - '
'3.83095660520737e+42*(-0.525599232981783*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) + '
'0.894055608820709*BCC_CR*BCC_FE*BCC_VA - 0.298657718120805*BCC_CR*BCC_VA + '
'BCC_FE*BCC_VA + 9.58772770853308e-13)^15/T^15 - '
'1.22565886734485e+72*(-0.525599232981783*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) + '
'0.894055608820709*BCC_CR*BCC_FE*BCC_VA - 0.298657718120805*BCC_CR*BCC_VA + '
'BCC_FE*BCC_VA + 9.58772770853308e-13)^25/T^25,if(T > '
'548.2*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) - 932.5*BCC_CR*BCC_FE*BCC_VA + '
'311.5*BCC_CR*BCC_VA - 1043.0*BCC_FE*BCC_VA & 548.2*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - '
'BCC_FE) - 932.5*BCC_CR*BCC_FE*BCC_VA + 311.5*BCC_CR*BCC_VA - 1043.0*BCC_FE*BCC_VA > '
'0,-79209031311018.7*(0.525599232981783*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) - '
'0.894055608820709*BCC_CR*BCC_FE*BCC_VA + 0.298657718120805*BCC_CR*BCC_VA - '
'BCC_FE*BCC_VA + 9.58772770853308e-13)^5/T^5 - '
'3.83095660520737e+42*(0.525599232981783*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) - '
'0.894055608820709*BCC_CR*BCC_FE*BCC_VA + 0.298657718120805*BCC_CR*BCC_VA - '
'BCC_FE*BCC_VA + 9.58772770853308e-13)^15/T^15 - '
'1.22565886734485e+72*(0.525599232981783*BCC_CR*BCC_FE*BCC_VA*(BCC_CR - BCC_FE) - '
'0.894055608820709*BCC_CR*BCC_FE*BCC_VA + 0.298657718120805*BCC_CR*BCC_VA - '
'BCC_FE*BCC_VA + 9.58772770853308e-13)^25/T^25,0))))*log((2.15*BCC_CR*BCC_FE*BCC_VA - '
'0.008*BCC_CR*BCC_VA + 2.22*BCC_FE*BCC_VA)*if(2.15*BCC_CR*BCC_FE*BCC_VA - '
'0.008*BCC_CR*BCC_VA + 2.22*BCC_FE*BCC_VA <= 0,-1.0,1.0) + 1)/(BCC_CR + BCC_FE) + '
'1.0*(BCC_CR*BCC_VA*if(T >= 298.15 & T < 2180.0,139250.0*1/T - 26.908*T*log(T) + '
'157.48*T + 0.00189435*T^2.0 - 1.47721e-6*T^3.0 - 8856.94,if(T >= 2180.0 & T < '
'6000.0,-2.88526e+32*T^(-9.0) - 50.0*T*log(T) + 344.18*T - 34869.344,0)) + '
'BCC_FE*BCC_VA*if(T >= 298.15 & T < 1811.0,77358.5*1/T - 23.5143*T*log(T) + 124.134*T '
'- 0.00439752*T^2.0 - 5.89269e-8*T^3.0 + 1225.7,if(T >= 1811.0 & T < '
'6000.0,2.2960305e+31*T^(-9.0) - 46.0*T*log(T) + 299.31255*T - 25383.581,0)))/(BCC_CR '
'+ BCC_FE) + 1.0*(BCC_CR*BCC_FE*BCC_VA*(500.0 - 1.5*T)*(BCC_CR - BCC_FE) + '
'BCC_CR*BCC_FE*BCC_VA*(24600.0 - 14.98*T) + BCC_CR*BCC_FE*BCC_VA*(9.15*T - '
'14000.0)*(BCC_CR - BCC_FE)^2)/(BCC_CR + BCC_FE); G/100000'
coupled_variables<<<{"description": "Vector of variables used in the parsed function"}>>> = 'BCC_CR'
constant_names<<<{"description": "Vector of constants used in the parsed function (use this for kB etc.)"}>>> = 'BCC_VA T eps'
constant_expressions<<<{"description": "Vector of values for the constants in constant_names (can be an FParser expression)"}>>> = '1 1000 0.01'
[]
[F_SIGMA]
type = DerivativeParsedMaterial<<<{"description": "Parsed Function Material with automatic derivatives.", "href": "../../../source/materials/DerivativeParsedMaterial.html"}>>>
property_name<<<{"description": "Name of the parsed material property"}>>> = F_SIGMA
outputs<<<{"description": "Vector of output names where you would like to restrict the output of variables(s) associated with this object"}>>> = exodus
output_properties<<<{"description": "List of material properties, from this material, to output (outputs must also be defined to an output type)"}>>> = F_SIGMA
expression<<<{"description": "Parsed function (see FParser) expression for the parsed material"}>>> = 'SIGMA_0FE := 1-SIGMA_0CR; SIGMA_1FE := 1-SIGMA_1CR; SIGMA_2FE := 1-SIGMA_2CR; G := '
'8.3145*T*(10.0*if(SIGMA_0CR > 1.0e-15,SIGMA_0CR*plog(SIGMA_0CR,eps),0) + '
'10.0*if(SIGMA_0FE > 1.0e-15,SIGMA_0FE*plog(SIGMA_0FE,eps),0) + 4.0*if(SIGMA_1CR > '
'1.0e-15,SIGMA_1CR*plog(SIGMA_1CR,eps),0) + 4.0*if(SIGMA_1FE > '
'1.0e-15,SIGMA_1FE*plog(SIGMA_1FE,eps),0) + 16.0*if(SIGMA_2CR > '
'1.0e-15,SIGMA_2CR*plog(SIGMA_2CR,eps),0) + 16.0*if(SIGMA_2FE > '
'1.0e-15,SIGMA_2FE*plog(SIGMA_2FE,eps),0))/(10.0*SIGMA_0CR + 10.0*SIGMA_0FE + '
'4.0*SIGMA_1CR + 4.0*SIGMA_1FE + 16.0*SIGMA_2CR + 16.0*SIGMA_2FE) + '
'(SIGMA_0FE*SIGMA_1CR*SIGMA_2CR*SIGMA_2FE*(-70.0*T - 170400.0) + '
'SIGMA_0FE*SIGMA_1FE*SIGMA_2CR*SIGMA_2FE*(-10.0*T - 330839.0))/(10.0*SIGMA_0CR + '
'10.0*SIGMA_0FE + 4.0*SIGMA_1CR + 4.0*SIGMA_1FE + 16.0*SIGMA_2CR + 16.0*SIGMA_2FE) + '
'(SIGMA_0CR*SIGMA_1CR*SIGMA_2CR*(30.0*if(T >= 298.15 & T < 2180.0,139250.0*1/T - '
'26.908*T*log(T) + 157.48*T + 0.00189435*T^2.0 - 1.47721e-6*T^3.0 - 8856.94,if(T >= '
'2180.0 & T < 6000.0,-2.88526e+32*T^(-9.0) - 50.0*T*log(T) + 344.18*T - 34869.344,0)) '
'+ 132000.0) + SIGMA_0CR*SIGMA_1CR*SIGMA_2FE*(-110.0*T + 16.0*if(T >= 298.15 & T < '
'1811.0,77358.5*1/T - 23.5143*T*log(T) + 124.134*T - 0.00439752*T^2.0 - '
'5.89269e-8*T^3.0 + 1225.7,if(T >= 1811.0 & T < 6000.0,2.2960305e+31*T^(-9.0) - '
'46.0*T*log(T) + 299.31255*T - 25383.581,0)) + 14.0*if(T >= 298.15 & T < '
'2180.0,139250.0*1/T - 26.908*T*log(T) + 157.48*T + 0.00189435*T^2.0 - '
'1.47721e-6*T^3.0 - 8856.94,if(T >= 2180.0 & T < 6000.0,-2.88526e+32*T^(-9.0) - '
'50.0*T*log(T) + 344.18*T - 34869.344,0)) + 123500.0) + '
'SIGMA_0CR*SIGMA_1FE*SIGMA_2CR*(4.0*if(T >= 298.15 & T < 1811.0,77358.5*1/T - '
'23.5143*T*log(T) + 124.134*T - 0.00439752*T^2.0 - 5.89269e-8*T^3.0 + 1225.7,if(T >= '
'1811.0 & T < 6000.0,2.2960305e+31*T^(-9.0) - 46.0*T*log(T) + 299.31255*T - '
'25383.581,0)) + 26.0*if(T >= 298.15 & T < 2180.0,139250.0*1/T - 26.908*T*log(T) + '
'157.48*T + 0.00189435*T^2.0 - 1.47721e-6*T^3.0 - 8856.94,if(T >= 2180.0 & T < '
'6000.0,-2.88526e+32*T^(-9.0) - 50.0*T*log(T) + 344.18*T - 34869.344,0)) + 140486.0) '
'+ SIGMA_0CR*SIGMA_1FE*SIGMA_2FE*(20.0*if(T >= 298.15 & T < 1811.0,77358.5*1/T - '
'23.5143*T*log(T) + 124.134*T - 0.00439752*T^2.0 - 5.89269e-8*T^3.0 + 1225.7,if(T >= '
'1811.0 & T < 6000.0,2.2960305e+31*T^(-9.0) - 46.0*T*log(T) + 299.31255*T - '
'25383.581,0)) + 10.0*if(T >= 298.15 & T < 2180.0,139250.0*1/T - 26.908*T*log(T) + '
'157.48*T + 0.00189435*T^2.0 - 1.47721e-6*T^3.0 - 8856.94,if(T >= 2180.0 & T < '
'6000.0,-2.88526e+32*T^(-9.0) - 50.0*T*log(T) + 344.18*T - 34869.344,0)) + 148800.0) '
'+ SIGMA_0FE*SIGMA_1CR*SIGMA_2CR*(10.0*if(T >= 298.15 & T < 1811.0,77358.5*1/T - '
'23.5143*T*log(T) + 124.134*T - 0.00439752*T^2.0 - 5.89269e-8*T^3.0 + 1225.7,if(T >= '
'1811.0 & T < 6000.0,2.2960305e+31*T^(-9.0) - 46.0*T*log(T) + 299.31255*T - '
'25383.581,0)) + 20.0*if(T >= 298.15 & T < 2180.0,139250.0*1/T - 26.908*T*log(T) + '
'157.48*T + 0.00189435*T^2.0 - 1.47721e-6*T^3.0 - 8856.94,if(T >= 2180.0 & T < '
'6000.0,-2.88526e+32*T^(-9.0) - 50.0*T*log(T) + 344.18*T - 34869.344,0)) + 56200.0) + '
'SIGMA_0FE*SIGMA_1CR*SIGMA_2FE*(26.0*if(T >= 298.15 & T < 1811.0,77358.5*1/T - '
'23.5143*T*log(T) + 124.134*T - 0.00439752*T^2.0 - 5.89269e-8*T^3.0 + 1225.7,if(T >= '
'1811.0 & T < 6000.0,2.2960305e+31*T^(-9.0) - 46.0*T*log(T) + 299.31255*T - '
'25383.581,0)) + 4.0*if(T >= 298.15 & T < 2180.0,139250.0*1/T - 26.908*T*log(T) + '
'157.48*T + 0.00189435*T^2.0 - 1.47721e-6*T^3.0 - 8856.94,if(T >= 2180.0 & T < '
'6000.0,-2.88526e+32*T^(-9.0) - 50.0*T*log(T) + 344.18*T - 34869.344,0)) + 152700.0) '
'+ SIGMA_0FE*SIGMA_1FE*SIGMA_2CR*(14.0*if(T >= 298.15 & T < 1811.0,77358.5*1/T - '
'23.5143*T*log(T) + 124.134*T - 0.00439752*T^2.0 - 5.89269e-8*T^3.0 + 1225.7,if(T >= '
'1811.0 & T < 6000.0,2.2960305e+31*T^(-9.0) - 46.0*T*log(T) + 299.31255*T - '
'25383.581,0)) + 16.0*if(T >= 298.15 & T < 2180.0,139250.0*1/T - 26.908*T*log(T) + '
'157.48*T + 0.00189435*T^2.0 - 1.47721e-6*T^3.0 - 8856.94,if(T >= 2180.0 & T < '
'6000.0,-2.88526e+32*T^(-9.0) - 50.0*T*log(T) + 344.18*T - 34869.344,0)) + 46200.0) + '
'SIGMA_0FE*SIGMA_1FE*SIGMA_2FE*(30.0*if(T >= 298.15 & T < 1811.0,77358.5*1/T - '
'23.5143*T*log(T) + 124.134*T - 0.00439752*T^2.0 - 5.89269e-8*T^3.0 + 1225.7,if(T >= '
'1811.0 & T < 6000.0,2.2960305e+31*T^(-9.0) - 46.0*T*log(T) + 299.31255*T - '
'25383.581,0)) + 173333.0))/(10.0*SIGMA_0CR + 10.0*SIGMA_0FE + 4.0*SIGMA_1CR + '
'4.0*SIGMA_1FE + 16.0*SIGMA_2CR + 16.0*SIGMA_2FE); G/100000'
coupled_variables<<<{"description": "Vector of variables used in the parsed function"}>>> = 'SIGMA_0CR SIGMA_1CR SIGMA_2CR'
constant_names<<<{"description": "Vector of constants used in the parsed function (use this for kB etc.)"}>>> = 'T eps'
constant_expressions<<<{"description": "Vector of values for the constants in constant_names (can be an FParser expression)"}>>> = '1000 0.01'
[]
# h(eta)
[h1]
type = SwitchingFunctionMaterial<<<{"description": "Helper material to provide $h(\\eta)$ and its derivative in one of two polynomial forms.\nSIMPLE: $3\\eta^2-2\\eta^3$\nHIGH: $\\eta^3(6\\eta^2-15\\eta+10)$", "href": "../../../source/materials/SwitchingFunctionMaterial.html"}>>>
function_name<<<{"description": "actual name for f(eta), i.e. 'h' or 'g'"}>>> = h1
h_order<<<{"description": "Polynomial order of the switching function h(eta)"}>>> = HIGH
eta<<<{"description": "Order parameter variable"}>>> = eta1
[]
[h2]
type = SwitchingFunctionMaterial<<<{"description": "Helper material to provide $h(\\eta)$ and its derivative in one of two polynomial forms.\nSIMPLE: $3\\eta^2-2\\eta^3$\nHIGH: $\\eta^3(6\\eta^2-15\\eta+10)$", "href": "../../../source/materials/SwitchingFunctionMaterial.html"}>>>
function_name<<<{"description": "actual name for f(eta), i.e. 'h' or 'g'"}>>> = h2
h_order<<<{"description": "Polynomial order of the switching function h(eta)"}>>> = HIGH
eta<<<{"description": "Order parameter variable"}>>> = eta2
[]
# g(eta)
[g1]
type = BarrierFunctionMaterial<<<{"description": "Helper material to provide $g(\\eta)$ and its derivative in a polynomial.\nSIMPLE: $\\eta^2(1-\\eta)^2$\nLOW: $\\eta(1-\\eta)$\nHIGH: $\\eta^2(1-\\eta^2)^2$", "href": "../../../source/materials/BarrierFunctionMaterial.html"}>>>
function_name<<<{"description": "actual name for f(eta), i.e. 'h' or 'g'"}>>> = g1
g_order<<<{"description": "Polynomial order of the barrier function g(eta)"}>>> = SIMPLE
eta<<<{"description": "Order parameter variable"}>>> = eta1
[]
[g2]
type = BarrierFunctionMaterial<<<{"description": "Helper material to provide $g(\\eta)$ and its derivative in a polynomial.\nSIMPLE: $\\eta^2(1-\\eta)^2$\nLOW: $\\eta(1-\\eta)$\nHIGH: $\\eta^2(1-\\eta^2)^2$", "href": "../../../source/materials/BarrierFunctionMaterial.html"}>>>
function_name<<<{"description": "actual name for f(eta), i.e. 'h' or 'g'"}>>> = g2
g_order<<<{"description": "Polynomial order of the barrier function g(eta)"}>>> = SIMPLE
eta<<<{"description": "Order parameter variable"}>>> = eta2
[]
# constant properties
[constants]
type = GenericConstantMaterial<<<{"description": "Declares material properties based on names and values prescribed by input parameters.", "href": "../../../source/materials/GenericConstantMaterial.html"}>>>
prop_names<<<{"description": "The names of the properties this material will have"}>>> = 'D L kappa'
prop_values<<<{"description": "The values associated with the named properties"}>>> = '10 1 0.1 '
[]
# Coefficients for diffusion equation
[Dh1]
type = DerivativeParsedMaterial<<<{"description": "Parsed Function Material with automatic derivatives.", "href": "../../../source/materials/DerivativeParsedMaterial.html"}>>>
material_property_names<<<{"description": "Vector of material properties used in the parsed function"}>>> = 'D h1(eta1)'
expression<<<{"description": "Parsed function (see FParser) expression for the parsed material"}>>> = D*h1
property_name<<<{"description": "Name of the parsed material property"}>>> = Dh1
coupled_variables<<<{"description": "Vector of variables used in the parsed function"}>>> = eta1
derivative_order<<<{"description": "Maximum order of derivatives taken"}>>> = 1
[]
[Dh2a]
type = DerivativeParsedMaterial<<<{"description": "Parsed Function Material with automatic derivatives.", "href": "../../../source/materials/DerivativeParsedMaterial.html"}>>>
material_property_names<<<{"description": "Vector of material properties used in the parsed function"}>>> = 'D h2(eta2)'
expression<<<{"description": "Parsed function (see FParser) expression for the parsed material"}>>> = D*h2*10/30
property_name<<<{"description": "Name of the parsed material property"}>>> = Dh2a
coupled_variables<<<{"description": "Vector of variables used in the parsed function"}>>> = eta2
derivative_order<<<{"description": "Maximum order of derivatives taken"}>>> = 1
[]
[Dh2b]
type = DerivativeParsedMaterial<<<{"description": "Parsed Function Material with automatic derivatives.", "href": "../../../source/materials/DerivativeParsedMaterial.html"}>>>
material_property_names<<<{"description": "Vector of material properties used in the parsed function"}>>> = 'D h2(eta2)'
expression<<<{"description": "Parsed function (see FParser) expression for the parsed material"}>>> = D*h2*4/30
property_name<<<{"description": "Name of the parsed material property"}>>> = Dh2b
coupled_variables<<<{"description": "Vector of variables used in the parsed function"}>>> = eta2
derivative_order<<<{"description": "Maximum order of derivatives taken"}>>> = 1
[]
[Dh2c]
type = DerivativeParsedMaterial<<<{"description": "Parsed Function Material with automatic derivatives.", "href": "../../../source/materials/DerivativeParsedMaterial.html"}>>>
material_property_names<<<{"description": "Vector of material properties used in the parsed function"}>>> = 'D h2(eta2)'
expression<<<{"description": "Parsed function (see FParser) expression for the parsed material"}>>> = D*h2*16/30
property_name<<<{"description": "Name of the parsed material property"}>>> = Dh2c
coupled_variables<<<{"description": "Vector of variables used in the parsed function"}>>> = eta2
derivative_order<<<{"description": "Maximum order of derivatives taken"}>>> = 1
[]
[]
[Kernels<<<{"href": "../../../syntax/Kernels/index.html"}>>>]
#Kernels for diffusion equation
[diff_time]
type = TimeDerivative<<<{"description": "The time derivative operator with the weak form of $(\\psi_i, \\frac{\\partial u_h}{\\partial t})$.", "href": "../../../source/kernels/TimeDerivative.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = cCr
[]
[diff_c1]
type = MatDiffusion<<<{"description": "Diffusion equation Kernel that takes an isotropic Diffusivity from a material property", "href": "../../../source/kernels/MatDiffusion.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = cCr
diffusivity<<<{"description": "The diffusivity value or material property"}>>> = Dh1
v<<<{"description": "Coupled concentration variable for kernel to operate on; if this is not specified, the kernel's nonlinear variable will be used as usual"}>>> = BCC_CR
args<<<{"description": "Optional vector of arguments for the diffusivity. If provided and diffusivity is a derivative parsed material, Jacobian contributions from the diffusivity will be automatically computed"}>>> = eta1
[]
[diff_c2a]
type = MatDiffusion<<<{"description": "Diffusion equation Kernel that takes an isotropic Diffusivity from a material property", "href": "../../../source/kernels/MatDiffusion.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = cCr
diffusivity<<<{"description": "The diffusivity value or material property"}>>> = Dh2a
v<<<{"description": "Coupled concentration variable for kernel to operate on; if this is not specified, the kernel's nonlinear variable will be used as usual"}>>> = SIGMA_0CR
args<<<{"description": "Optional vector of arguments for the diffusivity. If provided and diffusivity is a derivative parsed material, Jacobian contributions from the diffusivity will be automatically computed"}>>> = eta2
[]
[diff_c2b]
type = MatDiffusion<<<{"description": "Diffusion equation Kernel that takes an isotropic Diffusivity from a material property", "href": "../../../source/kernels/MatDiffusion.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = cCr
diffusivity<<<{"description": "The diffusivity value or material property"}>>> = Dh2b
v<<<{"description": "Coupled concentration variable for kernel to operate on; if this is not specified, the kernel's nonlinear variable will be used as usual"}>>> = SIGMA_1CR
args<<<{"description": "Optional vector of arguments for the diffusivity. If provided and diffusivity is a derivative parsed material, Jacobian contributions from the diffusivity will be automatically computed"}>>> = eta2
[]
[diff_c2c]
type = MatDiffusion<<<{"description": "Diffusion equation Kernel that takes an isotropic Diffusivity from a material property", "href": "../../../source/kernels/MatDiffusion.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = cCr
diffusivity<<<{"description": "The diffusivity value or material property"}>>> = Dh2c
v<<<{"description": "Coupled concentration variable for kernel to operate on; if this is not specified, the kernel's nonlinear variable will be used as usual"}>>> = SIGMA_2CR
args<<<{"description": "Optional vector of arguments for the diffusivity. If provided and diffusivity is a derivative parsed material, Jacobian contributions from the diffusivity will be automatically computed"}>>> = eta2
[]
# enforce pointwise equality of chemical potentials
[chempot1a2a]
# The BCC phase has only one sublattice
# we tie it to the first sublattice with site fraction 10/(10+4+16) in the sigma phase
type = KKSPhaseChemicalPotential<<<{"description": "KKS model kernel to enforce the pointwise equality of phase chemical potentials $dF_a/dc_a = dF_b/dc_b$. The non-linear variable of this kernel is $c_a$.", "href": "../../../source/kernels/KKSPhaseChemicalPotential.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = BCC_CR
cb<<<{"description": "Phase b concentration"}>>> = SIGMA_0CR
kb<<<{"description": "Site fraction for the cb variable (specify this if ca is a sublattice concentration, and make sure it is a true site fraction eg. 0.6666666) "}>>> = '${fparse 10/30}'
fa_name<<<{"description": "Base name of the free energy function Fa (f_name in the corresponding derivative function material)"}>>> = F_BCC_A2
fb_name<<<{"description": "Base name of the free energy function Fb (f_name in the corresponding derivative function material)"}>>> = F_SIGMA
args_b<<<{"description": "Vector of further parameters to Fb (optional, to add in second cross derivatives of Fb)"}>>> = 'SIGMA_1CR SIGMA_2CR'
[]
[chempot2a2b]
# This kernel ties the first two sublattices in the sigma phase together
type = SLKKSChemicalPotential<<<{"description": "SLKKS model kernel to enforce the pointwise equality of sublattice chemical potentials in the same phase.", "href": "../../../source/kernels/SLKKSChemicalPotential.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = SIGMA_0CR
a<<<{"description": "sublattice site fraction for the kernel variable"}>>> = 10
cs<<<{"description": "other sublattice concentration in the same phase"}>>> = SIGMA_1CR
as<<<{"description": "other sublattice site fraction in the same phase"}>>> = 4
F<<<{"description": "Base name of the free energy function"}>>> = F_SIGMA
coupled_variables<<<{"description": "Vector of variable arguments to the free energy function"}>>> = 'SIGMA_2CR'
[]
[chempot2b2c]
# This kernel ties the remaining two sublattices in the sigma phase together
type = SLKKSChemicalPotential<<<{"description": "SLKKS model kernel to enforce the pointwise equality of sublattice chemical potentials in the same phase.", "href": "../../../source/kernels/SLKKSChemicalPotential.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = SIGMA_1CR
a<<<{"description": "sublattice site fraction for the kernel variable"}>>> = 4
cs<<<{"description": "other sublattice concentration in the same phase"}>>> = SIGMA_2CR
as<<<{"description": "other sublattice site fraction in the same phase"}>>> = 16
F<<<{"description": "Base name of the free energy function"}>>> = F_SIGMA
coupled_variables<<<{"description": "Vector of variable arguments to the free energy function"}>>> = 'SIGMA_0CR'
[]
[phaseconcentration]
# This kernel ties the sum of the sublattice concentrations to the global concentration cCr
type = SLKKSMultiPhaseConcentration<<<{"description": "SLKKS multi-phase model kernel to enforce $c_i = \\sum_j h_j\\sum_k a_{jk} c_{ijk}$. The non-linear variable of this kernel is a phase's sublattice concentration", "href": "../../../source/kernels/SLKKSMultiPhaseConcentration.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = SIGMA_2CR
c<<<{"description": "Physical concentration"}>>> = cCr
ns<<<{"description": "Number of sublattices in each phase"}>>> = '1 3'
as<<<{"description": "Sublattice site fractions n all phases"}>>> = '1 10 4 16'
cs<<<{"description": "Array of sublattice concentrations phase for all phases"}>>> = 'BCC_CR SIGMA_0CR SIGMA_1CR SIGMA_2CR'
h_names<<<{"description": "Switching Function Materials for all phases"}>>> = 'h1 h2'
eta<<<{"description": "Order parameters for all phases"}>>> = 'eta1 eta2'
[]
# Kernels for Allen-Cahn equation for eta1
[deta1dt]
type = TimeDerivative<<<{"description": "The time derivative operator with the weak form of $(\\psi_i, \\frac{\\partial u_h}{\\partial t})$.", "href": "../../../source/kernels/TimeDerivative.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = eta1
[]
[ACBulkF1]
type = KKSMultiACBulkF<<<{"description": "KKS model kernel (part 1 of 2) for the Bulk Allen-Cahn. This includes all terms NOT dependent on chemical potential.", "href": "../../../source/kernels/KKSMultiACBulkF.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = eta1
Fj_names<<<{"description": "List of free energies for each phase. Place in same order as hj_names!"}>>> = 'F_BCC_A2 F_SIGMA'
hj_names<<<{"description": "Switching Function Materials that provide h. Place in same order as Fj_names!"}>>> = 'h1 h2'
gi_name<<<{"description": "Base name for the double well function g_i(eta_i)"}>>> = g1
eta_i<<<{"description": "Order parameter that derivatives are taken with respect to"}>>> = eta1
wi<<<{"description": "Double well height parameter"}>>> = 0.1
coupled_variables<<<{"description": "Vector of nonlinear variable arguments this object depends on"}>>> = 'BCC_CR SIGMA_0CR SIGMA_1CR SIGMA_2CR eta2'
[]
[ACBulkC1]
type = SLKKSMultiACBulkC<<<{"description": "Multi-phase SLKKS model kernel for the bulk Allen-Cahn. This includes all terms dependent on chemical potential.", "href": "../../../source/kernels/SLKKSMultiACBulkC.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = eta1
F<<<{"description": "Phase free energy function that is a function of 'c'"}>>> = F_BCC_A2
c<<<{"description": "Concentration variable F depends on"}>>> = BCC_CR
ns<<<{"description": "Number of sublattices in each phase"}>>> = '1 3'
as<<<{"description": "Sublattice site fractions n all phases"}>>> = '1 10 4 16'
cs<<<{"description": "Array of sublattice concentrations phase for all phases"}>>> = 'BCC_CR SIGMA_0CR SIGMA_1CR SIGMA_2CR'
h_names<<<{"description": "Switching Function Materials for all phases"}>>> = 'h1 h2'
eta<<<{"description": "Order parameters for all phases"}>>> = 'eta1 eta2'
[]
[ACInterface1]
type = ACInterface<<<{"description": "Gradient energy Allen-Cahn Kernel", "href": "../../../source/kernels/ACInterface.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = eta1
kappa_name<<<{"description": "The kappa used with the kernel"}>>> = kappa
[]
[lagrange1]
type = SwitchingFunctionConstraintEta<<<{"description": "Lagrange multiplier kernel to constrain the sum of all switching functions in a multiphase system. This kernel acts on a non-conserved order parameter eta_i.", "href": "../../../source/kernels/SwitchingFunctionConstraintEta.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = eta1
h_name<<<{"description": "Switching Function Materials that provides h(eta_i)"}>>> = h1
lambda<<<{"description": "Lagrange multiplier"}>>> = lambda
coupled_variables<<<{"description": "Vector of further variable arguments to the switching function"}>>> = 'eta2'
[]
# Kernels for Allen-Cahn equation for eta1
[deta2dt]
type = TimeDerivative<<<{"description": "The time derivative operator with the weak form of $(\\psi_i, \\frac{\\partial u_h}{\\partial t})$.", "href": "../../../source/kernels/TimeDerivative.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = eta2
[]
[ACBulkF2]
type = KKSMultiACBulkF<<<{"description": "KKS model kernel (part 1 of 2) for the Bulk Allen-Cahn. This includes all terms NOT dependent on chemical potential.", "href": "../../../source/kernels/KKSMultiACBulkF.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = eta2
Fj_names<<<{"description": "List of free energies for each phase. Place in same order as hj_names!"}>>> = 'F_BCC_A2 F_SIGMA'
hj_names<<<{"description": "Switching Function Materials that provide h. Place in same order as Fj_names!"}>>> = 'h1 h2'
gi_name<<<{"description": "Base name for the double well function g_i(eta_i)"}>>> = g2
eta_i<<<{"description": "Order parameter that derivatives are taken with respect to"}>>> = eta2
wi<<<{"description": "Double well height parameter"}>>> = 0.1
coupled_variables<<<{"description": "Vector of nonlinear variable arguments this object depends on"}>>> = 'BCC_CR SIGMA_0CR SIGMA_1CR SIGMA_2CR eta1'
[]
[ACBulkC2]
type = SLKKSMultiACBulkC<<<{"description": "Multi-phase SLKKS model kernel for the bulk Allen-Cahn. This includes all terms dependent on chemical potential.", "href": "../../../source/kernels/SLKKSMultiACBulkC.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = eta2
F<<<{"description": "Phase free energy function that is a function of 'c'"}>>> = F_BCC_A2
c<<<{"description": "Concentration variable F depends on"}>>> = BCC_CR
ns<<<{"description": "Number of sublattices in each phase"}>>> = '1 3'
as<<<{"description": "Sublattice site fractions n all phases"}>>> = '1 10 4 16'
cs<<<{"description": "Array of sublattice concentrations phase for all phases"}>>> = 'BCC_CR SIGMA_0CR SIGMA_1CR SIGMA_2CR'
h_names<<<{"description": "Switching Function Materials for all phases"}>>> = 'h1 h2'
eta<<<{"description": "Order parameters for all phases"}>>> = 'eta1 eta2'
[]
[ACInterface2]
type = ACInterface<<<{"description": "Gradient energy Allen-Cahn Kernel", "href": "../../../source/kernels/ACInterface.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = eta2
kappa_name<<<{"description": "The kappa used with the kernel"}>>> = kappa
[]
[lagrange2]
type = SwitchingFunctionConstraintEta<<<{"description": "Lagrange multiplier kernel to constrain the sum of all switching functions in a multiphase system. This kernel acts on a non-conserved order parameter eta_i.", "href": "../../../source/kernels/SwitchingFunctionConstraintEta.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = eta2
h_name<<<{"description": "Switching Function Materials that provides h(eta_i)"}>>> = h2
lambda<<<{"description": "Lagrange multiplier"}>>> = lambda
coupled_variables<<<{"description": "Vector of further variable arguments to the switching function"}>>> = 'eta1'
[]
# Lagrange-multiplier constraint kernel for lambda
[lagrange]
type = SwitchingFunctionConstraintLagrange<<<{"description": "Lagrange multiplier kernel to constrain the sum of all switching functions in a multiphase system. This kernel acts on the Lagrange multiplier variable.", "href": "../../../source/kernels/SwitchingFunctionConstraintLagrange.html"}>>>
variable<<<{"description": "The name of the variable that this residual object operates on"}>>> = lambda
h_names<<<{"description": "Switching function materials"}>>> = 'h1 h2'
etas<<<{"description": "eta order parameters"}>>> = 'eta1 eta2'
epsilon<<<{"description": "Shift factor to avoid a zero pivot"}>>> = 1e-6
[]
[]
[AuxKernels<<<{"href": "../../../syntax/AuxKernels/index.html"}>>>]
[GlobalFreeEnergy]
type = KKSMultiFreeEnergy<<<{"description": "Total free energy in multi-phase KKS system, including chemical, barrier and gradient terms", "href": "../../../source/auxkernels/KKSMultiFreeEnergy.html"}>>>
variable<<<{"description": "The name of the variable that this object applies to"}>>> = Fglobal
Fj_names<<<{"description": "List of free energies for each phase. Place in same order as hj_names and gj_names!"}>>> = 'F_BCC_A2 F_SIGMA'
hj_names<<<{"description": "Switching Function Materials that provide h. Place in same order as Fj_names and gj_names!"}>>> = 'h1 h2'
gj_names<<<{"description": "Barrier Function Materials that provide g. Place in same order as Fj_names and hj_names!"}>>> = 'g1 g2'
interfacial_vars<<<{"description": "Variable names that contribute to interfacial energy"}>>> = 'eta1 eta2'
kappa_names<<<{"description": "Vector of kappa names corresponding to each variable name in interfacial_vars in the same order."}>>> = 'kappa kappa'
w<<<{"description": "Double well height parameter"}>>> = 0.1
[]
[]
[Executioner<<<{"href": "../../../syntax/Executioner/index.html"}>>>]
type = Transient
solve_type = 'NEWTON'
line_search = none
petsc_options_iname = '-pc_type -sub_pc_type -sub_pc_factor_shift_type -ksp_gmres_restart'
petsc_options_value = 'asm lu nonzero 30'
l_max_its = 100
nl_max_its = 20
nl_abs_tol = 1e-10
end_time = 10000
[TimeStepper<<<{"href": "../../../syntax/Executioner/TimeStepper/index.html"}>>>]
type = IterationAdaptiveDT
optimal_iterations = 12
iteration_window = 2
growth_factor = 1.5
cutback_factor = 0.7
dt = 0.1
[]
[]
[VectorPostprocessors<<<{"href": "../../../syntax/VectorPostprocessors/index.html"}>>>]
[var]
type = LineValueSampler<<<{"description": "Samples variable(s) along a specified line", "href": "../../../source/vectorpostprocessors/LineValueSampler.html"}>>>
start_point<<<{"description": "The beginning of the line"}>>> = '-25 0 0'
end_point<<<{"description": "The ending of the line"}>>> = '25 0 0'
variable<<<{"description": "The names of the variables that this VectorPostprocessor operates on"}>>> = 'cCr eta1 eta2 SIGMA_0CR SIGMA_1CR SIGMA_2CR'
num_points<<<{"description": "The number of points to sample along the line"}>>> = 151
sort_by<<<{"description": "What to sort the samples by"}>>> = id
execute_on<<<{"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."}>>> = 'initial timestep_end'
[]
[mat]
type = LineMaterialRealSampler<<<{"description": "Samples real-valued material properties for all quadrature points in all elements that are intersected by a specified line", "href": "../../../source/vectorpostprocessors/LineMaterialRealSampler.html"}>>>
start<<<{"description": "The beginning of the line"}>>> = '-25 0 0'
end<<<{"description": "The end of the line"}>>> = '25 0 0'
property<<<{"description": "Name of the material property to be output along a line"}>>> = 'F_BCC_A2 F_SIGMA'
sort_by<<<{"description": "What to sort the samples by"}>>> = id
execute_on<<<{"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."}>>> = 'initial timestep_end'
[]
[]
[Postprocessors<<<{"href": "../../../syntax/Postprocessors/index.html"}>>>]
[F]
type = ElementIntegralVariablePostprocessor<<<{"description": "Computes a volume integral of the specified variable", "href": "../../../source/postprocessors/ElementIntegralVariablePostprocessor.html"}>>>
variable<<<{"description": "The name of the variable that this object operates on"}>>> = Fglobal
execute_on<<<{"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."}>>> = 'initial timestep_end'
[]
[cmin]
type = NodalExtremeValue<<<{"description": "Finds either the min or max elemental value of a variable over the domain.", "href": "../../../source/postprocessors/NodalExtremeValue.html"}>>>
value_type<<<{"description": "Type of extreme value to return. 'max' returns the maximum value. 'min' returns the minimum value. 'max_abs' returns the maximum of the absolute value."}>>> = min
variable<<<{"description": "The name of the variable that this postprocessor operates on"}>>> = cCr
execute_on<<<{"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."}>>> = 'initial timestep_end'
[]
[cmax]
type = NodalExtremeValue<<<{"description": "Finds either the min or max elemental value of a variable over the domain.", "href": "../../../source/postprocessors/NodalExtremeValue.html"}>>>
value_type<<<{"description": "Type of extreme value to return. 'max' returns the maximum value. 'min' returns the minimum value. 'max_abs' returns the maximum of the absolute value."}>>> = max
variable<<<{"description": "The name of the variable that this postprocessor operates on"}>>> = cCr
execute_on<<<{"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."}>>> = 'initial timestep_end'
[]
[ctotal]
type = ElementIntegralVariablePostprocessor<<<{"description": "Computes a volume integral of the specified variable", "href": "../../../source/postprocessors/ElementIntegralVariablePostprocessor.html"}>>>
variable<<<{"description": "The name of the variable that this object operates on"}>>> = cCr
execute_on<<<{"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."}>>> = 'initial timestep_end'
[]
[]
[Outputs<<<{"href": "../../../syntax/Outputs/index.html"}>>>]
exodus<<<{"description": "Output the results using the default settings for Exodus output."}>>> = true
print_linear_residuals<<<{"description": "Enable printing of linear residuals to the screen (Console)"}>>> = false
csv<<<{"description": "Output the scalar variable and postprocessors to a *.csv file using the default CSV output."}>>> = true
perf_graph<<<{"description": "Enable printing of the performance graph to the screen (Console)"}>>> = true
[]
(moose/modules/phase_field/examples/slkks/CrFe.i)References
- Aurélie Jacob, Erwin Povoden-Karadeniz, and Ernst Kozeschnik.
Revised thermodynamic description of the fe-cr system based on an improved sublattice model of the σ phase.
Calphad, 60:16–28, 2018.[BibTeX]
- D. Schwen, C. Jiang, and L.K. Aagesen.
A sublattice phase-field model for direct calphad database coupling.
Computational Materials Science, 195:110466, 2021.
URL: https://www.sciencedirect.com/science/article/pii/S0927025621001919, doi:https://doi.org/10.1016/j.commatsci.2021.110466.[BibTeX]