diff --git a/climada/hazard/test/test_trop_cyclone.py b/climada/hazard/test/test_trop_cyclone.py index 27a0851c9e..12ec1620b9 100644 --- a/climada/hazard/test/test_trop_cyclone.py +++ b/climada/hazard/test/test_trop_cyclone.py @@ -384,21 +384,21 @@ def test_holland_2010_pass(self): hol_x = _x_holland_2010(si_track, d_centr, close_centr) np.testing.assert_array_almost_equal(hol_x, [ [0.5, 0.500000, 0.485077, 0.476844, 0.457291], - [0.5, 0.500000, 0.410997, 0.289203, 0.000000], - [0.5, 0.497620, 0.131072, 0.000000, 0.000000], + [0.5, 0.500000, 0.410997, 0.400000, 0.000000], + [0.5, 0.497620, 0.400000, 0.400000, 0.400000], [0.5, 0.505022, 1.403952, 1.554611, 2.317948], ]) v_ang_norm = _stat_holland_2010(si_track, d_centr, close_centr, hol_x) - np.testing.assert_array_almost_equal(v_ang_norm, [ + np.testing.assert_allclose(v_ang_norm, [ # first column: converge to 0 when approaching storm eye # second column: vmax at RMW # fourth column: peripheral speed (17 m/s) at peripheral radius (unless x is clipped!) [0.0000000, 35.000000, 21.181497, 17.00000, 12.103461], - [1.2964800, 34.990037, 21.593755, 17.00000, 0.0000000], - [0.3219518, 15.997500, 13.585498, 16.00000, 16.000000], + [1.2964800, 34.990037, 21.593755, 12.89131, 0.0000000], + [0.3219518, 15.997500, 9.7120060, 8.087240, 6.2289690], [24.823469, 89.992938, 24.381965, 17.00000, 1.9292020], - ]) + ], atol=1e-5) def test_stat_holland_1980(self): """Test _stat_holland_1980 function. Compare to MATLAB reference.""" diff --git a/climada/hazard/trop_cyclone.py b/climada/hazard/trop_cyclone.py index b6d233760f..f54f7cedfa 100644 --- a/climada/hazard/trop_cyclone.py +++ b/climada/hazard/trop_cyclone.py @@ -1741,6 +1741,11 @@ def _x_holland_2010( r_max_norm = (r_max / r_n)**hol_b if vmax_in_brackets: x_n = np.log(v_n) / np.log(v_max_s**2 * r_max_norm * np.exp(1 - r_max_norm)) + + # Truncate to prevent maximum to be shifted too far away from RMW. The value 1.0 has been + # found to be reasonable by manual testing of thresholds. Note that this means that the + # peripheral wind speed v_n is not exactly attained in some cases. + x_n = np.fmin(x_n, 1.0) else: x_n = np.log(v_n / v_max_s) / np.log(r_max_norm * np.exp(1 - r_max_norm)) @@ -1748,11 +1753,10 @@ def _x_holland_2010( x_max = 0.5 hol_x[close_centr] = x_max + np.fmax(0, d_centr - r_max) * (x_n - x_max) / (r_n - r_max) - # Negative hol_x values appear when v_max_s is very close to or even lower than v_n (which - # should never happen in theory). In those cases, wind speeds might decrease outside of the eye - # wall and increase again towards the peripheral radius (which is actually unphysical). - # We clip hol_x to 0, otherwise wind speeds keep increasing indefinitely away from the eye: - hol_x[close_centr] = np.fmax(hol_x[close_centr], 0.0) + # Truncate to prevent wind speed from increasing again towards the peripheral radius (which is + # unphysical). A value of 0.4 has been found to be reasonable by manual testing of thresholds. + # Note that this means that the peripheral wind speed v_n is not exactly attained sometimes. + hol_x[close_centr] = np.fmax(hol_x[close_centr], 0.4) return hol_x