From 734695299083884c54768486a0359e728f7c727c Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Chi-Tu=C3=A2n=20Pham?= <chi-tuan.pham@edf.fr>
Date: Sun, 29 Dec 2024 19:47:35 +0100
Subject: [PATCH] [VnV][telemac3d] Update vent example with better choices for
 a few keywords - TREATMENT OF NEGATIVE DEPTHS = 2 (new default value since
 release 9.0) and TIDAL FLATS = YES (default value) to improve mass
 conservation of water - use of gradient conjugate + preconditioning 34 = 2*17
 to solve PPE, more efficient than GMRES for this example - accuracy for
 propagation decreased to 1.E-15 (improved differences between sequential and
 parallel runs) - switch to non-hydrostatic version

---
 documentation/data/biblio.bib        | 10 +++++
 examples/telemac3d/vent/doc/vent.tex | 63 ++++++++++++++++------------
 examples/telemac3d/vent/f3d_vent.slf |  2 +-
 examples/telemac3d/vent/t3d_vent.cas | 10 ++---
 examples/telemac3d/vent/vnv_vent.py  |  6 +--
 5 files changed, 54 insertions(+), 37 deletions(-)

diff --git a/documentation/data/biblio.bib b/documentation/data/biblio.bib
index 5a439af8a4..186d1dea90 100644
--- a/documentation/data/biblio.bib
+++ b/documentation/data/biblio.bib
@@ -40,6 +40,16 @@
   number = 	 {8},
   pages = 	 {1113-1134}
 }
+
+@Article{Wu1973,
+  author = 	 {Wu J.},
+  title = 	 {Prediction of near-surface drift currents from wind velocity},
+  journal = 	 {Journal of Hydraulic Division},
+  year = 	 {1973},
+  volume = 	 {99},
+  pages = 	 {1291-1302}
+}
+
 @Article{Smagorinsky1963,
   author = 	 {SMAGORINSKY J.},
   title = 	 {General simulation experiments with the primitive equations},
diff --git a/examples/telemac3d/vent/doc/vent.tex b/examples/telemac3d/vent/doc/vent.tex
index a9b25e2b9e..322b81ae0b 100644
--- a/examples/telemac3d/vent/doc/vent.tex
+++ b/examples/telemac3d/vent/doc/vent.tex
@@ -16,7 +16,7 @@ The geometry dimensions of rectangular channel is 500~m long and 100~m wide,
  with horizontal bed at depth -10~m. At the initial state, the channel is submitted
  to a constant 10~$\text{m}\cdot\text{s}^{-1}$ wind. \\
 The wind generates a slope of the free surface and a vertical two-dimensional circulation.
-Tsanis [1] has made an inventory of existing laboratory or in-situ
+Tsanis \cite{Tsanis1989} has made an inventory of existing laboratory or in-situ
 measurements and has plotted these values on a non-dimensional graph.
 He deduced a characteristic vertical velocity profile
 (Figure \ref{t3d:vent:fig:mixinglength}), which will be compared to
@@ -82,6 +82,9 @@ The mesh (Figures \ref{t3d:vent:fig:meshH}, \ref{t3d:vent:fig:meshV} and
 \ref{t3d:vent:fig:meshV_Zoom})  is composed of 543 triangular elements
 (315 nodes) with 15 planes  irregularly spaced on the vertical, to form prism elements.
 
+Note that the quality of this mesh can be improved by changing overconstrained
+triangles at each corner.
+
 \begin{figure}[!htbp]
  \centering
  \includegraphicsmaybe{[width=\textwidth]}{../img/Mesh.png}
@@ -107,27 +110,37 @@ The mesh (Figures \ref{t3d:vent:fig:meshH}, \ref{t3d:vent:fig:meshV} and
 The time step is 10~s for a simulated period of 20,000~s (5~h 33~min 20~s).
 
 \bigskip
-This case is computed with the hydrostatic pressure hypothesis. To solve advection,
-the LIPS scheme is used for the velocities
+This case is computed with the non-hydrostatic pressure hypothesis.
+To solve advection, the LIPS scheme is used for the velocities
 and Diagonal preconditioning is used for propagation (default option).
-The resolution accuracy is $10^{-8}$ for the diffusivity velocity
-and $10^{-8}$ for the propagation (default values). The implicitation coefficients
-for depth and velocities are equal to 0.6.
+
+The resolution accuracy is $10^{-8}$ for the diffusivity velocity (default
+value), so is for PPE whereas it is set to $10^{-15}$ for the propagation
+to have less differences between sequential and parallel runs.
+To solve the Poisson Pressure equation, conjugate gradient combined with
+a preconditioning equal to 34 are used as more efficient than GMRES.
+
+The implicitation coefficients for depth and velocities are equal to 0.6.
+\\
+
+Mass conservation is improved (up to machine precision) by using default value
+for the keyword \telkey{TREATMENT OF NEGATIVE DEPTHS} (= 2), instead of the old
+default value = 1 (smoothings) for this example (until release 9.0).
 
 \section{Results}
 \bigskip
 The mass balance is the following:
 \begin{lstlisting}[language=TelFortran]
-   INITIAL MASS                        :     500000.0
-   FINAL MASS                          :     500000.0
-   MASS LEAVING THE DOMAIN (OR SOURCE) :     0.000000
-   MASS LOSS                           :   -0.2517768E-03
+   INITIAL VOLUME                      :     500000.0
+   FINAL VOLUME                        :     500000.0
+   VOLUME EXITING (BOUNDARY OR SOURCE) :     0.000000
+   TOTAL VOLUME LOST                   :   -0.2910383E-09
 \end{lstlisting}
 
 \bigskip
-Thus, the mass balance is consistent with the accuracy asked
-(from 10$^{-8}$ to 10$^{-6}$ depending on the system to solve,
-for the diffusion of velocities or propagation).
+Thus, the mass balance is perfect (improved from release 9.0 by using
+default choice for \telkey{TREATMENT OF NEGATIVE DEPTHS} = 2).
+
 Figure \ref{t3d:vent:fig:vecVelo} shows the vertical circulation induced by the wind.
 \begin{figure}[!htbp]
  \centering
@@ -143,15 +156,18 @@ Figure \ref{t3d:vent:fig:vecVelo} shows the vertical circulation induced by the
 \end{figure}
 
 \bigskip
-The velocity at the surface is 0.18~$\text{m}\cdot \text{s}^{-1}$ and the return current reaches %0.177587
-a maximum value of 0.058~$\text{m}\cdot \text{s}^{-1}$ (Figure \ref{t3d:vent:fig:fieldvelo}). %0.05798
+The velocity at the surface is 0.18~$\text{m}\cdot \text{s}^{-1}$ and the return
+%Read USURF or MSURF in the 2D RESULT FILE: 0.177586
+current reaches a maximum value of 0.058~$\text{m}\cdot \text{s}^{-1}$ (Figure
+\ref{t3d:vent:fig:fieldvelo}). %0.05798
 However, it must be pointed out that the velocity at the surface depends
 on the refinement near the surface because the velocity gradient is
 very high in this area.
 A distance between the two first vertical points of 0.50~m instead of
 0.10~m induced a velocity at the surface of 0.12~$\text{m}\cdot \text{s}^{-1}$, but the velocity
 field below 1~m under the surface was only slightly modified.
-J. Wu [2] proposed for the velocity $\vec{u_s}$ us at the surface the expression:
+J. Wu \cite{Wu1973} proposed for the velocity $\vec{u_s}$ us at the surface the
+expression:
 \begin{equation*}
 \vec{u_s} = 0.55 \left(\frac{\vec{\tau}}{\rho_{air}}\right)^{1/2}.
 \end{equation*}
@@ -172,7 +188,8 @@ orientation as the wind, is close to measured profiles, but the lower
 part is smoothed.
 This could mean that the turbulence intensity is stronger in nature.
 The slope of the free surface presented in figure \ref{t3d:vent:fig:FreeSurface}
-is equal to $1.58\cdot 10^{-6}$ (between $x=100$~m and $x=400$~m). %1.57663
+is equal to $1.58\cdot 10^{-6}$ (between $x=100$~m and $x=400$~m). %1.5757212
+% = (0.00022398578+0.00024873059)/(400-100)
 The computation of the slope, assuming that the flow is homogeneous on the vertical,
 gives a slope equal to $1.66\cdot 10^{-6}$.
 This value is in agreement with the value given by \telemac{3D}.
@@ -184,20 +201,12 @@ This value is in agreement with the value given by \telemac{3D}.
  \label{t3d:vent:fig:FreeSurface}
 \end{figure}
 
-%\section{Conclusion}
+\section{Conclusion}
 
-\bigskip
+%\bigskip
 The velocity field produced by \telemac{3D} using the standard mixing length is correct.
 Near the surface, the quality of the results depends on the vertical resolution near the surface.
 The second level of the vertical mesh should be fixed 0.10~m below the surface for a good result.
 In the deeper part, the profile is a little bit more smoothed.
 Taking into account the effect of the wind directly into the turbulence model
 may improve this result.
-
-\section{Reference}
-
-[1] TSANIS I., Simulation of wind-induced water currents,
- Journal of Hydraulic Engineering, Vol.115, n 8, 1989, pp 1113-1134.\\
-
-[2] WU J., Prediction of near-surface drift currents from wind velocity,
- Journal of Hydraulic Division, 99, 1973, pp 1291-1302.
diff --git a/examples/telemac3d/vent/f3d_vent.slf b/examples/telemac3d/vent/f3d_vent.slf
index 5feb218af2..0c137e655b 100644
--- a/examples/telemac3d/vent/f3d_vent.slf
+++ b/examples/telemac3d/vent/f3d_vent.slf
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:8c3146133f1e41fbc0b84b3554faa6cb1a3f46395e24312843c32aa1a184a3e3
+oid sha256:5c49b42d3a29d1bfcadc81ed7bff722202739ab805fb6cedb1c5b9a248c29954
 size 334140
diff --git a/examples/telemac3d/vent/t3d_vent.cas b/examples/telemac3d/vent/t3d_vent.cas
index 2db78f69fd..cc65530498 100644
--- a/examples/telemac3d/vent/t3d_vent.cas
+++ b/examples/telemac3d/vent/t3d_vent.cas
@@ -21,7 +21,7 @@ MESH TRANSFORMATION : 2
 /
 TITLE : 'TEST CANAL Soumis a un vent'
 /
-VARIABLES FOR 2D GRAPHIC PRINTOUTS : 'U,V,H,B'
+VARIABLES FOR 2D GRAPHIC PRINTOUTS : 'U,V,H,B,S,USURF,VSURF,WSURF,MSURF'
 VARIABLES FOR 3D GRAPHIC PRINTOUTS : 'Z,U,V,W,NUZ'
 TIME STEP : 10.
 NUMBER OF TIME STEPS : 2000
@@ -54,10 +54,11 @@ COEFFICIENT OF WIND INFLUENCE : 1.625E-6
 WIND VELOCITY ALONG X         : 10.
 WIND VELOCITY ALONG Y         :  0.
 /----------------------------------------------------------------------
-ACCURACY FOR PROPAGATION : 1.E-8
+ACCURACY FOR PROPAGATION : 1.E-15
+SOLVER FOR PPE   : 1
+PRECONDITIONING FOR PPE : 34
 /
 / DEFAULT VALUES UNTIL V7P3 KEPT FOR NON REGRESSION
-NON-HYDROSTATIC VERSION : NO
 LAW OF BOTTOM FRICTION = 2
 /
 / DEFAULT VALUE UNTIL V8P0 KEPT FOR NON REGRESSION
@@ -65,6 +66,3 @@ FRICTION COEFFICIENT FOR THE BOTTOM = 60.
 /
 / DEFAULT VALUE UNTIL V8P1 KEPT FOR NON REGRESSION
 COEFFICIENT OF WIND INFLUENCE VARYING WITH WIND SPEED = NO
-/
-/ DEFAULT VALUE UNTIL V8P5 KEPT FOR NON REGRESSION
-TREATMENT OF NEGATIVE DEPTHS = 1
diff --git a/examples/telemac3d/vent/vnv_vent.py b/examples/telemac3d/vent/vnv_vent.py
index bad6b2e51c..ed45ba6ca0 100644
--- a/examples/telemac3d/vent/vnv_vent.py
+++ b/examples/telemac3d/vent/vnv_vent.py
@@ -55,17 +55,17 @@ class VnvStudy(AbstractVnvStudy):
         # Comparison with the last time frame of the reference file.
         self.check_epsilons('vnv_1:T3DRES',
                             'f3d_vent.slf',
-                            eps=[1.E-6, 1.E-7, 1.E-7, 1.E-8, 1.E-8])
+                            eps=[1.E-11])
 
         # Comparison with the last time frame of the reference file.
         self.check_epsilons('vnv_2:T3DRES',
                             'f3d_vent.slf',
-                            eps=[1.E-6, 1.E-7, 1.E-7, 1.E-8, 1.E-8])
+                            eps=[1.E-11])
 
         # Comparison between sequential and parallel run.
         self.check_epsilons('vnv_1:T3DRES',
                             'vnv_2:T3DRES',
-                            eps=[1.E-6, 1.E-7, 1.E-8, 1.E-8, 1.E-8])
+                            eps=[1.E-11])
 
 
     def _post(self):
-- 
GitLab