pyxtal.elasticity module
- pyxtal.elasticity.Voigt_6_to_full_3x3_strain(strain_vector)[source]
Form a 3x3 strain matrix from a 6 component vector in Voigt notation
- pyxtal.elasticity.Voigt_6_to_full_3x3_stress(stress_vector)[source]
Form a 3x3 stress matrix from a 6 component vector in Voigt notation
- pyxtal.elasticity.Voigt_6x6_to_cubic(C)[source]
Convert the Voigt 6x6 representation into the cubic elastic constants C11, C12 and C44.
- pyxtal.elasticity.Voigt_6x6_to_full_3x3x3x3(C)[source]
Convert from the Voigt representation of the stiffness matrix to the full 3x3x3x3 representation.
- Parameters:
C (array_like) – 6x6 stiffness matrix (Voigt notation).
- Returns:
C – 3x3x3x3 stiffness matrix.
- Return type:
array_like
- pyxtal.elasticity.elastic_moduli(C, l=array([1, 0, 0]), R=None, tol=1e-06)[source]
Calculate elastic moduli from 6x6 elastic constant matrix C_{ij}.
The elastic moduli calculated are: Young’s muduli, Poisson’s ratios, shear moduli, bulk mudulus and bulk mudulus tensor.
If a direction l is specified, the system is rotated to have it as its x direction (see Notes for details). If R is specified the system is rotated according to it.
- Parameters:
C (array_like) – 6x6 stiffness matrix (Voigt notation).
l (array_like, optional) – 3D direction vector for pull (the default is the x direction of the original system)
R (array_like, optional) – 3x3 rotation matrix.
tol (float, optional) – tolerance for checking validity of rotation and comparison of vectors.
- Returns:
E (array_like) – Young’s modulus for a stress in each of the three directions of the rotated system.
nu (array_like) – 3x3 matrix with Poisson’s ratios.
Gm (array_like) – 3x3 matrix with shear moduli.
B (float) – Bulk modulus.
K (array_like) – 3x3 matrix with bulk modulus tensor.
Notes
It works by rotating the elastic constant tensor to the desired direction, so it should be valid for arbitrary crystal structures. If only l is specified there is an infinite number of possible rotations. The chosen one is a rotation along the axis orthogonal to the plane defined by the vectors (1, 0, 0) and l.
Bulk modulus tensor as defined in O. Rand and V. Rovenski, “Analytical Methods in Anisotropic Elasticity”, Birkh”auser (2005), pp. 71.
- pyxtal.elasticity.elastic_properties(C)[source]
A quick summary of elastic properties from the 6*6 matrix
- pyxtal.elasticity.fit_elastic_constants(a, symmetry='triclinic', N_steps=5, delta=0.01, optimizer=None, verbose=True, GPa=True, tag='test', graphics=False, logfile=None, **kwargs)[source]
Compute elastic constants by linear regression of stress vs. strain
- Parameters:
a (ase.Atoms or list of ase.Atoms) – Either a single atomic configuration or a list of strained configurations. If a single configuration is given, it is passed
generate_strained_configs()
along with symmetry, N_steps, and delta to generate the set of strained configurations.symmetry (string) – Symmetry to use to determine which strain patterns are necessary. Default is ‘triclininc’, i.e. no symmetry.
N_steps (int) – Number of atomic configurations to generate for each strain pattern. Default is 5. Absolute strain values range from -delta*N_steps/2 to +delta*N_steps/2.
delta (float) – Strain increment for analytical derivatives of stresses. Default is 1e-2.
optimizer (ase.optimizer.*) – Optimizer to use for atomic positions (internal degrees of freedom) for each applied strain. Initial config a is not optimised, and should already have relaxed cell and atomic positions. Does not optimize atomic positions if set to None.
verbose (bool) – If True, print additional infomation about the quality of fit and summarise results of C_ij and estimated errors. Default True.
GPa (bool) – If True, convert to GPa
graphics (bool) – If True, use
matplotlib.pyplot
to plot the stress vs. strain curve for each C_ij component fitted. Default True.logfile (bool) – Log file to write optimizer output to. Default None (i.e. suppress output).
**kwargs (dict) – Additional arguments to pass to optimizer.run() method e.g. fmax.
- Returns:
C (array_like) – 6x6 matrix of the elastic constants in Voigt notation.
C_err (array_like) – If scipy.stats module is available then error estimates for each C_ij component are obtained from the accuracy of the linear regression. Otherwise an array of np.zeros((6,6)) is returned.
Notes
Code originally adapted from elastics.py script, available from http://github.com/djw/elastic-constants
- pyxtal.elasticity.full_3x3_to_Voigt_6_strain(strain_matrix)[source]
Form a 6 component strain vector in Voigt notation from a 3x3 matrix
- pyxtal.elasticity.full_3x3_to_Voigt_6_stress(stress_matrix)[source]
Form a 6 component stress vector in Voigt notation from a 3x3 matrix
- pyxtal.elasticity.full_3x3x3x3_to_Voigt_6x6(C)[source]
Convert from the full 3x3x3x3 representation of the stiffness matrix to the representation in Voigt notation. Checks symmetry in that process.
- pyxtal.elasticity.generate_strained_configs(at0, symmetry='triclinic', N_steps=5, delta=0.01)[source]
Generate a sequence of strained configurations
- Parameters:
a (ase.Atoms) – Bulk crystal configuration - both unit cell and atomic positions should be relaxed before calling this routine.
symmetry (string) – Symmetry to use to determine which strain patterns are necessary. Default is ‘triclininc’, i.e. no symmetry.
N_steps (int) – Number of atomic configurations to generate for each strain pattern. Default is 5. Absolute strain values range from -delta*N_steps/2 to +delta*N_steps/2.
delta (float) – Strain increment for analytical derivatives of stresses. Default 1e-2
- Returns:
Generator which yields a sequence of ase.Atoms instances corresponding
to the minima set of strained conifugurations required to evaluate the
full 6x6 C_ij matrix under the assumed symmetry.
- pyxtal.elasticity.invariants(s, syy=None, szz=None, syz=None, sxz=None, sxy=None, full_3x3_to_Voigt_6=<function full_3x3_to_Voigt_6_stress>)[source]
- pyxtal.elasticity.measure_triclinic_elastic_constants(a, delta=0.001, optimizer=None, logfile=None, **kwargs)[source]
Brute-force measurement of elastic constants for a triclinic (general) unit cell.
- Parameters:
a (ase.Atoms) – Atomic configuration.
optimizer (ase.optimizer.*) – Optimizer to use for atomic position. Does not optimize atomic position if set to None.
delta (float) – Strain increment for analytical derivatives of stresses.
- Returns:
C – 6x6 matrix of the elastic constants in Voigt notation.
- Return type:
array_like
- pyxtal.elasticity.poisson_ratio(C, l, m)[source]
Calculate approximate Poisson ratio
u_{lm} from 6x6 elastic constant matrix C_{ij}
This is the response in m direction to pulling in l direction. Result is dimensionless.
Formula is from W. Brantley, Calculated elastic constants for stress problems associated with semiconductor devices. J. Appl. Phys., 44, 534 (1973).
- pyxtal.elasticity.rotate_cubic_elastic_constants(C11, C12, C44, A, tol=1e-06)[source]
Return rotated elastic moduli for a cubic crystal given the elastic constant in standard C11, C12, C44 notation.
- Parameters:
C11 (float) – Cubic elastic moduli.
C12 (float) – Cubic elastic moduli.
C44 (float) – Cubic elastic moduli.
A (array_like) – 3x3 rotation matrix.
- Returns:
C – 6x6 matrix of rotated elastic constants (Voigt notation).
- Return type:
array
- pyxtal.elasticity.rotate_elastic_constants(C, A, tol=1e-06)[source]
Return rotated elastic moduli for a general crystal given the elastic constant in Voigt notation.
- Parameters:
C (array_like) – 6x6 matrix of elastic constants (Voigt notation).
A (array_like) – 3x3 rotation matrix.
- Returns:
C – 6x6 matrix of rotated elastic constants (Voigt notation).
- Return type:
array
- pyxtal.elasticity.youngs_modulus(C, l)[source]
Calculate approximate Youngs modulus E_l from 6x6 elastic constants matrix C_ij
This is the modulus for loading in the l direction. For the exact answer, taking into account elastic anisotropuy, rotate the C_ij matrix to the correct frame, compute the compliance matrix:
C = ... # 6x6 C_ij matrix in crystal frame A = ... # rotation matrix Cr = rotate_elastic_constants(C, A) S = np.inv(Cr) E_x = 1/S[0, 0] # Young's modulus for a pull in x direction E_y = 1/S[1, 1] # Young's modulus for a pull in y direction E_z = 1/S[0, 0] # Young's modulus for a pull in z direction
Notes
Formula is from W. Brantley, Calculated elastic constants for stress problems associated with semiconductor devices. J. Appl. Phys., 44, 534 (1973).