diff --git a/pyscses/site.py b/pyscses/site.py index c452e71..3346b7f 100644 --- a/pyscses/site.py +++ b/pyscses/site.py @@ -22,9 +22,21 @@ class Site: defects (list): List of DefectAtSite objects, containing the properties of all individual defects at the site. scaling (float): A scaling factor that can be applied in the charge calculation. valence (float): The charge of the defect present at the site (in atomic units). + saturation_parameter (float): Optional saturation parameter as described in + `Hendricks et al. Sol. Stat. Ionics (2002)`_ + and `Swift et al. Nature Comp. Sci. (2021)`_. + Setting `saturation_parameter` < `1.0` sets some proportion of excluded sites that are + unavailable for occupation by any explicit defects. + Default value is `1.0`, i.e., 100% of sites may be occupied. defects (list): List of Defect_Species objects for all defects present at the site. sites (list): List containing all x coordinates and corresponding defect segregation energies. + .. _Hendricks et al. Sol. Stat. Ionics (2002): + https://doi.org/10.1016/S0167-2738(02)00484-8 + + .. _Swift et al. Nature Comp. Sci. (2021): + https://doi.org/10.1038/s43588-021-00041-y + """ def __init__(self, @@ -33,7 +45,8 @@ def __init__(self, defect_species: List[DefectSpecies], defect_energies: List[float], scaling: Optional[np.ndarray] = None, - valence: float = 0.0) -> None: + valence: float = 0.0, + saturation_parameter: float = 1.0) -> None: """Initialise a Site object. Args: @@ -43,10 +56,21 @@ def __init__(self, defect_energies (list(float)): List of defect segregation energies for each defect species at this site. scaling (optional, list(float): Optional list of scaling factors for the net charge at this site. Default scaling for each defect species is 1.0. valence (optional, float): Optional formal valence for this site in the absence of any defects. Default is 0.0. + saturation_parameter (optional, float): Optional saturation parameter as described in + Hendricks et al. Sol. Stat. Ionics (2002) [#HendricksEtAl_SolStatIonics2002]_ + and Swift et al. Nature Comp. Sci. (2021). [#SwiftEtAl_NatureCompSci2021]_. + A saturation parameter < 1.0 introduces some proportion of excluded sites that are + unavailable for occupation by any explicit defects. Default is 1.0. Raises: ValueError if the number of DefectSpecies != the number of defect segregation energies != the number of scaling factors (if passed). + .. [HendricksEtAl_SolStatIonics2002]: + https://doi.org/10.1016/S0167-2738(02)00484-8 + + .. [SwiftEtAl_NatureCompSci2021]: + https://doi.org/10.1038/s43588-021-00041-y + """ if len(defect_species) != len(defect_energies): raise ValueError("len(defect_species) must be equal to len(defect_energies)") @@ -71,9 +95,10 @@ def __init__(self, self.scaling = np.ones_like(defect_energies, dtype=float) self.grid_point: Optional[GridPoint] = None self.valence = valence - self.fixed_defects = tuple(d for d in self.defects if d.fixed=True) - self.mobile_defects = tuple(d for d in self.defects if d.fixed=False) - self.alpha = 1.0 - sum((d.mole_fraction for d in fixed_defects)) + self.saturation_parameter = saturation_parameter + self.fixed_defects = tuple(d for d in self.defects if d.fixed) + self.mobile_defects = tuple(d for d in self.defects if not d.fixed) + self.alpha = self.saturation_parameter - sum((d.mole_fraction for d in self.fixed_defects)) def competing_defect_species(self) -> Dict[str, int]: """Returns a dictionary reporting the number of fixed and / or mobile defect species that can occupy this site. @@ -125,8 +150,8 @@ def average_local_energy(self, return self.grid_point.average_site_energy(method) else: raise ValueError("TODO") - - def probabilities(self: + + def probabilities(self, phi: float, temp: float) -> List[float]: """Calculates the probabilities of this site being occupied by each defect species. @@ -141,8 +166,8 @@ def probabilities(self: """ probabilities_dict = {} boltzmann_factors = {d.label: d.boltzmann_factor(phi, temp) for d in self.mobile_defects} - denominator = (self.alpha + - sum([d.mole_fraction * (boltzmann_factors[d.label] - 1.0) + denominator = (self.alpha + + sum([d.mole_fraction * (boltzmann_factors[d.label] - 1.0) for d in self.mobile_defects])) for defect in self.defects: if defect.fixed: @@ -150,7 +175,7 @@ def probabilities(self: else: numerator = self.alpha * defect.mole_fraction * boltzmann_factors[defect.label] probabilities_dict[defect.label] = numerator / denominator - return [probabilities_dict[d.label] for d in self.defects] + return [probabilities_dict[d.label] for d in self.defects] def defect_valences(self) -> np.ndarray: """Returns an array of valences for each defect from self.defects """ diff --git a/tests/test_site.py b/tests/test_site.py index 8c7f1b8..36ca3e6 100644 --- a/tests/test_site.py +++ b/tests/test_site.py @@ -22,11 +22,31 @@ def create_mock_defect_species(n): mock_defect_species.append(m) return mock_defect_species +def create_mock_defects_at_site(n): + labels = ['A', 'B', 'C', 'D', 'E'] + valence = [-2.0, -1.0, 0.0, 1.0, 2.0] + mole_fraction = [0.15, 0.25, 0.35, 0.45, 0.55] + mobility = [0.1, 0.2, 0.3, 0.4, 0.5] + energies = [-0.1, -0.2, -0.3, -0.4, -0.5] + mock_defects_at_site = [] + for i in range(n): + m = Mock(spec=DefectAtSite) + m.label = labels.pop() + m.valence = valence.pop() + m.mole_fraction = mole_fraction.pop() + m.mobility = mobility.pop() + m.energy = energies.pop() + m.fixed = False + mock_defects_at_site.append(m) + return mock_defects_at_site + class TestSiteInit(unittest.TestCase): def test_site_is_initialised(self): mock_defect_species = create_mock_defect_species(2) + mock_defects_at_site = create_mock_defects_at_site(2) with patch('pyscses.site.DefectAtSite', autospec=True) as mock_DefectAtSite: + mock_DefectAtSite.side_effect = mock_defects_at_site site = Site(label='A', x=1.5, defect_species=mock_defect_species, @@ -37,27 +57,58 @@ def test_site_is_initialised(self): self.assertEqual(site.defect_energies, [-0.2, +0.2]) np.testing.assert_equal(site.scaling, np.array([1.0, 1.0])) self.assertEqual(site.valence, 0.0) + self.assertEqual(site.saturation_parameter, 1.0) + self.assertEqual(site.fixed_defects, ()) + self.assertEqual(site.mobile_defects, tuple(mock_defects_at_site)) + self.assertEqual(site.alpha, 1.0) def test_site_is_initialised_with_optional_args(self): mock_defect_species = create_mock_defect_species(2) with patch('pyscses.site.DefectAtSite', autospec=True) as mock_DefectAtSite: + mock_DefectAtSite.side_effect = create_mock_defects_at_site(2) site = Site(label='B', x=1.5, defect_species=mock_defect_species, defect_energies=[-0.2, +0.2], scaling=[0.5, 0.4], - valence=-2.0) + valence=-2.0, + saturation_parameter=0.1) self.assertEqual(site.label, 'B') self.assertEqual(site.x, 1.5) self.assertEqual(site.defect_species, mock_defect_species) self.assertEqual(site.defect_energies, [-0.2, +0.2]) np.testing.assert_equal(site.scaling, np.array([0.5, 0.4])) self.assertEqual(site.valence, -2.0) + self.assertEqual(site.saturation_parameter, 0.1) + self.assertEqual(site.alpha, 0.1) + + def test_site_init_with_mixed_mobile_and_fixed_defects(self): + mock_defect_species = create_mock_defect_species(3) + mock_defects_at_site = create_mock_defects_at_site(3) + mock_defects_at_site[0].fixed = False + mock_defects_at_site[0].mole_fraction = 0.4 + mock_defects_at_site[1].fixed = True + mock_defects_at_site[1].mole_fraction = 0.3 + mock_defects_at_site[2].fixed = True + mock_defects_at_site[2].mole_fraction = 0.2 + with patch('pyscses.site.DefectAtSite', autospec=True) as mock_DefectAtSite: + mock_DefectAtSite.side_effect = mock_defects_at_site + site = Site(label='C', + x=1.5, + defect_species=mock_defect_species, + defect_energies=[-0.2, +0.2, 0.0]) + self.assertEqual(site.fixed_defects, (mock_defects_at_site[1], mock_defects_at_site[2])) + self.assertEqual(site.mobile_defects[0], mock_defects_at_site[0]) + self.assertEqual(site.alpha, 0.5) + # TODO: test that fixed defects are added to the site.fixed_defects tuple and mobile defects are added to the site.mobiile_defects tuple + + # TODO: test that site.alpha is computed correctly for various combinations of saturation_parameter and fixed defect species. def test_site_init_data_check_1(self): """Checks that initialising a Site object raises a ValueError if n(defect_species) != n(defect_energies)""" mock_defect_species = create_mock_defect_species(1) with patch('pyscses.site.DefectAtSite', autospec=True) as mock_DefectAtSite: + # mock_DefectAtSite.side_effect = create_mock_defects_at_site(1) with self.assertRaises(ValueError): site = Site(label='A', x=1.5, @@ -79,12 +130,13 @@ class TestSite(unittest.TestCase): def setUp(self): mock_defect_species = create_mock_defect_species(2) + mock_defects_at_site = create_mock_defects_at_site(2) with patch('pyscses.site.DefectAtSite', autospec=True) as mock_DefectAtSite: + mock_DefectAtSite.side_effect = mock_defects_at_site self.site = Site(label='A', x=1.5, defect_species=mock_defect_species, defect_energies=[-0.2, +0.2]) - self.site.defects = [Mock(spec=DefectAtSite), Mock(spec=DefectAtSite)] def test_defect_with_label(self):