A comparison of global ocean general circulation model solutions obtained with synchronous and accelerated integration methods

A comparison of global ocean general circulation model solutions obtained with synchronous and accelerated integration methods

Ocean Modelling 7 (2004) 323–341 www.elsevier.com/locate/ocemod A comparison of global ocean general circulation model solutions obtained with synchr...

2MB Sizes 38 Downloads 111 Views

Ocean Modelling 7 (2004) 323–341 www.elsevier.com/locate/ocemod

A comparison of global ocean general circulation model solutions obtained with synchronous and accelerated integration methods Gokhan Danabasoglu

*

National Center for Atmospheric Research, P.O. Box 3000, Boulder, CO 80307, USA Available online 13 November 2003

Abstract A 10 000-year synchronous control integration using a comprehensive ocean general circulation model (OGCM) subject to realistic, time-dependent forcing in a global domain is performed to quantitatively determine how well its solution is reproduced by two accelerated, but otherwise identical, equilibrium integrations. To our knowledge, this is the first long synchronous integration of such an OGCM. The two accelerated cases use tracer time steps increasing with depth and unequal momentum and tracer time steps, respectively. After accelerated equilibration, these cases are integrated further synchronously to achieve synchronous equilibration. The equilibration time is defined as the time when the annual- and global-mean potential temperature trend is below 105 C year1 . The synchronous control integration achieves equilibrium in about 3500 years. The accelerated equilibration times are 685 and 3000 surface years for the two cases, respectively. Their synchronous extensions require about 100 and 700 surface years. The accelerated solutions do differ from each other and from those of the control experiment. However, for practical purposes, many aspects of the accelerated and control solutions are more similar than otherwise. Because of severe non-conservation issues associated with tracer time-step variations with depth, this technique is not recommended. Instead, unequal momentum and tracer time steps should be used. Compared to the control case, this method of acceleration provides a factor of 2.5 reduction in the computational cost (accelerated + synchronous equilibration) with the possibility of further reductions. Any accelerated integration must be followed by a synchronous extension to recover the correct seasonal cycle and to eliminate any possible oscillatory behavior present in the accelerated phase. If the incremental gains with further synchronous extensions are judged to be minimal or of minor significance, shorter (<700 years) integrations may certainly suffice. Among others, the equilibration times are affected by the choice of surface forcing method.  2003 Elsevier Ltd. All rights reserved.

*

Tel.: +1-303-497-1604; fax: +1-303-497-1700. E-mail address: [email protected] (G. Danabasoglu).

1463-5003/$ - see front matter  2003 Elsevier Ltd. All rights reserved. doi:10.1016/j.ocemod.2003.10.001

324

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

Keywords: Synchronous equilibrium; Accelerated equilibrium; Global ocean general circulation model

1. Introduction The abyssal ocean diffusive processes with time scales of many thousand years ultimately control the approach to equilibrium of oceanic general circulation models (OGCMs). Despite ever increasing computational capabilities, particularly overall speed, obtaining equilibrium solutions is still practically unaffordable even for non-eddy-resolving models commonly used in climate system studies. Therefore, accelerated integration methods retain their wide usage in the oceanographic community to obtain equilibrium solutions (e.g. Large et al., 1997; Wood, 1998). Arguably, the most commonly applied technique is due to Bryan et al. (1975) in which the model tracer time step, also followed by the model calendar, is much longer than the momentum time step (order days vs. order hours/minutes). This approach was later extended by Bryan and Lewis (1979) to allow even longer tracer time steps with depth to exploit increased advective numerical stability limit due to diminishing current speeds in the deep ocean. Such a distorted physics technique has proven useful when an equilibrium baseline solution is sought without any interest in the transient behavior of the model. Detailed analyses of this method are presented in Bryan (1984) and Killworth et al. (1984), considering its effects on baroclinic instability, Rossby waves, and gravity wave dispersion and numerical stability, respectively. Clearly, the governing equations for the original and accelerated systems are identical only at steady state, and the primary assumption with time-independent forcing is the uniqueness of the final steady state solution. The usefulness of the distorted physics technique with unsteady forcing is explored by Danabasoglu et al. (1996). Without actually performing a companion synchronous integration (in which all model time steps are equal), the method is shown to work provided that the accelerated integration is followed by a short synchronous period to recover the proper seasonal cycle. This key issue is also recently addressed in Huang and Pedlosky (2002), where the lowest mode is shown to be no longer strictly barotropic when accelerated, leading to severely contaminated seasonal cycles. Using a series of OGCM integrations, including a companion synchronous one, Wang (2001) shows that although the choice of the acceleration parameters does affect the mean states of the accelerated equilibrium, these mean states are not too far from either each other or the synchronous mean state for practical purposes. The study also restresses the importance of a synchronous period following an accelerated integration. Although this work is the first in which a synchronous equilibrium solution is obtained for an OGCM, the model configuration is very simplistic, i.e. a box-type domain with simple sub-grid-scale parameterizations and without salinity and wind forcing. Furthermore, temperature forcing is of restoring type. In the present study, we revisit the accelerated integration technique one more time, using a state-of-the-art OGCM with the most recent parameterizations. A 10 000-year synchronous integration is performed subject to realistic, time-dependent forcing in a global domain. To our knowledge, this is the first, long synchronous integration with such a model. We compare two accelerated, but otherwise identical, integrations to this control case and to each other to determine if these acceleration methods do truly work. We also perform several multi-century

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

325

synchronous continuation integrations for the accelerated cases to assess their synchronous adjustments and to document any changes in their seasonal cycles and mean states. The equilibration times and the number of model time steps to get there are computed for each acceleratedsynchronous integration combination to measure computational gains. A brief description of the model along with some forcing details are given in Section 2. Section 3 describes our integration strategy and cases. The results are presented in Section 4, and some conclusions are given in Section 5.

2. Model description The oceanic model is the ocean component of the second generation of the National Center for Atmospheric Research (NCAR) Community Climate System Model (CCSM2). It is a state-ofthe-art, Bryan–Cox type (Bryan, 1969), level-coordinate model based on the Parallel Ocean Program (POP 1.4) of the Los Alamos National Laboratory (Smith et al., 1992). The model solves the primitive equations in general orthogonal coordinates in the horizontal (Smith et al., 1995) with the hydrostatic and Boussinesq approximations. An implicit free surface formulation is used for the barotropic equations (Dukowicz and Smith, 1994). Although the surface layer thickness is allowed to vary, we opt not to use the natural freshwater fluxes. Instead, these fluxes are converted to virtual salt fluxes using a reference salinity, and the globally integrated ocean volume remains constant. Only a summary of the ocean model physics and parameters is given here and further details can be found in Smith and Gent (2002). In order to complete this study in a reasonable time frame and to reduce the burden on our computational resources, the coarser, ·3 resolution version of the model is used. The domain is global. The Bering Strait and Hudson Bay are open, but the Gibraltar Strait and the Red Sea Outflow are closed. The grid is in spherical coordinates in the Southern Hemisphere. In the Northern Hemisphere, the North Pole is displaced into Greenland at 40W and 77N. The zonal resolution is constant at 3.6. The meridional resolution is variable with 0.86 at the Equator, monotonically increasing to 1.85 at 31S and staying constant further south. In the Northern Hemisphere high latitudes, the minimum and maximum meridional resolutions are about 1.3 and 2.2, occurring in the Atlantic and Pacific Oceans, respectively. This horizontal resolution results in a 100 (zonal) · 116 (meridional) grid. There are 25 vertical levels, monotonically increasing from 12 m near the surface to 450 m in the abyss. The minimum and maximum depths are 49 and 5000 m, respectively. The model tracer equations use the skew-flux form (Griffies, 1998) of the Gent and McWilliams (1990) isopycnal transport parameterization with a mixing coefficient of 800 m2 s1 . A general form of the anisotropic horizontal viscosity formulation due to Smith and McWilliams (2003) is implemented in the momentum equations. Its parameters are chosen following Large et al. (2001) with 1000 m2 s1 as the minimum (subject to numerical diffusive stability) background horizontal viscosity. The vertical mixing coefficients are determined from the KPP scheme of Large et al. (1994). In the ocean interior, the background internal wave mixing diffusivity varies in the vertical from 0.2 · 104 m2 s1 near the surface to 1.0 · 104 m2 s1 in the deep ocean. The increase in diffusivity occurs around 1000-m depth, as a crude representation of the enhanced deep vertical

326

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

mixing observed over rough topography (Ledwell et al., 2000). The background vertical viscosity obeys the same form, but with a turbulent Prandtl number of 10. Similar to other OGCMs, POP uses the Leap–Frog time-stepping method. The associated timesplitting error is eliminated via a time-averaging step in which the model calendar is advanced by one half the surface tracer time step. The model is driven by the net fluxes of heat, salt (freshwater), and momentum computed using the bulk forcing scheme described in Large et al. (1997). The same forcing method was also adopted in a recent inter-annual variability study by Doney et al. (2003). In the present work, a 3year cycle of atmospheric data sets covering the period 1991–1993 is used. Thus, the monthly mean International Satellite Cloud Climatology Project (ISCCP) solar radiation (Bishop and Rossow, 1991; Bishop et al., 1997) and cloud cover (Rossow and Schiffer, 1991) data sets include the most recent updates. To eliminate the low bias in the cloud cover data, it is further modified at high latitudes following the Hahn et al. (1987) observations. The monthly mean precipitation data represent a blending of the Microwave Sounding Unit (MSU) (Spencer, 1993) and Xie and Arkin (1996) data sets (see Doney et al., 2003 for blending details). The remaining atmospheric fields (winds, temperature, humidity, and density) are based on the 6-hourly NCEP/NCAR reanalysis data (Kalnay et al., 1996). A comparison of the NCEP/NCAR air temperatures to the available station data along the Antarctic coast shows an excessive cold bias in the NCEP/NCAR data (F. Bryan, personal communication). This bias is eliminated by limiting the minimum air temperatures south of 60S to roughly match the station data. When observed SSTs are used with these data sets, the net heat flux shows a positive bias. This is largely eliminated by subtracting uniformly 4 W m2 from the longwave-down heat flux component. The bulk forcing formulation uses the evolving model surface temperature. All the atmospheric data sets are on a nominal 2-grid and are interpolated to the ocean tracer points. The daily-means of the fluxes drive the OGCM. As described in Large and Yeager (2003), a climatological river runoff distribution from the continents contributes to the surface salt fluxes. Following Large et al. (1997), a precipitation correction factor is computed for each year based on the change in the global-mean salinity during that year. When the tracer time step increases with depth, the normalizing time for the tendency computation does accordingly increase. This global factor is then used to multiply precipitation and runoff fluxes to partially balance the evaporation for the next year. Due to the lack of any appreciable feedback between the model surface salinities and salt fluxes, any flux errors can lead to unbounded local salinity trends. Therefore, as a means for some local constraint, a local, weak salinity restoring is applied globally using a blending of the monthly mean salinity climatologies from Levitus et al. (1998) and Steele et al. (2001). The salinity enhancements along the Antarctic coast due to Doney and Hecht (2002) are also included in this blended climatology. The coefficient for this weak restoring is 11.5 mg m2 s1 psu1 , corresponding to a 1-year time scale over 12 m. Because the global-mean of this flux is subtracted at each time step, it does not contribute to the global salt budget. The surface salt fluxes over the land-locked, marginal seas are balanced over each marginal sea individually to prevent any unbounded, local salinity trends in these regions. Here, any deficit/excess is taken from/given to a designated open ocean region. This can be interpreted as an approximation to include some of the outflow effects of these seas which are otherwise closed. At high latitudes, the daily sea-ice concentration distributions from the Special Sensor Microwave/Imager (SSM/I) data set (Comiso, 1999) are used to determine the ice extent. The

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

327

usual air-sea fluxes are computed for open water and the wind stress passes directly through the sea-ice. No restoring boundary conditions are used under ice covered regions. Instead, we apply the freshwater fluxes diagnosed monthly from a companion coupled ocean–ice integration. In a few regions, there is no observed ice, but the ocean needs to form sea-ice to keep its surface temperature from falling below freezing. In such instances, the surface water is heated to the freezing point with a corresponding increase in salt. This heat flux is the only form of any additional heat flux applied under ice. This kind of ice is assumed to be 1-m thick for the purpose of computing a sea-ice concentration. It is not transported, and is melted locally before the surface temperature can rise above the freezing temperature of )1.8 C.

3. Acceleration details and cases Following Bryan (1984), the momentum and tracer equations for the distorted system are ou oðh; SÞ ¼ F; cðzÞ ¼ G: ot ot Here, u is the horizontal velocity vector, h is potential temperature, S is salinity, and t is time. F and G represent the rest of the terms in these equations. Finally, a and c, which can have a vertical, z, variation, are the acceleration factors. When a ¼ c ¼ 1, the original, synchronous equations are recovered. A list of the numerical experiments is given in Table 1. The 10 000-year synchronous control integration is designated as SYNC. In the first accelerated integration AðzÞ, a ¼ 1 but the tracer time step varies in z. Here, we use c ¼ 1 between 0- and 1099-m and c ¼ 0:02 between 2877- and 5000-m depth. In between, c has a linear variation. As discussed in Danabasoglu et al. (1996), these changes in the tracer time step result in tracer non-conservation. Therefore, the changes are restricted to intermediate and deep levels where the vertical fluxes are expected to be relatively small. In AC, the tracer time step is constant in z, but a ¼ 5. So, the momentum time step is five times smaller than the tracer one. This particular choice of a for AC is artificially dictated by the a

Table 1 List of numerical experiments Case

Initial condition (surface years)

Integration length (surface/ deep years)

Equilibration time (surface years)

Equilibration time-step (·106 )

DtT (surface/ deep in s)

DtU (s)

SYNC AðzÞ S_AðzÞ_450 S_AðzÞ_700 AC S_AC_3 S_AC_10

Levitus/rest Levitus/rest Year 450 of AðzÞ Year 700 of AðzÞ Levitus/rest Year 3000 of AC Year 10 000 of AC

10 000 700/35 000 400 500 10 000 800 800

3500 685 (450) 200 100 3000 700 750

16.6 6.5 (4.3) 0.9 0.5 3.3 3.3 3.5

6912 3456/172 800 6912 6912 34 560 6912 6912

6912 3456 6912 6912 6912 6912 6912

DtT and DtU are the tracer and momentum time steps, respectively. Equilibration time-step refers to the number of time steps to reach the equilibration time using the surface DtT . In AðzÞ, the numbers in parentheses represent alternative values based on the global-mean time. See text for the definitions.

328

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

CCSM2 requirement of at least daily coupling. This is further compounded by some POP model details, e.g. one averaging time step per coupling interval. Therefore, AC does not necessarily take full advantage of longer tracer time steps. SYNC, AðzÞ, and AC are all initialized with the January-mean Levitus et al. (1998) h and S and state of rest. Four additional, multi-century synchronous integrations are performed to study the synchronous adjustment following the accelerated integrations. S_AðzÞ_450 and S_AðzÞ_700 start from year 450 and 700 of AðzÞ, respectively. Similarly, S_AC_3 and S_AC_10 are the synchronous continuations of AC, starting at year 3000 and 10 000, respectively.

4. Results 4.1. Equilibration times t

The equilibration time (Table 1) is defined as the time when the annual- and global-mean hh i trend is below 105 C year1 , excluding any early transient behavior. This is a rather stringent criterion on the equilibration time for a global OGCM with realistic forcing. Indeed, even after a t fairly monotonic decrease below the given thresh-hold, hh i trend may occasionally be much higher. Our estimates, therefore, are somewhat subjective, but conservative. The trend is determined using t

t

hh iyear þ n  hh iyear ; nhti where hti ¼ 1 year is the global-mean time for all estimates except one estimate for AðzÞ. Here, we use hti ¼ 25:1 years to reflect the depth acceleration effects, i.e. 1 surface year equals 25.1 globalmean years. Due to the 3-year forcing cycle, we use n ¼ 3 in all cases except AC. In AC, we set n ¼ 66, t because a 66-year, regular oscillation with a hh i amplitude of 3 · 103 C develops after about 800 years (e.g. see Figs. 2 and 3). The associated anomalies start in the northern North Atlantic with local h and S amplitudes in excess of 1 C and 1 psu, respectively, in the upper ocean. The anomalies propagate southward along the American continent, decaying significantly when they reach the Southern Ocean. Therefore, the oscillations are more prominent in the Atlantic Ocean related fields. Unlike Wood (1998), there is no wave-like variability in the Southern Ocean associated with baroclinic instability (Bryan, 1984). Unfortunately, any further analysis of the present solutions does not easily indicate the culprit among many possibilities including surface forcing and the value of the acceleration factor (see also Wood, 1998 for other causes). We do not use S to determine the equilibration time. The surface forcing method described in Section 2 maintains the initial global-mean S value when the tracer time step has no z variation. Therefore, SYNC and all AC cases show no S trends (Fig. 1). For our purposes, the equilibration time-step count, representing the actual cost of an integration, is the most relevant information in Table 1. Here, because the model equations are integrated as usual, the time-averaging time steps are counted as full time steps. The table shows that the SYNC equilibrium is achieved in about 3500 years, taking 16.6 · 106 time steps. Based on

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

329

Fig. 1. Time series of the annual-mean potential temperature and salinity. Global and 3268-5000-m means are plotted. The lines are based on mean values for every third year. In all time series figures, the AðzÞ lines are plotted after multiplying the integration times by the global-mean time, 25.1. The black lines represent the synchronous extensions of the accelerated integrations.

surface years, AðzÞ and AC take 6.5 · 106 and 3.3 · 106 time steps, corresponding to 685 and 3000 years, respectively, for equilibration. Their synchronous extensions, S_AðzÞ_700 and S_AC_3, need additional 0.5 · 106 and 3.3 · 106 steps. Thus, the suite of AðzÞ and AC integration methods require approximately 7.0 · 106 and 6.6 · 106 steps, respectively. These represent about a factor of 2.5 reduction in the computational cost of an equilibrium solution compared to SYNC. An alternative count based on branching from the AðzÞ case earlier is computed as 5.2 · 106 using AðzÞ and S_AðzÞ_450. Finally, S_AC_10 is included to assess any effects of long, diffusive time scales. Even such, relatively long integrations are still much cheaper than a comparable length SYNC integration (14.5 · 106 vs. 47.4 · 106 time steps). The synchronous equilibration times in Table 1 for S_AC cases appear to be independent of the branching time. In contrast, S_AðzÞ cases require a shorter synchronous integration after a longer accelerated integration.

330

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

4.2. Time series We first consider the evolutions of h and S in Fig. 1. In addition to the global-mean time series, the volume-means for 3268–5000-m depth are included to represent the abyssal trends. Also, in all time series plots, AðzÞ line is plotted after multiplying its time axis by 25.1. The figure shows that t the global hh i values for all synchronous extensions following accelerated integrations are within t less than 0.01 C of the SYNC mean of 4.13 C. The AC global hh i time series (along with the t abyssal hS i below) definitely show that particularly this type of accelerated integration, in which it is more than 0.12 C warmer than in SYNC, must be followed by a synchronous integration. Both t t SYNC and AC reveal minor drifts in the deep ocean hh i. The AC deep hh i stays 0.05 C warmer than in SYNC even after the synchronous extensions. In contrast, the S_AðzÞ cases appear to t approach the SYNC deep-mean faster. However, neither S_AðzÞ nor S_AC deep hh is are as close to the SYNC deep-mean as their global-means to the SYNC one. t As noted earlier, if the accelerated integration technique is conservative, hS i changes very little due to our surface forcing method. The extent of non-conservation in AðzÞ, particularly obvious in S, is quite detrimental, showing a difference of more than 0.08 psu compared to either SYNC or AC. Similar magnitude differences are also present in the abyssal time series. The deep 0.04-psu fresh bias existing in AC vanishes upon switching to synchronous integrations. Fig. 2 shows the time series for the maximum (in [30–60N, 500–1200 m]) and minimum (in [60S–30N, 2000–5000 m]) global meridional overturning circulation computed using the Eulerian-mean velocity. The SYNC maximum and minimum transports are about 17.7 and 19.2 Sv, respectively. Both S_AðzÞ cases reproduce the same maximum transport as SYNC. The S_AC maximum transports are only slightly (0.3 Sv) higher than in SYNC. For the minimum transport, all of the synchronous extensions are within 0.2 Sv of the SYNC one. Note that the AC oscillation amplitude is about 1.7 and 0.2 Sv, respectively, for the maximum and minimum transports.

Fig. 2. Time series of the annual-mean global maximum and minimum meridional overturning circulations (MOCs). The maximum and minimum transports are searched in [30–60N, 500–1200 m] and [60S–30N, 2000–5000 m], respectively. The black lines represent the synchronous extensions of the accelerated integrations. Also, 1 Sv ” 106 m3 s1 .

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

331

Fig. 3. Annual-mean transport time series for the Antarctic circumpolar current (ACC) at Drake Passage, Indonesian Throughflow, Florida Strait, and Mozambique Channel. The black lines represent the synchronous extensions of the accelerated integrations.

The vertically integrated (barotropic) mass transports are given in Fig. 3 for the Antarctic circumpolar current (ACC) at Drake Passage, Indonesian Throughflow, Florida Strait, and Mozambique Channel. All transports show that the synchronous integrations succeeding the accelerated segments do produce the SYNC transports mostly within a few tenths of a Sv. The biggest difference is about 1 Sv for S_AðzÞ cases in the ACC. The amplitude of the AC oscillation is the largest (but still <1 Sv) in the Florida Strait transport. 4.3. Equilibrium solutions We first perform a quantitative assessment of how similar the solutions for synchronous equilibria following AðzÞ and AC are to those of SYNC. For this purpose, the global-mean rootmean-square (rms) differences for h and S (i.e. S_AðzÞ)SYNC and S_AC)SYNC) are computed

332

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

using 3-year means every 100 years (Fig. 4). The years 9997–9999 mean is used for SYNC. The rms h differences monotonically diminish with continued synchronous integrations. However, their trends become much smaller with time, suggesting that significantly longer integrations are needed for similar marginal reductions. The year 0 values actually show the rms differences using the last 3 years of the corresponding accelerated integrations, displaying that the AðzÞ method results in h distributions much closer to the SYNC ones than the AC method. At their respective equilibria, all synchronous extentions show 0.06–0.07 C rms differences. Note also that further accelerated integration with AðzÞ reduces both the accelerated rms difference and the length of the synchronous integration for comparable rms differences. In contrast with the very favorable rms h differences for AðzÞ cases, the corresponding rms S difference plots, again, reveal how significant the associated non-conservation issues are for this method of acceleration. Both S_AðzÞ cases have a difference of about 0.09 psu, unaffected by extended integrations. This is of course a significant concern for water mass properties. In contrast, the rms S differences diminish monotonically with time for S_AC cases with about 0.01-psu difference at equilibrium. Fig. 4 shows that half of the approach (half-life) towards the SYNC solution occurs in less than 200 years for all cases. This time scale is consistent with most of the time series given in Figs. 1–3, t following an initial, fast transient adjustment. The exceptions are the global hh i and ACC time series for AC where the half-life is a few decades. The zonal-mean h and S distributions for the entire globe and the Atlantic Ocean are presented in Figs. 5 and 6 in comparison to the SYNC mean for years 9997–9999. For all other cases except AC, the 3-year means are used at their respective equilibration times. In synchronous distributions, the largest global differences in h occur in the abyssal Arctic Ocean. Here, S_AðzÞ and S_AC

Fig. 4. Global-mean root-mean-square (rms) potential temperature and salinity difference time series for the synchronous extensions of the accelerated integrations. The differences are from the years 9997–9999 mean of SYNC and are computed every 100 years using 3-year means. Year 0 values correspond to the rms differences for the last 3 years of the corresponding accelerated integrations.

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

333

Fig. 5. Time- and zonal-mean global and Atlantic potential temperature. The top two panels are for SYNC and the contour interval is 2 C. The other panels show difference distributions between SYNC and the accelerated equilibrium and synchronous equilibrium of the accelerated integrations. The color key is for the difference plots with a contour interval of 0.1 C. The Atlantic distributions exclude the Hudson Bay, Labrador Sea, and the Arctic Ocean. A 66-year mean (years 3001–3066) is used for AC.

have comparable magnitude differences; >0.4 C colder and >0.3 C warmer, respectively. The cold bias is larger in S_AðzÞ_450 than in S_AðzÞ_700, suggesting that the equilibration time based

334

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

Fig. 6. Time- and zonal-mean global and Atlantic salinity. The top two panels are for SYNC and the contour interval is 0.4 psu. The other panels show difference distributions between SYNC and the accelerated equilibrium and synchronous equilibrium of the accelerated integrations. The color key is for the difference plots with a contour interval of 0.02 psu. The Atlantic distributions exclude the Hudson Bay, Labrador Sea, and the Arctic Ocean. A 66-year mean (years 3001–3066) is used for AC.

on the global-mean time may not be very reliable, i.e. the accelerated integration needs to be continued. This is also evident in AðzÞ–SYNC plots which show smaller differences than those of

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

335

S_AðzÞ_450)SYNC. In contrast, as indicated earlier, there is little gain in integrating AC beyond 3000 years. In the North Atlantic, S_AðzÞ cases are more than 0.1 C colder at intermediate depths. The differences are somewhat higher locally for S_AC cases in the northern North Atlantic. The S_AC abyssal Atlantic distributions are uniformly warmer by more than 0.1 C, associated with the changes in the meridional overturning circulation (see Fig. 7).

Fig. 7. Time-mean zonally integrated global and Atlantic meridional overturning circulation obtained with the Eulerian-mean velocity. The top two panels are for SYNC and the contour interval is 4 Sv. The other panels show difference distributions between SYNC and the synchronous equilibrium of the accelerated integrations with a contour interval of 0.1 Sv. The Atlantic distributions include the Hudson Bay, Labrador Sea, and the Arctic Ocean.

336

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

In S (Fig. 6), the fresh bias of AðzÞ and S_AðzÞ cases is quite evident. Indeed, this bias is more than 0.1 psu in the Arctic and North Atlantic. In contrast, S_AC cases reveal much reduced basin scale differences. Here, the largest differences are confined to the deep Arctic (0.04 psu) and near the surface in the North Atlantic (0.12 psu). The AC accelerated equilibrium solutions (Figs. 5 and 6) show non-negligible differences mostly confined to the upper 1000 m in both h and S. This further demonstrates that an AC accelerated equilibrium is simply not good enough and that a synchronous extension is absolutely necessary. The global and Atlantic meridional overturning circulation distributions due to the Eulerianmean velocity are given in Fig. 7. In S_AðzÞ, the biggest differences (0.3 Sv) occur in the Southern Ocean. Elsewhere, the comparisons show minor (0.1 Sv), relatively small scale differences. In contrast, S_AC differences are a little larger both in magnitude and their spatial extents. The Northern Hemisphere overturning circulation, dominated by the thermohaline circulation primarily occurring in the North Atlantic, penetrates slightly deeper in S_AC than in SYNC, producing 0.6 Sv difference. In the abyssal ocean, the circulation associated with the Antarctic Bottom Water appears to have somewhat higher transports into the Equatorial Northern Hemisphere. Note that there is little change in the maximum and minimum transports. These differences in the circulations are too small to result in any noticeable, significant changes in the northward heat transports given in Fig. 8. The peak differences are below 0.02 PW.

Fig. 8. Time-mean global and Atlantic northward heat transport due to the Eulerian-mean velocity. The top panels are for SYNC. The bottom panels show differences between SYNC and the synchronous equilibrium of the accelerated integrations. The Atlantic distributions include the Hudson Bay, Labrador Sea, and the Arctic Ocean.

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

337

4.4. Seasonal cycle The monthly mean transport time series presented in Fig. 9 show that AðzÞ and AC cycles are non-negligibly different from that of SYNC both in phase and amplitude. The phase errors disappear almost instantaneously upon switching to synchronous integrations. This is also true for the amplitude for some transports, e.g. Indonesian Throughflow and Mozambique Channel. In

Fig. 9. Monthly mean transport time series for the Antarctic circumpolar current (ACC) at Drake Passage, Indonesian Throughflow, Florida Strait, and Mozambique Channel for intra- and inter-annual cycle comparisons. The SYNC equilibrium is given by the red line. The AðzÞ and AC cycles at about 700 and 3000 years are represented by the solid blue and green lines, respectively. The remaining blue and green lines are for S_AðzÞ_700 and S_AC_3, respectively. The time series are shown for 3-year periods, corresponding to forcing years 1991–1993.

338

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

contrast, the S_AðzÞ_700 and S_AC_3 ACC time series are still far from the SYNC one even after 18 years of integration. This slow adjustment process is likely due to deep water formation processes off the Antarctic coast and associated diffusive time scales in the abyssal ocean. Gent et al. (2001) show the effects of this deep water formation on the ACC transport. After 18 years, the Florida Strait transport amplitude mismatches are also evident for S_AC_3. In general, early parts of the AðzÞ synchronous extensions appear to be closer to the SYNC cycle. These transport time series suggest that longer synchronous integrations are needed to fully recover the SYNC cycle, as evidenced by the equilibrium time series, in contrast with Danabasoglu et al. (1996) where no synchronous equilibrium solutions were obtained. In the ocean interior, away from the convective regions, there are no appreciable phase or amplitude differences between SYNC and all other synchronous or accelerated near-surface h time series (not shown) due to strong feedbacks between the surface h and surface heat fluxes. In contrast, likely as a result of lack of any significant feedbacks between the surface S and freshwater fluxes, the near-surface S time series (not shown) reveal some amplitude differences in comparison to the SYNC cycle. The only exception to this is the S_AC_3 equilibrium cycle which matches the SYNC one quite remarkably. The AC cycle also displays some phase differences. Further down in the deep ocean, the seasonal cycle amplitudes are very small. Nevertheless, there are no major phase and amplitude differences for h and S cycles for all cases except for the offset of the mean distributions. This offset is of course an important concern especially for S in all AðzÞ and S_AðzÞ cases. Finally, the AðzÞ seasonal cycle amplitudes are significantly larger than in other cases in the abyssal ocean, in agreement with Danabasoglu et al. (1996).

5. Discussion and conclusions We have performed a 10 000-year synchronous integration using a comprehensive OGCM subject to realistic, time-dependent forcing in a global domain. To our knowledge, this is the first such integration of an OGCM. Two accelerated, but otherwise identical, equilibrium integrations along with their multi-century equilibrated synchronous extensions have also been obtained to quantitatively determine how well they reproduce the synchronous control solution. The two accelerated cases use tracer time steps increasing with depth and unequal momentum and tracer time steps with no depth variations, respectively (Bryan, 1984). We define the equilibration time as the time when the annual- and global-mean potential temperature trend is below 105 C year1 . The synchronous control integration achieves equilibrium in about 3500 years. The accelerated equilibration times are 685 and 3000 surface years for the two cases, respectively. Their synchronous extensions require about 100 and 700 surface years. The accelerated solutions do indeed differ from each other and from those of the control experiment, in agreement with Wang (2001). However, for practical purposes, many aspects of the accelerated and control solutions are more similar to each other than they are different. One important exception here is the detrimental effects of non-conservation due to the depth acceleration method, particularly obvious in salinity distributions. Depending on the surface forcing type and details, the level of non-conservation may vary, and further synchronous integrations may not necessarily eliminate these errors. Therefore, we do not recommend this kind of acceleration for any meaningful integration.

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

339

Obviously, the desired equilibrium solution can only be obtained via a synchronous integration method. However, if one wishes to expedite this extended process, we do recommend using unequal momentum and tracer time steps instead. Even in our present unfavorable set-up where further increases in the tracer time step are limited by some technical issues, we get a factor of 2.5 reduction in the computational cost (accelerated + synchronous extension) compared to the synchronous equilibration. We speculate that these savings could approach order 10 when our technical limitations are eliminated. Any accelerated integration must definitely be followed by a synchronous segment. Among the reasons are the adjustments of the intra- and inter-annual cycles and any possible oscillatory behavior present in the accelerated phase. In our experience, these oscillations disappear promptly upon switching to synchronous integration. With our present definition of equilibrium, the synchronous equilibration times reported here are much longer than those of Danabasoglu et al. (1996) where a synchronous control case was not obtained. Having this control integration at hand now, one can decide how long to integrate synchronously after acceleration considering the root-mean-square difference time series. If the incremental gains with further integration are judged to be minimal or of minor significance, shorter integrations may certainly suffice. We note that the equilibration times will undoubtedly depend on the model physics and parameter choices, particularly diffusion coefficients in the abyssal ocean giving the longest time scales. Nevertheless, because most non-eddy-resolving, coarse resolution climate models, including the ·1 version of the CCSM2, use similar physics and parameter choices as in our present study, we do not anticipate significant deviations from our current equilibration times. However, the permissible values of the acceleration factors are likely to be lower at higher resolutions (Killworth et al., 1984), thus somewhat limiting the degree of payoff. Usage of any acceleration methods is of course not permissible in eddy-permitting or -resolving, high resolution models, because the presence of unstable meso-scale eddies will make their solutions highly time-dependent. Finally, we find that the equilibration times are strongly affected by the choice of surface boundary conditions. For example, replacing the present under-ice boundary conditions with restoring type boundary conditions (with order 10-day time scales over 12 m) produces shorter equilibration times in two sensitivity experiments with synchronous and depth accelerated time steps, respectively. We speculate that when surface properties are clamped as in restoring especially in these deep water formation regions, occurrence of more regular convective events leads to reduced diffusive time scales (via smaller length scales) in the abyssal ocean. Acknowledgements The author would like to thank Drs. P.R. Gent, W.G. Large, and F.O. Bryan for helpful discussions and suggestions. The National Center for Atmospheric Research is sponsored by the National Science Foundation. References Bishop, J.K.B., Rossow, W.B., 1991. Spatial and temporal variability of global surface solar irradiance. J. Geophys. Res. 96, 16839–16858.

340

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

Bishop, J.K.B., Rossow, W.B., Dutton, E.G., 1997. Surface solar irradiance from the International Satellite Cloud Climatology Project 1983–1991. J. Geophys. Res. 102, 6883–6910. Bryan, K., 1969. A numerical method for the study of the circulation of the world ocean. J. Comput. Phys. 4, 347–376. Bryan, K., 1984. Accelerating the convergence to equilibrium of ocean-climate models. J. Phys. Oceanogr. 14, 666–673. Bryan, K., Lewis, L.J., 1979. A water mass model of the world ocean. J. Geophys. Res. 84, 2503–2517. Bryan, K., Manabe, S., Pacanowski, R.C., 1975. A global ocean-atmosphere climate model. Part II: The oceanic circulation. J. Phys. Oceanogr. 5, 30–46. Comiso, J., 1999 (updated 2002). Bootstrap sea ice concentrations for NIMBUS-7 SMMR and DMSP SSM/I. National Snow and Ice Data Center, Boulder, CO, USA. Digital media. Available from . Danabasoglu, G., McWilliams, J.C., Large, W.G., 1996. Approach to equilibrium in accelerated global oceanic models. J. Climate 9, 1092–1110. Doney, S.C., Hecht, M.W., 2002. Antarctic bottom water formation and deep-water chlorofluorocarbon distributions in a global ocean climate model. J. Phys. Oceanogr. 32, 1642–1666. Doney, S.C., Yeager, S., Danabasoglu, G., Large, W.G., McWilliams, J.C., 2003. Modeling global oceanic inter-annual variability (1958–1997): simulation design and model-data evaluation. NCAR Tech. Note NCAR/TN-452 + STR, 48 pp (available from NCAR, P.O. Box 3000, Boulder, CO 80307, USA). Dukowicz, J.K., Smith, R.D., 1994. Implicit free-surface formulation of the Bryan–Cox–Semtner ocean model. J. Geophys. Res 99, 7991–8014. Gent, P.R., McWilliams, J.C., 1990. Isopycnal mixing in ocean circulation models. J. Phys. Oceanogr. 20, 150–155. Gent, P.R., Large, W.G., Bryan, F.O., 2001. What sets the mean transport through Drake Passage? J. Geophys. Res. 106, 2693–2712. Griffies, S.M., 1998. The Gent–McWilliams skew flux. J. Phys. Oceanogr. 28, 831–841. Hahn, C.J., Warren, S.G., London, J., Jenne, R.L., Chervin, R.M., 1987. Climatological data for clouds over the globe from surface observations. Rep. NDP-026, 57 pp (available from Carbon Dioxide Information Center, Oak Ridge, TN 37831-6050). Huang, R.X., Pedlosky, J., 2002. On aliasing Rossby waves induced by asynchronous time stepping. Ocean Modell. 5, 65–76. Kalnay, E., Kanamitsu, M., Kistler, R., Collins, W., Deaven, D., Gandin, L., Iredell, M., Saha, S., White, G., Woollen, J., Zhu, Y., Chelliah, M., Ebisuzaki, W., Higgins, W., Janowiak, J., Mo, K.C., Ropelewski, C., Leetmaa, A., Reynolds, R., Jenne, R., 1996. The NCEP/NCAR reanalysis project. Bull. Amer. Meteor. Soc. 77, 437–471. Killworth, P.D., Smith, J.M., Gill, A.E., 1984. Speeding up ocean circulation models. Ocean Modell. 56, 1–4. Large, W.G., Yeager, S.G., 2003. Forcing global ocean models with diurnal to decadal variability. NCAR Tech. Note (in press). Large, W.G., McWilliams, J.C., Doney, S.C., 1994. Oceanic vertical mixing: a review and a model with a nonlocal boundary layer parameterization. Rev. Geophys. 32, 363–403. Large, W.G., Danabasoglu, G., Doney, S.C., McWilliams, J.C., 1997. Sensitivity to surface forcing and boundary layer mixing in the NCAR CSM ocean model: annual-mean climatology. J. Phys. Oceanogr. 27, 2418–2447. Large, W.G., Danabasoglu, G., McWilliams, J.C., Gent, P.R., Bryan, F.O., 2001. Equatorial circulation of a global ocean climate model with anisotropic horizontal viscosity. J. Phys. Oceanogr. 31, 518–536. Ledwell, J.R., Montgomery, E., Polzin, K., Laurent, L.St., Schmitt, R., Toole, J., 2000. Evidence for enhanced mixing over rough topography in the abyssal ocean. Nature 403, 179–182. Levitus, S., Boyer, T., Conkright, M., Johnson, D., OÕBrien, T., Antonov, J., Stephens, C., Gelfeld, R., 1998. World Ocean Database 1998, vol. I: Introduction. NOAA Atlas NESDIS 18, US Government Printing Office, Washington, DC, 346 pp. Rossow, W.B., Schiffer, R.A., 1991. ISCCP cloud data products. Bull. Am. Meteor. Soc. 72, 2–20. Smith, R.D., Gent, P.R., 2002. Reference manual for the Parallel Ocean Program (POP), ocean component of the Community Climate System Model (CCSM2.0). Los Alamos National Laboratory Technical Report LA-UR-022484 (Available from LANL, Los Alamos, New Mexico 87545. Also available from ). Smith, R.D., McWilliams, J.C., 2003. Anisotropic horizontal viscosity for ocean models. Ocean Modell. 5, 129–156.

G. Danabasoglu / Ocean Modelling 7 (2004) 323–341

341

Smith, R.D., Dukowicz, J.K., Malone, R.C., 1992. Parallel ocean general circulation modeling. Physica D 60, 38–61. Smith, R.D., Kortas, S., Meltz, B., 1995. Curvilinear coordinates for global ocean models. Los Alamos National Laboratory Technical Report LA-UR-95-1146 (available from LANL, Los Alamos, New Mexico 87545). Spencer, R.W., 1993. Global oceanic precipitation from the MSU during 1979–91 and comparisons to other climatologies. J. Climate 6, 1301–1326. Steele, M., Morley, R., Ermold, W., 2001. PHC: a global ocean hydrography with a high quality Arctic Ocean. J. Climate 14, 2079–2087. Available from . Wang, D., 2001. A note on using the accelerated convergence method in climate models. Tellus 53A, 27–34. Wood, R.A., 1998. Time step sensitivity and accelerated spinup of an ocean GCM with a complex mixing scheme. J. Atm. Oceanic Tech. 15, 482–495. Xie, P., Arkin, P.A., 1996. Analysis of global monthly precipitation using gauge observations, satellite estimates, and numerical model predictions. J. Climate 9, 840–858.