Omega/add coupling fluxes

components/omega/configs/Default.yml

@@ -30,9 +30,9 @@ Omega:
     HorzTracerFluxOrder: 2
     VerticalTracerFluxLimiterEnable: true
     VerticalTracerFluxOrder: 3
-  WindStress:
+  SrfStress:
     InterpType: Isotropic
-  SurfaceRestoring:
+  SrfRestoring:
     TracersToRestore: [Temperature, Salinity]
     PistonVelocity: 1.585e-5
   PressureGrad:
@@ -47,10 +47,12 @@ Omega:
     VelHyperDiffTendencyEnable: true
     ViscDel4: 1.2e11
     DivFactor: 1.0
-    WindForcingTendencyEnable: false
-    SurfaceTracerRestoringEnable: false
+    SrfStressForcingTendencyEnable: false
+    SrfTracerRestoringEnable: false
     BottomDragTendencyEnable: false
     BottomDragCoeff: 0.0
+    SrfThicknessForcingTendencyEnable: false
+    SrfTracerForcingTendencyEnable: false
     TracerHorzAdvTendencyEnable: true
     TracerDiffTendencyEnable: true
     EddyDiff2: 10.0
@@ -125,6 +127,16 @@ Omega:
       EndTime: 99999-12-31_00:00:00
       Contents:
         - Restart
+    Forcing:
+      UsePointerFile: false
+      Filename: forcing.nc
+      Mode: read
+      Precision: double
+      Freq: 1
+      FreqUnits: days
+      UseStartEnd: false
+      Contents:
+        - Forcing
     RestartWrite:
       UsePointerFile: true
       PointerFilename: ocn.pointer

components/omega/doc/design/Tendencies.md

@@ -50,7 +50,7 @@ class Tendencies{
    SSHGradOnEdge SSHGrad;
    VelocityDiffusionOnEdge VelocityDiffusion;
    VelocityHyperDiffOnEdge VelocityHyperDiff;
-   WindForcingOnEdge WindForcing;
+   SrfStressForcingOnEdge SrfStressForcing;
    BottomDragOnEdge BottomDrag;
    TracerDiffOnCell TracerDiffusion;
    TracerHyperDiffOnCell TracerHyperDiff;

components/omega/doc/devGuide/AuxiliaryVariables.md

@@ -101,7 +101,7 @@ The following auxiliary variable groups are currently implemented:
 ||  Del2RelVortVertex ||
 | TracerAuxVars | HTracersEdge | Center or Upwind|
 || Del2TracersCell ||
-| WindForcingAuxVars |  ZonalStressCell ||
+| MomForcingAuxVars |  ZonalStressCell ||
 ||  MeridStressCell ||
 ||  NormalStressEdge ||
 | SurfTracerRestAuxVars | TracersMonthlySurfClimoCell ||

components/omega/doc/devGuide/EOS.md

@@ -70,6 +70,14 @@ volume arrays, do
 Eos.computeBruntVaisalaFreqSq(ConservTemp, AbsSalinity, Pressure, SpecVol);
 ```
 
+## Helper functions for conversion
+
+To provide the surface in-situ temperature that the coupler requires, use the helper function that converts conservative temperature (the state variable when using Teos-10) into potential temperature, given the absolute salinity and pressure arrays:
+
+```c++
+Eos.calcPtFromCt(ConsrvTemp, AbsSalinity, Pressure);
+```
+
 ## Removal of Eos
 
 To clear the Eos instance do:

components/omega/doc/devGuide/Forcing.md

@@ -5,67 +5,120 @@
 This page describes design and implementation details for forcing-related
 pathways in Omega, currently this includes:
 
-- Wind forcing
+- Surface momemtum forcing (e.g. stress from wind forcing)
+- Surface thickness and tracer forcing
 - Surface tracer restoring
 
-## Wind forcing design
+## Surface momentum forcing design
 
-### Wind forcing data flow
+### Surface stress forcing data flow
 
 1. External fields provide:
    - `WindStressZonal`
    - `WindStressMeridional`
 2. Auxiliary-state compute builds `NormalStressEdge` from those fields.
-3. Tendency term applies wind-stress forcing to edge-normal velocity tendency.
+3. Tendency term applies stress forcing to edge-normal velocity tendency.
 
-### Wind forcing key classes/components
+### Surface momemtum forcing key classes/components
 
-- `WindForcingAuxVars`
-  - Stores wind-stress cell fields and computed `NormalStressEdge`
+- `MomForcingAuxVars`
+  - Stores surface stress cell fields and computed `NormalStressEdge`
   - Applies configured interpolation choice (`InterpType`)
 - `AuxiliaryState::computeMomAux`
-  - Calls `WindForcingAuxVars::computeVarsOnEdge`
-- `WindForcingOnEdge` tendency term
+  - Calls `MomForcingAuxVars::computeVarsOnEdge`
+- `SrfStressForcingOnEdge` tendency term
   - Adds contribution proportional to normal stress and inverse layer
     thickness in the surface layer
 
-### Wind forcing config coupling
+### Surface stress forcing config coupling
 
-- `Omega.WindStress.InterpType`
+- `Omega.SrfStress.InterpType`
   - mapped to `InterpCellToEdgeOption`
-- `Omega.Tendencies.WindForcingTendencyEnable`
-  - gates execution of wind forcing tendency kernel
+- `Omega.Tendencies.SrfStressForcingTendencyEnable`
+  - gates execution of surface stress (e.g. wind forcing) tendency kernel
+
+## Surface flux forcing design
+
+### Surface flux forcing data flow
+
+**Thickness equation pathway:**
+
+1. External fields provide freshwater and salt flux components:
+   - `SnowFlux`, `RainFlux`, `EvaporationFlux`
+   - `SeaIceFreshWaterFlux`, `IceRunoffFlux`, `RiverRunoffFlux`
+   - `SeaIceSaltFlux`
+2. `Forcing` stores the flux fields in `TracerForcingAuxVars`
+3. The tendency term `SrfThicknessForcing`sums both the freshwater and salt mass fluxes, converted to be applied to
+the surface layer pseudo-thickness.
+
+**Tracer equation pathway:**
+
+1. External fields provide heat and salt flux components:
+   - `LatentHeatFlux`, `SensibleHeatFlux`
+   - `LongWaveHeatFluxUp`, `LongWaveHeatFluxDown`
+   - `SeaIceHeatFlux`, `ShortWaveHeatFlux`
+   - `SeaIceSaltFlux`, `SnowFlux`, `IceRunoffFlux`
+2. `Forcing` stores the flux fields in `TracerForcingAuxVars`
+3. The tendency terms for each tracer `SrfTracerForcing` converts external heat fluxes to tendencies in conservative temperature,
+ and external mass salt flux to salinity (g/kg) in the surface layer.
+
+### Coupled flux forcing key classes/components
+
+- `TracerForcingAuxVars`
+  - Stores 13 coupled flux cell-centered fields: 7 freshwater fluxes and 6 heat
+    fluxes, plus 1 salt flux component
+  - Fields initialized to zero and registered in `Forcing` field group
+- `SrfThicknessForcingOnCell` tendency term
+  - Computes freshwater flux contribution: $\sum (\text{SnowFlux} + \text{RainFlux} + \text{EvaporationFlux} + \text{SeaIceFreshWaterFlux} + \text{IceRunoffFlux} + \text{RiverRunoffFlux} + \text{SeaIceSaltFlux}) / \rho_{sw}$
+  - Applied only at surface layer (top active layer) using `MinLayerCell`
+- `SrfTracerForcingOnCell` tendency term
+  - For temperature: computes heat flux minus latent heat of fusion contribution: $(\sum \text{HeatFluxes} - (\text{SnowFlux} + \text{IceRunoffFlux}) L_i) \times H_{\text{FluxFac}}$
+  - For salinity: applies salt flux with unit conversion: $\text{SeaIceSaltFlux} \times S_{\text{FluxFac}}$
+  - Applied only at surface layer using `MinLayerCell`
+  - Uses tracer index validation to apply to specific tracers only
+- `Forcing`
+  - Manages `TracerForcingAuxVars` instance
+- `Tendencies`
+  - Calls `SrfThicknessForcingOnCell` in `computeThicknessTendenciesOnly`
+  - Calls `SrfTracerForcingOnCell` in `computeTracerTendenciesOnly` after surface tracer restoring
+
+### Coupled flux forcing config coupling
+
+- `Omega.Tendencies.SrfThicknessForcingTendencyEnable`
+  - gates execution of coupled flux thickness kernel
+  - controls freshwater and salt flux forcing on sea surface height
+- `Omega.Tendencies.SrfTracerForcingTendencyEnable`
+  - gates execution of coupled flux tracer kernel
+  - controls heat flux forcing on temperature and salt flux forcing on salinity
 
 ## Surface tracer restoring design
 
 ### Surface tracer restoring data flow
 
 1. External fields provide target values: `TracersMonthlySurfClimoCell` (values and units should match the state variables)
-2. Auxiliary-state compute forms restoring differences: `SurfTracerRestoringDiffsCell = target - tracer_surface`
-3. Tendency term applies restoring only at surface layer and only for tracers selected from `SurfaceRestoring.TracersToRestore`.
+2. `SurfTracerRestAuxVars` stores `TracersMonthlySurfClimoCell` for later restoring use.
+3. Tendency term applies restoring only at surface layer and only for tracers selected from `SrfRestoring.TracersToRestore`.
 
 ### Surface tracer restoring key classes/components
 
 - `SurfTracerRestAuxVars`
-  - Inputs: `TracersMonthlySurfClimoCell`, tracer state array
-  - Output: `SurfTracerRestoringDiffsCell`
-  - Uses `MinLayerCell` to select surface layer index
-- `SurfaceTracerRestoringOnCell` tendency term
-  - Applies `PistonVelocity * SurfTracerRestoringDiffsCell` at surface
+  - Stores `TracersMonthlySurfClimoCell`
+- `SrfTracerRestoringOnCell` tendency term
+  - Applies `PistonVelocity * (target - state)` at the surface layer
 - `Tendencies`
-  - Parses `SurfaceRestoring.TracersToRestore` and resolves tracer indices
+  - Parses `SrfRestoring.TracersToRestore` and resolves tracer indices
   - Builds `TracerIdsToRestore` and `NTracersToRestore`
   - Applies tracer-selection logic at call site in
     `computeTracerTendenciesOnly`
   - Aborts if restoring is enabled but no tracer IDs are available
 
 ### Surface tracer restoring config coupling
 
-- `Omega.SurfaceRestoring.PistonVelocity`
+- `Omega.SrfRestoring.PistonVelocity`
   - tendency scaling
-- `Omega.SurfaceRestoring.TracersToRestore`
+- `Omega.SrfRestoring.TracersToRestore`
   - tracer-level enable list used to build `TracerIdsToRestore`
-- `Omega.Tendencies.SurfaceTracerRestoringEnable`
+- `Omega.Tendencies.SrfTracerRestoringEnable`
   - gates restoring tendency execution
 
 ## Notes

components/omega/doc/devGuide/TendencyTerms.md

@@ -35,13 +35,13 @@ implemented:
 - `SSHGradOnEdge`
 - `VelocityDiffusionOnEdge`
 - `VelocityHyperDiffOnEdge`
-- `WindForcingOnEdge`
+- `SrfStressForcingOnEdge`
 - `BottomDragOnEdge`
 - `TracerHorzAdvOnCell`
 - `TracerHighOrderHorzAdvOnCell`
 - `TracerDiffOnCell`
 - `TracerHyperDiffOnCell`
-- `SurfaceTracerRestoringOnCell`
+- `SrfTracerRestoringOnCell`
 
 ## See Also
 

components/omega/doc/userGuide/AuxiliaryVariables.md

@@ -34,6 +34,7 @@ The following auxiliary variables are currently available:
 | NormalStressEdge | normal component of wind stress on edge
 | SurfTracerRestoringDiffsCell | surface tracer restoring differences on cells
 | TracersMonthlySurfClimoCell | monthly climatology values to restore to for surface tracer on cells
+| SurfInsituTemperature | insitu (potential) temperature at surface layer in Celsius
 
 ## See Also
 

components/omega/doc/userGuide/Forcing.md

@@ -5,6 +5,7 @@
 This page documents the user-facing configuration and behavior for current forcing in Omega:
 
 - Wind forcing
+- Coupled flux forcing
 - Surface tracer restoring
 
 ## Wind forcing
@@ -17,17 +18,17 @@ Wind forcing behavior is controlled by two configuration blocks:
 
 ```yaml
 Omega:
-  WindStress:
+  SrfStress:
     InterpType: Isotropic
 
   Tendencies:
-    WindForcingTendencyEnable: true
+    SrfStressForcingTendencyEnable: true
 ```
 
-- `WindStress.InterpType`
+- `SrfStress.InterpType`
   - `Isotropic`: isotropic cell-to-edge interpolation for wind stress
   - `Anisotropic`: anisotropic interpolation option
-- `Tendencies.WindForcingTendencyEnable`: switch to enable wind forcing tendency
+- `Tendencies.SrfStressForcingTendencyEnable`: switch to enable wind forcing tendency
 
 ### Required input fields
 
@@ -39,6 +40,71 @@ Wind forcing uses auxiliary wind-stress fields:
 These are used to form edge-normal stress (`NormalStressEdge`) that enters
 momentum tendencies.
 
+## Surface flux forcing
+
+Surface flux forcing applies ocean-atmosphere and ocean-sea ice fluxes from the other model
+components (atmosphere, sea ice) to the thickness and tracer equations. This enables
+the ocean to respond to heat, freshwater, and salt exchanges at the surface. These fluxes can be from data or (active) coupled components.
+
+### Surface flux forcing configuration
+
+Surface flux forcing is controlled by two configuration flags:
+
+```yaml
+Omega:
+  Tendencies:
+    SrfThicknessForcingTendencyEnable: false
+    SrfTracerForcingTendencyEnable: false
+```
+
+- `Tendencies.SrfThicknessForcingTendencyEnable`: enables coupled freshwater and salt flux forcing on thickness
+- `Tendencies.SrfTracerForcingTendencyEnable`: enables coupled heat and salt flux forcing on tracers
+
+### Required input fields
+
+Coupled flux forcing uses 13 auxiliary fields organized by type:
+
+**Freshwater mass fluxes (kg m⁻² s⁻¹):**
+- `SnowFlux`: precipitation from snow
+- `RainFlux`: precipitation from rain
+- `EvaporationFlux`: evaporative water loss
+- `SeaIceFreshWaterFlux`: freshwater input from sea-ice melt or formation
+- `IceRunoffFlux`: runoff from land ice
+- `RiverRunoffFlux`: runoff from rivers
+
+**Heat fluxes (W m⁻²):**
+- `LatentHeatFlux`: latent heat transfer
+- `SensibleHeatFlux`: sensible heat transfer
+- `LongWaveHeatFluxUp`: upward longwave radiation
+- `LongWaveHeatFluxDown`: downward longwave radiation
+- `SeaIceHeatFlux`: heat from sea-ice interaction
+- `ShortWaveHeatFlux`: shortwave (solar) radiation
+
+**Salt mass flux (kg m⁻² s⁻¹):**
+- `SeaIceSaltFlux`: salt flux from sea-ice formation/melt processes
+
+These fields are populated by external coupling components (typically atmosphere
+and ice models). Omega assumes the incoming values match the documented units.
+For now, there are assumed to come from a `forcing.nc` file, but later will be provided
+by the equivalent `ocn_comp_mct.F`.
+
+### Notes
+
+- Coupled fluxes are applied only at the surface layer (top active layer) for each cell.
+- Temperature tendency is computed as the sum of all heat fluxes minus latent heat
+  of fusion for snow and ice runoff, converted to temperature tendency via
+  $H_{\text{FluxFac}} = 1.0 / (\rho_{sw} c^0_{p,sw})$ where $c^0_{p,sw}$ is the reference
+  specific heat of seawater defined by TEOS-10. This allows the conversion to the
+  conservative temperature variable.
+- Salinity tendency from salt flux is scaled by $S_{\text{FluxFac}} = 1.0e3 / \rho_{sw}$
+  to account for unit conversion from kg/(m²·s) to salinity units (g/kg).
+- Fluxes are assumed to be in the documented units (i.e. net mass fluxes);
+  any unit conversion should be performed by the coupling component before providing flux
+  values to Omega.
+- The reference density used here ($\rho_{sw}$) is not a Boussinesq density, it is the
+  conversion factor from mass to pseudo-thickness.
+- No iceberg fluxes are included for now.
+
 ## Surface tracer restoring
 
 Surface tracer restoring applies a piston-velocity tendency, or damping, at the ocean
@@ -53,17 +119,17 @@ Surface tracer restoring is controlled by two configuration blocks:
 
 ```yaml
 Omega:
-  SurfaceRestoring:
+  SrfRestoring:
     TracersToRestore: [Temperature, Salinity]
     PistonVelocity: 1.585e-5
 
   Tendencies:
-    SurfaceTracerRestoringEnable: true
+    SrfTracerRestoringEnable: true
 ```
 
 - `TracersToRestore`: list of tracer names that restoring is applied to
 - `PistonVelocity`: restoring rate coefficient
-- `SurfaceTracerRestoringEnable`: switch to enable surface tracer restoring
+- `SrfTracerRestoringEnable`: switch to enable surface tracer restoring
 
 When restoring is enabled, Omega resolves `TracersToRestore` into an internal
 list of tracer IDs and applies restoring only to tracers in that list.
@@ -73,10 +139,9 @@ list of tracer IDs and applies restoring only to tracers in that list.
 Surface restoring uses auxiliary fields:
 
 - `TracersMonthlySurfClimoCell`: restoring target climatological values
-- `SurfTracerRestoringDiffsCell`: computed target-minus-state differences
 
-The restoring tendency is computed at the surface layer only and is limited by
-the configured `PistonVelocity` and target-minus-state difference.
+The restoring tendency `SrfTracerRestoring` is computed at the surface layer only using the
+configured `PistonVelocity` and the inline target-minus-state difference.
 
 ## Notes