initial surface flow models

wip model1

wip model1

wip bc

wip driver

wip river bc

wip river bc

added analytic function to test

both no river and analytic case 1 run

wip cleanup

wip maxwell v catchment

still wip maxwell tilted V

 added in wavespeed method for Rusanov flux

changed names, documentation, added tilted V artifact, kept old River stuff

soil tests seem to pass

overland seems to work

deleted old River.jl and test_no_river.jl files from remote

fixed unit test

removed old file, fixed doc

formatter doesnt do the right thing

move long running tilted v into its own test, run on cluster

fix boundry condition test to match new interface

gpu v catchment test fixes

format

minor doc stuff

more minor doc stuff

fixed precompute method

changed variable name from area to height, removed width

formatting

fixed doc to reflect new name

.buildkite/pipeline.yml

@@ -762,6 +762,15 @@ steps:
       queue: central
       slurm_ntasks: 1
 
+  - label: "cpu_land_overland_flow_vcatchment"
+    key: "cpu_land_overland_flow_vcatchment"
+    command:
+      - "mpiexec julia --color=yes --project test/Land/Model/test_overland_flow_vcatchment.jl"
+    agents:
+      config: cpu
+      queue: central
+      slurm_ntasks: 1
+
   - label: "gpu_thermodynamics"
     key: "gpu_thermodynamics"
     command:
@@ -2012,6 +2021,16 @@ steps:
       slurm_ntasks: 1
       slurm_gres: "gpu:1"
 
+  - label: "gpu_land_overland_flow_vcatchment"
+    key: "gpu_land_overland_flow_vcatchment"
+    command:
+      - "mpiexec julia --color=yes --project test/Land/Model/test_overland_flow_vcatchment.jl"
+    agents:
+      config: gpu
+      queue: central
+      slurm_ntasks: 1
+      slurm_gres: "gpu:1"
+
   - label: "gpu_unittests"
     key: "gpu_unittests"
     command:

docs/list_of_apis.jl

@@ -22,6 +22,9 @@ apis = [
             "APIs/Land/SoilWaterParameterizations.md",
         "Soil Heat Parameterizations" =>
             "APIs/Land/SoilHeatParameterizations.md",
+        "Surface Flow" => "APIs/Land/SurfaceFlow.md",
+        "Radiative Boundary Conditions" =>
+            "APIs/Land/RadiativeEnergyFlux.md",
     ],
     "Common" => [
         "Orientations" => "APIs/Common/Orientations.md",

docs/src/APIs/Land/SurfaceFlow.md

@@ -0,0 +1,17 @@
+# SurfaceFlow
+
+```@meta
+CurrentModule = ClimateMachine.Land.SurfaceFlow
+```
+
+## SurfaceFlow types and functions
+```@docs
+OverlandFlowModel
+NoSurfaceFlowModel
+surface_boundary_flux!
+surface_boundary_state!
+calculate_velocity
+Precip
+VolumeAdvection
+SurfaceWaterHeight
+```

src/Driver/driver_configs.jl

@@ -670,7 +670,8 @@ function MultiColumnLandModel(
     zmin = zero(FT),
     array_type = ClimateMachine.array_type(),
     mpicomm = MPI.COMM_WORLD,
-    boundary = ((3, 3), (3, 3), (1, 2)),
+    boundary = ((3, 4), (5, 6), (1, 2)),
+    solver_type = ExplicitSolverType(),
     periodicity = (false, false, false),
     meshwarp = (x...) -> identity(x),
     numerical_flux_first_order = CentralNumericalFluxFirstOrder(),

src/Land/Model/LandModel.jl

@@ -30,13 +30,14 @@ import ..BalanceLaws:
     compute_gradient_flux!,
     nodal_init_state_auxiliary!,
     init_state_prognostic!,
-    nodal_update_auxiliary_state!
+    nodal_update_auxiliary_state!,
+    wavespeed
 using ..DGMethods: LocalGeometry, DGModel
 export LandModel
 
 
 """
-    LandModel{PS, S, LBC, SRC, IS} <: BalanceLaw
+    LandModel{PS, S, SF, LBC, SRC, SRCDT, IS} <: BalanceLaw
 
 A BalanceLaw for land modeling.
 Users may over-ride prescribed default values for each field.
@@ -46,19 +47,23 @@ Users may over-ride prescribed default values for each field.
     LandModel(
         param_set,
         soil;
+        surface,
         boundary_conditions,
         source,
+        source_dt,
         init_state_prognostic
     )
 
 # Fields
 $(DocStringExtensions.FIELDS)
 """
-struct LandModel{PS, S, LBC, SRC, SRCDT, IS} <: BalanceLaw
+struct LandModel{PS, S, SF, LBC, SRC, SRCDT, IS} <: BalanceLaw
     "Parameter set"
     param_set::PS
     "Soil model"
     soil::S
+    "Surface Flow model"
+    surface::SF
     "struct of boundary conditions"
     boundary_conditions::LBC
     "Source Terms (Problem specific source terms)"
@@ -75,6 +80,7 @@ parameter_set(m::LandModel) = m.param_set
     LandModel(
         param_set::AbstractParameterSet,
         soil::BalanceLaw;
+        surface::BalanceLaw = NoSurfaceFlowModel(),
         boundary_conditions::LBC = (),
         source::SRC = (),
         init_state_prognostic::IS = nothing
@@ -85,6 +91,7 @@ Constructor for the LandModel structure.
 function LandModel(
     param_set::AbstractParameterSet,
     soil::BalanceLaw;
+    surface::BalanceLaw = NoSurfaceFlowModel(),
     boundary_conditions::LBC = LandDomainBC(),
     source::SRC = (),
     init_state_prognostic::IS = nothing,
@@ -94,6 +101,7 @@ function LandModel(
     land = (
         param_set,
         soil,
+        surface,
         boundary_conditions,
         source,
         source_dt,
@@ -106,26 +114,32 @@ end
 function vars_state(land::LandModel, st::Prognostic, FT)
     @vars begin
         soil::vars_state(land.soil, st, FT)
+        surface::vars_state(land.surface, st, FT)
     end
 end
 
 
 function vars_state(land::LandModel, st::Auxiliary, FT)
     @vars begin
+        x::FT
+        y::FT
         z::FT
         soil::vars_state(land.soil, st, FT)
+        surface::vars_state(land.surface, st, FT)
     end
 end
 
 function vars_state(land::LandModel, st::Gradient, FT)
     @vars begin
         soil::vars_state(land.soil, st, FT)
+        surface::vars_state(land.surface, st, FT)
     end
 end
 
 function vars_state(land::LandModel, st::GradientFlux, FT)
     @vars begin
         soil::vars_state(land.soil, st, FT)
+        surface::vars_state(land.surface, st, FT)
     end
 end
 
@@ -135,8 +149,11 @@ function nodal_init_state_auxiliary!(
     tmp::Vars,
     geom::LocalGeometry,
 )
+    aux.x = geom.coord[1]
+    aux.y = geom.coord[2]
     aux.z = geom.coord[3]
     land_init_aux!(land, land.soil, aux, geom)
+    land_init_aux!(land, land.surface, aux, geom)
 end
 
 function flux_first_order!(
@@ -145,8 +162,19 @@ function flux_first_order!(
     state::Vars,
     aux::Vars,
     t::Real,
-    directions,
-) end
+    direction,
+)
+    tend = Flux{FirstOrder}()
+    _args = (; state, aux, t, direction)
+    args = merge(_args, (precomputed = precompute(land, _args, tend),))
+
+    map(prognostic_vars(land)) do prog
+        var, name = get_prog_state(flux, prog)
+        val = Σfluxes(prog, eq_tends(prog, land, tend), land, args)
+        setproperty!(var, name, val)
+    end
+    nothing
+end
 
 
 function compute_gradient_argument!(
@@ -210,6 +238,7 @@ function nodal_update_auxiliary_state!(
     t::Real,
 )
     land_nodal_update_auxiliary_state!(land, land.soil, state, aux, t)
+    land_nodal_update_auxiliary_state!(land, land.surface, state, aux, t)
 end
 
 
@@ -247,7 +276,6 @@ function init_state_prognostic!(
 end
 
 include("prog_types.jl")
-
 include("RadiativeEnergyFlux.jl")
 using .RadiativeEnergyFlux
 include("SoilWaterParameterizations.jl")
@@ -260,11 +288,19 @@ include("soil_heat.jl")
 include("Runoff.jl")
 using .Runoff
 include("land_bc.jl")
+include("SurfaceFlow.jl")
+using .SurfaceFlow
 include("soil_bc.jl")
-
 include("prognostic_vars.jl")
-
-include("source.jl")
 include("land_tendencies.jl")
+include("source.jl")
+
+function wavespeed(m::LandModel, n⁻, state::Vars, aux::Vars, t::Real, direction)
+    FT = eltype(state)
+    h = max(state.surface.height, FT(0.0))
+    v = calculate_velocity(m.surface, aux.x, aux.y, h)
+    speed = norm(v)
+    return speed
+end
 
 end # Module