An operational implementation of the GHER model for the Black Sea , with SST and CTD data assimilation

An operational implementation of the GHER model for the Black Sea, with SST and CTD data assimilation L. Vandenbulcke, A. Capet, J. M. Beckers, M. Grégoire, and S. Besiktepe GeoHydrodynamics and Environment Research, Université de Liège, Belgium NATO Undersea Research Centre, La Spezia, Italy Received: 7 August 2009 – Accepted: 18 August 2009 – Published: 31 August 2009 Correspondence to: L. Vandenbulcke (luc.vandenbulcke@ulg.ac.be) Published by Copernicus Publications on behalf of the European Geosciences Union.


Introduction
The "Turkish Straits Systems 08" (TSS08) cruise took place in September 2008 and aimed at the process-oriented real-time characterization of a marine region, in this case the Dardanelles Strait, the Marmara Sea and the Western Black Sea.During the campaign, numerous in situ measurements were carried out onboard the research vessel (R/V) Alliance as well as some other small ships of the riparian countries.Processes of different scales were measured: large scale hydrographic surveys provided a synoptic view of the basin, but mesoscale and sub-mesoscale measurements were also carried out by means of towed yoyo-systems, CTD chains, gliders etc.Fixed moorings also provided measurements of currents, sea level, temperature and salinity.Furthermore, satellite observations were acquired and sent in near real-time to the R/V.Finally, different ocean circulation models were also run onboard the R/V.In particular, the HOPS (Robins et al., 1996;Lozano et al., 1996) and GHER (Beckers, 1991) models were used.
The present study focuses on the near real-time implementation of the latter model Introduction

Conclusions References
Tables Figures

Back Close
Full in the Black Sea with a 5 km horizontal resolution.We examine whether the model can be run fast enough for delivering timely forecasts of reasonable skill, and what requirements are needed.The model domain was pre-implemented on a bootable external USB hard drive, which was brought on the R/V on 19 September 2008; simulations were only carried out from that day on.First, some hindcasts starting on 16 August were run up to 15 September, and compared with measurements.Afterward, the model was run up to 19 September, and then daily in order to provide operational forecasts, up to the end of September.Section 2 describes the GHER hydrodynamic model and its data assimilation scheme.Section 3 examines the model results during the cruise.Conclusions are given in Sect. 4.

GHER model and data assimilation scheme
The GHER hydrodynamic model is available online (http://modb.oce.ulg.ac.be/viewsvn/) with a user-manual (http://modb.oce.ulg.ac.be/mediawiki/index.php/ Category:GHER3D).It is a free-surface primitive equation model developed in the early '90s.It forecasts the prognostic variables of temperature, salinity, see surface elevation, horizontal velocity and turbulent kinetic energy, using the hydrostatic, β-plane and Boussinesq approximations.Horizontally, it uses an Arakawa-C grid.In the vertical, it uses a double sigma coordinate, the limit between the two zones being at 170 m depth.Its integration scheme is conservative for tracers.Furthermore, the model uses mode splitting: for computational efficiency, the baroclinic timestep is much larger than the barotropic one.The vertical turbulence uses a k turbulent kinetic energy closure scheme described in Nihoul et al. (1989).Further information about the model can be found in Beckers (1991).The GHER model has already been implemented in various basins; for an implementation in the Black Sea with horizontal resolutions of 15 and It comprises a data assimilation procedure described below, a nesting procedure for rapid-relocatable nested grids (Barth et al., 2005), and a module to couple it with passive tracers e.g. for ecosystem models, see Grgoire and Beckers (2004).The latter two possibilities are not used in the present implementation.The model is accompanied by the Singular Evolutive Extended Kalman (SEEK) data assimilation (DA) filter, described in Pham et al. (1998); Brasseur et al. (1999).This DA scheme is a reduced-rank approximation of the Extended Kalman filter.In the latter filter, only the first-and second-order moments of the error statistics are retained.Thus, we suppose that the anomalies can be considered quasilinear, and the model error is approximately Gaussian.This is obviously not the case, but, for relatively short time forecasts, we still will use this widespread approach.The low-rank model error covariance matrix is built based upon empirical orthogonal functions (EOFs) representing the model error variability.The implementation allows to use full 3-D multivariate (temperature (T ), salinity (S) and sea surface elevation (η)) EOFs.Hence, whenever data is assimilated, this leads to corrections on the T , S, and η variables.These corrections are multiplied by a radial Gaussian function centered on the corresponding observation in order to limit the spatial extent of the correction that an observation can yield.In the present case, the Gaussian's extent is put to 100 km.Let's notice that if the model error covariance matrix would accurately represent the model error covariance, this step would not be necessary; unphysical long-range correlations would not be present in the computed statistics.Finally, in view of the large scales of interest, the geostrophic velocity correction (corresponding to the T , S and η corrections) is computed and applied to the model horizontal velocity variables U and V .Some specific aspects of the current implementation are given hereafter.Although the region of interest if the South-Western Black Sea, our model domain covers the whole basin in order to avoid (unknown) open sea boundary conditions.The model bathymetry is the Smith and Sandwell (1997) bathymetry with some smoothing.The horizontal resolution is about 5 km; this is smaller than the first internal radius of deformation which is 20 km in the Black Sea (Ozsoy and Unluata, 1997) 31 vertical levels.The barotropic and baroclinic timesteps are 15 s and 15 min respectively.With these settings, the model (using OpenMP parallelization) takes about 5 min to simulate one day on a modern computer with an Intel quad-core cpu.Let's note that another implementation with a horizontal resolution of about 2 km, and timesteps of 2.5 s and minutes respectively, requires about 2 h and 30 min for one day of simulation; it was not used operationally during this experiment.The sea-atmosphere fluxes are computed interactively with bulk formula; the atmospheric forcing fields forecasts (for air temperature, air dew temperature, wind speed and atmospheric pressure) are obtained in near real-time from the COAMPS model.The cloud coverage field was unavailable and obtained from the NCEP forecasts website.Finally, the Bosphorus and 6 rivers (Danube, Dnestr, Dniepr, Kisil, Rioni, and Sakarya) are represented in the model with imposed fluxes obtained from hydrological simulations for the year 2000 (Ludwig et al., 2009).It was decided to start the model on 16 August, 2 weeks before the actual TSS08 start, and a few days after the start of COAMPS forcing fields availability.Initial conditions were taken from the MODB/MEDAR4 monthly climatology (Brasseur et al., 1996).Data from the HYDROB-LACK91 campaign (Aubrey et al., 1992) could also have served to build an initial condition, but this was not more recent than the climatology and did not correspond particularly to the condition of the month August.

Hindcast with data assimilation
Output fields from the simulation staring on 16 August and ending on 15 September were compared with sea surface temperature (SST) images obtained during the cruise from the MODIS satellite (http://modis.gsfc.nasa.gov); it could immediately be seen that the forecast was not adequate (the root mean square error between forecast and the first SST image was 2.32 • C, see to an inappropriate initial condition, obtained from a climatology over 15 years old.Thus, it was decided to restart this simulation while assimilating both SST images and temperature and salinity profiles obtained from CTD measurements onboard the R/V.It was hoped that the joint assimilation of synoptic surface data and vertical profiles would consistently correct the 3-D field.However, the data assimilation scheme requires an estimation of the model error covariance matrix, which was not available.Therefore, it was decided to use the time variability of the model during the previous hindcast as a proxy for the error covariance.Fifteen multi-variate EOFs were thus obtained; each variable (T, S, η) was normalized according to the variable's variance and the corresponding grid cell size.We noted from the start that this relatively small rank (15) would probably not allow mesoscale corrections.Another approximation in the data assimilation procedure consisted of not updating the model error space but rather keeping these 15 EOFs constant in time.The observations error covariance matrix was chosen diagonal, with standard deviations of 50 • C for the MODIS SST data (to account for the representativity error due to the relatively high spatial resolution of the SST observation), and 1 • C and 0.7 psu for the CTD measurements of temperature and salinity respectively.8 SST maps were selected between 24 August and 9 September (each of them around 11:00 in the morning); the other available images were affected by many clouds.There were 182 CTD profiles between between 26 August and 11 September, but some of them were in the Straits: CTD profiles in the Black Sea were only available between 2 and 8 September.However, the data assimilation scheme allows out-of-grid observations to be included, and rejects them automatically.Thus, no preprocessing was applied to the CTD data, and all of them were simply listed in chronological order in the model's data assimilation parameters file.Finally, let's note that only the western half of the Black Sea was covered by the CTD measurements and the (available part of the) SST images; hence no direct corrections were applied to the eastern part.
Table 1 shows the root mean square (rms) difference between SST observations and the model forecast (i.e.independant data, just before assimilation) and model analysis Introduction

Conclusions References
Tables Figures

Back Close
Full (just after assimilation).It can be seen that the mean difference before the first assimilation cycle, on 24 August, is large (2.32 • C), but strongly reduced with DA.Then, in between assimilation cycles, the model slightly increases the rms error.In particular, after the third SST analysis, CTD profiles are also assimilated.Due to the strongly approximated model error covariance matrix, this might create instabilities or at least inconsistent updates.Hence, between the fourth and the seventh SST assimilation cycles, the model (and CTD assimilation cycles) relatively strongly increases the surface temperature rms error; but of course SST assimilation cycles in turn reduce the rms error.One can suppose that a part of the rms error remaining after the latest DA cycles (1.3 to about 1.1 • C) is due to this back-and-forth game of corrections due to SST and CTD observations.The remaining part is orthogonal to the (low rank) model error space, and could never be corrected with our setup.
Figures 1 to 4 show the effect of data assimilation on model surface temperature.As announced before, assimilation of the first SST image (covering from about 30 • E to 35 • E) strongly heats the whole western basin, by 1.5 • C up to 4 • C (Fig. 2).However, the realistic cold water mass along the western part (around 42 • N, 33 • E) of the southern coast is also heated away (Fig. 1).This cold water will be re-created by the model and by later assimilation cycles, particularly the SST assimilated on 6 September (see Fig. 3).Of course, the eastern part of the basin (and in the case of the first assimilation cycle, the western coast not covered by data) are not corrected.From the second assimilation cycle on, the effect of assimilation is still mostly heating, but not everywhere anymore: some parts of the domain are cooled.The heating trend during assimilation cycles may be caused by too cool water still flowing in the western half basin from the uncorrected eastern half basin.However, corrections tend to become smaller and smaller.Figure 4 shows the effect of assimilation of a temperature and salinity profile; the 100th CTD profile is chosen as an example.The effect is mainly to put back in place the cold water mass along the Turkish coast.

Operational forecasts
The hindcast simulation described before was continued daily, in forecasting mode, starting on 20 September.Operational forecasts were run as soon as the COAMPS atmospheric forcing fields were available.The model currently requires atmospheric forcing fields in its own binary format; the conversion takes just a couple of seconds with the proper scripts.Afterward, the simulation takes about 5 min for one day of forecast.Outputs are written directly in netCDF format (K. et al., 2006),which were then stored on the cruise's database.Generally speaking, the forecasts showed good agreement with the data (MODIS SST mages and in situ measurements).Comparison with HOPS, the other model run onboard the R/V, was difficult due to differences between the two setups.In particular, the HOPS model was implemented on a domain covering only the south-western part of the Black Sea, hence using open sea boundary conditions.Nevertheless, Figs. 5 and 6 show the SST forecasts of the GHER and HOPS models, and the closest (in time) MODIS image.With respect to this SST image, the GHER model forecast is good, except the cold water area along the Turkish coast which is really too cold (too strong correction during the latest assimiation cycle).The HOPS model forecast looks generally too cold by 2 or 3 degrees, but the cold area along the Turbkish coast is not over cooled like in the GHER model; it is rather too warm in the HOPS model.forecasts were bad (compared to SST images).Hence, the hindcast was started over, and this time 8 SST images and about 100 CTD temperature and salinity profiles were assimilated.The effect of data assimilation was very large and resulted in a strong heating of the model temperature during the first assimilation cycles; but the following cycles lead to more modest and local effects.Following the hindcast simulations, daily operational forecasts were carried out up to the end of September, with results in good agreement with observations.On a modern personal computer (with a cpu with 4 cores), one day of simulation required about 5 min, plus a few minuts to download and convert the atmospheric forcing fields into the desired format.Outputs were saved directly to netCDF format by the model, and stored on the campaign database.T the GHER model on 8 September, midday , we have shown an operational implementation of the GHER model in the Black Sea during the TSS08 campaign.The model codewas pre-configured and stored on a bootable external USB hard drive.The drive was brought onboard the research vessel Alliance during the campaign, and connected on a computer with an intel quadcore cpu.First, a month of hindcast simulations was carried out.Mostly because of relatively bad initial conditions coming from a 15-years old climatology, the model Introduction

Fig. 1 .
Fig. 1.Surface temperature [ • C] and elevation [m] of the forecast on 24 August at 10h15, before (top panel) and after (bottom panel) assimilation of SST data in the western half of the basin.figure8

Fig. 1 .
Fig. 1.Surface temperature [ • C] and elevation [m] of the forecast on 24 August at 10:15, before (top panel) and after (bottom panel) assimilation of SST data in the western half of the basin.
. Our model uses Introduction

Table 1 )
. The inadequacy was attributed mainly Introduction

Table 1 .
Root mean square difference between SST observations and model forecast and analysis [ • C].