add ice_concentration methods, cleanup

Author: juliasloan25

docs/src/fieldexchanger.md

@@ -12,7 +12,7 @@ The specific fields that are exchanged depend on the requirements of the compone
 The fields imported from the atmosphere to the coupler are specified by extending `FieldExchanger.import_atmos_fields!`
 The default `import_atmos_fields!` imports radiative fluxes, liquid precipitation, and snow precipitation.
 
-The fields of a component model that get updated by the coupler are specified by extending `FieldExchanger.update_sim!`
+The fields of a component model that get updated by the coupler are specified by extending `FieldExchanger.update_sim!`.
 The default `update_sim!` for an atmosphere model updates the direct and diffuse surface albedos,
 the surface temperature, and the turbulent fluxes.
 The default `update_sim!` for a surface model updates the air density, radiative fluxes,
@@ -21,6 +21,10 @@ These updates are done via the `update_field!` function, which must be extended
 particular variable and component model.
 If an `update_field!` function is not defined for a particular component model, it will be ignored.
 
+Note that turbulent fluxes are not updated in `update_sim!`, but rather via
+`FluxCalculator.update_turbulent_fluxes!`, where fluxes are computed between
+the atmosphere and each surface model.
+
 ## FieldExchanger API
 
 ```@docs

docs/src/fluxcalculator.md

@@ -11,24 +11,23 @@ turbulent fluxes and ancillary quantities, such as the Obukhov length, using
 function is called at the end of each coupling step.
 
 All the quantities computed in `turbulent_fluxes!` are calculated
-separately for each surface type using the
+separately for each surface model using the
 [`FluxCalculator.compute_surface_fluxes!`](@ref) function. This function can be
 extended by component models if they need specific type of flux calculation, and
 a default is provided for models that can use the standard flux calculation.
-The final result of `turbulent_fluxes!` is an area-weighted sum of
-all the contributions of the various surfaces.
 
 The default method of [`FluxCalculator.compute_surface_fluxes!`](@ref), in turn,
 calls [`FluxCalculator.get_surface_fluxes`](@ref). This function uses a thermal
 state obtained by using the model surface temperature, extrapolates atmospheric
 density adiabatically to the surface, and with the surface humidity (if
 available, if not, assuming a saturation specific humidity for liquid phase).
-`compute_surface_fluxes!` also updates the component internal fluxes fields with
+`compute_surface_fluxes!` also updates the component internal fluxes fields via
 [`FluxCalculator.update_turbulent_fluxes!`](@ref), and adds the area-weighted
 contribution from this component model to the `CoupledSimulation` fluxes fields.
-Any extension of this function for a particular surface model is also expected
-to update both the component models' internal fluxes and the CoupledSimulation
-object's fluxes fields.
+
+Any extension of `FluxCalculator.compute_surface_fluxes!` for a particular
+surface model is also expected to update both the component models' internal
+fluxes and the CoupledSimulation object's fluxes fields.
 
 [`FluxCalculator.compute_surface_fluxes!`](@ref) sets:
 - the flux of momenta, `F_turb_ρτxz`, `F_turb_ρτyz`;

docs/src/interfacer.md

@@ -223,6 +223,7 @@ for the following properties:
 | `surface_direct albedo`  | bulk direct surface albedo                                     |         |
 | `surface_diffuse albedo` | bulk diffuse surface albedo                                    |         |
 | `surface_temperature`    | surface temperature                                            | K       |
+| `ice_concentration`      | sea ice concentration (*sea ice models only*)                  |         |
 
 !!! note
     `area_fraction` is expected to be defined on the boundary space of the simulation,

experiments/ClimaEarth/components/ocean/clima_seaice.jl

@@ -82,33 +82,25 @@ function ClimaSeaIceSimulation(land_fraction, ocean; output_dir)
     # Allocate space for the sea ice-ocean (io) fluxes
     ocean_ice_fluxes = ocean.ocean_ice_fluxes
 
-    # TODO clean this up
     # Get the initial area fraction from the fractional ice concentration
-    area_fraction = ice.model.ice_concentration
-
-    # Overwrite ice fraction with the static land area fraction anywhere we have nonzero land area
-    #  max needed to avoid Float32 errors (see issue #271; Heisenbug on HPC)
     boundary_space = axes(ocean.area_fraction)
     FT = CC.Spaces.undertype(boundary_space)
+    area_fraction = Interfacer.remap(ice.model.ice_concentration, boundary_space)
 
-    area_fraction_boundary_space = Interfacer.remap(area_fraction, boundary_space)
-    @. area_fraction_boundary_space =
-        max(min(area_fraction_boundary_space, FT(1) - land_fraction), FT(0))
+    # Overwrite ice fraction with the static land area fraction anywhere we have nonzero land area
+    #  max needed to avoid Float32 errors (see issue #271; Heisenbug on HPC)
+    @. area_fraction = max(min(area_fraction, FT(1) - land_fraction), FT(0))
 
     sim = ClimaSeaIceSimulation(
         ice,
-        area_fraction_boundary_space,
+        area_fraction,
         melting_speed,
         remapping,
         ocean_ice_fluxes,
     )
     return sim
 end
 
-function update_sic!(area_fraction, ice)
-    # TODO
-end
-
 ###############################################################################
 ### Functions required by ClimaCoupler.jl for a SurfaceModelSimulation
 ###############################################################################
@@ -146,6 +138,8 @@ Interfacer.get_field(sim::ClimaSeaIceSimulation, ::Val{:surface_temperature}) =
 Update the turbulent fluxes in the simulation using the values stored in the coupler fields.
 These include latent heat flux, sensible heat flux, momentum fluxes, and moisture flux.
 
+The input `fields` are already area-weighted, so there's no need to weight them again.
+
 Note that currently the moisture flux has no effect on the sea ice model, which has
 constant salinity.
 
@@ -202,10 +196,7 @@ function FluxCalculator.update_turbulent_fluxes!(sim::ClimaSeaIceSimulation, fie
     # Update the sea ice only where the concentration is greater than zero.
     si_flux_heat = ice_sim.ice.model.external_heat_fluxes.top
     OC.interior(si_flux_heat, :, :, 1) .=
-        OC.interior(si_flux_heat, :, :, 1) .+ (
-            (OC.interior(sim.ice.model.ice_concentration, :, :, 1) .> 0) .*
-            (remapped_F_lh .+ remapped_F_sh)
-        )
+        OC.interior(si_flux_heat, :, :, 1) .+ (remapped_F_lh .+ remapped_F_sh)
 
     return nothing
 end
@@ -245,16 +236,14 @@ function FieldExchanger.update_sim!(sim::ClimaSeaIceSimulation, csf, area_fracti
 
     # Update only the part due to radiative fluxes. For the full update, the component due
     # to latent and sensible heat is missing and will be updated in update_turbulent_fluxes.
-    # Update the sea ice only where the concentration is greater than zero.
     si_flux_heat = sim.ice.model.external_heat_fluxes.top
-    OC.interior(si_flux_heat, :, :, 1) .=
-        (OC.interior(sim.ice.model.ice_concentration, :, :, 1) .> 0) .* remapped_F_radiative
+    OC.interior(si_flux_heat, :, :, 1) .= remapped_F_radiative
 
     return nothing
 end
 
 """
-    ocean_seaice_fluxes!(ocean_sim, ice_sim)
+    ocean_seaice_fluxes!(ocean_sim::OceananigansSimulation, ice_sim::ClimaSeaIceSimulation)
 
 Compute the fluxes between the ocean and sea ice, storing them in the `ocean_ice_fluxes`
 fields of the ocean and sea ice simulations.

experiments/ClimaEarth/components/ocean/oceananigans.jl

@@ -283,7 +283,7 @@ Interfacer.get_field(sim::OceananigansSimulation, ::Val{:surface_temperature}) =
 Update the turbulent fluxes in the simulation using the values computed at this time step.
 These include latent heat flux, sensible heat flux, momentum fluxes, and moisture flux.
 
-The input fields are already area-weighted, so there's no need to weight them again.
+The input `fields` are already area-weighted, so there's no need to weight them again.
 
 A note on sign conventions:
 SurfaceFluxes and Oceananigans both use the convention that a positive flux is an upward flux.
@@ -402,7 +402,6 @@ function FieldExchanger.update_sim!(sim::OceananigansSimulation, csf, area_fract
     # Update only the part due to radiative fluxes. For the full update, the component due
     # to latent and sensible heat is missing and will be updated in update_turbulent_fluxes.
     oc_flux_T = surface_flux(sim.ocean.model.tracers.T)
-    # TODO document ordering of flux updates somewhere
     OC.interior(oc_flux_T, :, :, 1) .=
         OC.interior(oc_flux_T, :, :, 1) .+
         remapped_F_radiative ./ (ocean_reference_density * ocean_heat_capacity)

experiments/ClimaEarth/components/ocean/prescr_seaice.jl

@@ -162,8 +162,10 @@ function PrescribedIceSimulation(
 end
 
 # extensions required by Interfacer
-Interfacer.get_field(sim::PrescribedIceSimulation, ::Val{:area_fraction}) =
-    sim.integrator.p.area_fraction
+Interfacer.get_field(
+    sim::PrescribedIceSimulation,
+    ::Union{Val{:area_fraction}, Val{:ice_concentration}},
+) = sim.integrator.p.area_fraction
 Interfacer.get_field(sim::PrescribedIceSimulation, ::Val{:roughness_buoyancy}) =
     sim.integrator.p.params.z0b
 Interfacer.get_field(sim::PrescribedIceSimulation, ::Val{:roughness_momentum}) =

experiments/ClimaEarth/setup_run.jl

@@ -601,12 +601,11 @@ function CoupledSimulation(config_dict::AbstractDict)
 
         # 3. Update any fields in the model caches that can only be filled after the initial exchange.
         FieldExchanger.set_caches!(cs)
-        # TODO add a step "update area fractions" to avoid ordering between steps 3/4
 
         # 4. Compute any ocean-sea ice fluxes
         FluxCalculator.ocean_seaice_fluxes!(cs)
 
-        # 4. Calculate and update turbulent fluxes for each surface model,
+        # 5. Calculate and update turbulent fluxes for each surface model,
         #  and save the weighted average in coupler fields
         FluxCalculator.turbulent_fluxes!(cs)
     end
@@ -748,9 +747,12 @@ function step!(cs::CoupledSimulation)
     ## update the surface fractions for surface models
     FieldExchanger.update_surface_fractions!(cs)
 
-    ## exchange all non-turbulent flux fields between models
+    ## exchange all non-turbulent flux fields between models, including radiative and precipitation fluxes
     FieldExchanger.exchange!(cs)
 
+    ## compute any ocean-sea ice fluxes
+    FluxCalculator.ocean_seaice_fluxes!(cs)
+
     ## calculate turbulent fluxes in the coupler and update the model simulations with them
     FluxCalculator.turbulent_fluxes!(cs)
 

src/FieldExchanger.jl

@@ -210,11 +210,13 @@ end
 """
     update_sim!(sim::SurfaceModelSimulation, csf, area_fraction)
 
-Updates the surface component model cache with the current coupler fields besides turbulent fluxes.
+Updates the surface component model cache with the current coupler fields
+*besides turbulent fluxes*, which are updated in `update_turbulent_fluxes`.
 
 # Arguments
 - `sim`: [Interfacer.SurfaceModelSimulation] containing a surface model simulation object.
 - `csf`: [NamedTuple] containing coupler fields.
+- `area_fraction`: [CC.Fields.Field] containing the area fraction of this surface model.
 """
 function update_sim!(sim::Interfacer.SurfaceModelSimulation, csf, area_fraction)
     FT = eltype(area_fraction)

src/Interfacer.jl

@@ -219,6 +219,10 @@ get_field(
     },
 ) = get_field_error(sim, val)
 
+# Sea ice models need to provide ice concentration
+get_field(sim::SeaIceModelSimulation, val::Val{:ice_concentration}) =
+    get_field_error(sim, val)
+
 """
     get_field(sim::ComponentModelSimulation, val::Val)