Using the SKAGEX dataset for evaluation of ocean model skills

Journal of Marine Systems 18 Ž1999. 313–331

Using the SKAGEX dataset for evaluation of ocean model skills Jarle Berntsen a

a,b,)

, Einar Svendsen

b,1

Department of Mathematics, UniÕersity of Bergen, Johs. Bruns gt. 12, N-5008 Bergen, Norway b Institute of Marine Research, Nordnesparken 2, N-5024 Bergen, Norway Received 16 December 1996; accepted 23 July 1997

Abstract Numerical ocean models are now being applied in numerous oceanographic studies. However, the qualities of the model results are often uncertain and there is a great need for standards and procedures for evaluation of the skills of numerical general circulation models. In this paper measurements from repeated hydrographical sections across Skagerrak taken in 1990, the SKAGEX dataset, are used to evaluate the skills of two s-coordinate ocean models and to study the sensitivity of these models to model parameters. A methodology for quantification of model skills based on observations from repeated hydrographical sections in general is suggested. Area averages of absolute differences are for Skagerrak completely dominated by the discrepancies in the upper few meters of the ocean and may not be used to assess models’ abilities to reproduce the fields in the larger and deeper part of the ocean. Therefore, discrepancies between average values in time from the observed fields and time averaged values from model outputs are related to the natural variability of the fields. The numbers produced with the suggested measure are relative numbers that will be specific for each section and for each series of observation. Ideally we would therefore like to see the measures computed for a number of sections for various models and choices of model parameters in order to assess model skills. The value of the SKAGEX dataset as a tool for model improvements is demonstrated. Evidence to support the importance of applying non-oscillatory, gradient preserving advection schemes in areas with sharp density fronts is given. The method is used to identify that the forcingrinitial valuesrboundary values for the temperature field are inferior to the corresponding values for the salinity field. With the present coarse resolution, 11 layers in the vertical, it is shown that it is far from obvious that the quality of the model results improve when replacing simple Richardson number formulations for vertical mixing processes with higher order turbulence closure in the Skagerrak area. q 1999 Elsevier Science B.V. All rights reserved. Keywords: general circulation models; evaluation; Skagerrak; sections; hydrodynamics

1. Introduction A number of numerical ocean circulation models have been under development over the last decades and applied in numerous oceanographical studies. ) Corresponding author. Fax: q47-55-584885; e-mail: [email protected]. 1 Fax: q47-55-238584; e-mail: [email protected].

Hydrodynamical models are coupled to chemical– biological models with different degrees of complexity. The methodology for measuring the quality of the vast amount of model results being produced is, however, poorly developed. One exception appears to be root mean square errors for the evaluation of outputs from tidal models. In baroclinic studies for instance rms errors may not be appropriate as quality measures. Small shifts of eddies, meanders or fila-

0924-7963r99r$ - see front matter q 1999 Elsevier Science B.V. All rights reserved. PII: S 0 9 2 4 - 7 9 6 3 Ž 9 7 . 0 0 1 1 1 - 5

314

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

ments may produce large rms errors even if the dynamics are well reproduced by the model. The use of rms errors to compare model outputs and observations in baroclinic studies would require very high sampling rates and also high quality initial and boundary values for the numerical model. Both requirements above are hard to meet in practice. Scientists working in this field of research are often comparing a sequence of plots andror numbers with corresponding sequences produced with different models or the same model with different model inputs or parameters. Then the question is: ‘Which sequence of model outputs is in best agreement with the observations?’ The reasoning behind the choices made are today often qualitative because it is difficult to measure the improvements when going from one set of model outputs to another. More generally it will therefore be difficult to identify improvements in the qualities of complex ocean models, and unless progress is made on systematic methods for comparing model results and observations it is easy to foresee a growing number of ocean models appearing on the market, all with some skills and they may all claim to be better than the competitors in some respect. In their closing remarks to the workshop on Quantitative Skill Assessment for Coastal Ocean Models, Lynch et al. Ž1995. point at the need for high-quality data sets for field based intercomparison exercises and for appropriate methods for comparing model results and data. In the present paper the connections between applications, observations and model evaluation are discussed in Section 2. In Section 3 an evaluation technique based on using data from hydrographical sections in general is described. In Section 4 measurements from the SKAGEX experiment Žsee Dybern et al. Ž1994.. are used to evaluate the skills of two sigma-coordinate ocean models.

2. Applications, observations and evaluation Quality measures for model results should be relevant for major application areas for numerical ocean models. In the management of the ocean environment and ocean resources it will become increasingly more important to quantify the transports of water masses,

pollutants, nutrients and fish eggs and larvae between different areas of the ocean Žsee e.g., Berntsen et al. Ž1994. and Skogen et al. Ž1995... Models describing the physical variables are coupled to chemical–biological models to describe the growth of algae and the effects of anthropogenic inputs of nutrients to the ocean. For such applications vertical exchange processes and transports in the upper layers of the ocean are often of major importance. The basic hydrodynamic equations for the ocean models are the same for different ocean areas, but the major geophysical processes may differ. Therefore, the ‘best’ ocean model for one area of the ocean may or may not be the appropriate choice for another area. A project aimed at revealing the properties of models intended for applications in a specific area like the Skagerrak Žsee Fig. 1. should therefore be based on observations from the chosen model area, or an area with similar characteristics. With the present resolution of ocean models there will typically be a number of potentially important unresolved processes. For many applications, however, this may not be of major concern as long as the mean fields and transports are realistically modelled and the skills of models to reproduce average fields should therefore be assessed. To produce good estimates of average fields from observations we need enough measurements to resolve the major time and space scales. The fastest relevant waves are the long surface gravity waves that propagate with a speed of c 0 s gH where g is the gravitation constant and H the bottom depth. For H s 100 m, c 0 becomes approximately 30 m sy1 . The velocities are to a large extent dominated by these long gravity waves, and the variations in the velocity fields may therefore be large even on small time and space scales. This is confirmed by SKAGEX currentmeter measurements which sometimes showed changes of the order of 100 cm sy1 in 1 day at the southern coast of Norway. From shipborn ADCP measurements horizontal gradients between the in and out flowing Skagerrak water horizontal current gradients were of the order of 1–2 = 10y5 sy1 Ž10–20 cm sy1 per 10 km. while the model ŽPOM. gave approximately 1.5 = 10y5 sy1 ŽSvendsen et al., 1996.. Estimates of residual volume fluxes through specific sections may be computed from moored current

'

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

315

Fig. 1. Bottom topography of the Skagerrak and the adjacent areas. A, B, C, D, E, F, G and H show the different sections with the positions of the hydrographical stations. Areas deeper than 500 m are hatched, and the 50 m and 200 m bottom contours are enhanced ŽDanielssen et al., 1996..

meters where the observations normally are taken frequently enough in time. However, the spatial scales are often not resolved. A methodology for computing residual transport estimates and corresponding reliable error estimates from current meter rig measurements has not to our knowledge been established. In 1990, from May 21 to May 27, seven current meter rigs covered a 100-km long section from Denmark to Norway and residual transport estimates were produced Žsee Svendsen et al. Ž1996... Such a coverage with current meter rigs is rare, but even for this case conservative estimates of the errors may be as large as the residual transports, and therefore it will be difficult to quantify the accuracy of residual

model transports based on integrated residual transports from current meter observations. Scalar fields like temperature, salinity, density, nitrate etc. are transported with the speed of the water masses. In and out of Skagerrak, typical velocities are 0.1–0.3 m sy1 . The velocities or the momentum are to a great extent transported with the speed of the surface gravity waves, which in Skagerrak may be up to 80 m sy1 . The time variability of the scalar fields is therefore smaller than the time variability of the velocity fields. This means that far less observations in time are needed to resolve the relevant time scales for scalar fields than for the velocity fields. We regard the system to be advection dominated and for such systems smaller time vari-

316

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

ability is reflected in smaller space variability. Thus we expect that fewer observations in space are needed to resolve the scalar fields than to resolve the velocity fields. Integrated measures from observations are to be used as approximations for the exact integrated measures for the real fields. Since the exact integrated measures are unknown, integrated measures from model results will be compared to integrated measures from observations. If the errors in the integrated measures from observations are large, we may find good agreement in our comparisons also in cases where the model outputs are far from the real fields. The value of integrated measures from observations as tools for model evaluation thus depend on small errors in the measures. Since scalar fields are more conservative, the errors in integrated measures from observations will be smaller than the errors in integrated measures from observed velocities. With smaller error bars on the measures from observations, it will be much easier to decide which of a sequence of model results is the best and to decide whether a change in a numerical algorithm, a model parameter or in the model forcing is an improvement or not. A good agreement across a section between model results and observed values of a scalar field like salinity does not directly mean that the transports, which may be our main concern, are correct. In ocean areas with relatively strong horizontal and vertical gradients in the scalar fields and weak local forcing of the fields it may, however, indicate that the transports in average have been well represented by the model. On the other hand for cases where the local forcing of the scalar fields is dominant, a good agreement between the model scalar fields and the observed scalar fields will not be an indication that the real residual currents are well represented. Scalar fields like temperature, salinity and nutrients show a yearly variation. In studies where the goal is to produce climatological mean fields that are representative for different periods of the year, one good evaluation strategy will be to apply data from hydrographical stations taken at the same time of the year for a sequence of years. In other studies we may try to reproduce the mean fields of different periods of specific years. An alternative would then be to focus on data taken frequently along fixed sections over a limited period of the specific year. Using data

from Skagerrak, we will follow the last approach in Section 4.

3. Evaluation methodology In ocean research and management some hydrographical sections are taken routinely. The spatial resolution is typically good enough to resolve the major length scales whereas the resolution in time may be a problem. The quality of an estimate for the mean value in time of a field variable measured at the same spatial position may be related to the natural variability of the variable over the period in time considered. We therefore suggest to relate discrepancies between time averaged values from measurements and time averaged values from model results to the natural variability of the fields. As a measure of the natural variability of a field over the time period considered we have applied the standard deviation of the field. Primarily we would like to compare time averages produced by model outputs with exact time averages for the real fields. From a limited number of observations in time we can only produce estimates of the exact time averages and the standard deviations may be used as measures on the uncertainties in these estimates. Time averages from model outputs may be very close to exact time averages even if the differences between the model output time averages and the estimated averages from observations are large in cases with large standard deviations. We want the quality measure to reflect this. A second reason for relating the discrepancies to the standard deviations is that we expect it to be more difficult to model average fields correctly for systems with large natural variation. Therefore, we expect a model with a given skill to produce model output time averages that are further from the real systems time averages when modelling systems with large fluctuations. We have tried to construct a model evaluation methodology that ideally would give the same quantitative skill measure whether we are modelling a system with small or large fluctuations. The numbers produced with the suggested measure are relative numbers that will be specific for each section and for each series of observation. Attempts have been made to apply measures based on absolute differences only. However, area averages

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

of such measures are at least for Skagerrak completely dominated by the discrepancies in the upper few meters of the ocean and models abilities to reproduce the fields in the major and deeper parts of Skagerrak were not assessed at all.

3.1. Measures from obserÕations Let a scalar field in position x i measured at time t j have the value Fi j. The average value in time of this field in position x i over NT observations will 1 NT j then be Fi s Ý F and the standard deviation: NT js1 i 1 2 NT j FiSD s 'NT y 1 Ý js1 Ž Fi y Fi . Let x be the horizontal coordinate along the hydrographical section and z the vertical coordinate. Let F Ž x, z . be the time average of F in the point Ž x, z .. An area average of F over the section will then be:

(

F˜ s

1 Area

HFd xd z

where the integral is taken over the area of the section. The time averages Fi and FiSD may be interpolated to a grid that is equidistant in x and equidistant in z. One way of approximating the above integral will then be: F˜ ;

1 NG

NG

Ý

Fi g

igs1

where NG is the number of grid points and Fi g are the interpolated values of Fi . The errors in the interpolation step and in the integration step may be estimated by taking the differences between approximations obtained with two different methods. Large discrepancies between model fields and observed fields may be due to constant offsets of the model fields. In such cases the gradients of the model fields may be correct even if the discrepancies are large. The gradients of the density field along a section are of particular interest because the internal

317

pressure is a major driving force. From values of the density we may use the thermal wind relation:

Ey yf

Ez

g Er s

r0 E x

to estimate the velocities, y , normal to the sections. f is the Coriolis parameter, g gravity, r 0 the reference density, r the density computed from salinity and temperature using the equation of state. To be able to integrate the above equation, we have to assume that we know the velocity at some vertical level, for instance y s 0 at the bottom. From observations of salinity and temperature we may thus compute for each observation time a geostrophic transport GT j. A time average GT with corresponding standard deviation GT SD that is an indicator for the variability in the density field along the section, may also be obtained. 3.2. Comparing model results and data From our model results we may pick from model points along the section model approximations to scalar fields at selected observation times at the same times as the physical observations. Let FMJ I be a model approximation to the scalar field F picked at model point x I along the section at time t J . From these values we may for each model point compute FM I and FMSDI which should be compared to the corresponding estimates from observations. These time averaged measures may be interpolated to the same equidistant grid as Fi and FiSD . The discrepancy between the time averaged model and time averaged observed scalar field in grid point x i g may then be measured in numbers of standard deviations: SD DFi g s < FM i g y Fi g
318

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

tiesrinabilities to reproduce average fields may therefore be assessed.; From the equidistant grid values area averages DF, F˜M and F˜MSD may be com˜ The above puted with the same technique as F. discrepancy measure was also applied in Berntsen et al. Ž1996.. Notice that since the nominator and the denominator in the expression for DFi g above have the same dimensions, the discrepancy measures will be dimensionless numbers. The numbers produced with the suggested measure are relative numbers that will be specific for each section and for each series of observation. We may thus not conclude that a model producing model outputs giving quality measure 0.1 for one section is better than another model that produced model outputs giving quality measure 1.0 for another section. If the numbers 0.1 and 1.0 were for the same section, we would favour the first model. However, we have not been able to construct a method that helps us to decide how much smaller than for instance 1.0 the quality measure must be before we may state that we found a significant improvement. Such a method must be based on assumptions on the distributions of measurements and the correlations between them that we found difficult to make. Ideally we would like to see ; measures like DF computed for a number of sections for various models and choices of model parameters in order to assess model skills. If the measures are consistent, that is if all measures for Model 1 are smaller than the corresponding measures for Model 2 for all fields and all sections, we are willing to accept Model 1 as the better model of the two. Thus, in our search for better models or model set ups, it will be the relative sizes of the measures from one experiment to another that are of interest more than the absolute values. Model geostrophic transports may be computed at times t J and time averages GTM should be compared to GT and related to GT SD .

4. Data from SKAGEX and model results 4.1. The SKAGEX dataset A description of the Skagerrak Experiment is given in Dybern et al. Ž1994.. In Fig. 1 the SKAGEX area with the repeated sections is shown.

Here we will focus on data from section H, the outer Skagerrak, and section F, the central Skagerrak, gathered during SKAGEX I in May–June, 1990. In Svendsen et al. Ž1996., it is shown that the general circulation and major processes may more qualitatively be reproduced by the sigma-coordinate ocean model due to Blumberg and Mellor Ž1987.. In the period from 24 May to 20 June 1990, CTD measurements along section F were taken 10 times at 13 positions and totally, 5846 measurements of salinity and temperature were made. Section F is 114.7 km long, the maximum depth is 654 m and the total area of the section is approximately 33.5 = 10 6 m2 . Time and space integrated measures for the whole section F are given in Table 1 and corresponding figures for the upper 50 m are given in Table 2. The time average salinity field is shown in Fig. 2Ž1.. The standard deviations for this field are shown in Fig. 2Ž2.. The salinity field for the upper 50 m of section F is shown in Fig. 3Ž1.. During SKAGEX I 1990 CTD measurements along section H were taken nine times at approximately 15 positions and totally 4530 measurements of salinity and temperature were made. Section H is 177.3 km long, the maximum depth is 441 m and the total area of the section is approximately 24.3 = 10 6 m2 . Time and space integrated measures for the whole section H are given in Table 3 and corresponding measures for the upper 50 m are given in Table 1 Area integrated measures for the whole section F-comparisons datarmodel results-BOM Ž1 Sv s1 Sverdrups1=10 6 m3 sy1 .

S˜ S

& SD

Data

RUNB1

RUNB2

34.56 p.s.u.

34.72 p.s.u.

34.45 p.s.u.

0.15 p.s.u.

0.13 p.s.u.

0.24 p.s.u.

7.688C

8.128C

8.208C

T SD GTout GTin SD GTout GTinSD

0.338C 0.721 Sv 0.591 Sv 0.167 Sv 0.148 Sv

0.198C 0.921 Sv 0.750 Sv 0.164 Sv 0.229 Sv

0.188C 1.089 Sv 0.825 Sv 0.167 Sv 0.180 Sv

DS

&

)

4.36

2.78

DT

&

)

10.53

11.97

DG out

&

)

1.20

2.20

DG in Model transports

) )

1.07 1.122 Sv

1.58 0.879 Sv

T˜

&

&

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331 Table 2 Area integrated measures for the upper 50 m of section F-comparisons datarmodel results-BOM Ž1 Sv s1 Sverdrups1=10 6 m3 sy1 . Data S˜

RUNB1

RUNB2

33.00 p.s.u.

33.43 p.s.u.

33.14 p.s.u.

S SD

0.75 p.s.u.

0.48 p.s.u.

0.68 p.s.u.

T˜

&

9.308C

8.308C

8.068C

T SD GTout GTin SD GTout GTinSD

0.988C 0.107 Sv 0.108 Sv 0.038 Sv 0.022 Sv

0.458C 0.127 Sv 0.106 Sv 0.031 Sv 0.036 Sv

0.358C 0.103 Sv 0.062 Sv 0.032 Sv 0.025 Sv

DS

&

)

4.00

1.80

DT

&

)

2.03

2.64

DG out

&

)

0.53

0.11

DG in

)

0.09

2.09

&

&

Table 4. The time average salinity field for section H is shown in Fig. 4Ž1. and the corresponding standard deviations in Fig. 4Ž2.. The salinity field for the upper 50 m of section H is shown in Fig. 5Ž1.. To produce the space integrated measures the time integrated values are interpolated to a 50 = 50 grid. An interpolated value at a grid point is obtained from values of the field at approximately the same depth by linear interpolation. As long as we interpolate from values at the same depth, the uncertainties in the interpolation step, and also in the integration step, become typically an order of magnitude smaller than the standard deviations. We separate between inflowing and outflowing water masses and integrate all velocities inrout of Skagerrak separately to compute geostrophic transports inrout of the sections.

319

available from Princeton University and will here be denoted as the POM ŽPrinceton Ocean Model.. A more recent model developed by Berntsen et al. Ž1996b. at the Institute of Marine Research, Bergen and the University of Bergen, Norway, was later on applied in this study. The model will here be denoted as the BOM ŽBergen Ocean Model.. The prognostic variables of these models are the three components of the velocity field, temperature, salinity and surface elevation. The governing equations are the momentum equations, the continuity equation and conservation equations for temperature and salinity. The governing equations together with their boundary conditions are approximated by finite difference techniques. In the vertical a s-coordinate representation is used in both models. In this representation the sea surface is mapped to 0 and the sea bottom to y1, thus the depth at each value of sigma is proportional to the bottom depth. The models are implemented for an extended North Sea with a 20-km horizontal resolution ŽFig. 6Ž1.. and for SkagerrakrKattegat with a 4-km horizontal resolution ŽFig. 6Ž2... Vertically 11 s-coordinate layers are used for both models. The layers follow the bottom topography and are chosen to give high resolution near the surface. At 100-m depth the layers are 0.5 m, 0.7 m, 1.3 m, 2.5 m, 5 m, 10 m, 20 m, 20 m, 20 m, 15 m and 5 m thick. The same initial values and boundary conditions are used for both models. To represent horizontal sub-grid scale processes, the models utilizes the Smagorinsky Ž1963. diffusion formulation in which the horizontal viscosity, A M , and diffusivity coefficients, A H , are modelled by:

Ž A M , A H . s Ž C M ,C H . D y Ž Ž E urE x .

2

2

4.2. The ocean models

q Ž E urE y q EyrE x . r2

The circulation of the North Sea and Skagerrakr Kattegat was first approximated by a three-dimensional, primitive equation, time-dependent circulation model due to Blumberg and Mellor Ž1987.. This model performed favorably in the model evaluation project MOMOP Žsee Røed et al. Ž1995. and Hackett et al. Ž1995.. and is now being used by a number of Norwegian scientists. The version of the model applied in this study is based on the version publicly

q Ž EyrE y .

2 1r2

.

where D x and D y are the horizontal resolution in x and y, respectively. C M and C H are dimensionless viscosity and diffusivity parameters. 4.3. Differences between the models The governing equations and the numerical grid are common to both methods. The differences are

320

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

found in the choices of numerics and the major differences are summarized below. The POM applies time splitting whereas the BOM applies implicit methods for the surface gravity waves and therefore no time splitting is necessary. The internal pressure

terms in the BOM model are treated in z-coordinates to avoid artificial velocities in areas with strong stratification and varying topography. In the POM the leapfrog scheme is applied for advection whereas the BOM applies a gradient preserving and mono-

Fig. 2. Ž1. Mean values of observed salinity in p.s.u. for section F. Ž2. Standard deviations of observed salinity in p.s.u. for section F. Ž3. Mean values of model salinity in p.s.u. for section F from RUNP1. Ž4. Mean values of model salinity in p.s.u. for section F from RUNP2. Ž5. Mean values of model salinity in p.s.u. for section F from RUNB1. Ž6. Mean values of model salinity in p.s.u. for section F from RUNB2.

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

321

Fig. 2 Žcontinued..

tonic TVD scheme due to Roe and Sweby Žsee for instance Sweby, 1984.. 4.4. Initial Õalues The North Sea models are run from 15r10-89 to 1r8-90 and climatological values of velocity, temperature, salinity and water elevation for October are used as initial values ŽMartinsen et al., 1992.. Output of these fields from the 20-km model at 15r3-90 are interpolated to a 4-km grid and used as initial values for the SkagerrakrKattegat models which are run to 1r8-90. 4.5. Boundary conditions At the lateral open boundaries, except at the boundary to the Baltic, a flow relaxation scheme ŽFRS. is implemented ŽMartinsen and Engedahl, 1987.. The FRS-zones for both models are seven grid-cells wide. For the North Sea models climatological values of velocity, temperature, salinity and water elevation for the respective months, in addition to one tidal constituent Ž M2 ., are used to specify the lateral boundary conditions, and output from this model of the same state variables is every hour interpolated to the FRS-zones of the 4 km models

and used as lateral boundary conditions for these models. The flow to and from the Baltic is implemented after an algorithm due to Stigebrandt Ž1980.. The flow is determined from the difference in modelled water level between the southern Kattegat and the Baltic, taking climatological freshwater input to the Baltic into account. The water entering Kattegat from the Baltic is given a salinity of 8.0 p.s.u. In the 20 km model all inflowroutflow is placed at Storebelt. In the 4-km model, the flow is shared between Storebelt and Øresund Žsee Fig. 6Ž1 and 2... The models are run with hindcast atmospheric forcing Žmomentum flux and surface pressure every 6 h. provided by the Norwegian Meteorological Institute ŽReistad and Iden, 1995.. In lack of information on surface heat fluxes, we relax the sea surface temperature towards climatology ŽCox and Bryan, 1984.. The surface flux of heat is specified by: KH

ET Ez

sg ŽT )yT .

where K H is the vertical diffusivity in the top layer, T ) the climatological value of temperature and g a time constant selected to be 1.735 = 10y5 m sy1 .

322

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

Fig. 3. Ž1. Mean values of observed salinity in p.s.u. for section F-upper 50 m. Ž2. Mean values of model salinity in p.s.u. for section F from RUNB1-upper 50 m. Ž3. Mean values of model salinity in p.s.u. for section F from RUNB2-upper 50 m.

This means that during weak forcing the temperature in the upper 15 m or so of the surface layer return to the climatological value on a time scale of 10 days. The North Sea models are run with monthly mean river runoff from the Rhine, Meuse, Scheldt, Ems, Weser, Elbe, Humber, Tyne and Tees. Daily river

runoff from the six largest Swedish rivers between Øresund and Norway is supplied. Fresh water runoff from the coast of Norway, including Glomma, is based on monthly mean fluxes and distributed on nine outlets along the Norwegian coast. The SkagerrakrKattegat models are run with the same fresh

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331 Table 3 Area integrated measures for the whole section H-comparisons datarmodel results-BOM Ž1 Sv s1 Sverdrups1=10 6 m3 sy1 . Data

RUNB1

RUNB2

34.55 p.s.u.

34.89 p.s.u.

34.81 p.s.u.

0.40 p.s.u.

0.14 p.s.u.

0.30 p.s.u.

8.208C

8.628C

8.578C

T SD GTout GTin SD GTout GTinSD

0.448C 0.679 Sv 0.177 Sv 0.208 Sv 0.095 Sv

0.268C 0.540 Sv 0.126 Sv 0.190 Sv 0.059 Sv

0.288C 0.479 Sv 0.180 Sv 0.116 Sv 0.028 Sv

DS

&

)

0.75

0.89

DT

&

)

10.76

10.69

DG out

&

)

0.67

0.96

DG in Model transports

) )

0.54 0.744 Sv

0.03 0.716 Sv

S˜ S T˜

& SD

&

&

water discharge from the Swedish and Norwegian coasts, but because of the finer horizontal resolution the discharge from the Norwegian coast is distributed on 11 coastal cells. 4.6. Skagerrak as an experimental mixing tank Skagerrak may be regarded as a large experimental tank where a number of water masses entering the area are mixed before leaving as the Norwegian coastal current. Skagerrak is therefore well suited as a test area for sub-grid scale parametrizations. From the plots of average values of salinity Žsee Fig. 2Ž1. and Fig. 4Ž1.., we may study how deep the surface fresh water is mixed, and from the plots of standard deviations Žsee Fig. 2Ž2. and Fig. 4Ž2.., we note that there are relatively small variations below 50–100 m depth. At the western boarder of Skagerrak, section H ŽFig. 1., the 35.0 p.s.u. contour of the average salinity field meets the Norwegian coast at 130 m depth. From the standard deviations we see that this depth typically varies between 90 m to 200 m. In the central Skagerrak Žsee Fig. 2Ž1 and 2.., the variability is smaller and the depth of the 35.0 p.s.u. contour will typically be 60 m " 20 m. The density structure and the variations in the density structure are reflected in the geostrophic transports Žsee Tables 1–4..

323

There are a number of processes causing the down-mixing of surface fresh water and with the present resolution horizontally and not least vertically we cannot expect to resolve these processes. The best we can hope for is to obtain the effects of the sub grid scale processes on the large scale fields by choosing ‘good’ parametrizations of the processes. In the following we therefore particularly study the effects of choosing different parametrizations of sub-grid scale processes. 4.7. Model results— POM In the following, the major results from the numerical experiments performed with the POM are given. A more detailed description is given in Berntsen et al. Ž1996.. The Blumberg and Mellor model applies a turbulence closure model due to Mellor and Yamada Ž1982. for the turbulence kinetic energy and the turbulence macroscale. The vertical eddy viscosity, K M , and vertical eddy diffusivity, K H , are computed from these turbulence parameters. However, the actual values used in the computations are K M q K MIN1 and K H q K MIN2 , respectively and the model results are often very sensitive to the choice of the background values K MIN1 and K MIN2 . The same applies to the horizontal dimensionless viscosityrdiffusivity parameters C M and C H which

Table 4 Area integrated measures for the upper 50 m of section H-comparisons datarmodel results-BOM Ž1 Sv s1 Sverdrups1=10 6 m3 sy1 . Data S˜

RUNB1

RUNB2

33.62 p.s.u.

34.51 p.s.u.

34.54 p.s.u.

S SD

0.83 p.s.u.

0.29 p.s.u.

0.40 p.s.u.

T˜

&

9.418C

8.678C

8.428C

T SD GTout GTin SD GTout GTinSD

0.978C 0.100 Sv 0.075 Sv 0.045 Sv 0.016 Sv

0.368C 0.073 Sv 0.020 Sv 0.020 Sv 0.012 Sv

0.338C 0.070 Sv 0.015 Sv 0.026 Sv 0.012 Sv

DS

&

)

1.16

1.52

DT

&

)

2.24

2.73

DG out

&

)

0.60

0.67

DG in

)

3.44

3.75

&

&

324

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

are chosen to be C M s C H s 0.1 in the experiments with POM presented below. In RUNP1 the parameter settings are as in the public domain version of the Blumberg and Mellor code, that is K MIN1 s K MIN2 s 2 = 10y5 m2 sy1 . In

RUNP2 K MIN2 is changed. In Tables 5 and 6, the area averages for sections F and H of mean salinity and standard deviations are presented together with ˜ The geostrophic area averages of the discrepancies S. transports given in Tables 5 and 6 are based on the

Fig. 4. Ž1. Mean values of observed salinity in p.s.u. for section H.Ž2. Standard deviations of observed salinity in p.s.u. for section H. Ž3. Mean values of model salinity in p.s.u. for section H from RUNP1. Ž4. Mean values of model salinity in p.s.u. for section H from RUNP2. Ž5. Mean values of model salinity in p.s.u. for section H from RUNB1. Ž6. Mean values of model salinity in p.s.u. for section H from RUNB2.

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

325

Fig. 4 Žcontinued..

average salinity and temperature fields for the period 24 May to 20 June 1990. In the averages from data, the measurements taken along section H on 5 June 1990 were not included due to probable errors. Average model transports inrout of the sections are also given. We have not tried to locate a selection of parameters for this model that minimize our error measures. The focus is on demonstrating the effects of varying the parameters. In a recent study, Skogen et al. Ž1996. demonstrate that it is possible to select parameter settings that reduce the error measures. However, they still obtain some high salinity deep water at section F. 4.7.1. Results for RUNP1 From Table 5 we note that the time average model salinity field is more than 12 standard deviations off the time average observed salinity field as an area average across section F. From Table 6 we see that the corresponding number for section H is more than 5. One major reason for a larger discrepancy across section;F is that the area average standard deviation, S SD , is approximately 0.15 p.s.u. across section F whereas it is 0.33 p.s.u. across section H, and larger standard deviations give smaller

discrepancies with the present error measure. Com; paring the average geostrophic transports, DGo ut and ; DGi n, given in Table 5, we note that despite the large discrepancy in the salinity field the integrated internal model pressure out of section F may be in good agreement with the corresponding observed internal pressure. The model geostrophic transport into the section is on the other hand much larger than the geostrophic transport estimated from observations. For section H we note from Table 6 that the model geostrophic transport out of section H is approximately five times larger than the corresponding value from observations. To understand the discrepancies noted above we compare the average model salinities for RUNP1 for sections F and H, Fig. 2Ž3.Fig. 4Ž3., with the corresponding observed fields, Fig. 2Ž1.Fig. 4Ž1.. For section F we find: Ža. The salinity of the deep model water is much too high, up to almost 37.0 p.s.u. at the bottom. Žb. The model surface water is too saline. Žc. The intermediate model water is too fresh. The model 35.0 p.s.u. contour is typically more than 100 m too deep. Žd. The vertical near surface salinity gradient is too weak.

326

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

Fig. 5. Ž1. Mean values of observed salinity in p.s.u. for section H-upper 50 m. Ž2. Mean values of model salinity in p.s.u. for section H from RUNB1-upper 50 m. Ž3. Mean values of model salinity in p.s.u. for section H from RUNB2-upper 50 m.

The bottom water of section F deserves a comment. From Fig. 2Ž1., we note that this water is less saline than the intermediate water masses. This bottom water is formed probably at the northern North Sea plateau during cold winter conditions ŽSvendsen et al., 1991.. This deep water is not represented in the initial values of the model and the processes generat-

ing it are not well represented in the model. Therefore, we cannot expect to reproduce this observed feature. Model salinity field water masses with salinity above 35.5 p.s.u. we associate with the advection scheme. The numerical model applies the leapfrog scheme for advection both horizontally and verti-

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

327

Fig. 6. Ž1. Bottom topography of the North Sea model area. Ž2. Bottom topography of the Skagerrak model area.

cally, and it is well known that this scheme is non-monotonic. Especially we have seen that during upwelling heavy model water is generated in the bottom cells near the coast. This water sinks out and fills up Skagerrak. For section H, the comments a, b, c and d for section F qualitatively apply. The salinity of the model bottom water is almost 36.0 p.s.u. In the model salinity the horizontal gradients are too large, due to excessive downmixing of fresh water near the Norwegian coast, giving an unrealistic strong baroclinic forcing along the coast. There are therefore in our opinion two major reasons for the large discrepancies between the observed and the model fields. First, the use of an oscillatory advection scheme and excessive downmixing of fresh water masses. The second reason is probably connected to the first. 4.7.2. Results for RUNP2 In RUNP2 we have increased the minimum value of the vertical diffusivity constant, K MIN2 , near the bottom and in the upper 5 m of the surface layer to 10y3 m2 sy1 . The minimum allowed vertical diffusivity is reduced to 10y7 m2 sy1 in the intermediate waters.

From Tables 5 and 6, we note that the changes in the model parameters have reduced the discrepancy measure from 12.90 to 6.61 at section F whereas for section H there is an increase from 5.74 to 7.60. For the estimated geostrophic transports the most noticeable feature is again the large model geostrophic transport out of section H. Thus, in average we found an improvement in the discrepancy measures even if for one section the quality of the model field apparently was reduced. That a change in model parameters, numerical techniques, initial andror boundary values causes improvements in some quality measures and degrada-

Table 5 Area integrated measures for the whole section F-comparisons datarmodel results-POM Ž1 Sv s1 Sverdrups1=10 6 m3 sy1 .

S˜

Data

RUNP1

RUNP2

34.56 p.s.u.

35.07 p.s.u.

34.71 p.s.u.

S SD GTout GTin

0.15 p.s.u. 0.899 Sv 0.399 Sv

0.19 p.s.u. 0.882 Sv 0.734 Sv

0.12 p.s.u. 1.208 Sv 0.370 Sv

DS Model transports

) )

12.90 3.051 Sv

6.61 2.644 Sv

&

&

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

328

Table 6 Area integrated measures for the whole section H-comparisons datarmodel results-POM Ž1 Sv s1 Sverdrups1=10 6 m3 sy1 . Data

RUNP1

RUNP2

34.59 p.s.u.

34.66 p.s.u.

34.62 p.s.u.

S GTout GTin

0.33 p.s.u. 0.639 Sv 0.255 Sv

0.24 p.s.u. 3.328 Sv 0.087 Sv

0.19 p.s.u. 3.038 Sv 0.086 Sv

DS Model transports

) )

5.74 1.765 Sv

7.60 1.709 Sv

S˜

& SD

&

tions of others is a very typical situation. We would primarily like to see quality measures for a wide range of sections. If the results for all sections are better for one choice of model set up, it will be easy to make a choice. When the results from one section to another are conflicting, as in the present case, it will be much more difficult to conclude. A reason for conflicting results may be that we in our model changes not yet have focused on the main source of model errors and that we cannot expect a general improvement across a range of sections. To understand the changes in the discrepancy measures we focus on plots of the average salinity fields. For section F the main reason for a large discrepancy measure in RUNP1 was the very high salinities in the deeper parts of the section. From Fig. 2Ž4., we note that the increased minimum vertical diffusivities stops the generation of heavy water efficiently. The maximum salinity of the bottom water is now slightly above 35.1 p.s.u., and the quality measure is clearly improved. For section H the main reason for discrepancies in RUNP1 was perhaps the excessive downmixing of fresh water along the Norwegian coast. Comparing Fig. 4Ž3 and 4. we note that the change in model parameters has further increased this downmixing bringing the average model salinity field further away from the corresponding observed field. 4.8. Model results— BOM The main reason for the discrepancies between the model results and the data seen in the previous section were considered to be the choice of inappropriate numerical techniques and mainly the use of the leapfrog scheme for advection. To follow this up

the code documented in Berntsen et al. Ž1996b. was developed. The evaluation strategy was also somewhat extended including discrepancies between model and observed temperature fields. To maintain density fronts as well as possible, the horizontal diffusivity coefficient C H is chosen to be 0.0. Non-zero values have been tried, but the quality of model results measured by our suggested quality measures slowly degrades as C H is increased. The horizontal eddy viscosity coefficient C M is chosen to be 0.5. Two options on how to represent vertical exchange processes have been tried. In RUNB1 an implementation of the Mellor and Yamada Ž1982. 2 1r2 level model with the modifications due to Galperin et al. Ž1988. to compute K M and K H is applied. At the surface the boundary condition for the turbulent macroscale, l, is related to the wave height through the wind speed, u10 , 10 m above the 2 Ž sea surface by the equation: l s 0.0246 u10 SWAMP group, 1985.. In RUNB2 an implementation of the Richardson number formulation due to Munk and Anderson Ž1948. to compute K M is applied. In addition scalar fields, salinity and temperature, are swapped vertically whenever E rrEz ) 0. The background vertical eddy viscosity K MIN1 is 2 = 10y5 m2 sy1 and the background vertical eddy diffusivity K MIN2 is 5 = 10y8 m2 sy1 in the experiments with the BOM. 4.8.1. General remarks— BOM From Tables 1, 3, 5 and 6 we note that for the salinity fields the discrepancy measures for the BOM are smaller than the corresponding numbers for the ; POM. For section F DS are 4.36 and 2.78 in the two experiments reported with the BOM and 12.90 and 6.61 in the two experiments with the POM. For section H, 0.75 and 0.89 for the BOM should be compared to 5.74 and 7.60 for the POM. This is to the present authors a clear message that numerics matters and that especially the choice of advection scheme is of great importance. We note that our discrepancy measures for salinity are in general much betterrsmaller than the measures for temperature. This is due to small standard deviations in the deeper waters for the temperature combined with uncertainties in the initial and boundary values for model temperature fields and lack of data on heat fluxes.

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

However, measured in degrees Celcius, the discrepancies are typically less than 18, and we do not believe that errors in the temperature field have a major effect on the density structure and the internal pressure field. For salinity we have at least more realistic forcing and this is clearly reflected in the discrepancy measures. The model geostrophic transports are generally in reasonable agreement with the geostrophic transports from the observed fields indicating reasonable density structures. One exception is the model geostrophic transport of the upper 50 m of section H into Skagerrak Žsee Table 4.. Here the model transports are much smaller than the geostrophic transports from the observed fields due to too little vertical mixing of surface model waters towards the coast of Denmark. In the average observed salinity field across section H a core of higher salinity water is found close to the Danish coast. This is a filament of North Sea water transported southward due to a relatively fine scale feature in the bottom topography Žsee Fig. 1 and Danielssen et al. Ž1996... This is missing in all corresponding model fields due to the lack of this feature in the model topography. For section F the discrepancy measures for salinity are in general much poorer than the corresponding measures for section H. One reason is that the standard deviations of the measured salinity field are typically much smaller across section F. The water masses of the central and deeper parts of section F may also have a much longer flushing time, and a 2-month spin up time is probably too short for this section. For section F, we thus note that the model is able to reproduce the inflow of 35.2 p.s.u. water to Skagerrak, but the large core of 35.2 p.s.u. water between 250 m and 400 m depth is not reproduced. 4.8.2. Results— RUNB1 In RUNB1, the Mellor–Yamada turbulence closure scheme was applied. For section F, the least saline surface model water is 31.0 p.s.u. ŽFig. 3Ž2.. while in the average observed field we find water with salinity 24.0 p.s.u. The minimas in average surface salinity we find off the Norwegian coast both in the observed field and in the model field, but the minimum in the model field is further from the coast than the minimum in the observed field. This indicate that the core of the Norwegian Coastal Current

329

at this part of Skagerrak is located off the Norwegian coast at least during SKAGEX I. By comparing Fig. 2Ž1,2 and 5., we note that the 35.0 p.s.u. and 35.1 p.s.u. contours of the model fields meet the Norwegian coast at reasonable depths. These contours, however, are in the central parts of the section found much deeper in the model fields than in the observations. At section H, we note that the freshest surface water of the average model field is 33.0 p.s.u. while the corresponding observed field is 30.0 p.s.u. This indicates excessive mixing of the freshwater added to the system in the model. 4.8.3. Results— RUNB2 In RUNB2, the Richardson number formulation is applied to produce K M , whereas K H is chosen to be constant and equal to 5 = 10y8 m2 sy1 . In addition we let two vertically adjacent cells exchange scalar properties whenever the density of the upper cell is greater than the density of the lower cell. Comparing the quality measures for the two experiments with the BOM for section F, we note that ; the integrated measures DS are reduced both for the whole section and the upper 50 m Žsee Tables 1 and 2. when replacing the Mellor–Yamada scheme with the simple Richardson number-based formulation. The message is, however, not clear. From Tables 3 and 4 we note that the integrated quality measures for section H are somewhat larger in RUNB2 than in RUNB1. The reason for apparently conflicting results may again be that there are other sources of errors that are more important than the choice of turbulence closure scheme in the vertical. Coarse spatial resolution especially in the vertical is one good candidate. From the quality measures for the two reported experiments with the BOM, we may conclude that with the present model inputs and the present model resolution it does not lead to an overall improvement of model results to replace a simple Richardson based formulation with the Mellor–Yamada scheme. And unless the more complicated and time consuming technique gives a clear improvement, one should stick to the simple, which in this case is the Munk and Anderson Ž1948. formulation. From Fig. 3Ž5.Fig. 5Ž5., we note that we find 20.0 water at the surface in the central Skagerrak at

330

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

section F, and 26.0 water at the surface near the coast of Denmark at section H. This shows that it is possible with the present coarse resolution to avoid the excessive mixing of the surface waters seen in the previous experiments by using this very simple formulation for vertical exchange processes. The reason is that the values of K M tend to be much larger with this scheme than with the Mellor and Yamada formulation. The momentum at the surface produced by wind stress, is therefore ‘diffused’ more rapidly to the deeper layers and we get less convergencerdivergence in the surface flow field and less vertical exchange due to convection.

5. Summary The measurements gathered along the repeated sections across Skagerrak during SKAGEX I have been used to produce average fields for these sections over the observation period. A technique for quantifying the effects of changes in the numerical schemes, model parameters andror model forcing is suggested. Area averages of absolute differences between model outputs and observations are for Skagerrak completely dominated by the discrepancies in the upper few meters of the ocean. Therefore, discrepancies between time averaged values from the observed fields and model outputs are related to the natural variability of the fields. The numbers produced with the suggested measure are relative numbers that will be specific for each section and for each series of observation. Ideally we would therefore like to see the measures computed for a number of sections for various models and choices of model parameters in order to assess model skills. To understand the discrepancies and relative changes when changing the model, plots of the time averaged fields are necessary. In our studies we have used the measures produced by this technique in combination with plots of average cross section fields in the search for better models for this area. The search will continue in experiments where better resolution both horizontally and vertically and improved initial values and boundary conditions hopefully will lead us towards a better agreement with the data. There are numerous reasons for lacking agreement between observed and modelled fields. In our

study we have demonstrated that at least for the present model area the choice of numerical algorithms strongly affects the quality of the model outputs. It may be argued that the algorithms and parameters tried for representing sub-grid scale effects are different for the two models and that this may explain the differences. The choices of parametrizations of sub-grid scale effects are, however, strongly related to the choices of numerical algorithms and the choices should and cannot be done independently. The differences between the algorithms for representing vertical subgrid scale effects are more fundamental in the experiments performed with the Bergen Ocean Model than in the experiments with the Princeton Ocean Model where only the background diffusivities are varied. Even so the model results from the BOM show less variations and we feel more comfortable with the model that is least sensitive to the representations of sub-grid scale effects. In our opinion the major source of errors when applying the POM for Skagerrak studies is the leapfrog advection scheme. The under-rover-shooting associated with it will easily cause oscillations in the pressure fields and instabilities near fronts. To apply gradient preserving, monotonic advection schemes in areas like the Skagerrak is therefore essential. Using the suggested quality measures we have been able to identify clearly that the quality of the model temperature field is inferior to the quality of the model salinity field in our model set up for Skagerrak. This points us towards putting more efforts in acquiring realistic heat fluxes and better initial values and boundary conditions for temperature. From the experiments with the BOM, it is far from obvious that the quality of the model results improves by going from simple Richardson numberbased techniques to higher order turbulence closure with the present vertical resolution.

Acknowledgements We like to thank the SKAGEX participants for making their valuable dataset available for this model evaluation. We also thank two anonymous referees for valuable remarks and suggestions.

J. Berntsen, E. SÕendsenr Journal of Marine Systems 18 (1999) 313–331

References Berntsen, J., Skagen, D.W., Svendsen, E., 1994. Modelling the transport of particles in the North Sea with reference to sandeel larvae. Fish. Oceanogr. 3, 81–91. Berntsen, J., Svendsen, E., Ostrowski, M., 1996. Validation and sensitivity study of a sigma-coordinate ocean model using the SKAGEX dataset. International Council for the Exploration of the Sea. CM 1996rC:5, 28 pp. Berntsen, J., Skogen, M.D., Espelid, 1996b. Description of a sigma-coordinate ocean model. Fisken og Havet 12. Institute of Marine Research. Postbox 1870, N-5024 Bergen, Norway. Blumberg, A.F., Mellor, G.L., 1987. A description of a three-dimensional coastal ocean circulation model. In: Heaps, N. ŽEd.., Three-Dimensional Coastal Ocean Models, Vol. 4. American Geophysical Union. Cox, M.D., Bryan, K., 1984. A numerical model of the ventilated thermocline. J. Phys. Oceanogr. 14, 674–687. Danielssen, D.S., Edler, S., Fonselius, S., Hernroth, L., Ostrowski, M., Svendsen, E., Talpsepp, L., 1996. Oceanographic variability in SkagerrakrNorthern Kattegat, May–June 1990. Submitted to ICES Journal of Marine Sciences. Dybern, B.I., Danielssen, D.S., Hernroth, L., Svendsen, E., 1994. The Skagerrak EXPERIMENT-SKAGEX report 1988–1994. Technical Report 635, TemaNord 1994:635. Nordic Council of Ministers, Store Strandstraede 18, DK-1255 Copenhagen K. Galperin, B., Kantha, L.H., Hassid, S., Rosati, A., 1988. A quasi-equilibrium turbulent energy model for geophysical flows. J. Atmos. Sci. 45, 55–62. Hackett, B., Røed, L.P., Gjevik, B., Martinsen, E.A., Eide, L.I., 1995. A Review of the Metocean Modeling Project ŽMOMOP.: Part 2. Model Validation Study. In: Lynch, D.R., Davies, A.M. ŽEds.., Quantitative Skill Assessment for Coastal Ocean Models. American Geophysical Union, Washington, DC, pp. 307–327. Lynch, D.R., Davies, A.M., Gerritsen, H., Mooers, C.N.K., 1995. Closure: quantitative skill assessment for coastal ocean models. In: Lynch, D.R., Davies, A.M. ŽEds.., Quantitative Skill Assessment for Coastal Ocean Models. American Geophysical Union, Washington, DC, pp. 501–506. Martinsen, E.A., Engedahl, H., 1987. Implementation and testing of a lateral boundary scheme as an open boundary condition in a barotropic ocean model. Coastal Eng. 11, 603–627. ˚ Martinsen, E.A., Engedahl, H., Ottersen, G., Adlandsvik, B.,

331

Loeng, H., Balino, B., 1992. Climatological and hydrographical data for hindcast of ocean currents. Tech. Rep. 100, The Norwegian Meteorological Institute. Mellor, G.L., Yamada, T., 1982. Development of a turbulence closure model for geophysical fluid problems. Rev. Geophys. Space Phys. 20, 851–875. Munk, W.H., Anderson, E.R., 1948. Notes on a theory of the thermocline. J. Mar. Res. 7, 276–295. Reistad, M., Iden, K., 1995. Updating, correction and evaluation of a hindcast database of air pressure, wind and waves for the North Sea, the Norwegian Sea and the Barents Sea. Research Rep. 9. The Norwegian Meteorological Institute, 35 pp. Røed, L.P., Hackett, B., Gjevik, B., Eide, L.I., 1995. A Review of the Metocean Modeling Project ŽMOMOP.: Part 1. Model Comparison Study. In: Lynch, D.R., Davies, A.M. ŽEds.., Quantitative Skill Assessment for Coastal Ocean Models. American Geophysical Union, Washington, DC, pp. 285–305. Skogen, M.D., Svendsen, E., Berntsen, J., Aksnes, D., Ulvestad, K., 1995. Modelling the primary production in the North Sea using a coupled three-dimensional physical–chemical–biological ocean model. Estuarine, Coastal and Shelf Sci. 41, 545– 565. Skogen, M.D., Svendsen, E., Ostrowski, M., 1996. Quantifying volume transports during SKAGEX with the Norwegian Ecological Model System. Submitted to Continental Shelf Research. Smagorinsky, J., 1963. General circulation experiments with the primitive equations: I. The basic experiment. Monthly Weather Rev. 91, 99–164. Stigebrandt, A., 1980. Barotropic and Baroclinic response of a semi-enclosed basin to barotropic forcing of the sea. In: Freeland, H.J. et al. ŽEds.., Fjord Oceanography. Plenum, New York, pp. 141–164. ˚ Svendsen, E., Berntsen, J., Skogen, M., Adlandsvik, B., Martinsen, E., 1996. Model simulation of the Skagerrak circulation and hydrography during SKAGEX. J. Mar. Syst. 8, 219–236. Svendsen, E., Saetre, R., Mork, M., 1991. Features of the northern North Sea circulation. Continental Shelf Res. 11, 493–508. SWAMP group, 1985. Sea wave modelling project ŽSWAMP.. An intercomparison study of wind wave predictions models: Part 1. Principal results and conclusions. In: Ocean Wave Modelling, Plenum, New York, p. 256. Sweby, P.K., 1984. High resolution schemes using flux limiters for hyperbolic conservation laws. SIAM J. Numer. Anal. 21, 995–1011.