Add zenith angle dependence to sea ice albedo calculations

dependencies.yaml

@@ -22,8 +22,8 @@ casim:
     ref: 2026.03.2
 
 jules:
-    source: git@github.com:MetOffice/jules.git
-    ref: 2026.03.2
+    source: git@github.com:DanCopsey/jules.git
+    ref: zenith_angle_seaice
 
 lfric_apps:
     source:

interfaces/jules_interface/build/extract.yaml

@@ -106,6 +106,7 @@ jules:
     - src/science/params/veg_param_mod.F90
     - src/science/radiation/Jin11_osa_mod.F90
     - src/science/radiation/albpft_jls_mod.F90
+    - src/science/radiation/albpond_mod.F90
     - src/science/radiation/albsnow_jls_mod.F90
     - src/science/radiation/albsnow_ts_jls_mod.F90
     - src/science/radiation/calc_direct_albsoil_mod.F90

interfaces/jules_interface/rose-meta/jules-lfric/.versions.py

@@ -0,0 +1,54 @@
+import re
+import sys
+
+from metomi.rose.upgrade import MacroUpgrade  # noqa: F401
+
+from .version30_31 import *
+
+
+class UpgradeError(Exception):
+    """Exception created when an upgrade fails."""
+
+    def __init__(self, msg):
+        self.msg = msg
+
+    def __repr__(self):
+        sys.tracebacklimit = 0
+        return self.msg
+
+    __str__ = __repr__
+
+
+"""
+Copy this template and complete to add your macro
+class vnXX_txxx(MacroUpgrade):
+    # Upgrade macro for <TICKET> by <Author>
+    BEFORE_TAG = "vnX.X"
+    AFTER_TAG = "vnX.X_txxx"
+    def upgrade(self, config, meta_config=None):
+        # Add settings
+        return config, self.reports
+"""
+class vn31_t400(MacroUpgrade):
+    """Upgrade macro for ticket #400 by Dan Copsey."""
+
+    BEFORE_TAG = "vn3.1"
+    AFTER_TAG = "vn3.1_t400"
+
+    def upgrade(self, config, meta_config=None):
+        # Commands From: rose-meta/jules-lfric
+        # Is this a coupled model? The easiest way to get this is from the number of sea ice categories.
+        nice = int(self.get_setting_value(config, ["namelist:jules_sea_seaice", "nice"])
+        is_coupled = nice > 1
+        # Add the new namelist settings
+        self.add_setting(config, ["namelist:jules_sea_seaice", "l_zenith_albedo"], ".false.")
+        if is_coupled:
+            self.add_setting(config, ["namelist:jules_sea_seaice", "meltpond_alb_vn"], "'cice'")
+        else:
+            self.add_setting(config, ["namelist:jules_sea_seaice", "meltpond_alb_vn"], "'none'")
+        self.add_setting(config, ["namelist:jules_sea_seaice", "snow_grain_size_max"], "100.0")
+        self.add_setting(config, ["namelist:jules_sea_seaice", "snow_grain_size_min"], "70.0")
+        self.add_setting(config, ["namelist:jules_sea_seaice", "snowpatch"], "0.02")
+
+
+        return config, self.reports

interfaces/jules_interface/rose-meta/jules-lfric/HEAD/rose-meta.conf

@@ -70,6 +70,83 @@ trigger=namelist:files=albedo_vis_ancil_path: .true. ;
 [namelist:jules_radiation=l_sea_alb_var_chl]
 trigger=namelist:files=sea_ancil_path: .true. ;
 
+[namelist:jules_sea_seaice=l_zenith_albedo]
+compulsory=true
+description=Use zenith angles in the calculation of sea ice
+           = (and snow on sea ice) albedos
+help=With this set to true then the sea ice albedo is adjusted in the same
+    =way as bare soil in order to adjust the albedo to take into account
+    =zenith angle.
+    =The snow on sea ice will also use the same code as snow on land in
+    =order to calculate the snow albedo which also takes into account
+    =the zenith angle.
+!kind=default
+ns=namelist/Science/JULES Surface/Sea and sea-ice
+sort-key=Panel-C12
+type=logical
+
+[namelist:jules_sea_seaice=snow_grain_size_min]
+compulsory=true
+description=Minimum snow grain size
+help=The snow albedo scheme includes a snow grain size temperature
+    =dependence. For snow on sea ice this goes from a minimum snow
+    =grain size at -30oC to a maximum snow grain size at melting point.
+    =The minimum snow grain size is set here. The larger the snow grain
+    =size the lower the albedo in the NIR part of the spectrum. The
+    =recommended value is 50 um.
+!kind=default
+ns=namelist/Science/JULES Surface/Sea and sea-ice
+sort-key=Panel-C12a
+type=real
+
+[namelist:jules_sea_seaice=snow_grain_size_max]
+compulsory=true
+description=Maximum snow grain size
+help=The snow albedo scheme includes a snow grain size temperature
+    =dependence. For snow on sea ice this goes from a minimum snow
+    =grain size at -30oC to a maximum snow grain size at melting point.
+    =The maximum snow grain size is set here. The larger the snow grain
+    =size the lower the albedo in the NIR part of the spectrum. The
+    =recommended value is 200 um.
+!kind=default
+ns=namelist/Science/JULES Surface/Sea and sea-ice
+sort-key=Panel-C12b
+type=real
+
+[namelist:jules_sea_seaice=snowpatch]
+compulsory=true
+description=snowpatch length scale (m)
+help=Length scale for parameterizing non uniform snow coverage (m).
+    =This is the depth of snow at which half of the sea ice is
+    =covered in snow. The smaller this number the more the sea ice
+    =is covered by snow.
+!kind=default
+ns=namelist/Science/JULES Surface/Sea and sea-ice
+sort-key=Panel-C12c
+type=real
+
+[namelist:jules_sea_seaice=meltpond_alb_vn]
+compulsory=true
+description=Melt pond albedo scheme
+!enumeration=true
+help=The scheme used to calculate the albedo of melt ponds
+    =none:  Use this if no melt ponds are used in the model which
+    =       instead adds a temperature dependence on albedo to
+    =       represent how melt ponds would have made the surface
+    =       darker if they existed.
+    =cice:  Use the CICE melt pond albedo scheme which uses a fixed value
+    =       (from albpondv_cice and albpondi_cice) which is linearly
+    =       ramped to with melt pond depth. No zenith angle dependence is
+    =       included.
+    =malinka: Use the Malinka et al melt pond albedo scheme which uses
+    =       Fresnel equations for reflection off the surface and multiple
+    =       scattering within the pond. This has inbuilt zenith angle
+    =       dependence and more realistic depth dependence.
+ns=namelist/Science/JULES Surface/Sea and sea-ice
+sort-key=Panel-C11
+value-titles="none","cice","malinka"
+values='none','cice','malinka'
+
 [namelist:jules_snow=can_clump]
 fail-if=len(this) != namelist:jules_surface_types=npft; # A value must be given for each PFT
        =all(this == 0) and any(namelist:jules_snow=cansnowpft ==  '.true.'); # Results in floating point exception if 0. Only used if can_model = 4 (JULES default), cansnowpft = TRUE on that tile and l_embedded_snow = TRUE (LFRic default).

interfaces/jules_interface/rose-meta/jules-lfric/versions.py

@@ -31,3 +31,28 @@ def upgrade(self, config, meta_config=None):
         # Add settings
         return config, self.reports
 """
+
+
+class vn31_t400(MacroUpgrade):
+    """Upgrade macro for ticket #400 by Dan Copsey."""
+
+    BEFORE_TAG = "vn3.1"
+    AFTER_TAG = "vn3.1_t400"
+
+    def upgrade(self, config, meta_config=None):
+
+        # Is this a coupled model? The easiest way to get this is from the number of sea ice categories.
+        nice = int(self.get_setting_value(config, ["namelist:jules_sea_seaice", "nice"]))
+        is_coupled = nice > 1
+
+        # Add the new namelist settings
+        self.add_setting(config, ["namelist:jules_sea_seaice", "l_zenith_albedo"], ".false.")
+        if is_coupled:
+            self.add_setting(config, ["namelist:jules_sea_seaice", "meltpond_alb_vn"], "'cice'")
+        else:
+            self.add_setting(config, ["namelist:jules_sea_seaice", "meltpond_alb_vn"], "'none'")
+        self.add_setting(config, ["namelist:jules_sea_seaice", "snow_grain_size_max"], "100.0")
+        self.add_setting(config, ["namelist:jules_sea_seaice", "snow_grain_size_min"], "70.0")
+        self.add_setting(config, ["namelist:jules_sea_seaice", "snowpatch"], "0.02")
+
+        return config, self.reports

interfaces/jules_interface/source/support/jules_physics_init_mod.f90

@@ -56,6 +56,14 @@ module jules_physics_init_mod
                               i_high_wind_drag_null,                           &
                               i_high_wind_drag_limited,                        &
                               i_high_wind_drag_reduced_v1,                     &
+                              meltpond_alb_vn,                                 &
+                              meltpond_alb_vn_none,                            &
+                              meltpond_alb_vn_cice,                            &
+                              meltpond_alb_vn_malinka,                         &
+                              l_zenith_albedo_in => l_zenith_albedo,           &
+                              snow_grain_size_min_in => snow_grain_size_min,   &
+                              snow_grain_size_max_in => snow_grain_size_max,   &
+                              snowpatch_in => snowpatch,                       &
                               cdn_hw_sea_in => cdn_hw_sea,                     &
                               cdn_max_sea_in => cdn_max_sea,                   &
                               u_cdn_hw_in => u_cdn_hw,                         &
@@ -212,7 +220,8 @@ subroutine jules_physics_init()
          ip_ss_coare_mq, a_chrn_coare, b_chrn_coare, u10_max_coare,         &
          l_10m_neut, alpham, dtice, l_iceformdrag_lupkes,                   &
          l_stability_lupkes, l_use_dtstar_sea, hcap_sea, beta_evap,         &
-         l_sice_meltponds, l_sice_meltponds_cice,                           &
+         l_sice_meltponds, i_meltpond_alb_vn, l_zenith_albedo,              &
+         snow_grain_size_min, snow_grain_size_max, snowpatch,               &
          l_cice_alb, l_saldep_freeze, l_sice_multilayers,                   &
          l_sice_scattering, l_sice_swpen, l_ssice_albedo,                   &
          pen_rad_frac_cice, sw_beta_cice,                                   &
@@ -371,12 +380,25 @@ subroutine jules_physics_init()
     z0h_z0m_miz          = 0.2_r_um
     z0h_z0m_sice         = 0.2_r_um
     z0sice               = 5.0e-4_r_um
+    l_zenith_albedo      = l_zenith_albedo_in
+    snow_grain_size_min  = snow_grain_size_min_in
+    snow_grain_size_max  = snow_grain_size_max_in
+    snowpatch            = snowpatch_in
+
+    ! Setup the melt pond albedo scheme
+    select case (meltpond_alb_vn)
+      case(meltpond_alb_vn_none)
+        i_meltpond_alb_vn = 0
+      case(meltpond_alb_vn_cice)
+        i_meltpond_alb_vn = 1
+      case(meltpond_alb_vn_malinka)
+        i_meltpond_alb_vn = 2
+    end select
 
     ! Setup switches that vary depending if the model is
     ! coupled to an ocean/sea-ice model or not.
     if (l_couple_sea_ice) then
       l_sice_meltponds      = .true.
-      l_sice_meltponds_cice = .true.
       l_tstar_sice_new      = .false.
       l_cice_alb            = .true.
       l_sice_multilayers    = .true.
@@ -387,7 +409,6 @@ subroutine jules_physics_init()
       sw_beta_cice          = 0.3_r_um
     else
       l_sice_meltponds      = .false.
-      l_sice_meltponds_cice = .false.
       l_tstar_sice_new      = .true.
       l_cice_alb            = .false.
       l_sice_multilayers    = .false.