# Dispersal of owenia fusiformis larvae by wind-driven currents: turbulence, swimming behaviour and mortality in a three-dimensional stochastic model

Katell Guizien, Timothée Brochier, Jean-Claude Duchêne, Byoung Seol Koh, Patrick Marsaleix

► To cite this version:

Katell Guizien, Timothée Brochier, Jean-Claude Duchêne, Byoung Seol Koh, Patrick Marsaleix. Dispersal of owenia fusiformis larvae by wind-driven currents: turbulence, swimming behaviour and mortality in a three-dimensional stochastic model. Marine Ecology Progress Series, 2006, 311, pp.47-66. 10.3354/meps311047. hal-00138672

HAL Id: hal-00138672

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

Submitted on 16 Dec 2020

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
# Dispersal of *Owenia fusiformis* larvae by wind-driven currents: turbulence, swimming behaviour and mortality in a three-dimensional stochastic model

K. Guizien<fnref n="1" />*, T. Brochier<fnref n="1" />, J.-C. Duchêne<fnref n="1" />, B.-S. Koh<fnref n="1" />,<fnref n="3" />, P. Marsaleix<fnref n="2" />

**ABSTRACT:** A Lagrangian stochastic larva tracking model based on a 3-dimensional (3D) high resolution wind-driven coastal circulation model is used to study the dispersal of benthic larvae. The larva tracking model includes 3D advection and turbulence, and a species-specific larval swimming behaviour that accounts for ontogenic changes, sensitivity to light exposure and inter-individual variability. Larval mortality can also be included. The dispersal model is applied to *Owenia fusiformis* larvae, whose swimming behaviour description is based on both existing data and new complementary measurements. Larval velocities (resulting from both settling and swimming behaviour) were measured with actographic equipment and ranged between –1 and 0.9 mm s⁻<fnref n="1" />. Measured swimming activity rates were lower than 50 %. The sensitivity study of larval dispersal in March–April 1999 showed that: (1) the dispersal of neutrally buoyant passive larvae is more sensitive to the physical forcings resolution, because of both advection and diffusion processes, than to the variability of spawning locations within neighbouring grid cells (up to 1 km apart) in Banyuls Bay (France, NW Mediterranean); and (2) a physical barrier, located at 20 m deep in Argelès (France, NW Mediterranean) and 30 m deep in Banyuls Bay, separated nearshore and offshore larval dispersal in 1999. The final positions and local retention of larvae released in Banyuls Bay and Argelès result from: (1) the balance between the 3D turbulence, larval settling velocity (~0.8 mm s⁻<fnref n="1" />) and swimming activity rate; and (2) natural mortality, although the effect is not proportional to survival rates. High resolution larvae dispersal patterns for *O. fusiformis* in Banyuls Bay suggest that self-recruitment was low in the Banyuls population during spring 1999 and confirm that post-settlement deposit patterns observed there in May 1999 were insignificant. In addition, interconnections between the Argelès and Banyuls populations can exist.

**KEY WORDS:** Larval dispersal · 3D numerical modeling · Lagrangian/Eulerian models · *Owenia fusiformis* · Wind-driven currents

Resale or republication not permitted without written consent of the publisher

## INTRODUCTION

Most benthic invertebrates (70%, Thorson 1946) show a pelagic larval stage in their life cycle that can last from a few hours to several weeks. This observation raises a recurrent question on which many works (e.g. Crisp 1976) are based: Is a pelagic larval stage the most effective development mode for sustainable population dynamics in benthic species? Larvae released into the flow at any given location may either benefit

from currents circulating inside a retention zone, or increase their coverage at the risk of entering inadequate habitats. There is no definite answer to what the most effective development mode for recruitment success and sustainable population dynamics is (Scheltema 1986). This led Scheltema (1986) to suggest posing a different question: How do benthic animals cope with such a common development mode to ensure sustainable population dynamics? And consequently, when dispersion by currents is coupled with larval

<footnote n="1">¹Laboratoire d'Océanographie Biologique de Banyuls (CNRS/UPMC), BP 44, 66651 Banyuls-sur-Mer, France</footnote>

<footnote n="2">²Laboratoire d'Aérologie de Toulouse (CNRS/UPS), 14 avenue E. Belin, 31400 Toulouse, France</footnote>

<footnote n="3">³National Fisheries Research and Development Institute, Tidal-flat Research Center, 573-882 Gunsan, Republic of Korea</footnote>
behaviour, what population dynamics emerge? Benthic population dynamics is a major issue in the management of nearshore regions. Some benthic species are fished directly and others are used to assess the biochemical quality of the environment since their sedentary adult stage integrates bio-chemical environmental changes throughout their lifetime. The use of some of these populations as bio-indicators is questionable if their dynamics are unstable due to larval dispersal.

Pioneering larval dispersal studies have tried to identify physical retention structures (which are generally highly predictable) relative to the duration of the animal's larval stage, e.g. tidal estuaries for a short larval stage (de Wolf 1974, Chen et al. 1997) or a general thermo-haline circulation loop for a longer larval stage (Hill et al. 1996). But, in most cases, a given species' recruitment patterns depend on the competition between dispersion by highly changeable circulation (wind-driven) and recruitment failure in hostile habitats. In this context, it is no longer possible to rely on the only correlation between an observed physical structure and a species' distribution area. It is necessary to simultaneously integrate the locomotory behaviour of larvae from the time and location of their release until their settlement with motion caused by the flow. To achieve this goal, numerical modelling of larval dispersal is a powerful tool. Modelling larval dispersal requires at a minimum solving for the advection-diffusion of a pool of larvae as driven by the current and turbulence in an appropriate physical model in addition to a species-specific biological model describing the evolving larval stage.

Larval dispersal models in the literature can be divided into 2 groups: (1) Eulerian models solved at each mesh point for the time evolution of a concentration representing the cohort, including advection, turbulent diffusion and sometimes certain population dynamics parameters (Carlotti 1996, Ellien et al. 2000, 2004, Barnay et al. 2003); and (2) Lagrangian models solved for individual pathways, accounting for advection (Johnson & Hess 1990, Werner et al. 2001), and sometimes including either specific larval locomotory behaviour (Verdier-Bonnet et al. 1997, Perdersen et al. 2003), turbulence dispersion (Oliver et al. 1992, Ommundsen 2002) or in one case, mortality (Porch 1998). To our knowledge, few studies have accounted for 3D advection (Verdier-Bonnet et al. 1997, Perdersen et al. 2003) and only one included the effect of vertical turbulent diffusion on dispersal (Werner et al. 2001). Ideally, a comprehensive larval dispersal model should account for 3D advection and diffusion, individual larval locomotory behaviour and some population dynamics parameters (spawning and mortality rates). It should be noted that accounting for individual locomotory behaviour is straightforward in a Lagrangian approach whereas these effects have not yet been included in any Eulerian studies.

In this paper, we present a 3D Lagrangian stochastic biophysical model for larval dispersal which accounts for turbulence, individual larval swimming behaviour, and mortality. It is coupled to a 3D wind-driven circulation model. The swimming behaviour for *Owenia fusiformis* larvae was studied experimentally and the results were used to derive a specific swimming behaviour model. The 3D Lagrangian larval dispersal model was then combined with the *O. fusiformis* larval swimming behaviour model to study *O. fusiformis* larval dispersal in Banyuls Bay (France, NW Mediterranean) during spring 1999. The larval dispersal model's sensitivity to both physical and biological processes is evaluated, and larval dispersal patterns are compared to available post-settlement observations in May 1999 in Banyuls Bay.

## MATERIALS AND METHODS

**The Lagrangian stochastic biophysical model.** In the following section, the Lagrangian stochastic biophysical model for larval dispersal is presented. It includes (1) 3D physical forcings in advection and diffusion, (2) individual stochastic larval velocity, and (3) a stochastic mortality procedure. In particular, the species-specific stochastic larval velocity includes an age-dependent settling velocity following ontogenic changes in larvae for their benthic metamorphosis, and a swimming ability description including individual variability in light exposure and age dependencies.

**Accounting for 3D turbulence:** Although the Lagrangian approach is widespread amongst larval dispersal studies (Johnson & Hess 1990, Chen et al. 1997, Verdier-Bonnet et al. 1997, Werner et al. 2001, Perdersen et al. 2003), few of these studies accounted for turbulence and, more often than not, only in 2 horizontal dimensions (Oliver et al. 1992, Porch 1998, Ommundsen 2002). Since the pioneering work of Oliver et al. (1992) on the dispersal of neutrally buoyant coral larvae, small-scale turbulence has been represented by a 2D random walk procedure (Monte-Carlo Method) based on an analogy to molecular diffusion. The turbulent horizontal diffusion coefficient was assumed to be constant, except in Porsch's study (1998).

However, turbulent dispersion can be accounted for in Lagrangian particle-tracking models by other procedures. In particular, relying on the concept of particle eddy encounters, turbulence is straightforwardly implemented by adding a 3D fluctuating velocity to the mean flow velocity. The reliability of this procedure has been shown, for example, in a turbulent atmo-
spherical boundary layer (Ohba et al. 1997) or in a tur-
bulent plane mixing layer (Oliveira et al. 2002). Our
paper uses such a Lagrangian stochastic modelling
procedure for turbulence. The Lagrangian stochastic
larval dispersal model reads:

$$
\frac{d\bar{X}}{dt} = \bar{U}(\bar{X},t) + \bar{u}'(\bar{X},t) + \bar{u}_{sp}(\bar{X},t) \quad (1)
$$

in which $\bar{X}$ is the 3D space location and $t$ is time,
$\bar{U}(\bar{X},t)$ is the flow velocity at $(\bar{X},t)$, $\bar{u}'(\bar{X},t)$ is the fluc-
tuating turbulent velocity at $(\bar{X},t)$, and $\bar{u}_{sp}(\bar{X},t)$ is the
individual larva's velocity, which will be defined in
the next section.

Each of the 3 components of the turbulent velocity $\bar{u}'(\bar{X},t)$ is obtained by randomly sampling a Gaussian distribution with standard deviation $\sqrt{2k/3}$, where $k$ is the Turbulent Kinetic Energy (TKE), and therefore the turbulence is assumed to be isotropic. In practice, in our model, the Box-Muller algorithm simultaneously generates 4 Gaussian random variables at each time step. Three of them are used to produce the components of the turbulent velocity, and the fourth remains unused. A 3D Eulerian wind-driven current model provides 3D time- and space-varying velocities and the TKE forcings.

Eq. (1) is integrated over time through a Euler
scheme to produce the successive positions of each
Lagrangian particle. The integration time-step for
Eq. (1) is constrained by 2 factors. First, it depends on
the ratio of the particle to fluid response times τp/τ,
namely the Stokes number. In our case, it can be esti-
mated that τ = Δx/u ~ 500 s with a horizontal flow
velocity u ~ 0.2 m s⁻¹, a horizontal mesh size Δx ~ 100 m
and τp < 1 s depending on the larva size. When
τp/τ << 1 the particles essentially follow the fluid and
according to Oliveira et al. (2002) the integration time-
step should be around 0.1τ, i.e. 50 s. Second, fluctuat-
ing velocities should be updated every time a larva
encounters a new eddy. This particle-eddy interaction
time can be defined as the minimum between the eddy
life time and the eddy transit time through a cell. We
estimate that the particle-eddy interaction time is the
transit time of fast moving surface waves through a cell.
In practice, the latter is close to the barotropic mode
time step of the 3D current model, i.e. roughly ~1 s.

Theoretically, the integration time-step for Eq. (1)
should be lower than this particle-eddy interaction
time. However, since the circulation flow velocity and
the larva's own velocity vary slowly, in practice, Eq. (1)
is integrated over the 3D current model's baroclinic
time step of ~30 s, which fulfils the first limitations.
However, the fluctuating turbulent velocities used in
Eq. (1) result from averaging the fluctuating velocities
(sampled at the 3D current model barotropic time step
of ~1 s) over the baroclinic time step.

Preliminary tests estimated that about 30 particles
released every 10 min are sufficient to simulate the
main turbulent dispersion features. Lagrangian advec-
tion-diffusion simulations were evaluated and com-
pared to Eulerian ones based on the turbulent dif-
fusivity concept (Appendix 1, available at: www.
int-res.com/articles/suppl/m311p047_app.pdf).

**Accounting for swimming behaviour:** In Eq. (1), a specific larval velocity $\bar{u}_{sp}(\bar{X},t)$ has to be specified. For passive particles, this velocity equals the settling velocity. For growing larvae, the specific velocity may vary with time and space and between individuals, as a result of both ontogenic settling changes and swimming activity.

In order to account for inter-individual variability, we
suggest that the larval velocity should be generated
stochastically from a velocity probability distribution
given by its cumulative frequency distribution. This cu-
mulative velocity frequency distribution $F_{S,L}$ may vary
with ontogenic changes defined by the larva stage $S =
t - t_{spawn}$ where $t_{spawn}$ is the spawn date and light condi-
tions are defined by the light quantity $L(z,t)$ at $z$ meters
below the sea surface and at time $t$ in days. Thus, $L(z,t) =
L_{max}\sin[2\pi(t - 0.25)]e^{-\alpha z}$ where $L_{max}$ is the maximum
light quantity at the sea surface at solar midday and $\alpha$ is
the light attenuation coefficient of sea water. Ideally,
the cumulative velocity frequency distribution $F_{S,L}$ for a
given species is derived from experimental studies.

Given a cumulative larval velocity frequency distrib-
ution $F_{S,L}$, a random larva velocity $\bar{u}_{sp}(\bar{X},t)$ is generated
independently:

$$
\bar{u}_{sp}(\bar{X}, t) = F_{S,L}^{-1}(\xi) \qquad (2)
$$

where $\xi$ is an equiprobable random variable giving
the cumulated frequency corresponding to the larval
velocity $\bar{u}_{sp}(\bar{X},t)$ as schematised in Fig. 1. This stochas-

<figure><img src="image_2.png" /><figcaption>Fig. 1. *Owenia fusiformis*. Diagram for stochastic generation of larvae velocity from a cumulative velocity frequency distribution *F*. Larvae velocity is *F*⁻¹ (ξ) with ξ an equiprobable random variable; *V*<sub>min</sub>: minimum velocity; *V*<sub>max</sub>: maximum velocity</figcaption></figure>
tic generator is used for each larva at each integration
time step depending on the larval behaviour.

**Accounting for mortality:** Natural mortality and/or predation is accounted for in the Lagrangian biophysical model. At each integration time step, the number of larvae to be eliminated is computed by applying a daily mortality rate *M* to the live larvae from the previous time step. In general, this is a real number. The value is rounded towards zero to get an integer number of larvae to be randomly eliminated at this time step. The remaining fractional part is added to the number of larvae to be eliminated at the following time step. The survival rate after *J* days equals (1 − *M*)^J.

**Specific swimming behaviour study for *Owenia fusiformis* larvae.** The tubicolous polychaete *O. fusiformis* Delle Chiaje is a polythelic species with a ben-tho-pelagic life cycle. It has a mitraria larva, character-istic of the Oweniidae family with a ciliated umbrella and long straight bristles (Fig. 2a). Four ontogenic stages can be distinguished (Wilson 1932). Koh (2002) estimated stage durations for the NW Mediterranean Sea (Fig. 2b) with the entire planktonic stage lasting 22 to 23 d.

Although Wilson (1932) already mentioned a moder-
ate upward swimming behaviour, previous modelling
studies (Thiébaut et al. 1992, Verdier-Bonnet et al.
1997, Barnay et al. 2003) neglected any active swim-
ming. These studies focused on estimating the larva's
settling velocity and on its ontogenic variations.
Thiébaut (1994) derived an age-dependent settling
velocity (also shown in Fig. 2b) from field observations
of ontogenic vertical migration. He suggested an
upward velocity of 2.5 mm s⁻¹ for eggs (1 d), a zero
velocity for mitraria larvae younger than 20 d, and a
downward velocity of 1.13 × 10⁻² mm s⁻¹ for mitraria
larvae older than 20 d, valid until the larvae settled
after 28 d with a fall velocity of 2.5 mm s⁻¹. Similarly,
Verdier-Bonnet et al. (1997) considered that the *Owe-
nia fusiformis* larva is a true passive tracer with a zero
velocity for 17 d, but their study also suggested a
higher settling velocity of 0.76 mm s⁻¹ based on mea-
surements for the oldest stages (also shown in Fig. 2).
Recent direct measurements on larvae showed no
ontogenic evolution of the settling velocity: about
0.8 mm s⁻¹ for all stages (Koh 2002). Rey (1995) mea-
sured swimming velocities for larvae in Stages 1 and 2

<figure><img src="image_3.png" /><figcaption>Fig. 2. *Owenia fusiformis*. (a) Larvae evolution stage durations defined for the Mediterranean Sea (Koh 2002) and (b) settling velocity according to various authors: Thiébaut 1994 (— —), Verdier-Bonnet et al. 1997 (···) and Koh 2002 (—)</figcaption></figure>
and reported upward velocities reached 0.7 ± 0.37 mm s⁻¹ (13 ind.). In our study, new measurements were carried out for larvae in Stages 2 to 4 to complement the description of *O. fusiformis* larval velocities $\bar{u}_{sp}(\bar{X}, t)$ following its ontogenic evolution.

Swimming activity measurements were made using an actographic system that can detect planktonic organisms as small as 100 µm (Duchêne & Queiroga 2001). Twenty-five live *Owenia fusiformis* larvae were extracted from fresh plankton net samples on 23 and 26 April 2004. The larvae were kept in filtered seawater at 12°C in darkness until the experiment began and were fed with *Nitzchia* sp. algae 1 d after collection. On 23 April 2004, most larvae were at Stage 2 and evolved towards Stage 3 over the following days.

Groups of 3 larvae at the same ontogenic stage were transferred in quartz glass cells (10 × 10 × 45 mm) for the measurements and left undisturbed. Larvae were mechanically stimulated just before the measurements by slowly turning the bubble-free box upside down for 1 min. This procedure inverted gravitational forces while minimizing the creation of any circulation inside the box. Larvae were then exposed alternately to Infra Red (IR) light peaking at 830 nm (<1.5 × 10⁻³ W m⁻²) and to Cold White Light (CWL, 2.992 W m⁻²) during a series of 10 min long successive sequences. This procedure was designed to study the sensitivity to light exposure as experienced at sea during nycthemeral and depth variations while minimizing convection due to the slight warming under IR light. According to field measurements, 2.992 W m⁻² is the light flux to Banyuls Bay at around 15 m depth. A diagram of the experimental set-up is shown in Fig. 3. Projections of larval positions on a vertical plane were recorded every second. The spatial resolution for position-tracking ranged from 64 to 76 µm vertically and from 78 to 93 µm horizontally, depending on the distance between the box and the camera.

This procedure allowed the swimming behaviour of *Owenia fusiformis* larvae to be measured with respect to both their ontogenic stages (Stages 2 to 4) and light exposure (nycthemeral and depth variations) over the duration given in Table 1. In all, 4 h of trajectory were recorded. An average swimming activity rate, defined

as the ratio of non-zero velocity duration to total record duration, was computed for each stage, light conditions and time lag after mechanical stimulation. Because of the limited experimental box size, motion activity is considered as equivalent to swimming activity. Within active motion periods, vertical and horizontal velocity components were derived from larval trajectories and filtered out by eliminating unrealistic acceleration (>1 mm s⁻²) and applying a 4 s moving average.

The general biophysical model combined with the observations of *Owenia fusiformis* swimming behaviour were used to describe *O. fusiformis* larval dispersal in a Mediterranean embayment off the Catalan coast (Banyuls Bay, France).

*Owenia fusiformis* larval dispersal in Banyuls Bay in 1999. *Ecology of Owenia fusiformis along the French Catalan coast: O. fusiformis* has several interesting features for testing our Lagrangian biophysical larval dispersal model. Indeed, *O. fusiformis* has a short reproductive period in spring (maximum 3 mo), allowing 3D simulations of wind-driven coastal circulation to be carried out for the entire larval dispersal period.

<figure><img src="image_4.png" /><figcaption>Fig. 3. Diagram of the actographic experimental set-up. Dashed line in the light exposure protocol indicates that some experiments were prolonged with a third sequence of Cold White Light exposure</figcaption></figure>

<table><caption>Table 1. *Owenia fusiformis*. Actographic record duration and measured swimming activity rate with respect to larvae ontogenic stages and light conditions. IR: Infra Red; CWL: Cold White Light</caption><thead><tr><th rowspan="2"></th><th colspan="2">Stage 2</th><th colspan="2">Stage 3</th><th colspan="2">Stage 4</th></tr><tr><th>Record length (s)</th><th>Swimming activity rate (%)</th><th>Record length (s)</th><th>Swimming activity rate (%)</th><th>Record length (s)</th><th>Swimming activity rate (%)</th></tr></thead><tbody><tr><td>IR, during first 10 min</td><td>8400</td><td>25</td><td>9600</td><td>33</td><td>600</td><td>2</td></tr><tr><td>IR, after 10 min</td><td>6000</td><td>5</td><td>12000</td><td>28</td><td>1200</td><td>2</td></tr><tr><td>CWL, after 10 min</td><td>8400</td><td>7</td><td>12000</td><td>8</td><td>1200</td><td>8</td></tr></tbody></table>
<figure><img src="image_5.png" /><figcaption>Fig. 4. NW Mediterranean Sea bathymetric map with depth contour in meters (depth contour spacing is 200 m below 200 m). (♦) Locations of known *Owenia fusiformis* populations (M. Marseille; A, Argelès; B, Banyuls; C, Creus; R, Rosas; and E, Estarlit). (●) Other possible locations of *O. fusiformis* populations. Boxes show the extension of the 4 computational domains; the 2 largest boxes are labeled in the upper left-hand corner. 1: 1500 m resolution grid; 2: 700 m resolution grid</figcaption></figure>

In addition, discontinuous populations with restricted coverages have been clearly identified from the Italian Thyrrenian coast to the Spanish Catalan coast (Guille 1971, Massé 1971, Maron-Ramos 1973, Amouroux 1974, Gambi & Giangrande 1986, Pinedo et al. 2000). These population locations are plotted in Fig. 4. Other potential populations inferred from point observations are also indicated (J. M. Amouroux & C. Labrune pers. comm.).

Such a limited number of populations and consequently of spawning locations constitutes an advantage for our Lagrangian approach over the Eulerian alternative with respect to computing capacity requirements. Eulerian simulations require a fixed storage capacity proportional to the number of mesh points. In contrast, storage capacity for Lagrangian simulations is proportional to the number of released particles, which in turn depends on the number of spawning locations. Although this constraint is less and less important with the current progress in data storage technology, this explains why the Eulerian approach was preferred in larval dispersal studies when widespread populations were studied (Ellien et al. 2000, 2004, Barnay et al. 2003).

In this paper, we focus on the Banyuls Bay population which was studied in 1999 by Koh (2002). Two 0.1 m² grab samples were taken in May 1999 at each station in

the bay (shown by a circle in Fig. 5) and the resulting *Owenia fusiformis* density distribution is also shown. Two zones with densities ranging from 10 to 100 ind. m⁻² were observed. One is located in the middle of the northern cove of Banyuls Bay (at 20 to 30 m depth) and the other in the middle of the southern cove of Banyuls Bay (at 10 to 30 m depth). *O. fusiformis* worms found in these zones were small and not sexually mature. Reproductive adults were only observed in a very restricted area at 8 to 12 m depth in the southern cove of Banyuls Bay, where densities were up to 2000 ind. m⁻² (solid diamond marked 'B' in Fig. 5). Grab sample observations were complemented by daily vertical plankton hauls at Location B collected with a ring mouth net (WP2 international standard, 0.57 m mouth diameter, 2.6 m long, 200 µm mesh size) from 15 January to 31 May 1999. *O. fusiformis* larvae were observed from 18 March until 12 April 1999 (Koh 2002). From the duration of the planktonic stage, the spawning period in 1999 was estimated to last from

<figure><img src="image_6.png" /><figcaption>Fig. 5. *Owenia fusiformis*. Observations in spring 1999 in the 2 coves of Banyuls Bay: reproductive adult zone (♦) in Banyuls Bay (B) and *O. fusiformis* density (as ind. m⁻²) in May 1999 are shown. (●) Grab sample locations in May 1999. Depth contours are also plotted with 10 m spacing; 30 m deep contour is labeled</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>Fig. 6. *Owenia fusiformis*. Time series of wind speed (black lines) and wind direction (direction from which the wind blows, grey lines) measured at Cap Béar together with larvae presence period in Banyuls Bay (dashed line) for 1999. Lifetime of larvae released from successive spawnings at midday, once a day, are shown (solid horizontal lines). Dates given as d/mo/yr</figcaption></figure>

18 March 1999 until 21 March 1999 (shown in Fig. 6) and also suggests that benthic observations in May 1999 reflect *O. fusiformis* post-settlements patterns for the spring 1999 recruitment.

**Physical forcings along the French Catalan coast: 3D wind-driven current modeling:** Coastal circulation in the NW Mediterranean, as a result of strong meteorological forcings, has a high variability and requires reliance on numerical modelling to describe it. Wind-driven currents in Banyuls Bay were modelled using the SYMPHONIE model (Estournel et al. 1997). This model has been validated at various spatial scales (from a few km up to hundred of km) in the NW Mediterranean Sea with a series of current observations and hydrographic surveys over the years 1998 to 2001 (Estournel et al. 1997, Dufau-Julliand et al. 2004, Petrenko et al. 2005, Ulses et al. 2005). SYMPHONIE is a 3D Eulerian primitive equations and sigma-coordinate model with a free surface explicitly calculated by splitting the model equations into barotropic (rapid) and baroclinic (slow) modes (Blumberg & Mellor 1987). Vertical turbulence is estimated using a second-order closure scheme (Gaspar et al. 1990) which solves for the TKE equation with an algebraic formulation of the mixing and dissipation lengths, and assumes that the turbulence is isotropic. The horizontal viscosity for the momentum equation is not derived from this turbulence closure and is instead set to 10 m² s⁻¹, in order to reduce numerical instabilities. Bottom boundary conditions specify that the momentum flux is equal to the bottom stress, assuming the velocity follows a logarithmic law and turbulence equilibrium is reached. The bottom roughness length is constant over the domain

(10⁻² m). External forcings for the circulation model also include open boundary flow conditions and large river inputs (the Grand Rhône, Petit Rhône, Hérault, Orb and Aude rivers). A 1-way nesting procedure provides the open-boundary flow conditions for a given area from the larger area enclosing it.

Simulations of Banyuls Bay circulation were run from 10 March until 15 April 1999 over 4 nested grids of increasing resolution and decreasing size as shown in Fig. 4 (indicated by 4 boxes, the 2 largest are labeled 1 and 2 in the upper left corner). Details of the numerical treatment of open boundary conditions and the nesting procedure with grid sizes and resolutions are given in Appendix 2 (available at: www.int-res.com/articles/suppl/m311p047_app.pdf). The largest area was initialised and driven at its open boundaries by a winter climate regime for the

Liguro-Provençal current which was maintained constant during the entire simulation (Estournel et al. 1997). It is worth noting that in Banyuls Bay in late winter-early spring, the sea is in general well mixed with a nearly constant salinity of around 38 and a temperature of around 12°C. Surface boundary conditions specified the downward and upward radiative fluxes of long and short wave, and also the fluxes of heat and momentum (Dufau-Julliand et al. 2004).

In practice, surface boundary conditions were derived from bulk formulae (Geernaert 1990) based on meteorological parameters (wind, pressure, air temperature and humidity) and the sea surface temperature as computed by SYMPHONIE. Meteorological forcings were supplied by METEO FRANCE's weather forecast model ALADIN with a 0.1 degree resolution and interpolated on the SYMPHONIE grids. For the 2 finest grids (Grid 3 is 250 m × 250 m, Grid 4 is 100 m × 100 m, indicated on Fig. 7), Cap Béar's strong acceleration effect on local wind was included by mixing computed wind values and hourly wind measurements collected at METEO FRANCE's Cap Béar station (Fig. 6). The wind time series recorded at Cap Béar shows that a strong, although not regular, Tramontane wind dominated throughout the study period. A smooth transition between measured and computed values was applied, which varied with distance from the cape, over a 5 km radius area (Fig. 7). The hydrodynamic outputs required to compute larval dispersal (velocities and TKE) were stored hourly.

**Larval dispersal simulations:** Offline simulations of larval dispersal were carried out using interpolated hydrodynamic values between each hour from 18 March
<figure><img src="image_8.png" /><figcaption>Fig. 7. *Owenia fusiformis*. Example of forcing spatial variations in wind intensity after mixing Cap Béar measurements and computations by METEO FRANCE's weather forecast model ALADIN. Boxes show the extension of the computational domains (42°N Mercator projection) labeled in the upper left-hand corner: 3 is the 250 m resolution grid and 4 is the 100 m resolution grid. Depth contours ranging from 0 to 80 m are plotted with 10 m spacing. (●) Spawning locations for neutrally buoyant passive larvae at increasing depths and distance from the coast. (◆) Spawning locations in *O. fusiformis* populations: 2 in Argelès (A) and 2 in Banyuls Bay (B)</figcaption></figure>

1999 at 0:00 h until 15 April 1999 at 15:00 h. A sequence of larval dispersal simulations, progressively adding turbulence, settling, swimming and mortality processes (Table 2) was performed using physical forcings with a resolution of 250 m × 250 m and 100 m × 100 m. Swimming behaviours based on the previously mentioned studies and incorporating our new measurements were used. Daily mortality rates of 1

<table><caption>Table 2. *Owenia fusiformis*. Biophysical model forcings. Swimming behaviour is based on actographic measurements of swimming speed (1) with the measured swimming activity as plotted on Fig. 10: (2) with a 100% swimming activity rate for Stage 1 and as measured for other stages; or (3) with a 100% swimming activity rate for all stages. T: test</caption><thead><tr><th>Condition name</th><th>Advection</th><th>Diffusion</th><th>Settling velocity (mm s<sup>-1</sup>)</th><th>Swimming behaviour</th><th>Daily mortality rate M (%)</th></tr></thead><tbody><tr><td>T1</td><td>Yes</td><td>No</td><td>0</td><td>No</td><td>No</td></tr><tr><td>T2</td><td>Yes</td><td>Yes</td><td>0</td><td>No</td><td>No</td></tr><tr><td>T3</td><td>Yes</td><td>Yes</td><td>0.8</td><td>No</td><td>No</td></tr><tr><td>T4</td><td>Yes</td><td>Yes</td><td>0.8</td><td>1</td><td>No</td></tr><tr><td>T5</td><td>Yes</td><td>Yes</td><td>0.8</td><td>2</td><td>No</td></tr><tr><td>T6</td><td>Yes</td><td>Yes</td><td>0.8</td><td>3</td><td>No</td></tr><tr><td>T7</td><td>Yes</td><td>Yes</td><td>0.8</td><td>1</td><td>1</td></tr><tr><td>T8</td><td>Yes</td><td>Yes</td><td>0.8</td><td>1</td><td>10</td></tr></tbody></table>

and 10% were tested. These values were chosen within the range estimated for marine invertebrate larvae (i.e. 1.6 and 35.7% d⁻¹, Rumrill 1990), in the absence of any specific mortality measurements concerning *Owenia fusiformis* larvae.

Continuous spawning was simulated by releasing 30 eggs every 30 min from 18 March 1999 at 0:00 h until 21 March 1999 at 23:30 h (5760 eggs per spawning location). Two *Owenia fusiformis* populations are included in the simulation domain, one in Argelès (A) and the other in Banyuls Bay (B). Two neighbouring spawning locations (Fig. 7) were considered for each population. In Banyuls Bay, the 2 spawning locations were only 500 m apart, reflecting the limited area where reproductive adults have been found (Koh 2002). The Argelès population coverage is larger than in Banyuls Bay (J. M. Amouroux pers. comm.) and thus the 2 spawning locations were set 1150 m apart. These effective spawning locations for *O. fusiformis* larvae are complemented by 6 additional spawning locations along the sites' bathymetric gradients, 3 in southern Banyuls Bay and 3 in Argelès (Fig. 7). Larvae positions were stored every 10 min during the simulations.

## RESULTS

### Swimming behaviour measurements and model for *Owenia fusiformis* larvae

#### Swimming activity rates measurements

Swimming activity rates recorded during the first 10 min (under IR conditions) and after 10 min exposure, depending on the light conditions, are reported in Table 1 for each larval stage.

For Stage 2 larvae, the swimming activity rate during the first 10 min is around 5 times higher than the swimming activity rate after 10 min, whatever the light exposure.

Firstly, the absence of motion recorded after 10 min for Stage 2 larvae confirms that convection due to heating by the light was limited, even under IR light. Secondly, it shows that Stage 2 larvae swimming activity does not depend on light exposure. Thirdly, we suggest the difference between the first 10 min and later observations may be explained by the absence of mechanical stimulation and is considered as an experimental bias since at sea, the water mass is never at rest. In summary, we consider that Stage 2 larvae should show a swimming activity rate of 25% on average, without sensitivity to light exposure.
In contrast, Stage 3 larvae show similar swimming activity rates under IR conditions in the first 10 min and in later measurements, but their swimming activity is much lower under CWL exposure. We have already shown that convection due to light heating was limited. Thus, the succession of decreases and increases in swimming activity under the alternating light conditions indicates that Stage 3 larvae have a photo-negative reaction.

Swimming activity rate for Stage 4 larvae was always very low and corresponds to free fall conditions. Overall, these measurements show a decrease in larval swimming activity amongst the ontogenic Stages 2 to 4.

## Swimming velocity measurements

In Fig. 8, a sample raw trajectory record is shown (Fig. 8a) together with the vertical and horizontal velocity components derived from it (Fig. 8b). Although it is clear that the directions of larval velocities are constrained by the quartz glass cell size, velocities magnitudes are not. Indeed, the quartz glass cell size may limit the magnitude of horizontal and vertical velocities to 2.5 and 10 mm s⁻¹ respectively, over 4 s, but these values were never reached. Hence, larval motion is treated as random behaviour. In view of the comparable swimming activity rates under different light exposure conditions, all Stage 2 velocity records were merged together. In contrast, the Stage 3 larvae velocity records were grouped by light exposure conditions.

For Stage 4 larvae, trajectory records always indicated a falling motion; hence, the swimming behaviour of Stage 4 larvae was neglected. For all stages, only vertical velocities were accounted for, because at sea, horizontal advection is 1 order of magnitude larger than vertical advection.

Empirical distributions for vertical velocities were constructed to describe the active motion of larvae depending on both their stage and light exposure by cumulating the time spent in 25 velocity classes ranging from -1.15 to 1.25 mm s⁻¹ (Fig. 9). The average

<figure><img src="image_9.png" /><figcaption>Fig. 8. *Owenia fusiformis*. (a) Example of a raw trajectory record and (b) horizontal and vertical velocities for a Stage 3 larva under infra red light</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Fig. 9. *Owenia fusiformis*. Empirical velocity distributions derived from vertical velocity measurements for Stage 2 and 3 larvae, under different light conditions. IR: Infra Red light; CWL: Cold White Light; *D*: $\chi^2$ distance to Normal law</figcaption></figure>
velocity was –0.06 mm s⁻¹ for Stage 2 larvae and –0.02 mm s⁻¹ for Stage 3 larvae, with standard deviations lying between 0.2 and 0.3 mm s⁻¹. A χ² test showed that none of these empirical distributions follows a Gaussian law. The similarity between empirical distributions was tested using a non-parametric Kolmogorov-Smirnov test variant. The similarity hypothesis between empirical distributions for Stage 2 and Stage 3 under IR and for Stage 3 under IR and CWL was clearly rejected (1% confidence level). In contrast, the similarity hypothesis between empirical distributions for Stage 2 and Stage 3 under CWL was accepted at the 5% confidence level. However, due to the difference in swimming activity rates between these 2 stages, the respective swimming behaviours for Stage 2 and Stage 3 were treated separately.

Finally, considering a fall velocity of –0.8 mm s⁻¹ for all stages (as measured on dead larvae by Koh 2002), the swimming velocity of *Owenia fusiformis* larvae ranges between 0 and 1.6 mm s⁻¹, a value in agreement with the maximum speed reported for many annelid larvae (Chia et al. 1984). A swimming speed of 1 mm s⁻¹ is a common value for a polychaete larvae of 200 µm with at least 1 ciliary ring (Konstantinova 1969). We derive a behaviour model for *O. fusiformis* larvae from these new measurements.

## Swimming behaviour model for *Owenia fusiformis* larvae

The swimming behaviour for *Owenia fusiformis* larvae is described by a cumulative velocity distribution that depends on the larval stage and light exposure when relevant. This cumulative velocity distribution is the sum of the velocity distribution (which describes the active motion period) weighted by the swimming activity rate and the fall velocity weighted by the swimming inactivity rate. In other words, in the sea, inactive larvae would fall. We apply a fall velocity of –0.8 mm s⁻¹ for all stages (Koh 2002), except for eggs, which are given an upward velocity of 2.5 mm s⁻¹ during a period of 12 h. For Stage 1 larvae, we assume a Gaussian cumulative swimming velocity distribution having a mean upward velocity of 0.7 mm s⁻¹ with a standard deviation of 0.35 mm s⁻¹ (Rey 1995) and the Stage 1 larvae swimming activity rate is set either to 100 or 40%. When active, the average velocity is –0.06 mm s⁻¹ for Stage 2 larvae and –0.02 mm s⁻¹ for Stage 3 larvae, with standard deviations between 0.2 and 0.3 mm s⁻¹ for both. Stage 2 larvae have an active swimming behaviour 25% of the time, both in darkness and daylight. Stage 3 larvae have an active swimming behaviour 30% of the time in darkness and 8% of the time in daylight. The threshold defining the switch

between darkness and daylight behaviour is set to 10% of midday light at the free surface $L_{max}$ = 50 W m⁻², and the light attenuation within the water column, $\alpha$, equals 0.25. Light parameter values are based on a vertical light profile measured in Banyuls Bay on 8 April 2004. Stage 4 larvae fall only. The larval motion behaviour model is summarised in Fig. 10.

<figure><img src="image_11.png" /><figcaption>Fig. 10. *Owenia fusiformis*. Stochastic model for larvae motion behaviour for different ontogenic stages and light conditions</figcaption></figure>

## Larval dispersal sensitivity study

### Resolution of physical forcings

Fig. 11a,b shows a close up in Banyuls Bay of samples of the velocity fields computed on 2 grids of increasing resolution: 250 m for Grid 3 (Fig. 11a) and
<figure><img src="image_12.png" /><figcaption>Fig. 11. (a,b) Sample projection of a 3D velocity field onto a horizontal plane at 5 m below the sea surface. Arrows show current direction and grey scale shows current velocity in cm s⁻¹. (c,d) Final positions (■) of neutrally buoyant passive larvae released at the same location in southern Banyuls Bay from 18 March 1999 at 0:00 h until 21 March 1999 at 23:30 h; computed with advection only. (a,c) and (b,d) correspond to 250 m resolution (Grid 3) and 100 m resolution (Grid 4) computations, respectively. (♦) Spawning locations in Banyuls Bay (B). Depth contours ranging from 0 to 80 m are plotted with 10 m spacing (42° N Mercator projection)</figcaption></figure>

100 m for Grid 4 (Fig. 11b). A higher resolution enabling a finer bathymetric description leads to higher spatial variations in the velocity field, especially near the coast. Fig. 11c,d shows final positions within a reference area corresponding to Grid 4 coverage and reached after 22 d by neutrally buoyant passive larvae (no settling, no swimming) that were released in 1 of the 2 spawning locations (B). These final positions result from dispersal computations by advection only at each spatial resolution. A local retention rate has been defined as the ratio of the number of larvae remaining in a reference area corresponding to Grid 4 coverage (given above each plot in Fig. 11c,d) to the total number of the 5760 larvae released. Hence, with advection only, increasing the resolution of physical forcings led to similar final positions but decreased the local retention rate by a third: from 66% (250 m resolution,

Fig. 11c) to 43% (100 m resolution, Fig. 11d). Larvae trajectories within Grid 3 (not shown) indicate that the higher local retention rate for 250 m resolution is not explained by larvae returning back into the reference area after leaving it. The meteorological forcings in March–April 1999 (strong Tramontane winds) induced a persistent southward offshore current along the Catalan Coast that prevented most larvae which were expulsed from the reference area from returning back into it over the 22 d. Moreover, whatever the grid resolution, larvae released in neighbouring locations (500 m apart) in Banyuls Bay reach similar final positions and their local retention rate is about the same (not shown). Thus, sensitivity to the resolution of physical forcings is greater than to the spawning location resolution when solely the advection process is concerned.
<figure><img src="image_13.png" /><figcaption>Fig. 12. *Owenia fusiformis*. (a,b) Sample projection of a 3D Turbulent Kinetic Energy (TKE) field onto a horizontal plane at 5 m below the sea surface with grey scale in m² s⁻². (c,d) Final positions (■) of neutrally buoyant passive larvae released at the same location in southern Banyuls Bay from 18 March 1999 at 0:00 h until 21 March 1999 at 23:30 h; computed with advection and diffusion. (a,c) and (b,d) correspond to 250 m resolution (Grid 3) and 100 m resolution (Grid 4) computations, respectively. (◆) Spawning locations in Banyuls Bay (B). Depth contours ranging from 0 to 80 m are plotted with 10 m spacing (42° N Mercator projection)</figcaption></figure>

Including the effect of turbulence (Fig. 12a,b shows samples of the TKE fields) significantly changes the final dispersal patterns at both resolutions (Fig. 12c,d) when new final positions are compared to the advection only case. As in the advection only case, final positions on Grids 3 (250 m resolution) and 4 (100 m resolution) are similar. Yet, the local retention rate is 3 times larger with the 250 m resolution (26%) than with the 100 m resolution (8%). This means that discrepancies in local retention rates between the 2 resolutions are enhanced when turbulence is included. This is due to 2 factors: (1) turbulent diffusion enables larval spreading to benefit from the higher spatial variability in the 100 m resolution velocity field; and (2) the TKE itself has a higher spatial variability in the 100 m resolution

computations than in the 250 m ones. In addition, final positions and local retention rates for larvae released in 2 neighbouring spawning locations are again very similar (not shown). Hence, sensitivity to the resolution of physical forcings is even greater when both advection and diffusion processes are included, and this sensitivity is higher than the sensitivity to the spawning location resolution.

## Spawning locations along bathymetric gradients

Fig. 13 shows the final positions after 22 d computed on the finest Grid 4 (100 m resolution) for neutrally buoyant passive larvae (no settling, no swimming)
<figure><img src="image_14.png" /><figcaption>Fig. 13. *Owenia fusiformis*. Computations on Grid 4 (100 m × 100 m) of final positions (■) of neutrally buoyant passive larvae released at increasing depth and distance from the coast (from left to right) from 18 March 1999 at 0:00 h until 21 March 1999 at 23:30 h in (a to d) southern Banyuls Bay and (e to h) in front of Argelès. (♦) Spawning locations in Banyuls Bay (B) or Argelès (A). Depth contours ranging from 0 to 80 m are plotted with 10 m spacing (42° N Mercator projection)</figcaption></figure>

released in 8 locations at increasing distance from the coast, and along 2 bathymetric gradients in front of Banyuls Bay (Fig. 13a–d) and Argelès (Fig. 13e–h). Note that the bathymetry is steeper in Banyuls Bay than in front of Argelès, which means that the spawning location at 30 m depth in front of Argelès is 3 km away from the coastline while in Banyuls Bay, this same depth is located at a distance of 1.5 km. However, in the following presentation, we refer to water depth only in describing the evolution of local retention rates as spawning locations move away from the coast, although we recognize that both water depth and distance to the coastline define the hydrodynamic separation between nearshore and offshore circulations. In addition, due to the particular meteorological forcings in 1999, most larvae which exit the smallest computational domain (Grid 4) will not return back into it within the 22 d. Hence, sensitivity to spawning location is described in terms of local retention rates inside the reference area corresponding to the smallest grid (Grid 4) coverage.

At both locations, local retention rates in March–April 1999 decreased as the distance to the coast and depth at spawning location increased. In front of Argelès, the local retention rate in March–April 1999 ranged from 3% for nearshore spawning locations to zero for spawning locations deeper than 20 m. In contrast, in Banyuls Bay, the local retention rate in March–April 1999 remained close to 8% for spawning locations down to 30 m deep and falls to zero for spawning loca-

tions deeper than 40 m. In addition, larvae released in front of Argelès at 20 m and remaining in the domain spread over Banyuls Bay (8 km southward) while those released in Banyuls Bay spread locally within a 3 km radius.

## Swimming behaviour and mortality

Fig. 14 shows the final positions computed on the finest Grid 4 (100 m resolution) for larvae released in the Banyuls Bay population and accounting for: (a) advection-diffusion-settling, (c) advection-diffusion-settling with a 100% swimming activity rate for all stages, (b) and (d) advection-diffusion-settling with measured swimming activity rates for Stages 2, 3 and 4 and a variable Stage 1 swimming activity rate from 40% (Fig. 14b) to 100% (Fig. 14d).

Using even a small settling velocity (0.8 mm s⁻¹) increases the average local retention rate in the computational domain, which rises from 10 to 40% (Fig. 14a), but it also results in a wider final dispersal pattern. Adding the observed *Owenia fusiformis* swimming behaviour creates patterns intermediate between those simulated with advection-diffusion-settling (Fig. 12c) and those with advection-diffusion-settling (Fig. 14a). In particular, the final positions of larvae in Fig. 14b & d show variability due to the uncertainty in the swimming activity rate, which ranges from 40 to 100%, in Stage 1. In the following presentation, we
<figure><img src="image_15.png" /><figcaption>Fig. 14. *Owenia fusiformis*. Computations on Grid 4 (100 m × 100 m) of final positions (■) of larvae released in the same location in southern Banyuls Bay from 18 March 1999 at 0:00 h until 21 March 1999 at 23:30 h depending on various settings of the larval dispersal model: (a) advection-diffusion-settling; (b) advection-diffusion-settling and swimming activity as shown in Fig. 10; (c) advection-diffusion-settling and 100% swimming activity; (d) advection-diffusion-settling and swimming activity as (b) except Stage 1 with 100% swimming activity. Depth contours ranging from 0 to 80 m are plotted with 10 m spacing (42°N Mercator projection). (◆) Spawning location in Banyuls Bay (B)</figcaption></figure>

have chosen to set the Stage 1 swimming activity rate to 40%. This selection was based on the assumption that larvae do not swim continuously. In addition, since morphological changes between Stages 1 and 2 are less marked than between later stages, we suggest that Stages 1 and 2 could reasonably be assumed to have similar swimming activity rates.

Fig. 15 shows final position patterns for larvae released in Banyuls Bay while accounting for a daily mortality rate set to *M* = 1% and *M* = 10%, respectively, and using the swimming behaviour from Fig. 10. In our simulations, mortality decreases the local retention rate of larvae (to 15 and 7%, respectively) and limits the spread of the final positions.

Finally, the validity of our swimming behaviour model for *Owenia fusiformis* was tested against observations. Fig. 16 shows final position patterns for larvae released in 4 locations (2 in the Banyuls Bay reproduc-

tive adult area, B, and 2 in the Argelès population, A) considering the measured swimming behaviour without mortality. The maximum coverage of the *O. fusiformis* observed post-settlement deposits in Banyuls Bay in May 1999 is also shown. According to Fig. 16, the southern post-settlement area could result from larvae originated in both the Banyuls Bay and Argelès populations.

## DISCUSSION

### *Owenia fusiformis* larval swimming behaviour and vertical ontogenic migration

Previous larval dispersal studies have described the biological model *Owenia fusiformis* in the Mediterranean (Verdier-Bonnet et al. 1997) or in the English
<figure><img src="image_16.png" /><figcaption>Fig. 15. *Owenia fusiformis*. Computations on Grid 4 (100 m × 100 m) of final positions (■) of larvae released in the same location in southern Banyuls Bay from 18 March 1999 at 0:00 h until 21 March 1999 at 23:30 h accounting for advection-diffusion-settling and swimming activity as shown in Fig. 10 with a daily mortality rate set to (a) 1% and (b) 10%. Depth contours ranging from 0 to 80 m are plotted with 10 m spacing (42° N Mercator projection). (♦) Spawning location in Banyuls Bay (B)</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>Fig. 16. *Owenia fusiformis*. Computations on Grid 4 (100 m × 100 m) of final positions (■) of 5760 larvae released in each of 4 locations (2 in the Banyuls Bay population; 2 in the Argelès population) from 18 March 1999 at 0:00 h until 21 March 1999 at 23:30 h with advection-diffusion-settling and swimming activity as shown in Fig. 10. Depth contours ranging from 0 to 80 m are plotted with 10 m spacing (42° N Mercator projection). (♦) Spawning location in Banyuls Bay (B) or Argelès (A)</figcaption></figure>
Channel (Thiébaut et al. 1992, Barnay et al. 2003). Those studies neglected active swimming by *O. fusiformis* and considered only ontogenic-related variations in settling velocities. Fig. 12d shows dispersal patterns corresponding to the larval behaviour suggested by those studies, that is, taking a zero settling velocity for larvae younger than 20 d. In particular, Thiébaut et al. (1992) derived a settling velocity evolution from ontogenic migration observed in the English Channel. Without vertical advection (2D hydrodynamics), such an ontogenic migration could result from either a passive settling or a decrease in swimming activity. Our experimental results confirm the larvae's ability to swim with almost no difference in range for Stages 2 and 3, but with a swimming activity rate always less than 50% for stages older than 8 d and continuing to decrease with age.

Estimating the importance of vertical migration resulting from ontogenic changes either in swimming behaviour or in settling velocity appears difficult to accomplish from observations at sea. Vertical motion of the flow can hide part of the larvae's downward migration or even simulate it. Fig. 17a shows the presence of a small amount of vertical motion in our simulations (about 0.1 mm s⁻¹). Such a 3D velocity field is a general feature of nearshore wind-driven circulation. In our simulations without either settling or swimming, the average vertical position evolves from 21 m below the sea surface for Stage 1 larvae, to 47 m below the surface for Stage 2 larvae, a return up to the surface of the sea for Stage 3 larvae and finally a descent down to 47 m below the surface for Stage 4 larvae. This means that in the absence of any ontogenic mechanisms, the water flow alone can, in some cases, simulate a vertical migration resembling an ontogenic migration (Fig. 17b). In the mega-tidal English Channel environment, the water flow is strongly 2D but a weak wind-driven circulation superimposed on it can introduce small vertical motions. Such vertical movements may have introduced a bias in indirect estimates of the *Owenia fusiformis* settling velocity at sea, especially when settling velocities are on the order of magnitude of 1 mm s⁻¹. This may explain the slower settling velocity found by Thiébaut et al. (1992) as compared to recent direct settling estimates (Koh 2002) and our experimental estimates of swimming velocities.

## Selection of physical forcings resolution and a physical barrier

The numerical study by Verdier-Bonnet et al. (1997) on *Owenia fusiformis* recruitment in Banyuls Bay drew conclusions on recruitment locations and success rate from 3D advection-only computations performed on a 200 m × 200 m × 4 m Cartesian grid, considering 72 spawning locations enclosed in 600 m × 800 m zones. Verdier-Bonnet et al. (1997) suggested that their computations could be improved by increasing mesh and wind-forcing resolutions as we have done. However, it is worth noting that our use of σ coordinates following the bathymetry instead of cartesian ones improved the vertical resolution near the coast. This explains why our 250 m × 250 m computations already show a higher spatial variability than their 200 m × 200 m computations. Nevertheless, increasing mesh resolution when only advection is accounted for has less effect than including turbulence. Turbulence was neglected in the Verdier-Bonnet et al.

<figure><img src="image_18.png" /><figcaption>Fig. 17. *Owenia fusiformis*. (a) Sample of projection of 3D velocity field onto a vertical plane (small arrows show vertical velocities and grey scale shows the velocity field intensity in cm s⁻¹) computed on Grid 4 (100 m × 100 m). (◆) shows 1 of the 2 spawning points. A small vertical motion is observed. (b) Vertical locations 11 d after spawning began for Stage 1 and Stage 2 larvae. Although settling was not included, vertical motion simulates a downward vertical migration with ontogenic evolution</figcaption></figure>
(1997) study although they did emphasise that it was computed by their circulation model.

Since the dispersal computations are more sensitive to the physical forcings resolution or to turbulence than to the precise locations of spawning within the reproductive adult area, we believe their conclusions should be revised. But what should the minimal 3D physical forcings resolution be to reach stable dispersal? First, the grid resolution should accurately describe the bathymetry in the computed domain. For example, in Banyuls Bay, a 100 m resolution accurately describes the bathymetry up to 0.5 m depth, which is also the limit of our computations. Secondly, the grid resolution should not filter the atmospheric forcings. In Banyuls Bay, such a high resolution (100 m) in atmospheric forcings is not currently available, although the strong surrounding orography may lead to small-scale spatial variability.

In summary, we believe that larval dispersal off the Catalan coast can be accurately studied using a 100 m × 100 m resolution on a σ grid in 3D physical forcings. Nevertheless, as the low local retention rate suggests, larval transfer over a larger area than our smallest computational domain (Grid 4) should be considered. Taking advantage of the nested procedure used for the hydrodynamic computations, we are already able to compute the dispersal of larvae after they leave the smallest high resolution domain (Grid 4) over a larger low resolution domain (Grid 3). However, ideally it is preferable to keep a high resolution, especially near the coast, while increasing the size of the dispersal area. Due to computer capacity limitations, this requires computing physical forcings (velocities and TKE) at a high resolution on the adjacent domains and linking them while computing the dispersal offline. This is the object of future work.

Considering various spawning locations for neutrally buoyant passive larvae (no settling, no swimming), computations for March–April 1999 showed a strong southward larval flux with no physical barrier between Argelès and Banyuls benthic invertebrates populations. However, local retention rates dropped rapidly to zero for larvae released at both increasing depths and distance from the coastline. This result does suggest the existence of a physical barrier between nearshore and offshore larval dispersal in March–April 1999. The barrier was located at around 20 m deep in Argelès and 30 m deep in Banyuls, and 3 km from the coastline at both sites.

## Larval settling, swimming behaviour and mortality

Including even a small settling velocity (0.8 mm s⁻¹) affected larval dispersal in our 3D simulations (Fig. 14a). While accounting for vertical migration will not affect the planar dispersal patterns when physical

forcings are 2D (due to the absence of shear on the vertical), it does when 3D physical forcings are considered. In this respect, the settling process becomes a key parameter when both advection and diffusion are 3D. However, in contrast with passive particles, the larval settling process results from the balance between passive settling and active swimming behaviours that may contribute to reduce it.

In our measurements of *Owenia fusiformis* larval swimming behaviour, swimming speed compensates the settling speed, on average and this balance averages 0.02 mm s⁻¹ with a standard deviation of 0.3 mm s⁻¹. Disregarding the measured swimming activity rate, i.e. using a 100% swimming activity rate, the final dispersal pattern is very close to that of advection-diffusion only (Fig. 14b). Accounting for a 100% swimming activity rate means adding a stochastic velocity ranging between -0.3 and 0.3 mm s⁻¹ to the vertical turbulent velocities. Since turbulent velocities (typically 10 mm s⁻¹) are at least 1 order of magnitude larger than swimming speed, if a 100% swimming activity rate compensates the settling process on average, the overall swimming-settling larval behaviour may be neglected.

In our experiments, measured swimming activity rates for Stages 2, 3 and 4 were less than 40%. Thus, settling process affects these stages. However, some uncertainty remains over the Stage 1 swimming activity rate which was set to either 40 or 100% in the computations. It appears that when the settling process is assumed not to affect larvae younger than 8 d (100% swimming activity rate for Stage 1, Fig. 14d), the dispersal pattern of up to 4 km from the spawning location and local retention rate (15%) are very similar to those without settling at any stages. This means that 85% of released larvae exited the smallest computational domain (Grid 4) before they are 8 d old and 15% of released larvae also reach their final position also before the age of 8 d. In contrast, when a lower swimming activity rate is given to Stage 1 larvae (40%), their final positions are very similar to those obtained with settling only, except that the local retention rate decreases from 40% (without swimming) to 20% (with swimming). This result is interpreted as a composite of the advection-diffusion patterns of final positions for larvae with settling (60%) and without settling (40%). Considering the uncertainty about the swimming activity of Stage 1 larvae, a swimming rate of 40% was used in the following sections.

Varying daily mortality rates from *M* = 1% to *M* = 10% led to a survival rate (after 22 d) ranging from 80 to 9.8%, respectively, of the total released larvae pool. Since the probability of death of the same proportion of larvae in 2 complementary subgroups of the released pool is lower than 100%, the mortality procedure may have a non-linear effect on local retention rates. Consequently, local retention rates with mortal-
ity (respectively of 15 and 7%) are not exactly equal to the local retention rate without mortality (20%) multiplied by their respective survival rates. It can also be expected that conditional mortality will affect final positions and local retention rates even more. Future versions of the model will be improved by including aspects of the trophic dynamics. Note that the turbulence effect in predator/prey encounters can be easily included in our model (Werner et al. 2001).

## Simulations versus observations in 1999

In our simulations, final position dispersal patterns for larvae released in Banyuls Bay suggest that self-recruitment, even though it remains possible, was not important in the Banyuls Bay population during spring 1999 (Fig. 16a,b). In addition, connections with other populations in both the south and the north were possible. One striking feature of the patterns is the potential interconnection, although not symmetrical, with the Argelès population. Few larvae released in Banyuls Bay reached the Argelès population while many larvae released in Argelès reached Banyuls Bay. However, such dispersal patterns would have to be modified by recruitment constraints (e.g. substrate suitability) to obtain an effective post-metamorphosed worm settlement area and conclusions about larval flux amongst populations. For instance, in our simulations, many larvae have final positions on the rocky shores, but it is known that no *Owenia fusiformis* can survive there since this worm needs sand to make its tube. Without going into details of recruitment constraints for this species, we argue that simulated dispersal patterns should at least include effective post-settlement areas. This means that both simulated spatial coverage and densities of deposits are expected to be larger than the observed ones. From Koh's (2002) observations, the maximum coverage of the *O. fusiformis* post-settlement area in Banyuls Bay in May 1999 was in the middle of each Banyuls Bay cove. Our computed larval dispersal patterns include Koh's (2002) observations (shown by the solid contour on Fig. 16a–d). With respect to densities, observed post-settlement densities in May 1999 (less than 100 ind. m⁻²) were low compared to the observed densities of reproductive adults (2000 ind. m⁻²). After spatial integration of the data in Fig. 5, these observed densities gave rough estimates for the total number of post-settlement worms (10⁸ ind.) and the total number of reproductive adults in the Banyuls Bay population (5 × 10⁸ ind.). Assuming a 1:~1 sex ratio with a 10% fertilization rate for 20 000 ovocytes per female (according to average fertility of 2 yr old *O. fusiformis* females estimated in the English Channel by Thiébaut & Dauvin 1991), we suggest that about 1000 larvae per

adult would be released during an entire spawning cycle. Thus, we estimate from observed densities of reproductive adults that 5 × 10¹¹ larvae were released by the Banyuls Bay population in 1999. Compared to the number of observed post-settlement worms, only 1 out of 5000 larvae released in Banyuls settled in one of the 2 coves of Banyuls Bay in 1999. This estimate, based on observations, agrees with our simulations, in which only 1 larva reaches the same exact area where settlement has been observed. After doubling the number of released larvae (11 520 larvae) in our simulations, final dispersal patterns are similar to those using 5760 larvae, although precise final locations of outliers are not exactly the same. No settlement occurs in the same exact area where settlement has been observed, but settlements occur in the vicinity.

Fig. 18 shows the final positions reached by 4488 larvae after they leave the smallest high resolution domain (Grid 4, Fig. 16b) and are dispersed over the larger lower resolution domain (Grid 3). Only 88 larvae of the larvae transported off the smallest domain (less than 2% of the total) return back into it, which confirms that: (1) our small domain size does not introduce a bias in *Owenia fusiformis* local retention rate esti-

<figure><img src="image_19.png" /><figcaption>Fig. 18. *Owenia fusiformis*. Computations on Grid 3 (250 m × 250 m) of final positions of the 4488 larvae released in 1 of the 2 spawning locations in Banyuls Bay and transported off the smallest high resolution domain (Grid 4) through advection-diffusion-settling and swimming activity as shown in Fig. 10. Number of larvae returning back into the small domain (Grid 4, outlined by box) is indicated above the plot. Depth contours ranging from 0 to 80 m are plotted with 10 m spacing (42° N Mercator projection)</figcaption></figure>
mates for Banyuls Bay in spring 1999, and (2) post-settlement observations in Banyuls Bay in May 1999 represent insignificant deposits.

## CONCLUSIONS

The population dynamics of most benthic invertebrate species are driven by their planktonic larval stage, in particular through larval dispersal into inadequate or adequate habitats by hydrodynamics (Crisp 1976). However, the effect of dispersal may vary among species and hydrodynamic environments (de Wolf 1974, Hill et al. 1996, Verdier-Bonnet et al. 1997). Numerical simulations can be a powerful tool for estimating specific meroplanktonic population drops due to dispersal. This paper details a Lagrangian stochastic larva-tracking model driven by a high resolution 3D wind-driven current model. The larva tracking model includes 3D advection and turbulence, and a species-specific larval swimming behaviour. Ontogenic changes, sensitivity to light exposure and inter-individual variability can be specified in the swimming behaviour. Larval mortality is also included. In this paper, the dispersal model is applied to *Owenia fusiformis* larvae, in the Mediterranean Sea. To this end, a swimming behaviour model for *O. fusiformis* larvae was developed and tested with larvae velocities, resulting from measurements of both settling and swimming behaviour, ranging between –1 and 0.9 mm s⁻¹ and swimming activity rates lower than 50%.

The sensitivity study of larval dispersal in March–April 1999 showed that: (1) the dispersal of neutrally buoyant passive larvae is more sensitive to the resolution of physical forcings through both advection and diffusion processes than to the variability of spawning locations within neighbouring grid cells in Banyuls Bay (up to 1 km apart); and (2) a physical barrier located at 20 m deep in Argelès and 30 m deep in Banyuls Bay separated nearshore and offshore larval dispersal in 1999. For *Owenia fusiformis* larvae released in Banyuls Bay and Argelès, we found that: (1) final positions and local retention result from the balance between the 3D turbulence, larval settling velocity (although it is small ~0.8 mm s⁻¹) and swimming activity rate, and (2) natural mortality also affects final larval position and local retention, although it is not proportional to survival rates.

Computations of high resolution dispersal patterns for *Owenia fusiformis* in Banyuls Bay suggest that self recruitment was low in the Banyuls Bay population during spring 1999 and agrees with the insignificant post-settlement deposit patterns observed in Banyuls Bay in May 1999. In addition, interconnections between Argelès and Banyuls populations can exist.

Finally, since settling velocity and maximum swimming speed for *Owenia fusiformis* are the same order of

magnitude, the amount of time the worm spends swimming can compensate settling and affects dispersal. As *O. fusiformis* settling velocity and swimming speed fall within the average potential of many annelids (Chia et al. 1984), we recommend that dispersal studies of annelid species are based on an accurate knowledge of passive settling and active swimming behaviours throughout their ontogenic changes and as a function of environmental factors (e.g. light, starvation). To this end, actographic equipment is a useful tool.

For *Owenia fusiformis* larval dispersal, our model needs further validation. Ongoing studies include comparing plankton field sampling results and dispersal computations at various locations and times in Banyuls Bay. After the validation is completed, this larval dispersal model will be used to study inter-annual variability of *O. fusiformis* self-recruitment and potential larval flux between distant populations.

*Acknowledgements.* This work was partly funded by the French Programme National d'Environnement Côtier under Research Action 2: Dynamique des Populations. We thank the crew of the 'Nereis II' for their help in larvae collection at sea, and J. Neveux and F. Lantine for providing *Nitzchia* sp. algae to feed *Owenia fusiformis* larvae. We also thank METEO FRANCE for supplying the atmospherical forcings and C. Labrune for providing unpublished macrofauna data. Many thanks to F. Charles and 2 anonymous reviewers for their helpful comments.

## LITERATURE CITED

* Amouroux JM (1974) Etude des peuplements infralittoraux de la côte du Roussillon. II: Faunistique et caractéristiques bionomiques. Vie Milieu 24(1B):209–222

* Auclair F, Casitas S, Marsaleix P (2000) Application of an inverse method to coastal modelling. J Atmos Oceanic Technol 17:1368–1391

* Barnay AS, Ellien C, Gentil F, Thiébaut E (2003) A model study on variations in larval supply: Are populations of the polychaete *Owenia fusiformis* opened or closed? Helgol Mar Res 56:229–237

* Blumberg AF, Mellor G (1987) A description of a three dimensional coastal circulation model. In: Heaps N (ed) Three dimensional coastal ocean model. American Geophysical Union, Washington, DC, p 1–16

* Carlotti F (1996) Simulations of the population dynamics of *Pectinaria koreni* in a 2-D vertical advection-diffusion model of the Seine estuary. Oceanolog Acta 19(3–4): 369–383

* Chen YH, Shaw PT, Wolcott TG (1997) Enhancing estuarine retention of planktonic larvae by tidal currents. Estuar Coast Shelf Sci 45:525–533

* Chia FS, Buckland-Nicks J, Young CM (1984) Locomotion of marine invertebrates larvae: a review. Can J Zool 62: 1205–1222

* Crisp DJ (1976) The role of pelagic larvae. In: Davies PS (ed) Perspectives in experimental biology, Vol 1. Zoology. Pergamon Press, Oxford, p 145–155

* de Wolf P (1974) On the retention of marine invertebrates in estuaries. Thalassia Jugosl 10(1/2):415–424
Duchêne JC, Queiroga H (2001) Use of an intelligent CCD camera for the study of endogenous vertical migration rhythms in first zoeae of the crab *Carcinus maenas*. Mar Biol 139(5):901–909

Dufau-Julliand C, Marsaleix P, Petrenko A, Dekeyser I (2004) Three-dimensional modeling of the Gulf of Lion's hydro-dynamics (northwest Mediterranean) during January 1999 (MOOGLI3 Experiment) and late winter 1999: West-ern Mediterranean Intermediate Water's (WIW's) forma-tion and its cascading over the shelf break. J Geophys Res 109(C11): doi:10.1029/2003JC002019

Ellien C, Thiébaut E, Barnay AS, Dauvin JC, Gentil F, Salomon JC (2000) The influence of variability in larval dispersal on the dynamics of a marine metapopulation in the eastern Channel. Oceanol Acta 23(4):423–442

Ellien C, Thiébaut E, Dumas F, Salomon JC, Nival P (2004) A modelling study of the respective role of hydrodynamic processes and larval mortality on larval dispersal and recruitment of benthic invertebrates: examples of *Pectinaria koreni* (Annelida: Polychaeta) in the Bay of Seine (English Channel). J Plankton Res 26(2):117–132

Estournel C, Kondrachoff V, Marsaleix P, Vehil R (1997) The plume of the Rhone: numerical simulation and remote sensing. Cont Shelf Res 17:899–924

Estournel C, Durrieu de Madron X, Marsaleix P, Auclair F, Julliand C, Vehil R (2003) Observation and modelisation of the winter coastal oceanic circulation in the Gulf of Lion under wind conditions influenced by the continental orography (FETCH experiment). J Geophys Res 108(C3): doi:10.1029/2001JC000825

Gambi MC, Giangrande A (1986) Distribution of soft-bottom polychaetes in two coastal areas of the Thyrrenian Sea (Italy): structural analysis. Estuar Coast Shell Sci 23:847–862

Gaspar P, Gregoris Y, Lefevre JM (1990) A simple eddy kinetic energy model for simulations of the oceanic vertical mixing: tests at station Papa and long-term upper ocean study site. J Geophys Res 95:179–193

Geernaert GL (1990) Bulk parameterizations for the wind stress and heat fluxes, in surface waves and fluxes, Vol 1. In: Geernaert GL, Plant WJ (eds) Current theory. Kluwer Academic, Norwell, MA, p 91–172

Guille A (1971) Bionomie benthique du plateau de la côte catalane française. II: les communautés de la macrofaune. Vie Milieu 21(1B):149–280

Hill AE, Brown J, Fernand L (1996) The western Irish Sea gyre: a retention system for Norway lobster (*Nephrops norvegicus*). Oceanol Acta 19(3-4):357–368

Johnson DF, Hess KW (1990) Numerical simulations of blue crab larval dispersal and recruitment. Bull Mar Sci 46(1): 195–213

Koh BS (2002) Le cycle de vie d'Owenia fusiformis (Annélide polychète) en Méditerranée: les processus locaux de recrutement et la distribution géographique de l'espèce à l'échelle mondiale. Doctorat de l'Université Paris VI, Paris

Konstantinova MI (1969) Movement of polychaete larvae. Dokl Akad Nauk SSSR 188:942–945

Maron-Ramos J (1973) Annélides polychètes du plateau continental de la côte espagnole. Systématique et écolgie. Doctorat de l'Université Paris VI, 2 tomes, Paris

Massé H (1971) Contribution à l'étude quantitative et dynamique de la macrofaune de peuplements des sables fins infralittoraux des côtes de Provence. Doctorat de l'Université d'Aix-Marseille, Marseille

Oey LY, Chen P (1992) A model simulation of circulation in the northeast Atlantic shelves and seas. J Geophys Res 97: 20087–20115

Ohba R, Shao Y, Kouchi A (1997) A wind tunnel and numeri-

cal investigation of turbulent dispersal in coastal atmospheric boundary layers. Boundary-Layer Meteorol 87(2): 255–273

Oliveira LA, Costa VAF, Baliga BR (2002) A Lagrangian-Eulerian model of particle dispersal in a turbulent plane mixing layer. Int J Numer Meth Fluids 40:639–653

Oliver JK, King BA, Willis BL, Babcock RC, Wolanski E (1992) Dispersal of coral larvae from a lagoonal reef—II. Comparisons between model predictions and observed concentrations. Cont Shelf Res 12(7/8):873–889

Ommundsen A (2002) Models of cross shelf transport introduced by the Lofoten Maelstrom. Cont Shelf Res 22:93–113

Pedersen OP, Aschan M, Rasmussen T, Tande KS, Slagstad D (2003) Larval dispersal and mother populations of *Pandalus borealis* investigated by a Lagrangian particle-tracking model. Fish Res 65:173–190

Petrenko A, Leredde Y, Marsaleix P (2005) Circulation in a stratified and wind-forced Gulf of Lions, NW Mediterranean Sea: in situ and modeling data. Cont Shelf Res 25: 7–27

Pinedo S, Sardà R, Rey C, Bhaud M (2000) Effect of sediment particle size on recruitment of *Owenia fusiformis* in the Bay of Blanes (NW Mediterranean Sea): an experimental approach to explain field distribution. Mar Ecol Prog Ser 203:205–213

Porch CE (1998) A numerical study of larval fish retention along the southeast Florida coast. Ecol Model 109(1):35–59

Rumrill SS (1990) Natural mortality of marine invertebrate larvae. Ophelia 32:163–198

Rey C (1995) Originalité du recrutement benthique chez *Owenia fusiformis* polychaeta: dissémination larvale et rapport avec le substrat. Mémoire du DEA Océanologie Biologique de l'Université Paris VI, Paris

Scheltema RS (1986) On dispersal and planktonic larvae of benthic invertebrates: an eclectic overview and summary of problems. Bull Mar Sci 39(2):290–322

Thiébaut E (1994) Dynamique du recrutement et dispersal larvale de deux annélides polychètes *Owenia fusiformis* and *Pectinaria koreni* en régime mégatidal (Baie de Seine orientale, Manche). Doctorat de l'Université Paris VI, Paris

Thiébaut E, Dauvin JC (1991) Etude de la fécondité d'Owenia fusiformis Delle Chiade (Annélide polychète). J Rech Océanogr 16(1,2):39–41

Thiébaut E, Dauvin JC, Lagadeuc Y (1992) Transport of *Owenia fusiformis* larvae (Annelida: Polychaeta) in the Bay of Seine. I. Vertical distribution in relation to water column stratification and ontogenic vertical migration. Mar Ecol Prog Ser 80:29–39

Thorson G (1946) Reproduction and larval development of Danish marine bottom invertebrates, with special reference to the planktonic larvae in the Sound (Øresund). Medd DanFisk-Havunders (1):1–523

Ulses C, Grenz C, Marsaleix P, Schaaf E, Estournel C, Meulé S, Pinazo C (2005) Circulation in a semi-enclosed bay under influence of strong fresh water input. J Mar Syst 56: 113–132

Verdier-Bonnet C, Carlotti F, Rey C, Bhaud M (1997) A model of larval dispersal coupling wind-driven currents and vertical larval behaviour: application to the recruitment of the annelid *Owenia fusiformis* in Banyuls Bay, France. Mar Ecol Prog Ser 160:217–231

Werner FE, Mackenzie BR, Perry RI, Lough RG, Naimie CE, Blanton BO, Quinlan JA (2001) Larval trophodynamics, turbulence, and drift on Georges Bank: a sensitivity analysis of cod and haddock. Sci Mar 65(Suppl1):99–115

Wilson DP (1932) On the Mitraria larva of *Owenia fusiformis* Delle Chiade. Phil Trans R Soc Ser B 221:231–334