# Secondary flows induced by wind forcing in the Rhône region of freshwater influence

Guillaume Reffray, Philippe Fraunié, Patrick Marsaleix

► To cite this version:

Guillaume Reffray, Philippe Fraunié, Patrick Marsaleix. Secondary flows induced by wind forcing in the Rhône region of freshwater influence. Ocean Dynamics, 2004, 54 (2), pp.179-196. 10.1007/s10236-003-0079-y. hal-00138730

HAL Id: hal-00138730

https://hal.science/hal-00138730v1

Submitted on 5 Oct 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Secondary flows induced by wind forcing in the Rhône region of freshwater influence

Guillaume Reffray · P. Fraunié · P. Marsaleix

**Abstract** Secondary flows induced by the blocking effect of a river plume on a transverse upwelling are investigated in a microtidal region of freshwater influence (ROFI). A nested version of the SYMPHONIE primitive-equation free-surface model for 3-D baroclinic coastal flows has been developed for the Rhône ROFI. The main characteristics of the model are a generalized sigma coordinate system in finite differences, using a time splitting for external and internal modes and high-order numerical advection schemes for density fields in combination with an modified turbulence closure scheme. The nesting system consists of two grids forced by the high-resolution ALADIN model atmospheric data. The coarse grid of 3 km resolution for the whole Gulf of Lions allows the forcing of the Liguro-Provençal large-scale current when the fine mesh of 1-km resolution is centred on the river mouth of the Grand Rhône. Documented field experiments from the Biodypar 3 field campaign performed during March 1999 are used for validation. Numerical results, CTD profiles and a SPOT TSM visible image are in good agreement concerning the shape and structure of the river plume. Other coastal flow features can be observed from satellite imagery. Computations of realistic situations recover these main secondary structures. Complementary process-oriented runs give an explanation of how the coastal upwelling induced by an inhomogeneous offshore wind is destabilized by the combination of the river plume and along-shelf current-blocking effects. In the end, a factor-separation analysis provides evidence that the locally

non-linear effects in momentum contribute to the occurrence of secondary vortices.

**Keywords** Numerical simulation · Symphonie model · Factor-separation study

## 1 Introduction

This paper deals with the higher-order numerical simulation of non-linear interactions in coastal flows, especially concerning the secondary structures occurring in microtidal sea when a supercritical river plume interacts with an upwelling and a background shelf circulation.

The river plume circulation has been extensively observed (Garvine and Monk 1974; Garvine, 1977) and numerically investigated for basic situations (Garvine, 1982, 1987; Chao and Boicourt 1986; Chao, 1988a,b, 1990; Oey and Mellor, 1993, Xing and Davies, 1999). These observations and numerical studies, using either a two-layer model (Garvine 1982) or a primitive-equation model (Chao and Boicourt, 1986), established the main features of river plume flow.

Firstly, in the absence of wind and tide, the brackish water is deflected to the right by the Coriolis force in the Northern Hemisphere, reaches the coast and abruptly, the anticyclonic bulge is extended by a coastal jet. This anticyclonic circulation is characterized by a rising motion when, in contrast, a significant downwelling effect is occurring in the region of transition between the bulge and the coastal jet. For supercritical river plumes, when the internal Froude number is more than 1, the plume body is separated from the ambient heavier water by a horizontal density front which induces a dynamical front.

This rapidly responding surface flow is very sensitive to external forcings at larger space and time scales, mainly the wind stress, the ambient circulation and tidal pulsation. Then, non-linear three-dimensional interactions occur in the complex coastal bathymetry. For an upwelling favourable wind, the river plume is swept

<footnote>G. Reffray (✉) · P. Fraunié
Laboratoire de Sondages Electromagnétiques de l'Environnement
Terrestre CNRS, Université de Toulon et du Var,
BP 132 83957 La Garde cedex
e-mail: Guillaume.Reffray@lseet.univ-tln.fr</footnote>

<footnote>P. Marsaleix
Laboratoire d'Aérologie CNRS,
14 Avenue Edouard Belin,
31400 Toulouse</footnote>
offshore and the vertical mixing, weak in the region close to the river mouth, suddenly increases downstream as the Richardson number decreases, and the fluid tends towards homogeneity. Furthermore, when the tide is no longer negligible, the plume pulses at the tidal frequency and a comparable mixing process occurs on ebb tide (Simpson 1997). For downwelling favourable winds, the coastal jet deepens and accelerates.

In addition, the forcing of the large-scale background circulation on the river plume spreading has very recently been investigated. Xing and Davies (2002) studied the interaction processes between a river plume, an along-shelf current and an offshore wind on an idealized site from the Ebro delta shelf in Balearic sea. Change of the shelf width and Mistral wind from the land were driving the interaction between the river plume and the along-shelf current, providing an offshore flow and an anticyclonic surface eddy to the south, i.e. the right side of the river plume.

Next, realistic plume modelling was investigated on different sites (Kourafalou et al. 1996; Turner and Rabalais 1994; Ruddick et al. 1994, 1995), but especially in situations where the tidal background turbulent mixing is important.

More recently, the whole region of freshwater influence (ROFI), as introduced by Simpson (1997), an area larger than the river plume where salinity-induced baroclinic circulation dominates the thermal one, has been investigated in order to better describe the mixing processes in the surface layer. In particular, the time and spatial scales for river inputs of suspended matter and nutrients in the coastal zone are of special interest for biological bloom prediction and interpretation of the chlorophyll satellite pictures.

The Rhône river plume flows into the Gulf of Lions (south of France) in the western Mediterranean microtidal sea where the main forcings are the mesoscale Liguro-Provençal Current (hereafter designated LPC), part of the northern Mediterranean current, mostly flowing southwestward along the shelf break and dominant winds, especially both the Mistral and the Tramontane strong landwinds (Millot 1990).

Several studies have been performed of this site, by both field campaigns and remote sensing (Demarcq and Wald 1984; Forget et al. 1990) together with modelling. Up to now, these studies mainly addressed the Rhône river plume, for both academic situations (Marsaleix et al. 1998) and more realistic situations (Estournel et al. 1997; Arnoux-Chiavassa 1998; Arnoux-Chiavassa 2003).

In agreement with the usual classification table (Garvine 1985, Chao 1988a), the Rhône river plume is of a supercritical type, in which internal Froude number is larger than 1. As a consequence, the occurrence of horizontal density fronts and the vertical layering effect become major processes which are difficult to reproduce in computations, so higher-order advection schemes in combination with adapted turbulence closures have to

be carefully chosen (Arnoux-Chiavassa 1998; Xing and Davies 1999; Arnoux-Chiavassa 2003).

In particular, the vertical mixing processes occurring for large Richardson numbers have been investigated (Arnoux-Chiavassa, in press) in the Rhône river plume using a SMART (Sharp and monotonic algorithm for realistic transport, Gaskell and Lau 1988) scheme for both density and velocity fields to capture the density fronts and permitting the use of optimized anisotropic turbulent closures (Lerredde et al. 1999).

The work described herein extends these previous studies to the influence of the large-scale circulation and the time and spatial variability of the wind which is specified as the outputs of the high-resolution (0.1°) Aladin meteorological model. Moreover, a nested modelling strategy was employed to account for the LPC and the background shelf circulation and thus determine the time-dependent open-sea boundary conditions in a fine-grid computational domain centred on the Rhône delta. These boundary conditions are interpolated from the 3-D coarse-grid outputs, where initial and boundary conditions are fixed using a quasi-geostrophic simplified model for the coarse-grid large domain (Estournel et al. 2003).

A realistic model computation has been made for a situation where satellite images and salinity profiles were available from a field campaign. In addition, basic situations removing the river discharge or the external circulation forcings have been investigated in order to obtain a better understanding of this complex coastal flow.

## 2 The improved Symphonie model

### 2.1 The basic model

The Symphonie model used in this study is described in Estournel et al. (1997) and Marsaleix et al. (1998). It is a primitive-equation, free-surface model which solves conservation equations of momentum, mass, temperature and salinity, applying the Boussinesq approximation and hydrostatic equilibrium. Horizontal diffusion is retained for numerical stability and is of Laplacian form with a constant coefficient taken at $10 \text{ m}^2 \text{ s}^{-1}$ for the fine mesh size of 1 km (and $30 \text{ m}^2 \text{ s}^{-1}$ for the coarse mesh size of 3 km). Equations are discretized using an Arakawa C grid (Arakawa et al, 1983). An upward diffusing scheme is used for the time-stepping, while a mode-splitting technique (Blumberg and Mellor 1987) is used to compute the free surface elevation. Time steps for the external flows on the fine grid are about 10 and 200 s for the internal flow.

Generalized sigma coordinates (Pietrzak et al. 2002) with 30 levels were used (Fig. 1) for refinement of the vertical discretization in the surface layer where buoyancy effects in the river plume and the halocline in the ROFI are dominant. Preliminary tests (not shown) have been performed for comparison with a cartesian grid
<figure><img src="image_2.png" /><figcaption>**Fig. 1** Vertical levels using generalized sigma coordinate system</figcaption></figure>

(Arnoux-Chiavassa 2003) in order to make sure the pressure gradients are weakly sensitive to the coordinate system. Indeed, the vertical grid leads to a quasi-Cartesian coordinate system in the baroclinic surface layer.

## 2.2 The advection scheme

In order to accurately capture the large density gradients and fronts, we used Dimitrov's (2000) TVD advection scheme for a C-type grid which induces a weak implicit numerical diffusion. This scheme estimates the scalar field at each face of the cell using a quadratic interpolation. As these kind of schemes introduce undershoots and overshoots near the strong gradients, we apply the ultimate limiter (Leonard 1991) in order to restore the monotonicity of the solution.

## 2.3 Turbulence

The vertical diffusion coefficient was computed from turbulent kinetic energy $k$ and an algebraic mixing length $l$. Comparisons on river plume flows (Luyten et al. 1996; Xing and Davies, 1999) with other models like $k-\epsilon$ (Rodi 1987), where $\epsilon$ represents the turbulent dissipation, or $k-kl$ (Mellor and Yamada 1982) have shown that this type of model provides comparable results for a lower cost. We chose here the $k-l$ turbulence model presented in the Xing and Davies' 1999, where the mixing length $l$ can be locally adapted to stratified conditions. A simple algebraic form of $l$ is given by:

$$l_0 = 1 / \left( \frac{1}{l_1} + \frac{1}{l_2} \right), \qquad (1)$$

where

$$l_1 = K(\sigma H + z_0) \exp(\beta_1 \sigma) \text{ and } l_2 = K(H - \sigma H + z_s), \qquad (2)$$

with $K = 0.4$ the von Karman's constant, $\beta_1$ an empirical coefficient fixed to 2, $H$ the total water depth and $\sigma$ the normalized vertical coordinate. The boundary condition at the surface and at the bottom are controlled by $z_s$ and $z_b$, respectively, the surface and the bed roughness lengths are both fixed to $10^{-3}$ m.

The stratification is taken into consideration using the limiting condition of Galperin et al. (1989) applied to $l$:

$$l \le \frac{k_1 q}{N}, \qquad (3)$$

with $k_1$ a specified dimensionless constant, taking a value of 0.53 (Luyten et al. 1996), $q$ represents the turbulent velocity fluctuation scale ($k = 2q^2$) and $N$ the Brunt Vaisala buoyancy frequency.

The eddy viscosity and the eddy diffusivity, $A_v$ and $K_v$ respectively, are computed from

$$A_v = lqS_M \text{ and } K_v = lqS_H, \qquad (4)$$

with $S_M$ and $S_H$ the stability functions used by Galperin et al (1988, 1989).

## 2.4 The nested computational domains

The Rhône ROFI is a part of the Gulf of Lions, itself included in a synoptic system: the western Mediterranean Sea. The large structures of this latter system exhibit a space scale of about 1000 km and a duration from several days to even several months. Examples are the anticyclonic gyres produced by the instabilities of the Algerian Current (Beckers 1995) in the Baleares sea. The large-scale influence of the western Mediterranean circulation on the Gulf of Lions occurs through both Kelvin-type waves and the LPC (Echevin et al. in press).

In order to accurately reproduce the river plume spreading in the Rhône ROFI, a fine mesh of the order of 1 km was shown to be necessary (Arnoux-Chiavassa 1998). Indeed, the observed structures on a series of visible imagery exhibit a space scale of only several kilometres, smaller than the internal Rossby radius. The time scale of the local advective processes ranges from several hours (response of the river plume to an unsteady wind) to several days (local upwelling spin-up). At this stage, the low-order nesting method consists in forcing the open-boundary conditions on the fine grid using the fields recorded from outputs of the coarse-grid large model. The limiting factor of this direct nesting method is that no feedback from the fine-grid model is possible. So, the smallest structures cannot influence the
forcing fields. This method is easy to implement and numerically stable.

### 2.4.1 The coarse-grid large model

The 3-km horizontal resolution coarse grid extends 150 km offshore and fully covers the Gulf of Lions and the LPC. The bathymetry (Fig. 2a) exhibits the continental shelf and a shelf break including several canyons. The freshwater discharge in the gulf is mainly provided by two river branches, the grand Rhône and the petit Rhône, with additional inputs from three smaller sources, the Herault, the Orb and the Aude Rivers. The small nested fine-grid domain (Fig. 2b) centred on the Grand Rhône river mouth is marked by a square.

### 2.4.2 The fine grid

The velocity fields, the vertically integrated barotropic transport and the scalar fields (temperature, salinity and surface level) are recorded to be imposed on the two open boundaries of the 1-km resolution fine mesh, centred on the Rhône mouth (Fig. 2b) and extending only 40 km offshore. The southern open sea boundary is fixed on the shelf break and is influenced by the possible LPC penetration, whereas the western open boundary is mainly forced by the background continental shelf circulation induced by the wind. The bathymetry smoothly

<figure><img src="image_3.png" /><figcaption>**Fig. 2a, b** The two computational domains. **a** Representation of the whole Gulf of Lions. The isobaths are drawn every 50 m between 50 and 500 m and every 500 m between 500 m and 2500 m. **b** Representation of a domain focused on the Rhône river plume. The isobaths are drawn every 10 m between 10 and 100 m and every 50 m between 100 and 500 m</figcaption></figure>

decreases across the continental slope and the cone formed by the sediment fluxes in front of the river mouth is represented.

## 2.5 Initial and boundary conditions

For the two domains, zero flux condition was applied in the usual way on the closed boundaries. Other boundary conditions are presented hereafter.

### 2.5.1 Initial conditions and open-sea boundaries for the large domain

In the absence of large-scale validated operational OGCM outputs, the first difficulty for numerical experiments in this site is to accurately initialize a good representation of the LPC.

To overcome this problem, we use the simplified method proposed by (Estournel et al. 2003), which considers the steady linear diagnostic equations for the barotropic mode assuming a rigid lid linked to a surface pressures and introducing a barotropic stream function $\Psi$ for the vertically integrated velocity:

$$-f \frac{\partial \Psi}{\partial x} = -\bar{R} \frac{\partial P_s}{\rho_0 \partial \bar{x}} + \frac{\tau_{sx} - \tau_{bx}}{\rho_0} - H \left[ \frac{\partial}{\partial x} K_h \frac{\partial}{\partial x} \left( \frac{1}{H} \frac{\partial \Psi}{\partial y} \right) + \frac{\partial}{\partial y} K_h \frac{\partial}{\partial y} \left( \frac{1}{H} \frac{\partial \Psi}{\partial y} \right) \right] \quad (5)$$

$$-f \frac{\partial \Psi}{\partial y} = -\bar{R} \frac{\partial P_s}{\rho_0 \partial \bar{y}} + \frac{\tau_{sy} - \tau_{by}}{\rho_0} - H \left[ \frac{\partial}{\partial x} K_h \frac{\partial}{\partial x} \left( \frac{1}{H} \frac{\partial \Psi}{\partial x} \right) + \frac{\partial}{\partial y} K_h \frac{\partial}{\partial y} \left( \frac{1}{H} \frac{\partial \Psi}{\partial x} \right) \right], \quad (6)$$
where $H$ is the bottom depth, $K_h$ is the horizontal turbulent viscosity coefficient, $f$ is the Coriolis parameter, $\rho_0$ is the reference density, $\Psi$ is the stream function, $\tau_s$ and $\tau_b$ are, respectively, the wind and the bottom stresses and $\bar{R} = \int_{-H}^{0} R dz$ with

$$R(z) = \frac{\rho(z) - \rho(\text{bottom})}{\rho(\text{surface}) - \rho(\text{bottom})} \quad (7)$$

The $R(z)$ profile is determined from the density profile $\rho(z)$ as measured from CTD. Note that if $\rho(z)$ is constant on the water column, the basic shallow-water equations are recovered from the system of Eqs. (5) and (6).

Next, the stream function can be obtained from the curl of Eqs. (5) and (6) following the system:

$$A\Psi = B \quad (8)$$

solved on the C grid. $\Psi$ is a column matrix containing the stream function at each grid point. The square matrix $A$ contains the linear operators of the stream function. Finally, the column matrix $B$ contains the forcing terms, namely the curl of $\frac{\tau_s - \tau_b}{\rho_0 R}$ and the lateral boundary conditions.

The stream function becomes zero at the coast, so the only difficulty concerning this steady-state problem is to propose a good approximation of the stream function at the open sea. The assumption that the LPC mainly follows the continental slope, especially in the case of an offshore wind in which we are interested here, led us to

prescribe a linear bottom depth relationship. Estournel et al (2003) verified this choice from on-board ADCP current data. Thus, using this assumption, the stream function can be estimated just like the components of the vertically integrated velocity field and the free-surface elevation from the surface pressure from the diagnostic Eqs. (5) or (6). Next, using the vertical density profile $R$, the components of the baroclinic current are deduced. To complete the initial state, a perturbation of the density field can be deduced and then the salinity and temperature fields.

We apply this method on the coarse grid without wind forcing. A picture of the initialization of the vertically integrated velocity field is represented in Fig. 3. The LPC barotropic signal precisely follows the shelf-break bathymetry including canyons (Fig. 2a)

The first numerical tests (Estournel et al. 2003) using this method in a coarse-grid (3 km) model dealt with an analysis of the circulation on the continental shelf and the results were shown to be numerically stable. In particular, large structures of several dozen of kilometers, as, for instance, the two gyres on the continental shelf due to Mistral and Tramontane offshore wind forcing, were well reproduced, but the finest structures observed in the vicinity of the Rhône river mouth could not be recovered with precision.

During the run, the open-boundary conditions (velocity components, vertically integrated circulation

<figure><img src="image_4.png" /><figcaption>Fig. 3 Initial conditions on the coarse grid. The vertically integrated circulation</figcaption></figure>
and the scalar fields) are of steady Dirichlet type as the fixed values remain the initial computed values.

### 2.5.2 Initial conditions and open sea boundaries for the nested domain

Concerning the fine grid, two different types of simulation were performed. The first was the direct nesting method, in which the fields at the open boundary are simply interpolated from the 3-km coarse-grid model. A smooth interpolation procedure is applied to the whole nested domain concerning initial conditions.

Some test cases in the absence of external circulation have been performed on the fine-grid domain alone with a Neumann-type zero-gradient boundary condition applied to velocity components, salinity and temperature. The barotropic transport across the open boundaries was calculated using a radiation condition (Orlanski 1976).

### 2.5.3 The sea-surface and bottom-boundary conditions

The river plume dynamics is sensitive to the wind forcing (Chao 1988b). In addition, some particular river plume shapes can be successfully computed only if the high-frequency wind variability is taken into account (Estournel et al. 2001). For this, the high resolution of the 0.1° Aladin meteorological model has been chosen for the sea-surface forcing, including physical processes like wind stress and heat fluxes (evaporation and warming). The meteorological structure considered here (mainly Mistral and Tramontane southernward offshore winds) exhibits a quite large space scale, so a linear interpolation on the two grids is expected to be physically correct.

At the surface and at the bottom, the momentum fluxes were calculated according to:

$$ \rho_0 K \left( \frac{\partial u}{\partial z}, \frac{\partial u}{\partial \bar{z}} \right) = \begin{cases} (\tau_{sx}, \tau_{sy}) & \text{at the surface} \\ (\tau_{bx}, \tau_{by}) & \text{at the bottom} \end{cases} \quad (9) $$

where $\tau_s$ is the wind stress issued to the meteorological forcing model Aladin and $\tau_b$ the bottom stress calculated by a quadratic law (Blumberg and Mellor, 1987).

The turbulent energy at the surface and at the bottom is given by:

$$ q^2 = B_1^{2/3} u_*^2, \quad (10) $$

where $u_*$ is the local friction velocity.

For the scalar fields (temperature and salinity), a zero-flux condition is applied at the bottom. At the surface, the heat flux includes the warming effect as decomposed into four parts. The first part is the sensible heat ($Q_c$), the second the latent heat ($Q_e$) and finally, the third and the fourth components of the surface flux of temperature are due to the radiation flux for respectively short ($Q_s$) and long waves ($Q_l$) and so:

$$ \left( K_v \frac{\partial T}{\partial \bar{z}} \right)_{z=\eta} = (1 - T_r) Q_s + Q_l + Q_c + Q_e . \quad (11) $$

The part of the radiation flux in the domain of the short waves which penetrate is directly added in the heat conservation law as:

$$ \frac{\partial R_F}{\partial \bar{z}} = \lambda Q_s T \exp(\lambda z) , \quad (12) $$

with the transmission coefficient $T_r = 0.31$ and the attenuation coefficient $\lambda = 0.042 \text{ ms}^{-1}$ (Jerlov 1976).

The surface boundary condition for the salinity is

$$ \left( K_v \frac{\partial S}{\partial \bar{z}} \right)_{z=\eta} = E_s * S_z = \eta , \quad (13) $$

where $E_s$ is the rate of evaporation. The precipitation rate is not considered here mainly because this parameter is not well estimated by most of the atmospheric models. $E_s$ is performed from the latent heat flux:

$$ E_s = \frac{Q_e}{L_v} , \quad (14) $$

where $L_v$ is the latent heat for condensation.

### 2.5.4 River discharge

The Rhône mouth width is about 1 km, so the discharge takes place at a single grid point where the water column depth at rest is 2 m. For high-water discharge rates, a salt wedge does not occur and salinity was fixed to zero and temperature to 10 °C all across the river mouth. The barotropic transport $U$ was calculated on the basis of freshwater discharge flux $F$ by the following relationship:

$$ U = \frac{F}{l} , \quad (15) $$

where $l$ is the horizontal grid spacing. The same input condition was applied to the little Rhône river. The other rivers (used in the large domain), have runoff in the order of $50 \text{ m}^3\text{s}^{-1}$, which is considered as a point source at the surface level.

## 3 Modelling a realistic situation

### 3.1 The experimental data

The modelling strategy described above was applied to the field campaign Biodypar 3 (for ‘Biologie et Dynamique du Panache du Rhône’) performed from March 22 to April 2, 1999, and especially to the situation encountered on March 29, 1999. A sample of wind intensity and direction (Fig. 4) measured at the Tour du Valat located 5 km from the Grand Rhône river mouth, shows that this situation is characterized by a typical Mistral offshore wind condition in combination with a high Rhône River discharge (more than $2000 \text{ m}^3\text{s}^{-1}$; Fig. 5). CTD profiles and a primary SPOT image of suspended matter (Fig. 6), taken at 10.00 h, were made available for a qualitative comparison with numerical results.
<figure><img src="image_5.png" /><figcaption>**Fig. 4** Wind conditions sampled at several kilometres of the Grand Rhône river mouth during the Biodypar 3 field campaign</figcaption></figure>

### 3.1.1 Runs for a typical landwind situation

The numerical model was run on the coarse grid for 7 days from March 22, using Estournel et al (2003) initial and open sea boundary conditions as described above. The vertically integrated circulation on the shelf (Fig. 7) contemporary with the SPOT image shows the main flow structures in the open sea boundaries of the small nested domain (delimited by a square). The southern boundary of this small domain is marked by a weak intrusion of the LPC on the continental shelf, probably due to a relaxation effect for a decreasing wind strength. The western boundary is dominated by an anticyclonic large-scale gyre, which is a response to the curl of the strong, non-uniform Mistral wind blowing since the previous day (Fig. 8); (Estournel et al. 2003).

After this first step, providing the initial and 3-D boundary conditions, the run on the fine grid started from March 22, 05.00 h. If we represent the vertically integrated circulation for the situation of March 29,

10.00 h, on the whole fine-grid domain (Fig. 9), we can verify that the nesting method is stable and consistent. Large-scale features like the LPC instability, the anticyclonic gyre due to the river discharge and especially the anticyclonic background circulation in the vicinity of the western boundary, are well preserved from the coarse to the fine grid when secondary structures are added in the fine-grid domain flows as, for example, the Petit Rhône influence. In the following analysis of the situation encountered on March 29, 1999, the graphic representation will be zoomed in the Grand Rhône river mouth. However, Fig. 9 will provide important information for the understanding of the near-field ROFI system.

### 3.1.2 The plume structure

The surface velocity and salinity fields computed on March 29, 10.00 h (Fig. 10) show that the structure of

<figure><img src="image_6.png" /><figcaption>**Fig. 5** Grand Rhône discharge during the Biodypar 3 field campaign</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>**Fig. 6** SPOT satellite image of March 29, 1999 at 10:00 h</figcaption></figure>

the river plume is in good agreement with the SPOT TSM satellite image (Fig. 6). In particular, the location and structure of the two lateral fronts and the river plume width are well reproduced (note that the SPOT

imaged is rotated about 10° anticlockwise with respect to the computational domain).

There are two opposite mechanisms occurring under Mistral-type wind conditions. On one hand, the river plume is separated from the coast and its offshore spreading is significantly increased (Chao 1988b) and, on the other hand, the mixing between the freshwater surface layer and the underlayer marine water is increased offshore due to the wind friction acting on a weakly stratified flow (Xing and Davies 1999). This behaviour is clearly observed on both the SPOT image and the model outputs. The TSM picture exhibits marked fronts associated with high concentrations (saturated light colour) up to 8 km offshore, which correspond to the computed salinity field. The eastern front is sharpened by an additional offshore density current induced by the strong salinity gradient (see the confined salinity isolines) with an increasing Froude number when, in contrast, the western front is decelerated and rapidly dissipated with a subcritical Froude number induced by the density current flowing to the coast. Thus, an asymmetrical transition region from the supercritical river plume to a subcritical one is clearly observed about 8 km offshore (salinity isoline 34 and lighter TSM concentration), evidencing a sudden vertical mixing.

This vertical mixing mechanism is clearly noticeable on the vertical profiles of salinity and horizontal velocity (Fig. 11). The vertical salinity profiles have been sam-

<figure><img src="image_8.png" /><figcaption>**Fig. 7** Vertically integrated circulation on the continental shelf (bounded by 500-m depth isobath)</figcaption></figure>
<figure><img src="image_9.png" /><figcaption>**Fig. 9** Vertically integrated circulation in the fine grid (bounded by 200-m depth isobath. The arrows are drawn every two points)</figcaption></figure>

<figure><img src="image_10.png" /></figure>
**Fig. 10** Surface currents and contours of surface salinity (ppt) for March 29, 1999. The isohalines have been drawn every 3 ppt. The *five black points* represent the location of the CTD stations

<figure><img src="image_11.png" /></figure>

**Fig. 11** CTD profiles for the March 29, 1999. Each line represents a station from 1 to 5. The left panel represents the temperature profiles, the centre panel the salinity profiles and the right panel the velocity profiles. The measurements are represented as solid lines whereas the model outputs are dashed lines

<figure><img src="image_12.png" /></figure>
pled on five stations marked by black points on Fig. 10, using a quasi-Lagrangian strategy (Broche et al. 2000) following the surface layer along a drifting buoy trajectory.

The measured vertical salinity profiles rapidly rearrange in a two-layer type structure deepening offshore up to a transition point where a sudden vertical mixing occurs. The numerical results (dashed lines) Fig. 10 reproduce this mechanism only qualitatively because of the well-known failure of TKE-based closure model to reproduce staircase density vertical profiles as encountered for Richardson numbers greater than 1. However, computations (stations 2 and 3, Fig. 11) qualitatively reproduce the different mixing processes for momentum and scalars where the scalar profiles tend to a two-layer system, although the horizontal velocity fields present a smooth shear-type profile. Indeed, this important feature is due to the fact that for high values of the Richardson number (Posmentier 1977; Phillips 1977) the turbulent diffusion is not the same for the momentum and the scalar (salt or temperature) fields, i.e. there are strong variations of the turbulent Schmidt number (Linden 1980; Galperin et al. 1988; Berrabaa et al. 2000). Moreover, in this case, the two-layer density system leads to a significant reverse undercurrent.

Downstream (station 5), the observed salinity profile points out a transition to a quasi-homogeneous water column when the Richardson number decreases. The computed data are not able to reproduce such a transi-

tion due to the continuous closure laws in the turbulence model. Previous studies (Arnoux-Chiavassa 2003) have shown that a carefully fitted conditional turbulence closure model was necessary to account for such a discontinuous behaviour of the vertical mixing, depending on the Richardson number. However, the smoothed mixing levels obtained here using a single turbulence closure seem to be sufficient for operational modelling.

### 3.1.3 *The secondary flows*

The TSM image (Fig. 6) exhibits two clearwater patches (dark on the picture) along the western river plume front, between the river plume itself (marked in white colour by the surface suspended matter) and the loaded (white on the picture) nearshore zone due to resuspended matter from the bottom.

In the very surface level (Fig. 10), the coastal flow dynamics is dominated by the imposed wind-friction boundary condition but at 5 m depth (Fig. 12) the computed velocity field shows that these patches correspond to two counterrotating vortices generating a dipole. These two gyres are confined between the western river plume, the coast and the background anticyclonic shelf circulation (Fig. 9) that seems to feed the larger cyclonic gyre. This large gyre should then induce the counterrotating smaller one close to the river mouth. The computed dipole (Fig. 12) occurred about 6 h

<figure><img src="image_13.png" /><figcaption>**Fig. 12** Velocity field at 5 m depth for the situation of March 29, 1999</figcaption></figure>
earlier and will vanish several hours later due to un-
steady wind conditions. The locations of the computed
vortices and observed patches centres are surprisingly
close each other although the computed gyres are larger
than the observed clearwater patches, as they extend
from the western river plume front up to the coast.
However, this may be explained by an entrainment effect
of turbid water inside the vortices. Such entrainment
effects are also clearly observable on the downstream
western part of the subcritical river plume around the
larger cyclonic gyre. On the other side, away from the
eastern river plume front, another cyclonic vortex
occurring in the open boundary of the Gulf of Fos seems
to be generated by entrainment from the anticyclonic
background shelf circulation.

All these structures are also visible at 10 m depth (Fig. 13) and in the vertically integrated horizontal velocity field (Fig. 14). This implies that these vortices are of a barotropic type except they are not manifest at the surface because of the wind-stress boundary condition without feedback. We infer that they are induced by a baroclinic effect, like a single interaction between the river plume fronts and the local coastal upwelling, although they would be comparable in size to the baroclinic internal Rossby radius (a few km). From the previously noted features, many questions concerning associated physical mechanisms arise. Firstly, how is the dipole generated ? Is this a generic mechanism ? In order to investigate these questions, supplementary numerical

tests have been run. It seems that the background shelf
circulation has a non-negligible influence on the ROFI
dynamics. So, by removing the river discharge, only the
interaction between the wind-induced circulation, the
LPC and the bathymetry will be analyzed. Does the
dipole still occur ? If not, how is the background circu-
lation rearranged ? Secondly, two complementary tests
will be run without the LPC and large-scale circulation,
i.e. without nesting. Firstly, by running without river
outputs, the local coastal upwelling effect will be char-
acterized, secondly, by reintroducing the river discharge,
the local interaction between the coastal upwelling and
the river plume will be analyzed. What are the relevant
forcing effects inducing the occurrence of the dipole ?
What is the real influence of the LPC forcing ?

We present here three unrealistic test runs useful to
point out the comparative influence of the surface den-
sity flow induced by the river discharge, the LPC forcing
and the large-scale background circulation on the ver-
tically integrated circulation. Postprocessing analysis has
been performed following the factor-separation method
proposed by Stein and Alpert (1993).

**4 Factor-separation investigation (numerical tests)**

To understand the influence of each forcing factor on
the circulation offshore the Grand Rhône river mouth, a
process-oriented study was performed on the real

<figure><img src="image_14.png" /><figcaption>**Fig. 13** Velocity field at 10 m depth for the situation of March 29, 1999</figcaption></figure>
<figure><img src="image_15.png" /><figcaption>**Fig. 14** Vertically integrated circulation for the situation of March 29, 1999</figcaption></figure>

bathymetry, keeping the local wind stress effect. Indeed, the interaction between the Grand Rhône river plume and the mesoscale circulation (LPC) without any meteorological forcing will not be investigated here, as several previous studies (Marsaleix et al. 1998; Arnoux-Chiavassa in press) have shown that the Grand Rhône river plume, even if of supercritical-diffusive type, never reaches the shelf break and never meets the LPC by itself.

Thus, we will consider only two factors and their influence on the vertically integrated circulation will be analyzed using the factor-separation method (Stein and Alpert 1993). To apply this method to these two factors, we need four run outputs. The first reference one (Sr), corresponding to a realistic reference situation, when all the factors are taken into consideration, was presented in the last section (Fig. 14). The second (S0) considers only the local upwelling effect of the non-uniform wind stress on the real bathymetry without any other factor. The two other runs are performed from the basic (S0) case by adding, respectively, the river discharges (S1) and the LPC and large scale circulation (S2). When the LPC and the large-scale circulation are removed (simulations S0 and S1), the open sea boundary conditions on the fine-grid domain are the ones presented in Section 2.5.2. Only the cases (Sr) and (S2) use the nesting approach. The description of these numerical tests are condensed in Table 1.

## 4.1 The direct simulations

The vertically integrated horizontal velocity fields in the vicinity of the Grand Rhône river mouth, corresponding to the simulations described in Table 1, are presented together in Fig. 15.

This first run (S0: the single local upwelling effect) without the river discharges and external circulation considers only the effect of the given, horizontally sheared wind forcing on the fine grid of the real bathymetry. Three main large structures occurring on the shelf are induced by a barotropic perturbation of the local upwelling effect (see the full small domain in Fig. 16). Firstly, an anticyclonic gyre (called in the subsequent sections western gyre) is located along the western open sea boundary and fills the third western part of the domain. A second anticyclonic large structure (hereafter called the eastern gyre) is centred 15 km

<table><caption>**Table 1** Description of the four numerical tests</caption><thead><tr><td>Simulation description</td><td>Results: matrix name</td></tr></thead><tbody><tr><td>All factors (Fig. 14)</td><td>Sr</td></tr><tr><td>Local upwelling effect</td><td>S0</td></tr><tr><td>Local upwelling effect adding river discharges</td><td>S1</td></tr><tr><td>Local upwelling effect adding LPC and large-scale circulation</td><td>S2</td></tr></tbody></table>
<figure><img src="image_16.png" /><figcaption>b</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>c</figcaption></figure>

<figure><img src="image_18.png" /><figcaption>d</figcaption></figure>

<figure><img src="image_19.png" /><figcaption>**Fig. 15a–d** Vertically integrated circulation for the simulation a Sr; b S0; c S1; d S2</figcaption></figure>

southward of the Gulf of Fos. This gyre is kept away from the southern open sea boundary by the shelf break, due to the conservation of potential vorticity (Xing and Davies 2002). Between these two anticyclonic gyres, a third large structure is located (called in the subsequent sections central gyre), which is constrained to turn in the cyclonic sense. Secondary structures naturally appear in shallow waters near the eastern boundary and in the vicinity of the Gulf of Fos. Other secondary flows also appear along the southern open sea boundary because of the free-boundary condition and because of the shelf break. The three main structures are more clearly visible

on the zoom near the Grand Rhône river mouth (Fig. 15b).

Instead of the dipole observed in the (Sr) run, a divergence zone is now observed, generated by the central and eastern gyres. A part of the upwelled water is advected by the central gyre while the other part is diverted into the Gulf of Fos by the secondary flow.

The second (S1) test run is the simulation (S0), but adding the river discharges. This configuration corresponds to previous process-oriented computations introducing high-order numerical advection schemes (Arnoux-Chiavassa 2003) when the external flow was not accounted for. The circulation in the whole fine-grid domain (not given) is complicated by the influence of the freshwater discharges. Figure 15c shows the vertically
<figure><img src="image_20.png" /><figcaption>**Fig. 16** Vertically integrated circulation in the fine grid induced by local wind without the river discharges (S0). The arrows are drawn every two points</figcaption></figure>

integrated horizontal circulation magnified in the vicinity of the Grand Rhône river mouth. The main features are that the eastern gyre is strongly strengthened by the river discharge while the western gyre is conserved and the central gyre seems to degenerate into several small cyclonic vortices (4 in Fig. 15c). Nearshore, two of them take place in the same location as the dipole observed in the (Sr) run. The dominant eastern gyre confines the nearshore secondary cyclonic vortex inside the Gulf of Fos.

The third and last complementary test run (S2) is the same as the first numerical simulation (S0), but now adding the large-scale circulation output from the large-domain modelling including the LPC (thus introducing the nesting method) but without the river discharges. The vertically integrated horizontal circulation on the full continental shelf (not given) is very comparable to the (Sr) analogous situation including the river discharges (Fig. 7), except for only one difference concerning the eastern anticyclonic gyre which is significantly reduced in this case. Figure 17, directly comparable to the Fig. 16, shows the three large structures modified by the surrounding circulation taken into account by the nesting method. The western gyre is recovered as in the previous computations, and does not extend eastward because of a bathymetric effect. Its extension is blocked by the constriction of the isobaths due to the delta structure of the Grand Rhône (Fig. 2b).

The eastern gyre, which is very comparable to the (S0) case, is again restrained by a bathymetric effect. This is the reason why, in the absence of the river discharge, the secondary flow appearing in the Gulf of Fos extends in front of the river mouth. However (Fig. 15d), we observe that this eastern gyre is displaced westward because of the LPC influence and, as a consequence, the central gyre is swept southward.

The dipole secondary structure only occurs in the realistic (Sr) situation, but seems to be induced by the interaction of two main structures occurring in any case: i.e. a western and an eastern anticyclonic gyre, induced each side of the river mouth by the non-uniform Mistral offshore wind stress. In the absence of a river discharge – (S0) and (S2) cases – a central cyclonic gyre occurs that degenerates into either four cyclonic secondary vortices without LPC forcing (S1), or only two contrarotary vortices (the dipole) in the realistic case (Sr) because of the LPC confinement of the two main anticyclonic western and eastern gyre.

To better analyze the non-linear influence of each factor, we will proceed in the next section to a factor-separation study.

## 4.2 Factor-separation study

This simple but consistent method leads to a better understanding of the influence of each factor on the
<figure><img src="image_21.png" /><figcaption>**Fig. 17** Vertically integrated horizontal circulation in the nested fine grid induced by wind-tress but without the river discharges (S2) (bounded by 200-m depth isobath)</figcaption></figure>

vertically integrated horizontal circulation and their non-linear interaction. The contribution matrix of each factor, called (f0), (f1), (f2) and (f1-2) are obtained from the matrices (S0), (S1), (S2) and (Sr) in the following linearized way:

f0 = S0 (the single local upwelling flow),

f1 = S1 - S0 (influence of the river discharge),

f2 = S2 - S0 (influence of the LPC and the large-scale background circulation),

f1 - 2 = Sr - (S1 + S2) + S0 (influence of the non-linear interaction between the river discharge and the large-scale circulation).

It is easy to verify that the sum of the contribution matrix (f0, f1, f2 and f1-2) is equal to the realistic simulation Sr.

These results are presented in Fig. 18 concerning a zoom of the small computational domain in the vicinity of the Rhône river mouth (cf. Fig. 15).

Several processes can be identified from this technique : (1) Fig. 18b confirms the Grand Rhône river plume influence on the anticyclonic eastern circulation as far as the Gulf of Fos. (2) The influence of the LPC and the large-scale background circulation (Fig. 18c) is weak around the river mouth but induces an anticyclonic circulation in the central region that will contribute to considerably reduce the cyclonic central gyre (S0).

These two factors are necessary to destabilize the local 3-D upwelling for the dipole formation. Indeed, the necessary momentum is provided by the upwelling effect (Fig. 18a), when the non-linear interaction between the river plume and the mesoscale circulation clearly induces the genesis of the dipole (Fig. 18d). Furthermore, the river plume not only participates in the dipole generation, but also tends to block the upwelling circulation along the western dynamical front (Fig. 18b).

## 5 Conclusion

Coastal secondary structures have been observed on each side of the Rhône River plume under an upwelling favourable Mistral offshore wind blowing to the southeast. Nested 3-D baroclinic computations have been run on this site introducing the main forcings: the along-shelf Liguro Provençal Current, a realistic river discharge and inhomogeneous wind stress from ALADIN meteorological model, contemporary to a field campaign. The river plume structure and associated secondary structures have been reproduced in good agreement with a SPOT TSM image and in situ measured CTD vertical profiles. In particular, two clearwater patches appearing (in dark) in the SPOT image are shown to be associated with two vortices forming a dipole (one cyclonic gyre and the
<figure><img src="image_22.png" /><figcaption>Fig. 18a–d Influence of the factors a f0, b f1, c f2 and d f1–2</figcaption></figure>

other anticyclonic) along the western river plume front. The centres of the computed vortices and the clearwater patches observed on the SPOT image coincide. The computed structures are larger than the observed ones because of an entrainment effect of the TSM towards the centre of the vortices. Another single cyclonic gyre appears close to the Grand Rhône river mouth in the Gulf of Fos.

A process-oriented study has been performed to better understand the formation of these secondary structures each side of the Grand Rhône river plume. Moreover, a factor-separation technique enabled us to describe separately the influence of each of the two main

forcings: the river discharge and the mesoscale back-
ground shelf circulation when keeping the realistic
bathymetry as well as the temporal sampling of the
horizontally sheared wind which introduces a local
vorticity.

Thus, it can be stated from this analysis that the dipole is generated by the inhomogenous wind stress inducing 3-D upwelling confined by the western river front. However, the non-linear interaction between the mesoscale background circulation and the river plume is shown to be a key mechanism in the perturbation of the 3-D upwelling to the dipole occurrence.

The cyclonic gyre observed in the Gulf of Fos is generated simply by the wind stress but pushed into the Gulf by the river plume.
**Acknowledgements** This research was supported by the French Programme National Ocean Atmosphere Multi Echelles. BIO-DYPAR 3 CTD data were made available through the French Programme National d'Environnement Côtier by courtesy of Philippe Forget, Joel Gaggelli and Jean Jacques Naudin. Referees are acknowledged for constructive detailed comments on the first version of the paper. John Simpson is gratefully thanked for valuable comments.

## References

Arakawa A, Lamb VR (1997) Computational design of the basic dynamical processes of the UCLA general circulation model. Meth Comput Phys 17: 173-265

Arnoux-Chiavassa S (1998) Modélisation d'écoulements côtiers stratifiés présentant des fronts: application au panache du Rhône. Thèse de doctorat de l'université de Toulon et du Var

Arnoux-Chiavassa S, Rey V, Fraunié P (2003) Modeling 3D Rhône river plume using a higher-order advection scheme. Ocean Acta 26: 299-309

Beckers JM (1995) La méditerranée occidentale: de la modélisation mathématique à la simulation numérique. Thèse de doctorat de l'université de Liège

Berrabaa S, Fraunié P, Crochet M (2000) 2D Large-eddy simulation of highly stratified flows: the stepwise structure effect. In: Minev, Lin, Wong (eds). Advances in computation: theory and practice. Nova Science Books and Journals

Blumberg AF, Mellor GL (1987) A description of a three-dimensional coastal ocean circulation model. In: Heaps NS (ed). Three-dimensional coastal ocean models. Coastal and estuarines sciences, vol 4. American Geophysical Union, Washington, DC, pp 1-16

Chao S-Y, Boicourt WC (1986) Onset of estuarines plumes. J Phys Océanogr 16: 2137-2149

Chao S-Y (1988a) River-forced estuarines plumes. J Phys Océanogr 18: 72-88

Chao SY (1988b) Wind-driven motion of estuarine plumes. J Phys Océanogr 18: 1144-1166

Chao S-Y (1990) Tidal modulation of estuarine plumes. J Phys Océanogr 20: 1115-1123

Demarcq H, Wald L (1984) La dynamique superficielle du panache du Rhône d'après l'imagerie infrarouge satellitaire. Ocean Acta 7: 159-162

Dimitrov G (1999) QUICK scheme. Rapport LSEET

Durand N, Fiandrino A, Fraunié P, Ouillon S, Forget P, Naudin JJ (2002) Suspended matter dispersion in the Ebro ROFI: an integrated approach. Continental Shelf Res 22: Special issue FANS (fluxes across a narrow shelf) 267-284

Echevin V, Mortier L, Crépon M Interaction of a coastal current with a gulf: application to the shelf circulation of the Gulf of Lions in the Mediterranean Sea. J Phys Oceanogr (In press)

Estournel C, Kondrachoff V, Marsaleix P, Vehil R (1997) The plume of the Rhône: numerical simulation and remote sensing. Continental Shelf Res 17: 899-924

Estournel C, Broche P, Marsaleix P, Devenon JL, Auclair F, Vehil R (2001) The Rhône river plume in unsteady conditions: numerical and experimental results. Estuar Coast Shelf Sci 53: 25-38

Estournel C, Durrieu de Madron X, Marsaleix P, Auclair F, Juliand C, Vehil R (2003) Observation and modelisation of the winter coastal oceanic circulation in the Gulf of Lions under wind conditions influenced by the continental orography (FETCH experiment). J Geophys Res

Forget P, Devenon JL, De Maistre JC, Broche P, Leveau M (1990) VHF remote sensing for mapping river plume circulation. Geophys Res Lett 17: 1097-1100

Galperin B, Kantha LH, Hassid S, Rosati A (1988) A quasi-equilibrium turbulent energy model for geophysical flows. J Atmos Sci 45: 55-62

Galperin B, Rosati A, Kantha LH, Mellor GL (1989) Modeling rotating stratified flows with application to oceanic mixed layers. J Phys Oceanogr 19: 901-916

Garvine RW (1977) Observations of the motion field of the Connecticut River plume. J Geophys Res 83: 441-454

Garvine RW (1982) A steady-state model for buoyant surface plume hydrodynamics in coastal waters. Tellus 34: 293-306

Garvine RW (1987) Estuary plumes and front in shelf waters: a layer model. J Phys Oceanogr 17: 1877-1896

Garvine RW (1995) A dynamical system for classifying buoyant coastal discharges. Continental Shelf Res 15: 1585-1596

Garvine RW, Monk JD (1974) Frontal structure of a river plume. J Geophys Res 79: 2251-2259

Gaskell PH, Lau AKC (1988) Curvature-compensated convective transport: SMART, a new boundedness-preserving transport algorithm. Int J Num Meth Fluids 8: 617-641

Kourafalou VH, Lee TN, Oey L-Y, Wang JD (1996) The fate of river discharge on the continental shelf, 2. Transport of coastal low-salinity waters under realistic wind and tidal forcing. J Geophys Res 101: 3435-3455

Leonard BP (1991) The ULTIMATE conservative difference scheme applied to unsteady one-dimensional advection. Comp Meth Appl Mech Eng 88: 17-74

Leredde Y, Devenon JL, Dekeyser I (1999) Turbulent viscosity optimized by data assimilation. Annal Geophys 17: 1463-1477

Linden LF (1980) Mixing across a density interface produced by grid turbulence. J Fluid Mech 100: 691-703

Luyten PJ, Deleersnijder EL, Ozer J, Ruddick KG (1996) Presentation of a family of turbulence closure models for stratified shallow waters flow and preliminary application to the Rhine outflow region. Continental Shelf Res 16: 101-130

Marsaleix P, Estournel C, Kondrachoff V, Vehil R (1998) A numerical study of the formation of the Rhône river plume. J Mar Syst 14: 99-115

Mellor GL, Yamada T (1982) Development of a turbulence closure model for geophysical fluid problems. Rev Geol Space Phys 20: 851-875

Millot C (1990) The Gulf of Lions' hydrodynamics. Continental Shelf Res 10: 885-894

Oey L-Y, Mellor GL (1993) Subtidal variability of estuarine outflow plume and coastal current: a model study. J Phys Oceanogr 23: 164-171

Orlanski I (1976) A simple boundary condition for unbounded hyperbolic flows. J Comput Phys 21: 255-261

Pietrzak J, Jacobson JB, Burchard H, Vested HJ, Petersen O (2002) A three-dimensional hydrostatic model for coastal and ocean modelling using a generalised topography following co-ordinate system. Ocean Modelling 4: 173-205

Posmentier ES (1977) The generation of salinity fine structure by vertical diffusion. J Phys Oceanogr 7: 298-300

Phillips OM (1977) The dynamics of the upper ocean, 2nd ed. Cambridge University Press, Princeton

Rodi W (1984) Turbulence models and their application in hydraulics, 2nd ed. Intern Assoc Hydraul Res Delft, Netherlands

Ruddick KG, Deleersnijder E, De Mulder T, Luyten PJ (1994) A model study of the Rhine discharge front and downwelling circulation. Tellus 46(A): 149-159

Ruddick KG, Deleersnijder E, Luyten PJ, Ozer J (1995) Haline stratification in the Rhine-Meuse freshwater plume: a three-dimensional model sensitivity analysis. Continental Shelf Res 15: 1597-1630

Simpson JH (1997) Physical processes in the ROFI regime. J Mar Syst 12: 3-15

Turner RE, Rabalais NN (1994) Coastal eutrophisation near the Mississipi River delta. Nature 368: 619-621

Xing J, Davies AM (1999) The effect of wind direction and mixing upon the spreading of a buoyant plume in a non-tidal regime. Continental Shelf Res 19: 1437-1483

Xing J, Davies AM (2002) Influence of topographic features and along shelf flow upon the Ebro plume. Continental Shelf Res 22: 199-227