Complete implementation of Cigre 207, extract common solar heating

linerate/equations/cigre207/convective_cooling.py

@@ -299,7 +299,6 @@ def _compute_horizontal_natural_nusselt_number(
     elif GrPr < 1e6:
         return 0.480 * GrPr**0.250
     else:
-        warnings.warn("GrPr out of bounds: Must be < 10^6.")
         # Outside table range, what should we do here?
         return 0.125 * GrPr**0.333
 

linerate/equations/cigre207/solar_heating.py

@@ -1,80 +1,41 @@
 import numpy as np
 from numpy import pi
 
-from ...units import Meter, Unitless, WattPerMeter, WattPerSquareMeter
+from .. import cigre601
+from ...units import Meter, Unitless, WattPerSquareMeter
 
 
 def compute_direct_solar_radiation(
     sin_solar_altitude: Unitless,
-    clearness_ratio: Unitless,
     height_above_sea_level: Meter,
 ) -> WattPerSquareMeter:
     r"""Compute the direct solar radiation.
 
-    Equation (10-11) on page 19 of :cite:p:`cigre601`. Equation (10) states that the direct solar
+    On page 38 of :cite:p:`cigre601`. Equation (10) states that the direct solar
     radiation on a surface normal to the solar beam at sea level, :math:`I_{B(0)}`, is given by
 
     .. math::
 
         N_s \frac{1280 \sin(H_s)}{\sin(H_s) + 0.314},
 
-    where :math:`N_s` is the clearness ratio which is used to adjust the amount of radiation
-    compared to what goes through a standard Indian atmosphere, and :math:`H_s` is the solar
-    altitude.
-
-    While the solar radiation model is based on :cite:p:`sharma1965interrelationships` and
-    therefore have parameters estimated for an Indian atmosphere, it gives comparable results to
-    the solar radiation model in the IEEE standard :cite:p:`ieee738`. It is therefore reasonable to
-    assume that the parameters work in other climates as well.
+    where :math:`H_s` is the solar altitude.
+    To correct for height above sea level, we use the Eq. 19 from Cigre 601, since no equation is provided
+    in Cigre 207.
 
     Parameters
     ----------
     sin_solar_altitude:
         :math:`\sin\left(H_s\right)`. The sine of the solar altitude.
-    clearness_ratio:
-        :math:`N_s`. The clearness ratio (or clearness number in
-        :cite:p:`sharma1965interrelationships,cigre207`).
     height_above_sea_level:
         :math:`y~\left[\text{m}\right]`. The conductor's altitude.
 
     Returns
     -------
     Union[float, float64, ndarray[Any, dtype[float64]]]
         :math:`I_B~\left[\text{W}~\text{m}^{-2}\right]`. The direct solar radiation.
-
-    Note
-    ----
-    The 1280 originates and 0.314 in the above equation originates from
-    :cite:p:`sharma1965interrelationships`, which is cited in :cite:p:`morgan1982thermal` (which is
-    listed as the reference in :cite:p:`cigre601`). In :cite:p:`sharma1965interrelationships` the
-    empirical relationship
-
-    .. math::
-
-        I_{B(0)} = \frac{1.842 \sin(H_s)}{\sin(H_s) + 0.3135}~\text{Ly}~\text{min}^{-1}
-
-    is introduced, and by converting from Langley per minute to :math:`\text{W}~\text{m}^{-2}`, we
-    obtain
-
-    .. math::
-
-        I_{B(0)} = N_s \frac{1284.488 \sin(H_s)}{\sin(H_s) + 0.3135}~\text{W}~\text{m}^{-2},
-
-    which is equal to the equation we use (with three significant digits).
     """
-    sin_H_s = sin_solar_altitude
-    N_s = clearness_ratio
-    y = height_above_sea_level
-
-    I_B_0 = N_s * 1280 * sin_H_s / (sin_H_s + 0.314)
-    # Equation 19 says that
-    # I_B = I_B_0 * (1 + 1.4e-4 * y * (1367/I_B_0 - 1))
-    # However, if I_B_0 = 0, this will divide by 0. To return NaN-values if and only
-    # if the input is NaN, we therefore reformulate it
-    scaled_y = 1.4e-4 * y
-    I_B = I_B_0 * (1 - scaled_y) + 1367 * scaled_y
-
-    return np.where(sin_H_s >= 0, I_B, 0 * I_B)
+    clearness_ratio = 1.0
+    return cigre601.solar_heating.compute_direct_solar_radiation(sin_solar_altitude, clearness_ratio, height_above_sea_level)
 
 
 def compute_diffuse_sky_radiation(
@@ -83,9 +44,9 @@ def compute_diffuse_sky_radiation(
 ) -> WattPerSquareMeter:
     r"""Compute the diffuse radiation (light scattered in the atmosphere).
 
-    Equation (13) on page 20 of :cite:p:`cigre601`.
+    On page 38 of :cite:p:`cigre207`.
 
-    This equation differ from :cite:p:`cigre207`, however the difference is small, and the
+    This equation differ from :cite:p:`cigre601`, however the difference is small, and the
     diffuse radiation is a small contributor to the overall solar radiation, so the total
     discrepancy between the models is small.
 
@@ -103,7 +64,7 @@ def compute_diffuse_sky_radiation(
     """
     sin_H_s = sin_solar_altitude
     I_B = direct_solar_radiation
-    return np.maximum(0, (430.5 - 0.3288 * I_B)) * np.maximum(0, sin_H_s)
+    return np.maximum(0, (570 - 0.47 * I_B)) * np.maximum(0, sin_H_s)**1.2
 
 
 def compute_global_radiation_intensity(
@@ -115,14 +76,14 @@ def compute_global_radiation_intensity(
 ) -> WattPerSquareMeter:
     r"""Compute the global radiation intensity experienced by the conductor.
 
-    Equation (9) on page 18 of :cite:p:`cigre601` state that the global radiation intensity,
+    Equation (47) on page 38 of :cite:p:`cigre207` state that the global radiation intensity,
     :math:`I_T`, is given by
 
     .. math::
 
         I_T =
             I_B \left(\sin(\eta) + 0.5 F \pi \sin(H_s)\right) +
-            I_d \left(1 + 0.5 F \pi\right),
+            0.5 I_d \pi \left(1 + F \right),
 
     where :math:`\eta` is the incidence angle of the sun on the line, :math:`H_s` is the solar
     altitude and :math:`F` is the ground albedo (amount of radiation diffusely reflected from the
@@ -181,36 +142,4 @@ def compute_global_radiation_intensity(
     sin_eta = sin_angle_of_sun_on_line
     F_pi_half = 0.5 * pi * F
 
-    return I_B * (sin_eta + F_pi_half * sin_H_s) + I_d * (1 + F_pi_half)  # type: ignore
-
-
-def compute_solar_heating(
-    absorptivity: Unitless,
-    global_radiation_intensity: WattPerSquareMeter,
-    conductor_diameter: Meter,
-) -> WattPerMeter:
-    r"""Compute the solar heating experienced by the conductor.
-
-    Equation (8) on page 18 of :cite:p:`cigre601`.
-
-    Parameters
-    ----------
-    absorptivity:
-        :math:`\alpha_s`. Material constant. According to :cite:p:`cigre601`, it starts at
-        approximately 0.2 for new cables and reaches a constant value of approximately 0.9
-        after about one year.
-    global_radiation_intensity:
-        :math:`I_T~\left[\text{W}~\text{m}^{-2}\right]`.The global radiation intensity.
-    conductor_diameter:
-        :math:`D~\left[\text{m}\right]`. Outer diameter of the conductor.
-
-    Returns
-    -------
-    Union[float, float64, ndarray[Any, dtype[float64]]]
-        :math:`P_S~\left[\text{W}~\text{m}^{-1}\right]`. The solar heating of the conductor
-    """
-    alpha_s = absorptivity
-    I_T = global_radiation_intensity
-    D = conductor_diameter
-
-    return alpha_s * I_T * D
+    return I_B * (sin_eta + F_pi_half * sin_H_s) + I_d * pi/2 * (1 + F)  # type: ignore

linerate/equations/cigre601/solar_heating.py

@@ -85,7 +85,7 @@ def compute_diffuse_sky_radiation(
 
     Equation (13) on page 20 of :cite:p:`cigre601`.
 
-    This equation differ from :cite:p:`cigre207`, however the difference is small, and the
+    This equation differs from :cite:p:`cigre207`, however the difference is small, and the
     diffuse radiation is a small contributor to the overall solar radiation, so the total
     discrepancy between the models is small.
 
@@ -182,35 +182,3 @@ def compute_global_radiation_intensity(
     F_pi_half = 0.5 * pi * F
 
     return I_B * (sin_eta + F_pi_half * sin_H_s) + I_d * (1 + F_pi_half)  # type: ignore
-
-
-def compute_solar_heating(
-    absorptivity: Unitless,
-    global_radiation_intensity: WattPerSquareMeter,
-    conductor_diameter: Meter,
-) -> WattPerMeter:
-    r"""Compute the solar heating experienced by the conductor.
-
-    Equation (8) on page 18 of :cite:p:`cigre601`.
-
-    Parameters
-    ----------
-    absorptivity:
-        :math:`\alpha_s`. Material constant. According to :cite:p:`cigre601`, it starts at
-        approximately 0.2 for new cables and reaches a constant value of approximately 0.9
-        after about one year.
-    global_radiation_intensity:
-        :math:`I_T~\left[\text{W}~\text{m}^{-2}\right]`.The global radiation intensity.
-    conductor_diameter:
-        :math:`D~\left[\text{m}\right]`. Outer diameter of the conductor.
-
-    Returns
-    -------
-    Union[float, float64, ndarray[Any, dtype[float64]]]
-        :math:`P_S~\left[\text{W}~\text{m}^{-1}\right]`. The solar heating of the conductor
-    """
-    alpha_s = absorptivity
-    I_T = global_radiation_intensity
-    D = conductor_diameter
-
-    return alpha_s * I_T * D

linerate/equations/solar_heating.py

@@ -0,0 +1,33 @@
+from linerate.units import WattPerSquareMeter, Unitless, Meter, WattPerMeter
+
+
+def compute_solar_heating(
+    absorptivity: Unitless,
+    global_radiation_intensity: WattPerSquareMeter,
+    conductor_diameter: Meter,
+) -> WattPerMeter:
+    r"""Compute the solar heating experienced by the conductor.
+
+    Equation (8) on page 18 of :cite:p:`cigre601` and (11) on page 4 in :cite:p:`cigre207`.
+
+    Parameters
+    ----------
+    absorptivity:
+        :math:`\alpha_s`. Material constant. According to :cite:p:`cigre601`, it starts at
+        approximately 0.2 for new cables and reaches a constant value of approximately 0.9
+        after about one year.
+    global_radiation_intensity:
+        :math:`I_T~\left[\text{W}~\text{m}^{-2}\right]`.The global radiation intensity.
+    conductor_diameter:
+        :math:`D~\left[\text{m}\right]`. Outer diameter of the conductor.
+
+    Returns
+    -------
+    Union[float, float64, ndarray[Any, dtype[float64]]]
+        :math:`P_S~\left[\text{W}~\text{m}^{-1}\right]`. The solar heating of the conductor
+    """
+    alpha_s = absorptivity
+    I_T = global_radiation_intensity
+    D = conductor_diameter
+
+    return alpha_s * I_T * D

linerate/model.py

@@ -19,7 +19,7 @@
     joule_heating,
     math,
     radiative_cooling,
-    solar_angles, convective_cooling, dimensionless,
+    solar_angles, convective_cooling, dimensionless, solar_heating,
 )
 from linerate.equations.math import switch_cos_sin
 from linerate.types import Span, Weather
@@ -355,7 +355,7 @@ def compute_solar_heating(
         I_T = cigre601.solar_heating.compute_global_radiation_intensity(
             I_B, I_d, F, sin_eta, sin_H_s
         )
-        return cigre601.solar_heating.compute_solar_heating(
+        return solar_heating.compute_solar_heating(
             alpha_s,
             I_T,
             D,

linerate/models/Cigre207.md

@@ -15,16 +15,22 @@ We do not implement skin effect, since this is included in the AC resistance val
 ## Magnetic effect
 
 The effect of magnetic heating seems to be intertwined with the skin effect in Cigre 207.
-To make things simpler and hopefully accurate enough, we use the same linearised model as our implementation of
-Cigre 601.
-**Maybe we have to come back to this one...**
+To make things simpler and hopefully accurate enough, we use the same linearised model for magnetic effects in ACSR
+lines as our implementation of Cigre 601.
 
 ## Solar heating
 
-Solar heating is the same as for Cigre 601, but with clearness ratio 1.
+Cigre 207 assumes that the diffuse radiation is uniformly directed, this differs from Cigre 601.
+Apart from this, solar heating is the same as for Cigre 601, with clearness ratio 1.
 
 ## Convective cooling
+
 There are a few inconsistencies in the Cigre 207 modelling of convective cooling.
 The definition of the Reynolds number is given as
 $Re = \rho_r V D/\nu_f$, where \rho_r is the air density relative to density at sea level.
-This factor is not present in the ordinary definition of Reynolds number, so we choose to omit it.
+This factor is not present in the ordinary definition of Reynolds number, but seems to indicate that the kinematic
+viscosity $\nu_f$ has to be corrected for the air relative density.
+
+Cigre 207 has a non-smooth transition to low wind speeds (<0.5m/s) at angles of attack less than 45 degrees.
+This leads to an increase in cooling as the wind speed drops below 0.5m/s, which seems non-physical.
+This discontinuity is removed in Cigre 601.

linerate/models/cigre207.py

@@ -1,7 +1,8 @@
 from abc import abstractmethod
 
 from linerate import ThermalModel, Span, Weather
-from linerate.equations import solar_angles, cigre601, math, cigre207, dimensionless, convective_cooling, joule_heating
+from linerate.equations import solar_angles, cigre601, math, cigre207, dimensionless, convective_cooling, joule_heating, \
+    solar_heating
 from linerate.equations.math import switch_cos_sin
 from linerate.model import _copy_method_docstring
 from linerate.units import Date, Celsius, Ampere, WattPerMeter, OhmPerMeter
@@ -53,7 +54,7 @@ def compute_solar_heating(
         I_T = cigre601.solar_heating.compute_global_radiation_intensity(
             I_B, I_d, F, sin_eta, sin_H_s
         )
-        return cigre601.solar_heating.compute_solar_heating(
+        return solar_heating.compute_solar_heating(
             alpha_s,
             I_T,
             D,