From d8409c3d08172ba7803448d9b5fb6142e46e4c8b Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Chi-Tu=C3=A2n=20Pham?= <chi-tuan.pham@edf.fr>
Date: Wed, 8 Jan 2025 12:22:38 +0100
Subject: [PATCH] [VnV][waqtel] Update heat_exchange example with better
 choices for some keywords - TREATMENT OF NEGATIVE DEPTHS = 2 (new default
 value since release 9.0) to improve mass conservation (water depth +
 temperature, water depth up to machine accuracy) - SCHEME FOR ADVECTION OF
 VELOCITIES = 5 (PSI) rather than 1 (characteristics) for better mass
 conservation + reduced differences between seq/par runs - use of gradient
 conjugate + preconditioning 34 = 2*17 to solve PPE, more efficient than GMRES
 for this example - accuracies for diffusion of velocities/tracers/k-epsilon +
 PPE set to 1.E-15 as easily and efficiently reached with this set-up

---
 .../heat_exchange/doc/heat_exchange.tex       | 84 +++++++++----------
 .../f3d_heat_exchange_rain_wind.slf           |  2 +-
 .../t3d_heat_exchange_rain_wind.cas           | 14 ++--
 .../waqtel/heat_exchange/vnv_heat_exchange.py |  6 +-
 4 files changed, 50 insertions(+), 56 deletions(-)

diff --git a/examples/waqtel/heat_exchange/doc/heat_exchange.tex b/examples/waqtel/heat_exchange/doc/heat_exchange.tex
index 42f8b412f4..734b55a4c0 100644
--- a/examples/waqtel/heat_exchange/doc/heat_exchange.tex
+++ b/examples/waqtel/heat_exchange/doc/heat_exchange.tex
@@ -1,24 +1,19 @@
-% case name
 \chapter{heat\_exchange}
-%
-% - Purpose & Description:
-%     These first two parts give reader short details about the test case,
-%     the physical phenomena involved, the geometry and specify how the numerical solution will be validated
-%
+
 \section{Purpose}
-%
+
 This case is an example of the use of the heat exchange with atmosphere module
 of \waqtel coupled with \telemac{3D}.
-%
+
 \section{Description}
-%
+
 A square basin at rest is considered (length and width = 200~m)
 with flat bathymetry and elevation at 0~m.
-%
+
 \section{Computational options}
-%
+
 \subsection{Mesh}
-%
+
 The triangular mesh is composed of 272 triangular elements and 159 nodes
 (see Figure \ref{fig:heat_exchange:mesh}).
 
@@ -39,14 +34,14 @@ seen on Figure \ref{fig:heat_exchange:mesh_section}.
 \caption{Vertical mesh}
  \label{fig:heat_exchange:mesh_section}
 \end{figure}
-%
+
 \subsection{Physical parameters}
-%
+
 The heat exchange module is activated by setting \telkey{WATER QUALITY PROCESS}
 = 11 in the \telemac{3D} \telkey{STEERING FILE}.\\
 
 The tracer used is temperature.
-%
+
 The density is considered as a function of the water temperature:
 \begin{equation}
   \rho(T) = \rho_0 \left(1 - \alpha (T-T_0)^2\right),
@@ -125,48 +120,53 @@ user manual).\\
 The $k$-$\epsilon$ turbulence model is used for the vertical direction whereas
 a constant viscosity is used for the horizontal directions with the default value
 of 10$^{-6}$~m$^2$/s.
-%
+
 \subsection{Initial and Boundary Conditions}
-%
+
 The initial water depth is 5~m with a fluid at rest.\\
 The initial temperature is 20~$^{\circ}$C in the whole domain.\\
 
 There are only closed lateral boundaries with free slip condition and
 Nikuradse law is used to model friction at the bottom with a coefficient
 0.001~m.
-%
+
 \subsection{General parameters}
-%
+
 The time step is 2~min = 120~s for a simulated period of 30~days
 (= 2,592,000~s).
-%
+
 \subsection{Numerical parameters}
-%
+
 The non-hydrostatic version of \telemac{3D} is used.
 It is mandatory to model stratifications correctly.
 
-To solve the advection steps, the method of characteristics is chosen for the
-velocities + $k$-$\epsilon$ and the PSI-type MURD scheme for the temperature
-for CPU time reasons.
+To solve the advection steps, the PSI-type MURD scheme is chosen for the
+temperature and velocities whereas the method of characteristics is chosen for
+$k$-$\epsilon$ for CPU time reasons mainly.
+For velocities, PSI-type MURD scheme is a better choice than the method of
+characteristics for mass conservation and also to reduce the differencies
+between sequential and parallel runs.
 
-To reach convergence when solving linear systems, the keywords
-\telkey{MASS-LUMPING FOR DIFFUSION} has to be set to 1. (default = 0.) whereas
-\telkey{MAXIMUM NUMBER OF ITERATIONS FOR PPE} has to be set to 200 (Poisson
-Pressure Equation).
+To reach convergence when solving linear systems, the keyword
+\telkey{MASS-LUMPING FOR DIFFUSION} has to be set to 1. (default = 0.).
 
 To accelerate solving, preconditioning 34 (= 2 $\times$ 17) is used for
 every conjugate gradient solver for diffusion step (velocities, tracers and
 $k$-$\epsilon$).
+Moreover, to accelerate the solving of Poisson Pressure equation, conjugate
+gradient and preconditioning 34 are used.
+As the solving of these 4 operations is easy and efficient, the accuracies
+asked can be set to 10$^{-15}$.
+
+To enable not too long computations, \telkey{ACCURACY FOR PROPAGATION} is
+relaxed to 10$^{-6}$ (default = 10$^{-8}$).
+
+\telkey{TREATMENT OF NEGATIVE DEPTHS} is set to default value 2
+(flux control) to improve mass conservation (water depth + temperature,
+up to machine accuracy for water depth).
 
-To enable not too long computations, the following parameters have been relaxed:
-\telkey{ACCURACY FOR PPE} = 10$^{-4}$ and \telkey{ACCURACY FOR PROPAGATION} = 10$^{-6}$
-(default = 10$^{-8}$ for both of them).
-%
-%     We comment in this part the numerical results against the reference ones,
-%     giving understanding keys and making assumptions when necessary.
-%
 \section{Results}
-%
+
 Figure \ref{fig:heat_exchange:res_evol} shows the temperature evolution along
 time at 3 locations (top, middle and bottom of the water column at the center of
 the square, same $x$ and $y$).
@@ -204,16 +204,8 @@ Figure \ref{fig:heat_exchange:res_temp_section} at final time.
 \caption{Vertical distribution of temperature at final time}
 \label{fig:heat_exchange:res_temp_section}
 \end{figure}
-%
+
 \section{Conclusion}
-%
+
 \waqtel is able to model heat exchange with atmosphere phenomena when coupled
 with \telemac{3D}.
-%
-% Here is an example of how to include the graph generated by validate_telemac.py
-% They should be in test_case/img
-%\begin{figure} [!h]
-%\centering
-%\includegraphics[scale=0.3]{../img/mygraph.png}
-% \caption{mycaption}\label{mylabel}
-%\end{figure}
diff --git a/examples/waqtel/heat_exchange/f3d_heat_exchange_rain_wind.slf b/examples/waqtel/heat_exchange/f3d_heat_exchange_rain_wind.slf
index 26b62417b2..7164d09635 100644
--- a/examples/waqtel/heat_exchange/f3d_heat_exchange_rain_wind.slf
+++ b/examples/waqtel/heat_exchange/f3d_heat_exchange_rain_wind.slf
@@ -1,3 +1,3 @@
 version https://git-lfs.github.com/spec/v1
-oid sha256:72d021c9d1a18cacb0a1e3ddf33f7ae761d0ddc2018159227869d72db6040299
+oid sha256:837dc86692d30ee34329474642ce962196cad0e2613e8bf32efd2ab4b27e2f41
 size 295980
diff --git a/examples/waqtel/heat_exchange/t3d_heat_exchange_rain_wind.cas b/examples/waqtel/heat_exchange/t3d_heat_exchange_rain_wind.cas
index add786f2ce..ea9acca2db 100644
--- a/examples/waqtel/heat_exchange/t3d_heat_exchange_rain_wind.cas
+++ b/examples/waqtel/heat_exchange/t3d_heat_exchange_rain_wind.cas
@@ -45,14 +45,19 @@ HORIZONTAL TURBULENCE MODEL           : 1        / viscosite constante
 /----------------------------
 /SOLVER FOR PROPAGATION = 2
 /
+SOLVER FOR PPE = 1
+PRECONDITIONING FOR PPE = 34
+/
 PRECONDITIONING FOR DIFFUSION OF VELOCITIES = 34
 PRECONDITIONING FOR DIFFUSION OF TRACERS    = 34
 PRECONDITIONING FOR DIFFUSION OF K-EPSILON  = 34
 /
-ACCURACY FOR PPE                     : 1.E-4
+ACCURACY FOR PPE                     : 1.E-15
 ACCURACY FOR PROPAGATION             : 1.E-6
+ACCURACY FOR DIFFUSION OF VELOCITIES : 1.E-15
+ACCURACY FOR DIFFUSION OF TRACERS    : 1.E-15
+ACCURACY FOR DIFFUSION OF K-EPSILON  : 1.E-15
 /
-MAXIMUM NUMBER OF ITERATIONS FOR PPE = 200
 /----------------------------
 /     SOURCE TERMS
 /----------------------------
@@ -84,12 +89,9 @@ COEFFICIENT OF WIND INFLUENCE VARYING WITH WIND SPEED = YES
 /----------------------------
 NON-HYDROSTATIC VERSION = YES
 / METHOD OF CHARACTERISTICS TO ACCELERATE THE COMPUTATION
-SCHEME FOR ADVECTION OF VELOCITIES = 1
+SCHEME FOR ADVECTION OF VELOCITIES = 5
 SCHEME FOR ADVECTION OF K-EPSILON  = 1
 SCHEME OPTION FOR ADVECTION OF TRACERS = 1 / PSI LIGHTLY FASTER THAN LIPS HERE
 / MORE STABILITY
 FREE SURFACE GRADIENT COMPATIBILITY : 0.9
 MASS-LUMPING FOR DIFFUSION   : 1.
-/
-/ DEFAULT VALUE UNTIL V8P5 KEPT FOR NON REGRESSION
-TREATMENT OF NEGATIVE DEPTHS = 1
diff --git a/examples/waqtel/heat_exchange/vnv_heat_exchange.py b/examples/waqtel/heat_exchange/vnv_heat_exchange.py
index 8c565673b9..cf0b416015 100644
--- a/examples/waqtel/heat_exchange/vnv_heat_exchange.py
+++ b/examples/waqtel/heat_exchange/vnv_heat_exchange.py
@@ -54,17 +54,17 @@ class VnvStudy(AbstractVnvStudy):
         # Comparison with the last time frame of the reference file.
         self.check_epsilons('vnv_3:T3DRES',
                             'f3d_heat_exchange_rain_wind.slf',
-                            eps=[6.E-5, 0.007, 0.014, 1.E-3, 0.22])
+                            eps=[4e-5, 5e-3, 4e-3, 2e-4, 5.9e-2])
 
         # Comparison with the last time frame of the reference file.
         self.check_epsilons('vnv_4:T3DRES',
                             'f3d_heat_exchange_rain_wind.slf',
-                            eps=[7.E-5, 0.013, 0.013, 1.E-3, 0.67])
+                            eps=[4e-5, 5e-3, 4e-3, 2e-4, 5.9e-2])
 
         # Comparison between sequential and parallel run.
         self.check_epsilons('vnv_3:T3DRES',
                             'vnv_4:T3DRES',
-                            eps=[6.E-5, 0.013, 0.015, 1.E-3, 0.66])
+                            eps=[4e-5, 5e-3, 5e-3, 2e-4, 5.9e-2])
 
 
     def _post(self):
-- 
GitLab