Merge pull request #16 from chrisroadmap/twolayer-year

change tropical day length and add options for two layer model

CHANGELOG.rst

@@ -16,6 +16,7 @@ master
 v0.2.0
 ------
 - added: `solar` module (`#14 <https://github.com/chrisroadmap/climateforcing/pull/14>`_)
+- fixed: `SECPERYEAR` for `twolayermodel` now uses tropical year length of 365.24219 days instead of 365 (`#16 <https://github.com/chrisroadmap/climateforcing/pull/16>`_)
 - changed: all docstrings converted to numpy style (`#15 <https://github.com/chrisroadmap/climateforcing/pull/15>`_)
 
 v0.1.1

src/climateforcing/constants.py

@@ -1,4 +1,4 @@
 """Module-level constants."""
 EARTHRADIUS = 6371000  # m
-SECPERYEAR = 60 * 60 * 24 * 365  # s   # should be 365.2425
+SECPERYEAR = 60 * 60 * 24 * 365.24219  # s
 STEFAN_BOLTZMANN = 5.670374419e-8  # W m-2 K-4

src/climateforcing/twolayermodel/__init__.py

@@ -314,6 +314,8 @@ class TwoLayer:  # pylint: disable=too-few-public-methods,too-many-instance-attr
         "cdeep": 109.0,
         "gamma": 0.67,
         "eff": 1.28,
+        "earthradius": EARTHRADIUS,
+        "secperyear": SECPERYEAR,
     }
 
     def __init__(self, **kwargs):  # pylint: disable=too-many-locals
@@ -332,12 +334,22 @@ def __init__(self, **kwargs):  # pylint: disable=too-many-locals
                 Heat exchange coefficient between upper and deep layer, W m-2 K-1
             eff : float
                 Efficacy of deep ocean.
+            earthradius : float, default=`climateforcing.constants.EARTHRADIUS`
+                Earth radius (m)
+            secperyear : float, default=`climateforcing.constants.SECPERYEAR`
+                Seconds in a year (s)
         """
         self.params = self.default.copy()
         self.params.update(kwargs)
 
         # care with unit handling
-        self.ntoa2joule = 4 * math.pi * EARTHRADIUS * EARTHRADIUS * SECPERYEAR * 1e-22
+        self.ntoa2joule = (
+            4
+            * math.pi
+            * self.params["earthradius"] ** 2
+            * self.params["secperyear"]
+            * 1e-22
+        )
 
         # Define derived constants
         cdeep_p = self.params["cdeep"] * self.params["eff"]

tests/test_twolayermodel.py

@@ -59,3 +59,8 @@ def test_twolayermodel():
     assert np.allclose(result.qtot, EXPECTED_RESULT.qtot)
     assert np.allclose(result.tg, EXPECTED_RESULT.tg)
     assert np.allclose(result.tlev, EXPECTED_RESULT.tlev)
+
+    # Change secperyear and verify OHC is now different
+    scm365 = TwoLayerModel(scm_in=result, secperyear=365 * 24 * 60 * 60)
+    result365 = scm365.run()
+    assert ~np.allclose(result365.ohc, EXPECTED_RESULT.ohc)