Add getters and setters for variables used to override surface radiative fluxes

HISTORY.md

@@ -1,6 +1,26 @@
 History
 =======
 
+Unreleased
+----------
+
+Minor changes:
+- Added getters and setters for `Statein%adjsfcdlw_override`,
+  `Statein%adjsfcdsw_override`, and `Statein%adjsfcnsw_override`.  These
+  correspond to the time adjusted total sky downward longwave radiative flux at
+  the surface, the time adjusted total sky downward shortwave radiative flux at
+  the surface, and the time adjusted total sky net shortwave radiative flux at
+  the surface, respectively.  Note they are only available if the
+  `gfs_physics_nml.override_surface_radiative_fluxes` namelist parameter is set
+  to `.true.`.
+- Added getter and setter for `Radtend%sfalb`, the surface diffused shortwave
+  albedo.
+- Added flags for the physics timestep, `dt_atmos`, and namelist flag for
+  overriding the surface radiative fluxes, `override_surface_radiative_fluxes`,
+  to the `Flags` class.
+- Fixed a bug in the implementation of boolean flags that prevented them from
+  working properly; to date the only flag this impacted was `do_adiabatic_init`.
+
 v0.6.0 (2021-01-27)
 -------------------
 

fill_templates.py

@@ -23,13 +23,38 @@
 )
 SETUP_DIR = os.path.dirname(os.path.abspath(__file__))
 PROPERTIES_DIR = os.path.join(SETUP_DIR, "fv3gfs/wrapper")
+FORTRAN_TO_C_AND_CYTHON_TYPES = {
+    "integer": {"type_c": "c_int", "type_cython": "int"},
+    "real": {"type_c": "c_double", "type_cython": "REAL_t"},
+    "logical": {"type_c": "c_int", "type_cython": "bint"},
+}
+OVERRIDES_FOR_SURFACE_RADIATIVE_FLUXES = [
+    "override_for_time_adjusted_total_sky_downward_longwave_flux_at_surface",
+    "override_for_time_adjusted_total_sky_downward_shortwave_flux_at_surface",
+    "override_for_time_adjusted_total_sky_net_shortwave_flux_at_surface",
+]
 
 
 def get_dim_range_string(dim_list):
     token_list = [dim_ranges[dim_name] for dim_name in dim_list]
     return ", ".join(reversed(token_list))  # Fortran order is opposite of Python
 
 
+def assign_types_to_flags(flag_data):
+    flag_properties = []
+    for flag in flag_data:
+        type_fortran = flag["type_fortran"]
+        if type_fortran in FORTRAN_TO_C_AND_CYTHON_TYPES:
+            flag.update(FORTRAN_TO_C_AND_CYTHON_TYPES[type_fortran])
+        else:
+            unexpected_type = flag["type_fortran"]
+            raise NotImplementedError(
+                f"unexpected value for type_fortran: {unexpected_type}"
+            )
+        flag_properties.append(flag)
+    return flag_properties
+
+
 if __name__ == "__main__":
     requested_templates = sys.argv[1:]
 
@@ -51,7 +76,6 @@ def get_dim_range_string(dim_list):
     physics_2d_properties = []
     physics_3d_properties = []
     dynamics_properties = []
-    flagstruct_properties = []
 
     for properties in physics_data:
         if len(properties["dims"]) == 2:
@@ -65,22 +89,7 @@ def get_dim_range_string(dim_list):
         properties["dim_colons"] = ", ".join(":" for dim in properties["dims"])
         dynamics_properties.append(properties)
 
-    for flag in flagstruct_data:
-        if flag["type_fortran"] == "integer":
-            flag["type_c"] = "c_int"
-            flag["type_cython"] = "int"
-        elif flag["type_fortran"] == "real":
-            flag["type_c"] = "c_double"
-            flag["type_cython"] = "REAL_t"
-        elif flag["type_fortran"] == "logical":
-            flag["type_c"] = "c_bool"
-            flag["type_cython"] = "bint"
-        else:
-            unexpected_type = flag["type_fortran"]
-            raise NotImplementedError(
-                f"unexpected value for type_fortran: {unexpected_type}"
-            )
-        flagstruct_properties.append(flag)
+    flagstruct_properties = assign_types_to_flags(flagstruct_data)
 
     if len(requested_templates) == 0:
         requested_templates = all_templates
@@ -94,6 +103,7 @@ def get_dim_range_string(dim_list):
             physics_3d_properties=physics_3d_properties,
             dynamics_properties=dynamics_properties,
             flagstruct_properties=flagstruct_properties,
+            overriding_fluxes=OVERRIDES_FOR_SURFACE_RADIATIVE_FLUXES,
         )
         with open(out_filename, "w") as f:
             f.write(result)

fv3gfs/wrapper/_properties.py

@@ -9,6 +9,15 @@
 with open(os.path.join(DIR, "physics_properties.json"), "r") as f:
     PHYSICS_PROPERTIES = json.load(f)
 
+with open(os.path.join(DIR, "flagstruct_properties.json"), "r") as f:
+    FLAGSTRUCT_PROPERTIES = json.load(f)
+
+OVERRIDES_FOR_SURFACE_RADIATIVE_FLUXES = [
+    "override_for_time_adjusted_total_sky_downward_longwave_flux_at_surface",
+    "override_for_time_adjusted_total_sky_downward_shortwave_flux_at_surface",
+    "override_for_time_adjusted_total_sky_net_shortwave_flux_at_surface",
+]
+
 DIM_NAMES = {
     properties["name"]: properties["dims"]
     for properties in DYNAMICS_PROPERTIES + PHYSICS_PROPERTIES

fv3gfs/wrapper/_restart/io.py

@@ -1,13 +1,17 @@
 from .._wrapper import get_tracer_metadata
-from .._properties import DYNAMICS_PROPERTIES, PHYSICS_PROPERTIES
+from .._properties import (
+    DYNAMICS_PROPERTIES,
+    PHYSICS_PROPERTIES,
+    OVERRIDES_FOR_SURFACE_RADIATIVE_FLUXES,
+)
 
 # these variables are found not to be needed for smooth restarts
 # later we could represent this as a key in the dynamics/physics properties
 RESTART_EXCLUDE_NAMES = [
     "convective_cloud_fraction",
     "convective_cloud_top_pressure",
     "convective_cloud_bottom_pressure",
-]
+] + OVERRIDES_FOR_SURFACE_RADIATIVE_FLUXES
 
 
 def get_restart_names():

fv3gfs/wrapper/flagstruct_properties.json

@@ -28,5 +28,11 @@
         "fortran_name" : "do_adiabatic_init",
         "location" : "do_adiabatic_init",
         "type_fortran": "logical"
+    },
+    {
+        "name": "override_surface_radiative_fluxes",
+        "fortran_name" : "override_surface_radiative_fluxes",
+        "location" : "IPD_Control",
+        "type_fortran": "logical"
     }
 ]

fv3gfs/wrapper/physics_properties.json

@@ -112,6 +112,13 @@
         "container": "Radtend",
         "dims": ["y", "x"]
     },
+    {
+        "name": "surface_diffused_shortwave_albedo",
+        "fortran_name": "sfalb",
+        "units": "",
+        "container": "Radtend",
+        "dims": ["y", "x"]
+    },
     {
         "name": "total_sky_downward_shortwave_flux_at_top_of_atmosphere",
         "fortran_name": "topfsw",
@@ -455,5 +462,26 @@
         "units": "unknown",
         "container": "Sfcprop",
         "dims": ["z_soil", "y", "x"]
+    },
+    {
+        "name": "override_for_time_adjusted_total_sky_downward_longwave_flux_at_surface",
+        "fortran_name": "adjsfcdlw_override",
+        "units": "W/m^2",
+        "container": "Statein",
+        "dims": ["y", "x"]
+    },
+    {
+        "name": "override_for_time_adjusted_total_sky_downward_shortwave_flux_at_surface",
+        "fortran_name": "adjsfcdsw_override",
+        "units": "W/m^2",
+        "container": "Statein",
+        "dims": ["y", "x"]
+    },
+    {
+        "name": "override_for_time_adjusted_total_sky_net_shortwave_flux_at_surface",
+        "fortran_name": "adjsfcnsw_override",
+        "units": "W/m^2",
+        "container": "Statein",
+        "dims": ["y", "x"]
     }
 ]

templates/_wrapper.pyx

@@ -199,8 +199,16 @@ cdef int set_2d_quantity(name, REAL_t[:, ::1] array) except -1:
     {% endif %}
 {% endfor %}
 {% for item in physics_2d_properties %}
+    {% if item.name in overriding_fluxes %}
+    elif name == '{{ item.name }}':
+        if flags.override_surface_radiative_fluxes:
+            set_{{ item.fortran_name }}{% if "fortran_subname" in item %}_{{ item.fortran_subname }}{% endif %}(&array[0, 0])
+        else:
+            raise fv3gfs.util.InvalidQuantityError('Overriding surface fluxes can only be set if gfs_physics_nml.override_surface_radiative_fluxes is set to .true.')
+    {% else %}
     elif name == '{{ item.name }}':
         set_{{ item.fortran_name }}{% if "fortran_subname" in item %}_{{ item.fortran_subname }}{% endif %}(&array[0, 0])
+    {% endif %}
 {% endfor %}
     else:
         raise ValueError(f'no setter available for {name}')
@@ -261,10 +269,20 @@ def get_state(names, dict state=None, allocator=None):
         state['time'] = get_time()
 
 {% for item in physics_2d_properties %}
+    {% if item.name in overriding_fluxes %}
+    if '{{ item.name }}' in input_names_set:
+        if flags.override_surface_radiative_fluxes:
+            quantity = _get_quantity(state, "{{ item.name }}", allocator, {{ item.dims | safe }}, "{{ item.units }}", dtype=real_type)
+            with fv3gfs.util.recv_buffer(quantity.np.empty, quantity.view[:]) as array_2d:
+                get_{{ item.fortran_name }}{% if "fortran_subname" in item %}_{{ item.fortran_subname }}{% endif %}(&array_2d[0, 0])
+        else:
+            raise fv3gfs.util.InvalidQuantityError('Overriding surface fluxes can only be accessed if gfs_physics_nml.override_surface_radiative_fluxes is set to .true.')
+    {% else %}
     if '{{ item.name }}' in input_names_set:
         quantity = _get_quantity(state, "{{ item.name }}", allocator, {{ item.dims | safe }}, "{{ item.units }}", dtype=real_type)
         with fv3gfs.util.recv_buffer(quantity.np.empty, quantity.view[:]) as array_2d:
             get_{{ item.fortran_name }}{% if "fortran_subname" in item %}_{{ item.fortran_subname }}{% endif %}(&array_2d[0, 0])
+    {% endif %}
 {% endfor %}
 
 {% for item in physics_3d_properties %}
@@ -372,6 +390,11 @@ class Flags:
         get_{{item.fortran_name}}(&{{item.name}})
         return {{item.name}}
 {% endfor %}
+    @property
+    def dt_atmos(self):
+        cdef int dt_atmos
+        get_physics_timestep_subroutine(&dt_atmos)
+        return dt_atmos
 
 flags = Flags()
 

templates/flagstruct_data.F90

@@ -1,6 +1,7 @@
 module flagstruct_data_mod
 
 use atmosphere_mod, only: Atm, mytile
+use atmos_model_mod, only: IPD_Control
 use fv_nwp_nudge_mod, only: do_adiabatic_init
 use iso_c_binding
 
@@ -11,7 +12,7 @@ module flagstruct_data_mod
 {% for item in flagstruct_properties %}
     {% if item.fortran_name == "do_adiabatic_init" %}
     subroutine get_do_adiabatic_init(do_adiabatic_init_out) bind(c)
-        logical(c_bool), intent(out) :: do_adiabatic_init_out
+        logical(c_int), intent(out) :: do_adiabatic_init_out
         do_adiabatic_init_out = do_adiabatic_init
     end subroutine get_do_adiabatic_init
     {% else %}
@@ -21,6 +22,8 @@ subroutine get_{{ item.fortran_name }}({{ item.fortran_name }}_out) bind(c)
         {{ item.fortran_name }}_out = Atm(mytile)%flagstruct%{{ item.fortran_name }}
         {% elif item.location == "Atm" %}
         {{ item.fortran_name }}_out = Atm(mytile)%{{ item.fortran_name }}
+        {% elif item.location == "IPD_Control" %}
+        {{ item.fortran_name }}_out = IPD_Control%{{ item.fortran_name }}
         {% endif %}
     end subroutine get_{{ item.fortran_name }}
     {% endif %}

tests/test_all.py

@@ -1,9 +1,7 @@
 import unittest
 from mpi4py import MPI
 from util import run_unittest_script
-import os
 
-base_dir = os.path.dirname(os.path.realpath(__file__))
 
 # The packages we import will import MPI, causing an MPI init, but we don't actually
 # want to use MPI under this script. We have to finalize so mpirun will work on
@@ -13,25 +11,34 @@
 
 class UsingMPITests(unittest.TestCase):
     def test_getters(self):
-        run_unittest_script(os.path.join(base_dir, "test_getters.py"))
+        run_unittest_script("test_getters.py")
 
-    def test_setters(self):
-        run_unittest_script(os.path.join(base_dir, "test_setters.py"))
+    def test_setters_default(self):
+        run_unittest_script("test_setters.py", "false")
+
+    def test_setters_while_overriding_surface_radiative_fluxes(self):
+        run_unittest_script("test_setters.py", "true")
+
+    def test_overrides_for_surface_radiative_fluxes_modify_diagnostics(self):
+        run_unittest_script("test_overrides_for_surface_radiative_fluxes.py")
 
     def test_diagnostics(self):
-        run_unittest_script(os.path.join(base_dir, "test_diagnostics.py"))
+        run_unittest_script("test_diagnostics.py")
 
     def test_tracer_metadata(self):
-        run_unittest_script(os.path.join(base_dir, "test_tracer_metadata.py"))
+        run_unittest_script("test_tracer_metadata.py")
 
     def test_get_time_julian(self):
-        run_unittest_script(os.path.join(base_dir, "test_get_time.py"), "julian")
+        run_unittest_script("test_get_time.py", "julian")
 
     def test_get_time_thirty_day(self):
-        run_unittest_script(os.path.join(base_dir, "test_get_time.py"), "thirty_day")
+        run_unittest_script("test_get_time.py", "thirty_day")
 
     def test_get_time_noleap(self):
-        run_unittest_script(os.path.join(base_dir, "test_get_time.py"), "noleap")
+        run_unittest_script("test_get_time.py", "noleap")
+
+    def test_flags(self):
+        run_unittest_script("test_flags.py")
 
 
 if __name__ == "__main__":

tests/test_diagnostics.py

@@ -4,7 +4,7 @@
 import fv3gfs.util
 from mpi4py import MPI
 
-from util import main
+from util import get_default_config, main
 
 test_dir = os.path.dirname(os.path.abspath(__file__))
 
@@ -48,4 +48,5 @@ def test_get_diagnostic_data(self):
 
 
 if __name__ == "__main__":
-    main(test_dir)
+    config = get_default_config()
+    main(test_dir, config)