# Extra Factors in the KWN Model #

The default options in the KWN model assumes bulk nucleation and a spherical precipitate. However, in real-life systems, these assumptions may not be true. Several options are present to model heterogenous nucleation and non-spherical precipitate shapes.

## Nucleation #

The options for the nucleation site density includes the following: ‘bulk’, ‘dislocations’, ‘grain_boundaries’, ‘grain_edges’, ‘grain_corners’

The nucleation site density ($N_0$) for bulk nucleation is determined by the number of solutes in the bulk lattice. For dislocation, $N_0$ depends on the dislocation density. $N_0$ for grain boundaries, edges and corners depends on the grain size [1]. For grain boundary nucleation, the change in surface energy accounts for both the creation of the precipitate/matrix interface and removal of grain boundary, for which the grain boundary energy must be defined [2]. By default, the grain boundary energy is set to 0.3 $J/m^2$.

While the KWNModel will automatically calculate the nucleation site densities for each site type, these values can be manually set.

from kawin.KWNBase import PrecipitateBase

model = PrecipitateBase(t0 = 0, tf = 100, steps = 2e3, linearTimeSpacing = True)

#Change nucleation site type to grain boundaries
model.setNucleationSite('grain_boundaries')

#Set grain boundary energy for nucleation on grain boundaries/edges/corners
model.setGrainBoundaryEnergy(0.3)

#Change dislocation density and grain size
model.setNucleationDensity(grainSize = 10, aspectRatio = 1, dislocationDensity = 5e12)

#Manually set nucleation site density for each site type
model.bulkN0 = 1e30          #Bulk nucleation site density
model.dislocationN0 = 1e30   #Site density on dislocations
model.GBareaN0 = 1e30        #Site density on grain boundaries
model.GBedgeN0 = 1e30        #Site density on grain edges
model.GBcornerN0 = 1e30      #Site density on grain corners


## Shape factors #

Currently, the KWN model has support for ellipsoidal (prolate/needle and oblate/plate) and cuboidal precipitates. For cuboidal precipitates the cubic factor currently set constant at $\sqrt{2}$ [3,4]. These shapes are defined in the KWN model by their equivalent spherical radius ($R_{sph}$) and an aspect ratio ($\alpha$), where the $\alpha$ can either be constant or as a function of $R_{sph}$.

The aspect ratio is defined as the ratio of the long axis over the short axis. Conversion between the radius along the short axis ($r$) and the equivalent spherical radius ($R_{sph}$) is given by:

Needle: $R_{sph} = \sqrt[3]{\alpha} r$

Plate: $R_{sph} = \sqrt[3]{\alpha^2} r$

Cuboidal: $R_{sph} = \sqrt[3]{\frac{3 \alpha}{4 \pi}} r$

Deviation from a spherical precipitate changes both the thermodynamics (Gibbs-Thomson effect) and kinetics (growth rate). The free energy contribution from the Gibbs-Thomson effect is given by:

$$\Delta G_{TH} = g(\alpha) \frac{2 \gamma V_M^\beta}{R_{sph}}$$

The changes in the growth rate is given by:

$$\frac{dR}{dt} = f(\alpha) \frac{dR_{sph}}{dt}$$

The functions of $g(\alpha)$ and $f(\alpha)$ are taken from Ref. 3 and 4.

#Change precipitate shape
model.setAspectRatioNeedle(ratio = 1.5)
model.setAspectRatioPlate(ratio = 1.5)
model.setAspectRatioCuboidal(ratio = 1.5)

#Remove aspect ratio and set to spherical shape
model.setSpherical()

ar = lambda r: 2.3 * (r/1e-9)**1.1
model.setAspectRatioNeedle(ratio = ar)


## Strain Energy #

Molar volume differences between the matrix and precipitate phase can induce strains, which reduces the driving force for nucleation. For spherical and cuboidal precipitates, the strain energy can be calculated by Khachaturyan’s approximation. For ellipsoidal precipitates, the strain energy can be calculated using Eshelby’s tensor [4,5].

Similar to the Thermodynamics and Surrogate modules, the strain energy is calculated using a module separated from KWNBase. Inserting the strain energy parameters requires creating and setting up a StrainEnergy object, then inserting it into the KWN model for a specified phase.

The StrainEnergy object requires the elastic constants and eigenstrains to be defined. External stresses can also be defined if applicable. The eigenstrains and external stress can be defined as a tensor (3x3), values along the three axes (array of length 3), or a single value to be applied on all 3 axes. The elastic constants can be defined using its 6x6 tensor, the three elastic constants ($c_{11}$, $c_{12}$ and $c_{44}$), or by at least two moduli (e.g. elastic modulus, poission ratio, shear modulus, bulk modulus, etc.)

from kawin.ElasticFactors import StrainEnergy

#Create StrainEnergy object
se = StrainEnergy()

#Set elastic tensor by its elastic modulus and possion ratio
se.setModuli(E = 100e9, nu = -0.3)

#Set eigenstrains
# [[0.01, 0.00, 0.00]
#  [0.00, 0.01, 0.00]
#  [0.00, 0.00, 0.02]]
se.setEigenstrain([0.01, 0.01, 0.02])

#Insert StrainEnergy object into KWN model
model.setStrainEnergy(se)

#Use strain energy to calculate aspect ratio (for plate- and needle-like precipitates)
#This will override the aspect ratio that was defined when setting the precipitate shape
model.setStrainEnergy(se, calculateAspectRatio=True)


## Strain Energy Example (Cu-Ti system) #

In the Cu-Ti system (dilute Ti), the needle-like $Cu_4Ti$ precipitate creates lattice strains in the Cu-matrix. The following parameters are applicable to this system (from K. Wu et al (2018)):

Eigenstrains of the $Cu_4Ti$ precipitate:

$\epsilon_{11} = 0.022$

$\epsilon_{22} = 0.022$

$\epsilon_{33} = 0.003$

Elastic constants for the Cu matrix

$c_{11} = 168.4 \quad GPa$

$c_{12} = 121.4 \quad GPa$

$c_{44} = 75.4 \quad GPa$

We can use these values to determine the strain energy of the $Cu_4Ti$ precipitate for any given aspect ratio. In this example, we’ll vary the aspect ratio from 1 to 2 and calculate the strain energy. The volume of the precipitate will be set constant to the volume of a sphere with a radius of 4 nm.

%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np

#By default, StrainEnergy outputs 0
#This is changed within the KWN model before the model is solved for
#However, we can manually change it. For this example, we need to set it to the calculate for ellipsoidal shapes
se.setEllipsoidal()

#Set elastic tensor by c11, c12 and c44 values
se.setElasticConstants(168.4e9, 121.4e9, 75.4e9)

#Set eigenstrains
se.setEigenstrain([0.022, 0.022, 0.003])

#Aspect ratio
aspect = np.linspace(1, 2, 100)

#Equivalent spherical radius of 4 nm
rSph = 4e-9 / np.cbrt(aspect)
r = np.array([rSph, rSph, aspect*rSph]).T

E = se.strainEnergy(r)

plt.plot(aspect, E)
plt.xlim([1, 2])
plt.ylim([1.0e-17, 1.5e-17])
plt.xlabel('Aspect Ratio')
plt.ylabel('Strain Energy (J)')
plt.show()


## Calculating Aspect Ratio from Strain Energy #

The aspect ratio for plate- and needle-like precipitates can be determined by minimizing the energy contributions from the strain and interfacial energy contributions.

$$\alpha = argmin\left( \frac{4}{3}\pi R_{sph}^{3} \Delta G_{el}(\alpha) + 4 \pi R_{sph}^{2} g(\alpha) \gamma \right)$$

Where $R_{sph}$ is the equivalent spherical radius. The strain energy module has two options for calculating the equilibrium aspect ratio: iterative or searching. The iterative method (StrainEnergy.eqAR_byGR) performs a Golden Section search to find the minimum. The search method (StrainEnergy.eqAR_bySearch) will calculate the net energy contribution for a number of aspect ratios and will return the aspect ratio that gives the minimum. By default, this method is accurate up to 2 significant digits. In addition, due to caching, this method is also faster than the iterative method for large number of calculations.

## Example ($\gamma''$ in IN718) #

The $\gamma''$ precipitate in IN718 are plate shape where the aspect ratio depends on the size of the precipitate. Using the elastic properties of IN718 (shear modulus of 57.1 GPa and Poisson’s ratio of 0.33) and the eigenstrain of the $\gamma''$ precipitate, the relationship between the aspect ratio and precipitate diameter (long axis) can be found.

from kawin.ShapeFactors import ShapeFactor

#Strain energy parameters
se = StrainEnergy()
se.setEigenstrain([6.67e-3, 6.67e-3, 2.86e-2])
se.setModuli(G=57.1e9, nu=0.33)
se.setEllipsoidal()

#Shape factor parameters (only the shape needs to be defined)
sf = ShapeFactor()
sf.setPlateShape()

#Calculate equilibrium aspect ratio
gamma = 0.02375
Rsph = np.linspace(1e-10, 40e-9, 100)
eqAR = se.eqAR_bySearch(Rsph, gamma, sf)

#Convert spherical radius to diameter of the plate
R = 2*Rsph / np.cbrt(eqAR**2)*eqAR

#Plot diameter vs. aspect ratio
plt.plot(R, eqAR)
plt.xlim([0, 40e-9])
plt.ylim([1, 9])
plt.xlabel('Diameter (m)')
plt.ylabel('Aspect Ratio')
plt.show()


## References #

1. Kozeschnik, Ernst et al. Precipitation Modeling, Momentum Press, 2012
2. P. J. Clemm and J. C. Fisher, "The Influence of Grain Boundaries on the Nucleation of Secondary Phases" Acta Metallurgica 3 (1955) p. 70
3. B. Holmedal, E. Osmundsen, Q. Du, "Precipitation of Non-Spherical Particles in Aluminum Alloys Part I: Generalization of the Kampmann-Wagner Numerical Model" Metallurgical and Materials Transactions A 47 (2016) p. 581
4. K. Wu, Q. Chen and P. Mason, "Simulation of Precipitation Kinetics with Non-Spherical Particles" J. Phase Equilib. Diffus. 39 (2018) p. 571
5. C. Weinberger, W. Cai and D. Barnett, ME340B Lecture Notes - Elasticity of Microscopic Structures, Standford University 2005. http://micro.standford.edu/~caiwei/me340b/content/me340b-notes_v01.pdf