elaston.linear_elasticity module

Bases: object

Linear elastic field class based on the 3x3x3x3 elastic tensor :math:C_{ijkl}:

where :math:sigma_{ij} is the ij-component of stress and :math:epsilon_{kl} is the kl-component of strain.

Examples I: Get bulk modulus from the elastic tensor from the voigt average:

>>> medium = LinearElasticity(C_11=211.0, C_12=130.0, C_44=82.0)
>>> medium_voigt = medium.get_voigt_average()
>>> parameters = medium_voigt.get_elastic_moduli()
>>> print(parameters['bulk_modulus'])

Example II: Get strain field around a point defect:

>>> import numpy as np
>>> medium = LinearElasticity(C_11=211.0, C_12=130.0, C_44=82.0)  # Fe
>>> random_positions = np.random.random((10, 3)) - 0.5
>>> dipole_tensor = np.eye(3)
>>> print(medium.get_point_defect_strain(random_positions, dipole_tensor))

Example III: Get stress field around a dislocation:

>>> import numpy as np
>>> medium = LinearElasticity(C_11=211.0, C_12=130.0, C_44=82.0)
>>> random_positions = np.random.random((10, 3))-0.5
>>> # Burgers vector of a screw dislocation in bcc Fe
>>> burgers_vector = np.array([0, 0, 2.86 * np.sqrt(3) / 2])
>>> print(medium.get_dislocation_stress(random_positions, burgers_vector))

Example IV: Estimate the distance between partial dislocations:

>>> medium = LinearElasticity(C_11=110.5, C_12=64.8, C_44=30.9)  # Al
>>> lattice_constant = 4.05
>>> partial_one = np.array([-0.5, 0, np.sqrt(3) / 2]) * lattice_constant
>>> partial_two = np.array([0.5, 0, np.sqrt(3) / 2]) * lattice_constant
>>> distance = 100
>>> stress_one = medium.get_dislocation_stress([0, distance, 0], partial_one)
>>> print('Choose `distance` in the way that the value below corresponds to SFE')
>>> medium.get_dislocation_force(stress_one, [0, 1, 0], partial_two)

get_compliance_tensor(rotate: bool = True, voigt: bool = False) → ndarray[source]

get_dislocation_displacement(positions: ndarray, burgers_vector: ndarray) → ndarray[source]

Displacement field around a dislocation according to anisotropic elasticity theory described by [Eshelby](https://doi.org/10.1016/0001-6160(53)90099-6).

Parameters:

positions ((n,2) or (n,3)-array) – Position around a dislocation. The third axis coincides with the dislocation line.
burgers_vector ((3,)-array) – Burgers vector

Returns:

Displacement field (z-axis coincides with the: dislocation line)

Return type:

((n, 3)-array)

get_dislocation_energy(burgers_vector: ndarray, r_min: float, r_max: float, mesh: int = 100) → float[source]

Energy per unit length along the dislocation line.

Parameters:

burgers_vector ((3,)-array) – Burgers vector
r_min (float) – Minimum distance from the dislocation core
r_max (float) – Maximum distance from the dislocation core
mesh (int) – Number of grid points for the numerical integration along the angle

Returns:

Energy of dislocation per unit length

Return type:

(float)

The energy is defined by the product of the stress and strain (i.e. energy density), which is integrated over the plane vertical to the dislocation line. The energy density \(w\) according to the linear elasticity is given by:

Therefore, the energy per unit length \(U\) is given by:

This implies \(r_min\) cannot be 0 as well as \(r_max\) cannot be infinity. This is the consequence of the fact that the linear elasticity cannot describe the core structure properly, and a real medium is not infinitely large. While \(r_max\) can be defined based on the real dislocation density, the choice of \(r_min\) should be done carefully.

get_dislocation_energy_density(positions: ndarray, burgers_vector: ndarray) → ndarray[source]

Energy density field around a dislocation (product of stress and strain, cf. corresponding methods)

Parameters:

positions ((n,2) or (n,3)-array) – Position around a dislocation. The third axis coincides with the dislocation line.
burgers_vector ((3,)-array) – Burgers vector

Returns:

Energy density field

Return type:

((n,)-array)

static get_dislocation_force(stress: ndarray, glide_plane: ndarray, burgers_vector: ndarray) → ndarray[source]

Force per unit length along the dislocation line. This method is useful for estaimting the distance between partial dislocations. At equilibrium, the force acting on the dislocation line corresponds to the stacking fault energy (SFE).

Parameters:

stress ((n, 3, 3)-array) – External stress field at the dislocation line
glide_plane ((3,)-array) – Glide plane
burgers_vector ((3,)-array) – Burgers vector

Returns:

Force per unit length acting on the dislocation.

Return type:

((3,)-array)

Here is a short list of the stacking fault energies of a few materials:

Ni: 90 mJ/m^2 = 5.62 meV/Å^2
Ag: 25 mJ/m^2 = 1.56 meV/Å^2
Au: 75 mJ/m^2 = 4.68 meV/Å^2
Cu: 70-78 mJ/m^2 = 4.36-4.87 meV/Å^2
Mg: 125 mJ/m^2 = 7.80 meV/Å^2
Al: 160-250 mJ/m^2 = 10.0-15.6 meV/Å^2

Source: https://en.wikipedia.org/wiki/Stacking-fault_energy

get_dislocation_strain(positions: ndarray, burgers_vector: ndarray) → ndarray[source]

Strain field around a dislocation according to anisotropic elasticity theory described by [Eshelby](https://doi.org/10.1016/0001-6160(53)90099-6).

Parameters:

positions ((n,2) or (n,3)-array) – Position around a dislocation. The third axis coincides with the dislocation line.
burgers_vector ((3,)-array) – Burgers vector

Returns:

Strain field (z-axis coincides with the dislocation line)

Return type:

((n, 3, 3)-array)

get_dislocation_stress(positions: ndarray, burgers_vector: ndarray) → ndarray[source]

Stress field around a dislocation according to anisotropic elasticity theory described by [Eshelby](https://doi.org/10.1016/0001-6160(53)90099-6).

Parameters:

positions ((n,2) or (n,3)-array) – Position around a dislocation. The third axis coincides with the dislocation line.
burgers_vector ((3,)-array) – Burgers vector

Returns:

Stress field (z-axis coincides with the dislocation line)

Return type:

((n, 3, 3)-array)

get_elastic_moduli() → dict[str, float][source]

get_elastic_tensor(voigt: bool = False, rotate: bool = True) → ndarray[source]

get_point_defect_displacement(positions: ndarray, dipole_tensor: ndarray, n_mesh: int = 100, optimize: bool = True, check_unique: bool = False) → ndarray[source]

Displacement field around a point defect

Parameters:

positions ((n,3)-array) – Positions in real space or reciprocal space (if fourier=True).
dipole_tensor ((3,3)-array) – Dipole tensor
n_mesh (int) – Number of mesh points in the radial integration in case if anisotropic Green’s function (ignored if isotropic=True or fourier=True)
optimize (bool) – cf. optimize in numpy.einsum
check_unique (bool) – Whether to check the unique positions

Returns:

Displacement field

Return type:

((n,3)-array)

According to the definition of the Green’s function (cf. docstring of get_greens_function):

where \(u_i(r)\) is the displacement field of component \(i\) at position \(r\) and \(f_j(a)\) is the force component \(j\) of the atom at position \(a\). By taking the polynomial development we obtain:

The first term disappears because the sum of the forces is zero. From the second term we define the dipole tensor \(P_{jk} = a_k f_j(a)\). Following the definition above, we can obtain the displacement field, strain field, stress field and energy density field if the dipole tensor and the elastic tensor are known.

The dipole tensor of a point defect is commonly obtained from the following equation:

where \(U\) is the potential energy, \(V\) is the volume and \(\varepsilon\) is the strain field. At equilibrium, the derivative of the potential energy with respect to the strain disappears:

With this in mind, we can calculate the dipole tensor of Ni in Al with the following lines:

>>> from pyiron_atomistics import Project
>>> pr = Project('dipole_tensor')
>>> job = pr.create.job.Lammps('dipole')
>>> n_repeat = 3
>>> job.structure = pr.create.structure.bulk('Al', cubic=True).repeat(n_repeat)
>>> job.structure[0] = 'Ni'
>>> job.calc_minimize()
>>> job.run()
>>> dipole_tensor = -job.structure.get_volume() * job['output/generic/pressures'][-1]

Instead of working with atomistic calculations, the dipole tensor can be calculated by the lambda tensor [1], which is defined as:

where \(c\) is the concentration of the defect, \(V\) is the volume and \(\varepsilon\) is the strain field. Then the dipole tensor is given by:

ref:

[1] Nowick, Arthur S. Anelastic relaxation in crystalline solids. Vol. 1. Elsevier, 2012.

get_point_defect_energy_density(positions: ndarray, dipole_tensor: ndarray, n_mesh: int = 100, optimize: bool = True) → ndarray[source]

Energy density field around a point defect using the Green’s function method

Parameters:

positions ((n,3)-array) – Positions in real space or reciprocal space (if fourier=True).
dipole_tensor ((3,3)-array) – Dipole tensor
n_mesh (int) – Number of mesh points in the radial integration in case if anisotropic Green’s function (ignored if isotropic=True or fourier=True)
optimize (bool) – cf. optimize in numpy.einsum

Returns:

Energy density field

Return type:

((n,)-array)

According to the definition of the Green’s function (cf. docstring of get_greens_function):

where \(u_i(r)\) is the displacement field of component \(i\) at position \(r\) and \(f_j(a)\) is the force component \(j\) of the atom at position \(a\). By taking the polynomial development we obtain:

The first term disappears because the sum of the forces is zero. From the second term we define the dipole tensor \(P_{jk} = a_k f_j(a)\). Following the definition above, we can obtain the displacement field, strain field, stress field and energy density field if the dipole tensor and the elastic tensor are known.

The dipole tensor of a point defect is commonly obtained from the following equation:

where \(U\) is the potential energy, \(V\) is the volume and \(\varepsilon\) is the strain field. At equilibrium, the derivative of the potential energy with respect to the strain disappears:

With this in mind, we can calculate the dipole tensor of Ni in Al with the following lines:

>>> from pyiron_atomistics import Project
>>> pr = Project('dipole_tensor')
>>> job = pr.create.job.Lammps('dipole')
>>> n_repeat = 3
>>> job.structure = pr.create.structure.bulk('Al', cubic=True).repeat(n_repeat)
>>> job.structure[0] = 'Ni'
>>> job.calc_minimize()
>>> job.run()
>>> dipole_tensor = -job.structure.get_volume() * job['output/generic/pressures'][-1]

Instead of working with atomistic calculations, the dipole tensor can be calculated by the lambda tensor [1], which is defined as:

where \(c\) is the concentration of the defect, \(V\) is the volume and \(\varepsilon\) is the strain field. Then the dipole tensor is given by:

ref:

[1] Nowick, Arthur S. Anelastic relaxation in crystalline solids. Vol. 1. Elsevier, 2012.

get_point_defect_strain(positions: ndarray, dipole_tensor: ndarray, n_mesh: int = 100, optimize: bool = True, check_unique: bool = False) → ndarray[source]

Strain field around a point defect using the Green’s function method

Parameters:

positions ((n,3)-array) – Positions in real space or reciprocal space (if fourier=True).
dipole_tensor ((3,3)-array) – Dipole tensor
n_mesh (int) – Number of mesh points in the radial integration in case if anisotropic Green’s function (ignored if isotropic=True or fourier=True)
optimize (bool) – cf. optimize in numpy.einsum
check_unique (bool) – Whether to check the unique positions

Returns:

Strain field

Return type:

((n,3,3)-array)

According to the definition of the Green’s function (cf. docstring of get_greens_function):

where \(u_i(r)\) is the displacement field of component \(i\) at position \(r\) and \(f_j(a)\) is the force component \(j\) of the atom at position \(a\). By taking the polynomial development we obtain:

The first term disappears because the sum of the forces is zero. From the second term we define the dipole tensor \(P_{jk} = a_k f_j(a)\). Following the definition above, we can obtain the displacement field, strain field, stress field and energy density field if the dipole tensor and the elastic tensor are known.

The dipole tensor of a point defect is commonly obtained from the following equation:

where \(U\) is the potential energy, \(V\) is the volume and \(\varepsilon\) is the strain field. At equilibrium, the derivative of the potential energy with respect to the strain disappears:

With this in mind, we can calculate the dipole tensor of Ni in Al with the following lines:

>>> from pyiron_atomistics import Project
>>> pr = Project('dipole_tensor')
>>> job = pr.create.job.Lammps('dipole')
>>> n_repeat = 3
>>> job.structure = pr.create.structure.bulk('Al', cubic=True).repeat(n_repeat)
>>> job.structure[0] = 'Ni'
>>> job.calc_minimize()
>>> job.run()
>>> dipole_tensor = -job.structure.get_volume() * job['output/generic/pressures'][-1]

Instead of working with atomistic calculations, the dipole tensor can be calculated by the lambda tensor [1], which is defined as:

where \(c\) is the concentration of the defect, \(V\) is the volume and \(\varepsilon\) is the strain field. Then the dipole tensor is given by:

ref:

[1] Nowick, Arthur S. Anelastic relaxation in crystalline solids. Vol. 1. Elsevier, 2012.

get_point_defect_stress(positions: ndarray, dipole_tensor: ndarray, n_mesh: int = 100, optimize: bool = True) → ndarray[source]

Stress field around a point defect using the Green’s function method

Parameters:

positions ((n,3)-array) – Positions in real space or reciprocal space (if fourier=True).
dipole_tensor ((3,3)-array) – Dipole tensor
n_mesh (int) – Number of mesh points in the radial integration in case if anisotropic Green’s function (ignored if isotropic=True or fourier=True)
optimize (bool) – cf. optimize in numpy.einsum

Returns:

Stress field

Return type:

((n,3,3)-array)

According to the definition of the Green’s function (cf. docstring of get_greens_function):

where \(u_i(r)\) is the displacement field of component \(i\) at position \(r\) and \(f_j(a)\) is the force component \(j\) of the atom at position \(a\). By taking the polynomial development we obtain:

The first term disappears because the sum of the forces is zero. From the second term we define the dipole tensor \(P_{jk} = a_k f_j(a)\). Following the definition above, we can obtain the displacement field, strain field, stress field and energy density field if the dipole tensor and the elastic tensor are known.

The dipole tensor of a point defect is commonly obtained from the following equation:

where \(U\) is the potential energy, \(V\) is the volume and \(\varepsilon\) is the strain field. At equilibrium, the derivative of the potential energy with respect to the strain disappears:

With this in mind, we can calculate the dipole tensor of Ni in Al with the following lines:

>>> from pyiron_atomistics import Project
>>> pr = Project('dipole_tensor')
>>> job = pr.create.job.Lammps('dipole')
>>> n_repeat = 3
>>> job.structure = pr.create.structure.bulk('Al', cubic=True).repeat(n_repeat)
>>> job.structure[0] = 'Ni'
>>> job.calc_minimize()
>>> job.run()
>>> dipole_tensor = -job.structure.get_volume() * job['output/generic/pressures'][-1]

Instead of working with atomistic calculations, the dipole tensor can be calculated by the lambda tensor [1], which is defined as:

where \(c\) is the concentration of the defect, \(V\) is the volume and \(\varepsilon\) is the strain field. Then the dipole tensor is given by:

ref:

[1] Nowick, Arthur S. Anelastic relaxation in crystalline solids. Vol. 1. Elsevier, 2012.

get_reuss_average() → LinearElasticity[source]: Reuss average of the elastic tensor. The Reuss average is obtained from the arithmetic mean of the compliance tensor. The Reuss average is isotropic and can be used to calculate the elastic moduli of the medium.

get_voigt_average() → LinearElasticity[source]: Voigt average of the elastic tensor. The Voigt average is defined as the arithmetic mean of the elastic tensor. The Voigt average is isotropic and can be used to calculate the elastic moduli of the medium.

get_zener_ratio() → float[source]: Zener ratio or the anisotropy index. If 1, the medium is isotropic. If isotropic, the analytical form of the Green’s function is used for the calculation of strain and displacement fields.

is_cubic() → bool[source]

is_isotropic() → bool[source]

property orientation: ndarray | None

Rotation matrix that defines the orientation of the system. If set, the elastic tensor will be rotated accordingly. For example a box with a dislocation should get:

>>> medium.orientation = np.array([[1, 1, 1], [1, 0, -1], [1, -2, 1]])

If a non-orthogonal orientation is set, the second vector is orthogonalized with the Gram Schmidt process. It is not necessary to specify the third axis as it is automatically calculated.