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

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

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

-revision 1.4 by heimbach,
Fri Oct  5 22:22:20 2001 UTC
+revision 1.23 by jmc,
Mon Aug 30 23:09:19 2010 UTC
 Line 1
  % $Header$
  % $Name$
+ Author: Patrick Heimbach
+ \label{ask_the_author:doc_ad_2}
  {\sf Automatic differentiation} (AD), also referred to as algorithmic
  (or, more loosely, computational) differentiation, involves
- automatically deriving code to calculate
+ automatically deriving code to calculate partial derivatives from an
- partial derivatives from an existing fully non-linear prognostic code.
+ existing fully non-linear prognostic code.  (see \cite{gri:00}).  A
- (see \cite{gri:00}).
+ software tool is used that parses and transforms source files
- A software tool is used that parses and transforms source files
+ according to a set of linguistic and mathematical rules.  AD tools are
- according to a set of linguistic and mathematical rules.
+ like source-to-source translators in that they parse a program code as
- AD tools are like source-to-source translators in that
+ input and produce a new program code as output
- they parse a program code as input and produce a new program code
+ (we restrict our discussion to source-to-source tools, ignoring
- as output.
+ operator-overloading tools).  However, unlike a
- However, unlike a pure source-to-source translation, the output program
+ pure source-to-source translation, the output program represents a new
- represents a new algorithm, such as the evaluation of the
+ algorithm, such as the evaluation of the Jacobian, the Hessian, or
- Jacobian, the Hessian, or higher derivative operators.
+ higher derivative operators.  In principle, a variety of derived
- In principle, a variety of derived algorithms
+ algorithms can be generated automatically in this way.
- can be generated automatically in this way.
+ MITgcm has been adapted for use with the Tangent linear and Adjoint
- The MITGCM has been adapted for use with the
+ Model Compiler (TAMC) and its successor TAF (Transformation of
- Tangent linear and Adjoint Model Compiler (TAMC) and its successor TAF
+ Algorithms in Fortran), developed by Ralf Giering (\cite{gie-kam:98},
- (Transformation of Algorithms in Fortran), developed
+ \cite{gie:99,gie:00}).  The first application of the adjoint of MITgcm
- by Ralf Giering (\cite{gie-kam:98}, \cite{gie:99,gie:00}).
+ for sensitivity studies has been published by \cite{maro-eta:99}.
- The first application of the adjoint of the MITGCM for senistivity
+ \cite{stam-etal:97,stam-etal:02} use MITgcm and its adjoint for ocean
- studies has been published by \cite{maro-eta:99}.
+ state estimation studies.  In the following we shall refer to TAMC and
- \cite{sta-eta:97,sta-eta:01} use the MITGCM and its adjoint
+ TAF synonymously, except were explicitly stated otherwise.
- for ocean state estimation studies.
- In the following we shall refer to TAMC and TAF synonymously,
+ As of mid-2007 we are also able to generate fairly efficient
- except were explicitly stated otherwise.
+ adjoint code of the MITgcm using a new, open-source AD tool,
+ called OpenAD (see \cite{naum-etal:06,utke-etal:08}.
- TAMC exploits the chain rule for computing the first
+ This enables us for the first time to compare adjoint models
- derivative of a function with
+ generated from different AD tools, providing an additional
- respect to a set of input variables.
+ accuracy check, complementary to finite-difference gradient checks.
- Treating a given forward code as a composition of operations --
+ OpenAD and its application to  MITgcm is described in detail
- each line representing a compositional element, the chain rule is
+ in section \ref{sec_ad_openad}.
- rigorously applied to the code, line by line. The resulting
- tangent linear or adjoint code,
+ The AD tool exploits the chain rule for computing the first derivative of a
- then, may be thought of as the composition in
+ function with respect to a set of input variables.  Treating a given
- forward or reverse order, respectively, of the
+ forward code as a composition of operations -- each line representing
- Jacobian matrices of the forward code's compositional elements.
+ a compositional element, the chain rule is rigorously applied to the
+ code, line by line. The resulting tangent linear or adjoint code,
+ then, may be thought of as the composition in forward or reverse
+ order, respectively, of the Jacobian matrices of the forward code's
+ compositional elements.
  %**********************************************************************
  \section{Some basic algebra}
  \label{sec_ad_algebra}
+ \begin{rawhtml}
+ <!-- CMIREDIR:sec_ad_algebra: -->
+ \end{rawhtml}
  %**********************************************************************
  Let $ \cal{M} $ be a general nonlinear, model, i.e. a
-Line 52 
 $\vec{u}=(u_1,\ldots,u_m)$
+Line 62 
 $\vec{u}=(u_1,\ldots,u_m)$
  such as forcing functions) to the $n$-dimensional space
  $V \subset I\!\!R^n$ of
  model output variable $\vec{v}=(v_1,\ldots,v_n)$
- (model state, model diagnostcs, objective function, ...)
+ (model state, model diagnostics, objective function, ...)
  under consideration,
  %
  \begin{equation}
- \begin{split}
+ \begin{aligned}
  {\cal M} \, : & \, U \,\, \longrightarrow \, V \\
  ~      & \, \vec{u} \,\, \longmapsto \, \vec{v} \, = \,
  {\cal M}(\vec{u})
  \label{fulloperator}
- \end{split}
+ \end{aligned}
  \end{equation}
  %
  The vectors $ \vec{u} \in U $ and $ v \in V $ may be represented w.r.t.
-Line 141 
 w.r.t. their corresponding inner product
+Line 151 
 w.r.t. their corresponding inner product
  $\left\langle \,\, , \,\, \right\rangle $
  %
  \begin{equation}
- \begin{split}
+ \begin{aligned}
  {\cal J} & = \,
  {\cal J} |_{\vec{u}^{(0)}} \, + \,
  \left\langle \, \nabla _{u}{\cal J}^T |_{\vec{u}^{(0)}} \, , \, \delta \vec{u} \, \right\rangle
-Line 150 
 $\left\langle \,\, , \,\, \right\rangle
+Line 160 
 $\left\langle \,\, , \,\, \right\rangle
  {\cal J} |_{\vec{v}^{(0)}} \, + \,
  \left\langle \, \nabla _{v}{\cal J}^T |_{\vec{v}^{(0)}} \, , \, \delta \vec{v} \, \right\rangle
  \, + \, O(\delta \vec{v}^2)
- \end{split}
+ \end{aligned}
  \label{deljidentity}
  \end{equation}
  %
-Line 191 
 the gradient $ \nabla _{u}{\cal J} $ can
+Line 201 
 the gradient $ \nabla _{u}{\cal J} $ can
  invoking the adjoint $ M^{\ast } $ of the tangent linear model $ M $
  %
  \begin{equation}
- \begin{split}
+ \begin{aligned}
  \nabla _{u}{\cal J}^T |_{\vec{u}} &
  = \, M^T |_{\vec{u}} \cdot \nabla _{v}{\cal J}^T |_{\vec{v}}  \\
  ~ & = \, M^T |_{\vec{u}} \cdot \delta \vec{v}^{\ast} \\
  ~ & = \, \delta \vec{u}^{\ast}
- \end{split}
+ \end{aligned}
  \label{adjoint}
  \end{equation}
  %
-Line 220 
 model integration,
+Line 230 
 model integration,
  starting at step 0 and moving up to step $\Lambda$, with intermediate
  ${\cal M}_{\lambda} (\vec{u}) = \vec{v}^{(\lambda+1)}$ and final
  ${\cal M}_{\Lambda} (\vec{u}) = \vec{v}^{(\Lambda+1)} = \vec{v}$.
- Let ${\cal J}$ be a cost funciton which explicitly depends on the
+ Let ${\cal J}$ be a cost function which explicitly depends on the
  final state $\vec{v}$ only
  (this restriction is for clarity reasons only).
  %
-Line 244 
 $ \langle \, \nabla _{v}{\cal J}^T \, ,
+Line 254 
 $ \langle \, \nabla _{v}{\cal J}^T \, ,
  = \nabla_v {\cal J} \cdot \delta \vec{v} $ )
  %
  \begin{equation}
- \begin{split}
+ \begin{aligned}
  \nabla_v {\cal J} (M(\delta \vec{u})) & = \,
  \nabla_v {\cal J} \cdot M_{\Lambda}
  \cdot ...... \cdot M_{\lambda} \cdot ...... \cdot
  M_{1} \cdot M_{0} \cdot \delta \vec{u} \\
  ~ & = \, \nabla_v {\cal J} \cdot \delta \vec{v} \\
- \end{split}
+ \end{aligned}
  \label{forward}
  \end{equation}
  %
-Line 258 
 whereas in reverse mode we have
+Line 268 
 whereas in reverse mode we have
  %
  \begin{equation}
  \boxed{
- \begin{split}
+ \begin{aligned}
  M^T ( \nabla_v {\cal J}^T) & = \,
  M_{0}^T \cdot M_{1}^T
  \cdot ...... \cdot M_{\lambda}^T \cdot ...... \cdot
-Line 267 
 M_{\Lambda}^T \cdot \nabla_v {\cal J}^T
+Line 277 
 M_{\Lambda}^T \cdot \nabla_v {\cal J}^T
  \cdot ...... \cdot
  \nabla_{v^{(\lambda)}} {\cal J}^T \\
  ~ & = \, \nabla_u {\cal J}^T
- \end{split}
+ \end{aligned}
  }
  \label{reverse}
  \end{equation}
-Line 286 
 $ \vec{v}^{(\lambda)} $ at each intermed
+Line 296 
 $ \vec{v}^{(\lambda)} $ at each intermed
  %
  \begin{equation}
  \boxed{
- \begin{split}
+ \begin{aligned}
  \nabla_{v^{(\lambda)}} {\cal J}^T |_{\vec{v}^{(\lambda)}}
  & = \,
  M_{\lambda}^T |_{\vec{v}^{(\lambda)}} \cdot ...... \cdot
  M_{\Lambda}^T |_{\vec{v}^{(\lambda)}} \cdot \delta \vec{v}^{\ast} \\
  ~ & = \, \delta \vec{v}^{(\lambda) \, \ast}
- \end{split}
+ \end{aligned}
  }
  \end{equation}
  %
-Line 301 
 We note in passing that that the $\delta
+Line 311 
 We note in passing that that the $\delta
  are the Lagrange multipliers of the model equations which determine
  $ \vec{v}^{(\lambda)}$.
- In coponents, eq. (\ref{adjoint}) reads as follows.
+ In components, eq. (\ref{adjoint}) reads as follows.
  Let
  \[
  \begin{array}{rclcrcl}
-Line 322 
 Let
+Line 332 
 Let
  \end{array}
  \]
  denote the perturbations in $\vec{u}$ and $\vec{v}$, respectively,
- and their adjoint varaiables;
+ and their adjoint variables;
  further
  \[
  M \, = \, \left(
-Line 409 
 and the shorthand notation for the adjoi
+Line 419 
 and the shorthand notation for the adjoi
  $ \delta v^{(\lambda) \, \ast}_{j} = \frac{\partial}{\partial v^{(\lambda)}_{j}}
  {\cal J}^T $, $ j = 1, \ldots , n_{\lambda} $,
  for intermediate components, yielding
+ {\small
  \begin{equation}
- \small
+ \begin{aligned}
- \begin{split}
  \left(
  \begin{array}{c}
  \delta v^{(\lambda) \, \ast}_1 \\
-Line 456 
 for intermediate components, yielding
+Line 466 
 for intermediate components, yielding
  \delta v^{\ast}_{n} \\
  \end{array}
  \right)
- \end{split}
+ \end{aligned}
  \end{equation}
+ }
  Eq. (\ref{forward}) and (\ref{reverse}) are perhaps clearest in
  showing the advantage of the reverse over the forward mode
-Line 468 
 variables $u$
+Line 479 
 variables $u$
  {\it all} intermediate states $ \vec{v}^{(\lambda)} $) are sought.
  In order to be able to solve for each component of the gradient
  $ \partial {\cal J} / \partial u_{i} $ in (\ref{forward})
- a forward calulation has to be performed for each component seperately,
+ a forward calculation has to be performed for each component separately,
  i.e. $ \delta \vec{u} = \delta u_{i} {\vec{e}_{i}} $
  for  the $i$-th forward calculation.
  Then, (\ref{forward}) represents the
-Line 487 
 M^T \left( \nabla_v {\cal J}^T \left(\de
+Line 498 
 M^T \left( \nabla_v {\cal J}^T \left(\de
  \nabla_u {\cal J}^T \cdot \delta \vec{J}
  \]
  where now $ \delta \vec{J} \in I\!\!R^l $ is a vector of
- dimenison $ l $.
+ dimension $ l $.
  In this case $ l $ reverse simulations have to be performed
  for each $ \delta J_{k}, \,\, k = 1, \ldots, l $.
  Then, the reverse mode is more efficient as long as
  $ l < n $, otherwise the forward mode is preferable.
- Stricly, the reverse mode is called adjoint mode only for
+ Strictly, the reverse mode is called adjoint mode only for
  $ l = 1 $.
  A detailed analysis of the underlying numerical operations
-Line 528 
 operator which maps the model state spac
+Line 539 
 operator which maps the model state spac
  Then, $ \nabla_v {\cal J} $ takes the form
  %
  \begin{equation*}
- \begin{split}
+ \begin{aligned}
  \nabla_v {\cal J}^T & = \, 2 \, \, H \cdot
  \left( \, {\cal H}(\vec{v}) - \vec{d} \, \right) \\
  ~          & = \, 2 \sum_{j} \left\{ \sum_k
  \frac{\partial {\cal H}_k}{\partial v_{j}}
  \left( {\cal H}_k (\vec{v}) - d_k \right)
  \right\} \, {\vec{f}_{j}} \\
- \end{split}
+ \end{aligned}
  \end{equation*}
  %
  where $H_{kj} = \partial {\cal H}_k / \partial v_{j} $ is the
-Line 557 
 Because of the local character of the de
+Line 568 
 Because of the local character of the de
  (a derivative is defined w.r.t. a point along the trajectory),
  the intermediate results of the model trajectory
  $\vec{v}^{(\lambda+1)}={\cal M}_{\lambda}(v^{(\lambda)})$
- are needed to evaluate the intermediate Jacobian
+ may be required to evaluate the intermediate Jacobian
  $M_{\lambda}|_{\vec{v}^{(\lambda)}} \, \delta \vec{v}^{(\lambda)} $.
+ This is the case e.g. for nonlinear expressions
+ (momentum advection, nonlinear equation of state), state-dependent
+ conditional statements (parameterization schemes).
  In the forward mode, the intermediate results are required
  in the same order as computed by the full forward model ${\cal M}$,
  but in the reverse mode they are required in the reverse order.
-Line 569 
 point of evaluation has to be recomputed
+Line 583 
 point of evaluation has to be recomputed
  A method to balance the amount of recomputations vs.
  storage requirements is called {\sf checkpointing}
- (e.g. \cite{res-eta:98}).
+ (e.g. \cite{gri:92}, \cite{res-eta:98}).
- It is depicted in \reffig{3levelcheck} for a 3-level checkpointing
+ It is depicted in \ref{fig:3levelcheck} for a 3-level checkpointing
  [as an example, we give explicit numbers for a 3-day
  integration with a 1-hourly timestep in square brackets].
  \begin{itemize}
-Line 580 
 In a first step, the model trajectory is
+Line 594 
 In a first step, the model trajectory is
  $ {n}^{lev3} $ subsections [$ {n}^{lev3} $=3 1-day intervals],
  with the label $lev3$ for this outermost loop.
  The model is then integrated along the full trajectory,
- 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
  [i.e. 3 times, at
  $ i = 0,1,2 $ corresponding to $ k_{i}^{lev3} = 0, 24, 48 $].
+ In addition, the cost function is computed, if needed.
  %
  \item [$lev2$]
  In a second step each subsection itself is divided into
- $ {n}^{lev2} $ sub-subsections
+ $ {n}^{lev2} $ subsections
  [$ {n}^{lev2} $=4 6-hour intervals per subsection].
  The model picks up at the last outermost dumped state
  $ v_{k_{n}^{lev3}} $ and is integrated forward in time along
  the last subsection, with the label $lev2$ for this
  intermediate loop.
- 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
  timestep
  [i.e. 4 times, at
  $ 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 615 
 $ i = 0,1,2,3 $ corresponding to $ k_{i}
  \item [$lev1$]
  Finally, the model picks up at the last intermediate dump state
  $ v_{k_{n}^{lev2}} $ and is integrated forward in time along
- the last sub-subsection, with the label $lev1$ for this
+ the last subsection, with the label $lev1$ for this
  intermediate loop.
- Within this sub-subsection only, the model state is stored
+ Within this sub-subsection only, parts of the model state is stored
- at every timestep
+ to memory at every timestep
  [i.e. every hour $ i=0,...,5$ corresponding to
  $ k_{i}^{lev1} = 66, 67, \ldots, 71 $].
- Thus, the  final state $ v_n = v_{k_{n}^{lev1}} $ is reached
+ The  final state $ v_n = v_{k_{n}^{lev1}} $ is reached
- and the model state of all peceeding timesteps along the last
+ and the model state of all preceding timesteps along the last
- sub-subsections are available, enabling integration backwards
+ innermost subsection are available, enabling integration backwards
- in time along the last sub-subsection.
+ in time along the last subsection.
- Thus, the adjoint can be computed along this last
+ The adjoint can thus be computed along this last
- sub-subsection $k_{n}^{lev2}$.
+ subsection $k_{n}^{lev2}$.
  %
  \end{itemize}
  %
  This procedure is repeated consecutively for each previous
- sub-subsection $k_{n-1}^{lev2}, \ldots, k_{1}^{lev2} $
+ subsection $k_{n-1}^{lev2}, \ldots, k_{1}^{lev2} $
  carrying the adjoint computation to the initial time
  of the subsection $k_{n}^{lev3}$.
  Then, the procedure is repeated for the previous subsection
-Line 627 
 $k_{1}^{lev3}$.
+Line 642 
 $k_{1}^{lev3}$.
  For the full model trajectory of
  $ n^{lev3} \cdot n^{lev2} \cdot n^{lev1} $ timesteps
  the required storing of the model state was significantly reduced to
- $ n^{lev1} + n^{lev2} + n^{lev3} $
+ $ n^{lev2} + n^{lev3} $ to disk and roughly $ n^{lev1} $ to memory
  [i.e. for the 3-day integration with a total oof 72 timesteps
- the model state was stored 13 times].
+ the model state was stored 7 times to disk and roughly 6 times
+ to memory].
  This saving in memory comes at a cost of a required
 full forward integrations of the model (one for each
  checkpointing level).
- The balance of storage vs. recomputation certainly depends
+ The optimal balance of storage vs. recomputation certainly depends
- on the computing resources available.
+ on the computing resources available and may be adjusted by
+ adjusting the partitioning among the
+ $ n^{lev3}, \,\, n^{lev2}, \,\, n^{lev1} $.
  \begin{figure}[t!]
- \centering
+ \begin{center}
  %\psdraft
- \psfrag{v_k1^lev3}{\mathinfigure{v_{k_{1}^{lev3}}}}
+ %\psfrag{v_k1^lev3}{\mathinfigure{v_{k_{1}^{lev3}}}}
- \psfrag{v_kn-1^lev3}{\mathinfigure{v_{k_{n-1}^{lev3}}}}
+ %\psfrag{v_kn-1^lev3}{\mathinfigure{v_{k_{n-1}^{lev3}}}}
- \psfrag{v_kn^lev3}{\mathinfigure{v_{k_{n}^{lev3}}}}
+ %\psfrag{v_kn^lev3}{\mathinfigure{v_{k_{n}^{lev3}}}}
- \psfrag{v_k1^lev2}{\mathinfigure{v_{k_{1}^{lev2}}}}
+ %\psfrag{v_k1^lev2}{\mathinfigure{v_{k_{1}^{lev2}}}}
- \psfrag{v_kn-1^lev2}{\mathinfigure{v_{k_{n-1}^{lev2}}}}
+ %\psfrag{v_kn-1^lev2}{\mathinfigure{v_{k_{n-1}^{lev2}}}}
- \psfrag{v_kn^lev2}{\mathinfigure{v_{k_{n}^{lev2}}}}
+ %\psfrag{v_kn^lev2}{\mathinfigure{v_{k_{n}^{lev2}}}}
- \psfrag{v_k1^lev1}{\mathinfigure{v_{k_{1}^{lev1}}}}
+ %\psfrag{v_k1^lev1}{\mathinfigure{v_{k_{1}^{lev1}}}}
- \psfrag{v_kn^lev1}{\mathinfigure{v_{k_{n}^{lev1}}}}
+ %\psfrag{v_kn^lev1}{\mathinfigure{v_{k_{n}^{lev1}}}}
- \mbox{\epsfig{file=part5/checkpointing.eps, width=0.8\textwidth}}
+ %\mbox{\epsfig{file=s_autodiff/figs/checkpointing.eps, width=0.8\textwidth}}
+ \resizebox{5.5in}{!}{\includegraphics{s_autodiff/figs/checkpointing.eps}}
  %\psfull
- \caption
+ \end{center}
- {Schematic view of intermediate dump and restart for
+ \caption{
+ Schematic view of intermediate dump and restart for
 -level checkpointing.}
  \label{fig:3levelcheck}
  \end{figure}
-Line 662 
 on the computing resources available.
+Line 682 
 on the computing resources available.
  % \subsection{Error covariance estimate and Hessian matrix}
  % \label{sec_hessian}
  \newpage
  %**********************************************************************
- \section{AD-specific setup by example: sensitivity of carbon sequestration}
+ \section{TLM and ADM generation in general}
- \label{sec_ad_setup_ex}
+ \label{sec_ad_setup_gen}
+ \begin{rawhtml}
+ <!-- CMIREDIR:sec_ad_setup_gen: -->
+ \end{rawhtml}
  %**********************************************************************
- The MITGCM has been adapted to enable AD using TAMC or TAF.
+ In this section we describe in a general fashion
- The present description, therefore, is specific to the
+ the parts of the code that are relevant for automatic
- use of TAMC or TAF as AD tool.
+ differentiation using the software tool TAF.
- The following sections describe the steps which are necessary to
+ Modifications to use OpenAD are described in \ref{sec_ad_openad}.
- generate a tangent linear or adjoint model of the MITGCM.
- We take as an example the sensitivity of carbon sequestration
+ \input{s_autodiff/text/doc_ad_the_model}
- in the ocean.
- The AD-relevant hooks in the code are sketched in
+ The basic flow is depicted in \ref{fig:adthemodel}.
- \reffig{adthemodel}, \reffig{adthemain}.
+ If CPP option \texttt{ALLOW\_AUTODIFF\_TAMC} is defined,
+ the driver routine
- \subsection{Overview of the experiment}
+ {\it the\_model\_main}, instead of calling {\it the\_main\_loop},
+ invokes the adjoint of this routine, {\it adthe\_main\_loop}
- We describe an adjoint sensitivity analysis of outgassing from
+ (case \texttt{\#define ALLOW\_ADJOINT\_RUN}), or
- the ocean into the atmosphere of a carbon-like tracer injected
+ the tangent linear of this routine {\it g\_the\_main\_loop}
- into the ocean interior (see \cite{hil-eta:01}).
+ (case \texttt{\#define ALLOW\_TANGENTLINEAR\_RUN}),
+ which are the toplevel routines in terms of automatic differentiation.
- \subsubsection{Passive tracer equation}
+ The routines {\it adthe\_main\_loop} or {\it g\_the\_main\_loop}
+ are generated by TAF.
- For this work the MITGCM was augmented with a thermodynamically
+ It contains both the forward integration of the full model, the
- inactive tracer, $C$. Tracer residing in the ocean
+ cost function calculation,
- model surface layer is outgassed according to a relaxation time scale,
+ any additional storing that is required for efficient checkpointing,
- $\mu$. Within the ocean interior, the tracer is passively advected
+ and the reverse integration of the adjoint model.
- by the ocean model currents. The full equation for the time evolution
- %
+ [DESCRIBE IN A SEPARATE SECTION THE WORKING OF THE TLM]
- \begin{equation}
- \label{carbon_ddt}
+ In Fig. \ref{fig:adthemodel}
- \frac{\partial C}{\partial t} \, = \,
+ the structure of {\it adthe\_main\_loop} has been strongly
- -U\cdot \nabla C \, - \, \mu C \, + \, \Gamma(C) \,+ \, S
+ simplified to focus on the essentials; in particular, no checkpointing
- \end{equation}
+ procedures are shown here.
- %
+ Prior to the call of {\it adthe\_main\_loop}, the routine
- also includes a source term $S$. This term
+ {\it ctrl\_unpack} is invoked to unpack the control vector
- represents interior sources of $C$ such as would arise due to
+ or initialise the control variables.
- direct injection.
+ Following the call of {\it adthe\_main\_loop},
- The velocity term, $U$, is the sum of the
+ the routine {\it ctrl\_pack}
- model Eulerian circulation and an eddy-induced velocity, the latter
+ is invoked to pack the control vector
- parameterized according to Gent/McWilliams
+ (cf. Section \ref{section_ctrl}).
- (\cite{gen-mcw:90, gen-eta:95}).
+ If gradient checks are to be performed, the option
- The convection function, $\Gamma$, mixes $C$ vertically wherever the
+ {\tt ALLOW\_GRADIENT\_CHECK} is defined. In this case
- fluid is locally statically unstable.
+ the driver routine {\it grdchk\_main} is called after
+ the gradient has been computed via the adjoint
- The outgassing time scale, $\mu$, in eqn. (\ref{carbon_ddt})
+ (cf. Section \ref{sec:ad_gradient_check}).
- 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
+ \subsection{General setup
- powerful tool for examining the impact of large-scale ocean circulation
+ \label{section_ad_setup}}
- on $ CO_2 $ outgassing due to interior injections.
- As source we choose a constant in time injection of
+ In order to configure AD-related setups the following packages need
- $ S = 1 \,\, {\rm mol / s}$.
+ to be enabled:
+ {\it
- \subsubsection{Model configuration}
+ \begin{table}[!ht]
+ \begin{tabular}{l}
- The model configuration employed has a constant
+ autodiff \\
- $4^\circ \times 4^\circ$ resolution horizontal grid and realistic
+ ctrl \\
- geography and bathymetry. Twenty vertical layers are used with
+ cost \\
- vertical spacing ranging
+ grdchk \\
- from 50 m near the surface to 815 m at depth.
+ \end{tabular}
- Driven to steady-state by climatalogical wind-stress, heat and
+ \end{table}
- fresh-water forcing the model reproduces well known large-scale
+ }
- features of the ocean general circulation.
+ The packages are enabled by adding them to your experiment-specific
+ configuration file
- \subsubsection{Outgassing cost function}
+ {\it packages.conf} (see Section ???).
- To quantify and understand outgassing due to injections of $C$
+ The following AD-specific CPP option files need to be customized:
- 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 ECCO\_CPPOPTIONS.h} \\
- %
+ This header file collects CPP options for the packages
- \item {\it COST\_CPPOPTIONS.h}
+ {\it autodiff, cost, ctrl} as well as AD-unrelated options for
- %
+ the external forcing package {\it exf}.
- \item {\it CPP\_EEOPTIONS.h}
+ \footnote{NOTE: These options are not set in their package-specific
- %
+ headers such as {\it COST\_CPPOPTIONS.h}, but are instead collected
- \item {\it CPP\_OPTIONS.h}
+ in the single header file {\it ECCO\_CPPOPTIONS.h}.
- %
+ The package-specific header files serve as simple
- \item {\it CTRL\_OPTIONS.h}
+ placeholders at this point.}
  %
- \item {\it ECCO\_OPTIONS.h}
+ \item {\it tamc.h} \\
- %
+ This header configures the splitting of the time stepping loop
- \item {\it SIZE.h}
+ w.r.t. the 3-level checkpointing (see section ???).
- %
- \item {\it adcommon.h}
- %
- \item {\it tamc.h}
  %
  \end{itemize}
+ %------------------------------------------------------------------
+ \subsection{Building the AD code using TAF
+ \label{section_ad_build}}
+ The build process of an AD code is very similar to building
+ the forward model. However, depending on which AD code one wishes
+ to generate, and on which AD tool is available (TAF or TAMC),
+ the following {\tt make} targets are available:
+ \begin{table}[!ht]
+ {\footnotesize
+ \begin{tabular}{|ccll|}
+ \hline
+ ~ & {\it AD-target} & {\it output} & {\it description} \\
+ \hline
+ \hline
+ (1) & {\tt <MODE><TOOL>only} & {\tt <MODE>\_<TOOL>\_output.f}  &
+ generates code for $<$MODE$>$ using $<$TOOL$>$ \\
+ ~ & ~ & ~ & no {\tt make} dependencies on {\tt .F .h} \\
+ ~ & ~ & ~ & useful for compiling on remote platforms \\
+ \hline
+ (2) & {\tt <MODE><TOOL>} & {\tt <MODE>\_<TOOL>\_output.f}  &
+ generates code for $<$MODE$>$ using $<$TOOL$>$ \\
+ ~ & ~ & ~ & includes {\tt make} dependencies on {\tt .F .h} \\
+ ~ & ~ & ~ & i.e. input for $<$TOOL$>$ may be re-generated \\
+ \hline
+ (3) & {\tt <MODE>all} & {\tt mitgcmuv\_<MODE>}  &
+ generates code for $<$MODE$>$ using $<$TOOL$>$ \\
+ ~ & ~ & ~ & and compiles all code \\
+ ~ & ~ & ~ & (use of TAF is set as default) \\
+ \hline
+ \end{tabular}
+ }
+ \end{table}
  %
- The runtime flag and parameters settings are contained in
+ Here, the following placeholders are used
- {\it verification/carbon/input/},
- together with the forcing fields and and restart files:
  %
  \begin{itemize}
  %
- \item {\it data}
+ \item $<$TOOL$>$
- %
- \item {\it data.cost}
- %
- \item {\it data.ctrl}
- %
- \item {\it data.gmredi}
  %
- \item {\it data.grdchk}
+ \begin{itemize}
- %
- \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*}
+ \item {\tt TAF}
+ \item {\tt TAMC}
  %
  \end{itemize}
  %
- Finally, the file to generate the adjoint code resides in
+ \item $<$MODE$>$
- $ adjoint/ $:
  %
  \begin{itemize}
  %
- \item {\it makefile}
+ \item {\tt ad} generates the adjoint model (ADM)
+ \item {\tt ftl} generates the tangent linear model (TLM)
+ \item {\tt svd} generates both ADM and TLM for \\
+ singular value decomposition (SVD) type calculations
  %
  \end{itemize}
  %
+ \end{itemize}
- Below we describe the customisations of this files which are
+ For example, to generate the adjoint model using TAF after routines ({\tt .F})
- specific to this experiment.
+ or headers ({\tt .h}) have been modified, but without compilation,
+ type {\tt make adtaf};
- \subsubsection{File {\it .genmakerc}}
+ or, to generate the tangent linear model using TAMC without
- This file overwrites default settings of {\it genmake}.
+ re-generating the input code, type {\tt make ftltamconly}.
- 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 ???.
- \subsubsection{File {\it ECCO\_OPTIONS.h}}
- The CPP options of several AD-related packages are grouped
- in this file:
- %
- \begin{itemize}
- %
- \item
- Adjoint support package: {\it pkg/autodiff/} \\
- This package contains hand-written adjoint code such as
- active file handling, flow directives for files which must not
- be differentiated, and TAMC-specific header files. \\
- \hspace*{4ex} {\tt \#define ALLOW\_AUTODIFF\_TAMC} \\
- defines TAMC-related features in the code. \\
- \hspace*{4ex} {\tt \#define ALLOW\_TAMC\_CHECKPOINTING} \\
- enables the checkpointing feature of TAMC
- (see Section \ref{???}).
- In the present example a 3-level checkpointing is implemented.
- The code contains the relevant store directives, common block
- and tape initialisations, storing key computation,
- 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}
- %
- \end{itemize}
- \subsubsection{File {\it SIZE.h}}
+ A typical full build process to generate the ADM via TAF would
+ look like follows:
+ \begin{verbatim}
+ % mkdir build
+ % cd build
+ % ../../../tools/genmake2 -mods=../code_ad
+ % make depend
+ % make adall
+ \end{verbatim}
- 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
+ \subsection{The AD build process in detail
- above header files of the forward code, the structure
+ \label{section_ad_build_detail}}
- 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}.
- \subsubsection{File {\it tamc.h}}
+ The {\tt make <MODE>all} target consists of the following procedures:
- This routine contains the dimensions for TAMC checkpointing.
+ \begin{enumerate}
  %
+ \item
+ A header file {\tt AD\_CONFIG.h} is generated which contains a CPP option
+ on which code ought to be generated. Depending on the {\tt make} target,
+ the contents is one of the following:
  \begin{itemize}
- %
+ \item
- \item {\tt \#ifdef ALLOW\_TAMC\_CHECKPOINTING} \\
+ {\tt \#define ALLOW\_ADJOINT\_RUN}
--level checkpointing is enabled, i.e. the timestepping
+ \item
- is divided into three different levels (see Section \ref{???}).
+ {\tt \#define ALLOW\_TANGENTLINEAR\_RUN}
- The model state of the outermost ({\tt nchklev\_3}) and the
+ \item
- intermediate ({\tt nchklev\_2}) timestepping loop are stored to file
+ {\tt \#define ALLOW\_ECCO\_OPTIMIZATION}
- (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}.
- %
  \end{itemize}
- The following parameters may be worth describing: \\
  %
- \hspace*{4ex} {\tt isbyte} \\
+ \item
- \hspace*{4ex} {\tt maxpass} \\
+ A single file {\tt <MODE>\_input\_code.f} is concatenated
- ~
+ consisting of all {\tt .f} files that are part of the list {\bf AD\_FILES}
+ and all {\tt .flow} files that are part of the list {\bf AD\_FLOW\_FILES}.
- \subsubsection{File {\it makefile}}
+ %
+ \item
- This file contains all relevant paramter flags and
+ The AD tool is invoked with the {\tt <MODE>\_<TOOL>\_FLAGS}.
- lists to run TAMC or TAF.
+ The default AD tool flags in {\tt genmake2} can be overrwritten by
- It is assumed that TAMC is available to you, either locally,
+ an {\tt adjoint\_options} file (similar to the platform-specific
- being installed on your network, or remotely through the 'TAMC Utility'.
+ {\tt build\_options}, see Section ???.
- TAMC is called with the command {\tt tamc} followed by a
+ The AD tool writes the resulting AD code into the file
- number of options. They are described in detail in the
+ {\tt <MODE>\_input\_code\_ad.f}
- TAMC manual \cite{gie:99}.
+ %
- Here we briefly discuss the main flags used in the {\it makefile}
+ \item
+ A short sed script {\tt adjoint\_sed} is applied to
+ {\tt <MODE>\_input\_code\_ad.f}
+ to reinstate {\bf myThid} into the CALL argument list of active file I/O.
+ The result is written to file {\tt <MODE>\_<TOOL>\_output.f}.
+ %
+ \item
+ All routines are compiled and an executable is generated
+ (see Table ???).
  %
- \begin{itemize}
+ \end{enumerate}
- \item [{\tt tamc}] {\tt
- -input <variable names>
+ \subsubsection{The list AD\_FILES and {\tt .list} files}
- -output <variable name> -r4 ... \\
- -toplevel <S/R name> -reverse <file names>
+ Not all routines are presented to the AD tool.
- }
+ Routines typically hidden are diagnostics routines which
- \end{itemize}
+ do not influence the cost function, but may create
+ artificial flow dependencies such as I/O of active variables.
+ {\tt genmake2} generates a list (or variable) {\bf AD\_FILES}
+ which contains all routines that are shown to the AD tool.
+ This list is put together from all files with suffix {\tt .list}
+ that {\tt genmake2} finds in its search directories.
+ The list file for the core MITgcm routines is in {\tt model/src/}
+ is called {\tt model\_ad\_diff.list}.
+ Note that no wrapper routine is shown to TAF. These are either
+ not visible at all to the AD code, or hand-written AD code
+ is available (see next section).
+ Each package directory contains its package-specific
+ list file {\tt <PKG>\_ad\_diff.list}. For example,
+ {\tt pkg/ptracers/} contains the file {\tt ptracers\_ad\_diff.list}.
+ Thus, enabling a package will automatically extend the
+ {\bf AD\_FILES} list of {\tt genmake2} to incorporate the
+ package-specific routines.
+ Note that you will need to regenerate the {\tt Makefile} if
+ you enable a package (e.g. by adding it to {\tt packages.conf})
+ and a {\tt Makefile} already exists.
+ \subsubsection{The list AD\_FLOW\_FILES and {\tt .flow} files}
+ TAMC and TAF can evaluate user-specified directives
+ that start with a specific syntax ({\tt CADJ}, {\tt C\$TAF}, {\tt !\$TAF}).
+ The main categories of directives are STORE directives and
+ FLOW directives. Here, we are concerned with flow directives,
+ store directives are treated elsewhere.
+ Flow directives enable the AD tool to evaluate how it should treat
+ routines that are 'hidden' by the user, i.e. routines which are
+ not contained in the {\bf AD\_FILES} list (see previous section),
+ but which are called in part of the code that the AD tool does see.
+ The flow directive tell the AD tool
  %
  \begin{itemize}
  %
- \item {\tt -toplevel <S/R name>} \\
+ \item which subroutine arguments are input/output
- Name of the toplevel routine, with respect to which the
+ \item which subroutine arguments are active
- control flow analysis is performed.
+ \item which subroutine arguments are required to compute the cost
- %
+ \item which subroutine arguments are dependent
- \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.
- %
- \item {\tt -r4} \\
- ~
  %
  \end{itemize}
+ %
+ The syntax for the flow directives can be found in the
+ AD tool manuals.
+ {\tt genmake2} generates a list (or variable) {\bf AD\_FLOW\_FILES}
+ which contains all files with suffix{\tt .flow} that it finds
+ in its search directories.
+ The flow directives for the core MITgcm routines of
+ {\tt eesupp/src/} and {\tt model/src/}
+ reside in {\tt pkg/autodiff/}.
+ This directory also contains hand-written adjoint code
+ for the MITgcm WRAPPER (section \ref{chap:sarch}).
+ Flow directives for package-specific routines are contained in
+ the corresponding package directories in the file
+ {\tt <PKG>\_ad.flow}, e.g. ptracers-specific directives are in
+ {\tt ptracers\_ad.flow}.
+ \subsubsection{Store directives for 3-level checkpointing}
+ The storing that is required at each period of the
+-level checkpointing is controled by three
+ top-level headers.
- \subsubsection{File {\it data}}
+ \begin{verbatim}
+ do ilev_3 = 1, nchklev_3
- \subsubsection{File {\it data.cost}}
+ #  include ``checkpoint_lev3.h''
+    do ilev_2 = 1, nchklev_2
- \subsubsection{File {\it data.ctrl}}
+ #     include ``checkpoint_lev2.h''
+       do ilev_1 = 1, nchklev_1
- \subsubsection{File {\it data.gmredi}}
+ #        include ``checkpoint_lev1.h''
- \subsubsection{File {\it data.grdchk}}
+ ...
- \subsubsection{File {\it data.optim}}
+       end do
+    end do
- \subsubsection{File {\it data.pkg}}
+ end do
+ \end{verbatim}
- \subsubsection{File {\it eedata}}
- \subsubsection{File {\it topog.bin}}
- \subsubsection{File {\it windx.bin, windy.bin}}
- \subsubsection{File {\it salt.bin, theta.bin}}
+ All files {\tt checkpoint\_lev?.h} are contained in directory
+ {\tt pkg/autodiff/}.
- \subsubsection{File {\it SSS.bin, SST.bin}}
- \subsubsection{File {\it pickup*}}
+ \subsubsection{Changing the default AD tool flags: ad\_options files}
- \subsection{Compiling the model and its adjoint}
- \newpage
+ \subsubsection{Hand-written adjoint code}
- %**********************************************************************
+ %------------------------------------------------------------------
- \section{TLM and ADM generation in general}
- \label{sec_ad_setup_gen}
- %**********************************************************************
- In this section we describe in a general fashion
- the parts of the code that are relevant for automatic
- differentiation using the software tool TAMC.
- \begin{figure}[b!]
- \input{part5/doc_ad_the_model}
- \caption{~}
- \label{fig:adthemodel}
- \end{figure}
- The basic flow is depicted in \reffig{adthemodel}.
- If the option {\tt ALLOW\_AUTODIFF\_TAMC} is defined, the driver routine
- {\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}).
  \subsection{The cost function (dependent variable)
  \label{section_cost}}
-Line 1165 
 the gradient has been computed via the a
+Line 1007 
 the gradient has been computed via the a
  The cost function $ {\cal J} $ is referred to as the {\sf dependent variable}.
  It is a function of the input variables $ \vec{u} $ via the composition
  $ {\cal J}(\vec{u}) \, = \, {\cal J}(M(\vec{u})) $.
- The input is referred to as the
+ The input are referred to as the
  {\sf independent variables} or {\sf control variables}.
  All aspects relevant to the treatment of the cost function $ {\cal J} $
- (parameter setting, initialisation, accumulation,
+ (parameter setting, initialization, accumulation,
  final evaluation), are controlled by the package {\it pkg/cost}.
+ The aspects relevant to the treatment of the independent variables
+ are controlled by the package {\it pkg/ctrl} and will be treated
+ in the next section.
- \begin{figure}[h!]
+ \input{s_autodiff/text/doc_cost_flow}
- \input{part5/doc_cost_flow}
- \caption{~}
+ \subsubsection{Enabling the package}
- \label{fig:costflow}
- \end{figure}
- \subsubsection{genmake and CPP options}
- %
- \begin{itemize}
- %
- \item
  \fbox{
  \begin{minipage}{12cm}
- {\it genmake}, {\it CPP\_OPTIONS.h}, {\it ECCO\_CPPOPTIONS.h}
+ {\it packages.conf}, {\it ECCO\_CPPOPTIONS.h}
  \end{minipage}
  }
- \end{itemize}
+ \begin{itemize}
  %
- The directory {\it pkg/cost} can be included to the
+ \item
- compile list in 3 different ways (cf. Section \ref{???}):
+ The package is enabled by adding {\it cost} to your file {\it packages.conf}
+ (see Section ???)
  %
- \begin{enumerate}
+ \item
- %
- \item {\it genmake}: \\
- Change the default settings in the file {\it genmake} by adding
+ \end{itemize}
- {\bf cost} to the {\bf enable} list (not recommended).
- %
- \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}.
  %
- \end{enumerate}
+ N.B.: In general the following packages ought to be enabled
+ simultaneously: {\it autodiff, cost, ctrl}.
  The basic CPP option to enable the cost function is {\bf ALLOW\_COST}.
  Each specific cost function contribution has its own option.
  For the present example the option is {\bf ALLOW\_COST\_TRACER}.
  All cost-specific options are set in {\it ECCO\_CPPOPTIONS.h}
  Since the cost function is usually used in conjunction with
  automatic differentiation, the CPP option
- {\bf ALLOW\_ADJOINT\_RUN} should be defined
+ {\bf ALLOW\_ADJOINT\_RUN} (file {\it CPP\_OPTIONS.h}) and
- (file {\it CPP\_OPTIONS.h}).
+ {\bf ALLOW\_AUTODIFF\_TAMC} (file {\it ECCO\_CPPOPTIONS.h})
+ should be defined.
- \subsubsection{Initialisation}
+ \subsubsection{Initialization}
  %
- The initialisation of the {\it cost} package is readily enabled
+ The initialization of the {\it cost} package is readily enabled
- as soon as the CPP option {\bf ALLOW\_ADJOINT\_RUN} is defined.
+ as soon as the CPP option {\bf ALLOW\_COST} is defined.
  %
  \begin{itemize}
  %
-Line 1250 
 Variables: {\it cost\_init}
+Line 1082 
 Variables: {\it cost\_init}
  }
  \\
  This S/R
- initialises the different cost function contributions.
+ initializes the different cost function contributions.
- The contribtion for the present example is {\bf objf\_tracer}
+ The contribution for the present example is {\bf objf\_tracer}
  which is defined on each tile (bi,bj).
  %
  \end{itemize}
-Line 1274 
 Within this 'driver' routine, S/R are ca
+Line 1106 
 Within this 'driver' routine, S/R are ca
  the chosen cost function contributions.
  In the present example ({\bf ALLOW\_COST\_TRACER}),
  S/R {\it cost\_tracer} is called.
- It accumulates {\bf objf\_tracer} according to eqn. (\ref{???}).
+ It accumulates {\bf objf\_tracer} according to eqn. (\ref{ask_the_author:doc_ad_2}).
  %
  \subsubsection{Finalize all contributions}
  %
-Line 1294 
 from each contribution and sums over all
+Line 1126 
 from each contribution and sums over all
  \begin{equation}
  {\cal J} \, = \,
  {\rm fc} \, = \,
- {\rm mult\_tracer} \sum_{bi,\,bj}^{nSx,\,nSy}
+ {\rm mult\_tracer} \sum_{\text{global sum}} \sum_{bi,\,bj}^{nSx,\,nSy}
  {\rm objf\_tracer}(bi,bj) \, + \, ...
  \end{equation}
  %
  The total cost function {\bf fc} will be the
- 'dependent' variable in the argument list for TAMC, i.e.
+ 'dependent' variable in the argument list for TAF, i.e.
  \begin{verbatim}
- tamc -output 'fc' ...
+ taf -output 'fc' ...
  \end{verbatim}
  %%%% \end{document}
- \begin{figure}
+ \input{s_autodiff/text/doc_ad_the_main}
- \input{part5/doc_ad_the_main}
- \caption{~}
- \label{fig:adthemain}
- \end{figure}
  \subsection{The control variables (independent variables)
  \label{section_ctrl}}
-Line 1327 
 as variable assignments. Therefore, file
+Line 1155 
 as variable assignments. Therefore, file
  active variables are written and from which active variables
  are read are called {\sf active files}.
  All aspects relevant to the treatment of the control variables
- (parameter setting, initialisation, perturbation)
+ (parameter setting, initialization, perturbation)
- are controled by the package {\it pkg/ctrl}.
+ are controlled by the package {\it pkg/ctrl}.
- \begin{figure}[h!]
+ \input{s_autodiff/text/doc_ctrl_flow}
- \input{part5/doc_ctrl_flow}
- \caption{~}
- \label{fig:ctrlflow}
- \end{figure}
  \subsubsection{genmake and CPP options}
  %
-Line 1350 
 are controled by the package {\it pkg/ct
+Line 1174 
 are controled by the package {\it pkg/ct
  %
  To enable the directory to be included to the compile list,
  {\bf ctrl} has to be added to the {\bf enable} list in
- {\it .genmakerc} (or {\it genmake} itself).
+ {\it .genmakerc} or in {\it genmake} itself (analogous to {\it cost}
+ package, cf. previous section).
  Each control variable is enabled via its own CPP option
  in {\it ECCO\_CPPOPTIONS.h}.
- \subsubsection{Initialisation}
+ \subsubsection{Initialization}
  %
  \begin{itemize}
  %
-Line 1391 
 and their gradients: {\it ctrl\_unpack}
+Line 1216 
 and their gradients: {\it ctrl\_unpack}
  \\
  %
  Two important issues related to the handling of the control
- variables in the MITGCM need to be addressed.
+ variables in MITgcm need to be addressed.
  First, in order to save memory, the control variable arrays
  are not kept in memory, but rather read from file and added
- to the initial fields during the model initialisation phase.
+ to the initial fields during the model initialization phase.
  Similarly, the corresponding adjoint fields which represent
  the gradient of the cost function w.r.t. the control variables
  are written to file at the end of the adjoint integration.
-Line 1423 
 and gradient are generated and initialis
+Line 1248 
 and gradient are generated and initialis
  %
  The dependency flow for differentiation w.r.t. the controls
  starts with adding a perturbation onto the input variable,
- thus defining the independent or control variables for TAMC.
+ thus defining the independent or control variables for TAF.
  Three types of controls may be considered:
  %
  \begin{itemize}
-Line 1444 
 temperature and salinity are initialised
+Line 1269 
 temperature and salinity are initialised
  a perturbation anomaly is added to the field in S/R
  {\it ctrl\_map\_ini}
  %
+ %\begin{eqnarray}
  \begin{equation}
- \begin{split}
+ \begin{aligned}
  u         & = \, u_{[0]} \, + \, \Delta u \\
  {\bf tr1}(...) & = \, {\bf tr1_{ini}}(...) \, + \, {\bf xx\_tr1}(...)
  \label{perturb}
- \end{split}
+ \end{aligned}
  \end{equation}
+ %\end{eqnarray}
  %
  {\bf xx\_tr1} is a 3-dim. global array
  holding the perturbation. In the case of a simple
  sensitivity study this array is identical to zero.
  However, it's specification is essential in the context
- of automatic differentiation since TAMC
+ of automatic differentiation since TAF
  treats the corresponding line in the code symbolically
  when determining the differentiation chain and its origin.
  Thus, the variable names are part of the argument list
- when calling TAMC:
+ when calling TAF:
  %
  \begin{verbatim}
- tamc -input 'xx_tr1 ...' ...
+ taf -input 'xx_tr1 ...' ...
  \end{verbatim}
  %
- Now, as mentioned above, the MITGCM avoids maintaining
+ Now, as mentioned above, MITgcm avoids maintaining
  an array for each control variable by reading the
  perturbation to a temporary array from file.
- To ensure the symbolic link to be recognized by TAMC, a scalar
+ To ensure the symbolic link to be recognized by TAF, a scalar
  dummy variable {\bf xx\_tr1\_dummy} is introduced
  and an 'active read' routine of the adjoint support
  package {\it pkg/autodiff} is invoked.
  The read-procedure is tagged with the variable
- {\bf xx\_tr1\_dummy} enabbling TAMC to recognize the
+ {\bf xx\_tr1\_dummy} enabling TAF to recognize the
- initialisation of the perturbation.
+ initialization of the perturbation.
- The modified call of TAMC thus reads
+ The modified call of TAF thus reads
  %
  \begin{verbatim}
- tamc -input 'xx_tr1_dummy ...' ...
+ taf -input 'xx_tr1_dummy ...' ...
  \end{verbatim}
  %
  and the modified operation to (\ref{perturb})
-Line 1494 
 in the code takes on the form
+Line 1321 
 in the code takes on the form
  %
  Note, that reading an active variable corresponds
  to a variable assignment. Its derivative corresponds
- to a write statement of the adjoint variable.
+ to a write statement of the adjoint variable, followed by
+ a reset.
  The 'active file' routines have been designed
  to support active read and corresponding adjoint active write
  operations (and vice versa).
-Line 1586 
 variables are written to {\bf adxx\_ ...
+Line 1414 
 variables are written to {\bf adxx\_ ...
  \begin{itemize}
  %
  \item {\bf vector\_ctrl}: the control vector \\
- At the very beginning of the model initialisation,
+ At the very beginning of the model initialization,
  the updated compressed control vector is read (or initialised)
  and distributed to 2-dim. and 3-dim. control variable fields.
  %
-Line 1611 
 at intermediate times can be written usi
+Line 1439 
 at intermediate times can be written usi
  {\it addummy\_in\_stepping}.
  This routine is part of the adjoint support package
  {\it pkg/autodiff} (cf.f. below).
+ The procedure is enabled using via the CPP-option
+ {\bf ALLOW\_AUTODIFF\_MONITOR} (file {\it ECCO\_CPPOPTIONS.h}).
  To be part of the adjoint code, the corresponding S/R
  {\it dummy\_in\_stepping} has to be called in the forward
  model (S/R {\it the\_main\_loop}) at the appropriate place.
+ The adjoint common blocks are extracted from the adjoint code
+ via the header file {\it adcommon.h}.
  {\it dummy\_in\_stepping} is essentially empty,
  the corresponding adjoint routine is hand-written rather
-Line 1640 
 the common blocks
+Line 1472 
 the common blocks
  {\bf /adtr1\_r/}, {\bf /adffields/},
  which have been extracted from the adjoint code to enable
  access to the adjoint variables.
+ {\bf WARNING:} If the structure of the common blocks
+ {\bf /dynvars\_r/}, {\bf /dynvars\_cd/}, etc., changes
+ similar changes will occur in the adjoint common blocks.
+ Therefore, consistency between the TAMC-generated common blocks
+ and those in {\it adcommon.h} have to be checked.
  %
  \end{itemize}
-Line 1654 
 The gradient $ \nabla _{u}{\cal J} |_{u_
+Line 1492 
 The gradient $ \nabla _{u}{\cal J} |_{u_
  with the value of the cost function itself $ {\cal J}(u_{[k]}) $
  at iteration step $ k $ serve
  as input to a minimization routine (e.g. quasi-Newton method,
- conjugate gradient, ... \cite{gil_lem:89})
+ conjugate gradient, ... \cite{gil-lem:89})
  to compute an update in the
  control variable for iteration step $k+1$
  \[
-Line 1665 
 u_{[k+1]} \, = \,  u_{[0]} \, + \, \Delt
+Line 1503 
 u_{[k+1]} \, = \,  u_{[0]} \, + \, \Delt
  $ u_{[k+1]} $ then serves as input for a forward/adjoint run
  to determine $ {\cal J} $ and $ \nabla _{u}{\cal J} $ at iteration step
  $ k+1 $.
- Tab. \ref{???} sketches the flow between forward/adjoint model
+ Tab. \ref{ask_the_author:doc_ad_2} sketches the flow between forward/adjoint model
  and the minimization routine.
+ {\scriptsize
  \begin{eqnarray*}
- \scriptsize
  \begin{array}{ccccc}
  u_{[0]} \,\, ,  \,\, \Delta u_{[k]}    & ~ & ~ & ~ & ~ \\
  {\Big\downarrow}
-Line 1720 
 ad \, v_{[k]} (\delta {\cal J}) =
+Line 1558 
 ad \, v_{[k]} (\delta {\cal J}) =
   ~ & ~ & ~ & ~ & \Delta u_{[k+1]} \\
  \end{array}
  \end{eqnarray*}
+ }
  The routines {\it ctrl\_unpack} and {\it ctrl\_pack} provide
  the link between the model and the minimization routine.
- As described in Section \ref{???}
+ As described in Section \ref{ask_the_author:doc_ad_2}
  the {\it unpack} and {\it pack} routines read and write
  control and gradient {\it vectors} which are compressed
  to contain only wet points, in addition to the full
-Line 1785 
 to {\it adxx\_...$<$k$>$}, again via the
+Line 1624 
 to {\it adxx\_...$<$k$>$}, again via the
  Finally, {\it ctrl\_pack} collects all adjoint files
  and writes them to the compressed vector file
  {\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.4
 


changed lines


 
Added in v.1.23
 Legend:



Removed from v.1.4
 


changed lines


 
Added in v.1.23
-Removed from v.1.4
+Added in v.1.23

	ViewVC Help
Powered by ViewVC 1.1.22