kawin.thermo (Thermodynamics)
kawin.kawin.thermo.Thermodynamics
ExtraFreeEnergyType(v.IndependentPotential)
class ExtraFreeEnergyType(self):
def ExtraFreeEnergyType.reduce(self):
ExtraGibbsModel(Model)
Subclass of pycalphad Model with extra variable GE
GE represents any extra contribution to the Gibbs free energy
such as the Gibbs-Thomson contribution
GeneralThermodynamics
class GeneralThermodynamics(self, database, elements, phases, drivingForceMethod = ‘tangent’, parameters = None):
Class for defining driving force and essential functions for
binary and multicomponent systems using pycalphad for equilibrium
calculations
Parameters
----------
database : Database or str
pycalphad Database or file name for database
elements : list
Elements to consider
Note: reference element must be the first index in the list
phases : list
Phases involved
Note: matrix phase must be first index in the list
drivingForceMethod : str (optional)
Method used to calculate driving force
Options are 'tangent' (default), 'approximate', 'sampling' and 'curvature' (not recommended)
parameters : list [str] or dict {str : float} or None
List of parameters to keep symbolic in the thermodynamic or mobility models
If None, then parameters are fixed
def GeneralThermodynamics.clearCache(self):
Removes any cached data
This is intended for surrogate training, where the cached data
will be removed
def GeneralThermodynamics._buildThermoModels(self):
Builds thermodynamic models for each phase
This assumes that the first phase is the parent phase and the rest of the phases are precipitate phases
For usage in a diffusion model, this won't affect anything
For each precipitate phase, it checks whether the phase has an order/disorder contribution
If so and if the disordered contribution is the parent phase (ex. gamma and gamma prime in Ni alloys),
then we add an new phase to the database that is only the disordered contribution as DIS_<phase name>
def GeneralThermodynamics._buildMobilityModels(self):
Builds mobility models for phases that have model parameters
def GeneralThermodynamics.updateParameters(self, parameters):
Update parameter dictionary with new values
Parameters
----------
parameters : dict {str : float}
Dictionary of parameters
NOTE: this does not have to be the full list and can also have other parameters in it
Only the parameters that are stored upon initialization will be changed
def GeneralThermodynamics._forceDisorder(self, phase):
For phases using an order/disorder model, pycalphad will neglect the disordered phase unless
it is the only phase set active, so the order and disordered portion of the phase will use the same model
For the Gibbs-Thomson effect to be applied, the ordered and disordered parts of the model will need to be kept separate
As a fix, a new phase is added to the database that uses only the disordered part of the model
We keep the original name of the disordered phase, but the disorder contribution will be the DIS_<phase name> phase
def GeneralThermodynamics.setDrivingForceMethod(self, drivingForceMethod):
Sets method for calculating driving force
Parameters
----------
drivingForceMethod - str
Options are ['approximate', 'sampling', 'curvature', 'tangent']
def GeneralThermodynamics.setDFSamplingDensity(self, density):
Sets sampling density for sampling method in driving
force calculations
Default upon initialization is 2000
Parameters
----------
density : int
Number of samples to take per degree of freedom in the phase
def GeneralThermodynamics.setEQSamplingDensity(self, density):
Sets sampling density for equilibrium calculations
Default upon initialization is 500
Parameters
----------
density : int
Number of samples to take per degree of freedom in the phase
def GeneralThermodynamics._generateTdependentFunction(self, func):
def GeneralThermodynamics.setMobility(self, mobility, phase, element = None):
Allows user to define mobility functions
Functions will assume that mobility is only a function of temperature
mobility : dict
Dictionary of functions for each element (including reference)
Each function takes in (v.GE, v.N, v.P, v.T, site fractions) and returns mobility
Optional - only required for multicomponent systems where
mobility terms are not defined in the TDB database
def GeneralThermodynamics.setDiffusivity(self, diffusivity, phase, element = None):
Allows user to define diffusivity functions
diffusivity : dict
Dictionary of functions for each element (including reference)
Each function takes in (v.GE, v.N, v.P, v.T, site fractions) and returns diffusivity
Optional - only required for multicomponent systems where
diffusivity terms are not defined in the TDB database
and if mobility terms are not defined
def GeneralThermodynamics.setMobilityCorrection(self, element, factor):
Factor to multiply mobility by for each element
Parameters
----------
element : str
Element to set factor for
If 'all', factor will be set to all elements
factor : float
Scaling factor
def GeneralThermodynamics._getConditions(self, x, T, gExtra = 0):
Creates dictionary of conditions from composition, temperature and gExtra
Parameters
----------
x : list
Composition (excluding reference element)
T : float
Temperature
gExtra : float (optional)
Gibbs free energy to add to phase. Defaults to 0
def GeneralThermodynamics._setupSubModels(self, precPhase = None):
def GeneralThermodynamics.getEq(self, x, T, gExtra = 0, precPhase = None):
Calculates equilibrium at specified x, T, gExtra
This is separated from the interfacial composition function so that this can be used for getting curvature for interfacial composition from mobility
Parameters
----------
x : float or array
Composition
Needs to be array for multicomponent systems
T : float
Temperature
gExtra : float
Gibbs-Thomson contribution (if applicable)
precPhase : str, int, list or None
Precipitate phase (default is first precipitate)
Options:
None - first precipitate phase in phase list
str - specific precipitate phase by name
list - all phases by name in list
-1 - no precipitate phase
Returns
-------
Workspace object
def GeneralThermodynamics.getLocalEq(self, x, T, gExtra = 0, precPhase = None, composition_sets = None):
Calculates local equilibrium at specified x, T, gExtra
Local equilibrium is defined as all phases have a single composition set
Miscibility gaps are ignored
In the case of a single phase with a miscibility gap, this returns the composition
set at x
Parameters
----------
x : float or array
Composition
Needs to be array for multicomponent systems
T : float
Temperature
gExtra : float
Gibbs-Thomson contribution (if applicable)
precPhase : str, int, list or None
Precipitate phase (default is first precipitate)
Options:
None - first precipitate phase in phase list
str - specific precipitate phase by name
list - all phases by name in list
-1 - no precipitate phase
composition_sets : list[CompositionSet] or None
Composition sets to use in equilibrium, this will override
the phase list and use only the supplied composition sets
If None, then composition sets are generated for each phase
Returns
-------
result - equilibrium convergence and chemical potentials
composition_sets - list of CompositionSet for phases in "equilibrium"
Note - "equilibrium" in terms of the matrix and singled out precipitate phase (or just matrix if precPhase is -1)
def GeneralThermodynamics.getInterdiffusivity(self, x, T, removeCache = True, phase = None):
Gets interdiffusivity at specified x and T
Requires TDB database to have mobility or diffusivity parameters
Parameters
----------
x : float, array or 2D array
Composition
Float or array for binary systems
Array or 2D array for multicomponent systems
T : float or array
Temperature
If array, must be same length as x
For multicomponent systems, must be same length as 0th axis
removeCache : boolean
Removes cached composition set if True, this will require a global
equilibrium calculation if called again
phase : str
Phase to compute diffusivity for (defaults to first or matrix phase)
This only needs to be used for multiphase diffusion simulations
Returns
-------
interdiffusivity - will return array if T is an array
For binary case - float or array of floats
For multicomponent - matrix or array of matrices
def GeneralThermodynamics._interdiffusivitySingle(self, x, T, removeCache = True, phase = None):
Gets interdiffusivity at unique composition and temperature
Parameters
----------
x : float or array
Composition
T : float
Temperature
removeCache : boolean
phase : str
Returns
-------
Interdiffusivity as a matrix (will return float in binary case)
def GeneralThermodynamics.getTracerDiffusivity(self, x, T, removeCache = True, phase = None):
Gets tracer diffusivity for element el at specified x and T
Requires TDB database to have mobility or diffusivity parameters
Parameters
----------
x : float, array or 2D array
Composition
Float or array for binary systems
Array or 2D array for multicomponent systems
T : float or array
Temperature
If array, must be same length as x
For multicomponent systems, must be same length as 0th axis
removeCache : boolean
Removes cached composition set if True, this will require a global
equilibrium calculation if called again
phase : str
Returns
-------
tracer diffusivity - will return array if T is an array
def GeneralThermodynamics._tracerDiffusivitySingle(self, x, T, removeCache = True, phase = None):
Gets tracer diffusivity at unique composition and temperature
Parameters
----------
x : float or array
Composition
T : float
Temperature
el : str
Element to calculate diffusivity
Returns
-------
Tracer diffusivity as a float
def GeneralThermodynamics.getDrivingForce(self, x, T, precPhase = None, removeCache = False, local_phase_sampling_conditions = None):
Gets driving force using method defined upon initialization
Parameters
----------
x : float, array or 2D array
Composition of minor element in bulk matrix phase
For binary system, use an array for multiple compositions
For multicomponent systems, use a 2D array for multiple compositions
Where 0th axis is for indices of each composition
T : float or array
Temperature in K
Must be same length as x if x is array or 2D array
precPhase : str (optional)
Precipitate phase to consider (default is first precipitate phase in list)
removeCache : bool (optional)
If True, this will not cache any equilibrium
This is used for training since training points may not be near each other
Returns
-------
(driving force, precipitate composition)
Driving force is positive if precipitate can form
Precipitate composition will be None if driving force is negative
def GeneralThermodynamics._resetDrivingForceCache(self, phase, removeCache):
def GeneralThermodynamics._getDrivingForceSampling(self, x, T, precPhase, removeCache = False, local_phase_sampling_conditions = None):
Gets driving force for nucleation by sampling
Steps
1. Compute local equilibrium at x and T of only the matrix phase
2. Sample precipitate phase
If ordered contribution to matrix phase, then sample ordering contribution
and remove points on the matrix free energy surface
3. Compute energy difference between precipitate samples and chemical potential hyperplane
4. Find sample that maximizes energy difference and return sample composition and driving force
Parameters
----------
x : float or array
Composition of minor element in bulk matrix phase
Use float for binary systems
Use array for multicomponent systems
T : float
Temperature in K
precPhase : str (optional)
Precipitate phase to consider (default is first precipitate phase in list)
removeCache : bool (optional)
If True, this will not cache any equilibrium
This is used for training since training points may not be near each other
Returns
-------
(driving force, precipitate composition)
Driving force is positive if precipitate can form
Precipitate composition will be None if driving force is negative
def GeneralThermodynamics._getDrivingForceApprox(self, x, T, precPhase, removeCache = False, local_phase_sampling_conditions = None):
Approximate method of driving force calculation
Assumes equilibrium composition of precipitate phase
Sampling method is used if driving force is negative
Steps:
1. Compute equilibrium and get composition sets for matrix and precipitate phase
2. Check for 2 phases and that one phase is the matrix and other phase is precipitate
If not, then resort to sampling method
3. Compute equilibrium at matrix composition and get chemical potential hyperplane
4. Driving force is the difference between the free energy of the precipitate (from step 1)
and the free energy on the chemical potential hyperplane (from step 3) at the precipitate composition
Parameters
----------
x : float or array
Composition of minor element in bulk matrix phase
Use float for binary systems
Use array for multicomponent systems
T : float
Temperature in K
precPhase : str (optional)
Precipitate phase to consider (default is first precipitate phase in list)
removeCache : bool (optional)
If True, this will not cache any equilibrium
This is used for training since training points may not be near each other
Returns
-------
(driving force, precipitate composition)
Driving force is positive if precipitate can form
Precipitate composition will be None if driving force is negative
def GeneralThermodynamics._getDrivingForceCurvature(self, x, T, precPhase, removeCache = False, local_phase_sampling_conditions = None):
Gets driving force from curvature of free energy function
Assumes small saturation
Steps:
1. Compute equilibrium and get composition sets for matrix and precipitate phase
2. Check for 2 phases and that one phase is the matrix and other phase is precipitate
If not, then resort to sampling method
3. Get dmu/dx (free energy curvature)
4. Compute (x_infty - x_matrix) * dmu/dx * (x_prec - x_matrix)^T
This does a first (or second?) order approximation of the driving force based off the curvature at x_infty
Sampling method is used if driving force is negative
Parameters
----------
x : float or array
Composition of minor element in bulk matrix phase
Use float for binary systems
Use array for multicomponent systems
T : float
Temperature in K
precPhase : str (optional)
Precipitate phase to consider (default is first precipitate phase in list)
removeCache : bool (optional)
If True, this will not cache any equilibrium
This is used for training since training points may not be near each other
Returns
-------
(driving force, precipitate composition)
Driving force is positive if precipitate can form
Precipitate composition will be None if driving force is negative
def GeneralThermodynamics._getDrivingForceTangent(self, x, T, precPhase, removeCache = False, local_phase_sampling_conditions = None):
Gets driving force from parallel tangent calculation
Steps
1. Compute equilibrium to get composition sets (or used previous cached CS)
2. Compute equilibrium of matrix phase at matrix composition
3. Remove composition and extra free energy from conditions
4. Add chemical potential for each component to conditions
5. Compute equilibrium of precipitate phase with new conditions
The calculated v.GE is the driving force
This will work for positive and negative driving forces
Parameters
----------
x : float or array
Composition of minor element in bulk matrix phase
Use float for binary systems
Use array for multicomponent systems
T : float
Temperature in K
precPhase : str (optional)
Precipitate phase to consider (default is first precipitate phase in list)
removeCache : bool (optional)
If True, this will not cache any equilibrium
This is used for training since training points may not be near each other
Returns
-------
(driving force, precipitate composition)
Driving force is positive if precipitate can form
Precipitate composition will be None if driving force is negative
def GeneralThermodynamics._getCompositionSetsForDF(self, x, T, precPhase):
Wrapper for getting composition set from x and T by either global equilibrium or local from a cached composition set
Parameters
----------
x : float or array
Composition of minor element in bulk matrix phase
Use float for binary systems
Use array for multicomponent systems
T : float
Temperature in K
precPhase : str (optional)
Precipitate phase to consider (default is first precipitate phase in list)
Returns
-------
chemical_potentials
cs_matrix - composition set of matrix phase
cs_precip - composition set of precipitate phase
def GeneralThermodynamics._getCompositionSetsEq(self, x, T, precPhase, cached_composition_sets = {}):
Gets composition set from x and T by global equilibrium
Steps
1. Compute equilibrium at x and T
If equilibrium did not converge or matrix phase is not stable, then return None
2. Get composition sets and add to cache
If precipitate is not stable, the return None
3. Resolve possible issues with miscibility gaps
4. Return values
Parameters
----------
x : float or array
Composition of minor element in bulk matrix phase
Use float for binary systems
Use array for multicomponent systems
T : float
Temperature in K
precPhase : str (optional)
Precipitate phase to consider (default is first precipitate phase in list)
Returns
-------
chemical_potentials
cs_matrix - composition set of matrix phase
cs_precip - composition set of precipitate phase
def GeneralThermodynamics._process_composition_sets(composition_sets):
def GeneralThermodynamics._update_composition_sets(composition_sets):
def GeneralThermodynamics._getPrecCompositionSetSamplingDF(self, x, T, matrix_chem_pot, precPhase, local_phase_sampling_conditions = None):
Gets samples for precipitate phase for use in sampling driving force method and returns driving force and precipitate composition
This is also use in tangent driving force method for when equilibrium is not (yet) cached
Steps
1. Compute local equilibrium at x and T of only the matrix phase
2. Sample precipitate phase
If ordered contribution to matrix phase, then sample ordering contribution
and remove points on the matrix free energy surface
3. Compute energy difference between precipitate samples and chemical potential hyperplane
4. Find sample that maximizes energy difference and return sample composition and driving force
Parameters
----------
x : float or array
Composition of minor element in bulk matrix phase
Use float for binary systems
Use array for multicomponent systems
T : float
Temperature in K
precPhase : str (optional)
Precipitate phase to consider (default is first precipitate phase in list)
Returns
-------
driving force - max free energy difference
precipitate composition set - corresponds to max driving force
kawin.kawin.thermo.BinTherm
BinaryThermodynamics (GeneralThermodynamics)
class BinaryThermodynamics (self, database, elements, phases, drivingForceMethod = ‘tangent’, interfacialCompMethod = ‘equilibrium’, parameters = None):
Class for defining driving force and interfacial composition functions
for a binary system using pyCalphad and thermodynamic databases
Parameters
----------
database : str
File name for database
elements : list
Elements to consider
Note: reference element must be the first index in the list
phases : list
Phases involved
Note: matrix phase must be first index in the list
drivingForceMethod : str (optional)
Method used to calculate driving force
Options are 'tangent' (default), 'approximate', 'sampling', and 'curvature' (not recommended)
interfacialCompMethod: str (optional)
Method used to calculate interfacial composition
Options are 'equilibrium' (default) and 'curvature' (not recommended)
parameters : list [str] or dict {str : float}
List of parameters to keep symbolic in the thermodynamic or mobility models
def BinaryThermodynamics .setInterfacialMethod(self, interfacialCompMethod):
Changes method for caluclating interfacial composition
Parameters
----------
interfacialCompMethod - str
Options are ['equilibrium', 'curvature']
def BinaryThermodynamics .setGuessComposition(self, conditions):
Sets initial composition when calculating equilibrium for interfacial energy
Parameters
----------
conditions : float, tuple or dict
Guess composition(s) to solve equilibrium for
This should encompass the region where a tieline can be found
between the matrix and precipitate phases
Options: float - will set to all precipitate phases
tuple - (min, max dx) will set to all precipitate phases
dictionary {phase name: scalar or tuple}
def BinaryThermodynamics .getInterfacialComposition(self, T, gExtra = 0, precPhase = None):
Gets interfacial composition accounting for Gibbs-Thomson effect
Parameters
----------
T : float or array
Temperature in K
gExtra : float or array (optional)
Extra contributions to the precipitate Gibbs free energy
Gibbs Thomson contribution defined as Vm * (2*gamma/R + g_Elastic)
Defaults to 0
precPhase : str
Precipitate phase to consider (default is first precipitate in list)
Note: for multiple conditions, only gExtra has to be an array
This will calculate compositions for multiple gExtra at the input Temperature
If T is also an array, then T and gExtra must be the same length
where each index will pertain to a single condition
Returns
-------
(parent composition, precipitate composition)
Both will be either float or array based off shape of gExtra
Will return (None, None) if precipitate is unstable
def BinaryThermodynamics ._interfacialCompositionFromEq(self, T, gExtra, precPhase):
Gets interfacial composition by calculating equilibrum with Gibbs-Thomson effect
Parameters
----------
T : float
Temperature in K
gExtra : float (optional)
Extra contributions to the precipitate Gibbs free energy
Gibbs Thomson contribution defined as Vm * (2*gamma/R + g_Elastic)
Defaults to 0
precPhase : str
Precipitate phase to consider (default is first precipitate in list)
Returns
-------
(parent composition, precipitate composition)
Both will be either float or array based off shape of gExtra
Will return (None, None) if precipitate is unstable
def BinaryThermodynamics ._interfacialCompositionFromCurvature(self, T, gExtra, precPhase):
Gets interfacial composition using free energy curvature
G''(x - xM)(xP-xM) = 2*y*V/R
Parameters
----------
T : float
Temperature in K
gExtra : float (optional)
Extra contributions to the precipitate Gibbs free energy
Gibbs Thomson contribution defined as Vm * (2*gamma/R + g_Elastic)
Defaults to 0
precPhase : str
Precipitate phase to consider (default is first precipitate in list)
Returns
-------
(parent composition, precipitate composition)
Both will be either float or array based off shape of gExtra
Will return (None, None) if precipitate is unstable
def BinaryThermodynamics .plotPhases(self, ax, T, gExtra = 0, plotGibbsOffset = False, *args, **kwargs):
Plots sampled points from the parent and precipitate phase
Parameters
----------
ax : Axis
T : float
Temperature in K
gExtra : float (optional)
Extra contributions to the Gibbs free energy of precipitate
Defaults to 0
plotGibbsOffset : bool (optional)
If True and gExtra is not 0, the sampled points of the
precipitate phase will be plotted twice with gExtra and
with no extra Gibbs free energy contributions
Defualts to False
kawin.kawin.thermo.MultiTherm
def _growthRateOutputFromCurvature(x, dG, R, gExtra, curvature: CurvatureOutput):
MulticomponentThermodynamics (GeneralThermodynamics)
class MulticomponentThermodynamics (self, database, elements, phases, drivingForceMethod = ‘tangent’, parameters = None):
Class for defining driving force and (possibly) interfacial composition functions
for a multicomponent system using pyCalphad and thermodynamic databases
Parameters
----------
database : str
File name for database
elements : list
Elements to consider
Note: reference element must be the first index in the list
phases : list
Phases involved
Note: matrix phase must be first index in the list
drivingForceMethod : str (optional)
Method used to calculate driving force
Options are 'tangent' (default), 'approximate', 'sampling' and 'curvature' (not recommended)
parameters : list [str] or dict {str : float}
List of parameters to keep symbolic in the thermodynamic or mobility models
def MulticomponentThermodynamics .clearCache(self):
def MulticomponentThermodynamics .getInterfacialComposition(self, x, T, gExtra = 0, precPhase = None):
Gets interfacial composition by calculating equilibrum with Gibbs-Thomson effect
Parameters
----------
T : float or array
Temperature in K
gExtra : float or array (optional)
Extra contributions to the precipitate Gibbs free energy
Gibbs Thomson contribution defined as Vm * (2*gamma/R + g_Elastic)
Defaults to 0
precPhase : str
Precipitate phase to consider (default is first precipitate in list)
Note: for multiple conditions, only gExtra has to be an array
This will calculate compositions for multiple gExtra at the input Temperature
If T is also an array, then T and gExtra must be the same length
where each index will pertain to a single condition
Returns
-------
(parent composition, precipitate composition)
Both will be either float or array based off shape of gExtra
Will return (None, None) if precipitate is unstable
def MulticomponentThermodynamics ._interfacialComposition(self, x, T, gExtra, precPhase):
Gets interfacial composition, will return None, None if composition is in single phase region
Parameters
----------
T : float
Temperature in K
gExtra : float (optional)
Extra contributions to the precipitate Gibbs free energy
Gibbs Thomson contribution defined as Vm * (2*gamma/R + g_Elastic)
Defaults to 0
precPhase : str
Precipitate phase to consider (default is first precipitate in list)
Returns
-------
(parent composition, precipitate composition)
Both will be either float or array based off shape of gExtra
Will return (None, None) if precipitate is unstable
def MulticomponentThermodynamics ._curvatureFactorFromEq(self, chemical_potentials, cs_matrix, cs_precip, precPhase):
Curvature factor (from Phillipes and Voorhees - 2013)
Steps
1. Check that there is 2 phases in equilibrium, one being the matrix and the other being precipitate
2. Get Dnkj, dmu/dx and inverse mobility term from composition set of matrix phase
3. Get dmu/dx of precipitate phase
4. Get difference in matrix and precipitate phase composition (we use a second order approximation to get precipitate composition as function of R)
5. Compute numerator, denominator, Gba and beta term
Denominator (X_bar^T * invMob * X_bar) is used for growth rate (eq 28)
Numerator (D^-1 * X_bar), denominator is used for matrix interfacial composition (eq 31)
Gba and matrix interfacial composition is used for precipitate interfacial composition (eq 36)
Gba here is (dmu/dx_beta)^-1 * dmu/dx_alpha
Note: these equations have a term X_bar^T * dmu/dx_alpha * X_bar_infty, but this is just the driving force so we don't need to calculate it here
Parameters
----------
chemical_potentials : 1-D float64 array
cs_matrix : CompositionSet
cs_precip : CompositionSet
precPhase : str (optional)
Precipitate phase (defaults to first precipitate in list)
Returns
-------
CurvatureOutput object
{D-1 dCbar / dCbar^T M-1 dCbar} - for calculating interfacial composition of matrix
{1 / dCbar^T M-1 dCbar} - for calculating growth rate
{Gb^-1 Ga} - for calculating precipitate composition
beta - Impingement rate
Ca - interfacial composition of matrix phase
Cb - interfacial composition of precipitate phase
def MulticomponentThermodynamics .curvatureFactor(self, x, T, precPhase = None, removeCache = False, searchDir = None, computeSearchDir = False):
Curvature factor (from Phillipes and Voorhees - 2013) from composition and temperature
This will find equilibrium and composition sets for x and T and call _curvatureFactorFromEq
Parameters
----------
x : array
Composition of solutes
T : float
Temperature
precPhase : str (optional)
Precipitate phase (defaults to first precipitate in list)
searchDir : None or array
If two-phase equilibrium is not present, then move x towards this composition to find two-phase equilibria
removeCache : bool (optional)
If True, this will not cache any equilibrium
This is used for training since training points may not be near each other
Returns
-------
CurvatureCache object
{D-1 dCbar / dCbar^T M-1 dCbar} - for calculating interfacial composition of matrix
{1 / dCbar^T M-1 dCbar} - for calculating growth rate
{Gb^-1 Ga} - for calculating precipitate composition
beta - Impingement rate
Ca - interfacial composition of matrix phase
Cb - interfacial composition of precipitate phase
Will return None if single phase
def MulticomponentThermodynamics ._process_invalid_eq(calc_source):
def MulticomponentThermodynamics ._searchForTwoPhaseEq(self, x, T, precPhase, removeCache = False, searchDir = None, computeSearchDir = False):
Given x and a search direction (which should correspond to the composition of precipitate from driving force calc)
Iteratively search between x and search direction until two phase equilibria is found
If two phase equilibria is not found, then return None
else return chemical_potential, cs_matrix, cs_precip
def MulticomponentThermodynamics .getGrowthAndInterfacialComposition(self, x, T, dG, R, gExtra, precPhase = None, removeCache = False, searchDir = None):
Returns growth rate and interfacial compostion given Gibbs-Thomson contribution
Parameters
----------
x : array
Composition of solutes
T : float
Temperature
dG : float
Driving force at given x and T
R : float or array
Precipitate radius
gExtra : float or array
Gibbs-Thomson contribution (must be same shape as R)
precPhase : str (optional)
Precipitate phase (defaults to first precipitate in list)
searchDir : None or array
If two-phase equilibrium is not present, then move x towards this composition to find a two-phase region
removeCache : bool (optional)
If True, this will not cache any equilibrium
This is used for training since training points may not be near each other
Returns
-------
GrowthRateOutput object
(growth rate, matrix composition, precipitate composition, equilibrium matrix comp, equilibrium precipitate comp)
growth rate will be float or array based off shape of R
matrix and precipitate composition will be array or 2D array based
off shape of R
def MulticomponentThermodynamics .impingementFactor(self, x, T, precPhase = None, removeCache = False, searchDir = None):
Returns impingement factor for nucleation rate calculations
Parameters
----------
x : array
Composition of solutes
T : float
Temperature
precPhase : str (optional)
Precipitate phase (defaults to first precipitate in list)
removeCache : bool (optional)
If True, this will not cache any equilibrium
This is used for training since training points may not be near each other
Returns
-------
beta - impingement factor