--- MITgcm_contrib/articles/ceaice/ceaice_adjoint.tex 2008/03/04 20:30:40 1.2 +++ 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,45 +346,13 @@ \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} -% -\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"