Wetting and Drying

Coastal flooding simulations require robust handling of moving shorelines where elements transition between wet and dry states [Medeiros2012] [Carrier1958]. Cocoa implements an element-based wetting/drying algorithm that maintains mass conservation and numerical stability.

Algorithm Overview

digraph wetdry { rankdir=TB; node [shape=box, style="filled,rounded", fontname="Helvetica", fontsize=10]; edge [color="#555555"]; bgcolor="transparent"; init [label="1. Initialize Cycle\nCopy committed -> working status", fillcolor="#E0E0E0", color="#616161"]; d1 [label="2. Nodal Drying (D1)\nH < H0 -> Dry", fillcolor="#FFEBEE", color="#E53935"]; w1 [label="3. Nodal Wetting (W1)\nWet neighbor + flow potential\n-> Wet", fillcolor="#E3F2FD", color="#1E88E5"]; barrier [label="4. Barrier Wetting\nForce-wet barrier neighbors", fillcolor="#E3F2FD", color="#1E88E5"]; de1 [label="5. Elemental Drying (DE1)\nShallow control section -> Inactive", fillcolor="#FFEBEE", color="#E53935"]; count [label="6. Active Element Count\nPer-node wet element tally", fillcolor="#FFF3E0", color="#FB8C00"]; d2 [label="7. Landlocked Drying (D2)\nZero active elements -> Dry", fillcolor="#FFEBEE", color="#E53935"]; commit [label="8. Status Reconciliation\nWorking -> committed", fillcolor="#E8F5E9", color="#43A047"]; alpha [label="9. Slope Limiter\nCompute pressure gradient\nlimiting factors", fillcolor="#F3E5F5", color="#8E24AA"]; init -> d1 -> w1 -> barrier -> de1 -> count -> d2 -> commit -> alpha; }

Fig. 32 Wet/dry algorithm phases executed after each GWCE solve

Cocoa uses a 9-phase wet/dry algorithm executed after each GWCE solve:

  1. Initialize Cycle: Copy committed node status to working status

  2. Nodal Drying D1: Mark nodes dry if depth falls below threshold

  3. Nodal Wetting W1: Mark dry nodes wet if neighbors can supply water

  4. Barrier Wetting: Force-wet dry neighbors of barrier (NIBNODECODE) nodes

  5. Elemental Drying DE1: Deactivate elements acting as shallow control sections (see Phase Details)

  6. Active Element Count: Count wet elements per node

  7. Landlocked Drying D2: Dry nodes with no active elements

  8. Status Reconciliation: Commit working status to final status

  9. Slope Limiter (ALPHA): Compute pressure gradient limiting factors

Node Status

Each node has a wet/dry status stored in WetDryData:

enum class NodeStatus : int {
  Dry = 0,
  Wet = 1
};

Two status arrays are maintained:

  • m_node_wet_status: Committed status used in computations

  • m_node_wet_working: Working status during algorithm phases

Element Status

Elements have an active/inactive status:

enum class ElementStatus : int {
  Inactive = 0,
  Active = 1
};

Inactive elements are excluded from GWCE and momentum assembly.

Depth Thresholds

The algorithm uses several depth thresholds, which can be configured via the physics section of the configuration file (see Configuration):

Parameter

Default

Config Key

Description

\(H_0\)

0.1 m

h0

Minimum depth for wet classification

\(H_{off}\)

0.12 m

(derived)

Wetting threshold (\(1.2 \times H_0\))

\(H_{abs}\)

0.08 m

(derived)

Absolute floor (\(0.8 \times H_0\))

\(V_{min}\)

0.01 m/s

velmin

Minimum velocity for wetting

Number line showing wet/dry depth thresholds and hysteresis

Fig. 33 Wet/dry depth thresholds. Nodes dry when \(H \le H_0\) and re-wet when \(H > H_{off}\), creating hysteresis that prevents oscillation at the shoreline.

Phase Details

Nodal Drying (D1)

A node becomes dry when its total depth falls below the threshold:

\[H_i = h_i + \zeta_i \leq H_0 \implies \text{Dry}\]

Nodal Wetting (W1)

A dry node becomes wet when:

  1. It has at least one wet neighbor

  2. Water can flow from wet to dry (potential gradient favorable)

  3. Estimated velocity exceeds minimum:

\[V = \sqrt{g H \frac{\Delta \zeta}{\Delta x}} > V_{min}\]

When a node wets, the TKM friction tensor is updated with TKNF scaling:

\[\tau = \min\left(\frac{TKNF \cdot C_d |\mathbf{V}|}{H}, \frac{2}{\Delta t}\right)\]

Elemental Drying (DE1)

DE1 deactivates elements that act as shallow control sections. Ordering the element’s nodal elevations \(\zeta_A \geq \zeta_B > \zeta_C\), the element is deactivated when either of the two higher-elevation nodes is shallower than the wetting threshold:

\[H_A < H_{off} \quad \text{or} \quad H_B < H_{off} \implies \text{Inactive}\]

Elements touching a node with a barrier contribution (NIBCNT > 0) are exempt, so weir overflow can keep its receiving elements active.

Active Element Count

For each node, count elements that are:

  • Active (not deactivated by DE1), AND

  • Have all three nodes wet

This determines if a node is “landlocked.”

Landlocked Drying (D2)

A node with zero active elements is marked dry, even if its depth suggests wet. This prevents isolated wet nodes from persisting.

Slope Limiter (ALPHA)

The slope limiter \(\alpha_L\) reduces pressure gradients in shallow water to prevent spurious flows. The limiter is controlled by the wet_dry_slope_limiter configuration parameter (ADCIRC’s SLIM), with a default of \(10^9\) which effectively disables limiting.

When limiting is enabled (e.g., wet_dry_slope_limiter: 0.0004), the scaling factor is computed based on the ratio of bathymetric gradient to water surface gradient:

(24)\[\alpha_L = \frac{\text{SLIM} \cdot |\nabla h|}{|\nabla \zeta|}\]

The limiter is only evaluated for active elements whose nodes are all wet and that contain shallow water (total depth at or below the shallow-water threshold). This factor multiplies the pressure gradient terms in both GWCE and momentum:

  • \(\alpha_L = 1\): Full pressure gradient (deep water or limiting disabled)

  • \(\alpha_L < 1\): Reduced pressure gradient (shallow water with limiting)

Data Structures

WetDryData

The WetDryData class stores all wet/dry state:

class WetDryData {
  NodeStatusView m_node_wet_status;      // Committed status
  NodeStatusView m_node_wet_working;     // Working during algorithm

  ElementStatusField m_element_active;   // TemporalField<ElementStatus, 2, 1>
  OrdinalView m_active_element_count;    // Per-node count

  FloatView m_slope_limiter;             // ALPHA values (float storage)
};

Node-to-Element Connectivity

For efficient iteration, DeviceMesh stores CSR (Compressed Sparse Row) format connectivity:

OrdinalView m_node_to_elem_offsets;  // [num_nodes + 1]
OrdinalView m_node_to_elem_list;     // [total_connections]

Implementation Files

  • WetDry.hpp: Algorithm orchestrator with static compute() method

  • WetDryKernels.hpp: Individual phase kernels

  • WetDryData.hpp: Data container class