diff --git a/elaston/linear_elasticity/green.py b/elaston/linear_elasticity/green.py index c52e233..1797a61 100644 --- a/elaston/linear_elasticity/green.py +++ b/elaston/linear_elasticity/green.py @@ -5,6 +5,7 @@ import numpy as np from elaston.linear_elasticity import tools from tqdm.auto import tqdm +from functools import cached_property __author__ = "Sam Waseda" __copyright__ = ( @@ -58,7 +59,7 @@ def _get_greens_function(self, r, derivative=0, fourier=False): ) def get_greens_function( - self, r, derivative=0, fourier=False, check_unique=False, save_memory=False + self, r, derivative=0, fourier=False, check_unique=False ): """ Args: @@ -68,8 +69,6 @@ def get_greens_function( check_unique (bool): If `True`, pyiron checks whether there are duplicate values and avoids calculating the Green's function multiple times for the same values - save_memory (bool): If `True`, for loop will be used instead of - vector calculation, which saves the memory but makes it slow. Returns: (numpy.array): Green's function values. If `derivative=0` or `fourier=True`, @@ -78,19 +77,9 @@ def get_greens_function( x = np.array(r) if check_unique: x, inv = np.unique(x.reshape(-1, 3), axis=0, return_inverse=True) - if save_memory: - g_tmp = np.array( - [ - self._get_greens_function( - r=xx, derivative=derivative, fourier=fourier - ) - for xx in tqdm(x) - ] - ) - else: - g_tmp = self._get_greens_function( - r=x, derivative=derivative, fourier=fourier - ) + g_tmp = self._get_greens_function( + r=x, derivative=derivative, fourier=fourier + ) if check_unique: g_tmp = g_tmp[inv].reshape(np.asarray(r).shape + (derivative + 1) * (3,)) return g_tmp @@ -120,15 +109,13 @@ def __init__(self, poissons_ratio, shear_modulus, min_distance=0, optimize=True) self.shear_modulus = shear_modulus self.min_dist = min_distance self.optimize = optimize - self._A = None self._B = None - @property + # Rewrite A using cached_property + @cached_property def A(self): """First coefficient of the Green's function. For more, cf. DocString in the class level.""" - if self._A is None: - self._A = (3 - 4 * self.poissons_ratio) * self.B - return self._A + return (3 - 4 * self.poissons_ratio) * self.B @property def B(self): diff --git a/elaston/linear_elasticity/linear_elasticity.py b/elaston/linear_elasticity/linear_elasticity.py index c7411ab..c993d9f 100644 --- a/elaston/linear_elasticity/linear_elasticity.py +++ b/elaston/linear_elasticity/linear_elasticity.py @@ -20,16 +20,6 @@ __date__ = "Aug 21, 2021" -def value_or_none(func): - def f(self): - if self.elastic_tensor is None: - return None - else: - return func(self) - - return f - - point_defect_explanation = """ According to the definition of the Green's function (cf. docstring of `get_greens_function`): @@ -73,6 +63,18 @@ def f(self): 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, which is defined as: + +.. math: + \\lambda_{ij} = \\frac{1]{V} \\frac{\\partial \\varepsilon_{ij}}{\\partial c} + +where :math:`c` is the concentration of the defect, :math:`V` is the volume +and :math:`\\varepsilon` is the strain field. Then the dipole tensor is given by: + +.. math: + P_{ij} = VC_{ijkl}\\lambda_{kl} """ @@ -197,7 +199,6 @@ def elastic_tensor(self, C): self._elastic_tensor = C @property - @value_or_none def elastic_tensor_voigt(self): """ Voigt notation of the elastic tensor, i.e. (i, j) = i, if i == j and @@ -206,7 +207,6 @@ def elastic_tensor_voigt(self): return tools.C_to_voigt(self.elastic_tensor) @property - @value_or_none def compliance_matrix(self): """Compliance matrix in Voigt notation.""" return np.linalg.inv(self.elastic_tensor_voigt) @@ -273,6 +273,17 @@ def youngs_modulus(self): """ Returns: ((3,)-array): xx-, yy-, zz-components of Young's modulus + + + Here is a list of Young's modulus of a few materials: + + - Al: 70.4 + - Cu: 170.0 + - Fe: 211.0 + - Mo: 442.0 + - Ni: 248.0 + - W: 630.0 + """ return 1 / self.compliance_matrix[:3, :3].diagonal() @@ -285,7 +296,6 @@ def get_greens_function( isotropic: bool = False, optimize: bool = True, check_unique: bool = False, - save_memory: bool = False, ): """ Green's function of the equilibrium condition: @@ -303,8 +313,6 @@ def get_greens_function( is isotropic, it will automatically be set to isotropic=True optimize (bool): cf. `optimize` in `numpy.einsum` check_unique (bool): Whether to check the unique positions - save_memory (bool): Whether to save memory by using a for loop - instead of `numpy.einsum` Returns: ((n,3,3)-array): Green's function values for the given positions @@ -320,7 +328,6 @@ def get_greens_function( derivative=derivative, fourier=fourier, check_unique=check_unique, - save_memory=save_memory, ) get_greens_function.__doc__ += Green.__doc__ @@ -333,7 +340,6 @@ def get_point_defect_displacement( isotropic: bool = False, optimize: bool = True, check_unique: bool = False, - save_memory: bool = False, ): """ Displacement field around a point defect @@ -347,8 +353,6 @@ def get_point_defect_displacement( is isotropic, it will automatically be set to isotropic=True optimize (bool): cf. `optimize` in `numpy.einsum` check_unique (bool): Whether to check the unique positions - save_memory (bool): Whether to save memory by using a for loop - instead of `numpy.einsum` Returns: ((n,3)-array): Displacement field @@ -361,7 +365,6 @@ def get_point_defect_displacement( isotropic=isotropic, optimize=optimize, check_unique=check_unique, - save_memory=save_memory, ) return -np.einsum("...ijk,...jk->...i", g_tmp, dipole_tensor) @@ -375,7 +378,6 @@ def get_point_defect_strain( isotropic: bool = False, optimize: bool = True, check_unique: bool = False, - save_memory: bool = False, ): """ Strain field around a point defect using the Green's function method @@ -389,8 +391,6 @@ def get_point_defect_strain( is isotropic, it will automatically be set to isotropic=True optimize (bool): cf. `optimize` in `numpy.einsum` check_unique (bool): Whether to check the unique positions - save_memory (bool): Whether to save memory by using a for loop - instead of `numpy.einsum` Returns: ((n,3,3)-array): Strain field @@ -403,7 +403,6 @@ def get_point_defect_strain( isotropic=isotropic, optimize=optimize, check_unique=check_unique, - save_memory=save_memory, ) v = -np.einsum("...ijkl,...kl->...ij", g_tmp, dipole_tensor) return 0.5 * (v + np.einsum("...ij->...ji", v)) diff --git a/tests/unit/elasticity/test_green.py b/tests/unit/elasticity/test_green.py index 4ad1e76..96cc8e9 100644 --- a/tests/unit/elasticity/test_green.py +++ b/tests/unit/elasticity/test_green.py @@ -59,14 +59,8 @@ def test_memory(self): aniso = Anisotropic(create_random_C()) positions = np.tile(np.random.random(3), 2).reshape(-1, 3) G_normal = aniso.get_greens_function(positions) - G_save = aniso.get_greens_function(positions, save_memory=True) - self.assertTrue(np.all(np.isclose(G_normal, G_save))) G_unique = aniso.get_greens_function(positions, check_unique=True) self.assertTrue(np.all(np.isclose(G_normal, G_unique))) - G = aniso.get_greens_function( - positions, check_unique=True, save_memory=True - ) - self.assertTrue(np.all(np.isclose(G_normal, G))) if __name__ == "__main__":