Skip to content

Erosion Calculation Process#

This document explains the bank erosion calculation process in the D-FAST Bank Erosion software, focusing on the compute_erosion_per_level method in the Erosion class.

Overview#

The compute_erosion_per_level method is a core component of the bank erosion calculation process. It calculates the bank erosion for a specific discharge level, processing both the left and right banks of the river. The method computes various parameters such as bank velocity, water depth, and erosion distances, and returns a SingleDischargeLevel object containing the results for both banks.

Data Model#

Before diving into the calculation process, it's important to understand the data model used in the bank erosion calculations:

classDiagram
    class DischargeLevels {
        +List~SingleDischargeLevel~ levels
        +__getitem__(int)
        +__iter__()
        +append(SingleDischargeLevel)
        +get_max_hfw_level()
        +total_erosion_volume()
    }

    class SingleDischargeLevel {
        +SingleCalculation left
        +SingleCalculation right
        +float hfw_max
    }

    class SingleCalculation {
        +array bank_velocity
        +array water_level
        +array chezy
        +array erosion_distance_tot
        +array erosion_volume_tot
        +array erosion_distance_shipping
        +array erosion_distance_flow
        +array erosion_distance_eq
        +array erosion_volume_eq
        +array volume_per_discharge
    }

    class SingleLevelParameters {
        +SingleParameters left
        +SingleParameters right
    }

    class SingleParameters {
        +float discharge
        +float probability
        +dict ship_parameters
    }

    DischargeLevels o-- SingleDischargeLevel : contains
    SingleDischargeLevel *-- SingleCalculation : has left/right
    SingleLevelParameters *-- SingleParameters : has left/right

The DischargeLevels class contains a list of SingleDischargeLevel objects, each representing the erosion calculation results for a specific discharge level. Each SingleDischargeLevel has left and right attributes of type SingleCalculation, which store the calculation results for the left and right banks, respectively.

Calculation Process#

The compute_erosion_per_level method follows these steps:

flowchart TD
    A[Initialize variables] --> B[Loop through each bank]
    B --> C[Calculate bank velocity]
    C --> D{Is first level?}
    D -- Yes --> E[Calculate bank height]
    D -- No --> F[Get fairway data]
    E --> F
    F --> F1[Extract water level, chezy, and water depth]
    F1 --> G{Is last level?}
    G -- Yes --> H[Calculate equilibrium erosion]
    G -- No --> I[Calculate bank erosion dynamics]
    H --> I
    I --> J[Accumulate eroded volumes per km]
    J --> K[Add calculation to results]
    K --> L{More banks?}
    L -- Yes --> B
    L -- No --> M[Create SingleDischargeLevel object]
    M --> N[Return results]

Detailed Steps#

  1. Initialize variables:
  2. Create an array to store eroded volumes per kilometer
  3. Initialize the maximum water depth to 0

  4. Create an empty list to store calculation results for each bank

  5. Loop through each bank (left and right):

  6. Create a new SingleCalculation object to store results

  7. Calculate bank velocity:

  8. Compute the velocity along the bank using the simulation data

  9. Calculate bank height (only for the first level):

  10. Determine the maximum bed elevation per cell along the bank

  11. Get fairway data:

  12. Get fairway face indices from the bank
  13. Retrieve fairway data (water level, chezy, water depth) in a single call
  14. Extract and assign the values to the SingleCalculation object

  15. Update the maximum water depth if necessary

  16. Calculate equilibrium erosion (only for the last level):

  17. Compute the equilibrium erosion distance and volume using the ErosionCalculator

  18. Store the results in the SingleCalculation object

  19. Calculate bank erosion dynamics:

  20. Compute the bank erosion during the current discharge level using the ErosionCalculator

  21. This includes calculating erosion due to both flow and shipping

  22. Accumulate eroded volumes per kilometer:

  23. Aggregate the eroded volumes by kilometer bins
  24. Store the results in the SingleCalculation object

  25. Add the volumes to the total eroded volumes array

  26. Output debug information (if debug mode is enabled):

  27. Write detailed information about the calculations to files

  28. Create a SingleDischargeLevel object:

    • Combine the calculation results for both banks into a single object
    • Include the maximum water depth in the fairway

  29. Return results:

    • Return the SingleDischargeLevel object and the eroded volumes array

Sequence Diagram#

The following sequence diagram illustrates the interaction between different components during the execution of the compute_erosion_per_level method:

sequenceDiagram
    participant Erosion
    participant BankData
    participant SimulationData
    participant ErosionCalculator
    participant SingleCalculation
    participant SingleDischargeLevel

    Erosion->>Erosion: compute_erosion_per_level(level_i, ...)
    loop for each bank (left and right)
        Erosion->>SingleCalculation: create new calculation object
        Erosion->>SimulationData: calculate_bank_velocity(bank_i, vel_dx)
        SimulationData-->>Erosion: bank_velocity
        Erosion->>SingleCalculation: set bank_velocity

        alt is first level
            Erosion->>SimulationData: calculate_bank_height(bank_i, zb_dx)
            SimulationData-->>Erosion: bank_height
        end

        Erosion->>SimulationData: get_fairway_data(fairway_face_indices)
        SimulationData-->>Erosion: data (containing water_level, chezy, water_depth)
        Erosion->>SingleCalculation: set water_level, chezy, water_depth from data

        alt is last level
            Erosion->>ErosionCalculator: comp_erosion_eq(...)
            ErosionCalculator-->>Erosion: erosion_distance_eq, erosion_volume_eq
            Erosion->>SingleCalculation: set erosion_distance_eq, erosion_volume_eq
        end

        Erosion->>ErosionCalculator: compute_bank_erosion_dynamics(...)
        ErosionCalculator-->>Erosion: updated SingleCalculation

        Erosion->>Erosion: get_km_eroded_volume(...)
        Erosion->>SingleCalculation: set volume_per_discharge
        Erosion->>Erosion: append to par_list

        alt debug mode enabled
            Erosion->>Erosion: _debug_output(...)
        end
    end

    Erosion->>SingleDischargeLevel: create with left=par_list[0], right=par_list[1]
    Erosion-->>Erosion: return level_calculation, dvol_bank

Data Flow#

The data flow in the compute_erosion_per_level method can be visualized as follows:

flowchart LR
    subgraph Inputs
        A[bank_data]
        B[simulation_data]
        C[fairway_data]
        D[single_parameters]
        E[erosion_inputs]
        F[km_bin]
    end

    subgraph Processing
        G[Calculate bank velocity]
        H[Calculate bank height]
        I1[Get fairway data]
        I2[Extract water level, chezy, water depth]
        J[Calculate equilibrium erosion]
        K[Calculate bank erosion dynamics]
        L[Accumulate eroded volumes]
    end

    subgraph Outputs
        M[SingleDischargeLevel]
        N[dvol_bank]
    end

    A --> G
    B --> G
    B --> H
    B --> I1
    A --> I1
    I1 --> I2

    C --> J
    D --> J
    E --> J
    A --> J
    I2 --> J

    G --> K
    H --> K
    I2 --> K
    C --> K
    D --> K
    E --> K

    K --> L
    F --> L
    A --> L

    L --> M
    L --> N
    J --> M

Key Calculations#

The compute_erosion_per_level method relies on two key calculations performed by the ErosionCalculator:

  1. Equilibrium Erosion Calculation (for the last discharge level):
  2. Calculates the equilibrium bank erosion distance and volume
  3. Takes into account bank height, water level, and ship-induced wave height

  4. Returns the maximum potential erosion distance and volume

  5. Bank Erosion Dynamics Calculation (for all discharge levels):

  6. Calculates the bank erosion during a specific discharge level
  7. Computes erosion due to both flow and shipping
  8. Takes into account bank velocity, water level, and ship parameters
  9. Returns the erosion distance and volume for the current level

Physical Processes and Mathematical Models#

The erosion calculation process is based on physical processes and mathematical models that describe the interaction between water flow, ship waves, and bank material. This section explains the key formulas and calculations used in the ErosionCalculator class.

Ship-Induced Wave Height Calculation#

The height of waves generated by passing ships is a key factor in bank erosion. The wave height at the bank is calculated using the following formula:

graph TD
    A[Calculate ship coefficient a1 based on ship type] --> B[Calculate Froude number]
    B --> C[Limit Froude number to 0.8]
    C --> D[Calculate wave attenuation factor A]
    D --> E[Calculate wave height h0]
  1. Ship coefficient (a1) depends on ship type:
  2. Multiple barge convoy: a1 = 0.5
  3. RHK ship / motor ship: a1 = 0.28 * ship_draught^1.25

  4. Towboat: a1 = 1

  5. Froude number is calculated as: 

    froude = ship_velocity / √(h * g)
    where:

  6. ship_velocity is the velocity of the ship [m/s]
  7. h is the water depth at the fairway [m]
  8. g is the gravitational acceleration [m/s²]

The Froude number is limited to 0.8 to prevent unrealistic wave heights.

  1. Wave attenuation factor (A) accounts for the reduction in wave height with distance from the fairway: 
    A = 0.5 * (1 + cos((bank_fairway_dist - fairway_wave_disappear_distance) / (fairway_wave_reduction_distance - fairway_wave_disappear_distance) * π))
    with the following conditions:
  2. If bank_fairway_dist < fairway_wave_disappear_distance, then A = 1

  3. If bank_fairway_dist > fairway_wave_reduction_distance, then A = 0

  4. Wave height (h0) is calculated as: 

    h0 = a1 * h * (bank_fairway_dist / h)^(-1/3) * froude^4 * A

Equilibrium Erosion Calculation#

The equilibrium erosion represents the maximum potential erosion that could occur over a very long time period. It is calculated as follows:

graph TD
    A[Calculate ship-induced wave height h0] --> B[Determine upper bound of erosion zone zup]
    B --> C[Determine lower bound of erosion zone zdo]
    C --> D[Calculate height of erosion zone ht]
    D --> E[Calculate height above erosion zone hs]
    E --> F[Calculate equilibrium erosion distance]
    F --> G[Calculate equilibrium erosion volume]

  1. Upper bound of erosion zone (zup): 

    zup = min(bank_height, water_level_fairway_ref + 2 * h0)

  2. Lower bound of erosion zone (zdo): 

    zdo = max(water_level_fairway_ref - 2 * h0, bank_protection_level)

  3. Height of erosion zone (ht): 

    ht = max(zup - zdo, 0)

  4. Height above erosion zone (hs): 

    hs = max(bank_height - water_level_fairway_ref + 2 * h0, 0)

  5. Equilibrium erosion distance (eq_erosion_distance): 

    eq_erosion_distance = ht / mu_slope
    where mu_slope is the bank slope [-]

  6. Equilibrium erosion volume (eq_erosion_volume): 

    eq_erosion_volume = (0.5 * ht + hs) * eq_erosion_distance * segment_length

Bank Erosion Dynamics Calculation#

The bank erosion dynamics calculation determines the erosion that occurs during a specific discharge level, taking into account both flow-induced and ship-induced erosion.

Flow-Induced Erosion#

Flow-induced erosion occurs when the flow velocity exceeds the critical velocity for the bank material:

graph TD
    A[Calculate erosion coefficient] --> B[Calculate critical velocity]
    B --> C[Calculate velocity ratio]
    C --> D[Calculate flow-induced erosion distance]
    D --> E[Calculate flow-induced erosion volume]

  1. Erosion coefficient (erosion_coef): 

    erosion_coef = 0.2 * √(tauc) * 10^-6
    where tauc is the critical shear stress [N/m²]

  2. Critical velocity (critical_velocity): 

    critical_velocity = √(tauc / (water_density * chezy² / g))
    where:

  3. water_density is the density of water [kg/m³]
  4. chezy is the Chezy coefficient [m^0.5/s]

  5. g is the gravitational acceleration [m/s²]

  6. Velocity ratio (crit_ratio): 

    crit_ratio = (vel / critical_velocity)²
    if vel > critical_velocity and water_level > bank_protection_level

  7. Flow-induced erosion distance (erosion_distance_flow): 

    erosion_distance_flow = erosion_coef * (crit_ratio - 1) * time_erosion * sec_year
    where:

  8. time_erosion is the erosion period [yr]

  9. sec_year is the number of seconds in a year (3600 * 24 * 365)

  10. Flow-induced erosion volume (dv_flow): 

    dv_flow = erosion_distance_flow * segment_length * h_line_flow
    where h_line_flow is the height of the bank subject to flow erosion

Ship-Induced Erosion#

Ship-induced erosion occurs due to waves generated by passing ships:

graph TD
    A[Calculate ship wave period] --> B[Calculate total wave exposure time]
    B --> C[Calculate strength coefficient]
    C --> D[Calculate ship-induced erosion distance]
    D --> E[Calculate ship-induced erosion volume]

  1. Ship wave period (ship_wave_period): 

    ship_wave_period = 0.51 * ship_velocity / g

  2. Total wave exposure time (ts): 

    ts = ship_wave_period * num_ship * num_waves_per_ship
    where:

  3. num_ship is the number of ships

  4. num_waves_per_ship is the number of waves per ship

  5. Strength coefficient (cE): 

    cE = 1.85 * 10^-4 / tauc

  6. Ship-induced erosion distance (erosion_distance_shipping): 

    erosion_distance_shipping = cE * wave_height² * ts * time_erosion
    if ship_wave_min < water_level_fairway_ref < ship_wave_max

  7. Ship-induced erosion volume (dv_ship): 

    dv_ship = erosion_distance_shipping * segment_length * h_line_ship
    where h_line_ship is the height of the bank subject to ship-induced erosion

Total Erosion#

The total erosion is the sum of flow-induced and ship-induced erosion:

  1. Total erosion distance (erosion_distance): 

    erosion_distance = erosion_distance_shipping + erosion_distance_flow

  2. Total erosion volume (erosion_volume): 

    erosion_volume = dv_ship + dv_flow

Usage Example#

The compute_erosion_per_level method is called within the _process_discharge_levels method of the Erosion class, which processes all discharge levels:

def _process_discharge_levels(self, ...):
    # ...
    for level_i in range(num_levels):
        # ...
        level_calculation, dvol_bank = self.compute_erosion_per_level(
            level_i,
            bank_data,
            simulation_data,
            fairway_data,
            single_parameters,
            erosion_inputs,
            km_bin,
            num_km,
            bank_height,
        )
        discharge_levels.append(level_calculation)
        # ...

The results are then used to calculate the total erosion volume and other parameters for the entire river.