\documentclass[a4paper]{article}

                                % Packages

\usepackage[utf8]{inputenc}
\usepackage{ifpdf}

\ifpdf
\else
\usepackage[active]{srcltx}
\fi

\usepackage{amsmath}

\usepackage[T1]{fontenc}
\usepackage{lmodern}

\usepackage{geometry}

\usepackage[np]{numprint}

\usepackage{hyperref}
\ifpdf
\hypersetup{pdftitle={Tutorial Tracers}, pdfauthor={LMDZ team},
  hypertexnames=false, pdfstartview=FitBH}
\fi

\author{LMDZ team}
\title{LDMZ tutorial: tracers}
\date{November 30th, 2017}


\geometry{top=2cm,bottom=2cm,left=2.5cm,right=2.5cm}


                                %---------------

\begin{document}

\maketitle

This tutorial focuses on using and adding tracers in LMDZ.

This document can be downloaded as a pdf file:
\begin{verbatim}
wget http://www.lmd.jussieu.fr/~lmdz/pub/Training/Tutorials/Tutorial_Tracers.pdf
\end{verbatim}
which should ease any copy/paste of command lines to issue.

\section{Prerequisits}

You should be familiar with setting up simulations, as described in tutorials
\#1 and \#2.

\section{Experimenting with tracers}

\subsubsection{Radon and lead}

In the file \verb+traceur.def+, you can see that two trace species,
named RN, for radon, and PB, for lead, are already defined in addition
to vapor and condensed water. Visualize these two tracers in the file
\verb+histhf.nc+. The NetCDF variables RN and PB are numbers of atoms
per unit mass, in kg$^{-1}$. Usually, activities are preferred for
comparison with observations. The activities per unit mass, in Bq
kg$^{-1}$, are:
\begin{align*}
  & A_\mathrm{Rn} = \frac{\mathrm{RN}}{\np{4.765E5}} \\
  & A_\mathrm{Pb} = \frac{\mathrm{PB}}{\np{1.028E9}}
\end{align*}
(The numerical values that appear in these fractions are the lifetimes
of radon and lead, in s.) You can take a look at
\href{http://lmdz.lmd.jussieu.fr/communication/lmdzinfo/lmdzinfo7.pdf/at_download/file}{LMDZ
  info number 7}, figure 19, or the plot in the
\href{http://lmdz.lmd.jussieu.fr/le-projet-lmdz/formation/2014/contents/general-introduction}{General
  introduction of the LMDZ model} (§ 4 ``Operating modes'') to check
that you have sensible results. A recent work on the modeling of
tracers with LMDZ is Pilon et al. (2015, QJRMS). If you want to see
activities per unit volume, you should divide by the mass density,
using NetCDF variables pres, temp and ovap in the file
\verb+histhf.nc+. (The Ferret color palette in figure 19 of LMDZ info
number 7 is saz2.)

\subsubsection{Inserting new tracers}

We will now explain how to modify the code in order to add idealized
tracers. As an example, we will add two tracers, and call them
\verb+NewTr1+ and \verb+NewTr2+. We will define a domain of the
horizontal grid in which the two tracers will initially have the same
constant value. The first tracer will be transported by the boundary
layer and convective sub-grid-scale motion, in addition to large-scale
advection. The second tracer will only be transported by large-scale
advection.

We will have to modify the Fortran program but let us first consider
the run-time parameters that we have to modify. At run-time, we decide
to include the tracers \verb+NewTr1+ and \verb+NewTr2+ in the
simulation by changing the file \verb+traceur.def+. Change the number
of tracers at the first line of \verb+traceur.def+ and append one line
for each tracer:
\begin{verbatim}
6
14 14 H2Ov
10 10 H2Ol
10 10 RN
10 10 PB
10 10 NewTr1
10 10 NewTr2
\end{verbatim}

Now let us turn to the Fortran program. The only file we need to
modify 
is

\begin{verbatim}
LMDZtesting/modipsl/modeles/LMDZ5/libf/phylmd/traclmdz_mod.F90
\end{verbatim}


Here are the changes you should make in that file:
\begin{itemize}
\item Declare two new module variables, \verb+id_NewTr1+ and
  \verb+id_NewTr2+, with type integer. These are the identifying
  numbers of the tracers in the program. You can take a previous
  declaration (\verb+id_pcsat+, \verb+id_pcocsat+ \dots) as a
  template. (Do not forget the OpenMP directives.)
\item All remaining changes will be made in the procedure
  \verb+traclmdz_init+, which is inside the module
  \verb+traclmdz_mod+. In order to define \verb+id_NewTr1+ and
  \verb+id_NewTr2+, the program will scan the file \verb+traceur.def+,
  looking for NewTr1 and NewTr2. So you should initialize
  \verb+id_NewTr1+ and \verb+id_NewTr2+ to 0 before the loop beginning
  at line 175. There is a comment just above, saying ``Recherche des
  traceurs connus'', which means ``looking for known tracers''.
\item In the body of this loop, set \verb+id_NewTr1+ to the value of
  index \verb+it+ if tname(iiq) equals NewTr1. You can add the test
  near line 265, for instance, after the test for pcq0. You can take
  another tracer as a template. Do the same for \verb+id_NewTr2+.
\item For NewTr2, just after setting \verb+id_NewTr2+, deactivate
  convective and boundary layer transport by setting
  \verb+conv_flg(it)+ and \verb+pbl_flg(it)+ to 0.
\item Finally, we will initialize the tracers. There is a loop on
  tracers, beginning at line 294, which tests whether the initial
  tracer field, read from the file \verb+start.nc+, is zero
  everywhere. At this point, the tracer field would also be zero if it
  was not found in \verb+start.nc+. There is a comment just above the
  test which says ``Initalize tracer that was not found in restart
  file'' (the typo ``Initalize'' is in the code!).
In the body of the test, for our two tracers, change the
  value at the surface in some horizontal domain. The value of tracers
  is in the variable \verb+tr_seri+. The first dimension of
  \verb+tr_seri+ is for the horizontal position, the second dimension
  is for the vertical level and the third dimension identifies the
  tracer. The index of the surface in the vertical dimension is 1. Use
  variables \verb+xlat+ (latitudes) and \verb+xlon+ (longitudes) to
  choose the horizontal domain.
\end{itemize}

In summary, after making those changes, svn diff should give you
something like this:
\begin{verbatim}
$ svn diff traclmdz_mod.F90
Index: traclmdz_mod.F90
===================================================================
--- traclmdz_mod.F90	(revision 2304)
+++ traclmdz_mod.F90	(working copy)
@@ -58,6 +58,8 @@
   LOGICAL,SAVE :: rnpb=.FALSE. ! Presence du couple Rn222, Pb210
 !$OMP THREADPRIVATE(rnpb)
 
+  INTEGER, SAVE:: id_newtr1, id_newtr2
+  !$OMP THREADPRIVATE(id_newtr1, id_newtr2)
 
 CONTAINS
 
@@ -172,6 +174,8 @@
 ! --------------------------------------------
     id_rn=0; id_pb=0; id_aga=0; id_be=0; id_o3=0
     id_pcsat=0; id_pcocsat=0; id_pcq=0; id_pcs0=0; id_pcos0=0; id_pcq0=0
+    id_newtr1 = 0
+    id_newtr2 = 0
     DO it=1,nbtr
        iiq=niadv(it+2)
        IF ( tname(iiq) == "RN" ) THEN
@@ -262,6 +266,12 @@
        ELSE IF ( tname(iiq) == "pcq0" .OR. tname(iiq) == "Pcq0" ) THEN
           id_pcq0=it
           conv_flg(it)=0 ! No transport by convection for this tracer
+       else if (tname(iiq) == "NewTr1") then
+          id_newtr1 = it
+       else if (tname(iiq) == "NewTr2") then
+          id_newtr2 = it
+          conv_flg(it) = 0
+          pbl_flg(it) = 0
        ELSE
           WRITE(lunout,*) 'This is an unknown tracer in LMDZ : ', trim(tname(iiq))
        END IF
@@ -325,6 +335,9 @@
                    tr_seri(i,:,it) = 100.
                 END IF
              END DO
+          else if (it == id_newtr1 .or. it == id_newtr2) then
+             where (xlat >= 40. .and. xlat <= 45. .and. xlon >=0. &
+                  .and. xlon <= 5.) tr_seri(:, 1, it) = 1.
           ELSE
              ! No specific initialization exist for this tracer
              tr_seri(:,:,it) = 0.
\end{verbatim}

Re-compile the program:
\begin{verbatim}
cd some_path/LMDZtesting/modipsl/modeles/LMDZ5/TUTO/SIMU1
rm gcm.e
cd ../..
\end{verbatim}
If you ran \verb+init.sh+ with veget=0:
\begin{verbatim}
./makelmdz -d 48x36x39 gcm
\end{verbatim}
If you set veget=1 in \verb+init.sh+:
\begin{verbatim}
./makelmdz -d 48x36x39 -v orchidee2.0 -cpp ORCHIDEE_NOZ0H gcm
\end{verbatim}
Rename restart files:
\begin{verbatim}
cd TUTO/SIMU1
mv restart.nc start.nc
mv restartphy.nc startphy.nc
mv sechiba_rest_out.nc sechiba_rest_in.nc
\end{verbatim}
Run the model:
\begin{verbatim}
../../gcm.e > listing
\end{verbatim}
(could take about 15 mn). Visualize the two new tracers in
\verb+histhf.nc+ and the difference between them.

\end{document}
