NEWS.txt

 Latest developments on main
 ===========================
 
+GAIA: fix units for CVSM serafin output file.
+
+GAIA: fix the CVSM output period for output of specified nodes in order to
+respect the value set in the steering file.
+
+Python: fix TelemacFile.meshz so that it actually returns the Z coordinates.
+
+TELEMAC-2D: add failure criteria for breach initiation.
+
+STBTEL: fix an issue when interpolating with 4 bathymetry files.
+
+GAIA: fix a coefficient in the bedload formulas of Engelund-Hansen
+according to the value given by the literature.
+
 Python: Telapy for TELEMAC-2D RESTART MODE
 
 TELEMAC-2D: Fortran API for RESTART MODE for perfect restart

configs/systel.edf.cfg

@@ -621,7 +621,7 @@ incs_mumps: -I%MUMPSHOME%/include
 flags_mumps: -DHAVE_MUMPS
 libs_mumps: -L%MUMPSHOME%/lib -ldmumps -lmumps_common -lpord
             -L%SCALAPACKHOME%/lib -lscalapack
-            -L%LAPACKHOME%/lib -llapack -lblas
+            -L%OPENBLASHOME%/lib -lopenblas
 #
 # AED info
 #
@@ -682,7 +682,7 @@ incs_mumps: -I%MUMPSHOME%/include
 flags_mumps: -DHAVE_MUMPS
 libs_mumps: -L%MUMPSHOME%/lib -ldmumps -lmumps_common -lpord
             -L%SCALAPACKHOME%/lib -lscalapack
-            -L%LAPACKHOME%/lib -llapack -lblas
+            -L%OPENBLASHOME%/lib -lopenblas
 #
 # AED info
 #
@@ -740,7 +740,7 @@ incs_mumps: -I%MUMPSHOME%/include
 flags_mumps: -DHAVE_MUMPS
 libs_mumps: -L%MUMPSHOME%/lib -ldmumps -lmumps_common -lpord
             -L%SCALAPACKHOME%/lib -lscalapack
-            -L%LAPACKHOME%/lib -llapack -lblas
+            -L%OPENBLASHOME%/lib -lopenblas
 #
 # AED info
 #
@@ -801,7 +801,7 @@ incs_mumps: -I%MUMPSHOME%/include
 flags_mumps: -DHAVE_MUMPS
 libs_mumps: -L%MUMPSHOME%/lib -ldmumps -lmumps_common -lpord
             -L%SCALAPACKHOME%/lib -lscalapack
-            -L%LAPACKHOME%/lib -llapack -lblas
+            -L%OPENBLASHOME%/lib -lopenblas
 #
 # AED info
 #

documentation/data/biblio.bib

@@ -4393,3 +4393,13 @@ doi = {10.1007/s11069-020-03862-8},
   pages    = {86-102},
   doi      = {https://doi.org/10.1016/j.coldregions.2008.03.004},
 }
+@Article{Goeury2022,
+  author  = {Goeury, Cédric and Bacchi, Vito and Zaoui, Fabrice and Bacchi, Sophie and Pavan, Sara and El kadi Abderrezzak, Kamal},
+  title   = {Uncertainty Assessment of Flood Hazard Due to Levee Breaching},
+  journal = {Water},
+  year    = {2022},
+  volume =   {14},
+  number =   {23},
+  doi      = {10.3390/w14233815},
+}
+

documentation/gaia/user/latex/3_1-SedimentTransportProcessesBedAndStratigraphy_bedload.tex

@@ -124,7 +124,7 @@ where the non-dimensional diameter $D_*=d[(\rho_s/\rho-1)g/\nu^2]^{1/3}$, with $
  \item \telkey{BED-LOAD TRANSPORT FORMULA FOR ALL SANDS = 3}
  \item The dimensionless current-induced sediment transport rate is given by:
   \begin{equation*}
-  \Phi_b=0.05 \frac{\theta_*^{\frac{5}{2}}}{c_f}
+  \Phi_b=0.1 \frac{\theta_*^{\frac{5}{2}}}{c_f}
   \end{equation*}
  with $c_f$ the adimensional friction coefficient and $\theta_*$ depending on the transport regime:
 \begin{equation*}
@@ -212,7 +212,7 @@ Further details on the WC-2003 formula and applications can be found in~\cite{Co
  \item \telkey{BED-LOAD TRANSPORT FORMULA FOR ALL SANDS = 30}
  \item The dimensionless current-induced sediment transport rate is given by:
   \begin{equation*}
-  \Phi_b=0.05 \frac{\theta^{\frac{5}{2}}}{c_f}
+  \Phi_b=0.1 \frac{\theta^{\frac{5}{2}}}{c_f}
   \end{equation*}
  with $c_f$ the adimensional friction coefficient and $\theta$ the Shields number without the correction factor for skin friction ($\theta=\frac{\tau_b}{(\rho_s-\rho)gd}$).
  \item Fortran subroutine {\ttfamily bedload\_engel\_gaia.f}.
@@ -468,7 +468,7 @@ For this case, the flag \texttt{LIEBOR} is set \texttt{= 5} and the flag \texttt
 4 5 5 0.0 0.0 0.0 0.0 (*@\color{PantoneRed}5@*) (*@\color{PantoneRed}0.0@*) 0.0 0.0 565 1
 \end{lstlisting}
 
-\subsubsection{Constant sediment discharge}
+\subsubsection{Constant sediment discharge in time}
 For this case, boundary condition files are needed for both \telemac{2D} and \gaia{}. In the \gaia{}'s boundary condition file, the flag \texttt{LIQBOR = 5} and \texttt{LIEBOR = 4}. The imposed solid discharge can be specified as follows:
 \begin{itemize}
 \item A value of the unit solid discharge [kg/(m~s)] in the column \texttt{Q2BOR} of the \gaia{}'s boundary condition file, as shown in the example below for an imposed unit discharge \texttt{Q2BOR=}$1.0$~kg/(m~s):
@@ -478,25 +478,37 @@ For this case, boundary condition files are needed for both \telemac{2D} and \ga
 
 Particular cases of \texttt{Q2BOR} can be programmed in the subroutine \texttt{conlit\_gaia.f}.
 
-\item A value of the total solid discharge (without pores) [kg/s] given through the keyword \telkey{PRESCRIBED SOLID DISCHARGES} (sequence of real values separated by semi-colons, one value per liquid boundary, no default value) in the steering file, as shown in the example below for an imposed total discharge equal to $1.0$~kg/s:
+\item A value of the total solid discharge (without pores) [kg/s] given through the keyword \telkey{PRESCRIBED SOLID DISCHARGES} (sequence of real values separated by semi-colons, one value per liquid boundary, no default value) in the steering file, as shown in the example below for an imposed total discharge equal to $10.0$~kg/s:
 \begin{lstlisting}[frame=trBL]
-4 (*@\color{PantoneRed}5@*) 5 0.0 0.0 0.0 0.0 (*@\color{PantoneRed}4@*) 0.0 0.0 0.0 565 1
+4 (*@\color{PantoneRed}5@*) 5 1.0 0.0 0.0 0.0 (*@\color{PantoneRed}4@*) 0.0 0.0 0.0 565 1
 \end{lstlisting}
-
 \begin{lstlisting}[frame=trBL]
-PRESCRIBED SOLID DISCHARGES : 1.0
+PRESCRIBED SOLID DISCHARGES : 10.0
 \end{lstlisting}
+When a value of solid discharge is given in the parameter file, \texttt{Q2BOR} is considered as only a profile (so it must be greater than 0.). For a constat profile, values equal to 1.0 (as here above) must be filled in.
 
 \end{itemize}
 
 \subsubsection{Time-series of sediment discharge}
 Time-series values of sediment discharge are specified in a file through the keyword \telkey{LIQUID BOUNDARIES FILE} (character type), declared in the hydrodynamic steering file. The \gaia{}'s boundary condition file must contain the flags as shown below:
 \begin{lstlisting}[frame=trBL]
-4 (*@\color{PantoneRed}5@*) 5 0.0 0.0 0.0 0.0 (*@\color{PantoneRed}4@*) 0.0 0.0 0.0 565 1
+4 (*@\color{PantoneRed}5@*) 5 1.0 0.0 0.0 0.0 (*@\color{PantoneRed}4@*) 0.0 0.0 0.0 565 1
 \end{lstlisting}
 
 The keyword \telkey{PRESCRIBED SOLID DISCHARGES} must be also included in the steering file, with an arbitrary value.
 
+\subsubsection{Repartition of imposed sediment discharge in case of graded
+distributions}
+In case of several classes of non-cohesive sediments, the sediment discharge can
+ be subdivided to the various classes through the keyword
+\telkey{CLASSES IMPOSED SOLID DISCHARGES DISTRIBUTION} (sequence of real values
+separated by
+semi-colons, one value per class whose sum is equal to 1.).
+
+In cases where this keyword is not used, the sediment discharge will be
+subdiveded among the various classes according to the sand ratio computed
+by \gaia.
+
 \subsection{Outflow boundary conditions}
 At the outflow boundary, bedload does not require any particular boundary condition.
 For this case, the flag \texttt{LIEBOR} is set \texttt{= 4} as shown in the example below:
@@ -505,10 +517,6 @@ For this case, the flag \texttt{LIEBOR} is set \texttt{= 4} as shown in the exam
 5 4 4 0.0 0.0 0.0 0.0 (*@\color{PantoneRed}4@*) 0.0 0.0 0.0 565 1
 \end{lstlisting}
 
-\begin{WarningBlock}{Note:}
-When the keyword \telkey{PRESCRIBED SOLID DISCHARGES} is used, the mass balance provided in the listing printouts information accounts for the pores $=Q_b/(1-\lambda)$, with $\lambda$ the porosity.
-\end{WarningBlock}
-
 %-------------------------------------------------------------------------------
 \subsection{Useful graphical printouts for bedload}
 %-------------------------------------------------------------------------------

documentation/gaia/user/latex/HowTos.tex

@@ -182,7 +182,7 @@ Values (real variables) can be specified as follows:
 \begin{itemize}
 \item \texttt{\textcolor{black}{EBOR:}} prescribed bed evolution
 \item \texttt{\textcolor{black}{CBOR:}} prescribed concentration
-\item \texttt{\textcolor{black}{Q2BOR:}} prescribed bedload discharge, expressed in m$^2/$s excluding voids.
+\item \texttt{\textcolor{black}{Q2BOR:}} prescribed bedload discharge, expressed in kg/(m~s).
 \end{itemize}
 
 For the particular case where a bedload solid discharge is imposed, an extra boundary condition file needs to be defined for \gaia{}. The treatment of boundary conditions for bedload and suspended sediment transport is given in \S\ref{sec:BedloadTransport} and \S\ref{chap:SuspendedSedimentTransport}, respectively.

documentation/telemac2d/user/latex/chap13.tex

@@ -196,8 +196,8 @@ Indeed the breach is represented in \telemac{2D} by a polygon (defined by users
 in the \telkey{BREACHES DATA FILE}) within the altitude of mesh nodes lowers
 with time according to the selected breach law.
 
-Current state-of-the-art on the breaching of earthen dykes due to overtopping flows
-shows that the dyke breaching expansion is progressive (non-instantaneous)
+Current state-of-the-art on the breaching of earthen dykes due to overtopping
+flows shows that the dyke breaching expansion is progressive (non-instantaneous)
 (\cite{Morris2009a,Risher2016,Wu2017,Rifai2017,Rifai2018,Rifai2019}, among
 others), following two main phases (referred to as breach formation
 and development period, respectively):
@@ -205,57 +205,107 @@ and development period, respectively):
 \item Phase 1 - Deepening and lateral widening: as the overtopping flow depth and
 velocity over the dyke increase, both breach deepening and widening are promoted
 with a shift of the breach centerline toward the channel downstream end.
-The breach sides collapse gradually. The breach expansion during this phase is fast,
+The breach sides collapse gradually. The breach expansion during this phase is
+fast,
 \item Phase 2 - Lateral widening: the main channel free surface decreases and the
 flow depth starts stabilizing at its minimum level (approaching the main channel
 critical flow depth). The breach development becomes slower, the upstream part of
 the breach stops evolving, and deepening becomes moderate tending to stabilize.
 The breach widens along the channel flow direction due to side slope failures.
 \end{itemize}
-The breach deepening (vertical incision of the breach) is faster than breach widening
-(\cite{Morris2009b,Wahl2017,Rifai2017,Rifai2018}). When the breach bottom
-reaches the foundation of the dyke or a non-erodible layer, no further deepening of the
-breach is possible and lateral widening is controlling the breach expansion until its
-stabilisation (i.e. fully formed breach, final width is reached and erosion is stopped).
-
-Selected empirical laws have been implemented in \telemac{2D} for simulating the time
-evolution of the breach expansion (widening and deepening). They are described here below. \\
+The breach deepening (vertical incision of the breach) is faster than breach
+widening (\cite{Morris2009b,Wahl2017,Rifai2017,Rifai2018}). When the breach bottom
+reaches the foundation of the dyke or a non-erodible layer, no further deepening
+of the breach is possible and lateral widening is controlling the breach
+expansion until its stabilisation (i.e. fully formed breach, final width is
+reached and erosion is stopped).
+
+Selected empirical laws have been implemented in \telemac{2D} for simulating
+the time evolution of the breach expansion (widening and deepening). They are
+described here below. \\
 As the flood period and inundation of the floodplain along with the dyke material
-characteristics impose certain limits on breach growth, the empirical laws are applied over
-a given duration. The breach expansion continues until the breach has expanded to its
-approximate maximum dimensions. Therefore the final (i.e. ultimate, maximum) breach
+characteristics impose certain limits on breach growth, the empirical laws are
+applied over a given duration.
+The breach expansion continues until the breach
+has expanded to its approximate maximum dimensions.
+
+\begin{WarningBlock}{Warning:}
+Therefore the final
+(i.e. ultimate, maximum) breach
 dimensions (width and bottom elevation), the duration of the breaching expansion
 (or duration of each phase), must be estimated outside
-of the \tel software by the user. These parameters are indeed mandatory in the
-\telkey{BREACHES DATA FILE} which contains the description and the characteristics of
+of the \tel software by the user.
+\end{WarningBlock}
+
+These parameters are indeed mandatory
+in the
+\telkey{BREACHES DATA FILE} which contains the description and the
+characteristics of
 the breaching process; it will be described here below.
-Except for the Froehlich model, the breach longitudinal shape is assumed rectangular.
+Except for the Froehlich model, the breach longitudinal shape is assumed
+rectangular.
 
-In addition to the breach expansion computation, it is important to define a criterion to start
-the breaching process. In the current release, 3 types of criteria are available:
+In addition to the breach expansion computation, it is important to define a
+criterion to start
+the breaching process. In the current release, 5 types of criteria are
+available:
 
 \begin{enumerate}
-\item at a given time. This option can be chosen filling the \telkey{BREACHES DATA FILE} with:
+\item at a given time. This option can be chosen filling the
+\telkey{BREACHES DATA FILE} with:
 \begin{lstlisting}[language=bash]
 # Option for breaching initiation
 1
 \end{lstlisting}
 
-\item when the water level above the dyke reaches a given value
-(in case of overflow). This option can be chosen filling the \telkey{BREACHES DATA FILE} with:
+\item when the water level above the dyke
+(the averaged value of water depth is computed for all points inside the
+polygone defining the dyke) reaches a given value.
+This option can be chosen filling the \telkey{BREACHES DATA FILE} with:
 \begin{lstlisting}[language=bash]
 # Option for breaching initiation
 2
 \end{lstlisting}
 
-\item when the water level at a given point reaches a certain value
-(in case of safety level). This option can be chosen filling the \telkey{BREACHES DATA FILE} with:
+\item when the water level at a given point reaches a certain value.
+This option can be chosen filling the \telkey{BREACHES DATA FILE} with:
 \begin{lstlisting}[language=bash]
 # Option for breaching initiation
 3
 \end{lstlisting}
-\end{enumerate}
 
+\item when the water level at a given point reaches a certain value and when
+the energy balance ($\Delta E$) defined as the difference between the upstream
+head (channel side, $E_{riv}$) and downstream head (floodplain side, $E_{pla}$)
+reaches a threshold value (see Figure \ref{fig:crit4+5}).
+According to the way used to compute the energy balance, two options can be chosen:
+\begin{itemize}
+ \item the hydraulic head is computed at points (mesh nodes) given by user (3
+   points are expected, see below).
+   This option can be chosen filling the \telkey{BREACHES DATA FILE} with:
+   \begin{lstlisting}[language=bash]
+   # Option for breaching initiation
+   4
+   \end{lstlisting}
+ \item the hydraulic head is automatically computed by \telemac{2D} taking an
+   averaged value of points which are located on the border of the polygon
+   used to define the breach
+   (it is assumed that the border should represent the foot of the dyke).
+   This option can be chosen filling the \telkey{BREACHES DATA FILE} with:
+   \begin{lstlisting}[language=bash]
+   # Option for breaching initiation
+   5
+   \end{lstlisting}
+\end{itemize}
+\end{enumerate}
+\begin{figure}
+\centering
+\includegraphics[width=\textwidth]{./graphics/goeury2022}
+\caption{Scheme illustrating variables used to identify the initiation of breach
+expansion on a profile across a levee (surmounted by an earth ridge). $v$ (m/s)
+is the mean flow velocity \cite{Goeury2022}.}
+\label{fig:crit4+5}
+\end{figure}
 Since release 7.0, it is possible to take into account a lateral growth
 of the breach (dyke opening by widening and deepening).
 Old breaching processes are not affected by this new feature: only the breach
@@ -366,7 +416,8 @@ meters, and $h_w$ = depth of water above the breach invert in meters.
 
 \item Verheij formula (2002) \\
 Verheij \cite{Verheij2002} provided a simple relationship between the breach
-width $B$ and time for sand and clay levees, based on field and laboratory data sets:
+width $B$ and time for sand and clay levees, based on field and laboratory data
+sets:
 \begin{itemize}
 \item for sand levees (i.e. non-cohesive dykes):
 \begin{equation}
@@ -391,22 +442,26 @@ This law can be selected filling the \telkey{BREACHES DATA FILE} with:
 8
 \end{lstlisting}
 \end{itemize}
-with $t$ = time in hours (after initiation of breaching), and $B$ = breach width in meters.
+with $t$ = time in hours (after initiation of breaching), and $B$ = breach width
+in meters.
 
 \item Verheij and Van der Knaap (2003) formula \\
 Verheij and Van der Knaap \cite{Verheij2003} improved the previous formulation
 by including the effect of the difference in water levels at both sides of the
 dyke at the breach location, and the critical flow velocity for the initiation
-erosion of the dyke material. The empirical equation in its integral form reads as:
+erosion of the dyke material. The empirical equation in its integral form reads
+as:
 \begin{equation}
 \label{eq:Verheij2003}
 \begin{array}{lc}
-B(t)=f_1\dfrac{g^{0.5}\Delta H^{1.5}}{u_c}\log\left(1+f_2\dfrac{g}{u_c}t\right) & \text{for}~t\leq T_f
+B(t)=f_1\dfrac{g^{0.5}\Delta H^{1.5}}{u_c}\log\left(1+f_2\dfrac{g}{u_c}t\right) &
+\text{for}~t\leq T_f
 \end{array}
 \end{equation}
 With $t$ = time in hours (after the initiation of breaching); $B$ = breach width
 in meters; $u_c$ = critical flow velocity for the initiation of erosion of dyke
-material (m/s); $f_1$ and $f_2$ = coefficients; $g$ = gravitational acceleration (m/s$^2$).
+material (m/s); $f_1$ and $f_2$ = coefficients; $g$ = gravitational acceleration
+(m/s$^2$).
 $\Delta H$ (m) denotes the difference in water level between the upstream and
 downstream sides of the breach. In \telemac{2D}, this difference is computed by
 considering the water head instead of water level, i.e. $\Delta H$ (m) =
@@ -418,12 +473,14 @@ integral one, which reads as:
 \begin{equation}
 \label{eq:Verheij2003:diff}
 \begin{array}{lc}
-B(t)=\dfrac{f_1 f_2}{ln 10}\dfrac{(g\Delta H)^{1.5}}{u_c^2}\dfrac{1}{1+\dfrac{f_2g}{u_c}t}\Delta t & \text{for}~t\leq T_f
+B(t)=\dfrac{f_1 f_2}{ln 10}\dfrac{(g\Delta H)^{1.5}}{u_c^2}
+\dfrac{1}{1+\dfrac{f_2g}{u_c}t}\Delta t & \text{for}~t\leq T_f
 \end{array}
 \end{equation}
 
 The suggested values and ranges have been proposed by Verheij and Van der Knaap
-\cite{Verheij2003} for coefficients $f_1$ and $f_2$ (see Table \ref{tab:coef:verheij2003})
+\cite{Verheij2003} for coefficients $f_1$ and $f_2$ (see Table
+\ref{tab:coef:verheij2003})
 (\cite{Curran2018,VanDamme2020}).
 Equation \eqref{eq:Verheij2003} contains the critical flow velocity $u_c$ for
 the surface erosion of dyke material.
@@ -444,7 +501,8 @@ $f_2$ & 0.04 & 0.01-1 \\
 \end{table}
 \begin{table}
 \centering
-\caption{Strength characteristics of various soil types \cite{Verheij2003}\cite{Verheij2009}}
+\caption{Strength characteristics of various soil types
+ \cite{Verheij2003}\cite{Verheij2009}}
 \begin{tabular}{ll}
 \hline
 Type of Soil & $u_c$ (m/s) \\
@@ -475,7 +533,8 @@ This law can be selected filling the \telkey{BREACHES DATA FILE} with:
 As mentioned previously, the breach deepening is evolving faster than the breach
 widening.
 In the \telemacsystem, the breach minimum level $Z_{b,min}$ (elevation
-of the dyke foundation, main channel bottom or of a rigid layer) is reached in a short period.
+of the dyke foundation, main channel bottom or of a rigid layer) is reached in a
+short period.
 Therefore, the time-evolution of the breach invert elevation is simulated
 according to the following linear-time progression law:
 \begin{equation}
@@ -517,7 +576,8 @@ The slope of the breach sides is time-varying.
 \begin{figure}
 \centering
 \includegraphics[width=\textwidth]{./graphics/bennani2016}
-\caption{Schematic representations of three empirical breach formation models (adapted from Bennani \cite{Bennani2016})}
+\caption{Schematic representations of three empirical breach formation models
+(adapted from Bennani \cite{Bennani2016})}
 \label{fig:froehlich}
 \end{figure}
 Froehlich \cite{Froehlich2008} used the concept expressed by Brunner
@@ -530,13 +590,16 @@ simulations as follows:
 \item the instantaneous top width of the breach is computed as:
 \begin{equation}
 \begin{array}{lll}
-B(t)=\beta(t)B_f & \text{for}~t\leq T_f &\text{with}~\beta(t)=\dfrac{1}{2}\left\{1+\sin\left[\pi\left(\dfrac{t}{T_f}-\dfrac{1}{2}\right)\right]\right\}
+B(t)=\beta(t)B_f & \text{for}~t\leq T_f &\text{with}~\beta(t)=\dfrac{1}{2}
+\left\{1+\sin\left[\pi\left(\dfrac{t}{T_f}-\dfrac{1}{2}\right)\right]\right\}
 \end{array}
 \end{equation}
 \item the instantaneous elevation of the breach bottom is calculated as follows:
 \begin{equation}
 \begin{array}{lll}
-Z_b(t)=Z_{b0}-\beta_1(t)(Z_{b0}-Z_{b,min}) & \text{for}~t\leq T_d &\text{with}~\beta_1(t)=\dfrac{1}{2}\left\{1+\sin\left[\pi\left(\dfrac{t}{T_d}-\dfrac{1}{2}\right)\right]\right\}
+Z_b(t)=Z_{b0}-\beta_1(t)(Z_{b0}-Z_{b,min}) & \text{for}~t\leq T_d &
+\text{with}~\beta_1(t)=\dfrac{1}{2}\left\{1+\sin\left[\pi\left(\dfrac{t}{T_d}
+-\dfrac{1}{2}\right)\right]\right\}
 \end{array}
 \end{equation}
 \end{itemize}
@@ -573,14 +636,19 @@ following way (the order of parameters cannot be changed):
 \begin{itemize}
 \item Number of breaches,
 \item Width of the polygon defining the breaches (in m),
-\item Option for the breaching initiation (from 1 to 3),
+\item Option for the breaching initiation (from 1 to 5),
 \item If the option of breaching initiation is at a given time (option 1):
 breach opening moment in seconds,
 \item Duration of the breaching process (in s),
 \item Option for lateral growth (from 1 to 10),
 \item Final bottom altitude of the breach (in m),
 \item If option of breaching initiation corresponds to controlling level of
-breach with global node (option 3): number of global node controlling the breach,
+breach with global node (option 3 and 5): number of global node controlling the breach,
+\item If option of breaching initiation is number 4: 3 mesh points have to be given.
+One to control the water depth, one to compute the hydraulic head on the channel and
+one to compute the hydraulic head on the floodplain,
+\item If option of breaching is number 4 or 5: threshold value for difference of hydraulic
+head (in m),
 \item If option of breaching initiation is not at given time: control level of
 the breach (in m),
 \item If the user defined formula with two growth rates is used for lateral

documentation/telemac3d/user/latex/general_setup.tex

@@ -974,7 +974,7 @@ In that example, the flow rate is prescribed at the first boundary and the free
 surface is prescribed at the second boundary.
 
 Since release 8.2, a time reference can be given:
-If a \#REFDAT with a date + hour in YYYY-MM-DD HH:MM:SS
+If a \#REFDATE with a date + hour in YYYY-MM-DD HH:MM:SS
 in year, month, day, hour, minute, second format is written in ASCII files
 related to time,
 the date+hour will be added to the times in these ASCII files

documentation/telemac3d/user/latex/physical_setup.tex

@@ -665,7 +665,8 @@ direction (resp. $x$ and $y$ components of wind velocity),
 \item TAIR for air temperature (in $^{\circ}$C),
 \item PATM for atmospheric pressure (in hPa = mbar),
 \item HREL for relative humidity (in \%),
-\item CLDC for cloud cover or nebulosity (in octas),
+\item CLDC for cloud cover or nebulosity (in octas for \waqtel or tenths for
+\khione),
 \item RAINI or RAINC for rain (in mm/s), see previous section for the meaning
 of the I and C additional charactere,
 \end{itemize}

documentation/waqtel/theory_guide/latex/module_thermic.tex

@@ -168,7 +168,7 @@ In 3D, the solar flux penetrating the water $RS$ is calculated with:
 \item the solar radiation reaching the surface with a clear sky. It depends on
 time and the location of the site. Perrin de Brichambaut's method
 \cite{perrin_ray_solaire_1963}, \cite{perrin_res_solaires_1975} is used,
-\item cloud coverage: a corrective term depending on the nebulosity is applied
+\item cloud cover: a corrective term depending on the nebulosity is applied
 (Berliand's formulae \cite{berliand_cloud_1952}, \cite{berliand_radiatssi_1960}),
 \item albedo that enables to compute the effective part of solar radiation
 penetrating the water.
@@ -312,7 +312,8 @@ with:
   \item $e_{air}$ = calibration coefficient of atmospheric radiation (1$^{\rm{st}}$ parameter of the model),
   \item $\sigma = $ Stefan-Boltzman constant, $\sigma = 5.67.10^{-8}$ Wm$^{-2}$K$^{-4}$,
   \item $T_{air}$ = air temperature ($^{\circ}$C),
-  \item $c$ = cloudiness (octa), given in the atmospheric data file
+  \item $c$ = cloudiness (octas), given in the atmospheric data file
+(WARNING: in \khione, it is given in tenths),
   \item $k$ = coefficient of cloud type, which depends on the type of clouds and their height
     (see Table \ref{tab_cloud_type}).
 \end{itemize}
@@ -366,19 +367,22 @@ To simplify calculations, an average value of $k$ = 0.2 is usually taken in 2D.\
 In 3D, clouds and albedo at
 the free surface determine the atmospheric radiation $RA$ penetrating the water:
 \begin{equation}
-RA = (1-alb_{lw}) e_{air}\sigma(T_{air}+273.15)^{4}(1+k . C^{2}),
+RA = (1-alb_{lw}) e_{air}\sigma(T_{air}+273.15)^{4}(1+k . \left( \frac{C}{8} \right)^{2}),
 \end{equation}
 where:
 \begin{itemize}
 \item $alb_{lw}$ = 0.03 is the water albedo for long radiative waves
   (common value used in the literature \cite{imerito_dyresm_2007},
   \cite{henderson-sellers_energy_balance_1986}),
+  (1-$alb_{lw}$) is a calibrating coefficient of atmospheric radiation,
 %\item $T_{air}$ ($^{\circ}$C) is the air temperature,
-\item $e_{air} = 0.937.10^{-5}(T_{air}+273.15)^{2}$ is the air emissivity,
+\item $e_{air}$ is the air emissivity (= $0.937.10^{-5}(T_{air}+273.15)^{2}$
+if using Swinbank formula, default option),
 \item $\sigma= 5.67.10^{-8}~\mathrm{{W.m^{-2}.K^{-4}}}$ is Stefan-Boltzmann's constant,
-\item $C$ is the nebulosity (tenths). Some meteorological services such as
+\item $C$ is the nebulosity (octas). Some meteorological services such as
 M\'{e}t\'{e}o France provide this data in octas, it needs to be converted
-into tenths,
+into tenths, hence the division by 8 in the formula
+(WARNING: in \khione, it is given in tenths),
 \item $k$ (dimensionless) is a parameter characterising the type of
 cloud. In practise, it is difficult to know the type of cloud during the
 period of simulation and a mean value of 0.17 is often used \cite{tva_heat_1972},
@@ -386,6 +390,15 @@ period of simulation and a mean value of 0.17 is often used \cite{tva_heat_1972}
 (see the table \ref{tab_cloud_type}).
 \end{itemize}
 
+When coupling \waqtel with \telemac{3D}, the formula to compute air emissivity
+$e_{air}$ times the factor depending on cloud cover $C$ can be chosen:
+\begin{itemize}
+\item 1: Idso and Jackson (1969),
+\item 2: Swinbank (1963),
+\item 3: Brutsaert (1975),
+\item 4: Yajima Tono Dam (2014).
+\end{itemize}
+
 The formulae in 2D and 3D are almost the same with few differences.
 
 %\begin{table}[ptbh]

documentation/waqtel/user/latex/module_thermic.tex

@@ -131,9 +131,11 @@ RA=e_{air}\sigma\left(T_{air}+273.15 \right)^4\left(1+k\left(\frac{c}{8}\right)^
 where:
 \begin{itemize}
 \item $e_{air}$ is a calibrating coefficient given by the keyword
-  \telkey{COEFFICIENTS FOR CALIBRATING ATMOSPHERIC RADIATION} (default 0.97),
+  \telkey{COEFFICIENTS FOR CALIBRATING ATMOSPHERIC RADIATION} (default = 0.97),
 \item $\sigma$ is the constant of Stefan-Boltzmann (= 5.67.10$^{-8}$ Wm$^{-2}$K$^{-4}$),
 \item $T_{air}$ is air temperature given in the \telkey{ASCII ATMOSPHERIC DATA FILE},
+\item $c$ = cloudiness (octas), given in the atmospheric data file
+(WARNING: in \khione, it is given in tenths),
 \item $k$ is the coefficient that represents the nature and elevation of clouds,
 it has a mean value of 0.2 and can be changed with the keyword
 \telkey{COEFFICIENT OF CLOUDING RATE}.
@@ -163,19 +165,24 @@ However, it varies like indicated in Table \ref{tab:kcloud}.
 In 3D, clouds and albedo at
 the free surface determine the atmospheric radiation $RA$ penetrating the water:
 \begin{equation*}
-RA = (1-alb_{lw}) e_{air}\sigma(T_{air}+273.15)^{4}(1+k . C^{2}),
+RA = (1-alb_{lw}) e_{air}\sigma(T_{air}+273.15)^{4}(1+k . \left( \frac{C}{8} \right)^{2}),
 \end{equation*}
 where:
 \begin{itemize}
 \item $alb_{lw}$ = 0.03 is the water albedo for long radiative waves
   (common value used in the literature \cite{imerito_dyresm_2007},
   \cite{henderson-sellers_energy_balance_1986}),
+  (1-$alb_{lw}$) is equal to the calibrating coefficient given by the keyword
+  \telkey{COEFFICIENTS FOR CALIBRATING ATMOSPHERIC RADIATION} (default = 0.97,
+  hence $alb_{lw}$ = 0.03 as default value),
 %\item $T_{air}$ ($^{\circ}$C) is the air temperature,
-\item $e_{air} = 0.937.10^{-5}(T_{air}+273.15)^{2}$ is the air emissivity,
+\item $e_{air}$ is the air emissivity (= $0.937.10^{-5}(T_{air}+273.15)^{2}$
+if using Swinbank formula, default option),
 \item $\sigma= 5.67.10^{-8}~\mathrm{{W.m^{-2}.K^{-4}}}$ is Stefan-Boltzmann's constant,
-\item $C$ is the nebulosity (tenths). Some meteorological services such as
+\item $C$ is the nebulosity (octas). Some meteorological services such as
 M\'{e}t\'{e}o France provide this data in octas, it needs to be converted
-into tenths,
+into tenths, hence the division by 8 in the formula,
+(WARNING: in \khione, it is given in tenths),
 \item $k$ (dimensionless) is a parameter characterising the type of
 cloud. In practise, it is difficult to know the type of cloud during the
 period of simulation and a mean value of 0.17 is often used \cite{tva_heat_1972},

documentation/waqtel/user/latex/practical_aspect.tex

@@ -73,7 +73,8 @@ be used in the headline of the \telkey{ASCII ATMOSPHERIC DATA FILE}:
 \item WINDX and WINDY: wind velocity components along $x$ and $y$ (in m/s),
 \item TAIR: air temperature (in $^{\circ}$C),
 \item PATM: atmospheric pressure (in mbar or hPa),
-\item CLDC: cloud coverage or nebulosity (in octas or tenths),
+\item CLDC: cloud cover or nebulosity (in octas for \waqtel or tenths for
+\khione),
 \item RAINI or RAINC: rainfall (last letter I or C depending on if it is
 an interpolated variable as other usual variables
 or if it is  given as cumulated variable),

examples/gaia/bosse-t2d/doc/bosse-t2d.tex

@@ -101,15 +101,15 @@ with $F(x)=c (z_b (x,t=0))$.
 
 We choose to expresse the solid discharge with the Engelund-Hansen formula:
   \begin{equation*}
-  q_s=0.05 \left[\left(s-1\right)gd^3\right]^{1/2}\frac{\theta^{\frac{5}{2}}}{c_f}
+  q_s=0.1 \left[\left(s-1\right)gd^3\right]^{1/2}\frac{\theta^{\frac{5}{2}}}{c_f}
   \end{equation*}
  with $c_f$ the adimensional friction coefficient, $\theta$ the Shields number ($\theta=\frac{\tau_b}{(\rho_s-\rho)gd}$) and $s$ the relative density of sediment.
 By explicitation of the Shields number, we obtain:
  \begin{equation*}
-  q_s=\frac{0.05}{c_f} \left[\left(s-1\right)gd^3\right]^{1/2}\left(\frac{1}{2}\frac{ c_f}{(s-1)gd} \right)^{5/2}||\vec{u}||^5
+  q_s=\frac{0.1}{c_f} \left[\left(s-1\right)gd^3\right]^{1/2}\left(\frac{1}{2}\frac{ c_f}{(s-1)gd} \right)^{5/2}||\vec{u}||^5
  \end{equation*}
 Then, the parameters to compute the bed celerity are: \\
-$A=\frac{0.05}{c_f}\left[\left( s-1\right)gd^3 \right]^{1/2}\left[\frac{c_f}{2g\left( s-1\right)d} \right]^{5/2}$ and $m=5$.
+$A=\frac{0.1}{c_f}\left[\left( s-1\right)gd^3 \right]^{1/2}\left[\frac{c_f}{2g\left( s-1\right)d} \right]^{5/2}$ and $m=5$.
 
 
 \subsection{Physical parameters}

examples/gaia/bosse-t2d/user_fortran/caracteristique.f

@@ -67,7 +67,7 @@
       D = 0.000150D0
       S = 2.65D0
       POROS = 0.375D0
-      CC_ENG = 0.05
+      CC_ENG = 0.1D0
       STRICKLER = 50.D0
       CFROT = 2.D0*GRAV/(STRICKLER**2*MOYENNE_H**(1.D0/3.D0))
       DEBIT = 0.25D0

examples/gaia/bosse-t2d/vnv_bosse-t2d.py

@@ -98,10 +98,10 @@ class VnvStudy(AbstractVnvStudy):
         # Comparison with exact solution
         self.check_epsilons('fe_seq:GAIRES','fe_seq:GAIRES',
                             var1='BOTTOM', var2='EXACT SOLUTION',
-                            eps=[0.09])
+                            eps=[0.098])
         self.check_epsilons('fv_seq:GAIRES','fv_seq:GAIRES',
                             var1='BOTTOM', var2='EXACT SOLUTION',
-                            eps=[0.07])
+                            eps=[0.092])
 
     def _post(self):
         """