--- MITgcm_contrib/articles/ceaice/ceaice_adjoint.tex	2008/02/26 19:27:26	1.1
+++ MITgcm_contrib/articles/ceaice/ceaice_adjoint.tex	2008/03/25 22:04:31	1.3
@@ -3,80 +3,55 @@
 
 \subsection{The adjoint of MITsim}
 
-The ability to generate tangent linear and adjoint model components
-of the MITsim has been a main design task.
+
+The ability to generate tangent linear and adjoint components
+of a coupled ocean sea-ice system was one of the main drivers
+behind the MITsim development.
 For the ocean the adjoint capability has proven to be an
-invaluable tool for sensitivity analysis as well as state estimation.
-In short, the adjoint enables very efficient computation of the gradient
-of scalar-valued model diagnostics (called cost function or objective function)
-with respect to many model "variables".
-These variables can be two- or three-dimensional fields of initial 
-conditions, model parameters such as mixing coefficients, or
-time-varying surface or lateral (open) boundary conditions.
+invaluable tool for sensitivity analysis as well as state estimation,
+as evidenced by various adjoint-based studies
+(for a recent summary, see \cite{heim:08}).
+
+The adjoint model operator (ADM) is the transpose of the tangent linear 
+model operator (TLM)
+of the full (in general nonlinear) forward model, i.e. the MITsim.
+It enables very efficient computation of gradients
+of scalar-valued model diagnostics 
+(so-called cost function or objective function)
+with respect to many model inputs (so-called independent or control variables).
+These inputs can be two- or three-dimensional fields of initial 
+conditions of the ocean or sea-ice state, model parameters such as 
+mixing coefficients, or time-varying surface or lateral (open) boundary conditions.
 When combined, these variables span a potentially high-dimensional
 (e.g. O(10$^8$)) so-called control space. Performing parameter perturbations
 to assess model sensitivities quickly becomes prohibitive at these scales.
-Alternatively, (time-varying) sensitivities of the objective function 
-to any element of the  control space can be computed very efficiently in 
-one single adjoint 
+Alternatively, transient sensitivities of the objective function 
+to any element of the  control and model state space can be computed 
+very efficiently in  one single adjoint 
 model integration, provided an efficient adjoint model is available.
 
-[REFERENCES]
-
-
-The adjoint operator (ADM) is the transpose of the tangent linear operator (TLM)
-of the full (in general nonlinear) forward model, i.e. the MITsim.
-The TLM maps perturbations of elements of the control space
-(e.g. initial ice thickness distribution)
-via the model Jacobian
-to a perturbation in the objective function 
-(e.g. sea-ice export at the end of the integration interval).
-\textit{Tangent} linearity ensures that the derivatives are evaluated
-with respect to the underlying model trajectory at each point in time.
-This is crucial for nonlinear trajectories and the presence of different
-regimes (e.g. effect of the seaice growth term at or away from the
-freezing point of the ocean surface).
-Ensuring tangent linearity can be easily achieved by integrating
-the full model in sync with the TLM to provide the underlying model state.
-Ensuring \textit{tangent} adjoints is equally crucial, but much more
-difficult to achieve because of the reverse nature of the integration:
-the adjoint accumulates sensitivities backward in time,
-starting from a unit perturbation of the objective function.
-The adjoint model requires the model state in reverse order.
-This presents one of the major complications in deriving an
-exact, i.e. \textit{tangent} adjoint model.
-
-Following closely the development and maintenance of TLM and ADM
-components of the MITgcm we have relied heavily on the
+Following closely the development and maintenance of the
+TLM and ADM components of the MITgcm we have relied heavily on the
 autmomatic differentiation (AD) tool
 "Transformation of Algorithms in Fortran" (TAF)
-developed by Fastopt (Giering and Kaminski, 1998)
-to derive TLM and ADM code of the MITsim.
+developed by Fastopt \citep{gier-kami:98}.
+to derive TLM and ADM code of the MITsim
+(for details see \cite{maro-etal:99}, \cite{heim-etal:05}).
 Briefly, the nonlinear parent model is fed to the AD tool which produces 
 derivative code for the specified control space and objective function.
-Following this approach has (apart from its evident success)
-several advantages:
-(1) the adjoint model is the exact adjoint operator of the parent model,
-(2) the adjoint model can be kept up to date with respect to ongoing
-development of the parent model, and adjustments to the parent model
-to extend the automatically generated adjoint are incremental changes
-only, rather than extensive re-developments,
-(3) the parallel structure of the parent model is preserved 
-by the adjoint model, ensuring efficient use in high performance
-computing environments.
-
-Some initial code adjustments are required to support dependency analysis
-of the flow reversal and certain language limitations which may lead
-to irreducible flow graphs (e.g. GOTO statements).
-The problem of providing the required model state in reverse order
-at the time of evaluating nonlinear or conditional
-derivatives is solved via balancing
-storing vs. recomputation of the model state in a multi-level
-checkpointing loop.
-Again, an initial code adjustment is required to support TAFs 
-checkpointing capability.
-The code adjustments are sufficiently simple so as not to cause
-major limitations to the full nonlinear parent model.
+Apart from its evident success, advantages of this approach have been 
+pointed out, e.g. by \cite{gier-kami:98}.
+
+Many issues underlying the efficient exact adjoint sea-ice code generation
+are similar to those arising for the ocean model's adjoint.
+Linearizing the model around the exact nonlinear model trajectory,
+as we do, is a crucial aspect in the presence of different
+regimes (e.g. effect of the seaice growth term at or away from the
+freezing point of the ocean surface).
+Adjusting the (parent) model code to support the AD tool in
+providing exact and efficient adjoint code is the main initial work.
+This may be substantial for legacy code, but fairly straightforward
+when coding with "AD application in mind".
 Once in place, an adjoint model of a new model configuration
 may be derived in about 10 minutes.
 
@@ -93,13 +68,25 @@
 * approximate adjoints
 
 
-\subsection{An example: sensitivities of sea-ice export through Fram Strait}
+\subsection{An example: sensitivities of sea-ice export through
+the Lancaster and Jones Sound}
 
 We demonstrate the power of the adjoint method
-in the context of investigating sea-ice export sensitivities through Fram Strait
-(for details of this study see Heimbach et al., 2007).
-%\citep[for details of this study see][]{heimbach07}. %Heimbach et al., 2007).
-The domain chosen is a coarsened version of the Arctic face of the
+in the context of investigating sea-ice export sensitivities through 
+Lancaster and Jones Sound. The rationale for doing so is to complement
+the analysis of sea-ice dynamics in the presence of narrow straits.
+Lancaster Sound is one of the main outflow paths of sea-ice flowing
+through the Canadian Arctic Archipelago (CAA). 
+Export sensitivities reflect dominant
+pathways through the CAA as resolved by the model.
+Sensitivity maps can shed a very detailed light on various quantities
+affecting the sea-ice export (and thus the underlying pathways).
+Note that while the dominant circulation through Lancaster Sound is
+toward the East, there is a small Westward flow to the North,
+hugging the coast of Devon Island [ARE WE RESOLVING THIS?],
+see e.g. \cite{mell:02, mich-etal:06,muen-etal:06}.
+
+The model domain is a coarsened version of the Arctic face of the
 high-resolution cubed-sphere configuration of the ECCO2 project
 \citep[see][]{menemenlis05}. It covers the entire Arctic,
 extends into the North Pacific such as to cover the entire
@@ -112,23 +99,221 @@
 (benchmarks have been performed both on an SGI Altix as well as an 
 IBM SP5 at NASA/ARC).
 
-Following a 1-year spinup, the model has been integrated for four
-years between 1992 and 1995. It is forced using realistic 6-hourly
+Following a 3-year spinup, the model has been integrated for four
+years and five months between January 1989 and May 1993.
+It is forced using realistic 6-hourly
 NCEP/NCAR atmospheric state variables. Over the open ocean these are
 converted into air-sea fluxes via the bulk formulae of
 \citet{large04}.  Derivation of air-sea fluxes in the presence of
 sea-ice is handled by the ice model as described in \refsec{model}.
-The objective function chosen is sea-ice export through Fram Strait
-computed for December 1995.  The adjoint model computes sensitivities
-to sea-ice export back in time from 1995 to 1992 along this
+The objective function chosen is 
+sea-ice export through 
+Lancaster Sound at XX$^{\circ}$W
+averaged over an 8-month period between October 1992 and May 1993.  
+
+The adjoint model computes sensitivities
+to sea-ice export back in time from 1993 to 1989 along this
 trajectory.  In principle all adjoint model variable (i.e., Lagrange
-multipliers) of the coupled ocean/sea-ice model are available to
-analyze the transient sensitivity behaviour of the ocean and sea-ice
-state.  Over the open ocean, the adjoint of the bulk formula scheme
+multipliers) of the coupled ocean/sea-ice model
+as well as the surface atmospheric state are available to
+analyze the transient sensitivity behaviour.  
+Over the open ocean, the adjoint of the bulk formula scheme
 computes sensitivities to the time-varying atmospheric state.  Over
 ice-covered parts, the sea-ice adjoint converts surface ocean
 sensitivities to atmospheric sensitivities.
 
+DISCUSS FORWARD STATE, INCLUDING SOME NUMBERS ON SEA-ICE EXPORT
+
+\subsection{Sensitivities to the sea-ice state}
+
+\paragraph{Sensitivities to the sea-ice thickness}
+
+The most readily interpretable ice-export sensitivity is that
+to ice thickness, $\partial J / \partial heff$.
+Fig. XXX depcits transient $\partial J / \partial heff$ using free-slip
+(left column) and no-slip (right column) boundary conditions.
+Sensitivity snapshots are depicted for (from top to bottom) 
+12, 24, 36, and 48 months prior to May 2003.
+The dominant features are in accordance with expectations:
+
+(*)
+Dominant pattern (for the free-slip run) is that of positive sensitivities, i.e.
+a unit increase in sea-ice thickness in most places upstream
+of Lancaster Sound will increase sea-ice export through Lancaster Sound.
+The dominant pathway follows (backward in time) through Barrow Strait
+into Viscount Melville Sound, and from there trough M'Clure Strait
+into the Arctic Ocean (the "Northwest Passage"). 
+Secondary paths are Northward from
+Viscount Melville Sound through Byam Martin Channel into
+Prince Gustav Adolf Sea and through Penny Strait into MacLean Strait.
+
+(*)
+As expected, at any given time the
+region of influence is larger for the free-slip than no-slip simulation.
+For the no-slip run, the region of influence is confined, after four years,
+to just West of Barrow Strait (North of Prince of Wales Island),
+and to the South of Penny Strait.
+In contrast, sensitivities of the free-slip run extend
+all the way to the Arctic interior both to the West
+(M'Clure St.) and to the North (Ballantyne St., Prince Gustav Adolf Sea,
+Massey Sound).
+
+(*)
+sensitivities seem to spread out in "pulses" (seasonal cycle)
+[PLOT A TIME SERIES OF ADJheff in Barrow Strait)
+
+(*)
+The sensitivity in Baffin Bay are more complex.
+The pattern evolves along the Western boundary, connecting
+the Lancaster Sound Polynya, the Coburg Island Polynya, and the
+North Water Polynya, and reaches into Nares Strait and the Kennedy Channel.
+The sign of sensitivities has an oscillatory character
+[AT FREQUENCY OF SEASONAL CYCLE?].
+First, we need to establish whether forward perturbation runs
+corroborate the oscillatory behaviour.
+Then, several possible explanations:
+(i) connection established through Nares Strait throughflow
+which extends into Western boundary current in Northern Baffin Bay.
+(ii) sea-ice concentration there is seasonal, i.e. partly
+ice-free during the year. Seasonal cycle in sensitivity likely
+connected to ice-free vs. ice-covered parts of the year.
+Negative sensitivities can potentially be attributed
+to blocking of Lancaster Sound ice export by Western boundary ice
+in Baffin Bay.
+(iii) Alternatively to (ii), flow reversal in Lancaster Sound is a possibility
+(in reality there's a Northern counter current hugging the coast of
+Devon Island which we probably don't resolve).
+
+Remote control of Kennedy Channel on Lancaster Sound ice export
+seems a nice test for appropriateness of free-slip vs. no-slip BCs.
+
+\paragraph{Sensitivities to the sea-ice area}
+
+Fig. XXX depcits transient sea-ice export sensitivities
+to changes in sea-ice concentration
+ $\partial J / \partial area$ using free-slip
+(left column) and no-slip (right column) boundary conditions.
+Sensitivity snapshots are depicted for (from top to bottom) 
+12, 24, 36, and 48 months prior to May 2003.
+Contrary to the steady patterns seen for thickness sensitivities,
+the ice-concentration sensitivities exhibit a strong seasonal cycle
+in large parts of the domain (but synchronized on large scale).
+The following discussion is w.r.t. free-slip run.
+
+(*)
+Months, during which sensitivities are negative:
+\\
+0 to 5   Db=N/A, Dr=5 (May-Jan) \\
+10 to 17 Db=7, Dr=5 (Jul-Jan) \\
+22 to 29 Db=7, Dr=5 (Jul-Jan) \\
+34 to 41 Db=7, Dr=5 (Jul-Jan) \\
+46 to 49 D=N/A \\
+%
+These negative sensitivities seem to be connected to months
+during which main parts of the CAA are essentially entirely ice-covered.
+This means that increase in ice concentration during this period
+will likely reduce ice export due to blocking
+[NEED TO EXPLAIN WHY THIS IS NOT THE CASE FOR dJ/dHEFF].
+Only during periods where substantial parts of the CAA are
+ice free (i.e. sea-ice concentration is less than one in larger parts of
+the CAA) will an increase in ice-concentration increase ice export.
+
+(*)
+Sensitivities peak about 2-3 months before sign reversal, i.e.
+max. negative sensitivities are expected end of July
+[DOUBLE CHECK THIS].
+
+(*) 
+Peaks/bursts of sensitivities for months
+14-17, 19-21, 27-29, 30-33, 38-40, 42-45
+
+(*) 
+Spatial "anti-correlation" (in sign) between main sensitivity branch
+(essentially Northwest Passage and immediate connecting channels),
+and remote places.
+For example: month 20, 28, 31.5, 40, 43.
+The timings of max. sensitivity extent are similar between
+free-slip and no-slip run; and patterns are similar within CAA,
+but differ in the Arctic Ocean interior.
+
+(*) 
+Interesting (but real?) patterns in Arctic Ocean interior.
+
+\paragraph{Sensitivities to the sea-ice velocity}
+
+(*)
+Patterns of ADJuice at almost any point in time are rather complicated
+(in particular with respect to spatial structure of signs).
+Might warrant perturbation tests.
+Patterns of ADJvice, on the other hand, are more spatially coherent,
+but still hard to interpret (or even counter-intuitive
+in many places).
+
+(*)
+"Growth in extent of sensitivities" goes in clear pulses:
+almost no change between months: 0-5, 10-20, 24-32, 36-44
+These essentially correspond to months of
+
+
+\subsection{Sensitivities to the oceanic state}
+
+\paragraph{Sensitivities to theta}
+
+\textit{Sensitivities at the surface (z = 5 m)}
+
+(*)
+mabye redo with caxmax=0.02 or even 0.05
+
+(*)
+Core of negative sensitivities spreading through the CAA as 
+one might expect [TEST]:
+Increase in SST will decrease ice thickness and therefore ice export.
+
+(*)
+What's maybe unexpected is patterns of positive sensitivities
+at the fringes of the "core", e.g. in the Southern channels
+(Bellot St., Peel Sound, M'Clintock Channel), and to the North
+(initially MacLean St., Prince Gustav Adolf Sea, Hazen St.,
+then shifting Northward into the Arctic interior).
+
+(*)
+Marked sensitivity from the Arctic interior roughly along 60$^{\circ}$W
+propagating into Lincoln Sea, then
+entering Nares Strait and Smith Sound, periodically
+warming or cooling[???] the Lancaster Sound exit.
+
+\textit{Sensitivities at depth (z = 200 m)}
+
+(*)
+Negative sensitivities almost everywhere, as might be expected.
+
+(*)
+Sensitivity patterns between free-slip and no-slip BCs
+are quite similar, except in Lincoln Sea (North of Nares St),
+where the sign is reversed (but pattern remains similar).
+
+\paragraph{Sensitivities to salt}
+
+T.B.D.
+
+\paragraph{Sensitivities to velocity}
+
+T.B.D.
+
+\subsection{Sensitivities to the atmospheric state}
+
+\begin{itemize}
+%
+\item
+plot of ATEMP for 12, 24, 36, 48 months
+%
+\item
+plot of HEFF for 12, 24, 36, 48 months
+%
+\end{itemize}
+
+
+
 \reffig{4yradjheff}(a--d) depict sensitivities of sea-ice export
 through Fram Strait in December 1995 to changes in sea-ice thickness
 12, 24, 36, 48 months back in time. Corresponding sensitivities to
@@ -161,41 +346,14 @@
 \subfigure[{\footnotesize -24 months}]
 {\includegraphics*[width=0.44\linewidth]{\fpath/run_4yr_ADJheff_arc_lev1_tim145_cmax2.0E+02.eps}}
 }
-
-\centerline{
-\subfigure[{\footnotesize
--36 months}]
-{\includegraphics*[width=0.44\linewidth]{\fpath/run_4yr_ADJheff_arc_lev1_tim218_cmax2.0E+02.eps}}
-%
-\subfigure[{\footnotesize
--48 months}]
-{\includegraphics*[width=0.44\linewidth]{\fpath/run_4yr_ADJheff_arc_lev1_tim292_cmax2.0E+02.eps}}
-}
+%
 \caption{Sensitivity of sea-ice export through Fram Strait in December 2005 to
 sea-ice thickness at various prior times.
 \label{fig:4yradjheff}}
 \end{figure}
 
 
-\begin{figure}[t!]
-\centerline{
-\subfigure[{\footnotesize -12 months}]
-{\includegraphics*[width=0.44\linewidth]{\fpath/run_4yr_ADJtheta_arc_lev1_tim072_cmax5.0E+01.eps}}
-%\includegraphics*[width=.3\textwidth]{H_c.bin_res_100_lev1.pdf}
-%
-\subfigure[{\footnotesize -24 months}]
-{\includegraphics*[width=0.44\linewidth]{\fpath/run_4yr_ADJtheta_arc_lev1_tim145_cmax5.0E+01.eps}}
-}
-
-\centerline{
-\subfigure[{\footnotesize
--36 months}] 
-{\includegraphics*[width=0.44\linewidth]{\fpath/run_4yr_ADJtheta_arc_lev1_tim218_cmax5.0E+01.eps}}
-%
-\subfigure[{\footnotesize
--48 months}]
-{\includegraphics*[width=0.44\linewidth]{\fpath/run_4yr_ADJtheta_arc_lev1_tim292_cmax5.0E+01.eps}}
-}
-\caption{Same as \reffig{4yradjheff} but for sea surface temperature
-\label{fig:4yradjthetalev1}}
-\end{figure}
+%%% Local Variables: 
+%%% mode: latex
+%%% TeX-master: "ceaice"
+%%% End: 
