/[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.1.1.1 by adcroft, Wed Aug 8 16:16:26 2001 UTC revision 1.7 by cnh, Thu Oct 25 18:36:55 2001 UTC
# Line 18  In principle, a variety of derived algor Line 18  In principle, a variety of derived algor
18  can be generated automatically in this way.  can be generated automatically in this way.
19    
20  The MITGCM has been adapted for use with the  The MITGCM has been adapted for use with the
21  Tangent linear and Adjoint Model Compiler (TAMC) and its succssor TAF  Tangent linear and Adjoint Model Compiler (TAMC) and its successor TAF
22  (Transformation of Algorithms in Fortran), developed  (Transformation of Algorithms in Fortran), developed
23  by Ralf Giering (\cite{gie-kam:98}, \cite{gie:99,gie:00}).  by Ralf Giering (\cite{gie-kam:98}, \cite{gie:99,gie:00}).
24  The first application of the adjoint of the MITGCM for senistivity  The first application of the adjoint of the MITGCM for sensitivity
25  studies has been published by \cite{maro-eta:99}.  studies has been published by \cite{maro-eta:99}.
26  \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
27  for ocean state estimation studies.  for ocean state estimation studies.
28    In the following we shall refer to TAMC and TAF synonymously,
29    except were explicitly stated otherwise.
30    
31  TAMC exploits the chain rule for computing the first  TAMC exploits the chain rule for computing the first
32  derivative of a function with  derivative of a function with
33  respect to a set of input variables.  respect to a set of input variables.
34  Treating a given forward code as a composition of operations --  Treating a given forward code as a composition of operations --
35  each line representing a compositional element -- the chain rule is  each line representing a compositional element, the chain rule is
36  rigorously applied to the code, line by line. The resulting  rigorously applied to the code, line by line. The resulting
37  tangent linear or adjoint code,  tangent linear or adjoint code,
38  then, may be thought of as the composition in  then, may be thought of as the composition in
39  forward or reverse order, respectively, of the  forward or reverse order, respectively, of the
40  Jacobian matrices of the forward code compositional elements.  Jacobian matrices of the forward code's compositional elements.
41    
42  %**********************************************************************  %**********************************************************************
43  \section{Some basic algebra}  \section{Some basic algebra}
# Line 50  $\vec{u}=(u_1,\ldots,u_m)$ Line 52  $\vec{u}=(u_1,\ldots,u_m)$
52  such as forcing functions) to the $n$-dimensional space  such as forcing functions) to the $n$-dimensional space
53  $V \subset I\!\!R^n$ of  $V \subset I\!\!R^n$ of
54  model output variable $\vec{v}=(v_1,\ldots,v_n)$  model output variable $\vec{v}=(v_1,\ldots,v_n)$
55  (model state, model diagnostcs, objective function, ...)  (model state, model diagnostics, objective function, ...)
56  under consideration,  under consideration,
57  %  %
58  \begin{equation}  \begin{equation}
# Line 105  In contrast to the full nonlinear model Line 107  In contrast to the full nonlinear model
107  $ M $ is just a matrix  $ M $ is just a matrix
108  which can readily be used to find the forward sensitivity of $\vec{v}$ to  which can readily be used to find the forward sensitivity of $\vec{v}$ to
109  perturbations in  $u$,  perturbations in  $u$,
110  but if there are very many input variables $(>>O(10^{6})$ for  but if there are very many input variables $(\gg O(10^{6})$ for
111  large-scale oceanographic application), it quickly becomes  large-scale oceanographic application), it quickly becomes
112  prohibitive to proceed directly as in (\ref{tangent_linear}),  prohibitive to proceed directly as in (\ref{tangent_linear}),
113  if the impact of each component $ {\bf e_{i}} $ is to be assessed.  if the impact of each component $ {\bf e_{i}} $ is to be assessed.
# Line 130  or a measure of some model-to-data misfi Line 132  or a measure of some model-to-data misfi
132  \label{compo}  \label{compo}
133  \end{eqnarray}  \end{eqnarray}
134  %  %
135  The linear approximation of $ {\cal J} $,  The perturbation of $ {\cal J} $ around a fixed point $ {\cal J}_0 $,
136  \[  \[
137  {\cal J} \, \approx \, {\cal J}_0 \, + \, \delta {\cal J}  {\cal J} \, = \, {\cal J}_0 \, + \, \delta {\cal J}
138  \]  \]
139  can be expressed in both bases of $ \vec{u} $ and $ \vec{v} $  can be expressed in both bases of $ \vec{u} $ and $ \vec{v} $
140  w.r.t. their corresponding inner product  w.r.t. their corresponding inner product
# Line 152  $\left\langle \,\, , \,\, \right\rangle Line 154  $\left\langle \,\, , \,\, \right\rangle
154  \label{deljidentity}  \label{deljidentity}
155  \end{equation}  \end{equation}
156  %  %
157  (note, that the gradient $ \nabla f $ is a pseudo-vector, therefore  (note, that the gradient $ \nabla f $ is a co-vector, therefore
158  its transpose is required in the above inner product).  its transpose is required in the above inner product).
159  Then, using the representation of  Then, using the representation of
160  $ \delta {\cal J} =  $ \delta {\cal J} =
# Line 168  transpose of $ A $, Line 170  transpose of $ A $,
170  \[  \[
171  A^{\ast} \, = \, A^T  A^{\ast} \, = \, A^T
172  \]  \]
173  and from eq. (\ref{tangent_linear}), we note that  and from eq. (\ref{tangent_linear}), (\ref{deljidentity}),
174    we note that
175  (omitting $|$'s):  (omitting $|$'s):
176  %  %
177  \begin{equation}  \begin{equation}
# Line 204  the adjoint variable of the model state Line 207  the adjoint variable of the model state
207  $ \delta \vec{u}^{\ast} $ the adjoint variable of the control variable $ \vec{u} $.  $ \delta \vec{u}^{\ast} $ the adjoint variable of the control variable $ \vec{u} $.
208    
209  The {\sf reverse} nature of the adjoint calculation can be readily  The {\sf reverse} nature of the adjoint calculation can be readily
210  seen as follows. Let us decompose ${\cal J}(u)$, thus:  seen as follows.
211    Consider a model integration which consists of $ \Lambda $
212    consecutive operations
213    $ {\cal M}_{\Lambda} (  {\cal M}_{\Lambda-1} (
214    ...... ( {\cal M}_{\lambda} (
215    ......
216    ( {\cal M}_{1} ( {\cal M}_{0}(\vec{u}) )))) $,
217    where the ${\cal M}$'s could be the elementary steps, i.e. single lines
218    in the code of the model, or successive time steps of the
219    model integration,
220    starting at step 0 and moving up to step $\Lambda$, with intermediate
221    ${\cal M}_{\lambda} (\vec{u}) = \vec{v}^{(\lambda+1)}$ and final
222    ${\cal M}_{\Lambda} (\vec{u}) = \vec{v}^{(\Lambda+1)} = \vec{v}$.
223    Let ${\cal J}$ be a cost function which explicitly depends on the
224    final state $\vec{v}$ only
225    (this restriction is for clarity reasons only).
226    %
227    ${\cal J}(u)$ may be decomposed according to:
228  %  %
229  \begin{equation}  \begin{equation}
230  {\cal J}({\cal M}(\vec{u})) \, = \,  {\cal J}({\cal M}(\vec{u})) \, = \,
# Line 215  seen as follows. Let us decompose ${\cal Line 235  seen as follows. Let us decompose ${\cal
235  \label{compos}  \label{compos}
236  \end{equation}  \end{equation}
237  %  %
238  where the ${\cal M}$'s could be the elementary steps, i.e. single lines  Then, according to the chain rule, the forward calculation reads,
239  in the code of the model,  in terms of the Jacobi matrices
 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}$  
 Then, according to the chain rule the forward calculation reads in  
 terms of the Jacobi matrices  
240  (we've omitted the $ | $'s which, nevertheless are important  (we've omitted the $ | $'s which, nevertheless are important
241  to the aspect of {\it tangent} linearity;  to the aspect of {\it tangent} linearity;
242  note also that per definition  note also that by definition
243  $ \langle \, \nabla _{v}{\cal J}^T \, , \, \delta \vec{v} \, \rangle  $ \langle \, \nabla _{v}{\cal J}^T \, , \, \delta \vec{v} \, \rangle
244  = \nabla_v {\cal J} \cdot \delta \vec{v} $ )  = \nabla_v {\cal J} \cdot \delta \vec{v} $ )
245  %  %
# Line 259  M_{\Lambda}^T \cdot \nabla_v {\cal J}^T Line 274  M_{\Lambda}^T \cdot \nabla_v {\cal J}^T
274  %  %
275  clearly expressing the reverse nature of the calculation.  clearly expressing the reverse nature of the calculation.
276  Eq. (\ref{reverse}) is at the heart of automatic adjoint compilers.  Eq. (\ref{reverse}) is at the heart of automatic adjoint compilers.
277  The intermediate steps $\lambda$ in  If the intermediate steps $\lambda$ in
278  eqn. (\ref{compos}) -- (\ref{reverse})  eqn. (\ref{compos}) -- (\ref{reverse})
279  could represent the model state (forward or adjoint) at each  represent the model state (forward or adjoint) at each
280  intermediate time step in which case  intermediate time step as noted above, then correspondingly,
281  $ {\cal M}(\vec{v}^{(\lambda)}) = \vec{v}^{(\lambda+1)} $, and correspondingly,  $ M^T (\delta \vec{v}^{(\lambda) \, \ast}) =
282  $ M^T (\delta \vec{v}^{(\lambda) \, \ast}) = \delta \vec{v}^{(\lambda-1) \, \ast} $,  \delta \vec{v}^{(\lambda-1) \, \ast} $ for the adjoint variables.
283  but they can also be viewed more generally as  It thus becomes evident that the adjoint calculation also
284  single lines of code in the numerical algorithm.  yields the adjoint of each model state component
285  In both cases it becomes evident that the adjoint calculation  $ \vec{v}^{(\lambda)} $ at each intermediate step $ \lambda $, namely
 yields at the same time the adjoint of each model state component  
 $ \vec{v}^{(\lambda)} $ at each intermediate step $ l $, namely  
286  %  %
287  \begin{equation}  \begin{equation}
288  \boxed{  \boxed{
# Line 285  M_{\Lambda}^T |_{\vec{v}^{(\lambda)}} \c Line 298  M_{\Lambda}^T |_{\vec{v}^{(\lambda)}} \c
298  %  %
299  in close analogy to eq. (\ref{adjoint})  in close analogy to eq. (\ref{adjoint})
300  We note in passing that that the $\delta \vec{v}^{(\lambda) \, \ast}$  We note in passing that that the $\delta \vec{v}^{(\lambda) \, \ast}$
301  are the Lagrange multipliers of the model state $ \vec{v}^{(\lambda)}$.  are the Lagrange multipliers of the model equations which determine
302    $ \vec{v}^{(\lambda)}$.
303    
304  In coponents, eq. (\ref{adjoint}) reads as follows.  In components, eq. (\ref{adjoint}) reads as follows.
305  Let  Let
306  \[  \[
307  \begin{array}{rclcrcl}  \begin{array}{rclcrcl}
# Line 308  Let Line 322  Let
322  \end{array}  \end{array}
323  \]  \]
324  denote the perturbations in $\vec{u}$ and $\vec{v}$, respectively,  denote the perturbations in $\vec{u}$ and $\vec{v}$, respectively,
325  and their adjoint varaiables;  and their adjoint variables;
326  further  further
327  \[  \[
328  M \, = \, \left(  M \, = \, \left(
# Line 395  and the shorthand notation for the adjoi Line 409  and the shorthand notation for the adjoi
409  $ \delta v^{(\lambda) \, \ast}_{j} = \frac{\partial}{\partial v^{(\lambda)}_{j}}  $ \delta v^{(\lambda) \, \ast}_{j} = \frac{\partial}{\partial v^{(\lambda)}_{j}}
410  {\cal J}^T $, $ j = 1, \ldots , n_{\lambda} $,  {\cal J}^T $, $ j = 1, \ldots , n_{\lambda} $,
411  for intermediate components, yielding  for intermediate components, yielding
412  \[  \begin{equation}
413  \footnotesize  \small
414    \begin{split}
415  \left(  \left(
416  \begin{array}{c}  \begin{array}{c}
417  \delta v^{(\lambda) \, \ast}_1 \\  \delta v^{(\lambda) \, \ast}_1 \\
# Line 404  for intermediate components, yielding Line 419  for intermediate components, yielding
419  \delta v^{(\lambda) \, \ast}_{n_{\lambda}} \\  \delta v^{(\lambda) \, \ast}_{n_{\lambda}} \\
420  \end{array}  \end{array}
421  \right)  \right)
422  \, = \,  \, = &
423  \left(  \left(
424  \begin{array}{ccc}  \begin{array}{ccc}
425  \frac{\partial ({\cal M}_{\lambda})_1}{\partial v^{(\lambda)}_1}  \frac{\partial ({\cal M}_{\lambda})_1}{\partial v^{(\lambda)}_1}
426  & \ldots &  & \ldots \,\, \ldots &
427  \frac{\partial ({\cal M}_{\lambda})_{n_{\lambda+1}}}{\partial v^{(\lambda)}_1} \\  \frac{\partial ({\cal M}_{\lambda})_{n_{\lambda+1}}}{\partial v^{(\lambda)}_1} \\
428  \vdots & ~ & \vdots \\  \vdots & ~ & \vdots \\
429  \frac{\partial ({\cal M}_{\lambda})_1}{\partial v^{(\lambda)}_{n_{\lambda}}}  \frac{\partial ({\cal M}_{\lambda})_1}{\partial v^{(\lambda)}_{n_{\lambda}}}
430  & \ldots  &  & \ldots \,\, \ldots  &
431  \frac{\partial ({\cal M}_{\lambda})_{n_{\lambda+1}}}{\partial v^{(\lambda)}_{n_{\lambda}}} \\  \frac{\partial ({\cal M}_{\lambda})_{n_{\lambda+1}}}{\partial v^{(\lambda)}_{n_{\lambda}}} \\
432  \end{array}  \end{array}
433  \right)  \right)
 %  
434  \cdot  \cdot
435  %  %
436    \\ ~ & ~
437    \\ ~ &
438    %
439  \left(  \left(
440  \begin{array}{ccc}  \begin{array}{ccc}
441  \frac{\partial ({\cal M}_{\lambda+1})_1}{\partial v^{(\lambda+1)}_1}  \frac{\partial ({\cal M}_{\lambda+1})_1}{\partial v^{(\lambda+1)}_1}
# Line 431  for intermediate components, yielding Line 448  for intermediate components, yielding
448  \frac{\partial ({\cal M}_{\lambda+1})_{n_{\lambda+2}}}{\partial v^{(\lambda+1)}_{n_{\lambda+1}}} \\  \frac{\partial ({\cal M}_{\lambda+1})_{n_{\lambda+2}}}{\partial v^{(\lambda+1)}_{n_{\lambda+1}}} \\
449  \end{array}  \end{array}
450  \right)  \right)
451  \cdot \ldots \ldots \cdot  \cdot \, \ldots \, \cdot
452  \left(  \left(
453  \begin{array}{c}  \begin{array}{c}
454  \delta v^{\ast}_1 \\  \delta v^{\ast}_1 \\
# Line 439  for intermediate components, yielding Line 456  for intermediate components, yielding
456  \delta v^{\ast}_{n} \\  \delta v^{\ast}_{n} \\
457  \end{array}  \end{array}
458  \right)  \right)
459  \]  \end{split}
460    \end{equation}
461    
462  Eq. (\ref{forward}) and (\ref{reverse}) are perhaps clearest in  Eq. (\ref{forward}) and (\ref{reverse}) are perhaps clearest in
463  showing the advantage of the reverse over the forward mode  showing the advantage of the reverse over the forward mode
# Line 450  variables $u$ Line 468  variables $u$
468  {\it all} intermediate states $ \vec{v}^{(\lambda)} $) are sought.  {\it all} intermediate states $ \vec{v}^{(\lambda)} $) are sought.
469  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
470  $ \partial {\cal J} / \partial u_{i} $ in (\ref{forward})  $ \partial {\cal J} / \partial u_{i} $ in (\ref{forward})
471  a forward calulation has to be performed for each component seperately,  a forward calculation has to be performed for each component separately,
472  i.e. $ \delta \vec{u} = \delta u_{i} {\vec{e}_{i}} $  i.e. $ \delta \vec{u} = \delta u_{i} {\vec{e}_{i}} $
473  for  the $i$-th forward calculation.  for  the $i$-th forward calculation.
474  Then, (\ref{forward}) represents the  Then, (\ref{forward}) represents the
# Line 460  In contrast, eq. (\ref{reverse}) yields Line 478  In contrast, eq. (\ref{reverse}) yields
478  gradient $\nabla _{u}{\cal J}$ (and all intermediate gradients  gradient $\nabla _{u}{\cal J}$ (and all intermediate gradients
479  $\nabla _{v^{(\lambda)}}{\cal J}$) within a single reverse calculation.  $\nabla _{v^{(\lambda)}}{\cal J}$) within a single reverse calculation.
480    
481  Note, that in case $ {\cal J} $ is a vector-valued function  Note, that if $ {\cal J} $ is a vector-valued function
482  of dimension $ l > 1 $,  of dimension $ l > 1 $,
483  eq. (\ref{reverse}) has to be modified according to  eq. (\ref{reverse}) has to be modified according to
484  \[  \[
# Line 468  M^T \left( \nabla_v {\cal J}^T \left(\de Line 486  M^T \left( \nabla_v {\cal J}^T \left(\de
486  \, = \,  \, = \,
487  \nabla_u {\cal J}^T \cdot \delta \vec{J}  \nabla_u {\cal J}^T \cdot \delta \vec{J}
488  \]  \]
489  where now $ \delta \vec{J} \in I\!\!R $ is a vector of dimenison $ l $.  where now $ \delta \vec{J} \in I\!\!R^l $ is a vector of
490    dimension $ l $.
491  In this case $ l $ reverse simulations have to be performed  In this case $ l $ reverse simulations have to be performed
492  for each $ \delta J_{k}, \,\, k = 1, \ldots, l $.  for each $ \delta J_{k}, \,\, k = 1, \ldots, l $.
493  Then, the reverse mode is more efficient as long as  Then, the reverse mode is more efficient as long as
494  $ l < n $, otherwise the forward mode is preferable.  $ l < n $, otherwise the forward mode is preferable.
495  Stricly, the reverse mode is called adjoint mode only for  Strictly, the reverse mode is called adjoint mode only for
496  $ l = 1 $.  $ l = 1 $.
497    
498  A detailed analysis of the underlying numerical operations  A detailed analysis of the underlying numerical operations
# Line 503  operator onto the $j$-th component ${\bf Line 522  operator onto the $j$-th component ${\bf
522  \paragraph{Example 2:  \paragraph{Example 2:
523  $ {\cal J} = \langle \, {\cal H}(\vec{v}) - \vec{d} \, ,  $ {\cal J} = \langle \, {\cal H}(\vec{v}) - \vec{d} \, ,
524   \, {\cal H}(\vec{v}) - \vec{d} \, \rangle $} ~ \\   \, {\cal H}(\vec{v}) - \vec{d} \, \rangle $} ~ \\
525  The cost function represents the quadratic model vs.data misfit.  The cost function represents the quadratic model vs. data misfit.
526  Here, $ \vec{d} $ is the data vector and $ {\cal H} $ represents the  Here, $ \vec{d} $ is the data vector and $ {\cal H} $ represents the
527  operator which maps the model state space onto the data space.  operator which maps the model state space onto the data space.
528  Then, $ \nabla_v {\cal J} $ takes the form  Then, $ \nabla_v {\cal J} $ takes the form
# Line 534  H \cdot \left( {\cal H}(\vec{v}) - \vec{ Line 553  H \cdot \left( {\cal H}(\vec{v}) - \vec{
553    
554  We note an important aspect of the forward vs. reverse  We note an important aspect of the forward vs. reverse
555  mode calculation.  mode calculation.
556  Because of the locality of the derivative,  Because of the local character of the derivative
557    (a derivative is defined w.r.t. a point along the trajectory),
558  the intermediate results of the model trajectory  the intermediate results of the model trajectory
559  $\vec{v}^{(\lambda+1)}={\cal M}_{\lambda}(v^{(\lambda)})$  $\vec{v}^{(\lambda+1)}={\cal M}_{\lambda}(v^{(\lambda)})$
560  are needed to evaluate the intermediate Jacobian  are needed to evaluate the intermediate Jacobian
561  $M_{\lambda}|_{\vec{v}^{(\lambda)}} \, \delta \vec{v}^{(\lambda)} $.  $M_{\lambda}|_{\vec{v}^{(\lambda)}} \, \delta \vec{v}^{(\lambda)} $.
562  In the forward mode, the intermediate results are required  In the forward mode, the intermediate results are required
563  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}$,
564  in the reverse mode they are required in the reverse order.  but in the reverse mode they are required in the reverse order.
565  Thus, in the reverse mode the trajectory of the forward model  Thus, in the reverse mode the trajectory of the forward model
566  integration ${\cal M}$ has to be stored to be available in the reverse  integration ${\cal M}$ has to be stored to be available in the reverse
567  calculation. Alternatively, the model state would have to be  calculation. Alternatively, the complete model state up to the
568  recomputed whenever its value is required.  point of evaluation has to be recomputed whenever its value is required.
569    
570  A method to balance the amount of recomputations vs.  A method to balance the amount of recomputations vs.
571  storage requirements is called {\sf checkpointing}  storage requirements is called {\sf checkpointing}
572  (e.g. \cite{res-eta:98}).  (e.g. \cite{res-eta:98}).
573  It is depicted in Fig. ... for a 3-level checkpointing  It is depicted in \ref{fig:3levelcheck} for a 3-level checkpointing
574  [as concrete example, we give explicit numbers for a 3-day  [as an example, we give explicit numbers for a 3-day
575  integration with a 1-hourly timestep in square brackets].  integration with a 1-hourly timestep in square brackets].
576  \begin{itemize}  \begin{itemize}
577  %  %
# Line 559  integration with a 1-hourly timestep in Line 579  integration with a 1-hourly timestep in
579  In a first step, the model trajectory is subdivided into  In a first step, the model trajectory is subdivided into
580  $ {n}^{lev3} $ subsections [$ {n}^{lev3} $=3 1-day intervals],  $ {n}^{lev3} $ subsections [$ {n}^{lev3} $=3 1-day intervals],
581  with the label $lev3$ for this outermost loop.  with the label $lev3$ for this outermost loop.
582  The model is then integrated over the full trajectory,  The model is then integrated along the full trajectory,
583  and the model state stored only at every $ k_{i}^{lev3} $-th timestep  and the model state stored only at every $ k_{i}^{lev3} $-th timestep
584  [i.e. 3 times, at  [i.e. 3 times, at
585  $ i = 0,1,2 $ corresponding to $ k_{i}^{lev3} = 0, 24, 48 $].  $ i = 0,1,2 $ corresponding to $ k_{i}^{lev3} = 0, 24, 48 $].
586  %  %
587  \item [$lev2$]  \item [$lev2$]
588  In a second step each subsection is itself divided into  In a second step each subsection itself is divided into
589  $ {n}^{lev2} $ subsubsections  $ {n}^{lev2} $ sub-subsections
590  [$ {n}^{lev2} $=4 6-hour intervals per subsection].  [$ {n}^{lev2} $=4 6-hour intervals per subsection].
591  The model picks up at the last outermost dumped state  The model picks up at the last outermost dumped state
592  $ v_{k_{n}^{lev3}} $ and is integrated forward in time over  $ v_{k_{n}^{lev3}} $ and is integrated forward in time along
593  the last subsection, with the label $lev2$ for this    the last subsection, with the label $lev2$ for this  
594  intermediate loop.  intermediate loop.
595  The model state is now stored only at every $ k_{i}^{lev2} $-th  The model state is now stored at every $ k_{i}^{lev2} $-th
596  timestep  timestep
597  [i.e. 4 times, at  [i.e. 4 times, at
598  $ 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 $].
599  %  %
600  \item [$lev1$]  \item [$lev1$]
601  Finally, the mode picks up at the last intermediate dump state  Finally, the model picks up at the last intermediate dump state
602  $ v_{k_{n}^{lev2}} $ and is integrated forward in time over  $ v_{k_{n}^{lev2}} $ and is integrated forward in time along
603  the last subsubsection, with the label $lev1$ for this    the last sub-subsection, with the label $lev1$ for this  
604  intermediate loop.  intermediate loop.
605  Within this subsubsection only, the model state is stored  Within this sub-subsection only, the model state is stored
606  at every timestep  at every timestep
607  [i.e. every hour $ i=0,...,5$ corresponding to  [i.e. every hour $ i=0,...,5$ corresponding to
608  $ k_{i}^{lev1} = 66, 67, \ldots, 71 $].  $ k_{i}^{lev1} = 66, 67, \ldots, 71 $].
609  Thus, the  final state $ v_n = v_{k_{n}^{lev1}} $ is reached  Thus, the  final state $ v_n = v_{k_{n}^{lev1}} $ is reached
610  and the model state of all peceeding timesteps over the last  and the model state of all  proceeding timesteps along the last
611  subsubsections are available, enabling integration backwards  sub-subsections are available, enabling integration backwards
612  in time over the last subsubsection.  in time along the last sub-subsection.
613  Thus, the adjoint can be computed over this last  Thus, the adjoint can be computed along this last
614  subsubsection $k_{n}^{lev2}$.  sub-subsection $k_{n}^{lev2}$.
615  %  %
616  \end{itemize}  \end{itemize}
617  %  %
618  This procedure is repeated consecutively for each previous  This procedure is repeated consecutively for each previous
619  subsubsection $k_{n-1}^{lev2}, \ldots, k_{1}^{lev2} $  sub-subsection $k_{n-1}^{lev2}, \ldots, k_{1}^{lev2} $
620  carrying the adjoint computation to the initial time  carrying the adjoint computation to the initial time
621  of the subsection $k_{n}^{lev3}$.  of the subsection $k_{n}^{lev3}$.
622  Then, the procedure is repeated for the previous subsection  Then, the procedure is repeated for the previous subsection
# Line 617  The balance of storage vs. recomputation Line 637  The balance of storage vs. recomputation
637  on the computing resources available.  on the computing resources available.
638    
639  \begin{figure}[t!]  \begin{figure}[t!]
640  \centering  \begin{center}
641  %\psdraft  %\psdraft
642  \psfrag{v_k1^lev3}{\mathinfigure{v_{k_{1}^{lev3}}}}  %\psfrag{v_k1^lev3}{\mathinfigure{v_{k_{1}^{lev3}}}}
643  \psfrag{v_kn-1^lev3}{\mathinfigure{v_{k_{n-1}^{lev3}}}}  %\psfrag{v_kn-1^lev3}{\mathinfigure{v_{k_{n-1}^{lev3}}}}
644  \psfrag{v_kn^lev3}{\mathinfigure{v_{k_{n}^{lev3}}}}  %\psfrag{v_kn^lev3}{\mathinfigure{v_{k_{n}^{lev3}}}}
645  \psfrag{v_k1^lev2}{\mathinfigure{v_{k_{1}^{lev2}}}}  %\psfrag{v_k1^lev2}{\mathinfigure{v_{k_{1}^{lev2}}}}
646  \psfrag{v_kn-1^lev2}{\mathinfigure{v_{k_{n-1}^{lev2}}}}  %\psfrag{v_kn-1^lev2}{\mathinfigure{v_{k_{n-1}^{lev2}}}}
647  \psfrag{v_kn^lev2}{\mathinfigure{v_{k_{n}^{lev2}}}}  %\psfrag{v_kn^lev2}{\mathinfigure{v_{k_{n}^{lev2}}}}
648  \psfrag{v_k1^lev1}{\mathinfigure{v_{k_{1}^{lev1}}}}  %\psfrag{v_k1^lev1}{\mathinfigure{v_{k_{1}^{lev1}}}}
649  \psfrag{v_kn^lev1}{\mathinfigure{v_{k_{n}^{lev1}}}}  %\psfrag{v_kn^lev1}{\mathinfigure{v_{k_{n}^{lev1}}}}
650  \mbox{\epsfig{file=part5/checkpointing.eps, width=0.8\textwidth}}  %\mbox{\epsfig{file=part5/checkpointing.eps, width=0.8\textwidth}}
651    \resizebox{5.5in}{!}{\includegraphics{part5/checkpointing.eps}}
652  %\psfull  %\psfull
653  \caption  \end{center}
654  {Schematic view of intermediate dump and restart for  \caption{
655    Schematic view of intermediate dump and restart for
656  3-level checkpointing.}  3-level checkpointing.}
657  \label{fig:erswns}  \label{fig:3levelcheck}
658  \end{figure}  \end{figure}
659    
660  \subsection{Optimal perturbations}  % \subsection{Optimal perturbations}
661  \label{optpert}  % \label{sec_optpert}
662    
663    
664  \subsection{Error covariance estimate and Hessian matrix}  % \subsection{Error covariance estimate and Hessian matrix}
665  \label{sec_hessian}  % \label{sec_hessian}
666    
667  \newpage  \newpage
668    
# Line 649  on the computing resources available. Line 671  on the computing resources available.
671  \label{sec_ad_setup_ex}  \label{sec_ad_setup_ex}
672  %**********************************************************************  %**********************************************************************
673    
674  The MITGCM has been adapted to enable AD using TAMC or TAF  The MITGCM has been adapted to enable AD using TAMC or TAF.
 (we'll refer to TAMC and TAF interchangeably, except where  
 distinctions are explicitly mentioned).  
675  The present description, therefore, is specific to the  The present description, therefore, is specific to the
676  use of TAMC as AD tool.  use of TAMC or TAF as AD tool.
677  The following sections describe the steps which are necessary to  The following sections describe the steps which are necessary to
678  generate a tangent linear or adjoint model of the MITGCM.  generate a tangent linear or adjoint model of the MITGCM.
679  We take as an example the sensitivity of carbon sequestration  We take as an example the sensitivity of carbon sequestration
680  in the ocean.  in the ocean.
681  The AD-relevant hooks in the code are sketched in  The AD-relevant hooks in the code are sketched in
682  \reffig{adthemodel}, \reffig{adthemain}.  \ref{fig:adthemodel}, \ref{fig:adthemain}.
683    
684  \subsection{Overview of the experiment}  \subsection{Overview of the experiment}
685    
686  We describe an adjoint sensitivity analysis of outgassing from  We describe an adjoint sensitivity analysis of out-gassing from
687  the ocean into the atmosphere of a carbon like tracer injected  the ocean into the atmosphere of a carbon-like tracer injected
688  into the ocean interior (see \cite{hil-eta:01}).  into the ocean interior (see \cite{hil-eta:01}).
689    
690  \subsubsection{Passive tracer equation}  \subsubsection{Passive tracer equation}
691    
692  For this work the MITGCM was augmented with a thermodynamically  For this work the MITGCM was augmented with a thermodynamically
693  inactive tracer, $C$. Tracer residing in the ocean  inactive tracer, $C$. Tracer residing in the ocean
694  model surface layer is outgassed according to a relaxation time scale,  model surface layer is out-gassed according to a relaxation time scale,
695  $\mu$. Within the ocean interior, the tracer is passively advected  $\mu$. Within the ocean interior, the tracer is passively advected
696  by the ocean model currents. The full equation for the time evolution  by the ocean model currents. The full equation for the time evolution
697  %  %
# Line 686  represents interior sources of $C$ such Line 706  represents interior sources of $C$ such
706  direct injection.  direct injection.
707  The velocity term, $U$, is the sum of the  The velocity term, $U$, is the sum of the
708  model Eulerian circulation and an eddy-induced velocity, the latter  model Eulerian circulation and an eddy-induced velocity, the latter
709  parameterized according to Gent/McWilliams (\cite{gen:90, dan:95}).  parameterized according to Gent/McWilliams
710    (\cite{gen-mcw:90, gen-eta:95}).
711  The convection function, $\Gamma$, mixes $C$ vertically wherever the  The convection function, $\Gamma$, mixes $C$ vertically wherever the
712  fluid is locally statically unstable.  fluid is locally statically unstable.
713    
714  The outgassing time scale, $\mu$, in eqn. (\ref{carbon_ddt})  The out-gassing time scale, $\mu$, in eqn. (\ref{carbon_ddt})
715  is set so that \( 1/\mu \sim 1 \ \mathrm{year} \) for the surface  is set so that \( 1/\mu \sim 1 \ \mathrm{year} \) for the surface
716  ocean and $\mu=0$ elsewhere. With this value, eqn. (\ref{carbon_ddt})  ocean and $\mu=0$ elsewhere. With this value, eqn. (\ref{carbon_ddt})
717  is valid as a prognostic equation for small perturbations in oceanic  is valid as a prognostic equation for small perturbations in oceanic
718  carbon concentrations. This configuration provides a  carbon concentrations. This configuration provides a
719  powerful tool for examining the impact of large-scale ocean circulation  powerful tool for examining the impact of large-scale ocean circulation
720  on $ CO_2 $ outgassing due to interior injections.  on $ CO_2 $ out-gassing due to interior injections.
721  As source we choose a constant in time injection of  As source we choose a constant in time injection of
722  $ S = 1 \,\, {\rm mol / s}$.  $ S = 1 \,\, {\rm mol / s}$.
723    
# Line 707  $4^\circ \times 4^\circ$ resolution hori Line 728  $4^\circ \times 4^\circ$ resolution hori
728  geography and bathymetry. Twenty vertical layers are used with  geography and bathymetry. Twenty vertical layers are used with
729  vertical spacing ranging  vertical spacing ranging
730  from 50 m near the surface to 815 m at depth.  from 50 m near the surface to 815 m at depth.
731  Driven to steady-state by climatalogical wind-stress, heat and  Driven to steady-state by climatological wind-stress, heat and
732  fresh-water forcing the model reproduces well known large-scale  fresh-water forcing the model reproduces well known large-scale
733  features of the ocean general circulation.  features of the ocean general circulation.
734    
735  \subsubsection{Outgassing cost function}  \subsubsection{Out-gassing cost function}
736    
737  To quantify and understand outgassing due to injections of $C$  To quantify and understand out-gassing due to injections of $C$
738  in eqn. (\ref{carbon_ddt}),  in eqn. (\ref{carbon_ddt}),
739  we define a cost function $ {\cal J} $ that measures the total amount of  we define a cost function $ {\cal J} $ that measures the total amount of
740  tracer outgassed at each timestep:  tracer out-gassed at each timestep:
741  %  %
742  \begin{equation}  \begin{equation}
743  \label{cost_tracer}  \label{cost_tracer}
744  {\cal J}(t=T)=\int_{t=0}^{t=T}\int_{A} \mu C \, dA \, dt  {\cal J}(t=T)=\int_{t=0}^{t=T}\int_{A} \mu C \, dA \, dt
745  \end{equation}  \end{equation}
746  %  %
747  Equation(\ref{cost_tracer}) integrates the outgassing term, $\mu C$,  Equation(\ref{cost_tracer}) integrates the out-gassing term, $\mu C$,
748  from (\ref{carbon_ddt})  from (\ref{carbon_ddt})
749  over the entire ocean surface area, $A$, and accumulates it  over the entire ocean surface area, $A$, and accumulates it
750  up to time $T$.  up to time $T$.
751  Physically, ${\cal J}$ can be thought of as representing the amount of  Physically, ${\cal J}$ can be thought of as representing the amount of
752  $CO_2$ that our model predicts would be outgassed following an  $CO_2$ that our model predicts would be out-gassed following an
753  injection at rate $S$.  injection at rate $S$.
754  The sensitivity of ${\cal J}$ to the spatial location of $S$,  The sensitivity of ${\cal J}$ to the spatial location of $S$,
755  $\frac{\partial {\cal J}}{\partial S}$,  $\frac{\partial {\cal J}}{\partial S}$,
756  can be used to identify regions from which circulation  can be used to identify regions from which circulation
757  would cause $CO_2$ to rapidly outgas following injection  would cause $CO_2$ to rapidly out-gas following injection
758  and regions in which $CO_2$ injections would remain effectively  and regions in which $CO_2$ injections would remain effectively
759  sequesterd within the ocean.  sequestered within the ocean.
760    
761  \subsection{Code configuration}  \subsection{Code configuration}
762    
763  The model configuration for this experiment resides under the  The model configuration for this experiment resides under the
764  directory {\it verification/carbon/}.  directory {\it verification/carbon/}.
765  The code customisation routines are in {\it verification/carbon/code/}:  The code customization routines are in {\it verification/carbon/code/}:
766  %  %
767  \begin{itemize}  \begin{itemize}
768  %  %
# Line 777  together with the forcing fields and and Line 798  together with the forcing fields and and
798  %  %
799  \item {\it data.ctrl}  \item {\it data.ctrl}
800  %  %
801    \item {\it data.gmredi}
802    %
803    \item {\it data.grdchk}
804    %
805    \item {\it data.optim}
806    %
807  \item {\it data.pkg}  \item {\it data.pkg}
808  %  %
809  \item {\it eedata}  \item {\it eedata}
# Line 803  $ adjoint/ $: Line 830  $ adjoint/ $:
830  \end{itemize}  \end{itemize}
831  %  %
832    
833  Below we describe the customisations of this files which are  Below we describe the customizations of this files which are
834  specific to this experiment.  specific to this experiment.
835    
836  \subsubsection{File {\it .genmakerc}}  \subsubsection{File {\it .genmakerc}}
837  This file overwites default settings of {\it genmake}.  This file overwrites default settings of {\it genmake}.
838  In the present example it is used to switch on the following  In the present example it is used to switch on the following
839  packages which are related to automatic differentiation  packages which are related to automatic differentiation
840  and are disabled by default: \\  and are disabled by default: \\
841  \hspace*{4ex} {\tt set ENABLE=( autodiff cost ctrl ecco )}  \\  \hspace*{4ex} {\tt set ENABLE=( autodiff cost ctrl ecco gmredi grdchk kpp )}  \\
842  Other packages which are not needed are switched off: \\  Other packages which are not needed are switched off: \\
843  \hspace*{4ex} {\tt set DISABLE=( aim obcs zonal\_filt shap\_filt cal exf )}  \hspace*{4ex} {\tt set DISABLE=( aim obcs zonal\_filt shap\_filt cal exf )}
844    
# Line 828  the standard include of the {\it CPP\_OP Line 855  the standard include of the {\it CPP\_OP
855    
856  This file contains 'wrapper'-specific CPP options.  This file contains 'wrapper'-specific CPP options.
857  It only needs to be changed if the code is to be run  It only needs to be changed if the code is to be run
858  in  parallel environment (see Section \ref{???}).  in a parallel environment (see Section \ref{???}).
859    
860  \subsubsection{File {\it CPP\_OPTIONS.h}}  \subsubsection{File {\it CPP\_OPTIONS.h}}
861    
# Line 837  This file contains model-specific CPP op Line 864  This file contains model-specific CPP op
864  Most options are related to the forward model setup.  Most options are related to the forward model setup.
865  They are identical to the global steady circulation setup of  They are identical to the global steady circulation setup of
866  {\it verification/exp2/}.  {\it verification/exp2/}.
867  The option specific to this experiment is \\  The three options specific to this experiment are \\
868    \hspace*{4ex} {\tt \#define ALLOW\_PASSIVE\_TRACER} \\
869    This flag enables the code to carry through the
870    advection/diffusion of a passive tracer along the
871    model integration. \\
872  \hspace*{4ex} {\tt \#define ALLOW\_MIT\_ADJOINT\_RUN} \\  \hspace*{4ex} {\tt \#define ALLOW\_MIT\_ADJOINT\_RUN} \\
873  This flag enables the inclusion of some AD-related fields  This flag enables the inclusion of some AD-related fields
874  concerning initialisation, link between control variables  concerning initialization, link between control variables
875  and forward model variables, and the call to the top-level  and forward model variables, and the call to the top-level
876  forward/adjoint subroutine {\it adthe\_main\_loop}  forward/adjoint subroutine {\it adthe\_main\_loop}
877  instead of {\it the\_main\_loop}.  instead of {\it the\_main\_loop}. \\
878    \hspace*{4ex} {\tt \#define ALLOW\_GRADIENT\_CHECK} \\
879    This flag enables the gradient check package.
880    After computing the unperturbed cost function and its gradient,
881    a series of computations are performed for which \\
882    $\bullet$ an element of the control vector is perturbed \\
883    $\bullet$ the cost function w.r.t. the perturbed element is
884    computed \\
885    $\bullet$ the difference between the perturbed and unperturbed
886    cost function is computed to compute the finite difference gradient \\
887    $\bullet$ the finite difference gradient is compared with the
888    adjoint-generated gradient.
889    The gradient check package is further described in Section ???.
890    
891  \subsubsection{File {\it ECCO\_OPTIONS.h}}  \subsubsection{File {\it ECCO\_OPTIONS.h}}
892    
# Line 864  enables the checkpointing feature of TAM Line 907  enables the checkpointing feature of TAM
907  (see Section \ref{???}).  (see Section \ref{???}).
908  In the present example a 3-level checkpointing is implemented.  In the present example a 3-level checkpointing is implemented.
909  The code contains the relevant store directives, common block  The code contains the relevant store directives, common block
910  and tape initialisations, storing key computation,  and tape initializations, storing key computation,
911  and loop index handling.  and loop index handling.
912  The checkpointing length at each level is defined in  The checkpointing length at each level is defined in
913  file {\it tamc.h}, cf. below.  file {\it tamc.h}, cf. below.
914  %  %
915  \item Cost function package: {\it pkg/cost/} \\  \item Cost function package: {\it pkg/cost/} \\
916  This package contains all relevant routines for  This package contains all relevant routines for
917  initialising, accumulating and finalizing the cost function  initializing, accumulating and finalizing the cost function
918  (see Section \ref{???}). \\  (see Section \ref{???}). \\
919  \hspace*{4ex} {\tt \#define ALLOW\_COST} \\  \hspace*{4ex} {\tt \#define ALLOW\_COST} \\
920  enables all general aspects of the cost function handling,  enables all general aspects of the cost function handling,
921  in particular the hooks in the foorward code for  in particular the hooks in the forward code for
922  initialising, accumulating and finalizing the cost function. \\  initializing, accumulating and finalizing the cost function. \\
923  \hspace*{4ex} {\tt \#define ALLOW\_COST\_TRACER} \\  \hspace*{4ex} {\tt \#define ALLOW\_COST\_TRACER} \\
924  includes the subroutine with the cost function for this  includes the call to the cost function for this
925  particular experiment, eqn. (\ref{cost_tracer}).  particular experiment, eqn. (\ref{cost_tracer}).
926  %  %
927  \item Control variable package: {\it pkg/ctrl/} \\  \item Control variable package: {\it pkg/ctrl/} \\
# Line 900  meridional wind stress \\ Line 943  meridional wind stress \\
943  freshwater flux \\  freshwater flux \\
944  \hspace*{2ex} {\tt \#define ALLOW\_HFLUX0\_CONTROL} &  \hspace*{2ex} {\tt \#define ALLOW\_HFLUX0\_CONTROL} &
945  heat flux \\  heat flux \\
946  \hspace*{2ex} {\tt \#undef ALLOW\_DIFFKR\_CONTROL} &  \hspace*{2ex} {\tt \#define ALLOW\_DIFFKR\_CONTROL} &
947  diapycnal diffusivity \\  diapycnal diffusivity \\
948  \hspace*{2ex} {\tt \#undef ALLOW\_KAPPAGM\_CONTROL} &  \hspace*{2ex} {\tt \#undef ALLOW\_KAPPAGM\_CONTROL} &
949  isopycnal diffusivity \\  isopycnal diffusivity \\
# Line 915  model. It is identical to the {\it verif Line 958  model. It is identical to the {\it verif
958  \hspace*{4ex} {\tt sNx = 90} \\  \hspace*{4ex} {\tt sNx = 90} \\
959  \hspace*{4ex} {\tt sNy = 40} \\  \hspace*{4ex} {\tt sNy = 40} \\
960  \hspace*{4ex} {\tt Nr = 20} \\  \hspace*{4ex} {\tt Nr = 20} \\
961  It correpsponds to a single-tile/single-processor setup:  It corresponds to a single-tile/single-processor setup:
962  {\tt nSx = nSy = 1, nPx = nPy = 1},  {\tt nSx = nSy = 1, nPx = nPy = 1},
963  with standard overlap dimensioning  with standard overlap dimensioning
964  {\tt OLx = OLy = 3}.  {\tt OLx = OLy = 3}.
# Line 932  The common blocks are used by the adjoin Line 975  The common blocks are used by the adjoin
975  \hspace*{4ex} is related to {\it DYNVARS.h} \\  \hspace*{4ex} is related to {\it DYNVARS.h} \\
976  \hspace*{4ex} {\tt common /addynvars\_cd/} &  \hspace*{4ex} {\tt common /addynvars\_cd/} &
977  \hspace*{4ex} is related to {\it DYNVARS.h} \\  \hspace*{4ex} is related to {\it DYNVARS.h} \\
978    \hspace*{4ex} {\tt common /addynvars\_diffkr/} &
979    \hspace*{4ex} is related to {\it DYNVARS.h} \\
980    \hspace*{4ex} {\tt common /addynvars\_kapgm/} &
981    \hspace*{4ex} is related to {\it DYNVARS.h} \\
982  \hspace*{4ex} {\tt common /adtr1\_r/} &  \hspace*{4ex} {\tt common /adtr1\_r/} &
983  \hspace*{4ex} is related to {\it TR1.h} \\  \hspace*{4ex} is related to {\it TR1.h} \\
984  \hspace*{4ex} {\tt common /adffields/} &  \hspace*{4ex} {\tt common /adffields/} &
# Line 956  This routine contains the dimensions for Line 1003  This routine contains the dimensions for
1003  3-level checkpointing is enabled, i.e. the timestepping  3-level checkpointing is enabled, i.e. the timestepping
1004  is divided into three different levels (see Section \ref{???}).  is divided into three different levels (see Section \ref{???}).
1005  The model state of the outermost ({\tt nchklev\_3}) and the  The model state of the outermost ({\tt nchklev\_3}) and the
1006  itermediate ({\tt nchklev\_2}) timestepping loop are stored to file  intermediate ({\tt nchklev\_2}) timestepping loop are stored to file
1007  (handled in {\it the\_main\_loop}).  (handled in {\it the\_main\_loop}).
1008  The innermost loop ({\tt nchklev\_1})  The innermost loop ({\tt nchklev\_1})
1009  avoids I/O by storing all required variables  avoids I/O by storing all required variables
# Line 968  In the present example the dimensions ar Line 1015  In the present example the dimensions ar
1015  \hspace*{4ex} {\tt nchklev\_2      =  30 } \\  \hspace*{4ex} {\tt nchklev\_2      =  30 } \\
1016  \hspace*{4ex} {\tt nchklev\_3      =  60 } \\  \hspace*{4ex} {\tt nchklev\_3      =  60 } \\
1017  To guarantee that the checkpointing intervals span the entire  To guarantee that the checkpointing intervals span the entire
1018  integration period the relation \\  integration period the following relation must be satisfied: \\
1019  \hspace*{4ex} {\tt nchklev\_1*nchklev\_2*nchklev\_3 $ \ge $ nTimeSteps} \\  \hspace*{4ex} {\tt nchklev\_1*nchklev\_2*nchklev\_3 $ \ge $ nTimeSteps} \\
1020  where {\tt nTimeSteps} is either specified in {\it data}  where {\tt nTimeSteps} is either specified in {\it data}
1021  or computed via \\  or computed via \\
# Line 982  Similar to above, the following relation Line 1029  Similar to above, the following relation
1029  %  %
1030  \end{itemize}  \end{itemize}
1031    
1032    The following parameters may be worth describing: \\
1033    %
1034    \hspace*{4ex} {\tt isbyte} \\
1035    \hspace*{4ex} {\tt maxpass} \\
1036    ~
1037    
1038  \subsubsection{File {\it makefile}}  \subsubsection{File {\it makefile}}
1039    
1040  This file contains all relevant paramter flags and  This file contains all relevant parameter flags and
1041  lists to run TAMC.  lists to run TAMC or TAF.
1042  It is assumed that TAMC is available to you, either locally,  It is assumed that TAMC is available to you, either locally,
1043  being installed on your network, or remotely through the 'TAMC Utility'.  being installed on your network, or remotely through the 'TAMC Utility'.
1044  TAMC is called with the command {\tt tamc} followed by a  TAMC is called with the command {\tt tamc} followed by a
# Line 996  Here we briefly discuss the main flags u Line 1049  Here we briefly discuss the main flags u
1049  \begin{itemize}  \begin{itemize}
1050  \item [{\tt tamc}] {\tt  \item [{\tt tamc}] {\tt
1051  -input <variable names>  -input <variable names>
1052  -output <variable name> ... \\  -output <variable name> -r4 ... \\
1053  -toplevel <S/R name> -reverse <file names>  -toplevel <S/R name> -reverse <file names>
1054  }  }
1055  \end{itemize}  \end{itemize}
# Line 1017  Dependent variable $ J $  which is to be Line 1070  Dependent variable $ J $  which is to be
1070  \item {\tt -reverse <file names>} \\  \item {\tt -reverse <file names>} \\
1071  Adjoint code is generated to compute the sensitivity of an  Adjoint code is generated to compute the sensitivity of an
1072  independent variable w.r.t.  many dependent variables.  independent variable w.r.t.  many dependent variables.
1073  The generated adjoint top-level routine computes the product  In the discussion of Section ???
1074    the generated adjoint top-level routine computes the product
1075  of the transposed Jacobian matrix $ M^T $ times  of the transposed Jacobian matrix $ M^T $ times
1076  the gradient vector $ \nabla_v J $.  the gradient vector $ \nabla_v J $.
1077  \\  \\
1078  {\tt <file names>} refers to the list of files {\it .f} which are to be  {\tt <file names>} refers to the list of files {\it .f} which are to be
1079  analyzed by TAMC. This list is generally smaller than the full list  analyzed by TAMC. This list is generally smaller than the full list
1080  of code to be compiled. The files not contained are either  of code to be compiled. The files not contained are either
1081  above the top-level routine (some initialisations), or are  above the top-level routine (some initializations), or are
1082  deliberately hidden from TAMC, either because hand-written  deliberately hidden from TAMC, either because hand-written
1083  adjoint routines exist, or the routines must not (or don't have to)  adjoint routines exist, or the routines must not (or don't have to)
1084  be differentiated. For each routine which is part of the flow tree  be differentiated. For each routine which is part of the flow tree
1085  of the top-level routine, but deliberately hidden from TAMC,  of the top-level routine, but deliberately hidden from TAMC
1086    (or for each package which contains such routines),
1087  a corresponding file {\it .flow} exists containing flow directives  a corresponding file {\it .flow} exists containing flow directives
1088  for TAMC.  for TAMC.
1089  %  %
1090    \item {\tt -r4} \\
1091    ~
1092    %
1093  \end{itemize}  \end{itemize}
1094    
1095    
1096  \subsubsection{File {\it data}}  \subsubsection{The input parameter files}
1097    
1098    \paragraph{File {\it data}}
1099    
1100    \paragraph{File {\it data.cost}}
1101    
1102    \paragraph{File {\it data.ctrl}}
1103    
1104  \subsubsection{File {\it data.cost}}  \paragraph{File {\it data.gmredi}}
1105    
1106  \subsubsection{File {\it data.ctrl}}  \paragraph{File {\it data.grdchk}}
1107    
1108  \subsubsection{File {\it data.pkg}}  \paragraph{File {\it data.optim}}
1109    
1110  \subsubsection{File {\it eedata}}  \paragraph{File {\it data.pkg}}
1111    
1112  \subsubsection{File {\it topog.bin}}  \paragraph{File {\it eedata}}
1113    
1114  \subsubsection{File {\it windx.bin, windy.bin}}  \paragraph{File {\it topog.bin}}
1115    
1116  \subsubsection{File {\it salt.bin, theta.bin}}  \paragraph{File {\it windx.bin, windy.bin}}
1117    
1118  \subsubsection{File {\it SSS.bin, SST.bin}}  \paragraph{File {\it salt.bin, theta.bin}}
1119    
1120  \subsubsection{File {\it pickup*}}  \paragraph{File {\it SSS.bin, SST.bin}}
1121    
1122  \subsection{Compiling the model and its adjoint}  \paragraph{File {\it pickup*}}
1123    
1124    \subsection{Compiling the model and its adjoint}
1125    
1126    The built process of the adjoint model is slightly more
1127    complex than that of compiling the forward code.
1128    The main reason is that the adjoint code generation requires
1129    a specific list of routines that are to be differentiated
1130    (as opposed to the automatic generation of a list of
1131    files to be compiled by genmake).
1132    This list excludes routines that don't have to be or must not be
1133    differentiated. For some of the latter routines flow directives
1134    may be necessary, a list of which has to be given as well.
1135    For this reason, a separate {\it makefile} is currently
1136    maintained in the directory {\tt adjoint/}. This
1137    makefile is responsible for the adjoint code generation.
1138    
1139    In the following we describe the build process step by step,
1140    assuming you are in the directory {\tt bin/}.
1141    A summary of steps to follow is given at the end.
1142    
1143    \paragraph{Adjoint code generation and compilation -- step by step}
1144    
1145    \begin{enumerate}
1146    %
1147    \item
1148    {\tt ln -s ../verification/???/code/.genmakerc .} \\
1149    {\tt ln -s ../verification/???/code/*.[Fh] .} \\
1150    Link your customized genmake options, header files,
1151    and modified code to the compile directory.
1152    %
1153    \item
1154    {\tt ../tools/genmake -makefile} \\
1155    Generate your Makefile (cf. Section ???).
1156    %
1157    \item
1158    {\tt make depend} \\
1159    Dependency analysis for the CPP pre-compiler (cf. Section ???).
1160    %
1161    \item
1162    {\tt make small\_f} \\
1163    This is the first difference between forward code compilation
1164    and adjoint code generation and compilation.
1165    Instead of going through the entire compilation process
1166    (CPP precompiling -- {\tt .f}, object code generation -- {\tt .o},
1167    linking of object files and libraries to generate executable),
1168    only the CPP compiler is invoked at this stage to generate
1169    the {\tt .f} files.
1170    %
1171    \item
1172    {\tt cd ../adjoint} \\
1173    {\tt make adtaf} or {\tt make adtamc} \\
1174    Depending on whether you have TAF or TAMC at your disposal,
1175    you'll choose {\tt adtaf} or {\tt adtamc} as your
1176    make target for the {\it makefile} in the directory {\tt adjoint/}.
1177    Several things happen at this stage.
1178    %
1179    \begin{enumerate}
1180    %
1181    \item
1182    The initial template file {\it adjoint\_model.F} which is part
1183    of the compiling list created by {\it genmake} is restored.
1184    %
1185    \item
1186    All Fortran routines {\tt *.f} in {\tt bin/} are
1187    concatenated into a single file (it's current name is
1188    {\it tamc\_code.f}).
1189    %
1190    \item
1191    Adjoint code is generated by TAMC or TAF.
1192    The adjoint code is written to the file {\it tamc\_code\_ad.f}.
1193    It contains all adjoint routines of the forward routines
1194    concatenated in {\it tamc\_code.f}.
1195    For a given forward routines {\tt subroutine routinename}
1196    the adjoint routine is named {\tt adsubroutine routinename}
1197    by default (that default can be changed via the flag
1198    {\tt -admark <markname>}).
1199    Furthermore, it may contain modified code which
1200    incorporates the translation of adjoint store directives
1201    into specific Fortran code.
1202    For a given forward routines {\tt subroutine routinename}
1203    the modified routine is named {\tt mdsubroutine routinename}.
1204    TAMC or TAF info is written to file
1205    {\it tamc\_code.prot} or {\it taf.log}, respectively.
1206    %
1207    \end{enumerate}
1208    %
1209    \item
1210    {\tt make adchange} \\
1211    The multi-threading capability of the MITGCM requires a slight
1212    change in the parameter list of some routines that are related to
1213    to active file handling.
1214    This post-processing invokes the sed script {\it adjoint\_ecco\_sed.com}
1215    to insert the threading counter {\bf myThId} into the parameter list
1216    of those subroutines.
1217    The resulting code is written to file {\it tamc\_code\_sed\_ad.f}
1218    and appended to the file {\it adjoint\_model.F}.
1219    This concludes the adjoint code generation.
1220    %
1221    \item
1222    {\tt cd ../bin} \\
1223    {\tt make} \\
1224    The file {\it adjoint\_model.F} now contains the full adjoint code.
1225    All routines are now compiled.
1226    %
1227    \end{enumerate}
1228    
1229    \paragraph{Adjoint code generation and compilation -- summary}
1230    ~ \\
1231    
1232    \[
1233    \boxed{
1234    \begin{split}
1235     ~ & \mbox{\tt cd bin} \\
1236     ~ & \mbox{\tt ln -s ../verification/my\_experiment/code/.genmakerc .} \\
1237     ~ & \mbox{\tt ln -s ../verification/my\_experiment/code/*.[Fh] .} \\
1238     ~ & \mbox{\tt ../tools/genmake -makefile} \\
1239     ~ & \mbox{\tt make depend} \\
1240     ~ & \mbox{\tt make small\_f} \\
1241     ~ & \mbox{\tt cd ../adjoint} \\
1242     ~ & \mbox{\tt make adtaf <OR: make adtamc>} \\
1243     ~ & \mbox{\tt make adchange} \\
1244     ~ & \mbox{\tt cd ../bin} \\
1245     ~ & \mbox{\tt make} \\
1246    \end{split}
1247    }
1248    \]
1249    
1250  \newpage  \newpage
1251    
1252  %**********************************************************************  %**********************************************************************
1253  \section{TLM and ADM code generation in general}  \section{TLM and ADM generation in general}
1254  \label{sec_ad_setup_gen}  \label{sec_ad_setup_gen}
1255  %**********************************************************************  %**********************************************************************
1256    
# Line 1068  In this section we describe in a general Line 1258  In this section we describe in a general
1258  the parts of the code that are relevant for automatic  the parts of the code that are relevant for automatic
1259  differentiation using the software tool TAMC.  differentiation using the software tool TAMC.
1260    
1261  \subsection{The cost function (dependent variable)}  \begin{figure}[b!]
1262    \input{part5/doc_ad_the_model}
1263    \caption{~}
1264    \label{fig:adthemodel}
1265    \end{figure}
1266    
1267    The basic flow is depicted in \ref{fig:adthemodel}.
1268    If the option {\tt ALLOW\_AUTODIFF\_TAMC} is defined, the driver routine
1269    {\it the\_model\_main}, instead of calling {\it the\_main\_loop},
1270    invokes the adjoint of this routine, {\it adthe\_main\_loop},
1271    which is the toplevel routine in terms of reverse mode computation.
1272    The routine {\it adthe\_main\_loop} has been generated using TAMC.
1273    It contains both the forward integration of the full model,
1274    any additional storing that is required for efficient checkpointing,
1275    and the reverse integration of the adjoint model.
1276    The structure of {\it adthe\_main\_loop} has been strongly
1277    simplified for clarification; in particular, no checkpointing
1278    procedures are shown here.
1279    Prior to the call of {\it adthe\_main\_loop}, the routine
1280    {\it ctrl\_unpack} is invoked to unpack the control vector,
1281    and following that call, the routine {\it ctrl\_pack}
1282    is invoked to pack the control vector
1283    (cf. Section \ref{section_ctrl}).
1284    If gradient checks are to be performed, the option
1285    {\tt ALLOW\_GRADIENT\_CHECK} is defined. In this case
1286    the driver routine {\it grdchk\_main} is called after
1287    the gradient has been computed via the adjoint
1288    (cf. Section \ref{section_grdchk}).
1289    
1290    \subsection{The cost function (dependent variable)
1291    \label{section_cost}}
1292    
1293  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}.
1294  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
# Line 1076  $ {\cal J}(\vec{u}) \, = \, {\cal J}(M(\ Line 1296  $ {\cal J}(\vec{u}) \, = \, {\cal J}(M(\
1296  The input is referred to as the  The input is referred to as the
1297  {\sf independent variables} or {\sf control variables}.  {\sf independent variables} or {\sf control variables}.
1298  All aspects relevant to the treatment of the cost function $ {\cal J} $  All aspects relevant to the treatment of the cost function $ {\cal J} $
1299  (parameter setting, initialisation, incrementation,  (parameter setting, initialization, accumulation,
1300  final evaluation), are controled by the package {\it pkg/cost}.  final evaluation), are controlled by the package {\it pkg/cost}.
1301    
1302    \begin{figure}[h!]
1303    \input{part5/doc_cost_flow}
1304    \caption{~}
1305    \label{fig:costflow}
1306    \end{figure}
1307    
1308  \subsubsection{genmake and CPP options}  \subsubsection{genmake and CPP options}
1309  %  %
# Line 1097  compile list in 3 different ways (cf. Se Line 1323  compile list in 3 different ways (cf. Se
1323  \begin{enumerate}  \begin{enumerate}
1324  %  %
1325  \item {\it genmake}: \\  \item {\it genmake}: \\
1326  Change the default settngs in the file {\it genmake} by adding  Change the default settings in the file {\it genmake} by adding
1327  {\bf cost} to the {\bf enable} list (not recommended).  {\bf cost} to the {\bf enable} list (not recommended).
1328  %  %
1329  \item {\it .genmakerc}: \\  \item {\it .genmakerc}: \\
# Line 1110  Call {\it genmake} with the option Line 1336  Call {\it genmake} with the option
1336  {\tt genmake -enable=cost}.  {\tt genmake -enable=cost}.
1337  %  %
1338  \end{enumerate}  \end{enumerate}
 Since the cost function is usually used in conjunction with  
 automatic differentiation, the CPP option  
 {\bf ALLOW\_ADJOINT\_RUN} should be defined  
 (file {\it CPP\_OPTIONS.h}).  
1339  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}.
1340  Each specific cost function contribution has its own option.  Each specific cost function contribution has its own option.
1341  For the present example the option is {\bf ALLOW\_COST\_TRACER}.  For the present example the option is {\bf ALLOW\_COST\_TRACER}.
1342  All cost-specific options are set in {\it ECCO\_CPPOPTIONS.h}  All cost-specific options are set in {\it ECCO\_CPPOPTIONS.h}
1343    Since the cost function is usually used in conjunction with
1344    automatic differentiation, the CPP option
1345    {\bf ALLOW\_ADJOINT\_RUN} should be defined
1346    (file {\it CPP\_OPTIONS.h}).
1347    
1348  \subsubsection{Initialisation}  \subsubsection{Initialization}
1349  %  %
1350  The initialisation of the {\it cost} package is readily enabled  The initialization of the {\it cost} package is readily enabled
1351  as soon as the CPP option {\bf ALLOW\_ADJOINT\_RUN} is defined.  as soon as the CPP option {\bf ALLOW\_ADJOINT\_RUN} is defined.
1352  %  %
1353  \begin{itemize}  \begin{itemize}
# Line 1152  Variables: {\it cost\_init} Line 1378  Variables: {\it cost\_init}
1378  }  }
1379  \\  \\
1380  This S/R  This S/R
1381  initialises the different cost function contributions.  initializes the different cost function contributions.
1382  The contribtion for the present example is {\bf objf\_tracer}  The contribution for the present example is {\bf objf\_tracer}
1383  which is defined on each tile (bi,bj).  which is defined on each tile (bi,bj).
1384  %  %
1385  \end{itemize}  \end{itemize}
1386  %  %
1387  \subsubsection{Incrementation}  \subsubsection{Accumulation}
1388  %  %
1389  \begin{itemize}  \begin{itemize}
1390  %  %
# Line 1206  The total cost function {\bf fc} will be Line 1432  The total cost function {\bf fc} will be
1432  tamc -output 'fc' ...  tamc -output 'fc' ...
1433  \end{verbatim}  \end{verbatim}
1434    
1435  \begin{figure}[t!]  %%%% \end{document}
 \input{part5/doc_ad_the_model}  
 \label{fig:adthemodel}  
 \caption{~}  
 \end{figure}  
1436    
1437  \begin{figure}  \begin{figure}
1438  \input{part5/doc_ad_the_main}  \input{part5/doc_ad_the_main}
 \label{fig:adthemain}  
1439  \caption{~}  \caption{~}
1440    \label{fig:adthemain}
1441  \end{figure}  \end{figure}
1442    
1443  \subsection{The control variables (independent variables)}  \subsection{The control variables (independent variables)
1444    \label{section_ctrl}}
1445    
1446  The control variables are a subset of the model input  The control variables are a subset of the model input
1447  (initial conditions, boundary conditions, model parameters).  (initial conditions, boundary conditions, model parameters).
1448  Here we identify them with the variable $ \vec{u} $.  Here we identify them with the variable $ \vec{u} $.
1449  All intermediate variables whose derivative w.r.t. control  All intermediate variables whose derivative w.r.t. control
1450  variables don't vanish are called {\sf active variables}.  variables do not vanish are called {\sf active variables}.
1451  All subroutines whose derivative w.r.t. the control variables  All subroutines whose derivative w.r.t. the control variables
1452  don't vanish are called {\sf active routines}.  don't vanish are called {\sf active routines}.
1453  Read and write operations from and to file can be viewed  Read and write operations from and to file can be viewed
# Line 1232  as variable assignments. Therefore, file Line 1455  as variable assignments. Therefore, file
1455  active variables are written and from which active variables  active variables are written and from which active variables
1456  are read are called {\sf active files}.  are read are called {\sf active files}.
1457  All aspects relevant to the treatment of the control variables  All aspects relevant to the treatment of the control variables
1458  (parameter setting, initialisation, perturbation)  (parameter setting, initialization, perturbation)
1459  are controled by the package {\it pkg/ctrl}.  are controlled by the package {\it pkg/ctrl}.
1460    
1461    \begin{figure}[h!]
1462    \input{part5/doc_ctrl_flow}
1463    \caption{~}
1464    \label{fig:ctrlflow}
1465    \end{figure}
1466    
1467  \subsubsection{genmake and CPP options}  \subsubsection{genmake and CPP options}
1468  %  %
# Line 1253  To enable the directory to be included t Line 1482  To enable the directory to be included t
1482  Each control variable is enabled via its own CPP option  Each control variable is enabled via its own CPP option
1483  in {\it ECCO\_CPPOPTIONS.h}.  in {\it ECCO\_CPPOPTIONS.h}.
1484    
1485  \subsubsection{Initialisation}  \subsubsection{Initialization}
1486  %  %
1487  \begin{itemize}  \begin{itemize}
1488  %  %
# Line 1293  Two important issues related to the hand Line 1522  Two important issues related to the hand
1522  variables in the MITGCM need to be addressed.  variables in the MITGCM need to be addressed.
1523  First, in order to save memory, the control variable arrays  First, in order to save memory, the control variable arrays
1524  are not kept in memory, but rather read from file and added  are not kept in memory, but rather read from file and added
1525  to the initial (or first guess) fields.  to the initial fields during the model initialization phase.
1526  Similarly, the corresponding adjoint fields which represent  Similarly, the corresponding adjoint fields which represent
1527  the gradient of the cost function w.r.t. the control variables  the gradient of the cost function w.r.t. the control variables
1528  are written to to file.  are written to file at the end of the adjoint integration.
1529  Second, in addition to the files holding the 2-dim. and 3-dim.  Second, in addition to the files holding the 2-dim. and 3-dim.
1530  control variables and the gradient, a 1-dim. {\sf control vector}  control variables and the corresponding cost gradients,
1531    a 1-dim. {\sf control vector}
1532  and {\sf gradient vector} are written to file. They contain  and {\sf gradient vector} are written to file. They contain
1533  only the wet points of the control variables and the corresponding  only the wet points of the control variables and the corresponding
1534  gradient.  gradient.
1535  This leads to a significant data compression.  This leads to a significant data compression.
1536  Furthermore, the control and the gradient vector can be passed to a  Furthermore, an option is available
1537    ({\tt ALLOW\_NONDIMENSIONAL\_CONTROL\_IO}) to
1538    non-dimensionalise the control and gradient vector,
1539    which otherwise would contain different pieces of different
1540    magnitudes and units.
1541    Finally, the control and gradient vector can be passed to a
1542  minimization routine if an update of the control variables  minimization routine if an update of the control variables
1543  is sought as part of a minimization exercise.  is sought as part of a minimization exercise.
1544    
# Line 1314  and gradient are generated and initialis Line 1549  and gradient are generated and initialis
1549    
1550  \subsubsection{Perturbation of the independent variables}  \subsubsection{Perturbation of the independent variables}
1551  %  %
1552  The dependency chain for differentiation starts  The dependency flow for differentiation w.r.t. the controls
1553  with adding a perturbation onto the the input variable,  starts with adding a perturbation onto the input variable,
1554  thus defining the independent or control variables for TAMC.  thus defining the independent or control variables for TAMC.
1555  Three classes of controls may be considered:  Three types of controls may be considered:
1556  %  %
1557  \begin{itemize}  \begin{itemize}
1558  %  %
# Line 1332  Three classes of controls may be conside Line 1567  Three classes of controls may be conside
1567  Consider as an example the initial tracer distribution  Consider as an example the initial tracer distribution
1568  {\bf tr1} as control variable.  {\bf tr1} as control variable.
1569  After {\bf tr1} has been initialised in  After {\bf tr1} has been initialised in
1570  {\it ini\_tr1} (dynamical variables including  {\it ini\_tr1} (dynamical variables such as
1571  temperature and salinity are initialised in {\it ini\_fields}),  temperature and salinity are initialised in {\it ini\_fields}),
1572  a perturbation anomaly is added to the field in S/R  a perturbation anomaly is added to the field in S/R
1573  {\it ctrl\_map\_ini}  {\it ctrl\_map\_ini}
# Line 1345  u         & = \, u_{[0]} \, + \, \Delta Line 1580  u         & = \, u_{[0]} \, + \, \Delta
1580  \end{split}  \end{split}
1581  \end{equation}  \end{equation}
1582  %  %
1583  In principle {\bf xx\_tr1} is a 3-dim. global array  {\bf xx\_tr1} is a 3-dim. global array
1584  holding the perturbation. In the case of a simple  holding the perturbation. In the case of a simple
1585  sensitivity study this array is identical to zero.  sensitivity study this array is identical to zero.
1586  However, it's specification is essential since TAMC  However, it's specification is essential in the context
1587    of automatic differentiation since TAMC
1588  treats the corresponding line in the code symbolically  treats the corresponding line in the code symbolically
1589  when determining the differentiation chain and its origin.  when determining the differentiation chain and its origin.
1590  Thus, the variable names are part of the argument list  Thus, the variable names are part of the argument list
# Line 1366  dummy variable {\bf xx\_tr1\_dummy} is i Line 1602  dummy variable {\bf xx\_tr1\_dummy} is i
1602  and an 'active read' routine of the adjoint support  and an 'active read' routine of the adjoint support
1603  package {\it pkg/autodiff} is invoked.  package {\it pkg/autodiff} is invoked.
1604  The read-procedure is tagged with the variable  The read-procedure is tagged with the variable
1605  {\bf xx\_tr1\_dummy} enabbling TAMC to recognize the  {\bf xx\_tr1\_dummy} enabling TAMC to recognize the
1606  initialisation of the perturbation.  initialization of the perturbation.
1607  The modified call of TAMC thus reads  The modified call of TAMC thus reads
1608  %  %
1609  \begin{verbatim}  \begin{verbatim}
# Line 1388  Note, that reading an active variable co Line 1624  Note, that reading an active variable co
1624  to a variable assignment. Its derivative corresponds  to a variable assignment. Its derivative corresponds
1625  to a write statement of the adjoint variable.  to a write statement of the adjoint variable.
1626  The 'active file' routines have been designed  The 'active file' routines have been designed
1627  to support active read and corresponding active write  to support active read and corresponding adjoint active write
1628  operations.  operations (and vice versa).
1629  %  %
1630  \item  \item
1631  \fbox{  \fbox{
# Line 1406  with the symbolic perturbation taking pl Line 1642  with the symbolic perturbation taking pl
1642  Note however an important difference:  Note however an important difference:
1643  Since the boundary values are time dependent with a new  Since the boundary values are time dependent with a new
1644  forcing field applied at each time steps,  forcing field applied at each time steps,
1645  the general problem may be be thought of as  the general problem may be thought of as
1646  a new control variable at each time step, i.e.  a new control variable at each time step
1647    (or, if the perturbation is averaged over a certain period,
1648    at each $ N $ timesteps), i.e.
1649  \[  \[
1650  u_{\rm forcing} \, = \,  u_{\rm forcing} \, = \,
1651  \{ \, u_{\rm forcing} ( t_n ) \, \}_{  \{ \, u_{\rm forcing} ( t_n ) \, \}_{
# Line 1432  calendar ({\it cal}~) and external forci Line 1670  calendar ({\it cal}~) and external forci
1670  %  %
1671  This routine is not yet implemented, but would proceed  This routine is not yet implemented, but would proceed
1672  proceed along the same lines as the initial value sensitivity.  proceed along the same lines as the initial value sensitivity.
1673    The mixing parameters {\bf diffkr} and {\bf kapgm}
1674    are currently added as controls in {\it ctrl\_map\_ini.F}.
1675  %  %
1676  \end{itemize}  \end{itemize}
1677  %  %
1678    
1679  \subsubsection{Output of adjoint variables and gradient}  \subsubsection{Output of adjoint variables and gradient}
1680  %  %
1681  Two ways exist to generate output of adjoint fields.  Several ways exist to generate output of adjoint fields.
1682  %  %
1683  \begin{itemize}  \begin{itemize}
1684  %  %
1685  \item  \item
1686  \fbox{  \fbox{
1687  \begin{minipage}{12cm}  \begin{minipage}{12cm}
1688  {\it ctrl\_pack}:  {\it ctrl\_map\_ini, ctrl\_map\_forcing}:
1689  \end{minipage}  \end{minipage}
1690  }  }
1691  \\  \\
 At the end of the forward/adjoint integration, the S/R  
 {\it ctrl\_pack} is called which mirrors S/R {\it ctrl\_unpack}.  
 It writes the following files:  
 %  
1692  \begin{itemize}  \begin{itemize}
1693  %  %
1694  \item {\bf xx\_...}: the control variable fields  \item {\bf xx\_...}: the control variable fields \\
1695    Before the forward integration, the control
1696    variables are read from file {\bf xx\_ ...} and added to
1697    the model field.
1698  %  %
1699  \item {\bf adxx\_...}: the adjoint variable fields, i.e. the gradient  \item {\bf adxx\_...}: the adjoint variable fields, i.e. the gradient
1700  $ \nabla _{u}{\cal J} $ for each control variable,  $ \nabla _{u}{\cal J} $ for each control variable \\
1701    After the adjoint integration the corresponding adjoint
1702    variables are written to {\bf adxx\_ ...}.
1703  %  %
1704  \item {\bf vector\_ctrl}: the control vector  \end{itemize}
1705  %  %
1706  \item {\bf vector\_grad}: the gradient vector  \item
1707    \fbox{
1708    \begin{minipage}{12cm}
1709    {\it ctrl\_unpack, ctrl\_pack}:
1710    \end{minipage}
1711    }
1712    \\
1713    %
1714    \begin{itemize}
1715    %
1716    \item {\bf vector\_ctrl}: the control vector \\
1717    At the very beginning of the model initialization,
1718    the updated compressed control vector is read (or initialised)
1719    and distributed to 2-dim. and 3-dim. control variable fields.
1720    %
1721    \item {\bf vector\_grad}: the gradient vector \\
1722    At the very end of the adjoint integration,
1723    the 2-dim. and 3-dim. adjoint variables are read,
1724    compressed to a single vector and written to file.
1725  %  %
1726  \end{itemize}  \end{itemize}
1727  %  %
# Line 1474  $ \nabla _{u}{\cal J} $ for each control Line 1733  $ \nabla _{u}{\cal J} $ for each control
1733  }  }
1734  \\  \\
1735  In addition to writing the gradient at the end of the  In addition to writing the gradient at the end of the
1736  forward/adjoint integration, many more adjoint variables,  forward/adjoint integration, many more adjoint variables
1737  representing the Lagrange multipliers of the model state  of the model state
1738  w.r.t. the model state  at intermediate times can be written using S/R
 at different times can be written using S/R  
1739  {\it addummy\_in\_stepping}.  {\it addummy\_in\_stepping}.
1740  This routine is part of the adjoint support package  This routine is part of the adjoint support package
1741  {\it pkg/autodiff} (cf.f. below).  {\it pkg/autodiff} (cf.f. below).
# Line 1491  than generated automatically. Line 1749  than generated automatically.
1749  Appropriate flow directives ({\it dummy\_in\_stepping.flow})  Appropriate flow directives ({\it dummy\_in\_stepping.flow})
1750  ensure that TAMC does not automatically  ensure that TAMC does not automatically
1751  generate {\it addummy\_in\_stepping} by trying to differentiate  generate {\it addummy\_in\_stepping} by trying to differentiate
1752  {\it dummy\_in\_stepping}, but rather takes the hand-written routine.  {\it dummy\_in\_stepping}, but instead refers to
1753    the hand-written routine.
1754    
1755  {\it dummy\_in\_stepping} is called in the forward code  {\it dummy\_in\_stepping} is called in the forward code
1756  at the beginning of each  at the beginning of each
# Line 1501  each timestep in the adjoint calculation Line 1760  each timestep in the adjoint calculation
1760  {\it addynamics}.  {\it addynamics}.
1761    
1762  {\it addummy\_in\_stepping} includes the header files  {\it addummy\_in\_stepping} includes the header files
1763  {\it adffields.h, addynamics.h, adtr1.h}.  {\it adcommon.h}.
1764  These header files are also hand-written. They contain  This header file is also hand-written. It contains
1765  the common blocks {\bf /addynvars\_r/}, {\bf /addynvars\_cd/},  the common blocks
1766    {\bf /addynvars\_r/}, {\bf /addynvars\_cd/},
1767    {\bf /addynvars\_diffkr/}, {\bf /addynvars\_kapgm/},
1768  {\bf /adtr1\_r/}, {\bf /adffields/},  {\bf /adtr1\_r/}, {\bf /adffields/},
1769  which have been extracted from the adjoint code to enable  which have been extracted from the adjoint code to enable
1770  access to the adjoint variables.  access to the adjoint variables.
# Line 1521  The gradient $ \nabla _{u}{\cal J} |_{u_ Line 1782  The gradient $ \nabla _{u}{\cal J} |_{u_
1782  with the value of the cost function itself $ {\cal J}(u_{[k]}) $  with the value of the cost function itself $ {\cal J}(u_{[k]}) $
1783  at iteration step $ k $ serve  at iteration step $ k $ serve
1784  as input to a minimization routine (e.g. quasi-Newton method,  as input to a minimization routine (e.g. quasi-Newton method,
1785  conjugate gradient, ...) to compute an update in the  conjugate gradient, ... \cite{gil_lem:89})
1786    to compute an update in the
1787  control variable for iteration step $k+1$  control variable for iteration step $k+1$
1788  \[  \[
1789  u_{[k+1]} \, = \,  u_{[0]} \, + \, \Delta u_{[k+1]}  u_{[k+1]} \, = \,  u_{[0]} \, + \, \Delta u_{[k+1]}
# Line 1535  Tab. \ref{???} sketches the flow between Line 1797  Tab. \ref{???} sketches the flow between
1797  and the minimization routine.  and the minimization routine.
1798    
1799  \begin{eqnarray*}  \begin{eqnarray*}
1800  \footnotesize  \scriptsize
1801  \begin{array}{ccccc}  \begin{array}{ccccc}
1802  u_{[0]} \,\, ,  \,\, \Delta u_{[k]}    & ~ & ~ & ~ & ~ \\  u_{[0]} \,\, ,  \,\, \Delta u_{[k]}    & ~ & ~ & ~ & ~ \\
1803  {\Big\downarrow}  {\Big\downarrow}
# Line 1552  v_{[k]} = M \left( u_{[k]} \right) & Line 1814  v_{[k]} = M \left( u_{[k]} \right) &
1814  {\cal J}_{[k]} = {\cal J} \left( M \left( u_{[k]} \right) \right)} \\  {\cal J}_{[k]} = {\cal J} \left( M \left( u_{[k]} \right) \right)} \\
1815  \multicolumn{1}{|c}{~} & ~ & ~ & ~ & \multicolumn{1}{c|}{~} \\  \multicolumn{1}{|c}{~} & ~ & ~ & ~ & \multicolumn{1}{c|}{~} \\
1816  \hline  \hline
1817    \multicolumn{1}{|c}{~} & ~ & ~ & ~ & \multicolumn{1}{c|}{~}  \\
1818    \multicolumn{1}{|c}{~} & ~ & ~ & ~ & \multicolumn{1}{c|}{{\Big\downarrow}} \\
1819    \multicolumn{1}{|c}{~} & ~ & ~ & ~ & \multicolumn{1}{c|}{~}  \\
1820  \hline  \hline
1821  \multicolumn{1}{|c}{~} & ~ & ~ & ~ & \multicolumn{1}{c|}{~} \\  \multicolumn{1}{|c}{~} & ~ & ~ & ~ & \multicolumn{1}{c|}{~} \\
1822  \multicolumn{1}{|c}{  \multicolumn{1}{|c}{
1823  \nabla_u {\cal J}_{[k]} (\delta {\cal J}) =  \nabla_u {\cal J}_{[k]} (\delta {\cal J}) =
1824  T\!\!^{\ast} \cdot \nabla_v {\cal J} |_{v_{[k]}} (\delta {\cal J})} &  T^{\ast} \cdot \nabla_v {\cal J} |_{v_{[k]}} (\delta {\cal J})} &
1825  \stackrel{\bf adjoint}{\mathbf \longleftarrow} &  \stackrel{\bf adjoint}{\mathbf \longleftarrow} &
1826  ad \, v_{[k]} (\delta {\cal J}) =  ad \, v_{[k]} (\delta {\cal J}) =
1827  \nabla_v {\cal J} |_{v_{[k]}} (\delta {\cal J}) &  \nabla_v {\cal J} |_{v_{[k]}} (\delta {\cal J}) &
# Line 1565  ad \, v_{[k]} (\delta {\cal J}) = Line 1830  ad \, v_{[k]} (\delta {\cal J}) =
1830  \multicolumn{1}{|c}{~} & ~ & ~ & ~ & \multicolumn{1}{c|}{~} \\  \multicolumn{1}{|c}{~} & ~ & ~ & ~ & \multicolumn{1}{c|}{~} \\
1831  \hline  \hline
1832   ~ & ~ & ~ & ~ & ~ \\   ~ & ~ & ~ & ~ & ~ \\
1833  ~ & ~ &  \hspace*{15ex}{\Bigg\downarrow}  
1834  {\cal J}_{[k]} \qquad {\Bigg\downarrow}  \qquad \nabla_u {\cal J}_{[k]}  \quad {\cal J}_{[k]}, \quad \nabla_u {\cal J}_{[k]}
1835   & ~ & ~ \\   & ~ & ~ & ~ & ~ \\
1836   ~ & ~ & ~ & ~ & ~ \\   ~ & ~ & ~ & ~ & ~ \\
1837  \hline  \hline
1838  \multicolumn{1}{|c}{~} & ~ & ~ & ~ & \multicolumn{1}{c|}{~} \\  \multicolumn{1}{|c}{~} & ~ & ~ & ~ & \multicolumn{1}{c|}{~} \\
# Line 1595  The corresponding I/O flow looks as foll Line 1860  The corresponding I/O flow looks as foll
1860    
1861  \vspace*{0.5cm}  \vspace*{0.5cm}
1862    
1863    {\scriptsize
1864  \begin{tabular}{ccccc}  \begin{tabular}{ccccc}
1865  {\bf vector\_ctrl\_$<$k$>$ } & ~ & ~ & ~ & ~ \\  {\bf vector\_ctrl\_$<$k$>$ } & ~ & ~ & ~ & ~ \\
1866  {\big\downarrow}  & ~ & ~ & ~ & ~ \\  {\big\downarrow}  & ~ & ~ & ~ & ~ \\
# Line 1605  The corresponding I/O flow looks as foll Line 1871  The corresponding I/O flow looks as foll
1871  \cline{3-3}  \cline{3-3}
1872  \multicolumn{1}{l}{\bf xx\_theta0...$<$k$>$} & ~ &  \multicolumn{1}{l}{\bf xx\_theta0...$<$k$>$} & ~ &
1873  \multicolumn{1}{|c|}{~} & ~ & ~ \\  \multicolumn{1}{|c|}{~} & ~ & ~ \\
1874  \multicolumn{1}{l}{\bf xx\_salt0...$<$k$>$} & $\longrightarrow$ &  \multicolumn{1}{l}{\bf xx\_salt0...$<$k$>$} &
1875    $\stackrel{\mbox{read}}{\longrightarrow}$ &
1876  \multicolumn{1}{|c|}{forward integration} & ~ & ~ \\  \multicolumn{1}{|c|}{forward integration} & ~ & ~ \\
1877  \multicolumn{1}{l}{\bf \vdots} & ~ & \multicolumn{1}{|c|}{~}    \multicolumn{1}{l}{\bf \vdots} & ~ & \multicolumn{1}{|c|}{~}  
1878  & ~ & ~ \\  & ~ & ~ \\
1879  \cline{3-3}  \cline{3-3}
1880  ~ & ~ & ~ & ~ & ~ \\  ~ & ~ & $\downarrow$ & ~ & ~ \\
1881  \cline{3-3}  \cline{3-3}
1882  ~ & ~ &  ~ & ~ &
1883  \multicolumn{1}{|c|}{~} & ~ &  \multicolumn{1}{|c|}{~} & ~ &
1884  \multicolumn{1}{l}{\bf adxx\_theta0...$<$k$>$}  \\  \multicolumn{1}{l}{\bf adxx\_theta0...$<$k$>$}  \\
1885  ~ & ~ & \multicolumn{1}{|c|}{adjoint integration} &  ~ & ~ & \multicolumn{1}{|c|}{adjoint integration} &
1886  $\longrightarrow$ &  $\stackrel{\mbox{write}}{\longrightarrow}$ &
1887  \multicolumn{1}{l}{\bf adxx\_salt0...$<$k$>$} \\  \multicolumn{1}{l}{\bf adxx\_salt0...$<$k$>$} \\
1888  ~ & ~ & \multicolumn{1}{|c|}{~}    ~ & ~ & \multicolumn{1}{|c|}{~}  
1889  & ~ & \multicolumn{1}{l}{\bf \vdots} \\  & ~ & \multicolumn{1}{l}{\bf \vdots} \\
# Line 1628  $\longrightarrow$ & Line 1895  $\longrightarrow$ &
1895  ~ & ~ & ~ & ~ &  {\big\downarrow} \\  ~ & ~ & ~ & ~ &  {\big\downarrow} \\
1896  ~ & ~ & ~ & ~ &  {\bf vector\_grad\_$<$k$>$ } \\  ~ & ~ & ~ & ~ &  {\bf vector\_grad\_$<$k$>$ } \\
1897  \end{tabular}  \end{tabular}
1898    }
1899    
1900  \vspace*{0.5cm}  \vspace*{0.5cm}
1901    
1902    
1903  {\it ctrl\_unpack} reads in the updated control vector  {\it ctrl\_unpack} reads the updated control vector
1904  {\bf vector\_ctrl\_$<$k$>$}.  {\bf vector\_ctrl\_$<$k$>$}.
1905  It distributes the different control variables to  It distributes the different control variables to
1906  2-dim. and 3-dim. files {\it xx\_...$<$k$>$}.  2-dim. and 3-dim. files {\it xx\_...$<$k$>$}.
1907  During the forward integration the control variables  At the start of the forward integration the control variables
1908  are read from {\it xx\_...$<$k$>$}.  are read from {\it xx\_...$<$k$>$} and added to the
1909  Correspondingly, the adjoint fields are written  field.
1910    Correspondingly, at the end of the adjoint integration
1911    the adjoint fields are written
1912  to {\it adxx\_...$<$k$>$}, again via the active file routines.  to {\it adxx\_...$<$k$>$}, again via the active file routines.
1913  Finally, {\it ctrl\_pack} collects all adjoint field files  Finally, {\it ctrl\_pack} collects all adjoint files
1914  and writes them to the compressed vector file  and writes them to the compressed vector file
1915  {\bf vector\_grad\_$<$k$>$}.  {\bf vector\_grad\_$<$k$>$}.
1916    
# Line 1648  and writes them to the compressed vector Line 1918  and writes them to the compressed vector
1918    
1919    
1920    
1921  \subsection{Flow directives and adjoint support routines}  \subsection{Flow directives and adjoint support routines \label{section_flowdir}}
1922    
1923  \subsection{Store directives and checkpointing}  \subsection{Store directives and checkpointing \label{section_checkpointing}}
1924    
1925  \subsection{Gradient checks}  \subsection{Gradient checks \label{section_grdchk}}
1926    
1927  \subsection{Second derivative generation via TAMC}  \subsection{Second derivative generation via TAMC}
1928    

Legend:
Removed from v.1.1.1.1  
changed lines
  Added in v.1.7

  ViewVC Help
Powered by ViewVC 1.1.22