The MIT/SIO  OceanStateEstimation

|  Global State Estimation  |   The North Atlantic    |   The Indian Ocean   |

People

Publications

Report of  GODAE/WOCE Workshop on Ocean State  Estimation

This project is part of the NPACI ESS Enterprise .

Summary:

This project is directed at the exploitation of ocean data for the purpose of determining the general circulation of the ocean, on a global continuous basis. The central vehicle for this exploitation is a global adjoint which we have developed locally into a (nearly) complete system for global data assimilation and analysis. Our present focus is on the time-evolving global circulation as it emerges from WOCE hydrography and current data and from the TOPEX/POSEIDON altimeter observations globally and in regional higher-resolution models. Results will be used to understand, in particular, the oceanic heat and fresh water fluxes, their divergences, their dynamical causes and effects, as well as a variety of related issues connecting the oceanic circulation to climate variability. The ocean state estimation methodology is also intended for use with gravity data from the upcoming GRACE mission and when combined with that data should come close to providing the optimum marine gravity field estimates. Moreover, the model results are presently used by our collegues Mick Fellows and John Marshall to calculate oceanic biological fluxes and their relation on ocean color, as well as tracer and carbon cycles.

A primary objective of the World Ocean Circulation Experiment (WOCE) is to obtain an understanding of the absolute time-varying large-scale circulation of the world ocean, and its impact on climate. To meet this goal for a global turbulent fluid and the associated very large number of degrees of freedom one must exploit all available data in a quantitative fashion along with the information contained in the known dynamics. The latter consists of Newton's laws of motion plus thermodynamic relations, and are practically constituted in ocean general circulation models (OGCMs). If an OGCM can be forced to full consistency with a variety of global data sets, the resulting circulation estimate can be employed to study the consequences of the circulation and its temporal variability. Simultaneously one ought to be able to make quantitative estimates of the uncertainty of the results and their sensitivity to observational strategies -- important elements in determining what is known about climate change and of utmost importance in designing future observational programs to reduce the remaining uncertainty. From another point of view, the quantitative combination of an OGCM with observations represents an initialization of the model--an essential step should one wish to undertake climate forecasting, probably in conjunction with a coupled atmospheric model.

Our present focus is on the evolving global circulation as it emerges as a base-line solution from altimetric measurements alone. Precise and accurate TOPEX/POSEIDON (T/P) sea surface height observations are now available on a routine basis for 5+ years covering the period September 1992 through the present. Combined with a general circulation model, they carry unprecedented information about the large-scale circulation (see Wunsch and Stammer (1998) for a recent review of satellite altimetry). Information from WOCE in situ data such as XBTs, floats and the WOCE hydrography is being used first as independent information to test the results. These data will subsequently be included in the estimation process.

As a complementary effort, regional experiments are being run with enhanced resolution. These experiments are nested into the global estimates and use the T,S and horizontal flow field at the open boundaries as part of the control vector. They are being run in the Indian Ocean (PI: J. Marotzke) and in the Atlantic Ocean north of 30$^{\circ }$S as part of the regional WOCE AIMS efforts. Another effort being addressed at MIT is the generation of uncertainty estimates of results and their sensitivity to specific data (R. Giering and C. Wunsch).
 
 

Figure: The ocean circulation and overlaying atmospheric surface forcing fields: The upper part shows monthly mean fields of the net sea surface heat flux (positive into the ocean) and the surface wind stress (white vector field, in N/m²) as provided by the National Center for Environmental Prediction (NCEP) re-analysis for January 1993. The contour increment for the heat flux is 30 W/m² and the range covers -670 W/m² over the Kuroshio and about 215 W/m² in the southern hemisphere. 

 

Figure: The same fields as the above figure, but displayed over the Atlantic Ocean.

A formal discussion of the adjoint formalism is given in a report of a first global aplication of our procedure, using one year of TOPEX/POSEIDON data. The report can be found at http://puddle.mit.edu/~detlef/OSE/report_0/index.html . Here we can only give a brief summary and an status report on present work.

The Adjoint Optimization Formalism

In general terms, the adjoint model is employed to bring, in an extremely efficient way, the model into consistence with the observations, which the forward model by itself predicts. From the resulting misfit between the model and the observations the adjoint of the original model determines which uncertain model parameter must be adjusted and by how much so as to bring the forward computation into accord with the data. Changes in those fields -- often referred to as ``control'' terms -- are determined as a best-fit in an least-squares sense over the full data period. Control terms are usually the model initial conditions and the surface forcing fields, but can likewise include internal model parameters such as friction coefficients and eddy tracer transfers. It is essential to realize that all adjustment of control parameters are enforced to stay in a dynamically acceptable range and that otherwise the model equations are preserved (see the report for mode details).

Because of the large dimensions of geophysical applications, the adjoint model is usually solved iteratively. To this end, a first guess time history of the model state is obtained from a first guess estimate of the unknown controls, initial state, and prescribed boundary forcing by integrating forward model forward in time. The adjoint model is then integrated backward in time to produce estimates of the sensitivities (Lagrange multipliers). A standard optimizing search algorithm (here a limited memory quasi-Newton BFGS) can be used with this information to improve the previous estimates of the control state vector. Another forward-backward iteration is then made, the cycle continuing until convergence is obtained

The Forward and Adjoint Models

We work with two models: the forward model and its adjoint; both are required for the estimation procedure.

The forward component is a Ocean General Circulation Model which is based on the incompressible Navier-Stokes equations on a sphere under the Boussinesq approximation. It was developed recently by John Marshall and his group (in collaboration with Arvid of the Laboratory for Computer Science at MIT) as a tool to study the ocean general circulation over a broad range of scales and physical processes, and is specifically coded for optimal use of modern computer architectures [Shaw et al. 1997; Hill and Shaw 1997]. For the present purpose, we use the hydrostatic, free-surface version of the model which consists of conservation equations for horizontal and vertical momentum, volume, heat, and salt, and an equation of state which are solved using finite-volume techniques.

Generally the coding of the ``adjoint'' of a complex numerical model is extremely time consuming and difficult, comparable in effort to developing the forward code itself. However, the modern computer code of the forward model rendered it possible to obtain the adjoint component from the forward code in a semi-automatic way by using the Tangent Linear and Adjoint Model Compiler (TAMC) which was written by Ralf Giering [Giering and Kaminsky, 1998] while at the Max-Planck Institut für Meteorologie, Hamburg. In practice, this system of automatic adjoint code generation has proven to be extremely flexible because it permits easy regeneration of the adjoint code whenever a change in the forward code, including the objective function, is made. It specifically provides flexibility in adding observations and related additional constrains in the estimation procedure. Equally important, it provides an ``adjoint'' component which has the same high-level coding structure and computational optimality present in the forward code.

First Results

In a preliminary attempt, configured to test a real ocean state estimation system, the model was run with 2$^{\circ }$ horizontal resolution and 20 levels in the vertical. The model was then constrained by absolute T/P sea surface height data from the year period covering 1993, relative to the most recent geoid estimate (Lamoine et al., 1997) (the forward model was constrained separately in the mean and time-dependent components of surface pressure, thus separating the geoid error from the distinctly different error in the time-evolving components) as well as NCEP surface forcing fields (wind stress, heat and fresh water fluxes) and the mean Levitus et al (1994) hydrography. In our present examples, the cotrol fields include the initial potential temperature ($\theta$) and salinity (S) fields are modified, as well as the surface forcing fields from one year. In that configuration, there are 1.5 $\times 10^6$ elements in the control vector.

Formally, all this led to the following cost function:

$$\frac{1}{2}[({\overline{\mbox{\boldmath$\eta$ } }}-{\overline{\mb......({\mbox{\boldmath$\eta$ } ^{\prime }}-{\mbox{\boldmath$\eta$ }^{\prime }}_{tp})$$ J = (9)

$$+({\mbox{\boldmath$\tau$ } _{x}}-{\mbox{\boldmath$\tau$ } _{x}}_{......ldmath$\tau$ } _{y}}-{\mbox{\boldmath$\mbox{\boldmath$\tau$ }$ } _{y}}_{ncep})$$

$$+ ({{\mathbf{H}}_{Q}}-{{\mathbf{H}}_{Q}}_{ncep})^T{\mathbf{W}}_{......)^T{\mathbf{W}}_{H_{F}}({{\mathbf{H}}_{F}}-{{\mathbf{H}}_{F}}_{ncep}) \nonumber$$

$$+({\overline{{\mathbf{T}}}}-{\overline{{\mathbf{T}}}}_{lev})^T{\m......T{\mathbf{W}}_{S}({\overline{{\mathbf{S}}}}-{\overline{{\mathbf{S}}}}_{lev})],$$

in which $\overline{\left(.\right)}$ indicates averages over the one year estimation period, and all other terms represent 10-day averages of the model state and observations. At the present time, the model error is being set to zero except at its boundary grid points and in its initial $\theta$ and S conditions.

Only results from this preliminary run are shown here (all results obtained after the cost function was reduced close to the desired minimum value) which illustrate a time-evolving model state which differs from both the forward model and the data. A full presentation of the output can be seen at the web-site http://puddle.mit.edu/$\sim$detlef/global.html.

Figure 1 shows the optimized monthly mean fields of the surface elevation together with the velocity at 120m depth for October 1993. Apart from the missing eddy variability, major elements of the general circulation are present in both fields. In contrast to the pure T/P observations, various spurious elments in the data due to errors in the geoid are absent in the combined estimate. This is especially clear in the tropical Pacific and the North Atlantic where previously known inconsistencies in the geoid led to unacceptable features in the inferred ocean circulation (Stammer and Wunsch, 1994).

Figure 1: Estimates of a monthly mean sea surface height field and velocity field at 120 m depth from October 1993 as it results after optimization. Contour interval is 20 cm for the surface elevation, and a reference velocity of 20 cm/s is shown below the figure. Velocities smaller than 1 cm/s in length were omitted from the plots.

Changes in the estimated model state relative to the unconstrained model are illustrated in Fig. 2 from variations in the velocity and temperature fields taken from the same month at 120 m and 610 m depth, respectively. Various spatially coherent temperature changes on eddy-to-basin scales can be seen. Maximum amplitudes are associated with western boundary currents in the northern hemisphere and along the Antarctic Circumpolar Current. As an example, the path of the North Atlantic Current is characterized by increased temperatures relative to the unconstrained solution. Note also the wave-like pattern in the velocity corrections of the tropical Pacific and Indian Ocean and the enhanced boundary currents at many locations of the world ocean. T/P surface observations demand changes in the model state over the entire water column.

Figure 2: Change of monthly mean temperature and velocity fields in the estimated state relative to the unconstrained model. The difference fields are taken from October 1993 at (a) 120 m and (b) 610 m depth, respectively. Contour intervals are 0.5 $^{\circ }$C, and a reference velocity of 7 cm/s is shown below the figure.

The impact of the T/P data on estimates of oceanic transports are illustrated in Fig. 3 in terms of meridional heat transports from the Atlantic Ocean (Fig. 3a), the Pacific Ocean (Fig. 3b) and the Indian Ocean (Fig. 3c). Estimates are shown of the annual-average transports from the unconstrained run (blue lines) and the constrained estimate (red lines). In the North Atlantic, the T/P data inquire an increase in the model-alone estimate of northward heat transport from about 0.8 PW to almost 1.2 PW at 25$^{\circ }$N. Note that the estimates from the constrained state are consistent with previous results from Macdonald and Wunsch (1996) and their errors estimates which are shown as open green circles and error bars. An exception can be found at 10$^{\circ }$N, where the model estimate is significantly lower than the former result. The basis for this ``adjustment'' is a strengthening of the vertical overturning flow field with an up to 1 cm/s stronger and about 1$^{\circ }$C warmer poleward flow in the Gulf Stream, and an enhanced and colder return flow at depth. It should be noted that changes in the model state due to T/P surface observations are enforced over the entire water column.

Figure 3: Time-mean meridional heat transport (in 10¹⁵ W) estimated from the last 6 month of the constrained (solid lines) and unconstrained model (dashed lines) for the Atlantic Ocean, the Pacific, and the Indian Ocean, respectively. Red bars represent the rms variability of the transports estimated over individual 10-day periods. Also shown as open circle and error bars are the heat flux estimates and their uncertainties obtained by Macdonald and Wunsch (1996) in the Atlantic from their box-inversion of the global ocean.

To indicate the degree to which the external time-dependent forcing fields are modified during the assimilation, Fig. 4 shows a typical example of the adjustment in the 10-day averages of heat and fresh water fluxes representing the T/P repeat cycle 21 (early September 1993). Variations in the NCEP heat fluxes are of the order of

Figure 4: Changes in (a) the 10-day averaged fresh water flux field and (b) the flux relative to the NCEP fields as they emerge from the optimization for the T/P repeat cycle 21 (early September, 1993). Contour intervals in (a) and (b) are 1 cm/year and 2 W m^-2, respectively.

$\pm$ 20 W m^-2 with maximum changes along the boundary currents in the norther hemisphere. Note the pronounced differences in the adjustment of the fresh water fluxes which are larges in the tropical oceans. Changes of similar amplitudes occur in all remaining repeat cycle. The wind stress fields are likewise adjusted to best fit the model to the observed TOPEX/POSEIDON absolute sea surface height observations. Overall, the amplitudes of changes in the external forcing fields are consistent with accepted uncertainties in the meteorological analyses.

Work in Progress:

Computations now underway are directed at greatly extending and improving the preliminary computation. Improvements include an increased model resolution to 1$^{\circ }$, an extended estimation period to 6 years (1992 through 1997), as well as more acurate representations of straits and sills in the topography. Equally important, a complete mixed layer model (Large et al., 1994) and an eddy parameterization (Gent and McWilliams, 1990) have been incorporated into the adjoint model. Moreover, a full (non-diagonal) geoid error covariance matrix is being used as well as other improvements in the remaining weights.

Data now include the absolute and time-varying T/P data from October 1992 through December 1997, SSH anomalies from the ERS-1 and ERS-2 satellites, monthly mean SST data (Reynolds and Smith, 1994), time-varying NCEP Re-Analysis fluxes of momentum, heat and freshwater, and NSCAT estimates of wind stress errors are being employed. Monthly means of the model state are required to remain within assigned bounds of the monthly mean Levitus et al. (1994) climatology. Note that the control vector now contains 8 million elements and the costfunction has the form:

$$\begin{eqnarray*}J = \frac{1}{2} [ ({\bf\overline{\zeta}} - {\bf\overline{\z...... SST})^{T} W_{SST} ({\bf T} - {\bf SST}). \\ \nonumber\nonumber\end{eqnarray*}$$ The set of computations now underway is still regarded as preliminary. The model resolution will eventually be greatly increased and as much as possible of the entire WOCE data set will be used. The combination of the OGCM with data, each appropriately weighted by its uncertainty estimates, should provide a useful basic description of the ocean circulation and its variability. Among many possible applications, one expects the production of much better estimates of property flux divergences, and the study of oceanic biogeochemical cycles which are dependent upon the circulation.

First examples of the GCM coupled to chemical and biological models are presented by M. Follows.
 
 
 

References

1.

Gent, P. R., and J. C. McWilliams, Isopycnal mixing in ocean models.

J. Phys. Oceanogr., 20, 150-155, 1990.

2.

Giering, R., and T. Kaminsky

Recipes for adjoint code construction
ACM Trans. Math. Software, in press, 1997.

3.

Large, W. G., J.C. McWilliams, and S.C. Doney, Oceanic vertical mixing: a review and a model with nonlocal boundary layer parameterization, Rev. Geophys., 32, 363-403, 1994.

4.

Lemoine, F., and 17 others, The development of the NASA GSFC and NIMA Joint Geopotential Model. In Gravity, Geoid and Marine Geodesy, International Association of Geodesy Symposia, 117, ed. Segawa et al. Springer-Verlag., Berlin Heidelberg, 1997.

5.

Levitus, S., R. Burgett, and T. Boyer,

World Ocean Atlas 1994, vol. 3, Salinity, and vol. 4, Temperature,
NOAA Atlas NESDIS 3 & 4, U.S. Dep. of Comm., Washington, D.C., 1994.

6.

Macdonald, A., and C. Wunsch, The global ocean circulation and heat flux.

Nature, 382, 436-439, 1995.

7.

Marotzke, J., Q. K. Zhang, R. Giering, D. Stammer, C. N. Hill, and T. Lee,

The linearization and adjoint of the MIT general circulation model,
in preparation, 1998.

8.

Marshall, J., A. Adcroft, C. Hill, L. Perelman, and C. Heisey,

A finite-volume, incompressible navier-stokes model for studies of the ocean on parallel computers.
J. Geophys. Res., 5753-5766, 1997a.

9.

Marshall, J., C. Hill, L. Perelman, and A. Adcroft,

Hydrostatic, quasi-hydrostatic and non-hydrostatic ocean modeling.
J. Geophys. Res., 5733-5752, 1997b.

10.

Reynolds, R. W., and T. M. Smith, Improved global sea surface temperature analyses using optimum interpolation, J. of Climate, 7, 929-948, 1994.

11.

Stammer, D. and C. Wunsch, Preliminary assessment of the accuracy and precision of TOPEX/Poseidon altimeter data with respect to the large scale ocean circulation, J. Geophys. Res., 99, 24000 - 25500, 1994.

12.

Stammer, D., C. Wunsch, R. Giering, Q. Zhang, J. Marotzke, J. Marshall, and C.N. Hill, The global ocean circulation estimated from TOPEX/POSEIDON altimetry and the MIT general circulation model, MIT Center for Global Change Science, Report 49, 1997.

13.

Wunsch, C., and D. Stammer, Satellite Altimetery, the Marine Geoid and the Oceanic General Circulation, Annual Reviews of Earth and Planetary Sciences, in press, 1998.

1998-05-08