/[MITgcm]/manual/s_autodiff/text/doc_ad_2.tex
ViewVC logotype

Diff of /manual/s_autodiff/text/doc_ad_2.tex

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.6 by adcroft, Thu Oct 11 19:37:39 2001 UTC revision 1.18 by heimbach, Mon Aug 1 22:31:36 2005 UTC
# Line 1  Line 1 
1  % $Header$  % $Header$
2  % $Name$  % $Name$
3    
4    Author: Patrick Heimbach
5    
6  {\sf Automatic differentiation} (AD), also referred to as algorithmic  {\sf Automatic differentiation} (AD), also referred to as algorithmic
7  (or, more loosely, computational) differentiation, involves  (or, more loosely, computational) differentiation, involves
8  automatically deriving code to calculate  automatically deriving code to calculate
# Line 21  The MITGCM has been adapted for use with Line 23  The MITGCM has been adapted for use with
23  Tangent linear and Adjoint Model Compiler (TAMC) and its successor TAF  Tangent linear and Adjoint Model Compiler (TAMC) and its successor TAF
24  (Transformation of Algorithms in Fortran), developed  (Transformation of Algorithms in Fortran), developed
25  by Ralf Giering (\cite{gie-kam:98}, \cite{gie:99,gie:00}).  by Ralf Giering (\cite{gie-kam:98}, \cite{gie:99,gie:00}).
26  The first application of the adjoint of the MITGCM for senistivity  The first application of the adjoint of the MITGCM for sensitivity
27  studies has been published by \cite{maro-eta:99}.  studies has been published by \cite{maro-eta:99}.
28  \cite{sta-eta:97,sta-eta:01} use the MITGCM and its adjoint  \cite{sta-eta:97,sta-eta:01} use the MITGCM and its adjoint
29  for ocean state estimation studies.  for ocean state estimation studies.
# Line 42  Jacobian matrices of the forward code's Line 44  Jacobian matrices of the forward code's
44  %**********************************************************************  %**********************************************************************
45  \section{Some basic algebra}  \section{Some basic algebra}
46  \label{sec_ad_algebra}  \label{sec_ad_algebra}
47    \begin{rawhtml}
48    <!-- CMIREDIR:sec_ad_algebra: -->
49    \end{rawhtml}
50  %**********************************************************************  %**********************************************************************
51    
52  Let $ \cal{M} $ be a general nonlinear, model, i.e. a  Let $ \cal{M} $ be a general nonlinear, model, i.e. a
# Line 52  $\vec{u}=(u_1,\ldots,u_m)$ Line 57  $\vec{u}=(u_1,\ldots,u_m)$
57  such as forcing functions) to the $n$-dimensional space  such as forcing functions) to the $n$-dimensional space
58  $V \subset I\!\!R^n$ of  $V \subset I\!\!R^n$ of
59  model output variable $\vec{v}=(v_1,\ldots,v_n)$  model output variable $\vec{v}=(v_1,\ldots,v_n)$
60  (model state, model diagnostcs, objective function, ...)  (model state, model diagnostics, objective function, ...)
61  under consideration,  under consideration,
62  %  %
63  \begin{equation}  \begin{equation}
# Line 220  model integration, Line 225  model integration,
225  starting at step 0 and moving up to step $\Lambda$, with intermediate  starting at step 0 and moving up to step $\Lambda$, with intermediate
226  ${\cal M}_{\lambda} (\vec{u}) = \vec{v}^{(\lambda+1)}$ and final  ${\cal M}_{\lambda} (\vec{u}) = \vec{v}^{(\lambda+1)}$ and final
227  ${\cal M}_{\Lambda} (\vec{u}) = \vec{v}^{(\Lambda+1)} = \vec{v}$.  ${\cal M}_{\Lambda} (\vec{u}) = \vec{v}^{(\Lambda+1)} = \vec{v}$.
228  Let ${\cal J}$ be a cost funciton which explicitly depends on the  Let ${\cal J}$ be a cost function which explicitly depends on the
229  final state $\vec{v}$ only  final state $\vec{v}$ only
230  (this restriction is for clarity reasons only).  (this restriction is for clarity reasons only).
231  %  %
# Line 301  We note in passing that that the $\delta Line 306  We note in passing that that the $\delta
306  are the Lagrange multipliers of the model equations which determine  are the Lagrange multipliers of the model equations which determine
307  $ \vec{v}^{(\lambda)}$.  $ \vec{v}^{(\lambda)}$.
308    
309  In coponents, eq. (\ref{adjoint}) reads as follows.  In components, eq. (\ref{adjoint}) reads as follows.
310  Let  Let
311  \[  \[
312  \begin{array}{rclcrcl}  \begin{array}{rclcrcl}
# Line 322  Let Line 327  Let
327  \end{array}  \end{array}
328  \]  \]
329  denote the perturbations in $\vec{u}$ and $\vec{v}$, respectively,  denote the perturbations in $\vec{u}$ and $\vec{v}$, respectively,
330  and their adjoint varaiables;  and their adjoint variables;
331  further  further
332  \[  \[
333  M \, = \, \left(  M \, = \, \left(
# Line 468  variables $u$ Line 473  variables $u$
473  {\it all} intermediate states $ \vec{v}^{(\lambda)} $) are sought.  {\it all} intermediate states $ \vec{v}^{(\lambda)} $) are sought.
474  In order to be able to solve for each component of the gradient  In order to be able to solve for each component of the gradient
475  $ \partial {\cal J} / \partial u_{i} $ in (\ref{forward})  $ \partial {\cal J} / \partial u_{i} $ in (\ref{forward})
476  a forward calulation has to be performed for each component seperately,  a forward calculation has to be performed for each component separately,
477  i.e. $ \delta \vec{u} = \delta u_{i} {\vec{e}_{i}} $  i.e. $ \delta \vec{u} = \delta u_{i} {\vec{e}_{i}} $
478  for  the $i$-th forward calculation.  for  the $i$-th forward calculation.
479  Then, (\ref{forward}) represents the  Then, (\ref{forward}) represents the
# Line 487  M^T \left( \nabla_v {\cal J}^T \left(\de Line 492  M^T \left( \nabla_v {\cal J}^T \left(\de
492  \nabla_u {\cal J}^T \cdot \delta \vec{J}  \nabla_u {\cal J}^T \cdot \delta \vec{J}
493  \]  \]
494  where now $ \delta \vec{J} \in I\!\!R^l $ is a vector of  where now $ \delta \vec{J} \in I\!\!R^l $ is a vector of
495  dimenison $ l $.  dimension $ l $.
496  In this case $ l $ reverse simulations have to be performed  In this case $ l $ reverse simulations have to be performed
497  for each $ \delta J_{k}, \,\, k = 1, \ldots, l $.  for each $ \delta J_{k}, \,\, k = 1, \ldots, l $.
498  Then, the reverse mode is more efficient as long as  Then, the reverse mode is more efficient as long as
499  $ l < n $, otherwise the forward mode is preferable.  $ l < n $, otherwise the forward mode is preferable.
500  Stricly, the reverse mode is called adjoint mode only for  Strictly, the reverse mode is called adjoint mode only for
501  $ l = 1 $.  $ l = 1 $.
502    
503  A detailed analysis of the underlying numerical operations  A detailed analysis of the underlying numerical operations
# Line 557  Because of the local character of the de Line 562  Because of the local character of the de
562  (a derivative is defined w.r.t. a point along the trajectory),  (a derivative is defined w.r.t. a point along the trajectory),
563  the intermediate results of the model trajectory  the intermediate results of the model trajectory
564  $\vec{v}^{(\lambda+1)}={\cal M}_{\lambda}(v^{(\lambda)})$  $\vec{v}^{(\lambda+1)}={\cal M}_{\lambda}(v^{(\lambda)})$
565  are needed to evaluate the intermediate Jacobian  may be required to evaluate the intermediate Jacobian
566  $M_{\lambda}|_{\vec{v}^{(\lambda)}} \, \delta \vec{v}^{(\lambda)} $.  $M_{\lambda}|_{\vec{v}^{(\lambda)}} \, \delta \vec{v}^{(\lambda)} $.
567    This is the case e.g. for nonlinear expressions
568    (momentum advection, nonlinear equation of state), state-dependent
569    conditional statements (parameterization schemes).
570  In the forward mode, the intermediate results are required  In the forward mode, the intermediate results are required
571  in the same order as computed by the full forward model ${\cal M}$,  in the same order as computed by the full forward model ${\cal M}$,
572  but in the reverse mode they are required in the reverse order.  but in the reverse mode they are required in the reverse order.
# Line 569  point of evaluation has to be recomputed Line 577  point of evaluation has to be recomputed
577    
578  A method to balance the amount of recomputations vs.  A method to balance the amount of recomputations vs.
579  storage requirements is called {\sf checkpointing}  storage requirements is called {\sf checkpointing}
580  (e.g. \cite{res-eta:98}).  (e.g. \cite{gri:92}, \cite{res-eta:98}).
581  It is depicted in \ref{fig:3levelcheck} for a 3-level checkpointing  It is depicted in \ref{fig:3levelcheck} for a 3-level checkpointing
582  [as an example, we give explicit numbers for a 3-day  [as an example, we give explicit numbers for a 3-day
583  integration with a 1-hourly timestep in square brackets].  integration with a 1-hourly timestep in square brackets].
# Line 580  In a first step, the model trajectory is Line 588  In a first step, the model trajectory is
588  $ {n}^{lev3} $ subsections [$ {n}^{lev3} $=3 1-day intervals],  $ {n}^{lev3} $ subsections [$ {n}^{lev3} $=3 1-day intervals],
589  with the label $lev3$ for this outermost loop.  with the label $lev3$ for this outermost loop.
590  The model is then integrated along the full trajectory,  The model is then integrated along the full trajectory,
591  and the model state stored only at every $ k_{i}^{lev3} $-th timestep  and the model state stored to disk only at every $ k_{i}^{lev3} $-th timestep
592  [i.e. 3 times, at  [i.e. 3 times, at
593  $ i = 0,1,2 $ corresponding to $ k_{i}^{lev3} = 0, 24, 48 $].  $ i = 0,1,2 $ corresponding to $ k_{i}^{lev3} = 0, 24, 48 $].
594    In addition, the cost function is computed, if needed.
595  %  %
596  \item [$lev2$]  \item [$lev2$]
597  In a second step each subsection itself is divided into  In a second step each subsection itself is divided into
598  $ {n}^{lev2} $ sub-subsections  $ {n}^{lev2} $ subsections
599  [$ {n}^{lev2} $=4 6-hour intervals per subsection].  [$ {n}^{lev2} $=4 6-hour intervals per subsection].
600  The model picks up at the last outermost dumped state  The model picks up at the last outermost dumped state
601  $ v_{k_{n}^{lev3}} $ and is integrated forward in time along  $ v_{k_{n}^{lev3}} $ and is integrated forward in time along
602  the last subsection, with the label $lev2$ for this    the last subsection, with the label $lev2$ for this  
603  intermediate loop.  intermediate loop.
604  The model state is now stored at every $ k_{i}^{lev2} $-th  The model state is now stored to disk at every $ k_{i}^{lev2} $-th
605  timestep  timestep
606  [i.e. 4 times, at  [i.e. 4 times, at
607  $ i = 0,1,2,3 $ corresponding to $ k_{i}^{lev2} = 48, 54, 60, 66 $].  $ i = 0,1,2,3 $ corresponding to $ k_{i}^{lev2} = 48, 54, 60, 66 $].
# Line 600  $ i = 0,1,2,3 $ corresponding to $ k_{i} Line 609  $ i = 0,1,2,3 $ corresponding to $ k_{i}
609  \item [$lev1$]  \item [$lev1$]
610  Finally, the model picks up at the last intermediate dump state  Finally, the model picks up at the last intermediate dump state
611  $ v_{k_{n}^{lev2}} $ and is integrated forward in time along  $ v_{k_{n}^{lev2}} $ and is integrated forward in time along
612  the last sub-subsection, with the label $lev1$ for this    the last subsection, with the label $lev1$ for this  
613  intermediate loop.  intermediate loop.
614  Within this sub-subsection only, the model state is stored  Within this sub-subsection only, parts of the model state is stored
615  at every timestep  to memory at every timestep
616  [i.e. every hour $ i=0,...,5$ corresponding to  [i.e. every hour $ i=0,...,5$ corresponding to
617  $ k_{i}^{lev1} = 66, 67, \ldots, 71 $].  $ k_{i}^{lev1} = 66, 67, \ldots, 71 $].
618  Thus, the  final state $ v_n = v_{k_{n}^{lev1}} $ is reached  The  final state $ v_n = v_{k_{n}^{lev1}} $ is reached
619  and the model state of all peceeding timesteps along the last  and the model state of all preceding timesteps along the last
620  sub-subsections are available, enabling integration backwards  innermost subsection are available, enabling integration backwards
621  in time along the last sub-subsection.  in time along the last subsection.
622  Thus, the adjoint can be computed along this last  The adjoint can thus be computed along this last
623  sub-subsection $k_{n}^{lev2}$.  subsection $k_{n}^{lev2}$.
624  %  %
625  \end{itemize}  \end{itemize}
626  %  %
627  This procedure is repeated consecutively for each previous  This procedure is repeated consecutively for each previous
628  sub-subsection $k_{n-1}^{lev2}, \ldots, k_{1}^{lev2} $  subsection $k_{n-1}^{lev2}, \ldots, k_{1}^{lev2} $
629  carrying the adjoint computation to the initial time  carrying the adjoint computation to the initial time
630  of the subsection $k_{n}^{lev3}$.  of the subsection $k_{n}^{lev3}$.
631  Then, the procedure is repeated for the previous subsection  Then, the procedure is repeated for the previous subsection
# Line 627  $k_{1}^{lev3}$. Line 636  $k_{1}^{lev3}$.
636  For the full model trajectory of  For the full model trajectory of
637  $ n^{lev3} \cdot n^{lev2} \cdot n^{lev1} $ timesteps  $ n^{lev3} \cdot n^{lev2} \cdot n^{lev1} $ timesteps
638  the required storing of the model state was significantly reduced to  the required storing of the model state was significantly reduced to
639  $ n^{lev1} + n^{lev2} + n^{lev3} $  $ n^{lev2} + n^{lev3} $ to disk and roughly $ n^{lev1} $ to memory
640  [i.e. for the 3-day integration with a total oof 72 timesteps  [i.e. for the 3-day integration with a total oof 72 timesteps
641  the model state was stored 13 times].  the model state was stored 7 times to disk and roughly 6 times
642    to memory].
643  This saving in memory comes at a cost of a required  This saving in memory comes at a cost of a required
644  3 full forward integrations of the model (one for each  3 full forward integrations of the model (one for each
645  checkpointing level).  checkpointing level).
646  The balance of storage vs. recomputation certainly depends  The optimal balance of storage vs. recomputation certainly depends
647  on the computing resources available.  on the computing resources available and may be adjusted by
648    adjusting the partitioning among the
649    $ n^{lev3}, \,\, n^{lev2}, \,\, n^{lev1} $.
650    
651  \begin{figure}[t!]  \begin{figure}[t!]
652  \begin{center}  \begin{center}
# Line 664  Schematic view of intermediate dump and Line 676  Schematic view of intermediate dump and
676  % \subsection{Error covariance estimate and Hessian matrix}  % \subsection{Error covariance estimate and Hessian matrix}
677  % \label{sec_hessian}  % \label{sec_hessian}
678    
679  \newpage  \newpage
680    
681  %**********************************************************************  %**********************************************************************
682  \section{AD-specific setup by example: sensitivity of carbon sequestration}  \section{TLM and ADM generation in general}
683  \label{sec_ad_setup_ex}  \label{sec_ad_setup_gen}
684    \begin{rawhtml}
685    <!-- CMIREDIR:sec_ad_setup_gen: -->
686    \end{rawhtml}
687  %**********************************************************************  %**********************************************************************
688    
689  The MITGCM has been adapted to enable AD using TAMC or TAF.  In this section we describe in a general fashion
690  The present description, therefore, is specific to the  the parts of the code that are relevant for automatic
691  use of TAMC or TAF as AD tool.  differentiation using the software tool TAF.
692  The following sections describe the steps which are necessary to  
693  generate a tangent linear or adjoint model of the MITGCM.  \input{part5/doc_ad_the_model}
694  We take as an example the sensitivity of carbon sequestration  
695  in the ocean.  The basic flow is depicted in \ref{fig:adthemodel}.
696  The AD-relevant hooks in the code are sketched in  If CPP option {\tt ALLOW\_AUTODIFF\_TAMC} is defined, the driver routine
697  \ref{fig:adthemodel}, \ref{fig:adthemain}.  {\it the\_model\_main}, instead of calling {\it the\_main\_loop},
698    invokes the adjoint of this routine, {\it adthe\_main\_loop},
699  \subsection{Overview of the experiment}  which is the toplevel routine in terms of automatic differentiation.
700    The routine {\it adthe\_main\_loop} has been generated by TAF.
701  We describe an adjoint sensitivity analysis of outgassing from  It contains both the forward integration of the full model, the
702  the ocean into the atmosphere of a carbon-like tracer injected  cost function calculation,
703  into the ocean interior (see \cite{hil-eta:01}).  any additional storing that is required for efficient checkpointing,
704    and the reverse integration of the adjoint model.
705  \subsubsection{Passive tracer equation}  
706    [DESCRIBE IN A SEPARATE SECTION THE WORKING OF THE TLM]
 For this work the MITGCM was augmented with a thermodynamically  
 inactive tracer, $C$. Tracer residing in the ocean  
 model surface layer is outgassed according to a relaxation time scale,  
 $\mu$. Within the ocean interior, the tracer is passively advected  
 by the ocean model currents. The full equation for the time evolution  
 %  
 \begin{equation}  
 \label{carbon_ddt}  
 \frac{\partial C}{\partial t} \, = \,  
 -U\cdot \nabla C \, - \, \mu C \, + \, \Gamma(C) \,+ \, S  
 \end{equation}  
 %  
 also includes a source term $S$. This term  
 represents interior sources of $C$ such as would arise due to  
 direct injection.  
 The velocity term, $U$, is the sum of the  
 model Eulerian circulation and an eddy-induced velocity, the latter  
 parameterized according to Gent/McWilliams  
 (\cite{gen-mcw:90, gen-eta:95}).  
 The convection function, $\Gamma$, mixes $C$ vertically wherever the  
 fluid is locally statically unstable.  
   
 The outgassing time scale, $\mu$, in eqn. (\ref{carbon_ddt})  
 is set so that \( 1/\mu \sim 1 \ \mathrm{year} \) for the surface  
 ocean and $\mu=0$ elsewhere. With this value, eqn. (\ref{carbon_ddt})  
 is valid as a prognostic equation for small perturbations in oceanic  
 carbon concentrations. This configuration provides a  
 powerful tool for examining the impact of large-scale ocean circulation  
 on $ CO_2 $ outgassing due to interior injections.  
 As source we choose a constant in time injection of  
 $ S = 1 \,\, {\rm mol / s}$.  
   
 \subsubsection{Model configuration}  
   
 The model configuration employed has a constant  
 $4^\circ \times 4^\circ$ resolution horizontal grid and realistic  
 geography and bathymetry. Twenty vertical layers are used with  
 vertical spacing ranging  
 from 50 m near the surface to 815 m at depth.  
 Driven to steady-state by climatalogical wind-stress, heat and  
 fresh-water forcing the model reproduces well known large-scale  
 features of the ocean general circulation.  
   
 \subsubsection{Outgassing cost function}  
   
 To quantify and understand outgassing due to injections of $C$  
 in eqn. (\ref{carbon_ddt}),  
 we define a cost function $ {\cal J} $ that measures the total amount of  
 tracer outgassed at each timestep:  
 %  
 \begin{equation}  
 \label{cost_tracer}  
 {\cal J}(t=T)=\int_{t=0}^{t=T}\int_{A} \mu C \, dA \, dt  
 \end{equation}  
 %  
 Equation(\ref{cost_tracer}) integrates the outgassing term, $\mu C$,  
 from (\ref{carbon_ddt})  
 over the entire ocean surface area, $A$, and accumulates it  
 up to time $T$.  
 Physically, ${\cal J}$ can be thought of as representing the amount of  
 $CO_2$ that our model predicts would be outgassed following an  
 injection at rate $S$.  
 The sensitivity of ${\cal J}$ to the spatial location of $S$,  
 $\frac{\partial {\cal J}}{\partial S}$,  
 can be used to identify regions from which circulation  
 would cause $CO_2$ to rapidly outgas following injection  
 and regions in which $CO_2$ injections would remain effectively  
 sequesterd within the ocean.  
   
 \subsection{Code configuration}  
   
 The model configuration for this experiment resides under the  
 directory {\it verification/carbon/}.  
 The code customisation routines are in {\it verification/carbon/code/}:  
 %  
 \begin{itemize}  
 %  
 \item {\it .genmakerc}  
 %  
 \item {\it COST\_CPPOPTIONS.h}  
 %  
 \item {\it CPP\_EEOPTIONS.h}  
 %  
 \item {\it CPP\_OPTIONS.h}  
 %  
 \item {\it CTRL\_OPTIONS.h}  
 %  
 \item {\it ECCO\_OPTIONS.h}  
 %  
 \item {\it SIZE.h}  
 %  
 \item {\it adcommon.h}  
 %  
 \item {\it tamc.h}  
 %  
 \end{itemize}  
 %  
 The runtime flag and parameters settings are contained in  
 {\it verification/carbon/input/},  
 together with the forcing fields and and restart files:  
 %  
 \begin{itemize}  
 %  
 \item {\it data}  
 %  
 \item {\it data.cost}  
 %  
 \item {\it data.ctrl}  
 %  
 \item {\it data.gmredi}  
 %  
 \item {\it data.grdchk}  
 %  
 \item {\it data.optim}  
 %  
 \item {\it data.pkg}  
 %  
 \item {\it eedata}  
 %  
 \item {\it topog.bin}  
 %  
 \item {\it windx.bin, windy.bin}  
 %  
 \item {\it salt.bin, theta.bin}  
 %  
 \item {\it SSS.bin, SST.bin}  
 %  
 \item {\it pickup*}  
 %  
 \end{itemize}  
 %  
 Finally, the file to generate the adjoint code resides in  
 $ adjoint/ $:  
 %  
 \begin{itemize}  
 %  
 \item {\it makefile}  
 %  
 \end{itemize}  
 %  
707    
708  Below we describe the customisations of this files which are  In Fig. \ref{fig:adthemodel}
709  specific to this experiment.  the structure of {\it adthe\_main\_loop} has been strongly
710    simplified to focus on the essentials; in particular, no checkpointing
711    procedures are shown here.
712    Prior to the call of {\it adthe\_main\_loop}, the routine
713    {\it ctrl\_unpack} is invoked to unpack the control vector
714    or initialise the control variables.
715    Following the call of {\it adthe\_main\_loop},
716    the routine {\it ctrl\_pack}
717    is invoked to pack the control vector
718    (cf. Section \ref{section_ctrl}).
719    If gradient checks are to be performed, the option
720    {\tt ALLOW\_GRADIENT\_CHECK} is defined. In this case
721    the driver routine {\it grdchk\_main} is called after
722    the gradient has been computed via the adjoint
723    (cf. Section \ref{section_grdchk}).
724    
725  \subsubsection{File {\it .genmakerc}}  %------------------------------------------------------------------
 This file overwrites default settings of {\it genmake}.  
 In the present example it is used to switch on the following  
 packages which are related to automatic differentiation  
 and are disabled by default: \\  
 \hspace*{4ex} {\tt set ENABLE=( autodiff cost ctrl ecco gmredi grdchk kpp )}  \\  
 Other packages which are not needed are switched off: \\  
 \hspace*{4ex} {\tt set DISABLE=( aim obcs zonal\_filt shap\_filt cal exf )}  
   
 \subsubsection{File {\it COST\_CPPOPTIONS.h,  CTRL\_OPTIONS.h}}  
   
 These files used to contain package-specific CPP-options  
 (see Section \ref{???}).  
 For technical reasons those options have been grouped together  
 in the file {\it ECCO\_OPTIONS.h}.  
 To retain the modularity, the files have been kept and contain  
 the standard include of the {\it CPP\_OPTIONS.h} file.  
   
 \subsubsection{File {\it CPP\_EEOPTIONS.h}}  
   
 This file contains 'wrapper'-specific CPP options.  
 It only needs to be changed if the code is to be run  
 in a parallel environment (see Section \ref{???}).  
   
 \subsubsection{File {\it CPP\_OPTIONS.h}}  
   
 This file contains model-specific CPP options  
 (see Section \ref{???}).  
 Most options are related to the forward model setup.  
 They are identical to the global steady circulation setup of  
 {\it verification/exp2/}.  
 The three options specific to this experiment are \\  
 \hspace*{4ex} {\tt \#define ALLOW\_PASSIVE\_TRACER} \\  
 This flag enables the code to carry through the  
 advection/diffusion of a passive tracer along the  
 model integration. \\  
 \hspace*{4ex} {\tt \#define ALLOW\_MIT\_ADJOINT\_RUN} \\  
 This flag enables the inclusion of some AD-related fields  
 concerning initialisation, link between control variables  
 and forward model variables, and the call to the top-level  
 forward/adjoint subroutine {\it adthe\_main\_loop}  
 instead of {\it the\_main\_loop}. \\  
 \hspace*{4ex} {\tt \#define ALLOW\_GRADIENT\_CHECK} \\  
 This flag enables the gradient check package.  
 After computing the unperturbed cost function and its gradient,  
 a series of computations are performed for which \\  
 $\bullet$ an element of the control vector is perturbed \\  
 $\bullet$ the cost function w.r.t. the perturbed element is  
 computed \\  
 $\bullet$ the difference between the perturbed and unperturbed  
 cost function is computed to compute the finite difference gradient \\  
 $\bullet$ the finite difference gradient is compared with the  
 adjoint-generated gradient.  
 The gradient check package is further described in Section ???.  
726    
727  \subsubsection{File {\it ECCO\_OPTIONS.h}}  \subsection{General setup
728    \label{section_ad_setup}}
729    
730  The CPP options of several AD-related packages are grouped  In order to configure AD-related setups the following packages need
731  in this file:  to be enabled:
732    {\it
733    \begin{table}[h!]
734    \begin{tabular}{l}
735    autodiff \\
736    ctrl \\
737    cost \\
738    grdchk \\
739    \end{tabular}
740    \end{table}
741    }
742    The packages are enabled by adding them to your experiment-specific
743    configuration file
744    {\it packages.conf} (see Section ???).
745    
746    The following AD-specific CPP option files need to be customized:
747  %  %
748  \begin{itemize}  \begin{itemize}
749  %  %
750  \item  \item {\it ECCO\_CPPOPTIONS.h} \\
751  Adjoint support package: {\it pkg/autodiff/} \\  This header file collects CPP options for the packages
752  This package contains hand-written adjoint code such as  {\it autodiff, cost, ctrl} as well as AD-unrelated options for
753  active file handling, flow directives for files which must not  the external forcing package {\it exf}.
754  be differentiated, and TAMC-specific header files. \\  \footnote{NOTE: These options are not set in their package-specific
755  \hspace*{4ex} {\tt \#define ALLOW\_AUTODIFF\_TAMC} \\  headers such as {\it COST\_CPPOPTIONS.h}, but are instead collected
756  defines TAMC-related features in the code. \\  in the single header file {\it ECCO\_CPPOPTIONS.h}.
757  \hspace*{4ex} {\tt \#define ALLOW\_TAMC\_CHECKPOINTING} \\  The package-specific header files serve as simple
758  enables the checkpointing feature of TAMC  placeholders at this point.}
759  (see Section \ref{???}).  %
760  In the present example a 3-level checkpointing is implemented.  \item {\it tamc.h} \\
761  The code contains the relevant store directives, common block  This header configures the splitting of the time stepping loop
762  and tape initialisations, storing key computation,  w.r.t. the 3-level checkpointing (see section ???).
763  and loop index handling.  
 The checkpointing length at each level is defined in  
 file {\it tamc.h}, cf. below.  
 %  
 \item Cost function package: {\it pkg/cost/} \\  
 This package contains all relevant routines for  
 initialising, accumulating and finalizing the cost function  
 (see Section \ref{???}). \\  
 \hspace*{4ex} {\tt \#define ALLOW\_COST} \\  
 enables all general aspects of the cost function handling,  
 in particular the hooks in the foorward code for  
 initialising, accumulating and finalizing the cost function. \\  
 \hspace*{4ex} {\tt \#define ALLOW\_COST\_TRACER} \\  
 includes the call to the cost function for this  
 particular experiment, eqn. (\ref{cost_tracer}).  
 %  
 \item Control variable package: {\it pkg/ctrl/} \\  
 This package contains all relevant routines for  
 the handling of the control vector.  
 Each control variable can be enabled/disabled with its own flag: \\  
 \begin{tabular}{ll}  
 \hspace*{2ex} {\tt \#define ALLOW\_THETA0\_CONTROL} &  
 initial temperature \\  
 \hspace*{2ex} {\tt \#define ALLOW\_SALT0\_CONTROL} &  
 initial salinity \\  
 \hspace*{2ex} {\tt \#define ALLOW\_TR0\_CONTROL} &  
 initial passive tracer concentration \\  
 \hspace*{2ex} {\tt \#define ALLOW\_TAUU0\_CONTROL} &  
 zonal wind stress \\  
 \hspace*{2ex} {\tt \#define ALLOW\_TAUV0\_CONTROL} &  
 meridional wind stress \\  
 \hspace*{2ex} {\tt \#define ALLOW\_SFLUX0\_CONTROL} &  
 freshwater flux \\  
 \hspace*{2ex} {\tt \#define ALLOW\_HFLUX0\_CONTROL} &  
 heat flux \\  
 \hspace*{2ex} {\tt \#define ALLOW\_DIFFKR\_CONTROL} &  
 diapycnal diffusivity \\  
 \hspace*{2ex} {\tt \#undef ALLOW\_KAPPAGM\_CONTROL} &  
 isopycnal diffusivity \\  
 \end{tabular}  
764  %  %
765  \end{itemize}  \end{itemize}
766    
767  \subsubsection{File {\it SIZE.h}}  %------------------------------------------------------------------
   
 The file contains the grid point dimensions of the forward  
 model. It is identical to the {\it verification/exp2/}: \\  
 \hspace*{4ex} {\tt sNx = 90} \\  
 \hspace*{4ex} {\tt sNy = 40} \\  
 \hspace*{4ex} {\tt Nr = 20} \\  
 It correpsponds to a single-tile/single-processor setup:  
 {\tt nSx = nSy = 1, nPx = nPy = 1},  
 with standard overlap dimensioning  
 {\tt OLx = OLy = 3}.  
   
 \subsubsection{File {\it adcommon.h}}  
   
 This file contains common blocks of some adjoint variables  
 that are generated by TAMC.  
 The common blocks are used by the adjoint support routine  
 {\it addummy\_in\_stepping} which needs to access those variables:  
   
 \begin{tabular}{ll}  
 \hspace*{4ex} {\tt common /addynvars\_r/} &  
 \hspace*{4ex} is related to {\it DYNVARS.h} \\  
 \hspace*{4ex} {\tt common /addynvars\_cd/} &  
 \hspace*{4ex} is related to {\it DYNVARS.h} \\  
 \hspace*{4ex} {\tt common /addynvars\_diffkr/} &  
 \hspace*{4ex} is related to {\it DYNVARS.h} \\  
 \hspace*{4ex} {\tt common /addynvars\_kapgm/} &  
 \hspace*{4ex} is related to {\it DYNVARS.h} \\  
 \hspace*{4ex} {\tt common /adtr1\_r/} &  
 \hspace*{4ex} is related to {\it TR1.h} \\  
 \hspace*{4ex} {\tt common /adffields/} &  
 \hspace*{4ex} is related to {\it FFIELDS.h}\\  
 \end{tabular}  
   
 Note that if the structure of the common block changes in the  
 above header files of the forward code, the structure  
 of the adjoint common blocks will change accordingly.  
 Thus, it has to be made sure that the structure of the  
 adjoint common block in the hand-written file {\it adcommon.h}  
 complies with the automatically generated adjoint common blocks  
 in {\it adjoint\_model.F}.  
768    
769  \subsubsection{File {\it tamc.h}}  \subsection{Building the AD code
770    \label{section_ad_build}}
771    
772  This routine contains the dimensions for TAMC checkpointing.  The build process of an AD code is very similar to building
773    the forward model. However, depending on which AD code one wishes
774    to generate, and on which AD tool is available (TAF or TAMC),
775    the following {\tt make} targets are available:
776    
777    \begin{table}[h!]
778    {\footnotesize
779    \begin{tabular}{ccll}
780    ~ & {\it AD-target} & {\it output} & {\it description} \\
781    \hline
782    \hline
783    (1) & {\tt <MODE><TOOL>only} & {\tt <MODE>\_<TOOL>\_output.f}  &
784    generates code for $<$MODE$>$ using $<$TOOL$>$ \\
785    ~ & ~ & ~ & no {\tt make} dependencies on {\tt .F .h} \\
786    ~ & ~ & ~ & useful for compiling on remote platforms \\
787    \hline
788    (2) & {\tt <MODE><TOOL>} & {\tt <MODE>\_<TOOL>\_output.f}  &
789    generates code for $<$MODE$>$ using $<$TOOL$>$ \\
790    ~ & ~ & ~ & includes {\tt make} dependencies on {\tt .F .h} \\
791    ~ & ~ & ~ & i.e. input for $<$TOOL$>$ may be re-generated \\
792    \hline
793    (3) & {\tt <MODE>all} & {\tt mitgcmuv\_<MODE>}  &
794    generates code for $<$MODE$>$ using $<$TOOL$>$ \\
795    ~ & ~ & ~ & and compiles all code \\
796    ~ & ~ & ~ & (use of TAF is set as default) \\
797    \hline
798    \hline
799    \end{tabular}
800    }
801    \end{table}
802    %
803    Here, the following placeholders are used
804  %  %
805  \begin{itemize}  \begin{itemize}
806  %  %
807  \item {\tt \#ifdef ALLOW\_TAMC\_CHECKPOINTING} \\  \item [$<$TOOL$>$]
 3-level checkpointing is enabled, i.e. the timestepping  
 is divided into three different levels (see Section \ref{???}).  
 The model state of the outermost ({\tt nchklev\_3}) and the  
 intermediate ({\tt nchklev\_2}) timestepping loop are stored to file  
 (handled in {\it the\_main\_loop}).  
 The innermost loop ({\tt nchklev\_1})  
 avoids I/O by storing all required variables  
 to common blocks. This storing may also be necessary if  
 no checkpointing is chosen  
 (nonlinear functions, if-statements, iterative loops, ...).  
 In the present example the dimensions are chosen as follows: \\  
 \hspace*{4ex} {\tt nchklev\_1      =  36 } \\  
 \hspace*{4ex} {\tt nchklev\_2      =  30 } \\  
 \hspace*{4ex} {\tt nchklev\_3      =  60 } \\  
 To guarantee that the checkpointing intervals span the entire  
 integration period the following relation must be satisfied: \\  
 \hspace*{4ex} {\tt nchklev\_1*nchklev\_2*nchklev\_3 $ \ge $ nTimeSteps} \\  
 where {\tt nTimeSteps} is either specified in {\it data}  
 or computed via \\  
 \hspace*{4ex} {\tt nTimeSteps = (endTime-startTime)/deltaTClock }.  
 %  
 \item {\tt \#undef ALLOW\_TAMC\_CHECKPOINTING} \\  
 No checkpointing is enabled.  
 In this case the relevant counter is {\tt nchklev\_0}.  
 Similar to above, the following relation has to be satisfied \\  
 \hspace*{4ex} {\tt nchklev\_0 $ \ge $ nTimeSteps}.  
808  %  %
809  \end{itemize}  \begin{itemize}
   
 The following parameters may be worth describing: \\  
810  %  %
811  \hspace*{4ex} {\tt isbyte} \\  \item {\tt TAF}
812  \hspace*{4ex} {\tt maxpass} \\  \item {\tt TAMC}
 ~  
   
 \subsubsection{File {\it makefile}}  
   
 This file contains all relevant paramter flags and  
 lists to run TAMC or TAF.  
 It is assumed that TAMC is available to you, either locally,  
 being installed on your network, or remotely through the 'TAMC Utility'.  
 TAMC is called with the command {\tt tamc} followed by a  
 number of options. They are described in detail in the  
 TAMC manual \cite{gie:99}.  
 Here we briefly discuss the main flags used in the {\it makefile}  
813  %  %
 \begin{itemize}  
 \item [{\tt tamc}] {\tt  
 -input <variable names>  
 -output <variable name> -r4 ... \\  
 -toplevel <S/R name> -reverse <file names>  
 }  
814  \end{itemize}  \end{itemize}
815  %  %
816    \item [$<$MODE$>$]
817    %
818  \begin{itemize}  \begin{itemize}
819  %  %
820  \item {\tt -toplevel <S/R name>} \\  \item {\tt ad} generates the adjoint model (ADM)
821  Name of the toplevel routine, with respect to which the  \item {\tt ftl} generates the tangent linear model (TLM)
822  control flow analysis is performed.  \item {\tt svd} generates both ADM and TLM for \\
823  %  singular value decomposition (SVD) type calculations
 \item {\tt -input <variable names>} \\  
 List of independent variables $ u $ with respect to which the  
 dependent variable $ J $ is differentiated.  
 %  
 \item {\tt -output <variable name>} \\  
 Dependent variable $ J $  which is to be differentiated.  
 %  
 \item {\tt -reverse <file names>} \\  
 Adjoint code is generated to compute the sensitivity of an  
 independent variable w.r.t.  many dependent variables.  
 In the discussion of Section ???  
 the generated adjoint top-level routine computes the product  
 of the transposed Jacobian matrix $ M^T $ times  
 the gradient vector $ \nabla_v J $.  
 \\  
 {\tt <file names>} refers to the list of files {\it .f} which are to be  
 analyzed by TAMC. This list is generally smaller than the full list  
 of code to be compiled. The files not contained are either  
 above the top-level routine (some initialisations), or are  
 deliberately hidden from TAMC, either because hand-written  
 adjoint routines exist, or the routines must not (or don't have to)  
 be differentiated. For each routine which is part of the flow tree  
 of the top-level routine, but deliberately hidden from TAMC  
 (or for each package which contains such routines),  
 a corresponding file {\it .flow} exists containing flow directives  
 for TAMC.  
824  %  %
825  \item {\tt -r4} \\  \end{itemize}
 ~  
826  %  %
827  \end{itemize}  \end{itemize}
828    
829    For example, to generate the adjoint model using TAF after routines ({\tt .F})
830    or headers ({\tt .h}) have been modified, but without compilation,
831    type {\tt make adtaf};
832    or, to generate the tangent linear model using TAMC without
833    re-generating the input code, type {\tt make ftltamconly}.
834    
 \subsubsection{The input parameter files}  
   
 \paragraph{File {\it data}}  
   
 \paragraph{File {\it data.cost}}  
   
 \paragraph{File {\it data.ctrl}}  
   
 \paragraph{File {\it data.gmredi}}  
   
 \paragraph{File {\it data.grdchk}}  
   
 \paragraph{File {\it data.optim}}  
   
 \paragraph{File {\it data.pkg}}  
   
 \paragraph{File {\it eedata}}  
   
 \paragraph{File {\it topog.bin}}  
   
 \paragraph{File {\it windx.bin, windy.bin}}  
835    
836  \paragraph{File {\it salt.bin, theta.bin}}  A typical full build process to generate the ADM via TAF would
837    look like follows:
838  \paragraph{File {\it SSS.bin, SST.bin}}  \begin{verbatim}
839    % mkdir build
840  \paragraph{File {\it pickup*}}  % cd build
841    % ../../../tools/genmake2 -mods=../code_ad
842  \subsection{Compiling the model and its adjoint}  % make depend
843    % make adall
844    \end{verbatim}
845    
846  The built process of the adjoint model is slightly more  %------------------------------------------------------------------
 complex than that of compiling the forward code.  
 The main reason is that the adjoint code generation requires  
 a specific list of routines that are to be differentiated  
 (as opposed to the automatic generation of a list of  
 files to be compiled by genmake).  
 This list excludes routines that don't have to be or must not be  
 differentiated. For some of the latter routines flow directives  
 may be necessary, a list of which has to be given as well.  
 For this reason, a separate {\it makefile} is currently  
 maintained in the directory {\tt adjoint/}. This  
 makefile is responsible for the adjoint code generation.  
847    
848  In the following we describe the build process step by step,  \subsection{The AD build process in detail
849  assuming you are in the directory {\tt bin/}.  \label{section_ad_build_detail}}
 A summary of steps to follow is given at the end.  
850    
851  \paragraph{Adjoint code generation and compilation -- step by step}  The {\tt make <MODE>all} target consists of the following procedures:
852    
853  \begin{enumerate}  \begin{enumerate}
854  %  %
855  \item  \item
856  {\tt ln -s ../verification/???/code/.genmakerc .} \\  A header file {\tt AD\_CONFIG.h} is generated which contains a CPP option
857  {\tt ln -s ../verification/???/code/*.[Fh] .} \\  on which code ought to be generated. Depending on the {\tt make} target,
858  Link your customized genmake options, header files,  the contents is
859  and modified code to the compile directory.  \begin{itemize}
 %  
860  \item  \item
861  {\tt ../tools/genmake -makefile} \\  {\tt \#define ALLOW\_ADJOINT\_RUN}
 Generate your Makefile (cf. Section ???).  
 %  
862  \item  \item
863  {\tt make depend} \\  {\tt \#define ALLOW\_TANGENTLINEAR\_RUN}
 Dependency analysis for the CPP pre-compiler (cf. Section ???).  
 %  
864  \item  \item
865  {\tt make small\_f} \\  {\tt \#define ALLOW\_ECCO\_OPTIMIZATION}
866  This is the first difference between forward code compilation  \end{itemize}
 and adjoint code generation and compilation.  
 Instead of going through the entire compilation process  
 (CPP precompiling -- {\tt .f}, object code generation -- {\tt .o},  
 linking of object files and libraries to generate executable),  
 only the CPP compiler is invoked at this stage to generate  
 the {\tt .f} files.  
867  %  %
868  \item  \item
869  {\tt cd ../adjoint} \\  A single file {\tt <MODE>\_input\_code.f} is concatenated
870  {\tt make adtaf} or {\tt make adtamc} \\  consisting of all {\tt .f} files that are part of the list {\bf AD\_FILES}
871  Depending on whether you have TAF or TAMC at your disposal,  and all {\tt .flow} files that are part of the list {\bf AD\_FLOW\_FILES}.
 you'll choose {\tt adtaf} or {\tt adtamc} as your  
 make target for the {\it makefile} in the directory {\tt adjoint/}.  
 Several things happen at this stage.  
 %  
 \begin{enumerate}  
872  %  %
873  \item  \item
874  The initial template file {\it adjoint\_model.F} which is part  The AD tool is invoked with the {\bf <MODE>\_<TOOL>\_FLAGS}.
875  of the compiling list created by {\it genmake} is restored.  The default AD tool flags in {\tt genmake2} can be overrwritten by
876    an {\tt adjoint\_options} file (similar to the platform-specific
877    {\tt build\_options}, see Section ???.
878    The AD tool writes the resulting AD code into the file
879    {\tt <MODE>\_input\_code\_ad.f}
880  %  %
881  \item  \item
882  All Fortran routines {\tt *.f} in {\tt bin/} are  A short sed script {\tt adjoint\_sed} is applied to
883  concatenated into a single file (it's current name is  {\tt <MODE>\_input\_code\_ad.f}
884  {\it tamc\_code.f}).  to reinstate {\bf myThid} into the CALL argument list of active file I/O.
885    The result is written to file {\tt <MODE>\_<TOOL>\_output.f}.
886  %  %
887  \item  \item
888  Adjoint code is generated by TAMC or TAF.  All routines are compiled and an executable is generated
889  The adjoint code is written to the file {\it tamc\_code\_ad.f}.  (see Table ???).
 It contains all adjoint routines of the forward routines  
 concatenated in {\it tamc\_code.f}.  
 For a given forward routines {\tt subroutine routinename}  
 the adjoint routine is named {\tt adsubroutine routinename}  
 by default (that default can be changed via the flag  
 {\tt -admark <markname>}).  
 Furthermore, it may contain modified code which  
 incorporates the translation of adjoint store directives  
 into specific Fortran code.  
 For a given forward routines {\tt subroutine routinename}  
 the modified routine is named {\tt mdsubroutine routinename}.  
 TAMC or TAF info is written to file  
 {\it tamc\_code.prot} or {\it taf.log}, respectively.  
890  %  %
891  \end{enumerate}  \end{enumerate}
892    
893    \subsubsection{The list AD\_FILES and {\tt .list} files}
894    
895    Not all routines are presented to the AD tool.
896    Routines typically hidden are diagnostics routines which
897    do not influence the cost function, but may create
898    artificial flow dependencies such as I/O of active variables.
899    
900    {\tt genmake2} generates a list (or variable) {\bf AD\_FILES}
901    which contains all routines that are shown to the AD tool.
902    This list is put together from all files with suffix {\tt .list}
903    that {\tt genmake2} finds in its search directories.
904    The list file for the core MITgcm routines is in {\tt model/src/}
905    is called {\tt model\_ad\_diff.list}.
906    Note that no wrapper routine is shown to TAF. These are either
907    not visible at all to the AD code, or hand-written AD code
908    is available (see next section).
909    
910    Each package directory contains its package-specific
911    list file {\tt <PKG>\_ad\_diff.list}. For example,
912    {\tt pkg/ptracers/} contains the file {\tt ptracers\_ad\_diff.list}.
913    Thus, enabling a package will automatically extend the
914    {\bf AD\_FILES} list of {\tt genmake2} to incorporate the
915    package-specific routines.
916    Note that you will need to regenerate the {\tt Makefile} if
917    you enable a package (e.g. by adding it to {\tt packages.conf})
918    and a {\tt Makefile} already exists.
919    
920    \subsubsection{The list AD\_FLOW\_FILES and {\tt .flow} files}
921    
922    TAMC and TAF can evaluate user-specified directives
923    that start with a specific syntax ({\tt CADJ}, {\tt C\$TAF}, {\tt !\$TAF}).
924    The main categories of directives are STORE directives and
925    FLOW directives. Here, we are concerned with flow directives,
926    store directives are treated elsewhere.
927    
928    Flow directives enable the AD tool to evaluate how it should treat
929    routines that are 'hidden' by the user, i.e. routines which are
930    not contained in the {\bf AD\_FILES} list (see previous section),
931    but which are called in part of the code that the AD tool does see.
932    The flow directive tell the AD tool
933  %  %
934  \item  \begin{itemize}
 {\tt make adchange} \\  
 The multi-threading capability of the MITGCM requires a slight  
 change in the parameter list of some routines that are related to  
 to active file handling.  
 This postprocessing invokes the sed script {\it adjoint\_ecco\_sed.com}  
 to insert the threading counter {\bf myThId} into the parameter list  
 of those subroutines.  
 The resulting code is written to file {\it tamc\_code\_sed\_ad.f}  
 and appended to the file {\it adjoint\_model.F}.  
 This concludes the adjoint codel generation.  
935  %  %
936  \item  \item which subroutine arguments are input/output
937  {\tt cd ../bin} \\  \item which subroutine arguments are active
938  {\tt make} \\  \item which subroutine arguments are required to compute the cost
939  The file {\it adjoint\_model.F} now contains the full adjoint code.  \item which subroutine arguments are dependent
 All routines are now compiled.  
940  %  %
941  \end{enumerate}  \end{itemize}
942    %
943    The syntax for the flow directives can be found in the
944    AD tool manuals.
945    
946  \paragraph{Adjoint code generation and compilation -- summary}  {\tt genmake2} generates a list (or variable) {\bf AD\_FLOW\_FILES}
947  ~ \\  which contains all files with suffix{\tt .flow} that it finds
948    in its search directories.
949    The flow directives for the core MITgcm routines of
950    {\tt eesupp/src/} and {\tt model/src/}
951    reside in {\tt pkg/autodiff/}.
952    This directory also contains hand-written adjoint code
953    for the MITgcm WRAPPER (see Section ???).
954    
955    Flow directives for package-specific routines are contained in
956    the corresponding package directories in the file
957    {\tt <PKG>\_ad.flow}, e.g. ptracers-specific directives are in
958    {\tt ptracers\_ad.flow}.
959    
960    \subsubsection{Store directives for 3-level checkpointing}
961    
962    The storing that is required at each period of the
963    3-level checkpointing is controled by three
964    top-level headers.
965    
966  \[  \begin{verbatim}
967  \boxed{  do ilev_3 = 1, nchklev_3
968  \begin{split}  #  include ``checkpoint_lev3.h''
969   ~ & \mbox{\tt cd bin} \\     do ilev_2 = 1, nchklev_2
970   ~ & \mbox{\tt ln -s ../verification/my\_experiment/code/.genmakerc .} \\  #     include ``checkpoint_lev2.h''
971   ~ & \mbox{\tt ln -s ../verification/my\_experiment/code/*.[Fh] .} \\        do ilev_1 = 1, nchklev_1
972   ~ & \mbox{\tt ../tools/genmake -makefile} \\  #        include ``checkpoint_lev1.h''
973   ~ & \mbox{\tt make depend} \\  
974   ~ & \mbox{\tt make small\_f} \\  ...
975   ~ & \mbox{\tt cd ../adjoint} \\  
976   ~ & \mbox{\tt make adtaf <OR: make adtamc>} \\        end do
977   ~ & \mbox{\tt make adchange} \\     end do
978   ~ & \mbox{\tt cd ../bin} \\  end do
979   ~ & \mbox{\tt make} \\  \end{verbatim}
 \end{split}  
 }  
 \]  
980    
981  \newpage  All files {\tt checkpoint\_lev?.h} are contained in directory
982    {\tt pkg/autodiff/}.
983    
 %**********************************************************************  
 \section{TLM and ADM generation in general}  
 \label{sec_ad_setup_gen}  
 %**********************************************************************  
984    
985  In this section we describe in a general fashion  \subsubsection{Changing the default AD tool flags: ad\_options files}
 the parts of the code that are relevant for automatic  
 differentiation using the software tool TAMC.  
986    
 \begin{figure}[b!]  
 \input{part5/doc_ad_the_model}  
 \caption{~}  
 \label{fig:adthemodel}  
 \end{figure}  
987    
988  The basic flow is depicted in \ref{fig:adthemodel}.  \subsubsection{Hand-written adjoint code}
989  If the option {\tt ALLOW\_AUTODIFF\_TAMC} is defined, the driver routine  
990  {\it the\_model\_main}, instead of calling {\it the\_main\_loop},  %------------------------------------------------------------------
 invokes the adjoint of this routine, {\it adthe\_main\_loop},  
 which is the toplevel routine in terms of reverse mode computation.  
 The routine {\it adthe\_main\_loop} has been generated using TAMC.  
 It contains both the forward integration of the full model,  
 any additional storing that is required for efficient checkpointing,  
 and the reverse integration of the adjoint model.  
 The structure of {\it adthe\_main\_loop} has been strongly  
 simplified for clarification; in particular, no checkpointing  
 procedures are shown here.  
 Prior to the call of {\it adthe\_main\_loop}, the routine  
 {\it ctrl\_unpack} is invoked to unpack the control vector,  
 and following that call, the routine {\it ctrl\_pack}  
 is invoked to pack the control vector  
 (cf. Section \ref{section_ctrl}).  
 If gradient checks are to be performed, the option  
 {\tt ALLOW\_GRADIENT\_CHECK} is defined. In this case  
 the driver routine {\it grdchk\_main} is called after  
 the gradient has been computed via the adjoint  
 (cf. Section \ref{section_grdchk}).  
991    
992  \subsection{The cost function (dependent variable)  \subsection{The cost function (dependent variable)
993  \label{section_cost}}  \label{section_cost}}
# Line 1293  the gradient has been computed via the a Line 995  the gradient has been computed via the a
995  The cost function $ {\cal J} $ is referred to as the {\sf dependent variable}.  The cost function $ {\cal J} $ is referred to as the {\sf dependent variable}.
996  It is a function of the input variables $ \vec{u} $ via the composition  It is a function of the input variables $ \vec{u} $ via the composition
997  $ {\cal J}(\vec{u}) \, = \, {\cal J}(M(\vec{u})) $.  $ {\cal J}(\vec{u}) \, = \, {\cal J}(M(\vec{u})) $.
998  The input is referred to as the  The input are referred to as the
999  {\sf independent variables} or {\sf control variables}.  {\sf independent variables} or {\sf control variables}.
1000  All aspects relevant to the treatment of the cost function $ {\cal J} $  All aspects relevant to the treatment of the cost function $ {\cal J} $
1001  (parameter setting, initialisation, accumulation,  (parameter setting, initialization, accumulation,
1002  final evaluation), are controlled by the package {\it pkg/cost}.  final evaluation), are controlled by the package {\it pkg/cost}.
1003    The aspects relevant to the treatment of the independent variables
1004    are controlled by the package {\it pkg/ctrl} and will be treated
1005    in the next section.
1006    
 \begin{figure}[h!]  
1007  \input{part5/doc_cost_flow}  \input{part5/doc_cost_flow}
 \caption{~}  
 \label{fig:costflow}  
 \end{figure}  
1008    
1009  \subsubsection{genmake and CPP options}  \subsubsection{Enabling the package}
1010  %  
 \begin{itemize}  
 %  
 \item  
1011  \fbox{  \fbox{
1012  \begin{minipage}{12cm}  \begin{minipage}{12cm}
1013  {\it genmake}, {\it CPP\_OPTIONS.h}, {\it ECCO\_CPPOPTIONS.h}  {\it packages.conf}, {\it ECCO\_CPPOPTIONS.h}
1014  \end{minipage}  \end{minipage}
1015  }  }
1016  \end{itemize}  \begin{itemize}
 %  
 The directory {\it pkg/cost} can be included to the  
 compile list in 3 different ways (cf. Section \ref{???}):  
1017  %  %
1018  \begin{enumerate}  \item
1019    The package is enabled by adding {\it cost} to your file {\it packages.conf}
1020    (see Section ???)
1021  %  %
1022  \item {\it genmake}: \\  \item
1023  Change the default settings in the file {\it genmake} by adding  
1024  {\bf cost} to the {\bf enable} list (not recommended).  
1025  %  \end{itemize}
 \item {\it .genmakerc}: \\  
 Customize the settings of {\bf enable}, {\bf disable} which are  
 appropriate for your experiment in the file {\it .genmakerc}  
 and add the file to your compile directory.  
 %  
 \item genmake-options: \\  
 Call {\it genmake} with the option  
 {\tt genmake -enable=cost}.  
1026  %  %
1027  \end{enumerate}  
1028    N.B.: In general the following packages ought to be enabled
1029    simultaneously: {\it autodiff, cost, ctrl}.
1030  The basic CPP option to enable the cost function is {\bf ALLOW\_COST}.  The basic CPP option to enable the cost function is {\bf ALLOW\_COST}.
1031  Each specific cost function contribution has its own option.  Each specific cost function contribution has its own option.
1032  For the present example the option is {\bf ALLOW\_COST\_TRACER}.  For the present example the option is {\bf ALLOW\_COST\_TRACER}.
1033  All cost-specific options are set in {\it ECCO\_CPPOPTIONS.h}  All cost-specific options are set in {\it ECCO\_CPPOPTIONS.h}
1034  Since the cost function is usually used in conjunction with  Since the cost function is usually used in conjunction with
1035  automatic differentiation, the CPP option  automatic differentiation, the CPP option
1036  {\bf ALLOW\_ADJOINT\_RUN} should be defined  {\bf ALLOW\_ADJOINT\_RUN} (file {\it CPP\_OPTIONS.h}) and
1037  (file {\it CPP\_OPTIONS.h}).  {\bf ALLOW\_AUTODIFF\_TAMC} (file {\it ECCO\_CPPOPTIONS.h})
1038    should be defined.
1039    
1040  \subsubsection{Initialisation}  \subsubsection{Initialization}
1041  %  %
1042  The initialisation of the {\it cost} package is readily enabled  The initialization of the {\it cost} package is readily enabled
1043  as soon as the CPP option {\bf ALLOW\_ADJOINT\_RUN} is defined.  as soon as the CPP option {\bf ALLOW\_COST} is defined.
1044  %  %
1045  \begin{itemize}  \begin{itemize}
1046  %  %
# Line 1378  Variables: {\it cost\_init} Line 1070  Variables: {\it cost\_init}
1070  }  }
1071  \\  \\
1072  This S/R  This S/R
1073  initialises the different cost function contributions.  initializes the different cost function contributions.
1074  The contribtion for the present example is {\bf objf\_tracer}  The contribution for the present example is {\bf objf\_tracer}
1075  which is defined on each tile (bi,bj).  which is defined on each tile (bi,bj).
1076  %  %
1077  \end{itemize}  \end{itemize}
# Line 1422  from each contribution and sums over all Line 1114  from each contribution and sums over all
1114  \begin{equation}  \begin{equation}
1115  {\cal J} \, = \,  {\cal J} \, = \,
1116  {\rm fc} \, = \,  {\rm fc} \, = \,
1117  {\rm mult\_tracer} \sum_{bi,\,bj}^{nSx,\,nSy}  {\rm mult\_tracer} \sum_{\text{global sum}} \sum_{bi,\,bj}^{nSx,\,nSy}
1118  {\rm objf\_tracer}(bi,bj) \, + \, ...  {\rm objf\_tracer}(bi,bj) \, + \, ...
1119  \end{equation}  \end{equation}
1120  %  %
# Line 1434  tamc -output 'fc' ... Line 1126  tamc -output 'fc' ...
1126    
1127  %%%% \end{document}  %%%% \end{document}
1128    
 \begin{figure}  
1129  \input{part5/doc_ad_the_main}  \input{part5/doc_ad_the_main}
 \caption{~}  
 \label{fig:adthemain}  
 \end{figure}  
1130    
1131  \subsection{The control variables (independent variables)  \subsection{The control variables (independent variables)
1132  \label{section_ctrl}}  \label{section_ctrl}}
# Line 1455  as variable assignments. Therefore, file Line 1143  as variable assignments. Therefore, file
1143  active variables are written and from which active variables  active variables are written and from which active variables
1144  are read are called {\sf active files}.  are read are called {\sf active files}.
1145  All aspects relevant to the treatment of the control variables  All aspects relevant to the treatment of the control variables
1146  (parameter setting, initialisation, perturbation)  (parameter setting, initialization, perturbation)
1147  are controled by the package {\it pkg/ctrl}.  are controlled by the package {\it pkg/ctrl}.
1148    
 \begin{figure}[h!]  
1149  \input{part5/doc_ctrl_flow}  \input{part5/doc_ctrl_flow}
 \caption{~}  
 \label{fig:ctrlflow}  
 \end{figure}  
1150    
1151  \subsubsection{genmake and CPP options}  \subsubsection{genmake and CPP options}
1152  %  %
# Line 1478  are controled by the package {\it pkg/ct Line 1162  are controled by the package {\it pkg/ct
1162  %  %
1163  To enable the directory to be included to the compile list,  To enable the directory to be included to the compile list,
1164  {\bf ctrl} has to be added to the {\bf enable} list in  {\bf ctrl} has to be added to the {\bf enable} list in
1165  {\it .genmakerc} (or {\it genmake} itself).  {\it .genmakerc} or in {\it genmake} itself (analogous to {\it cost}
1166    package, cf. previous section).
1167  Each control variable is enabled via its own CPP option  Each control variable is enabled via its own CPP option
1168  in {\it ECCO\_CPPOPTIONS.h}.  in {\it ECCO\_CPPOPTIONS.h}.
1169    
1170  \subsubsection{Initialisation}  \subsubsection{Initialization}
1171  %  %
1172  \begin{itemize}  \begin{itemize}
1173  %  %
# Line 1522  Two important issues related to the hand Line 1207  Two important issues related to the hand
1207  variables in the MITGCM need to be addressed.  variables in the MITGCM need to be addressed.
1208  First, in order to save memory, the control variable arrays  First, in order to save memory, the control variable arrays
1209  are not kept in memory, but rather read from file and added  are not kept in memory, but rather read from file and added
1210  to the initial fields during the model initialisation phase.  to the initial fields during the model initialization phase.
1211  Similarly, the corresponding adjoint fields which represent  Similarly, the corresponding adjoint fields which represent
1212  the gradient of the cost function w.r.t. the control variables  the gradient of the cost function w.r.t. the control variables
1213  are written to file at the end of the adjoint integration.  are written to file at the end of the adjoint integration.
# Line 1602  dummy variable {\bf xx\_tr1\_dummy} is i Line 1287  dummy variable {\bf xx\_tr1\_dummy} is i
1287  and an 'active read' routine of the adjoint support  and an 'active read' routine of the adjoint support
1288  package {\it pkg/autodiff} is invoked.  package {\it pkg/autodiff} is invoked.
1289  The read-procedure is tagged with the variable  The read-procedure is tagged with the variable
1290  {\bf xx\_tr1\_dummy} enabbling TAMC to recognize the  {\bf xx\_tr1\_dummy} enabling TAMC to recognize the
1291  initialisation of the perturbation.  initialization of the perturbation.
1292  The modified call of TAMC thus reads  The modified call of TAMC thus reads
1293  %  %
1294  \begin{verbatim}  \begin{verbatim}
# Line 1622  in the code takes on the form Line 1307  in the code takes on the form
1307  %  %
1308  Note, that reading an active variable corresponds  Note, that reading an active variable corresponds
1309  to a variable assignment. Its derivative corresponds  to a variable assignment. Its derivative corresponds
1310  to a write statement of the adjoint variable.  to a write statement of the adjoint variable, followed by
1311    a reset.
1312  The 'active file' routines have been designed  The 'active file' routines have been designed
1313  to support active read and corresponding adjoint active write  to support active read and corresponding adjoint active write
1314  operations (and vice versa).  operations (and vice versa).
# Line 1714  variables are written to {\bf adxx\_ ... Line 1400  variables are written to {\bf adxx\_ ...
1400  \begin{itemize}  \begin{itemize}
1401  %  %
1402  \item {\bf vector\_ctrl}: the control vector \\  \item {\bf vector\_ctrl}: the control vector \\
1403  At the very beginning of the model initialisation,  At the very beginning of the model initialization,
1404  the updated compressed control vector is read (or initialised)  the updated compressed control vector is read (or initialised)
1405  and distributed to 2-dim. and 3-dim. control variable fields.  and distributed to 2-dim. and 3-dim. control variable fields.
1406  %  %
# Line 1739  at intermediate times can be written usi Line 1425  at intermediate times can be written usi
1425  {\it addummy\_in\_stepping}.  {\it addummy\_in\_stepping}.
1426  This routine is part of the adjoint support package  This routine is part of the adjoint support package
1427  {\it pkg/autodiff} (cf.f. below).  {\it pkg/autodiff} (cf.f. below).
1428    The procedure is enabled using via the CPP-option
1429    {\bf ALLOW\_AUTODIFF\_MONITOR} (file {\it ECCO\_CPPOPTIONS.h}).
1430  To be part of the adjoint code, the corresponding S/R  To be part of the adjoint code, the corresponding S/R
1431  {\it dummy\_in\_stepping} has to be called in the forward  {\it dummy\_in\_stepping} has to be called in the forward
1432  model (S/R {\it the\_main\_loop}) at the appropriate place.  model (S/R {\it the\_main\_loop}) at the appropriate place.
1433    The adjoint common blocks are extracted from the adjoint code
1434    via the header file {\it adcommon.h}.
1435    
1436  {\it dummy\_in\_stepping} is essentially empty,  {\it dummy\_in\_stepping} is essentially empty,
1437  the corresponding adjoint routine is hand-written rather  the corresponding adjoint routine is hand-written rather
# Line 1768  the common blocks Line 1458  the common blocks
1458  {\bf /adtr1\_r/}, {\bf /adffields/},  {\bf /adtr1\_r/}, {\bf /adffields/},
1459  which have been extracted from the adjoint code to enable  which have been extracted from the adjoint code to enable
1460  access to the adjoint variables.  access to the adjoint variables.
1461    
1462    {\bf WARNING:} If the structure of the common blocks
1463    {\bf /dynvars\_r/}, {\bf /dynvars\_cd/}, etc., changes
1464    similar changes will occur in the adjoint common blocks.
1465    Therefore, consistency between the TAMC-generated common blocks
1466    and those in {\it adcommon.h} have to be checked.
1467  %  %
1468  \end{itemize}  \end{itemize}
1469    
# Line 1782  The gradient $ \nabla _{u}{\cal J} |_{u_ Line 1478  The gradient $ \nabla _{u}{\cal J} |_{u_
1478  with the value of the cost function itself $ {\cal J}(u_{[k]}) $  with the value of the cost function itself $ {\cal J}(u_{[k]}) $
1479  at iteration step $ k $ serve  at iteration step $ k $ serve
1480  as input to a minimization routine (e.g. quasi-Newton method,  as input to a minimization routine (e.g. quasi-Newton method,
1481  conjugate gradient, ... \cite{gil_lem:89})  conjugate gradient, ... \cite{gil-lem:89})
1482  to compute an update in the  to compute an update in the
1483  control variable for iteration step $k+1$  control variable for iteration step $k+1$
1484  \[  \[
# Line 1913  to {\it adxx\_...$<$k$>$}, again via the Line 1609  to {\it adxx\_...$<$k$>$}, again via the
1609  Finally, {\it ctrl\_pack} collects all adjoint files  Finally, {\it ctrl\_pack} collects all adjoint files
1610  and writes them to the compressed vector file  and writes them to the compressed vector file
1611  {\bf vector\_grad\_$<$k$>$}.  {\bf vector\_grad\_$<$k$>$}.
   
 \subsection{TLM and ADM generation via TAMC}  
   
   
   
 \subsection{Flow directives and adjoint support routines \label{section_flowdir}}  
   
 \subsection{Store directives and checkpointing \label{section_checkpointing}}  
   
 \subsection{Gradient checks \label{section_grdchk}}  
   
 \subsection{Second derivative generation via TAMC}  
   
 \section{Example of adjoint code}  

Legend:
Removed from v.1.6  
changed lines
  Added in v.1.18

  ViewVC Help
Powered by ViewVC 1.1.22