\documentclass[a4paper,12pt]{article}

\usepackage[utf8]{inputenc}
\usepackage{graphicx}
%%\def\eqref#1{Eq.~\ref{#1}}
\usepackage{subcaption}
\usepackage{a4wide}
\usepackage{natbib}
\bibliographystyle{apalike}
%\bibliographystyle{jas}
\usepackage{amsmath}
\usepackage{tabularx}
\usepackage{color}
\def\fig#1{Fig.~\ref{#1}}

\begin{document}

\def\pupper{h_{m}}  % pupper = altitude maximale de subsidence des masses d'air sec
\def\ptop{h_{wk}}   % ptop = lesommet de la poche
\def\sigmawk{k_{twk}}
\def\remarque#1{{\bf\textcolor{green}{#1}}}
\newcommand{\afaire}[1]{\textbf{\textcolor{red}{#1}}}


\title{Wakes tbd}
\maketitle

\section{Abstract}
\remarque{Pas de citations dans l'abstract}
Cold pools, formed under clouds by the evaporation of precipitation, play a central role in maintaining and organizing atmospheric convection. Their absence in climate models leads to significant errors in the representation of convection, such as premature convection onset.
The introduction of a cold pools parameterization proposed by into the LMDZ climate model has significantly improved the representation of convection, in particular with regard to its diurnal cycle. However, this parameterization had never been finely evaluated before in terms of representation of the cold pools themselves. This work provides for the first time such an evaluation based on Large Eddy Simulation (LES). We evaluate the underlying physics of the model, its internal variables as well as those used in the coupling with the deep convection scheme. The LES analyses demonstrate the relevance of the assumptions underlying the parameterization. We also show that its initial version represents well, at first order, the characteristics of the cold pools, although some biases were identified. These were corrected thanks to substantial modifications made to the cold pools scheme and a readjustment of some free parameters. Persistent defects can be corrected by adding thermal-related mixing in the cold pools and by considering  evolution of their density in a more physical way.

\newpage
\section{Introduction}


     During thunderstorms, a significant amount of precipitation evaporates before reaching the ground, generating cold air masses in the layers below the clouds. This cooled air, denser than its surroundings, collapses and then spreads horizontally across the surface, forming what are called cold pools. These are often associated with a gust front, capable of lifting the surrounding warm air and thus promoting the development of new convective cells. In organized propagative systems such as squall lines, convective columns are permanently generated by cold pool fronts at the front of the system \citep{rotunno1988,weisman2004}. When the cold pools is accompanied by a gust front, it is called a density current. These density currents are fueled by precipitating downdrafts, which is their main dynamic driver. Although present over both continents and oceans, density currents are generally deeper, colder, and propagate more rapidly over continents. They play a key role in the self-aggregation of tropical convection \citep{jeevanjee2013}, as well as in the transition between shallow and deep convection \citep{khairoutdinov2006,boing2012}.


In atmospheric Global Circulation Models (GCMs), as those used for climate change studies, convection has to parameterized due to the coarse horizontal resolution (30 to 300~km).
Simulating convective rainfall with parameterized physics is challenging \citep{randall2003}.
GCMs often underestimate rainfall rates \citep{kendon2012,pantillon2015,tan2018} and produce peak precipitation at noon, in phase with insolation, while the maximum precipitation is generally observed in late afternoon or during night \citep{randall2003,guichard2004,stephens2010,dirmeyer2012}.
Density currents probably play a key role in this timing, by self-maintaining convection \citep{pantillon2015,grant2018}.
One of the first attempts to parameterize density currents was proposed by \cite{qian1998}.
Later on, \cite{GL10I} proposed a parameterization based on a population of identical circular density currents that are cooled by convective precipitation. The coupling of the \cite{emanuel1991} parameterization of deep convection with this cold pool parameterization and with the thermal plume
model of \cite{rio2008} in the LMDZ climate model significantly improved the simulation of the diurnal cycle of precipitation in the tropics \citep{rio2009}, shifting its maximum from noon to mid afternoon. A further improvement was brought by the introduction of the stochastic triggering of deep convection \cite{rochetin2014} which made the simulated convection more intermittent.
Despite this success, and the use of the cold pool model in the standard version of the LMDZ atmospheric and IPSL (Institut Pierre Siméon Laplace) coupled models \citep{Hourdin2020,boucher2020}, it was not evaluated in details so far. This is explained not only by a lack of observational data but also by the fact that the internal variables of parameterizations are not directly accessible from observations.


     Large Eddy Simulations (LES) are a useful complement to observations. Their fine horizontal resolution enables them to simulate explicitly turbulent and convective motions in the boundary layer \citep{brown2002,siebesma2003}. One advantage of LES compared to observations is that they provide full three-dimensional information. They have been used extensively to develop and evaluate boundary layer and convection parameterizations \citep{rio2010,dorrestijn2013,strauss2019,legay2025}. LES have been used also to simulate, understand and develop parameterizations of cold pool parameterizations \citep{tompkins2001,khairoutdinov2006,couvreux2012,feng2015,kurowski2018}. However, their use for a cold pool parameterization assessment remains unexplored.


     Here we propose to use LES to evaluate in details the parameterization of cold pools of LMDZ \citep{GL10I,GL10II}. We first use LES to evaluate some of the fundamental relationships between large scale state variables (for LES, the horizontal average over the domain) and internal variables which are at the basis of the parameterization. We then propose improvements which are further assessed in simulations with a Single-Column-Model (SCM) version of LMDZ against LES.
In such simulations, the parameterization interact with all the other parameterizations and depend on the values of a number of free parameters. To explore the sensitivity of the results to those free parameters and retune the model after improvement of its physical consent, we use a tool for automatic calibration,  High-Tune-Explorer, developed recently \citep{couvreux2021,hourdin2021}. This tool, based on history matching, can be used to characterize the subspace of parameter values for which the model is in agreement with LES, given a series of target metrics and associated tolerance to error \citep{couvreux2021}. It is used here to explore the sensitivity of the agreement between SCM simulations and LES to the model free parameters. 


%\afaire{Ne pas oublier de reprendre en fonction du plan final}
     The paper starts by presenting in the section 1 the tools used : the LMDZ model, the cold pool parameterization by \cite{GL10I} (referred to as the GL10 hereafter), and the LES used for evaluation. The presentation of the tuning tool (largely published) and the setup of its use is let to an appendix to concentrate on model physics and improvement in the core of the paper. In section 2, we detail the cold pool sampling in LES, designed to assess the physical laws internal to the cold pool parameterization and its coupling with deep convection. Section 3 is devoted to a comparison of cold pool model variables simulated by LMDZ in SCM mode and those calculated in LES, in order to identify the model's limitations. These results will then be discussed, and proposed improvements will be detailed in section 4. Finally, we conclude with a synthesis and discussion of prospects in section 5.
     

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\section{Tools and methods}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\subsection{LMDZ and its single-column version}

     LMDZ is the General Circulation Model (GCM) used in this work. Developed in the 1970s at the Laboratoire de Meteorologie Dynamique \citep{sadourny1984,hourdin2006}, the “Z” in LMDZ refers to the model's ability to refine its horizontal grid over a specific region. This climate model is based on simplified Navier-Stokes equations for fluid mechanics, as well as transport equations. It represents the second generation \citep{hourdin2013} of a climate model initially described by Sadourny and Laval (1984). LMDZ is the atmospheric component of the IPSL coupled model. The latter is one of around twenty coupled models taking part in major international model intercomparison exercises, such as those of the CMIP (Coupled Model Intercomparison Project), the results of which are used in IPCC (Intergovernmental Panel on Climate Change) reports.

     LMDZ consists of two main parts, from a physical, mathematical and computational point of view. The first part, called “the dynamics”, concerns the numerical resolution of the atmospheric general circulation equations. This component manages horizontal exchanges between the model's grid cells. The second part, called “physics”, calculates the impact of radiation, small-scale processes (subgrid) and phase changes of water on dynamic variables via “physical parameterizations”. This “physical” part is made up of juxtaposed atmospheric columns, which do not interact with each other. Within each column, the variables are assumed to be statistically homogeneous in the horizontal plane.
\remarque{Le 1D/LES n'est pas une obligation et n'est vraiment utilisé que pour les paramétrisations de la convection et des nuages}

    The SCM version of LMDZ is built by extracting an atmospheric column from the GCM, incorporating all subgrid-scale parameterizations, and running it in a large-scale constrained environment. This approach has become central in the development and tuning of parameterizations of convection and associated clouds in several climate modeling groups. Parameterizations are often developed and evaluated within this single-column framework by comparing them with LES of the same atmospheric column. The SCM/LES approach was promoted in particular by GCSS (GEWEX Cloud Systems Study), a program aimed at improving the parameterization of cloud systems in climate models ($\afaire(ref)$). A major advantage of the SCM is its low computational cost, which allows a large number of simulations, even on a laptop, making it particularly useful in the development phase, where extensive testing is required.
      
\remarque{Pas possible de partir à ce point bille en tête je trouve, sans avoir un peu dressé le portrait des paramétrisations de la convection dans le modèle}.

The parameterization of turbulence, convection and clouds in LMDZ is based on a multi-scale, or object view.
The small scale turbulence, mainly active near the surface, is accounted for following \cite{yamada1983} scheme, with an eddy diffusive approach in which the eddy diffusivity relies on a prognostic equation for the turbulent kinetic energy.
A specific mass flux parameterization, the thermal plume model, accounts for the vertical transport by organised thermal plume, cells or rolls in the convective boundary layer \cite{hourdin2002,rio2008}. It was later on coupled to a bigaussian statistical representation of the subgrid distribution of humidity, leading to a strong improvement in the representation of cumulus and stratocumulus clouds \cite{jam2013,hourdin2019}.
Deep convection is represented with a modified version of the \cite{emanuel1991} scheme.

\remarque{En mettre suffisamment ici pour introduire déjà les thermes Q1 et Q2, unsat n'co. Pour les avoir sous la main plus loin}

The main modification of the deep convection scheme concerns the triggering and closure ...
\remarque{Compléter}
\remarque{Les phrases ci-dessous sont reprises de l'introduction. Il faut peut être les enlever plus haut ?}
proposed a parameterization based on a population of identical circular density currents that are cooled by convective precipitation.
Later on, \cite{GL10I}
model of \cite{rio2008} in the LMDZ climate model significantly improved the simulation of the diurnal cycle of precipitation in the tropics \citep{rio2009}, shifting its maximum from noon to mid afternoon. A further improvement was brought by the introduction of the stochastic triggering of deep convection \cite{rochetin2014} which made the simulated convection more intermittent.



\subsection{The cold pool model}

    The cold pool model represents a population of identical circular cold pools (the wakes) over an infinite plane containing the grid cell. All the wakes have the same height, radius, and vertical profiles of thermodynamic variables. Their centers are statistically distributed with a uniform density $D_{wk}$. Cold pools divide the space into two parts : (i) the interior of cold pools (w) is where convective precipitating downdrafts are located; in these downdrafts, the re- evaporation of precipitation generates intense cooling and strong negative buoyancy; (ii) the exterior of cold pools (x) contains the warm air that fuels the saturated convective currents (\fig{fig:poche}). The top height of the cold pool ($\ptop$) is defined as the altitude where the temperature difference between (w) and (x) becomes zero. Below this level cold pools are cooler than their exterior: they collapse and spread out as they are denser than their surrounding. The boundary between the cold pool and the environment is considered to be infinitely thin, and at each point on this boundary, the cold pool spreads at a rate $C$. $C$ is considered to be a random variable whose mean $C_{*}$ will give the rate at which the cold pool spreads. In the GL10 model, $C_{*}$ scales with the square root of the potential energy available in the cold pools, i.e the cold pool's collapse energy, $WAPE$ (Wake Available Potential Energy), given by:
    
    \begin{equation}
       WAPE = g\int \frac{\delta \rho}{\overline{\rho}} = -g\int_{0}^{\ptop}\frac{\delta{\theta_{v}}}{\overline{\theta_{v}}}dz
        \label{eq:fp_wape}
   \end{equation}

so that:

    \begin{equation}
        C_{*} = k\sqrt{2WAPE}
        \label{eq:fp_cstar}
   \end{equation}

 where $\rho$ is the air density; $\theta_{v}$ is the virtual potential temperature.
   
For any variable $X$, $\delta X = X_{w} - X_{x}$ is the difference of its mean value in the two subdomain and $\overline{X}$ the average over the horizontal domain.

  \begin{figure}
 \centering
     \includegraphics[width=7cm,height=7cm]{figures/poche_env.png}
	  \caption{Conceptual diagram of a density current \citep{GL10I}.}
       % \label{figure_poche}
	 \label{fig:poche}
 \end{figure}
  
\remarque{Il faut faire attention. Ici on décrit le model. Dans le modèle, k ne dépend pas de la structure des poches. Il est fixé. C'est important de faire la distinction pour bien que le lecteur comprenne le modèle. On impose une proportionalité entre deux choses. On ne peut même pas dire que ce coefficient dépend de quelque chose. C'est un modèle. Dans des développement ultérieur, on peut imaginer ...}
So far, $k$ in equation \eqref{eq:fp_cstar} is fixed, to a value between 0 and 1. This coefficient should probably depend on the structure of cold pools. Based on 3D CRM (Cloud Resolving Models) simulations, Lafore (2000) (oral communication) estimated this coefficient at 0.33 in the case of a linear structure such as squall grain. This is the value used in the GL10 model.

The spread rate of cold pools $C_{*}$ is deduced from the following relationship :
     \begin{equation}
  \partial_{t} \sigma_{wk} = 2\pi rC_{*} D_{wk} = 2C_{*} \sqrt{\pi D_{wk} \sigma_{wk}}
        \label{taux_exp}
\end{equation}

where $\sigma_{wk}$ is the surface fraction covered by cold pools ($\sigma_{wk}$ = $D_{wk} \pi r^{2}$).
    %\\
    %\begin{equation}
%\partial_{t} \sigma_{wk} = 2C_{*} \sqrt{\pi D_{wk} \sigma_{wk}}
        %\label{taux_exp}
    %\end{equation}
    %\\
    Due to the complex life cycle of cold pools (including birth, death, collisions and mergers), calculating their number density requires an other parameterization. In GL10 model, so far in its LMDZ implementation, and in this study, the value of the cold pool number density is imposed. In standard configurations of LMDZ, \remarque{C'est dans LMDZ ou vraiment dans GL10 ?} this density is fixed to a different value over ocean (10.10$^{-10}$ $m^{-2}$, i.e. 10 cold pools over 100~km$\times$100~km) and over the continent (8.10$^{-12}$ $m^{-2}$, i.e. around 8 cold pools over 1000~km$\times$1000~km). In the GL10 model, cold pools initially appear with a surface fraction of 2$\%$ and evolve over time according to equation \ref{taux_exp}. The evolution of $\sigma_{wk}$ is arbitrarily limited to a maximum of 40$\%$ ($\sigma_{wk} \leq 0.4$).

\remarque{les 2\% et 40\% sont des paramètres ajustables fondamentalement.}
    
    It is assumed that below the top of cold pool ($\ptop$), the vertical velocity profile associated with the subsidence of the cold pool results solely from the spreading at the surface, without lateral entrainment ($e_{w}$) or detrainment ($d_{w}$) between the cold pool and its environment.
Above this level, the subsidence is driven both by the need to fed the subsidence in the cold pool (air mass continuity equation) and by the additional reevaporation of rainfall bellow stratiform clouds.
The shape of the vertical profile of the velocity difference $\delta\omega$ is imposed as a piecewise linear function of pressure:
 $\delta{\omega}$ increases linearly from zero at the surface up to a maximum value at $\ptop$ and then decreases linearly between $\ptop$ and a maximum height at which it vanishes : $\pupper$.
The vertical subsidence which thus increases downward between $\pupper$ and $\ptop$ is fed by lateral entrainment ($e_w>0$) without detrainment.
This lateral entrainment accounts for the horizontal component of the meso-scale circulation known to entrain air from low- or mid- tropospheric air into the cold pool.
\remarque{Pour moi il est très importants qu'on se mette d'accord sur le paragraphe ci dessus, et que Lamine comprenne l'importance des ces discussions}

At $\pupper$, the top of the cold pool model, $\delta_X$, cancels for all cold pool state variables.

In GL10 model, $\pupper$ was set to 600 hPa and there was also a nonzero velocity difference ($\delta{\omega}^{cv}$) at $\pupper$, accounting for the difference of the convective mass fluxes between (w) and (x). In the version use on the paper, this difference is now zero ($\delta{\omega}^{cv} = 0$) above this level. 

%Above $\pupper$ is considered the altitude above which the thermodynamic differences between the cold pools and their exterior become negligible, except within convective currents (saturated currents in (x), unsaturated descending currents in (w)).


 
 The evolution of the potential temperature difference ($\delta \theta$) between (w) and (x) is controlled by differential heating ($\delta Q_{1}^{cv}$, $\delta Q_{1}^{wk}$) due to deep convection and cold pools, as well as by damping due to gravity waves ($\tau_{gw}$). The humidity difference ($\delta q$) follows a similar pattern, but without the damping effect of gravity waves. Heat sources are replaced by moisture sources ($\delta Q_{2}^{cv}$ for convection and $\delta Q_{2}^{wk}$ for cold pools).
    
     %  \[
   \begin{equation}
          \begin{cases}  
            \partial_{t} \delta \theta = - \overline{\omega} \partial_{p} \delta \theta + \frac{\delta Q_{1}^{wk} + \delta Q_{1}^{cv}}{C_{p}} - \frac{K_{gw}}{\tau_{gw}} \delta \theta, \\
            \\
            \partial_{t} \delta q = - \overline{\omega} \partial_{p} \delta q + \frac{\delta Q_{2}^{wk} + \delta Q_{2}^{cv}}{L_{v}}.   
          \end{cases}
          \label{couplage}
   \end{equation}

   where 
\begin{equation}
\tau_{gw} = \frac{\sqrt{\sqrt{\sigma_{wk}} - (1  - \sqrt{\sigma_{wk}})}}{4Nz \sqrt{D_{wk}}}
\end{equation}
is estimated as the time required for a wave with speed $Nz$ to travel a distance equal to the geometric mean of the cold pool size and the interval between cold pools. $C_{p}$ is the heat capacity of dry air, $N$ is the Brunt-Väisälä frequency, and $z$ is altitude. $K_{gw}$ is an efficiency of gravity waves. Finally, $L_{v}$ is the latent heat of vaporization of water.
   
   $\delta Q_{1}^{wk}$ (respectively $\delta Q_{2}^{wk}$) depend on the entrainment ($e_{w}$) of dry air, the differential advection of $\overline{\theta}$ (respectively $\overline{q}$) and $\delta \theta$ (respectively $\delta q$):
\begin{equation}
    \begin{cases} 
    \frac{\delta Q_{1}^{wk}}{C_{p}} = \frac{e_{w}}{\sigma_{wk}} \delta \theta - \delta \omega \partial_{p} \overline{\theta} - (1 - 2\sigma_{wk}) \delta \omega \partial_{p} \delta \theta \\
    \frac{\delta Q_{2}^{wk}}{L_{v}} = \frac{e_{w}}{\sigma_{wk}} \delta q - \delta \omega \partial_{p} \overline{q} - (1 - 2\sigma_{wk}) \delta \omega \partial_{p} \delta q \\
    \end{cases}
\end{equation}


 Similarly, $\delta Q_{1}^{cv}$ (respectively $\delta Q_{2}^{cv}$) are computed from the by heating tendencies associated with unsaturated currents ($Q_{cv}^{1, unsat}$, or $Q_{cv}^{2,unsat}$ for humidity) and saturated currents ($Q_{cv}^{1,sat}$, or $Q_{cv}^{2,sat}$ for humidity):
\begin{equation}
    \begin{cases} 
    \delta Q_{1}^{cv} = \frac{Q_{cv}^{1,unsat}}{\sigma_{wk}} - \frac{Q_{cv}^{1,sat}}{1- \sigma_{wk}} \\
    \\
    \delta Q_{2}^{cv} = \frac{Q_{cv}^{2,unsat}}{\sigma_{wk}} - \frac{Q_{cv}^{2,sat}}{1- \sigma_{wk}}
    \end{cases}
\end{equation}

$Q_{cv}^{x,unsat}$ and $Q_{cv}^{x,sat}$ ($x=1,2$) are given by the deep convection scheme but we will not go into details here.

   Entrainment is determined from the vertical gradient of $\delta \omega$ and the cold pool spreading rate, according to the following relationship:
   
\remarque{Bien vérifier la formule ci dessous. Elle ne me saute pas aux yeux}

   \begin{equation}
  e_{w} = \sigma_{wk}(1 - \sigma_{wk}) \partial_{p} \delta \omega + \partial_{t} \sigma_{wk}
        %\label{flux_masse}
   \end{equation}
   
   Equation \ref{couplage}, via the variables $\delta Q_{1}^{cv}$ and $\delta Q_{2}^{cv}$, describes the impact of deep convection on cold pools which results in their cooling due to precipitating descents, as discussed above.
   
   The cold pool model is now fully described. It includes:
  \begin{itemize}
        \item 
        three prognostic variables, derived directly from the model equations: the profiles of $\delta{\theta}$ and $\delta{q}$ and $\sigma_{wk}$.
        \item 
        two diagnostic variables, evaluated from the profile of $\delta{\theta}$: $\ptop$, $C_{*}$ and $WAPE$
        \item 
        three main free parameters: the coefficient $k$, the density $D_{wk}$ and $\tau_{gw}$ \remarque{Est-ce qu'on ne veut pas mettre les surfaces min et max ? Et d'autres peut être ?}
        \\
  \end{itemize} 

  \subsection{Coupling with deep convection model}
   In LMDZ, the deep convection scheme used is based on that of \cite{emanuel1991}, with modifications made by \cite{grandpeix2004}. The main modifications concern the conditions under which deep convection is triggered and closed.
   
   Deep convection is triggered when the Available Lifting Energy ($ALE$) exceeds the convective inhibition  (CIN) threshold. This can be caused either by uplift energy from the convective boundary layer ($ALE_{bl}$), provided by the thermals model \citep{rio2008}, or by energy generated by cold pools ($ALE_{wk}$):
   \begin{equation}
           \max(ALE_{bl},ALE_{wk}) > \left| CIN \right|
           %\label{flux_masse}
   \end{equation}
  
For deep convection triggering by thermal plumes, an additional condition is required: at least one cumulus cloud in a mesh exceeds a given size, specified by $S_{trig}$. For this, a probability of non-triggering is estimated, based on the characteristics of the spectrum of type 2 thermals ($N_{2}$, $S_{2}$) derived from the original thermal plume model \citep{rio2008}$ and S_{trig}$, according to the relationship proposed by \cite{rochetin2014}.

   \begin{equation}
           P_{\Delta t} = ([1 - exp(\frac{-S_{trig}}{S_{2}})]^{N_{2}})^{\frac{\Delta t}{\tau}}
           \label{prob}
   \end{equation}

where $\Delta t$ is the model time step; $\tau$ is the decorrelation time between cloud scenes.
$S{2}$, representing the mean effective cross-sectional area of cloud base thermal currents, is determined by the following relation:

   \begin{equation}
	   S_{2} = [a(\overline{Z_{top}} - \overline{Z_{lcl}}) + b \overline{Z_{lcl}}]^{2}
           \label{S2}
   \end{equation}

where $Z_{top}$ is the mean cloud depth; $Z_{lcl}$ is the mean cloud base altitude; a and b are free parameters

$N{2}$ is the corresponding thermal population in the mesh and is calculated from the following relation:

   \begin{equation}
	   N_{2} = \frac{(1 - \epsilon) \alpha_{tot} S_{d}}{S_{2}}
           \label{S2}
   \end{equation}

where $\alpha_{tot}$ is the surface covered by thermals; $S_{d}$ is the surface of the domain



In this framework, deep convection is triggered by thermals whenever a uniform random number R, between 0 and 1, exceeds the non-trigger probability.


   \begin{equation}
           R > P_{\Delta t}
           \label{rrrr}
   \end{equation}
  

The intensity of the convection depends on the mass flux ($M_{b}$) at the cloud base, determined by $ALP$, provided by thermals ($ALP_{bl}$) and cold pools ($ALP_{wk}$).

   \begin{equation}
           M_{b} = k \frac{ALP_{bl} + ALP_{wk}}{(2 w_{b}^{2} + \left| CIN \right|)}
           %\label{flux_masse}
   \end{equation}

where $k$ and $w_{b}$ are free parameters.

The two variables $ALE_{wk}$ and $ALP_{wk}$, to take account of the effect of cold pools on convection, have been introduced into the cold pools model by \cite{GL10II}.

   To calculate $ALE_{wk}$, the model assumes that the maximum speed ($C_{max}$) on the cold pool contour will trigger convection. This is assumed to be proportional to the square root of $WAPE$, with a higher coefficient of proportionality here (arbitrarily estimated at 1), leading to the following relationship: 
   
    \begin{equation}
C_{max} = k'\sqrt{2WAPE}
        \label{eq:fp_cstar_max}
\end{equation}
    
    where $k'$ = 1
    
\remarque{Il faut peut etre mieux dire ce qu'on a en tête en distinguant C* et Cmax. Et on aura un joli modèle de ca avec la distribution sous maille du vent. J'avais pas complètement réalisé ... Cmax sort tout seul du modèle de vent.}

     The Available Lifting Energy associated with cold pools is thus expressed by the following relationship : 
     \begin{equation}
ALE_{wk} = \frac{1}{2} C_{max}^{2}
        \label{eq:fp_alewk}
\end{equation}
    
    %$C_{max}$ est la vitesse maximale de la variable aléatoire C
    %\\
    %\\
    Combining equations \eqref{eq:fp_alewk} and \eqref{eq:fp_cstar_max} gives the expression for $ALE_{wk}$ below:
    
    \begin{equation}
ALE_{wk} = k'^{2}WAPE
        \label{eq:fp_ale_wk}
\end{equation}
    
    With $k'=1$, this equation says that, in the cold pool model, the Available Lifting Energy associated with cold pools is equal to the collapse energy.
    
   $ALP_{wk}$ is calculated assuming that cold pools exert a horizontal power on the surrounding air during its spreading. This horizontal power is then converted into vertical power. During this conversion, the model assumes that a large part of the horizontal power is dissipated, and that only 25$\%$ contributes to increasing the intensity of convection. 
   
   Each cold pool generates its own lifting power, depending on its spreading speed ($C_{*}$), height ($\ptop$) and the length ($L_{g}$) of its gust front. The total power ($ALP_{wk}$) of the cold pools is the product of the power supplied by each pool times the cold pool number density ($D_{wk}$).
   
    \begin{equation}
ALP_{wk} = \epsilon \frac{1}{2}\rho C_{*}^{3} \ptop L_{g} D_{wk}
        \label{eq:fp_alpwk}
\end{equation}
    
     where $\epsilon$ = 0.25 is the lifting efficiency with 
     
     \begin{equation}
L_{g} = 2\pi r
        \label{eq:long_rafale}
   \end{equation}

     \begin{equation}
   \sigma_{wk} = D_{wk} \pi r^{2}
          \label{eq:sigmawk}
   \end{equation} 
    
   Then, the lifting power $ALP_{wk}$ reads :
    
    \begin{equation}
ALP_{wk} = \epsilon \rho C_{*}^{3} \ptop \sqrt{\sigma_{wk} D_{wk} \pi}
        \label{eq:fp_alp_wk}
\end{equation}
    \subsection{Large Eddy Simulations}

    Large Eddy Simulations (LES) are numerical tools for simulating atmospheric phenomena with a horizontal resolution of tens to hundreds of meters. They are particularly well suited to the study of the thermodynamic structure of the boundary layer, as they resolve the eddies that form there. They offer an explicit and detailed representation of turbulent and convective movements within the boundary layer and associated clouds \citep{brown2002,siebesma2003}. While they are able to reproduce atmospheric thermodynamics and structure satisfactorily, the representation of cloud characteristics remains more delicate. They enable fairly direct simulation of turbulent and convective movements. In the presence of water phase changes, however, these simulations can become highly dependent their on the microphysical schemes used. One of the major strengths of LES lies in its ability to provide three-dimensional information not available from observations, making then an indispensable complement to the latter for understanding processes. In addition, LES can be used to validate the internal variables of parameterizations, enabling their physical realism to be assessed. They have been used to evaluate boundary layer and convection parameterizations \citep{rio2010}. In recent years, they are increasingly used to document the characteristics of cold pools and guide their parameterization \citep{couvreux2012,feng2015}.

    In this study, we use the outputs of two oceanic LES and one continental LES. 
    
    Both oceanic LES were carried out in Radiative-Convective Equilibrium (RCE) mode. RCE is a concept in which equilibrium is achieved between convective heating and radiative cooling of the atmosphere. A detailed description of RCE simulation protocols is provided in \cite{daleu2015}. In the RCE simulations used here, radiative computation is replaced by a constant cooling of -1.5~K per day, while the surface temperature is imposed. The destabilization leads to convection. The associated heating rate, largely corresponding to the release of latent heat, compensates for the cooling once quasi-equilibrium has been reached. Two oceanic LES of this RCE are used here, one is performed with the SAM model \citep{khairoutdinov2003} and the other one with MesoNH \citep{lac2018}. Both simulations cover an oceanic domain of 200~km$\times$200~km with horizontal resolution of 25~m. The lateral boundary conditions are cyclic for both models. The sea surface temperature is set at 300 K. These two RCE simulations run for 44 days, with equilibrium reached on simulation day 40. Output frequency for LES SAM is set to every 3 hours, while that for LES MesoNH is set to every 24 hours.
    
    The continental LES is based on the AMMA (African Monsoon Multidisciplinary Analysis) case. This case is derived from observations made on July 10, 2006 during the AMMA field campaign \citep{redelsperger2006}, during which a relatively small, short-lived convective system formed over Niamey \citep{couvreux2012,lothon2011}. This system, with a lifetime of around 6 hours, was observed by various instruments (radar and atmospheric soundings), supplemented by satellite data. This case study represents a typical example of deep convection in the Sahel regions. LES for this continental case is carried out with the MesoNH model over a 100~km $\times$ 100~km domain, with a horizontal resolution of 200~m. Lateral boundary conditions are cyclic and surface fluxes are imposed. Outputs are generated at a frequency of 10 minutes.
   
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    \section{Assessment of the cold pool model internal equations from LES}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    \subsection{Sampling of cold pools} \label{LES_LES}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
In order to use LES for the assessment of the cold pool parameterization, the first challenge is to separate cold pools from their environment. Indeed, there is no a priori established framework for objectively identifying cold pools in observations and numerical models \citep{rochetin2021}, and choices may depend in part on the physical picture one has of cold pools, and also, for the purpose at hand, on the picture underlying the parameterization. The first method for identifying cold pools proposed by \cite{young1995} was based on surface precipitation rates. In more recent studies, such as those by \cite{provod2016,zuidema2017,vogel2021,rochetin2021,touze2022}, the detection of cold pools is closer to a density current oriented detection, in which variations in temperature, pressure and wind are taken into account.

In the present study, the aim is not to isolate individual “cold pools objects”, but only to know whether a grid box is inside or outside a cold pool.
Also the boundary conditions are idealized targeting the statistical homogeneity assumption that is at the basis of the Reynolds decomposition between dynamical core and physics parameterizations. In this idealized case with uniform surface temperature, cold pools can be identified fairly immediately using a threshold on the anomaly of temperature at 10~m above surface ($T_{10m}$, first model layer).

    The map of divergence of wind at 10~m, smoothed on a 3.25~km$\times$3.25~km box, enables us to visually identify centers and gust fronts of cold pools, represented respectively by the maximum and minimum of divergence values (\fig{divergence_vent_rce} and \fig{divergence_vent_amma}).
Maxima of divergence of surface wind indicate the center of cold pools where cold air masses collapse. Precipitation is generally co-located with these divergence maxima. The fairly strong wind convergence observed around cold pools centers corresponds to the strong lift of air masses created upstream of the gust front at the cold pool's periphery.    

    %\begin{figure}
 %\centering
 %\includegraphics[width=11cm,height=11cm]{figure_chap3/carte_divergence_l1_paper_RCE.png}
 %\includegraphics[width=11cm,height=11cm]{figure_chap3/carte_divergence_l4_paper_RCE.png}
    % \caption{Carte de la divergence du vent à 10 m (en $s^{-1}$) multipliée par 1000 et lissée horizontalement sur 3.25 km$\times$3.25 km représentée sur deux instants (a et b) de la LES SAM réalisée sur le cas RCE océanique et superposés avec les contours d'anomalies de températures à 10 m à -0.4 K (vert), -0.2 K (rouge) et 0 K (noir)}
    %    \label{divergence_vent_rce}
 %\end{figure} 
   \begin{figure}
    \centering
    % Premiere sous-figure
    \begin{subfigure}[t]{\textwidth}
        \centering
        \includegraphics[width=11.5cm,height=11.5cm]{figures/carte_divergence_l1_paper_RCE.png}
       % \caption{Divergence du vent à 10 m - Instant (a).}
        \label{fig:suba}
    \end{subfigure}
    
    \vspace{-1.5cm} % Ajuste l'espace entre les sous-figures (facultatif)
    %\vspace{0pt} % Ajuste l'espace entre les sous-figures (facultatif)

    % Deuxieme sous-figure
    \begin{subfigure}[t]{\textwidth}
        \centering
        \includegraphics[width=11.5cm,height=11.5cm]{figures/carte_divergence_l4_paper_RCE.png}
        %\caption{Divergence du vent à 10 m - Instant (b).}
        \label{fig:subb}
    \end{subfigure}

    \vspace{-0.8cm} % Ajuste l'espace entre les sous-figures (facultatif)
    
    \caption{Map of divergence of wind at 10~m (in $s^{-1}$) multiplied by 1000 and smoothed horizontally over 3.25~km $\times$ 3.25~km represented on two instants (a and b) of the LES SAM carried out on the oceanic RCE case and superimposed with the contours of temperature anomalies at 10~m at -0.4~K (green), -0.2~K (red) and 0~K (black).}
    \label{divergence_vent_rce}
 \end{figure}
 \begin{figure}
    \centering
    % Premiere sous-figure
    \begin{subfigure}[t]{\textwidth}
        \centering
        \includegraphics[width=11.5cm,height=11.5cm]{figures/carte_divergence_l2_paper_AMMA.png}
       % \caption{Divergence du vent à 10 m - Instant (a).}
        \label{fig:sub_a}
    \end{subfigure}
    
    \vspace{-1.5cm} % Ajuste l'espace entre les sous-figures (facultatif)
    %\vspace{0pt} % Ajuste l'espace entre les sous-figures (facultatif)

    % Deuxieme sous-figure
    \begin{subfigure}[t]{\textwidth}
        \centering
        \includegraphics[width=11.5cm,height=11.5cm]{figures/carte_divergence_l7_paper_AMMA.png}
        %\caption{Divergence du vent à 10 m - Instant (b).}
        \label{fig:sub_b}
    \end{subfigure}

    \vspace{-0.8cm} % Ajuste l'espace entre les sous-figures (facultatif)
    
    \caption{ Map of divergence of wind at 10~m (in $s^{-1}$) multiplied by 1000 and smoothed horizontally over 3.25~km $\times$ 3.25~km represented on the 17:10 (a) and 18:00 (b) instants of the LES MESONH carried out on the AMMA case and superimposed with the contours of anomalies of temperature at 10~m at -1~K (green), -0.5~K (red) and 0~K (black).}
    \label{divergence_vent_amma}
   \end{figure}


   Both the LES in RCE and the LES in the AMMA case show cold pools groupings (or very close cold pools centers) forming a common gust front. This can be explained by the fact that, during propagation, cold pools can merge to create a single, larger cold pools. We can also observe that wind convergence is generally more intense between the centers of grouped cold pools, indicating that updrafts of air masses associated with gust fronts is more pronounced when these cold pools meet. This is in line with some studies that indicate that convection initiation on gust fronts is more efficient when two or more cold pools collide \citep{meyer2020,torri2019,haerter2018,feng2015}.

We superimpose on this map $T_{10m}$ anomaly contours with different values to determine an optimal threshold for this anomaly.
In the RCE case, the $T_{10m}$ anomaly at 0 K sometimes includes regions without cold pools centers, where divergence of surface wind is low (\fig{divergence_vent_rce}a and \fig{divergence_vent_rce}b) while anomaly contours -0.2~K and -0.4~K surround the centers of cold pools quite well. In the AMMA case, figure \ref{divergence_vent_amma}a clearly shows that the 0~K threshold is too high to identify cold pools. Figure \ref{divergence_vent_amma}b, on the other hand, shows that the -1~K threshold follows gust fronts of cold pools better than the -0.5~K threshold. On the basis of these analyses, we retain the $T_{10m}$ anomaly thresholds at -0.2~K and -1~K to identify cold pools in the RCE and AMMA cases respectively.

After selecting values for the $T_{10m}$ anomaly to separate the inside and outside of cold pools in the RCE and AMMA cases, we carry out sampling to calculate certain variables of cold pool model in the LES. We first determine the vertical profiles of temperature ($\delta{T}$), humidity ($\delta{q}$) and vertical velocity ($\delta{w}$) differences between (w) and (x). To do this, we apply the mask to the entire column to determine the vertical profiles. This vision of vertical cylinders is obviously open to question. But it does seem to apply, at least to the idealized cases studied here.

    \subsubsection*{Calculation of the spreading speed, $C_{*}$} 

It is assumed in the parameterization that cold pools are identical disks of the same radius (r). This assumption makes it easy to determine $C_{*}$ by of the divergence theorem. 
    
 \begin{equation}
 \int\int div(\overrightarrow{V_{10}})dS_{wk} =  C_{*}L_{g}  
        \label{eq:theoreme_divergence}
 \end{equation}
 
 \begin{equation}
 C_{*} = \frac{\overline{div(\overrightarrow{V_{10}})}{S_{wk}}}{L_{g}} 
        \label{eq:e_cstar}
 \end{equation} 
 
 where $S_{wk}$ is the surface of cold pools
 
 \begin{equation}
 S_{wk} = \pi r^{2} 
        \label{surf_wk}
 \end{equation} 
 
    Equations \ref{eq:long_rafale}, \ref{eq:sigmawk} and \ref{surf_wk} allow us to express $C_{*}$ as a function of the mean divergence of wind at 10 m, the surface fraction ($\sigma_{wk}$) and the density ($D_{wk}$) of cold pools by the relation : 
 
    \begin{equation}
 C_{*} = \frac{1}{2} \overline{div(\overrightarrow{V_{10m}})} \sqrt{\frac{\sigma_{wk}}{D_{wk}\pi}}
        \label{eq:e_cstar_f}
 \end{equation}
   
  To apply this calculation of $C_{*}$ in the LES, we take the horizontal average of the surface wind divergence inside cold pools. The surface fraction ($\sigma_{wk}$) of cold pools calculated in the LES is 0.12 for the AMMA case and 0.25 for the RCE case. To determine $D_{wk}$, we manually counted the centers of cold pools visible on the surface wind divergence maps (Fig. \ref{divergence_vent_rce} and \ref{divergence_vent_amma}), as we did not use automated detection methods in this study that could generate their number automatically. We find an approximate density, $D_{wk}$, of 5 cold pools per 100~km $\times$ 100~km for both the RCE and AMMA cases.

   \subsubsection*{Calculation of collapse energy}

 We finally calculate the collapse energy ($WAPE$) of cold pools in the LES using equation (\ref{eq:fp_wape}) proposed by \cite{GL10II}. The task consists of determining $\overline{\theta_{v}}$, as well as the profiles of $\delta{\theta_{v}}$ and $\ptop$ in the LES. To do this, we first computed $\delta{T}$ in the LES, then derived $\overline{\theta_{v}}$ and the profile of $\delta{\theta_{v}}$. Regarding the determination of $\ptop$, as suggested by \cite{GL10II}, we take this height at the altitude where the $\delta{T}$ profile cancels out. This altitude is around 950~hPa (approximately 600~m) in the oceanic RCE case and around 800~hPa (approximately 2~km) in the AMMA case (\fig{profil_deltat}).
   \begin{figure}
 \centering
  \includegraphics[width=9cm,height=11cm]{figures/profil_deltat.png}
 \caption{Vertical profiles of the temperature difference between the inside and the outside of cold pools calculated at an instant of the LES (SAM and MESONH) of the RCE case and an instant of the LES MESONH of the AMMA case.}
        \label{profil_deltat}
 \end{figure}
 
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   \subsection{Computing ALP and ALE form gust front vertical velocities}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

Here we compute the variables of Available Lifting Energy ($ALE_{wk}$) and Power ($ALP_{wk}$) associated with cold pools in the LES. To do this, we proceed in several steps: 

   \begin{enumerate}
       \item 
       We first determine an average cloud-base height at which we extract vertical velocities $w_{b}(x,y)$. This height corresponds to the altitude at which the average profile of condensed water reaches its first non-zero value. It is observed at around 950 hPa on the two oceanic LES (SAM and MesoNH) and at around 750 hPa on the LES for the AMMA case (MESONH) (\fig{profil_nuage}).
       \item 
       We then separate the updrafts on gust fronts from those associated with thermal plumes. Since the updrafts on gust fronts are both stronger and more coherent horizontally than the thermals observed in the environment of cold pools, we defined a mask based on a threshold an horizontal smoothed $w_{b}$ over 1.25~km$\times$1.25~km (RCE) et 2~km$\times$2~km (AMMA). The smoothed $w_{b}$ values are denoted as $\tilde{w}_{b}(x,y)$ in the rest of the text. After several analyses, we selected $\tilde{w}_{b}(x,y)$ thresholds of 0.6 m/s for the RCE case and 2 m/s for the AMMA case to identify gust fronts.
   \end{enumerate} 

      Figure \ref{gust_thermiques} presents maps of $T_{10m}$ anomaly, smoothed horizontally on a 2.5~km$\times$2.5~km grid, for the RCE and AMMA cases. On these maps, we have overlaid the contours of the $T_{10m}$ anomalies used to identify cold pools (-0.2~K for RCE and -1~K for AMMA), as well as the updrafts on gust fronts (in red) and thermals (in green). Visually, the gust fronts computed with $\tilde{w}_{b}(x,y)$ thresholds of 0.6~m/s (RCE) and 2~m/s (AMMA) align well with the contours of cold pools identified using these $T_{10m}$ anomaly thresholds. It also appears that most thermals are located in the environment of cold pools for both the RCE and AMMA cases (\fig{gust_thermiques}). This retrospectively validates a choice made in version 6A of the model, where the effect of thermals was only computed outside cold pools.
   Finally, to determine $ALE_{wk}$, we take the maximum kinetic energy in the domain, considering only $w_{b}(x,y)$ in the gust fronts mask ($w_{bgust}(x,y)$), as it is the maximum vertical velocity on the gust front that triggers convection. As for $ALP_{wk}$, which represents the average updrafts power provided by all cold pools in the domain, it is calculated from the horizontal average of the cube of $w_{bgust}$, weighted by the surface fraction ($\sigma_{gust}$) covered by gust fronts. The mask applied to gust fronts was used to calculate $\sigma_{gust}$, which is 0.017 for the RCE case and 0.067 for the AMMA case, for the times shown in figure. \ref{gust_thermiques}.
   \begin{equation}
 ALE_{wk} = max (\frac{1}{2}w_{bgust}^{2})
        \label{eq:e_ale_wk}
 \end{equation}
    
    \begin{equation}
 ALP_{wk} = \sigma_{gust}\frac{1}{2}\rho \overline{{w_{bgust}^{3}}}
        \label{eq:e_alp_wk}
 \end{equation}. 
   
    \begin{figure}
 \centering
 \includegraphics[width=9cm,height=11cm]{figures/profil_nuage_bis.png}
 \caption{Vertical profile of condensed water averaged horizontally on the LES in oceanic RCE carried out with the SAM and MésoNH models and the continental LES of the AMMA case carried out with MésoNH.}
        \label{profil_nuage}
   \end{figure}  
  \begin{figure}
    \centering
    % Premiere sous-figure
    \begin{subfigure}[t]{\textwidth}
        \centering
        \includegraphics[width=10.5cm,height=10.5cm]{figures/fig_poche_thermique_gust_RCE.png}
       % \caption{Divergence du vent à 10 m - Instant (a).}
        \label{fig:sub_a}
    \end{subfigure}
    
    \vspace{-1.5cm} % Ajuste l'espace entre les sous-figures (facultatif)
    %\vspace{0pt} % Ajuste l'espace entre les sous-figures (facultatif)

    % Deuxieme sous-figure
    \begin{subfigure}[t]{\textwidth}
        \centering
        \includegraphics[width=10.5cm,height=10.5cm]{figures/fig_poche_thermique_gust_AMMA.png}
        %\caption{Divergence du vent à 10 m - Instant (b).}
        \label{fig:sub_b}
    \end{subfigure}

    \vspace{-0.8cm} % Ajuste l'espace entre les sous-figures (facultatif)
    
    \caption{Maps of anomaly of temperature at 10~m, smoothed horizontally over 2.5~km$\times$2.5~km, represented on an instant of the LES SAM of the RCE case (a) and on the instant 6:00 PM of the LES of the AMMA case with black contours indicating thresholds of temperature at 10~m anomaly of -0.2~K (RCE) and -1~K (AMMA). The red color indicates the updrafts on the gust fronts given by the vertical velocities at cloud base ($w_{b}$) in the gust fronts mask, which is determined by the $w_{b}$ smoothed horizontally over 1.25~km$\times$1.25~km and exceeding 0.6~m/s (RCE) and over 2~km$\times$2~km with a value greater than 2~m/s (AMMA). The green dots represent thermals, defined by $w_{b}$ outside the gust front mask.}
    \label{gust_thermiques}
\end{figure}
   %\begin{figure}
 %    \centering      
    %    \includegraphics[width=9.5cm,height=9.5cm]{figure_chap3/fig_poche_thermique_gust_RCE.png}
    %    \includegraphics[width=9.5cm,height=9.5cm]{figure_chap3/fig_poche_thermique_gust_AMMA.png} 
    %    \caption{Cartes d'anomalie de température à 10 m, lissées horizontalement sur 2.5 km $\times$ 2.5 km, représentées sur un instant de la LES SAM du cas RCE (a) et sur l'instant 18h de la LES du cas AMMA avec contours noirs indiquant les seuils d'anomalies de température à 10 m de -0.2 K (RCE) et -1 K (AMMA). La couleur rouge indique les ascendances sur les fronts de rafales données par les vitesses verticales à la base des nuages ($w_{b}$) dans le masque des fronts de rafales qui est déterminé par les $w_{b}$ lissées horizontalement sur 1.25 km $\times$ 1.25 km et dépassant 0.6 m/s (RCE) et sur 2 km $\times$ 2 km avec une valeur supérieure à 2 m/s (AMMA). La points verts représentent les thermiques, définis par les $w_{b}$ en dehors du masque des fronts de rafales.}
    %    \label{gust_thermiques}
  %\end{figure}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    \subsection{Validation of Phenomenological Laws} \label{formule_param}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
 Physical parameterizations are defined by sets of mathematical equations intended to represent the subgrid process within a column of the model. The formulation of these equations is based on both a phenomenological understanding of the processes concerned and fundamental principles of physics. These parameterizations can be assessed in bulk, or piecewise, by isolating certain equations or relations between internal variables, or between internal variables and state variables of the GCM. LES offer the possibility of performing a priori validation and adjustment of these laws.

    In the cold pool model, variables $ALE_{wk}$, $ALP_{wk}$ and $C_{*}$ are determined from the collapse energy ($WAPE$) (see equations \eqref{eq:fp_cstar}, \eqref{eq:fp_ale_wk} and \eqref{eq:fp_alp_wk}), based on assumptions derived from physical laws. 

Here, we compare the values of $ALE_{wk}$, $ALP_{wk}$ and $C_{*}$ computed from the $WAPE$ with those obtained directly from the vertical speed at cloud base near the gust fronts ($w_{bgust}$) for $ALE_{wk}$ and $ALP_{wk}$, and from the mean divergence of wind at 10~m in cold pools for $C_{*}$.

    Table \ref{LES_LMDZ_CTRL} shows such comparison for the three LES available. The values of $ALE_{wk}$ calculated from $w_{bgust}$ and $WAPE$ from the $\delta{\theta_{v}}$ profiles are very close. In the LES of the RCE case made with SAM, $ALE_{wk}$ calculated from $w_{bgust}$ is slightly higher than the $WAPE$ from the $\delta{\theta_{v}}$ profile (table \ref{LES_LMDZ_CTRL}). However, even in this case, $ALE_{wk}$ determined from $w_{bgust}$ remains comparable to the $WAPE$ derived from the $\delta{\theta_{v}}$ profile. These results for the three LES are compatible with the hypothesis of equality between $ALE_{wk}$ and $WAPE$, as estimated by the parameterization.

    \begin{table}[ht]
    %\begin{center}
     \caption{Comparison of the variables of $WAPE$, $ALE_{wk}$, $C_{*}$ and $ALP_{wk}$ calculated in the samplings (E) and those calculated with the formulas of the parameterization (FP) for the coefficient $k$ = 0.33 and $k$ = 0.66 in the oceanic LES in RCE carried out with SAM and MESONH and in the continental LES of the AMMA case carried out with MESONH}
     \label{LES_LMDZ_CTRL}
     \centering
    % \begin{adjustbox}{max width=\textwidth} 
     \small
     \begin{tabularx}{\textwidth}{|X|X|X|X|X|X|X|X|X|}
     %\begin{tabularx}{\textwidth}{|l|l|l|l|l|l|l|l|l|} (\cite{de2021subseasonal})
     \hline
	     \textbf{} & \textbf{$WAPE$ }(J/Kg) & \textbf{$ALE_{wk}$ }(J/kg) (E) & \textbf{$C_{*}$ }(m/s) (FP) k=0.33 & \textbf{$C_{*}$ }(m/s) (E) & \textbf{$C_{*}$ }(m/s) (FP) k=0.66 & \textbf{$ALP_{wk}$ }(J/kg) (FP) k=0.33 & \textbf{$ALP_{wk}$ }($W/m^{2}$) (E) & \textbf{$ALP_{wk}$ }($W/m^{2}$) (FP) k=0.66 \\ \hline
    RCE SAM & 7.962 & 10.460 & 1.315 & 2.228 & 2.630 & 0.008 & 0.054 & 0.071 \\ \hline
    RCE MESO & 7.912 & 6.965 & 1.313 & 2.264 & 2.625 & 0.008 & 0.020 & 0.071 \\ \hline
    AMMA MESO & 34.250 & 33.480 & 2.727 & 4.939 & 5.454 & 0.104 & 0.982 & 0.831 \\ \hline
    \end{tabularx}
   % \end{adjustbox}
   %\end{center}
   \end{table}

	    Table \ref{LES_LMDZ_CTRL} shows that, $C_{*}$ values calculated from the $WAPE$ are systematically lower than those coming from the mean divergence of wind at 10 m in cold pools. This difference could be due to an underestimation of the coefficient $k$, imposed here at 0.33. By setting $k$ to 0.66, the calculation of $C_{*}$ based on the $WAPE$ becomes comparable to those obtained from the mean divergence of wind at 10~m in cold pools, notably for the RCE and AMMA cases (table \ref{LES_LMDZ_CTRL}). As discussed above, the value of 0.33 was retained following an oral communication by Lafore (2000). But other studies propose different values: \cite{lafore1989} estimate $k$ at 0.68 based on CRM simulations of 2D squall grain, \cite{bryan2005} estimate it at 0.5 from observations of cold pools during the BAMEX experiment in the American Great Plains. These results are thus compatible with the hypothesis of the model which postulates that the kinetic energy of cold pools results from the transformation of $WAPE$ into kinetic energy with a coefficient $k$ compatible with the published estimates.

  Table~\ref{LES_LMDZ_CTRL} also shows that, for the three LES cases, the values of $ALP_{wk}$ calculated with $C*$ from $WAPE$ are at least three times lower than those obtained from $w_{bgust}$. Two coefficients are involved in the calculation of $ALP_{wk}$ with the parameterization formula: the coefficient $k$ and the lifting efficiency $\epsilon$, imposed respectively at 0.33 and 0.25. Using $k$=0.66 however in the calculation of $C_{*}$, and keeping $\epsilon$ at is nominal value of 0.25 allows to reconcile the various estimates.
This is compatible with the hypothesis of the parameterization according to which 25$\%$ of the horizontal power provided by the cold pools during its propagation would be used to reinforce the intensity of the convection while a large part dissipates.
   
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   \section{Comparison between LES and standard LMDZ}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   \subsection{Vertical profiles of $\delta{T}$, $\delta{q}$ and $\delta{w}$}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

  In this section, we evaluate the profiles of $\delta{T}$, $\delta{q}$ and $\delta{w}$ computed by LMDZ versus LES ones. In SCM mode, the comparison is more demanding than those previously discussed, since all parameterizations interact with each other to arrive at the simulated values several hours (AMMA) or days (RCE) after initialization. 

  For the RCE case, we represent the profiles once a quasi-steady state has been reached.
Regarding the AMMA case, intermediate analyses show that cold pools appear in the afternoon around 5:00 PM with relatively low temperatures and develop during the day. Due to the variations in cold pools characteristics at different times on the continent, we average the results over the 7 available times between 5:00 PM and 6:00 PM to simplify our analyses. To compare with the model, we perform a single-column LMDZ control simulation (LMDZ CTRL) for the RCE and AMMA cases. These LMDZ simulations are performed with exactly the same initial and boundary conditions as the corresponding LES. For the RCE case, we perform a 44-day LMDZ CTRL simulation to reach a quasi-equilibrium. For the AMMA case, the LMDZ CTRL simulations are performed over the day of July 10, 2006, from 6:00 AM to midnight.


	    For the AMMA case, the cloud size threshold ($S_{trig}$) controlling the triggering of deep convection is adjusted so that convection triggers at the same time as in the LES in order to allow a precise comparison. Indeed, convection triggers before 2:00 PM in the AMMA case with the standard LMDZ configuration, while in the LES, it appears around 5:00 PM. To obtain a triggering simultaneous with that of the LES, we performed tests by modifying the value of $S_{trig}$. These tests made it possible to obtain the triggering of convection in the LMDZ simulation of the AMMA case at 4:50 PM by setting $S_{trig}$ at 24 $km{^2}$. In order to facilitate comparisons between LMDZ and LES, we also impose in the LMDZ simulations the density of cold pools estimated in the LES. We thus set a density of 5 cold pools per 100~km$\times$100~km, both for the RCE and AMMA cases. To represent the profiles of $\delta{T}$, $\delta{q}$ and $\delta{w}$ in LMDZ CTRL for the RCE case, we perform a time average between the 41st and 43rd day of simulation, in order to compare with the LES at the same times. For the AMMA case, we perform an average between 5:00 PM and 6:00 PM, as in the LES.

     \begin{figure}
     %\centering
     \includegraphics[width=16cm,height=8cm]{figures/fig_deltaT_Q_W_LES_SAM24_MESONH2_SCM_STAND.png}
     \includegraphics[width=16cm,height=8cm]{figures/fig_deltaT_Q_W_LES_LMDZ_AMMA_ptop0_SCM_STD.png}
        \caption{Vertical profiles of $\delta{T}$, $\delta{q}$ and $\delta{w}$ calculated in the LES and simulated by LMDZ control (LMDZ CTRL) on the RCE case (a, b, c) and on the AMMA case (d, e, f).}
        \label{profil_LES_CTRL}
 \end{figure}


	    The analysis of the $\delta{T}$ profiles in the LES confirms that cold pools are colder at the surface with temperatures increasing towards the top for the three LES. The cold pools are about three times deeper in AMMA (\fig{profil_LES_CTRL}a) than for the RCE case (\fig{profil_LES_CTRL}d). In the LES, we observe that the temperatures of cold pools for the AMMA case (around –2~K) are relatively close to those of the RCE case (around –1.2~K). However, observations indicate that cold pools are significantly more intense over continents than over oceans. The results obtained with the LES for the AMMA case could be explained by the fact that the analysis is carried out from the first moments following the formation of the cold pools. However, for this same case, observations reveal a temperature drop of around –5~K during the passage of cold pools \citep{lothon2011}. It should also be noted that the AMMA case corresponds to a particularly weak and atypical episode of continental convection. The $\delta{q}$ profiles indicate that at the surface, cold pools are wetter than their surroundings in the RCE case (\fig{profil_LES_CTRL}b), while in the AMMA case, they are slightly drier (\fig{profil_LES_CTRL}e). In both cases, the humidity within the cold pools decreases with altitude until they reach their summit, where they are dried by the subsidence of dry air masses into cold pools (\fig{profil_LES_CTRL}c and \ref{profil_LES_CTRL}f). On the RCE case, this subsidence vanishes below 800 hPa (\fig{profil_LES_CTRL}c), while for the AMMA case, it vanishes at a higher level, around 600 hPa (\fig{profil_LES_CTRL}f).	    


    The $\delta{T}$ profiles simulated with LMDZ CTRL are qualitatively consistent with LES, with a cold pool top (where $\delta T$ cancelled) at about the right altitude. Cold pools simulated with LMDZ are however warmer than LES for the RCE case (\fig{profil_LES_CTRL}a), and slightly colder at the surface than the LES for the AMMA case (\fig{profil_LES_CTRL}d). Consistently with LES, cold pools are also wetter at the surface and drier close to their top top (\fig{profil_LES_CTRL}b and \fig{profil_LES_CTRL}e). However the variations of $\delta{q}$ are much larger in LMDZ than LES. In particular, the cold pools are much too dry at their top in LMDZ. In both cases, cold pools are associated with subsidence. The height at which the subsidence of air masses in cold pools begins, fixed at 600~hPa in LMDZ CTRL, is too high compared to LES for the RCE case \ref{profil_LES_CTRL}e). 

The comparisons also reveal that the model simulates wetter cold pools at the surface than those in the LES in both cases, with a more pronounced difference for the RCE case.
    
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    \subsection{WAPE, ALE and ALP}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

   \begin{table}[ht]
   %\begin{center}
    \caption{Comparison of the $WAPE$, $ALE_{wk}$, $C_{*}$ and $ALP_{wk}$ computed from sampling of the LES and by LMDZ control (LMDZ CTRL) for the RCE AMMA cases. 
\remarque{La phrase ci dessous a typiquement sa place dans la caption qui doit décrire précisemenent ce qu'on montre}
The results are averaged over the days following the achievement of equilibrium (days 41, 42 and 43) for RCE and over the available instants between 5:00 PM and 6:00 PM for AMMA.  \label{LES_LMDZ_CTRL2}}
    \centering
    %\begin{adjustbox}{max width=\textwidth} 
    \small
    \begin{tabularx}{\textwidth}{|X|X|X|X|X|}
    \hline
	    \textbf{} & \textbf{$WAPE$ }(J/Kg) & \textbf{$ALE_{wk}$ }(J/kg) & \textbf{$C_{*}$} (m/s) &  \textbf{$ALP_{wk}$ }($W/m^{2}$) \\ \hline
    \multicolumn{5}{|c|}{RCE} \\ \hline
   LES SAM & 7.962 & 10.460 & 2.228 & 0.054 \\ \hline
   LES MESONH & 7.912 & 6.965 & 2.264 & 0.020 \\ \hline
   LMDZ CTRL & 2.957  & 2.957  & 0.802 & 0.001  \\ \hline
   \multicolumn{5}{|c|}{AMMA} \\ \hline
   LES MESONH & 34.250 & 33.480 & 4.939 & 0.982 \\ \hline
   LMDZ CTRL & 30.430 & 30.430 & 2.574 & 0.042  \\ \hline
   \end{tabularx}
  \end{table}

    For the RCE case, the $WAPE$ is significanlty smaller in LMDZ CTRL than in the LES, with a difference of at least a factor of 2
(Table \ref{LES_LMDZ_CTRL2}).
These low values of $WAPE$ in LMDZ CTRL also translate into low $ALE_{wk}$ values compared to LES (table \ref{LES_LMDZ_CTRL2}).

Indeed,  $ALE_{wk}$ in the RCE case it at least twice as low in LMDZ CTRL as in the LES. On the other hand, for the AMMA case, the $WAPE$ simulated by the model are globally in agreement with the values calculated in the LES (table \ref{LES_LMDZ_CTRL2}), which allows the model to obtain $ALE_{wk}$ comparable to those of the LES for this case (table \ref{LES_LMDZ_CTRL2}).

The value of $C_{*}$ simulated by LMDZ CTRL is also at least three times smaller than in the LES for all the cases (table \ref{LES_LMDZ_CTRL2}). This leads to a $ALP_{wk}$ value ten times smaller than in the LES.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    \section{Improvements of cold pool model}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

	   Here, we start by correcting the observed discrepancies between the LES and the model concerning the value of the coefficient $k$ and the altitude $\pupper$, and by assessing the impact of these changes on the temperature and humidity difference between the cold pools and their environment, before exploring other avenues for improvement.  
    
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    \subsection{Coefficient $k$}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    We present here the impact of increasing the coefficient $k$ from 0.33 to 0.66 (LMDZ V1 simulation) on the profiles of $\delta{T}$, $\delta{q}$, $\delta{w}$ as well as on the variables $C_{*}$, $WAPE$, $ALP_{wk}$ and $ALE_{wk}$.
	   In the RCE case, this modification significantly improves the profile of $\delta{w}$ below $\ptop$ (\fig{profil_LES_LMDZ_MODIF}c). This improvement is directly linked to an increase in $C_{*}$ (Table \ref{LES_LMDZ_MODIF}), since the profile of $\delta{w}$ below $\ptop$ depends on the spreading of cold pools. The increase in $C_{*}$ could be associated with a stronger air mass subsidence in the cold pool, which would contribute to a slight drying near the surface (\fig{profil_LES_LMDZ_MODIF}b). Similar effects, although less marked, are observed for the profiles of $\delta{w}$ and $\delta{q}$ in the AMMA case (\fig{profil_LES_LMDZ_MODIF}e,f). Nevertheless, the increase of $k$ also leads to an improvement of $C_{*}$ for this case (Table \ref{LES_LMDZ_MODIF}). The improvement of $C_{*}$ in both cases, RCE and AMMA, is also at the origin of a better representation of $ALP_{wk}$ (increase by a factor of 6 for RCE and by a factor of 10 for AMMA), even if this variable remains underestimated (Table \ref{LES_LMDZ_MODIF}). However, despite this modification of $k$, cold pools remain too dry at their top and wetter at the surface, in both cases (\fig{profil_LES_LMDZ_MODIF}b,e). In the RCE case, they also remain less cold in LMDZ V1 than in the LES (\fig{profil_LES_LMDZ_MODIF}a). For the AMMA case, a slight flattening of the $\delta{T}$ profile is observed, which nevertheless remains globally consistent with the LES (\fig{profil_LES_LMDZ_MODIF}d). The impact on the $\delta{T}$ profiles in the AMMA and RCE cases is at the origin of the decrease in the values of $WAPE$ and $ALE_{wk}$ for these two cases (Table \ref{LES_LMDZ_MODIF}).
 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    \subsection{Altitude $h_{m}$ } \label{modif2}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    \begin{figure}
     %\centering
     \includegraphics[width=16cm,height=8cm]{figures/fig_deltaT_Q_W_LES_SAM24_MESONH2_LMDZ_MODIF1234_NEW.png}
    \includegraphics[width=16cm,height=8cm]{figures/fig_deltaT_Q_W_LES_LMDZ_AMMA_LMDZ_MODIF1234_NEW.png}
        \caption{Vertical profiles of $\delta{T}$, $\delta{q}$ and $\delta{w}$ calculated in the LES and simulated in the control LMDZ (CTRL), LMDZ with the adjustment of the coefficient k to 0.66 (V1), LMDZ with the drop in altitude ($\pupper$) at which the subsidence of the air masses in cold pools is zero (V2) and LMDZ with the activation of thermals in the entire domain (V3) on the RCE case (a, b, c) and on the AMMA case (d, e, f).}
        \label{profil_LES_LMDZ_MODIF}
 \end{figure}

   In the previous sections, we found that the altitude at which the subsidence of dry air above cold pools initiate is observed in LES below 800 hPa for the RCE case and below 600 hPa for the AMMA case, while in LMDZ, this altitude $\pupper$ was arbitrarily set to 600 hPa in the original version of the parameterization. In version V2, in addition to the change of the value of $k$ from 0.33 to 0.66, we compute $\pupper$ as $\alpha \ptop$ with $\alpha = 3$ ($\alpha$ is considered as a new free parameter in the following section). A slight adjustment of $\ptop$ was also made thanks to the new numerical scheme proposed for its calculation, although the details are not discussed here. This adjustment, however, has no impact on the vertical profiles of $\delta{T}$, $\delta{q}$.

Comparisons between LMDZ V2 simulations and LES show a better representation of the $\delta{q}$ profiles at the top of cold pools in both the RCE and AMMA cases (\fig{profil_LES_LMDZ_MODIF}b and \fig{profil_LES_LMDZ_MODIF}e). These results show that the dry bias at the top of the cold pool in the original version was due to advection of dry air from too high an altitude. The fact that the model is able to reproduce consistently the vertical profile of vertical wind and the humidity at the top of cold pools when adjusting $\pupper$ validates the physics implemented in this cold pools scheme. This reveals that a limitation of this scheme lays in the choice of the value of a parameter (the height above 600~hPa) rather than in the formulation themselves. This modification also reduces slightly the humidity at the surface of cold pools in the RCE case. We however observe that cold pools always remain more humid at the surface in LMDZ V2 than in the LES. Concerning the $\delta{T}$ profiles, \fig{profil_LES_LMDZ_MODIF}a and \fig{profil_LES_LMDZ_MODIF}d indicate that this modification has a very limited impact on the $\delta{T}$ profiles in the AMMA and RCE cases. 

Table \ref{LES_LMDZ_MODIF} shows that the decrease in $\pupper$ weakly affects the variables $WAPE$, $C_{*}$, $ALE_{wk}$ and $ALP_{wk}$ for these two cases.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    \subsection{Activation of thermals throughout the domain}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
 To understand the origin of the wet bias in surface  of cold pools in LMDZ, we test to activate thermals throughout the domain. In the standard LMDZ configuration, thermals only interact with temperature and humidity profiles outside cold pools. This choice was originally made to account for the fact that the atmosphere is more stable inside cold pools, which would inhibit convection in these regions. In the version V3, we interact thermals with the grid-averaged temperature and humidity profiles, starting from the V2 configuration. For the RCE case, the LMDZ V3 simulations show a clear decrease in the surface humidity of cold pools, corresponding better to the results obtained with the LES (\fig{profil_LES_LMDZ_MODIF}b). This result is expected because the vertical transport by thermals systematically dries the surface \citep{diallo2017}. We also observe a slight modification of the profile at the top of cold pools. In the AMMA case, this modification has almost no effect (Fig. \ref{profil_LES_LMDZ_MODIF}e), probably because the simulation duration is too short. In the 30-day 3D test, cold pools also dry out significantly on the continent with this modification.These results suggest a key role for thermals in regulating surface humidity, via the mixing of moist air with dry air above. To represent this effect in the model, one could incorporate shallow, cloud-free thermals that primarily serve to mix the air.

 Intermediate tests have allow to assess the impact of surface evaporation flux on cold pool moisture by activating splitting, which differentiates this flux between (w) and (x). In the standard configuration, this flux is treated uniformly for both regions. The tests showed a limited effect of this flux on cold pools moisture for RCE. This test was not carried out for AMMA, as LMDZ does not yet allow it on the continent. However, it would be relevant to explore it.

 LMDZ V3 simulations show a cooling of cold pools by thermals for the RCE case, although cold pools remain less cold compared to the LES (Fig. \ref{profil_LES_LMDZ_MODIF}a). This cooling is at the cause of an increase in the variables $WAPE$, $C_{*}$, $ALE_{wk}$ and $ALP_{wk}$ in this case (Table \ref{LES_LMDZ_MODIF}). For the AMMA case, the impact of this modification on the $\delta{T}$ profile as well as on the associated variables remains very limited (\fig{profil_LES_LMDZ_MODIF}e,Table \ref{LES_LMDZ_MODIF}). But 3D tests also indicate a cooling of cold pools on the continent when thermals are activated in the whole domain.

    \begin{table}[ht]
    %\begin{center}
     \caption{Description of simulations performed with LMDZ in the standard configuration and with various modifications}
     \label{param_tun}
     \centering
     %\begin{adjustbox}{max width=\textwidth} 
     \small
     %\begin{tabularx}{\textwidth}{|X|X|}
     % \begin{tabularx}{|p{2cm}|X|}
     \begin{tabularx}{\textwidth}{|p{3.5cm}|X|}
     %\begin{tabularx}{\textwidth}{|l|l|l|l|l|l|l|l|l|} (\cite{de2021subseasonal})
     \hline
     \textbf{Simulations} & \textbf{Protocols} \\ \hline
	     LMDZ CTRL & simulation of LMDZ with the standard configuration by imposing $D_{wk}$ to 5$10^{-10}$  \\ \hline
     LMDZ V1 & LMDZ CTRL + change of $k$ to 0.66 \\ \hline
     LMDZ V2 & LMDZ V1 + drop of $\pupper$  \\ \hline
     LMDZ V3 & LMDZ V2 + activation of thermals throughout the domain  \\ \hline
    \end{tabularx}
    %\end{adjustbox}
   %\end{center}
    \end{table}
        

     \begin{table}[ht]
    %\begin{center}
     \caption{Comparison of the variables $WAPE$, $ALE_{wk}$, $C_{*}$ and $ALP_{wk}$ calculated from the samplings in the LES, with those simulated in LMDZ control (LMDZ CTRL), LMDZ with the adjustment of the coefficient k to 0.66 (V1), LMDZ with the drop in altitude ($\pupper$) at which the subsidence of the air masses in cold pools is zero (V2) and LMDZ with the activation of thermals in the entire domain (V3) on the RCE case and on the AMMA case.}
     \label{LES_LMDZ_MODIF}     
     %\centering
    % \begin{adjustbox}{max width=\textwidth} 
     \small
     \begin{tabularx}{\textwidth}{|X|X|X|X|X|}
     %\begin{tabularx}{\textwidth}{|C|C|C|C|C|}
     %\begin{tabularx}{\textwidth}{|l|l|l|l|l|l|l|l|l|} (\cite{de2021subseasonal})
     \hline
	     \textbf{} & \textbf{$WAPE$ }(J/Kg) & \textbf{$ALE_{wk}$ }(J/kg) & \textbf{$C_{*}$} (m/s) &  \textbf{$ALP_{wk}$ }($W/m^{2}$) \\ \hline
    \multicolumn{5}{|c|}{RCE} \\ \hline
    LES SAM & 7.962 & 10.460 & 2.228 & 0.054 \\ \hline
    LES MESONH & 7.912 & 6.965 & 2.264 & 0.020 \\ \hline
    LMDZ CTRL & 2.957 & 2.957 & 0.802 & 0.001  \\ \hline
    LMDZ V1 & 2.528 & 2.528 & 1.484 & 0.006  \\ \hline
    LMDZ V2 & 2.465 & 2.465 & 1.465 & 0.006  \\ \hline
    LMDZ V3 & 3.408 & 3.408 & 1.723 & 0.009  \\ \hline
    \multicolumn{5}{|c|}{AMMA} \\ \hline
    LES MESONH & 34.250 & 33.480 & 4.939 & 0.982 \\ \hline
    LMDZ CTRL & 30.430 & 30.430 & 2.574 & 0.042  \\ \hline
    LMDZ V1 & 22.020 & 22.020 & 4.380 & 0.479 \\ \hline
    LMDZ V2 & 20.580 & 20.580 & 4.234 & 0.399 \\ \hline
    LMDZ V3 & 20.640 & 20.640 & 4.240 & 0.404 \\ \hline
    \end{tabularx}
   % \end{adjustbox}
   %\end{center}
   %\footnotes
   \end{table}    


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    \subsection{Effect of changes on large-scale variables}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

    Although the modifications presented above have improved the representation of cold pools in the model, it is also essential to examine their impact on large-scale variables. In this section, we analyze the effect of these adjustments on variables such as potential temperature ($\theta$) and specific humidity ($qv$) profiles. For this, the same profiles are recalculated in the LES for the RCE and AMMA cases and then compared with those obtained in each modified version of the LMDZ model. The profiles of $\theta$ and $qv$ in the LES are calculated by a horizontal average of these variables over the domain.
    \\
    \begin{figure}
     %\centering
     \includegraphics[width=16cm,height=8cm]{figures/fig_profilsmoyens_RCE.png}
     \includegraphics[width=16cm,height=8cm]{figures/fig_profilsmoyens_AMMA.png}
        \caption{Vertical profiles of potential temperature ($\theta$) and specific humidity ($qv$) calculated in the LES and simulated in control LMDZ (LMDZ CTRL), LMDZ with the adjustment of the coefficient k to 0.66 (V1), LMDZ with the drop in altitude ($\pupper$) at which the subsidence of the air masses in cold pools is zero (V2) and LMDZ with the activation of thermals in the entire domain (V3) on the RCE case (a, b, c) and on the AMMA case (d, e, f).}
        \label{profilmoyen_LES_LMDZ}
 \end{figure}
    
	    \fig{profilmoyen_LES_LMDZ} shows that, in the RCE case, the modifications introduced in versions V1, V2 and V3 have a low impact on the $\theta$ profiles (\fig{profilmoyen_LES_LMDZ}a), while in the AMMA case, their influence remains negligible (\fig{profilmoyen_LES_LMDZ}c). Overall, all simulations (LMDZ CTRL, V1, V2 and V3) manage to reproduce the $\theta$ profiles well, both for RCE and for AMMA, although a slight warm bias is observed around 200~hPa in the RCE case. In contrast, versions V2 and V3 lead to a drying of the mid-troposphere in the RCE case (\fig{profilmoyen_LES_LMDZ}b), while for AMMA, the three modified versions have little effect on the humidity profiles (\fig{profilmoyen_LES_LMDZ}d). However, it appears that the LMDZ model insufficiently reproduces the humidity profiles in all simulations, both for RCE and AMMA. In the RCE case, a dry bias is clearly visible in the boundary layer, as well as between 800 and 400~hPa, for versions CTRL, V1, V2 and V3. For the AMMA case, a wet bias is observed in the boundary layer and above 600~hPa, while a dry bias is present between 700 and 600~hPa.
To correct these biases observed in the model, as well as the simulated hot cold pools in RCE case, a tuning of parameter is performed


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    \section{Tuning of free parameters} \label{modif3}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

	    Since the cold pool model is coupled to the deep convection model, the $\delta{T}$ profile could be influenced by the latter. Thus, more than twenty free parameters were selected, including those related to the cold pool model as well as to the convection schemes. The metrics retained are the $\delta{T}$, $q_v$ and $\theta$ profiles, evaluated through vertical averages at different levels and time averages (between 5 P.M. and 6 P.M. for the AMMA case, and between days 41 and 43 for RCE). In these tuning tests, only modifications affecting the $k$ and $\pupper$ coefficients were integrated. Adjustments related to thermals are not taken into account here, because, in LMDZ v3, our objective is only to highlight their role in the drying of surface cold pools. We suggest that further parameterization work address this point.

	    The results obtained after tuning show a clear improvement in the representation of the $\delta{T}$ profile in the RCE case (\fig{profilpoches_tun}a). The $\delta{q}$ and $\delta{w}$ profiles also remain well represented for this case (\fig{profilpoches_tun}b,c). \fig{profilmoyens_tun}b and \fig{profilmoyens_tun}d also illustrate a significant improvement in the humidity profiles for the RCE and AMMA cases. However, a humid bias persists in the boundary layer for the AMMA case (\fig{profilmoyens_tun}d), accompanied by a slight cooling in this same layer (\fig{profilmoyens_tun}c). A slight increase in the temperature of the cold pool and a drying at their top are also observed in AMMA after tuning (\fig{profilpoches_tun}d,c). However, the fact to obtain a better representation of the humidity profiles in both cases, while improving the representation of the cold pools, particularly in the RCE case, constitutes a satisfactory result. The values of the optimized parameters resulting from this tuning can now be used in the 3D version of the LMDZ model.
   \\
    \begin{figure}
     %\centering
     \includegraphics[width=16cm,height=8cm]{figures/fig_profilspoches_V123TUN_RCE.png}
     \includegraphics[width=16cm,height=8cm]{figures/fig_profilspoches_V123TUN_AMMA.png}
	    \caption{Vertical profiles of $\delta{T}$, $\delta{q}$ and $\delta{w}$ calculated in the LES and simulated in LMDZ TUNING, LMDZ with the adjustment of the coefficient k to 0.66 (V1), LMDZ with the drop in altitude ($\pupper$) at which the subsidence of the air masses in cold pools is zero (V2) and LMDZ with the activation of thermals in the entire domain (V3) on the RCE case (a, b, c) and on the AMMA case (d, e, f).}
        \label{profilpoches_tun}
 \end{figure}
    \\
    \begin{figure}
     %\centering
     \includegraphics[width=16cm,height=8cm]{figures/fig_profilsmoyens_V123TUN_RCE.png}
     \includegraphics[width=16cm,height=8cm]{figures/fig_profilsmoyens_V123TUN_AMMA.png}
	    \caption{Vertical profiles of potential temperature ($\theta$) and specific humidity ($qv$) calculated in the LES and simulated in LMDZ TUNING, with the adjustment of the coefficient $k$ to 0.66 (V1), LMDZ with the drop in altitude ($\pupper$) at which the subsidence of the air masses in cold pools is zero (V2) and LMDZ with the activation of thermals in the entire domain (V3) on the RCE case (a, b, c) and on the AMMA case (d, e, f).}
        \label{profilmoyens_tun}
 \end{figure}


    



% \begin{figure}
     %\centering
    % \includegraphics[width=16cm,height=8cm]{FIG_PAPIER/fig_deltaT_Q_W_LES_SAM24_MESONH2_LMDZ_MODIF1234.png}
      %  \includegraphics[width=16cm,height=8cm]{FIG_PAPIER/fig_deltaT_Q_W_LES_LMDZ_AMMA_LMDZ_MODIF1234.png}
     %   \caption{Profils verticaux de $\delta{T}$, $\delta{q}$ et $\delta{w}$ calculés dans les LES et simulées dans LMDZ contrôle (LMDZ CTRL), LMDZ avec l'ajustement du coefficient k à 0.66 (LMDZ MODIF1), LMDZ avec la baisse de l'altitude maximale de subsidence des masses d'air (LMDZ MODIF2), LMDZ avec l'activation des thermiques dans les poches (LMDZ MODIF3) et LMDZ avec l'augmentation de l'évaporation des précipitations (LMDZ MODIF4) sur l'océan (a, b, c) et sur le continent (d, e, f).}
      %  \label{profil_LES_LMDZ_TUNING}
 %\end{figure}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    \section{Conclusions}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
   
   Although the cold pool model proposed by \cite{GL10I} has improved the representation of convection in the LMDZ climate model \cite{rio2009}, its internal variables and physical properties have never been evaluated. This work proposes, for the first time, a detailed evaluation of the cold pool model, based on explicit simulations called LES. We evaluate both the physics of the model, its internal variables and those involved in the coupling with deep convection, based on two oceanic LES in the RCE regime and a continental LES of the AMMA case. For this, we first performed sampling in the LES, separating the interior and exterior of cold pools on the RCE and AMMA cases by surface temperature anomalies lower than -0.2~K and -1~K respectively, in order to calculate the targeted variables. The internal variables analyzed include the profiles of temperature ($\delta{T}$), humidity ($\delta{q}$) and vertical velocity ($\delta{\omega}$) differences between the inside and outside of cold pools, the collaspe energy ($WAPE$), the spreading speed ($C_{*}$), as well as the Available Liftting Energy ($ALE_{wk}$) and Power ($ALP_{wk}$) variables related to cold pools for the coupling with deep convection. 

    We first validated the physics of the cold pool model for calculations of $ALE_{wk}$, $C_{*}$ and $ALP_{wk}$ based on the $WAPE$. For this, these three variables were recalculated in the LES using the $WAPE$, derived from the sampled $\delta{\theta_{v}}$ profiles, according to the parameterization. The values obtained were then compared to those calculated from the divergence of wind at 10 m inside cold pools (for $C_{*}$) and the vertical velocities ($w_{bgust}$) at the cloud base at the gust fronts (for $ALE_{wk}$ and $ALP_{wk}$), also sampled in the same LES. The results show that the $ALE_{wk}$ calculated from the $WAPE$ is comparable to that estimated from the $w_{bgust}$. This result is consistent with the model hypothesis, which estimates an equality between $ALE_{wk}$ and $WAPE$. The spreading speed ($C_{*}$), determined from the mean of divergence of wind at 10 m inside cold pools, is consistent with the estimate based on the square root of $WAPE$. The proportionality coefficient $k$, evaluated here at 0.66, is consistent with the work of Lafore and Moncrieff (1998), and differs from the initially assumed value of 0.33 in the model. $ALP_{wk}$, calculated using $C_{*}$ from the $WAPE$ (with $k = 0.66$), is close to the estimate derived directly from the $w_{bgust}$. This result is compatible with the model hypothesis according to which $ALP_{wk}$ translates a transformation of the horizontal power of the pockets into a vertical power, with a conversion coefficient of 25$\%$. All of these results show the overall consistency of the model hypotheses with the three LES (RCE and AMMA) used in this study.

   We then compared the variables simulated by the model to those calculated in the LES by performing a simulation with a version single-column of LMDZ for the RCE and AMMA cases, using the same initial and boundary conditions as the LES. The results show that the initial version of the parameterization represents the cold pool properties well to first order, even if some biases could be identified. The dry bias observed at the top of the pockets in the model is attributed to a maximum subsidence altitude ($\pupper$) imposed at a too high level in the initial parameterization (above 600 hPa, both on the ocean and on the continent). In the LES, $\pupper$ is observed at lower altitudes: below 800~hPa for the RCE case, and below 600~hPa for the AMMA case. By making $\pupper$ dependent on $\ptop$ to account for its regional variation, and lowering it to a level consistent with that observed in LES, we significantly improve the simulated humidity at the top of cold pools for both cases. These results highlight the significant impact of descending air masses in cold pools on the vertical humidity profile. They also confirm the relevance of the physical model, which postulates impermeability of cold pools below the top and penetration of dry air only above this level, canceling out at a certain altitude. The differences observed in the profiles of $\delta_w$ below $\ptop$ in the RCE case between the model and the LES are corrected when we impose, in the model, a value of the coefficient $k$ equal to 0.66 instead of 0.33, as suggested by the LES. This change significantly improves the representation of $C_{*}$, which directly influences the profile of $\delta_{w}$, depending of this speed below $\ptop$. This modification of $k$ also significantly improved the representation of $ALP_{wk}$ in the model. A wet bias is also noted at the surface of cold pools in the RCE and AMMA cases. Our analyses show that this bias is linked to the absence, in the model, of the effect of thermals on the variation of humidity at the surface of cold pools. The evaporation flux plays a weak role in this variation, which seems to be mainly controlled by thermals. To account for the effect of thermals on humidity variation in cold pools, we propose to introduce a parameterization of shallow thermals (not producing clouds) inside cold pools. We also find that, in the RCE case, cold pools are warmer in the model than those in the LES. Given the coupling between the cold pools model and the deep convection model, it is likely that the cold pools temperature is influenced by convective processes. To correct for this warm bias, we conducted a calibration experiment using the HighTune tool to jointly adjust the free parameters of the cold pools and deep convection models. This tuning procedure also aimed to correct the dry and wet biases still present in the potential temperature and specific humidity profiles, despite the improvements made to the cold pools model. These adjustments led to a significant improvement in the representation of cold pool temperature in the RCE case, as well as specific humidity for the RCE and AMMA cases, even if a humid bias persists in the boundary layer for the AMMA case. Changes to the cold pocket model are now incorporated into the 3D version of LMDZ. Calibrated values of the free parameters for cold pools and convection can also be incorporated into the 3D version of LMDZ. 

    Although significant progress has been made in recent years in modeling cold pools, due to their important role in convection, challenges remain. For example, the life cycle of cold pools, including their birth, death or collisions, needs to be addressed. After highlighting the impact of thermals on humidity variations within cold pools, we encourage the development of a parameterization of thermals capable of taking into account their influence without leading to cloud formation. The issue of the propagation of cold pools from grid cell to grid cell needs to be also integrated into GCMs, as well as the wind gusts associated with their spreading.

\bibliography{bibliographie}



\appendix
\section{Tungin experiment}
    
    \subsection{High-Tune Explorer (HTExplo) tool}

    General circulation models, used for global warming projections, are essentially based on a separation between the dynamical core, which manages large-scale air movements, and the physical parameterizations, enabling the impact of subgrid processes on the large scale to be represented. Progress in improving these models has been slow in recent years, not only because of the difficulties of integrating these processes into the parameterizations, but also because of the complex tuning of the many free parameters involved in their formulation. This is the background to the development of the High-Tune Explorer (HTExplo) tool.

    HTExplo has been developed in collaboration between the LMD (Paris), the Centre National de Recherche Météorologiques (CNRM/Météo-France) and the University of Exeter (UK). It is an automatic calibration tool for free parameters, based on machine learning techniques from the uncertainty quantification community \citep{williamson2013}. This approach proposes a new calibration paradigm: instead of optimizing parameter values, it aims to identify the subset of parameters that enables the model to reproduce certain observables to a certain accuracy. The main steps involved in using the tool, as well as its mathematical foundations, are well described in \cite{couvreux2021}. The HTExplo tool was used for the first time in a SCM/LES comparison on several boundary layer cases of the LMDZ model, in order to characterize the subspace of free parameter values for which SCM simulations are consistent with LES for certain metrics and a given tolerance \citep{couvreux2021}. This information was then used by \cite{hourdin2021} to calibrate the 3D configuration. These authors demonstrated how reducing the parameter space using this method significantly saves computing and human resources. They also pointed out that this approach eases the burden on the modeler, enabling him or her to concentrate more on understanding and improving the physical parameterizations of the model.

\end{document}
