diff --git a/.gitignore b/.gitignore
index 3da10150..fb57579a 100644
--- a/.gitignore
+++ b/.gitignore
@@ -94,6 +94,9 @@ celerybeat-schedule
 # dotenv
 .env
 
+# direnv
+.envrc
+
 # virtualenv
 .venv
 venv/
diff --git a/iodata/formats/extxyz.py b/iodata/formats/extxyz.py
index 7f6b2423..1db534ae 100644
--- a/iodata/formats/extxyz.py
+++ b/iodata/formats/extxyz.py
@@ -27,7 +27,6 @@
 
 """
 
-from distutils.util import strtobool
 import shlex
 from typing import Iterator
 
@@ -35,7 +34,7 @@
 
 from ..docstrings import document_load_one, document_load_many
 from ..periodic import sym2num, num2sym
-from ..utils import angstrom, amu, LineIterator
+from ..utils import angstrom, amu, LineIterator, strtobool
 
 from .xyz import load_one as load_one_xyz
 
diff --git a/iodata/formats/qchemlog.py b/iodata/formats/qchemlog.py
index 55d09d51..58af4cdb 100644
--- a/iodata/formats/qchemlog.py
+++ b/iodata/formats/qchemlog.py
@@ -22,14 +22,13 @@
 """
 
 from typing import Tuple
-from distutils.util import strtobool
 
 import numpy as np
 
 from ..docstrings import document_load_one
 from ..orbitals import MolecularOrbitals
 from ..periodic import sym2num
-from ..utils import LineIterator, angstrom, kcalmol, calmol, kjmol
+from ..utils import LineIterator, angstrom, kcalmol, calmol, kjmol, strtobool
 
 __all__ = []
 
diff --git a/iodata/test/test_utils.py b/iodata/test/test_utils.py
index 42b582b5..d3b55536 100644
--- a/iodata/test/test_utils.py
+++ b/iodata/test/test_utils.py
@@ -18,9 +18,18 @@
 # --
 """Unit tests for iodata.utils."""
 
+import pytest
 
-from ..utils import amu
+from ..utils import amu, strtobool
 
 
 def test_amu():
     assert abs(amu * 1.008 - 1837.47) < 1e-1
+
+
+def test_strtobool():
+    assert strtobool("T") is True
+    assert strtobool("false") is False
+    assert strtobool("y") is True
+    with pytest.raises(ValueError):
+        strtobool("whatever")
diff --git a/iodata/utils.py b/iodata/utils.py
index facb8e53..6035d952 100644
--- a/iodata/utils.py
+++ b/iodata/utils.py
@@ -31,7 +31,7 @@
 
 
 __all__ = ['LineIterator', 'Cube', 'set_four_index_element', 'volume',
-           'derive_naturals', 'check_dm']
+           'derive_naturals', 'check_dm', 'strtobool']
 
 
 # The unit conversion factors below can be used as follows:
@@ -263,3 +263,27 @@ def check_dm(dm: np.ndarray, overlap: np.ndarray, eps: float = 1e-4, occ_max: fl
     if occupations.max() > occ_max + eps:
         raise ValueError('The density matrix has eigenvalues considerably larger than '
                          'max. error=%e' % (occupations.max() - 1))
+
+
+STRTOBOOL = {
+    'y': True,
+    'yes': True,
+    't': True,
+    'true': True,
+    'on': True,
+    '1': True,
+    'n': False,
+    'no': False,
+    'f': False,
+    'false': False,
+    'off': False,
+    '0': False
+}
+
+
+def strtobool(value: str) -> bool:
+    """Interpret string as a boolean."""
+    result = STRTOBOOL.get(value.lower())
+    if result is None:
+        raise ValueError(f"'{value}' cannot be converted to boolean")
+    return result