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