# HG changeset patch
# User vlabarre <vincent.labarre@oca.eu>
# Date 1680282475 -7200
#      Fri Mar 31 19:07:55 2023 +0200
# Node ID 934f6e88fe042f9e7a8096de69e5d5003c5e4b42
# Parent  81dc1d2a5616255b9cffd683db164dcdfe4e215b
section mixing in progress

diff --git a/2022strat_turb_polo/input/results.tex b/2022strat_turb_polo/input/results.tex
--- a/2022strat_turb_polo/input/results.tex
+++ b/2022strat_turb_polo/input/results.tex
@@ -22,7 +22,7 @@
 where $\Idiss=1$, i.e. $\epsK_z = \epsK/3$, corresponds to a small-scale isotropy
 where the kinetic energy is dissipated equally along each direction. Similar coefficients have been used as diagnostics of 2D numerical simulations of stratified turbulence \cite{Linares2020numerical}. \\
 
-\begin{figure}[H]
+\begin{figure}
 	\centerline{
 		\includegraphics[width=0.48\textwidth]{%
 			../tmp/fig_isotropy_velo_vs_Fh}
@@ -54,7 +54,7 @@
 two-dimensional simulations of \cite{Linares2020numerical}. The limit $F_h=1$ defines the transition
 to another regime in which the buoyancy behaves as a passive scalar. Then, there is a range $F_{hc} \simeq 0.13 \leq F_h \leq 1$ for a weakly stratified regime. Finally, viscosity is strongly important at all scales typically when $Re \lesssim 200$, such that corresponding simulations should belong to other flow regimes. 
 
-\begin{figure}[H]
+\begin{figure}
 \centerline{
 \includegraphics[width=0.7\textwidth]{%
 ../tmp/fig_isotropy_coef_vs_FhR}
@@ -63,7 +63,7 @@
 \end{figure}
 
 
-% I think that the following part is not useful
+% I think that the following part is not useful since it leads to the same conclusions
 \iffalse
 Another way to quantify the anisotropy of the flow is to look at the integral scales
 which are defined as
@@ -81,7 +81,7 @@
 \end{equation}
 
 
-\begin{figure}[H]
+\begin{figure}
 \centerline{
 \includegraphics[width=0.48\textwidth]{%
 ../tmp/fig_ratio_velo_vs_Fh}
@@ -107,50 +107,59 @@
 \label{sec:res-mixing}
 
 The mixing coefficient can be defined in several ways. In geophysical applications, it
-is defined as the ratio between the buoyancy flux and the kinetic energy dissipation
+is defined as the ratio between the buoyancy flux and the kinetic energy dissipation rates
 \cite{Gregg-etal2018mixing}. Doing so, the mixing coefficient is simply related to the
 turbulent diffusion that is needed in parametrization in climate models. Yet, if we want to quantify only the irreversible
 part of the mixing, it is better to define the mixing coefficient as the ratio of
-dissipations \cite{Maffioli2016mixing, Caulfield2020review}:
+dissipation rates \cite{Maffioli2016mixing, Caulfield2020review}:
 \begin{equation}
 \Gamma \equiv \frac{\epsA}{\epsK},
 \end{equation}
 
-which is a positive quantity. Figure~\ref{fig:mixing-coefficients-vs-FhR} represents the mixing ratio vs $(F_h, \R)$ in our simulations
-We first observe that the mixing coefficient becomes independent of $\R$ when $\R \gg 1$. This is consistent with the numerical simulations of stratified sheared flows \cite{Portwood-BruynKops-Caulfield2019} where no significant variation of $\Gamma$ up to $\R \simeq 10^3$ was observed. Secondly, for high $\R$, $\Gamma$ stays almost constant in the strongly stratified regime ($F_h \ll 1$), before increasing to an optimal value in the weakly stratified regime (around $F_h \simeq F_{hc} = 0.13$), and then decreasing for higher $F_h$. 
+\noindent which is a positive quantity. \\
 
-\begin{figure}[H]
+\begin{figure}
 \centerline{
 \includegraphics[width=0.7\textwidth]{%
 ../tmp/fig_mixing_coef_vs_FhR}
 }
-\caption{Mixing coefficient. Red letters correspond to
+\caption{Mixing coefficient vs $(F_h, \R)$. Red letters correspond to
 	simulations of table~\ref{table-simuls-regimes} analyzed in
-	subsection~\ref{spectra-seb-regimes}. The dashed blue lines are the same as in figure \ref{fig:isotropy-coefficients}. \label{fig:mixing-coefficients-vs-FhR}}
+	subsection~\ref{spectra-seb-regimes}. The dashed blue lines and the black line are the same as in figure \ref{fig:isotropy-coefficients}. \label{fig:mixing-coefficients-vs-FhR}}
+\end{figure}
+
+Figure~\ref{fig:mixing-coefficients-vs-FhR} represents the mixing coefficient vs $(F_h, \R)$ in our simulations
+We first observe that the mixing coefficient becomes independent of $\R$ when $\R \gg 1$. This is consistent with the numerical simulations of stratified sheared flows \cite{Portwood-BruynKops-Caulfield2019} where no significant variation of $\Gamma$ up to $\R \simeq 10^3$ was observed. Secondly, for high $\R$, $\Gamma$ stays almost constant in the strongly stratified regime ($F_h \ll 1$), before increasing to an optimal value in the weakly stratified regime (around $F_h \simeq F_{hc} = 0.13$), and then decreasing to zero for higher $F_h$. Interestingly, the largest values of the mixing ratio are obtained for viscous stratified flows ($F_h, \R \ll 1$).\\
+
+
+\begin{figure}
+	\centerline{
+	\includegraphics[width=0.48\textwidth]{%
+		../tmp/fig_mixing_coef_vs_Fh}
+	\includegraphics[width=0.48\textwidth]{%
+		../tmp/fig_ratio_Etoro_E_vs_Fh} 
+	}
+	\centerline{
+		\includegraphics[width=0.48\textwidth]{%
+			../tmp/fig_ratio_EA_E_vs_Fh}
+		\includegraphics[width=0.48\textwidth]{%
+			../tmp/fig_ratio_Epolo_E_vs_Fh}
+	}
+	\caption{Mixing coefficient (a), toroidal energy ratio (b), potential energy ratio (c), and poloidal energy ratio (d) versus the horizontal Froude number. The vertical dashed lines represent $F_h = F_{hc}  = 0.13$, for which we observe the optimal mixing when $\R \gg 1$, and $F_h =1$. Only simulations for which $\R > 10$ are shown. }
+	\label{fig:mixing-coefficients-vs-Fh}
 \end{figure}
 
 
-
-\begin{figure}[H]
-\centerline{
-\includegraphics[width=0.48\textwidth]{%
-	../tmp/fig_mixing_coef_vs_Fh}
-\includegraphics[width=0.48\textwidth]{%
-	../tmp/fig_ratio_Etoro_E_vs_Fh} 
-}
-\centerline{
-	\includegraphics[width=0.48\textwidth]{%
-		../tmp/fig_ratio_Epolo_E_vs_Fh}
-	\includegraphics[width=0.48\textwidth]{%
-		../tmp/fig_ratio_EA_E_vs_Fh}
-}
-\caption{Mixing coefficient (a), ratio of potential energy and total energy (b), ratio of poloidal energy and total energy (c), and ratio of toroidal energy and total energy (d) versus the horizontal Froude number. The vertical dashed lines represent $F_h = F_{hc}  = 0.13$, for which we observe the optimal mixing when $\R \gg 1$, and $F_h =1$. Only simulations for which $\R > 10$ are shown. }
-\label{fig:mixing-coefficients-vs-Fh}
-\end{figure}
-
-
-
-In figure~\ref{fig:mixing-coefficients-vs-Fh}, we plot the mixing coefficient and energies ratios as a function of $F_h$ for simulations with $\R > 10$. The evolution of the mixing coefficient as a function of $F_h$ in the strongly stratified and passive scalar regimes was described in \cite{Maffioli2016mixing, Garanaik-Venayagamoorthy2019mixing}. The authors showed that the mixing efficiency coincides with an optimal value of $E_{\rm A}/E_{\rm K}$. We refer the reader to section 2.1 of \cite{Garanaik-Venayagamoorthy2019mixing} for the physical explanations of the different scaling laws. Our simulations are fully consistent with these predictions (\ref{eq:Maffioli-Garanaik-Venayagamoorthy}). Interestingly, we observe that limit $\lim\limits_{F_h \rightarrow 0} \Gamma \simeq 0.42$ which is different from the value measured in \cite{Maffioli2016mixing} with a vertical vorticity, uncorrelated forcing. Also, the position of the optimal mixing is slightly different ($F_h \simeq 0.13$ in our simulations vs $F_h \simeq 0.3$ in \cite{Maffioli2016mixing}). We observe that the poloidal energy ratio $E_{\rm A}/E$ increases monotonically with $F_h$. The toroidal energy represents almost all the total energy ($E_{toro}/E \rightarrow 1$) when $F_h \ll 1$. $E_{toro}/E$ first decreases with $F_h$ to a minimum whose position coincides with the optimal mixing. Interestingly, the mixing efficiency also coincide with the threshold $F_{hc} = 0.13$ that was defined by the condition $\Ivelo(F_{hc}, \R \gg 1) \simeq 0.5$ (wich also coincides with the optimum of potential energy $E_{\rm A}/E$ \cite{Maffioli2016mixing}). Then, $E_{toro}/E$ increases with $F_h$ to reach a plateau when $F_h \gg 1$. The fact that the optimal mixing coincides with a minimum of toroidal energy is not surprising since the toroidal velocity field dissipates kinetic energy and cannot be converted into potential energy ``directly'', while the poloidal kinetic energy can. In the passive scalar regime, we observe that $\lim\limits_{F_h \gg 1} E_{polo}/E \simeq 0.6 > 1/2 > \lim\limits_{F_h \gg 1} E_{toro}/E \simeq 0.4$. This is due to the fact that we force the poloidal velocity. In the case where the toroidal velocity is forced, we rather obtain $\lim\limits_{F_h \gg 1} E_{toro}/E > 1/2 > \lim\limits_{F_h \gg 1} E_{polo}/E$ \todo{cite the manuscript for toroidal dateset}.
+\todo{continue here}
+In figure~\ref{fig:mixing-coefficients-vs-Fh}, we plot the mixing coefficient and energy ratios as a function of $F_h$ for simulations with $\R > 10$. 
+We observe that the mixing efficiency (i.e. the maximum of $\Gamma$) coincides with both a maximum of $E_{\rm A}/E$ \cite{Maffioli2016mixing} and a minimum of $E_{\rm toro}/E$. 
+We refer the reader to section 2.1 of \cite{Garanaik-Venayagamoorthy2019mixing} for the physical explanations of the different scaling laws. Our simulations are fully consistent with these predictions (\ref{eq:Maffioli-Garanaik-Venayagamoorthy}). 
+Interestingly, we observe that limit $\lim\limits_{F_h \rightarrow 0} \Gamma \simeq 0.42$ which is different from the value measured in \cite{Maffioli2016mixing} with a vertical vorticity, uncorrelated forcing. 
+Also, the position of the optimal mixing is slightly different ($F_h \simeq 0.13$ in our simulations vs $F_h \simeq 0.3$ in \cite{Maffioli2016mixing}). 
+We observe that the poloidal energy ratio $E_{\rm A}/E$ increases monotonically with $F_h$. The toroidal energy represents almost all the total energy ($E_{toro}/E \rightarrow 1$) when $F_h \ll 1$. 
+$E_{toro}/E$ first decreases with $F_h$ to a minimum whose position coincides with the optimal mixing. Interestingly, the mixing efficiency also coincide with the threshold $F_{hc} = 0.13$ that was defined by the condition $\Ivelo(F_{hc}, \R \gg 1) \simeq 0.5$ (wich also coincides with the optimum of potential energy $E_{\rm A}/E$ \cite{Maffioli2016mixing}). 
+Then, $E_{toro}/E$ increases with $F_h$ to reach a plateau when $F_h \gg 1$. The fact that the optimal mixing coincides with a minimum of toroidal energy is not surprising since the toroidal velocity field dissipates kinetic energy and cannot be converted into potential energy ``directly'', while the poloidal kinetic energy can. In the passive scalar regime, we observe that $\lim\limits_{F_h \gg 1} E_{polo}/E \simeq 0.6 > 1/2 > \lim\limits_{F_h \gg 1} E_{toro}/E \simeq 0.4$. 
+This is due to the fact that we force the poloidal velocity. In the case where the toroidal velocity is forced, we rather obtain $\lim\limits_{F_h \gg 1} E_{toro}/E > 1/2 > \lim\limits_{F_h \gg 1} E_{polo}/E$ \todo{cite the manuscript for toroidal dateset}.
 
 
 
@@ -204,7 +213,7 @@
 \subsubsection{One dimensional integrated spectra vs $(F_h, \R)$}
 
 
-\begin{figure}[H]
+\begin{figure}
 	\centerline{
 		\includegraphics[width=0.98\textwidth]{%
 			../tmp/fig_spectra_1strat}
@@ -215,7 +224,7 @@
 Figure~\ref{fig:spectra-1strat} represents the compensated one-dimensional spectra 
 \begin{align}
 	\label{eq:1dspectraA}
-	&E_{\rm A}(k_h) = \delta k_z ~ \sum\limits_{k_z} ~ E_{\rm A}(k_h, k_z) , ~~~  E_{\rm pot}(k_z) = \delta k_h ~ \sum\limits_{k_h} ~ E_{\rm A}(k_h, k_z) \\
+	&E_{\rm A}(k_h) = \delta k_z ~ \sum\limits_{k_z} ~ E_{\rm A}(k_h, k_z) , ~~~  E_{\rm A}(k_z) = \delta k_h ~ \sum\limits_{k_h} ~ E_{\rm A}(k_h, k_z) \\
 	\label{eq:1dspectrapolo}
 	&E_{polo}(k_h) = \delta k_z ~ \sum\limits_{k_z} ~ E_{\rm A}(k_h, k_z) , ~~~  E_{polo}(k_z) = \delta k_h ~ \sum\limits_{k_h} ~ E_{polo}(k_h, k_z) \\
 	\label{eq:1dspectratoro}
@@ -223,7 +232,7 @@
 \end{align}
 for $N=20$ and various $\R$. For all the spectra, we naturally observe that increasing $\R$ allows the energy transfer towards smaller and smaller horizontal and vertical scales. We observe a pic in both potential and poloidal energy spectra at the large forced scales. This pic is not present on toroidal spectra, showing that the injected poloidal kinetic energy is efficiently converted into potential energy. When $\R \simeq 10$, the dissipative range starts to come after the Ozmidov wavevector, as observed in previous studies (e.g. \cite{Kimura-Herring2012spectra}). The ratio $k_b/k_O$ slightly changes with $\R$, but converge to a constant value when $\R$ increases. For a same $\R$, poloidal, toroidal and potential energy spectra are remarquably similar above $k_O$, indicating that this range of scales should not be deeply affected by stratification. However, the spectra do present differences at $k_h, k_z < k_O$: (i) We observe a bounce in the buoyancy vertical energy spectra around $k_b$ that is not present in the poloidal and toroidal spectra. (ii) The compensated poloidal and toroidal horizontal energy spectra are steeper than the potential energy spectra before $k_b$.
 
-\begin{figure}[H]
+\begin{figure}
 	\centerline{
 		\includegraphics[width=0.98\textwidth]{%
 			../tmp/fig_spectra_1R}
@@ -250,7 +259,7 @@
 
 
 % regime P
-\begin{figure}[H]
+\begin{figure}
 	\centerline{
 		\includegraphics[width=0.98\textwidth]{%
 			../tmp/fig_seb_regimes_new_P}
@@ -259,7 +268,7 @@
 		table~\ref{table-simuls-regimes}). Ratios of energies (left) and transfer terms (\ref{eq:transferK}- \ref{eq:convK2toA}) (right). The dashed orange lines correspond to the maximal horizontal and vertical wavenumbers of the forcing region. The dashed green line corresponds to the dissipative scale. The dotted and dashed black lines correspond respectively to $k_b$ and $k_O$ (not visible for this simulation). The continuous black line indicates $k_z = k_h$. the dashed gray line correspond to $\chi_{\kk} = 3$. \label{fig:seb-P}}
 \end{figure}
 
-\begin{figure}[H]
+\begin{figure}
 	\centerline{
 		\includegraphics[width=0.55\textwidth]{%
 			../tmp/fig_spectra_slices_regime_P_3D}
@@ -273,7 +282,7 @@
 
 
 % regime W
-\begin{figure}[H]
+\begin{figure}
 \centerline{
 \includegraphics[width=0.98\textwidth]{%
 ../tmp/fig_seb_regimes_new_W}
@@ -281,7 +290,7 @@
 \caption{$(k_h, k_z)$ spectral energy budget for simulation W. See figure~\ref{fig:seb-P}. \label{fig:seb-W}}
 \end{figure}
 
-\begin{figure}[H]
+\begin{figure}
 	\centerline{
 		\includegraphics[width=0.55\textwidth]{%
 			../tmp/fig_spectra_slices_regime_W_kfkb}
@@ -296,7 +305,7 @@
 
 
 % regime O
-\begin{figure}[H]
+\begin{figure}
 \centerline{
 \includegraphics[width=0.98\textwidth]{%
 ../tmp/fig_seb_regimes_new_O}
@@ -304,7 +313,7 @@
 \caption{$(k_h, k_z)$ spectral energy budget for simulation O. See figure~\ref{fig:seb-P}. \label{fig:seb-O}}
 \end{figure}
 
-\begin{figure}[H]
+\begin{figure}
 \centerline{
 \includegraphics[width=0.55\textwidth]{%
 ../tmp/fig_spectra_slices_regime_O_kfkb}
@@ -319,7 +328,7 @@
  
 
 % regime L
-\begin{figure}[H]
+\begin{figure}
 \centerline{
 \includegraphics[width=0.98\textwidth]{%
 ../tmp/fig_seb_regimes_new_L}
@@ -327,7 +336,7 @@
 \caption{$(k_h, k_z)$ spectral energy budget for simulation L. See figure~\ref{fig:seb-L}. \label{fig:seb-L}}
 \end{figure}
 
-\begin{figure}[H]
+\begin{figure}
 \centerline{
 \includegraphics[width=0.55\textwidth]{%
 ../tmp/fig_spectra_slices_regime_L_kfkb}
@@ -341,7 +350,7 @@
 Figure~\ref{fig:seb-L} represents the spectral energy budget for simulation L. It is very similar to the spectral energy budget of simulation O (figure~\ref{fig:seb-O}), so we will not enter into the details. The main differences with the optimal regime are: (i) the forcing region fully lies in the buoyancy range. (ii) The potential energy is less dominant at large $k_z$. (iii) Toroidal energy is now dominant at very large horizontal scales $k_h \ll k_b$, supporting the fact that the toroidal energy tends to accumulate at small $k_h$ when we increases stratification \cite{Smith2002generation, Herbert2016cascade}. Figure~\ref{fig:spectra-L} represents slices of $(k_h, k_z)$ spectra for simulation L. When compared to simulation O, the simulation L naturally has larger buoyancy and Ozmidov ranges since it has a lower $F_h$. Horizontal and vertical spectra can be dissociated clearly. In the buoyancy range, spectra are close to $k_h^2~k_z^0$ away from the forcing region, with almost an equipartition between potential, poloidal and toroidal energies. Slopes of vertical slices of the spectra decreases between $k_b$ and $k_O$. For $k_z > k_b$, horizontal slices goes roughly as $k_h^{-5/3}$. At small-scales ($k_h$ or $k_z \gg k_O$), spectra are consistent with an equipartition of energy, i.e. $~ k_h^1 ~k_z^0$, as in the weak, and optimal regimes (not shown). In the LAST regime, the dissipation occurs at wavevectors larger than $k_O$. Yet, our resolution does not allow to attain sufficiently large $\R$ to observe an isotropic inertial range after $k_O$. 
 
 % regime D
-\begin{figure}[H]
+\begin{figure}
 	\centerline{
 		\includegraphics[width=0.98\textwidth]{%
 			../tmp/fig_seb_regimes_new_D}
@@ -349,7 +358,7 @@
 	\caption{$(k_h, k_z)$ spectral energy budget for simulation D. See figure~\ref{fig:seb-P}. \label{fig:seb-D}}
 \end{figure}
 
-\begin{figure}[H]
+\begin{figure}
 \centerline{
 \includegraphics[width=0.55\textwidth]{%
 ../tmp/fig_spectra_slices_regime_D_kfkb}
diff --git a/2022strat_turb_polo/py/save_E_vs_time_N40_Ri20.py b/2022strat_turb_polo/py/save_E_vs_time_N40_Ri20.py
--- a/2022strat_turb_polo/py/save_E_vs_time_N40_Ri20.py
+++ b/2022strat_turb_polo/py/save_E_vs_time_N40_Ri20.py
@@ -52,7 +52,7 @@
 norm = 1
 EKp = data_hist["EKz"] + data_hist["EKhd"]
 
-ax.plot(times, data_hist["EA"] / norm, "b", label=r"$E_{\rm pot}$")
+ax.plot(times, data_hist["EA"] / norm, "b", label=r"$E_{\rm A}$")
 ax.plot(times, data_hist["EKhr"] / norm, "r-", label=r"$E_{\rm toro}$")
 ax.plot(times, EKp / norm, "g", label=r"$E_{\rm polo}$", alpha=0.7)
 ax.plot(times, data_hist["E"] / norm, "k", label="$E$")
diff --git a/2022strat_turb_polo/py/save_mixing_coef_vs_Fh.py b/2022strat_turb_polo/py/save_mixing_coef_vs_Fh.py
--- a/2022strat_turb_polo/py/save_mixing_coef_vs_Fh.py
+++ b/2022strat_turb_polo/py/save_mixing_coef_vs_Fh.py
@@ -1,11 +1,12 @@
 import numpy as np
 import matplotlib.pyplot as plt
+plt.rcParams['text.usetex'] = True
 
-from util import save_fig, plot, params_simuls_regimes, Fh_limit
+from util import save_fig, plot, params_simuls_regimes, Fh_limit, R2_limit
 
 from util_dataframe import df
 
-tmp = df[df.R2 > 10].copy()
+tmp = df[df.R2 > R2_limit].copy()
 
 ax = plot(
     tmp,
@@ -41,8 +42,8 @@
 
     #ax.text(coef_x * Fh, coef_y * Gamma, letter, color="r")
 
-ax.set_xlabel("$F_h$", fontsize=12)
-ax.set_ylabel(r"$\Gamma=\epsilon_A / \epsilon_K$", fontsize=12)
+ax.set_xlabel("$F_h$", fontsize=14)
+ax.set_ylabel(r"$\Gamma=\varepsilon_{\rm A} / \varepsilon_{\rm K}$", fontsize=14)
 
 
 ax.set_xlim((1e-3, 20))
@@ -52,21 +53,21 @@
 
 
 ax.plot([1e-3, 1e-2], [0.42, 0.42], "b-")
-ax.text(2e-3, 0.2, r"$0.42 \times F_h^0$")
+ax.text(2e-3, 0.2, r"$0.42 \times F_h^0$", fontsize=12)
 
 xs = np.linspace(0.3, 1, 4)
 ax.plot(xs, 5e-2 * xs**-1, "g-")
-ax.text(0.3, 0.03, "${F_h}^{-1}$")
+ax.text(0.3, 0.03, r"${F_h}^{-1}$", fontsize=12)
 
 xs = np.linspace(1, 10, 4)
 ax.plot(xs, 5e-2 * xs**-2, "r-")
-ax.text(1.5, 0.001, "${F_h}^{-2}$")
+ax.text(1.5, 0.001, r"${F_h}^{-2}$", fontsize=12)
 
-
+ax.set_title(r"$(a)$", fontsize=14)
 
 fig = ax.figure
 fig.tight_layout()
-fig.text(0.84, 0.07, r"$\log_{10}(\mathcal{R})$", fontsize=12)
+fig.text(0.84, 0.07, r"$\log_{10}(\mathcal{R})$", fontsize=14)
 
 save_fig(fig, "fig_mixing_coef_vs_Fh.png")
 
diff --git a/2022strat_turb_polo/py/save_mixing_coef_vs_FhR.py b/2022strat_turb_polo/py/save_mixing_coef_vs_FhR.py
--- a/2022strat_turb_polo/py/save_mixing_coef_vs_FhR.py
+++ b/2022strat_turb_polo/py/save_mixing_coef_vs_FhR.py
@@ -1,4 +1,6 @@
+import numpy as np
 import matplotlib.pyplot as plt
+plt.rcParams['text.usetex'] = True
 
 from util import save_fig, plot, params_simuls_regimes, Fh_limit, R2_limit
 
@@ -33,12 +35,30 @@
 Fh_min, Fh_max = ax.get_xlim()
 ax.plot([Fh_min, Fh_limit], [R2_limit, R2_limit], linestyle=":")
 
-ax.set_xlabel(r"$F_h$", fontsize=12)
-ax.set_ylabel(r"$\mathcal{R}$", fontsize=12)
+ax.set_xlabel(r"$F_h$", fontsize=14)
+ax.set_ylabel(r"$\mathcal{R}$", fontsize=14)
+
+Fh = np.array([1e-3, 2e1])
+ax.plot(Fh, 200*Fh**2, "k-")
+
+ax.text(2e-3, 1e3, r"Strongly", fontsize=12, alpha=0.5)
+ax.text(2e-3, 3e2, r"stratified", fontsize=12, alpha=0.5)
+
+ax.text(1.5e-1, 1e4, r"Weakly", fontsize=12, alpha=0.5)
+ax.text(1.5e-1, 3e3, r"stratified", fontsize=12, alpha=0.5)
+
+ax.text(1.5e0, 3e4, r"Passive", fontsize=12, alpha=0.5)
+ax.text(1.5e0, 9e3, r"scalar", fontsize=12, alpha=0.5)
+
+ax.text(2e-3, 1e0, r"Viscosity", fontsize=12, alpha=0.5)
+ax.text(2e-3, 3e-1, r"affected", fontsize=12, alpha=0.5)
+
+ax.text(8e-1, 1e1, r"Viscous flows", fontsize=12, alpha=0.5)
+ax.text(8e-1, 3e0, r"$Re < 200$", fontsize=12, alpha=0.5)
 
 fig = ax.figure
 fig.tight_layout()
-fig.text(0.85, 0.07, r"$\Gamma=\frac{\epsilon_A}{\epsilon_K}$", fontsize=12)
+fig.text(0.85, 0.07, r"$\Gamma=\frac{\varepsilon_{\rm A}}{\varepsilon_{\rm K}}$", fontsize=14)
 
 save_fig(fig, "fig_mixing_coef_vs_FhR.png")
 
diff --git a/2022strat_turb_polo/py/save_ratio_Ekeys_E_vs_Fh.py b/2022strat_turb_polo/py/save_ratio_Ekeys_E_vs_Fh.py
--- a/2022strat_turb_polo/py/save_ratio_Ekeys_E_vs_Fh.py
+++ b/2022strat_turb_polo/py/save_ratio_Ekeys_E_vs_Fh.py
@@ -1,16 +1,17 @@
 from curses import keyname
 import numpy as np
 import matplotlib.pyplot as plt
+plt.rcParams['text.usetex'] = True
 
-from util import save_fig, plot, Fh_limit
+from util import save_fig, plot, Fh_limit, R2_limit
 
 from util_dataframe import df
 
 key = ["Epolo", "Etoro", "EA"]
-ylabel = [r"$E_{polo}/E$", r"$E_{toro}/E$", r"$E_{\rm pot}/E$"]
+ylabel = [r"$E_{\rm polo}/E$", r"$E_{\rm toro}/E$", r"$E_{\rm A}/E$"]
 
 for n in range(3):
-    tmp = df[df.R2 > 10].copy()
+    tmp = df[df.R2 > R2_limit].copy()
     #tmp = df.copy()
     tmp["ratioE"] = tmp[key[n]] / tmp["E"]
 
@@ -18,8 +19,8 @@
         tmp, "Fh", "ratioE", c=np.log10(tmp["R2"]), vmin=1, vmax=4, logy=True
     )
 
-    ax.set_xlabel("$F_h$", fontsize=12)
-    ax.set_ylabel(ylabel[n], fontsize=12)
+    ax.set_xlabel("$F_h$", fontsize=14)
+    ax.set_ylabel(ylabel[n], fontsize=14)
 
     ax.set_xlim((1e-3, 20))
 
@@ -32,20 +33,23 @@
 
         xs = np.linspace(0.3, 1, 4)
         ax.plot(xs, 5e-2 * xs**-1, "g-")
-        ax.text(0.3, 0.03, r"${F_h}^{-1}$")
+        ax.text(0.3, 0.03, r"${F_h}^{-1}$", fontsize=12)
 
         xs = np.linspace(1, 10, 4)
         ax.plot(xs, 5e-2 * xs**-2, "r-")
-        ax.text(1., 0.003, r"${F_h}^{-2}$")
+        ax.text(1., 0.003, r"${F_h}^{-2}$", fontsize=12)
 
         ax.set_ylim((1e-3, 1e0))
+        ax.set_title(r"$(c)$", fontsize=14)
+
     elif key[n] == "Etoro":
         xs = np.linspace(1, 10, 4)
         ax.plot(xs, 5e-1 * xs**0, "r-")
-        ax.text(3., 0.55, r"$0.5$")
+        ax.text(3., 0.55, r"$0.5$", fontsize=12)
 
-        #ax.set_ylim((0.99e-1, 1e0))
+        ax.set_ylim((0.99e-1, 1e0))
         ax.set_yticks([1e-1, 1e0])
+        ax.set_title(r"$(b)$", fontsize=14)
     
     elif key[n] == "Epolo":
         #xs = np.linspace(1e-3, 1e-2, 4)
@@ -54,16 +58,17 @@
 
         xs = np.linspace(1, 10, 4)
         ax.plot(xs, 5e-1 * xs**0, "r-")
-        ax.text(3., 0.4, r"$0.5$")
+        ax.text(3., 0.4, r"$0.5$", fontsize=12)
 
         ax.set_ylim((0.99e-1, 1e0))
         ax.set_yticks([1e-1, 1e0])
+        ax.set_title(r"$(d)$", fontsize=14)
         
         
 
     fig = ax.figure
 
-    fig.text(0.84, 0.07, r"$\log_{10}(\mathcal{R})$", fontsize=12)
+    fig.text(0.84, 0.07, r"$\log_{10}(\mathcal{R})$", fontsize=14)
 
     fig.tight_layout()
     save_fig(fig, f"fig_ratio_{key[n]}_E_vs_Fh.png")
diff --git a/2022strat_turb_polo/py/save_ratio_Ekeys_E_vs_FhR.py b/2022strat_turb_polo/py/save_ratio_Ekeys_E_vs_FhR.py
--- a/2022strat_turb_polo/py/save_ratio_Ekeys_E_vs_FhR.py
+++ b/2022strat_turb_polo/py/save_ratio_Ekeys_E_vs_FhR.py
@@ -7,7 +7,7 @@
 
 
 key = ["Epolo", "Etoro", "EA"]
-ylabel = [r"$E_{polo}/E$", r"$E_{toro}/E$", r"$E_{\rm pot}/E$"]
+ylabel = [r"$E_{polo}/E$", r"$E_{toro}/E$", r"$E_{\rm A}/E$"]
 
 for n in range(3):
 
diff --git a/2022strat_turb_polo/py/save_ratio_spectra_regimes.py b/2022strat_turb_polo/py/save_ratio_spectra_regimes.py
--- a/2022strat_turb_polo/py/save_ratio_spectra_regimes.py
+++ b/2022strat_turb_polo/py/save_ratio_spectra_regimes.py
@@ -49,7 +49,7 @@
     kx,
     ratio_A / mean_values["Gamma"],
     label=(
-        r"$(E_{\rm pot} \varepsilon_{\rm K})/(E_{\rm K} \varepsilon_{\rm A})$, with "
+        r"$(E_{\rm A} \varepsilon_{\rm K})/(E_{\rm K} \varepsilon_{\rm A})$, with "
         rf"$\varepsilon_{\rm A} / \varepsilon_{\rm K} = {mean_values['Gamma']:.3f}$"
     ),
 )
diff --git a/2022strat_turb_polo/py/save_seb_regimes_new.py b/2022strat_turb_polo/py/save_seb_regimes_new.py
--- a/2022strat_turb_polo/py/save_seb_regimes_new.py
+++ b/2022strat_turb_polo/py/save_seb_regimes_new.py
@@ -70,7 +70,7 @@
 
 # EA / E (kh, kz)
 cs = ax0.contourf(KH, KZ, EA / E, cmap=cm.jet, levels=levels)
-ax0.set_title(r"$E_{\rm pot}/E$", fontsize=12)
+ax0.set_title(r"$E_{\rm A}/E$", fontsize=12)
 # cbar = fig.colorbar(cs)
 # cbar.set_ticks([0.0, 0.2, 0.4, 0.6, 0.8, 1.0])
 # cbar.set_ticks([0.0, 200, 400, 600, 800, 1000])
diff --git a/2022strat_turb_polo/py/save_spectra_kh_kz.py b/2022strat_turb_polo/py/save_spectra_kh_kz.py
--- a/2022strat_turb_polo/py/save_spectra_kh_kz.py
+++ b/2022strat_turb_polo/py/save_spectra_kh_kz.py
@@ -83,13 +83,13 @@
 
 if key == "potential":
     spectrum = EA  # Potential energy
-    key_tex = "E_{\rm pot}"
+    key_tex = "E_{\rm A}"
 elif key == "poloidal":
     spectrum = Epolo
     key_tex = "E_{\rm K}"
 elif key == "toroidal":
     spectrum = Etoro
-    key_tex = "E_{\rm pot} + E_{\rm K}"
+    key_tex = "E_{\rm A} + E_{\rm K}"
 else:
     print(f"Do not know key {key} \n")
     # continue
@@ -185,9 +185,9 @@
 plt.xlabel(xlabel, fontsize=12)
 
 if coef_compensate == 0.0:
-    plt.ylabel(r"$E_{\rm pot}(k_h, k_z)$", fontsize=12)
+    plt.ylabel(r"$E_{\rm A}(k_h, k_z)$", fontsize=12)
 else:
-    plt.ylabel(rf"$E_{\rm pot}(k_h, k_z)$" + xlabel + r"$^{5/3}$", fontsize=12)
+    plt.ylabel(rf"$E_{\rm A}(k_h, k_z)$" + xlabel + r"$^{5/3}$", fontsize=12)
 
 
 plt.scatter(
diff --git a/2022strat_turb_polo/py/save_spectra_regimes_new.py b/2022strat_turb_polo/py/save_spectra_regimes_new.py
--- a/2022strat_turb_polo/py/save_spectra_regimes_new.py
+++ b/2022strat_turb_polo/py/save_spectra_regimes_new.py
@@ -75,7 +75,7 @@
 
 # EA / E (kh, kz)
 cs = ax0.contourf(KH, KZ, EA / E, cmap=cm.jet, levels=levels)
-ax0.set_title(r"$E_{\rm pot}/E$", fontsize=12)
+ax0.set_title(r"$E_{\rm A}/E$", fontsize=12)
 # cbar = fig.colorbar(cs)
 # cbar.set_ticks([0.0, 0.2, 0.4, 0.6, 0.8, 1.0])
 # cbar.set_ticks([0.0, 200, 400, 600, 800, 1000])
diff --git a/2022strat_turb_polo/py/save_spectra_slices_regimes.py b/2022strat_turb_polo/py/save_spectra_slices_regimes.py
--- a/2022strat_turb_polo/py/save_spectra_slices_regimes.py
+++ b/2022strat_turb_polo/py/save_spectra_slices_regimes.py
@@ -198,7 +198,7 @@
             kh,
             EA_kh * kh**coef_compensate,
             "b-",
-            label=(r"$E_{\rm pot}$"),
+            label=(r"$E_{\rm A}$"),
         )
         ax.plot(
             kh,
diff --git a/2022strat_turb_polo/py/save_spectra_slices_vs_omega_kh_kz.py b/2022strat_turb_polo/py/save_spectra_slices_vs_omega_kh_kz.py
--- a/2022strat_turb_polo/py/save_spectra_slices_vs_omega_kh_kz.py
+++ b/2022strat_turb_polo/py/save_spectra_slices_vs_omega_kh_kz.py
@@ -84,13 +84,13 @@
         omegas = f['omegas'][:]
         if key == "potential":
             spectrum = f["spectrum_A"][:]  # Potential energy
-            key_tex = "E_{\rm pot}"
+            key_tex = "E_{\rm A}"
         elif key == "kinetic":
             spectrum = f["spectrum_K"][:]
             key_tex = "E_{\rm K}"
         elif key == "total":
             spectrum = f["spectrum_A"][:] + f["spectrum_K"][:]
-            key_tex = "E_{\rm pot} + E_{\rm K}"
+            key_tex = "E_{\rm A} + E_{\rm K}"
         else:
             print(f"Do not know key {key} \n")
             #continue
@@ -110,7 +110,7 @@
             plt.xscale('log')
             plt.yscale('log')
             plt.xlabel(r'$\omega / \omega_{k}$', fontsize=12)
-            plt.ylabel(r'$E_{\rm pot}(k_h, k_z, \omega) / \max_{\omega} \left\{ E_{\rm pot}(k_h, k_z, \omega) \right\}$', fontsize=12)
+            plt.ylabel(r'$E_{\rm A}(k_h, k_z, \omega) / \max_{\omega} \left\{ E_{\rm A}(k_h, k_z, \omega) \right\}$', fontsize=12)
             plt.ylim([1e-9,5e-4])
             plt.ylim([1e-3,1e0])
             plt.xlim([7e-2,1e1])
diff --git a/2022strat_turb_polo/py/save_spectra_spatiotemporal_1couple_1nh.py b/2022strat_turb_polo/py/save_spectra_spatiotemporal_1couple_1nh.py
--- a/2022strat_turb_polo/py/save_spectra_spatiotemporal_1couple_1nh.py
+++ b/2022strat_turb_polo/py/save_spectra_spatiotemporal_1couple_1nh.py
@@ -159,13 +159,13 @@
 
     if key == "potential":
         spectrum = f["spectrum_A"][:]  # Potential energy
-        key_tex = "E_{\rm pot}"
+        key_tex = "E_{\rm A}"
     elif key == "kinetic":
         spectrum = f["spectrum_K"][:]
         key_tex = "E_{\rm K}"
     elif key == "total":
         spectrum = f["spectrum_A"][:] + f["spectrum_K"][:]
-        key_tex = "E_{\rm pot} + E_{\rm K}"
+        key_tex = "E_{\rm A} + E_{\rm K}"
     else:
         print(f"Do not know key {key} \n")
         #continue
@@ -225,7 +225,7 @@
         #rf"$proj={proj}  ~~~  N={N}  ~~~  {key_tex} ~~~ k_h / \delta k_h={ikh}$"
         ax.set_title(rf"${key_tex} ~~~ k_h =$" + kh_leg_txt, fontsize =12)
         cbar.set_label(
-            r"$\log \left( E_{\rm pot}(k_z, \omega) / \int E_{\rm pot}(k_z, \omega) \mathrm{d}\omega \right)$", fontsize =12
+            r"$\log \left( E_{\rm A}(k_z, \omega) / \int E_{\rm A}(k_z, \omega) \mathrm{d}\omega \right)$", fontsize =12
         )
         fig.tight_layout()
         # cbar.set_label(r"$\log S(k_z, \omega)$")
@@ -263,7 +263,7 @@
 
         ax.set_title(rf"${key_tex} ~~~ k_z =$" + kz_leg_txt, fontsize =12)
         cbar.set_label(
-            r"$\log \left( E_{\rm pot}(k_h, \omega) / \int E_{\rm pot}(k_h, \omega) \mathrm{d}\omega \right)$", fontsize =12
+            r"$\log \left( E_{\rm A}(k_h, \omega) / \int E_{\rm A}(k_h, \omega) \mathrm{d}\omega \right)$", fontsize =12
         )
         fig.tight_layout()
         # cbar.set_label(r"$\log S(k_h, \omega)$")
diff --git a/2022strat_turb_polo/py/save_spectra_spatiotemporal_kh_kz.py b/2022strat_turb_polo/py/save_spectra_spatiotemporal_kh_kz.py
--- a/2022strat_turb_polo/py/save_spectra_spatiotemporal_kh_kz.py
+++ b/2022strat_turb_polo/py/save_spectra_spatiotemporal_kh_kz.py
@@ -92,13 +92,13 @@
 
     if key == "potential":
         spectrum = f["spectrum_A"][:]  # Potential energy
-        key_tex = "E_{\rm pot}"
+        key_tex = "E_{\rm A}"
     elif key == "kinetic":
         spectrum = f["spectrum_K"][:]
         key_tex = "E_{\rm K}"
     elif key == "total":
         spectrum = f["spectrum_A"][:] + f["spectrum_K"][:]
-        key_tex = "E_{\rm pot} + E_{\rm K}"
+        key_tex = "E_{\rm A} + E_{\rm K}"
     else:
         print(f"Do not know key {key} \n")
         #continue
@@ -139,7 +139,7 @@
     #rf"$proj={proj}  ~~~  N={N}  ~~~  {key_tex} ~~~ k_h / \delta k_h={ikh}$"
     ax.set_title(rf"{projtxt}         $F_h={Fh:.3f}$        " + r"$\mathcal{R} =$" + rf"{R2:.2f}" +  "\n" + rf"$\omega=$" + rf"${ratio_omega_slice_N}$" + r" $N$" , fontsize =12)
     cbar.set_label(
-        r"$\log E_{\rm pot}(k_h, k_z, \omega)$", fontsize =12
+        r"$\log E_{\rm A}(k_h, k_z, \omega)$", fontsize =12
     )
     fig.tight_layout()
     # cbar.set_label(r"$\log S(k_z, \omega)$")