Merge pull request #30 from bjmorgan/warning_handling

add warning handler for overflows

changes.md

@@ -0,0 +1,6 @@
+# Change log
+
+## V1.1.0
+- new warning handling for overflow errors
+
+## V1.0.0 full as-reviewed release

docs/source/tutorial.ipynb

@@ -1,6 +1,7 @@
 {
  "cells": [
   {
+   "attachments": {},
    "cell_type": "markdown",
    "metadata": {},
    "source": [
@@ -35,6 +36,7 @@
    ]
   },
   {
+   "attachments": {},
    "cell_type": "markdown",
    "metadata": {},
    "source": [
@@ -64,6 +66,7 @@
    ]
   },
   {
+   "attachments": {},
    "cell_type": "markdown",
    "metadata": {},
    "source": [
@@ -96,6 +99,7 @@
    ]
   },
   {
+   "attachments": {},
    "cell_type": "markdown",
    "metadata": {},
    "source": [
@@ -119,6 +123,7 @@
    ]
   },
   {
+   "attachments": {},
    "cell_type": "markdown",
    "metadata": {},
    "source": [
@@ -135,6 +140,15 @@
      "name": "stdout",
      "output_type": "stream",
      "text": [
+      "DOSOverflowWarning: An overflow occurred during computation of\n",
+      "                          electron and hole concentrations. This is likely a natural result of the use of\n",
+      "                          a numerical solver for the Fermi energy search. This can likely be ignored\n",
+      "                          though you should always check the final results are reasonable.\n",
+      "DefectOverflowWarning: An overflow occurred during computation of\n",
+      "                          defect concentrations. This is likely a natural result of the use of\n",
+      "                          a numerical solver for the Fermi energy search. This can likely be ignored\n",
+      "                          though you should always check the final results are reasonable.\n",
+      "Temperature :      500  (K)\n",
       "SC Fermi level :      1.5000000000000888  (eV)\n",
       "Concentrations:\n",
       "n (electrons)  : 3.9275293083444446e-18 cm^-3\n",
@@ -153,25 +167,14 @@
       "           :  1  1.979342e+17          100.00 \n",
       "\n"
      ]
-    },
-    {
-     "name": "stderr",
-     "output_type": "stream",
-     "text": [
-      "/Users/alex/work/py-sc-fermi/py_sc_fermi/dos.py:231: RuntimeWarning: overflow encountered in exp\n",
-      "  + np.exp(\n",
-      "/Users/alex/work/py-sc-fermi/py_sc_fermi/dos.py:240: RuntimeWarning: overflow encountered in exp\n",
-      "  + np.exp((self.edos[self._n0_index() :] - e_fermi) / (kboltz * temperature))\n",
-      "/Users/alex/work/py-sc-fermi/py_sc_fermi/defect_charge_state.py:177: RuntimeWarning: overflow encountered in exp\n",
-      "  concentration = self.degeneracy * np.exp(expfac)\n"
-     ]
     }
    ],
    "source": [
     "defect_system.report()"
    ]
   },
   {
+   "attachments": {},
    "cell_type": "markdown",
    "metadata": {},
    "source": [
@@ -229,6 +232,7 @@
    ]
   },
   {
+   "attachments": {},
    "cell_type": "markdown",
    "metadata": {},
    "source": [
@@ -277,6 +281,7 @@
    ]
   },
   {
+   "attachments": {},
    "cell_type": "markdown",
    "metadata": {},
    "source": [
@@ -336,6 +341,7 @@
    ]
   },
   {
+   "attachments": {},
    "cell_type": "markdown",
    "metadata": {},
    "source": [
@@ -400,6 +406,7 @@
    ]
   },
   {
+   "attachments": {},
    "cell_type": "markdown",
    "metadata": {},
    "source": [
@@ -455,6 +462,7 @@
    ]
   },
   {
+   "attachments": {},
    "cell_type": "markdown",
    "metadata": {},
    "source": [
@@ -530,7 +538,7 @@
    "name": "python",
    "nbconvert_exporter": "python",
    "pygments_lexer": "ipython3",
-   "version": "3.10.8"
+   "version": "3.10.9"
   },
   "orig_nbformat": 4,
   "vscode": {

py_sc_fermi/__init__.py

@@ -1 +1 @@
-__version__ = "1.0.0"
+__version__ = "1.1.0"

py_sc_fermi/defect_system.py

@@ -3,6 +3,42 @@
 from py_sc_fermi.defect_species import DefectSpecies
 from py_sc_fermi.inputs import InputSet
 import numpy as np
+import warnings
+
+
+
+class CustomWarningManager:
+    def __init__(self):
+        self.dos_overflow_warning_issued = False
+        self.defect_overflow_warning_issued = False
+
+    def custom_warning(self, message, category, filename, lineno, file=None, line=None):
+        if category == RuntimeWarning:
+            if "dos" in str(filename):
+                if not self.dos_overflow_warning_issued:
+                    print(
+                        """DOSOverflowWarning: An overflow occurred during computation of
+                        electron and hole concentrations. This is likely a natural result of the use of
+                        a numerical solver for the Fermi energy search. This can likely be ignored
+                        though you should always check the final results are reasonable."""
+                    )
+                    self.dos_overflow_warning_issued = True
+            elif "defect" in str(filename):
+                if not self.defect_overflow_warning_issued:
+                    print(
+                        """DefectOverflowWarning: An overflow occurred during computation of
+                        defect concentrations. This is likely a natural result of the use of
+                        a numerical solver for the Fermi energy search. This can likely be ignored
+                        though you should always check the final results are reasonable."""
+                    )
+                    self.defect_overflow_warning_issued = True
+            else:
+                warnings.warn(message, category, filename, lineno, file, line)
+
+
+# Create a CustomWarningManager and set the custom_warning method as the warning handler
+warning_manager = CustomWarningManager()
+warnings.showwarning = warning_manager.custom_warning
 
 
 class DefectSystem(object):
@@ -153,29 +189,35 @@ def get_sc_fermi(self) -> Tuple[float, float]:
         reached_e_max = False
 
         # loop until convergence or max number of steps reached
-        for i in range(self.n_trial_steps):
-            q_tot = self.q_tot(e_fermi=e_fermi)
-            if e_fermi > emax:
-                if reached_e_min or reached_e_max:
-                    raise RuntimeError(f"No solution found between {emin} and {emax}")
-                reached_e_max = True
-                direction = -1.0
-            if e_fermi < emin:
-                if reached_e_max or reached_e_min:
-                    raise RuntimeError(f"No solution found between {emin} and {emax}")
-                reached_e_min = True
-                direction = +1.0
-            if abs(q_tot) < self.convergence_tolerance:
-                break
-            if q_tot > 0.0:
-                if direction == +1.0:
-                    step *= 0.25
+        with warnings.catch_warnings():
+            warnings.filterwarnings("once")
+            for i in range(self.n_trial_steps):
+                q_tot = self.q_tot(e_fermi=e_fermi)
+                if e_fermi > emax:
+                    if reached_e_min or reached_e_max:
+                        raise RuntimeError(
+                            f"No solution found between {emin} and {emax}"
+                        )
+                    reached_e_max = True
                     direction = -1.0
-            elif q_tot < 0.0:
-                if direction == -1.0:
-                    step *= 0.25
+                if e_fermi < emin:
+                    if reached_e_max or reached_e_min:
+                        raise RuntimeError(
+                            f"No solution found between {emin} and {emax}"
+                        )
+                    reached_e_min = True
                     direction = +1.0
-            e_fermi += step * direction
+                if abs(q_tot) < self.convergence_tolerance:
+                    break
+                if q_tot > 0.0:
+                    if direction == +1.0:
+                        step *= 0.25
+                        direction = -1.0
+                elif q_tot < 0.0:
+                    if direction == -1.0:
+                        step *= 0.25
+                        direction = +1.0
+                e_fermi += step * direction
 
         # return results
         residual = abs(q_tot)
@@ -284,7 +326,9 @@ def get_transition_levels(self) -> Dict[str, List[List]]:
         return transition_levels
 
     def as_dict(
-        self, decomposed: bool = False, per_volume: bool = True,
+        self,
+        decomposed: bool = False,
+        per_volume: bool = True,
     ) -> Dict[str, Any]:
         """Returns a dictionary of the properties of the ``DefectSystem`` object
         after solving for the self-consistent Fermi energy.
@@ -336,10 +380,10 @@ def as_dict(
             return {**run_stats, **decomp_concs}
 
     def site_percentages(
-        self, 
+        self,
     ) -> Dict[str, float]:
         """Returns a dictionary of the DefectSpecies in the DefectSystem which
-        giving the percentage of the sites in the structure that will host that 
+        giving the percentage of the sites in the structure that will host that
         defect.
 
         Returns:
@@ -350,9 +394,9 @@ def site_percentages(
         e_fermi = self.get_sc_fermi()[0]
 
         sum_concs = {
-                str(ds.name): float(
-                    (ds.get_concentration(e_fermi, self.temperature) / ds.nsites) * 100
-                )
-                for ds in self.defect_species
-            }
+            str(ds.name): float(
+                (ds.get_concentration(e_fermi, self.temperature) / ds.nsites) * 100
+            )
+            for ds in self.defect_species
+        }
         return sum_concs

py_sc_fermi/inputs.py

@@ -6,7 +6,7 @@
 from py_sc_fermi.dos import DOS
 from pymatgen.core import Structure
 from typing import Optional, List
-import yaml
+import yaml # type: ignore
 import os
 
 InputFermiData = namedtuple(

tests/test_defect_system.py

@@ -1,10 +1,12 @@
 import unittest
-from unittest.mock import Mock, PropertyMock, patch
+from unittest.mock import Mock, patch
+from io import StringIO
 
 import os
+import textwrap
 from py_sc_fermi.defect_species import DefectSpecies
 from py_sc_fermi.dos import DOS
-from py_sc_fermi.defect_system import DefectSystem
+from py_sc_fermi.defect_system import DefectSystem, CustomWarningManager
 from py_sc_fermi.defect_charge_state import DefectChargeState
 
 
@@ -24,6 +26,36 @@
 )
 
 
+class TestCustomWarningManager(unittest.TestCase):
+    def setUp(self):
+        self.warning_manager = CustomWarningManager()
+
+    @patch('sys.stdout', new_callable=StringIO)
+    def test_dos_overflow_warning(self, mock_stdout):
+        self.warning_manager.custom_warning('overflow', RuntimeWarning, 'dos_file.py', 42)
+        expected_warning = textwrap.dedent(
+                        """DOSOverflowWarning: An overflow occurred during computation of
+                        electron and hole concentrations. This is likely a natural result of the use of
+                        a numerical solver for the Fermi energy search. This can likely be ignored
+                        though you should always check the final results are reasonable.""")
+        self.assertEqual(mock_stdout.getvalue().strip(), expected_warning.strip())
+
+    @patch('sys.stdout', new_callable=StringIO)
+    def test_defect_overflow_warning(self, mock_stdout):
+        self.warning_manager.custom_warning('overflow', RuntimeWarning, 'defect_file.py', 42)
+        expected_warning = textwrap.dedent(
+                        """DefectOverflowWarning: An overflow occurred during computation of
+                        defect concentrations. This is likely a natural result of the use of
+                        a numerical solver for the Fermi energy search. This can likely be ignored
+                        though you should always check the final results are reasonable.""")
+        self.assertEqual(mock_stdout.getvalue().strip(), expected_warning.strip())
+
+    @patch('warnings.warn')
+    def test_other_warning(self, mock_warn):
+        self.warning_manager.custom_warning('other warning', RuntimeWarning, 'other_file.py', 42, None, None)
+        mock_warn.assert_called_once_with('other warning', RuntimeWarning, 'other_file.py', 42, None, None)
+
+
 class TestDefectSystemInit(unittest.TestCase):
     def test_defect_system_is_initialised(self):
         volume = 100