Erosion Processes

Resuspension Erosion flux E : (Partheniades)

\(E=E_0*\alpha*[\frac{\tau-\tau_{ce}}{\tau_{ce}}]^n\)
  • \(E_0\) : erodability
  • \(\alpha\) : correction factor
  • \(\tau\) : shear stress
  • \(\tau_{ce}\) : critical shear stress for erosion
  • \(n\) : power coefficient applied to excess erosion

Erosion Processes : Version 1

  • Erosion flux is applied in proportion to the classes present, by modulating parameters depending on the nature of the superficial sediment

  • the critical erosion shear stress of erosion \(\tau_{ce}\) depends on the content of sand / mud and on sediment consolidation.

  • \(\alpha\) is a correction factor
    • \(\alpha=1.+corfluer1*\max(0,corfluer2-Csed_{total-surf})\) (corfluer1 et corfluer2 given in paraMUSTANGV1.txt)
  • \(E_0\) and \(n\) depend on the mud fraction in sediment
    • \(fr_{mud} < frv1\) : sandy sediment (frv1=0.3 typically = frmudcr1 calculated from coef_frmudcr1 in paraMUSTANGV1.txt)
    • \(fr_{mud} > frv2\) : muddy sediment (frv2=0.7 typically = frmudcr2 given in paraMUSTANGV1.txt)
  • For muddy sediment :
    • \(E_0=E0_{mud}\) (given in paraMUSTANGV1.txt)
    • \(\tau_{ce}=x1toce_{mud}*Crel_{mud}^{x2toce_{mud}}\) (parameters given in paraMUSTANGV1.txt) with \(Crel_{mud}\) which is the relative mud concentration between sand grains
  • For sandy sediment :

    • \(E_0=E0_{sand}*E0_{sand-para}\)

    • \(E0_{sand}\) given in paraMUSTANGV1.txt or estimated by formulations depending on the option chosen by E0_sand_option

      • E0_sand_option = 0 ==> E0_sand= E0_sand_cst read in this namelist paraMUSTANGV1.txt

      • E0_sand_option = 1 ==> E0_sand evaluated with Van Rijn (1984) formulation

      • E0_sand_option = 2 ==> E0_sand evaluated with erodimetry formulation

        \(E0_{sand}=\min(0.27,1000*d50-0.01)*(\tau-\tau_{ce})^{n_{eros-sand}}\)

      • E0_sand_para is a multiplicative factor to correct erosion flux

    • \(\tau_{ce}\) estimated from sand characteristics

  • For sand/mud mixing : Erosion flux depends on proportions of the mixture

    • The option choice is given by ero_option in paraMUSTANGV1.txt
      • ero_option = 0 : pure mud behavior (for all particles and whatever the mixture)
      • ero_option = 1 : linear interpolation between sand and mud behavior, depend on proportions of the mixture
      • ero_option = 2 : formulation derived from that of J. Vareilles (2013)
      • ero_option = 3 : formulations proposed by B. Mengual (2015) with exponential coefficients depend on proportions of the mixture

Note

Several options are available but all are questionable and should be tested and used carefully

Erosion and Bedload Processes : Version 2

  • The new model allows a sharper distinction of erosion regimes, non-cohesive and cohesive, by offering the possibility of adding the transport process by bedload for gravel and sand and a flux of resuspension independent of the different classes in non-cohesive regime.

  • Just like the version 1, the transition from one regime to another depends on the nature of the sediment through the concept of a critical mud fraction frmudcr1 (calculated as in Le Hir et al (2011) as a function of average diameter of the sands). From the composition of the surface layer and the frmudcr1 associated, it is possible to know if the sediment in the superficial layer is non-cohesive or cohesive.

  • During a time step, different iterations can occur by alternating between the two regime. Whatever the number of iterations, the flows by bedload and resuspension are gradually integrated during the time step for each class of sediment and in each mesh.

  • Basic principles for calculating bedload :
    • The horizontal fluxes (kg/m/s) are projected parallel to the shear stress and evaluated at the center of the meshes.
    • The two components of the flux are taken entirely during time step, and treated as erosion. The bedload contributions of the adjacent meshes are then treated as an additional deposit.
    • The flux taken is limited according to the quantity available in the active layer
    • Each component is assigned to the “downstream” mesh edge: the flow is “decentered upstream”

  • Initialization of erosion parameters for non cohesive sediments (gravel and sand) :

    • initialization of the erosion parameters is done at the beginning of the simulation according to the options chosen in paraMUSTANGV2.txt

    • \(E_0\) depending on the boolean *E0_sand_option*
      • E0_sand_option = 0 ==> E0_sand= E0_sand_cst read in this namelist paraMUSTANGV2.txt
      • E0_sand_option = 1 ==> E0_sand evaluated with Van Rijn (1984) formulation, modified by J. Vareilles (2013)
      • E0_sand_option = 2 ==> E0_sand evaluated with erodimetry formulation

      \(E0_{sand}=\min(0.27,1000*d50-0.01)*(\tau-\tau_{ce})^{n_{eros-sand}}\)

      • E0_sand_option = 3 ==> formulation of Wu and Lin (2014): E0 is deduced from the equilibrium relation between the flux of erosion and the product of the settling rate by the reference concentration C_{ref} at a given height relative to the bottom (with a C_{ref} formulation specific to Wu and Lin, 2014).
      • E0 is then multiplicated by E0_sand_para (corrective factor)
    • \(\tau_{ce}\) depending on the boolean *tau_cri_option*
      • tau_cri_option = 0 ==> \(\tau_{ce}\) calculated from the Shields formula, modified by Soulsby (1997) for small diameters.
      • tau_cri_option = 1 ==> \(\tau_{ce}\) calculated according to the formulation of Wu and Lin (2014) ignoring the factor associated with masking / exposure effects. This formulation has the particularity of using a Shields mobility criterion constant of 0.03. The effects of masking / exposure can be added thanks to the activation of Booleans to inform in paraMUSTANGV2.txt.

    • \(n=n_{eros-sand}\) given in paraMUSTANGV2.txt, applicated to shear stress excess for all non-cohesives classes.

  • Estimation of erosion flux :

    • the surface sediment is characterized by its mud fraction which is compared with frmudcr1 calculated as a function of average diameter of the sands

    • if frmud > frmudcr1 : the sediment is treated in cohesif regime, and non cohesif if not.

    • There are two ways to deal with erosion of non-cohesive sediment, and the choice is made with the boolean l_eroindep_noncoh :
      • if l_eroindep_noncoh = False :
        • the erosion is treated very similar to version V1, without integrating bedload.
        • All classes erode homogeneously and the erosion parameters are modulated according to the nature of the sediment.
        • The procedure is the same as in a cohesive regime.
      • if l_eroindep_noncoh = True :
        • the erosion is treated independently for different classes of sand and silt;
        • bedload is taken into account

    • workflow to manage erosion fluxes and reworking of sediment layers is shown in erosion routine diagram - Version V2

    • non-cohesive regime, i.e. the superficial mud fraction is less than frmudcr1. erosion is treated class by class independently (if l_eroindep_noncoh = True*)

      • the different classes of sand and mud are independant (the gravel being transported only by bedload).
      • If the critical stress representative of the sedimentary mixture in the surface layer is greater than the shear stress near the bottom (waves + current), then an active layer is created.
      • The thickness of the active layer is calculated from the formulation of Harris and Wiberg (1997). A theoretical active layer thickness is calculated based on the excess shear near the bottom and the representative diameter of the whole non-cohesive classes (gravel and sand). We use two coefficients k1HW97 and k2HW97, whose value is given in paraMUSTANGV2.txt (default 0.07 and 6, respectively, according to Harris and Wiberg, 1997).
        • thanks to an iterative procedure, the layers of sediment from the water / sediment interface are merged until one of the following conditions is validated:
            1. the thickness of fused sediment reaches the theoretical thickness of the active layer,
            1. a cohesive layer is encountered (fm > frmudcr1),
            1. there is only one layer left.
        • this iterative process is modulated by using a parameter fusion_para_activlayer (given in paraMUSTANGV2.txt) :
          • if fusion_para_activlayer =0 ==> no fusion : superficial layer = active layer
          • if fusion_para_activlayer =1 ==> criterion (ii) used ( fm > frmudcr1 )
          • if 0 < fusion_para_activlayer < frmudcr2 ==> criterion (ii) modulated (\(fm > fusion_{para_activlayer}*frmudcr1\) )
      • If key_MUSTANG_bedload is used, bedload fluxes are evaluated for all non cohesive classes whode boolean l_bedload is equal to .TRUE. (given in variable.dat). These developments come from A. Rivier et al (2017)
        • The formulation chosen for calculating the fluxes carried is that of Wu and Lin (2014). This calculation depends in particular excess shear relative to the critical erosion stress.
        • By activating the Boolean l_peph_bedload in paraMUSTANGV2.txt, the critical stress of erosion relative to a given class will account of a masking / exposure effect related to the presence of other non-cohesive classes (see Wu and Lin, 2014).
        • These fluxes (in kg/m/s) are then projected in (x, y) according to the direction of the bottom shear stress exerted by the current.
        • The slope effects likely to modulate the bedload fluxes can be taken into account by informing l_slope_effect_bedload to .TRUE. in paraMUSTANGV2.txt, (from Lesser’s theory, 2004).
        • Finally, the fluxes carried in (x, y) are set to 0 if the neighboring mesh is exonde.
      • resupension fluxes are then evaluated for each sediment variable independently (Partheniades-Ariathurai) with different parameters (\(E_0\), \(\tau_{ce}\) and \(n\) ) specific to each variable.
        • For non-cohesive classes, the activation of a boolean l_peph_suspension (given in paraMUSTANGV2.txt) allows to take into account the effects of masking / exposure related to the presence of other noncohesive classes (see Wu and Lin, 2014) on the critical erosion stress estimation relative to a given class.
        • Still in the case of non-cohesive classes if the way of calculating the parameter of erodibility \(E_0\) does not follow the formalism of Wu and Lin (2014) (i.e., E0_sand_option ≠ 3) and that class is also carried by bedload, so it is possible to correct the resuspension flow of the class to satisfy the total transport and the bedload / resuspension ratios described by Wu and Lin (2014) in response to a constraint forcing. This is possible thanks to the activation of the boolean l_fsusp in the paraMUSTANGV2.txt, which leads to the application of a factor to the \(E_0\) of the class concerned (fsusp) equal to the ratio of transport in suspension and total transport associated with the constraint regime.
        • For cohesive classes, only transported in suspension, different options have been introduced concerning \(E_0\) and \(\tau_{ce}\).
          • for \(E_0\) , it is possible to apply a factor ( E0_mud_para_noncoh ) to the parameter erosion of the E0_mud (defined in paraMUSTANGV2.txt) to make the mud more or less erodible in a non-cohesive regime.
          • for \(\tau_{ce}\) , there are several options with *tau_cri_mud_option_eroindep* (given in paraMUSTANGV2.txt)
            • tau_cri_mud_option_eroindep = 0 ==> evaluated from the relative mud concentration betwwen non cohesives grains
            • tau_cri_mud_option_eroindep = 1 ==> minimum value between the average critical stress representative of sands within the mixture and the critical stress calculated from the relative mud concentration in the non-cohesive matrix.
            • tau_cri_mud_option_eroindep = 2 ==> minimum value of products between critical constraints and fractions for all sands present within the mixture.
            • tau_cri_mud_option_eroindep = 3 ==> minimum value between the critical stress of the finest sand and that calculated from the relative mud concentration between the non-cohesive grains.
      • sediment erosion in the active layer . The total erosion flux including both resuspension and bedload is calculated class by class. This flux multiplied by the time allotted to erosion is compared to the sediment mass available in the active layer.
        • If the mass to erode is greater than or equal to the mass of the available sediment class, then the active layer is completely emptied of the concerned class.
        • If not, only part of the concerned sediment class contained in the active layer is eroded
        • The update of the characteristics post-erosion of the surface layer (thickness, concentrations) is done via a new porosity calculation taking into account the different classes remaining.
        • If the new layer thickness is less than dzsmin, then the layer is merged with that from below to prevent cases where one would obtain infinitely small concentrations.
        • The resuspension and bedload fluxes related to the iteration are added to those from previous iterations made during the same time step (including iterations in non-cohesive and cohesive regimes).

    • cohesive regime

      • In cohesive regime, i.e., when the superficial mud fraction (fm) is greater than frmudcr1, the physical behavior of the matrix gravel-sandy is influenced by muddy sediments. This influence grows as the mud fraction increases within the mixture. In this case the erosion of the mixture applies homogeneous.

      • \(E_0\), \(\tau_{ce}\) and \(n=n\) are calculated according to the respective fractions of the different classes within the mixture. The erosion parameters in a cohesive regime tend (more or less rapidly) towards those associated with a purely muddy behavior (when fm increases). Different options have been introduced to manage the type and speed of transition into a cohesive regime. Depending on the value of the ero_option (to be informed in paraMUSTANGV2.txt), the parameters of erosion in cohesive regime are computed in different ways
        • ero_option = 0
          : pure mud behavior
          • \(E_0=E0_{mud}\) (given in paraMUSTANGV2.txt)
          • \(\tau_{ce}=x1toce_{mud}*Crel_{mud}^{x2toce_{mud}}\) (parameters given in paraMUSTANGV2.txt) with \(Crel_{mud}\) which is the relative mud concentration between non cohesives grains
          • \(n=n_{eros-mud}\) given in paraMUSTANGV2.txt, applicated to shear stress excess for all mud variables.
        • ero_option = 1 : linear interpolation between sand and mud behavior, depend on proportions of the mixture. The transition regime applies between frmudcr1 and a second critical fraction frmudcr2, to be informed in the paraMUSTANGV2.txt (default: 0.5).
        • ero_option = 2 : same logic as case 1 except that the transition follows a power law (e.g. Carniello et al, 2012)
        • ero_option = 3 : same logic as case 1 except that the transition varies exponentially and more or less fast (formulations proposed by B. Mengual (2017) and modified in 2018). The coefficient controlling the speed of the exponential transition is either equal to that prescribed by the user via the parameter xexp_ero (in paraMUSTANGV2.txt; default: 10), either, if l_xexp_ero_cst =.FALSE. , increases linearly depending on frmudcr1. So the more the average diameter of the sands increases, the more the transition between non-cohesive and cohesive regime intervenes for a large fraction of mud but the faster this transition is.
        ../_images/toce_ero_option3.jpg

      • if \(\tau > \tau_{ce}\) , erosion flux is then evaluated (Partheniades)

      • the effective erosion of the first layer of sediment is carried out. In the case where quantity to be eroded is smaller than or equal to the sandy-muddy stock in the superficial layer, then the latter is eroded in proportion to the classes. The “post-erosion” characteristics (thickness, concentrations) of the layer are updated ensuring that its thickness remains strictly greater than dzsmin. If this condition is not met, then the remaining quantities of sediment are placed in the layer underneath.

      • If the amount of sandy-mud sediment in the superficial layer is not enough to satisfy erosion, then all classes of sands and muds are eroded.
        • If this layer is devoid of gravel then the erosion of the latter is total.
        • If gravel is present, a new porosity calculation makes it possible to update the thickness of the layer and a paving effect is put in place.

      • If erosion time has not been fully consumed and if the sediment has more than one layer, then a new iteration will be initiated, but this time in a non-cohesive regime.

Lateral erosion

  • A lateral erosion process of a wet or dry cell was introduced into the model (only for suspensions).
    • Parameters coef_erolat, l_erolat_wet_cell and coef_tenfon_lat are used for this option., given in paraMUSTANGV1.txt or paraMUSTANGV2.txt
    dry cell :
    \(Ero_{lat}=(\alpha_{lat} * U_{neighb}^2-\tau_{ce})*H_{neighb}\)

    wet cell (if l_erolat_wet_cell = .TRUE.)
    \(Ero_{lat}=(\beta_{lat} * \alpha_{lat} * U_{neighb}^2-\tau_{ce})*\delta_H\)

    with \(\alpha_{lat}\) =coef_tenfon_lat
    with \(\beta_{lat}\) =coef_erolat
    with \(U_{neighb}\) =water current in the neigboring cell
    with \(H_{neighb}\) =depth in the neighboring cell
    with \(\delta_H\) =depth difference between eroded cell and the neighboring
    		 ../_images/lateral_erosion1.jpg

erosion routine diagram - Version V2

../_images/routine_erosion_MUSTANGV2.jpg