soil moisture and vegetation

src/dynamics/diagnostic_variables.jl

@@ -138,6 +138,7 @@ Base.@kwdef struct SurfaceVariables{NF<:AbstractFloat,Grid<:AbstractGrid{NF}}
 
     precip_large_scale::Grid = zeros(Grid,nlat_half)    # large scale precipitation (for output)
     precip_convection::Grid = zeros(Grid,nlat_half)     # convective precipitation (for output)
+    soil_moisture_availability::Grid = zeros(Grid,nlat_half)
 end
 
 # generator function based on a SpectralGrid

src/dynamics/models.jl

@@ -153,7 +153,7 @@ Base.@kwdef struct PrimitiveDryModel{NF<:AbstractFloat, D<:AbstractDevice} <: Pr
     # BOUNDARY CONDITIONS
     land_sea_mask::AbstractLandSeaMask{NF} = LandSeaMask(spectral_grid)
     ocean::AbstractOcean{NF} = SeasonalOceanClimatology(spectral_grid)
-    land::AbstractLand{NF} = SeasonalLandClimatology(spectral_grid)
+    land::AbstractLand{NF} = SeasonalLandTemperature(spectral_grid)
 
     # PHYSICS/PARAMETERIZATIONS
     physics::Bool = true
@@ -237,7 +237,9 @@ Base.@kwdef struct PrimitiveWetModel{NF<:AbstractFloat, D<:AbstractDevice} <: Pr
     # BOUNDARY CONDITIONS
     land_sea_mask::AbstractLandSeaMask{NF} = LandSeaMask(spectral_grid)
     ocean::AbstractOcean{NF} = SeasonalOceanClimatology(spectral_grid)
-    land::AbstractLand{NF} = SeasonalLandClimatology(spectral_grid)
+    land::AbstractLand{NF} = SeasonalLandTemperature(spectral_grid)
+    soil::AbstractSoil{NF} = SeasonalSoilMoisture(spectral_grid)
+    vegetation::AbstractVegetation{NF} = VegetationClimatology(spectral_grid)
 
     # PHYSICS/PARAMETERIZATIONS
     physics::Bool = true
@@ -287,6 +289,8 @@ function initialize!(model::PrimitiveWet)
     initialize!(model.land_sea_mask)
     initialize!(model.ocean)
     initialize!(model.land)
+    initialize!(model.soil)
+    initialize!(model.vegetation)
 
     # parameterizations
     initialize!(model.boundary_layer_drag,model)

src/dynamics/ocean.jl

@@ -21,6 +21,9 @@ Base.@kwdef struct SeasonalOceanClimatology{NF,Grid<:AbstractGrid{NF}} <: Abstra
     "filename of sea surface temperatures"
     file::String = "sea_surface_temperature.nc"
 
+    "Variable name in netcdf file"
+    varname::String = "sst"
+
     "Grid the sea surface temperature file comes on"
     file_Grid::Type{<:AbstractGrid} = FullGaussianGrid
 
@@ -39,33 +42,42 @@ function SeasonalOceanClimatology(SG::SpectralGrid;kwargs...)
 end
 
 function initialize!(ocean::SeasonalOceanClimatology{NF,Grid}) where {NF,Grid}
+    load_monthly_climatology!(ocean.monthly_temperature,ocean)
+end
+
+function load_monthly_climatology!( 
+    monthly::Vector{Grid},
+    scheme;
+    varname::String = scheme.varname
+) where {Grid<:AbstractGrid}
+
     # LOAD NETCDF FILE
-    if ocean.path == "SpeedyWeather.jl/input_data"
-        path = joinpath(@__DIR__,"../../input_data",ocean.file)
+    if scheme.path == "SpeedyWeather.jl/input_data"
+        path = joinpath(@__DIR__,"../../input_data",scheme.file)
     else
-        path = joinpath(ocean.path,ocean.file)
+        path = joinpath(scheme.path,scheme.file)
     end
     ncfile = NCDataset(path)
 
     # create interpolator from grid in file to grid used in model
     nx, ny = ncfile.dim["lon"], ncfile.dim["lat"]
     npoints = nx*ny
-    NF_file = typeof(ncfile["sst"].attrib["_FillValue"])
-    sst = ocean.file_Grid(zeros(NF_file,npoints))
-    interp = RingGrids.interpolator(NF,ocean.monthly_temperature[1],sst)
+    NF_file = typeof(ncfile[varname].attrib["_FillValue"])
+    grid = scheme.file_Grid(zeros(NF_file,npoints))
+    interp = RingGrids.interpolator(Float32,monthly[1],grid)
 
-    # interpolate and store in ocean
+    # interpolate and store in monthly
     for month in 1:12
-        sst_this_month = ncfile["sst"][:,:,month]
+        this_month = ncfile[varname][:,:,month]
         ij = 0
         for j in 1:ny
             for i in 1:nx
                 ij += 1
-                x = sst_this_month[i,j]
-                sst[ij] = ismissing(x) ? ocean.missing_value : x
+                x = this_month[i,j]
+                grid[ij] = ismissing(x) ? scheme.missing_value : x
             end
         end
-        interpolate!(ocean.monthly_temperature[month],sst,interp)
+        interpolate!(monthly[month],grid,interp)
     end
     return nothing
 end

src/dynamics/prognostic_variables.jl

@@ -136,11 +136,11 @@ Base.@kwdef mutable struct PrognosticVariablesLand{NF<:AbstractFloat,Grid<:Abstr
     "Snow depth [m]"
     const snow_depth::Grid = zeros(Grid,nlat_half)
 
-    "Soil wetness layer 1, volume fraction [1]"
-    const soil_wetness_layer1::Grid = zeros(Grid,nlat_half)
+    "Soil moisture layer 1, volume fraction [1]"
+    const soil_moisture_layer1::Grid = zeros(Grid,nlat_half)
 
-    "Soil wetness layer 2, volume fraction [1]"
-    const soil_wetness_layer2::Grid = zeros(Grid,nlat_half)
+    "Soil moisture layer 2, volume fraction [1]"
+    const soil_moisture_layer2::Grid = zeros(Grid,nlat_half)
 end
 
 # generator function based on a SpectralGrid

src/dynamics/time_integration.jl

@@ -221,12 +221,17 @@ function timestep!( progn::PrognosticVariables{NF}, # all prognostic variables
     model.feedback.nars_detected && return nothing  # exit immediately if NaRs already present
     (;time) = progn.clock                           # current time
 
-    ocean_timestep!(progn.ocean,time,model)
-    land_timestep!(progn.land,time,model)
-
-    # switch on/off all physics
-    model.physics && parameterization_tendencies!(diagn,progn,time,model)
-    model.physics || zero_tendencies!(diagn)        # set tendencies to zero otherwise
+    if model.physics                                # switch on/off all physics parameterizations
+        # time step ocean (temperature and TODO sea ice) and land (temperature and soil moisture)
+        ocean_timestep!(progn.ocean,time,model)
+        land_timestep!(progn.land,time,model)
+        soil_moisture_availability!(diagn.surface,progn.land,model)
+
+        # calculate all parameterizations
+        parameterization_tendencies!(diagn,progn,time,model)
+    else                                            # set tendencies to zero otherwise for accumulators
+        zero_tendencies!(diagn)
+    end       
 
     dynamics_tendencies!(diagn,progn,model,lf2)         # dynamical core
     implicit_correction!(diagn,model.implicit,progn)    # semi-implicit time stepping corrections

src/physics/column_variables.jl

@@ -40,6 +40,7 @@ function get_column!(   C::ColumnVariables,
     # TODO skin = surface approximation for now
     C.skin_temperature_sea = P.ocean.sea_surface_temperature[ij]
     C.skin_temperature_land = P.land.land_surface_temperature[ij]
+    C.soil_moisture_availability = D.surface.soil_moisture_availability[ij]
 end
 
 """Recalculate ring index if not provided."""

src/physics/define_column.jl

@@ -53,6 +53,7 @@ Base.@kwdef mutable struct ColumnVariables{NF<:AbstractFloat} <: AbstractColumnV
     surface_wind_speed::NF = 0
     skin_temperature_sea::NF = 0
     skin_temperature_land::NF = 0
+    soil_moisture_availability::NF = 0
 
     # THERMODYNAMICS
     surface_air_density::NF = 0