Numerical modelling of phytoplankton bloom in the upwelling ecosystem of Cabo Frio (Brazil)

Ecological Modelling 116 (1999) 135 – 148

Numerical modelling of phytoplankton bloom in the upwelling ecosystem of Cabo Frio (Brazil) Carlos A.A. Carbonel H. a,*,1, Jean Louis Valentin b a

Laborato´rio Nacional de Computac¸a˜o Cientı´fica, A6. Getulio Vargas 233, Cep: 25651 -070 Petro´polis, RJ, Brazil b Instituto de Biologia, Uni6ersidade Federal do Rio de Janeiro, Cidade Uni6ersitaria— Ilha do Funda˜o, Cep: 21949 -900 Rio de Janeiro, RJ, Brazil Received 20 January 1998; accepted 19 August 1998

Abstract The phenomenon of phytoplankton bloom during upwelling events in the coastal region of Cabo Frio (Brazil), is investigated using a 112 ocean model including the physics of coastal waters and the biological changes of the marine primary biomass. The vertical structure of the coupled, physical – biological model is described by an active layer overlaying a deep inert layer where the pressure gradient is set to zero. The physical model describes the changes of momentum, mass and heat in the dynamic layer and is forced by wind acting at the surface. The biological model describes the changes of its three components (nutrients, phytoplankton and zooplankton) and is forced by the primary biomass concentration injected into the upper layer as a consequence of upwelling favorable winds and by the light producing the photosynthesis of the phytoplankton. Coupling is established by the horizontal-velocity, layer-thickness and upwelling velocity fields from the physical solution. The equations are solved numerically in space and time by the finite difference method. Model results forced by transient winds are arranged initially for a standard run causing coastal upwelling. During the spin-up of winds, the nutrient, phytoplankton and zooplankton concentrations are injected into the upper layer by upwelling. Consequently a maximum of nutrients takes place (generated by the nutrient injection), decreasing from the increasing assimilation by phytoplankton and by reduction of nutrients supply (spin-down stage of winds). One day after the maximum of nutrients takes place, the phytoplankton bloom occurs with amplitude and duration similar to that observed in this region. This bloom is very short due the increase of grazing by zooplankton. The results of the standard run, reproduce observed features of physical and biological fields during upwelling events; in particular, the time dependent pattern of the short-lived phytoplankton bloom, suggesting that the physical–biological model presented here, contains much of this important time-dependent phenomenon. Experiments indicate that the change of transient wind duration affects the time response of the biological components, changing the instant of the maximum of nutrients concentration, but keeping unchanged the

* Corresponding author. Tel.: +55-024-2336226; fax: +55-024-2315595; e-mail: carbonel(@lncc.br. 1 ´ rsula, Instituto de Cie`ncias Biolo´gicas e Ambientais, Rua Jornalista Orlando Dantas 59, Present address: Universidade Santa U Botafogo, 22231-010 Rio de Janeiro, Brazil. Tel.: + 55-021-5525422; fax: +55-021-5516446. 0304-3800/99/$ - see front matter © 1999 Elsevier Science B.V. All rights reserved. PII: S 0 3 0 4 - 3 8 0 0 ( 9 8 ) 0 0 2 0 1 - 4

136

C.A.A. Carbonel H., J.L. Valentin / Ecological Modelling 116 (1999) 135–148

instant of maximum of phytoplankton. The reduction of the constant of optimum photosynthesis rate affects the bloom of phytoplankton delaying it and reducing its amplitude. Experiments indicate that the maximum grazing rate is an important factor in limiting the persistence and amplitude of the phytoplankton bloom. © 1999 Elsevier Science B.V. All rights reserved. Keywords: Numerical model; Upwelling; Primary biomass; Phytoplankton bloom; Cabo Frio (Brazil)

1. Introduction Coastal upwelling is a mesoscale phenomenon (  100 km) which is forced by winds blowing along the continental margins. Due to the Coriolis effect, a surface flow component in offshore direction is generated. This offshore flow, called the Ekman drift is compensated by the raising of deep-ocean waters to the surface bringing the cold deep nutrient-rich waters into the euphotic surface layer, where the photosynthesis takes place as nutrients become available to phytoplankton. The resulting consequence is high primary productivity and high fish yields. The most important upwelling zones are located off the west coast of South America, West Africa and the west coast of North America. In Brazil, the coastal region of Cabo Frio is characterized by the occurrence of very rapid upwelling events and is considered an unusual system of the tropical Brazilian coast. In this region, over the coastal band of 25 nautical miles (42°30%W – 42°06%W), the subtropical shelf waters (12– 15°C and 35.1 – 35.5‰) constitute the principal water mass of the upwelling process (Mascarenhas et al., 1971). The occurrence of phytoplankton blooms during the summer is reported by Valentin (1984) and Valentin et al. (1985). The duration of such blooms is less than 24 h, during upwelling events in the coastal region of Cabo Frio. The blooms occurred generally, after 3 or 4 days of upwelling favorable NE winds with velocities greater as 10 knots, bringing cold deep water with nutrients concentrations to the surface and resulting in biological production. The sea surface temperature (SST) decreases to 15°C, the nutrients concentration increases up to 10–12 mM NO3 – N (less than 0.5 mg l − 1 chlorophyll a) and the local phytoplankton

bloom occurs near shore with concentrations of up to 6 mg l − 1 of chlorophyll as presented in Fig. 1b (observations at fixed station, 22°59%06%%S and 42°13%0%%W). When the NE wind velocities decrease, or by the presence of polar fronts the direction of the wind changes, generating downwelling conditions, the nutrient concentration decreases to 0.4 mM and the SST increases to 23°C. Coastal ocean models range in dynamical complexity from simple, linear surface–layer models to sophisticated, nonlinear, continuously stratified systems. To represent ocean processes and at the same time keep their simplicity in the vertical structure, layered models must allow for lateral inhomogeneities. The simplest possibility correspond to vertically averaging in each layer, velocity, layer-thickness and concentrations fields. The physical inhomogeneous layers models were introduced more than 2 decades ago by Lavoie (1972) in a study of the effects of lakes on the lower atmosphere circulation and later Schop and Cane (1983) for the study of equatorial dynamics. Since then, many works which use the same idea were published (e.g. Anderson, 1984; Anderson and Mc Creary, 1985; Mc Creary and Kundu, 1988; Ripa, 1993; Carbonel, 1998). Since the classic prey–predator equation of Lotka–Volterra, marine scientists have attempted to describe the processes of biological production in the ocean by deducing equations which represent the interaction of biological elements with their environment (Patten, 1968). The experience in mathematical modeling simulation of planktonic systems was systematized in different upwelling regions: Peru (Walsh, 1975), Oregon and California (Wroblewski, 1980), Northwest (Riley, 1965) and Southwest (Jones and Henderson, 1987) Africa coast and West Florida (O’Brien and Wroblewski,

C.A.A. Carbonel H., J.L. Valentin / Ecological Modelling 116 (1999) 135–148

137

Fig. 1. (a) The coastal region of Cabo Frio (23°S latitude) located at the central-southern Brazilian littoral characterized by a coastal upwelling due to prevailing E–NE winds; and (b) Phytoplankton bloom observed at a fixed offshore station located at the west side of Cabo Frio (from 17 to 21 of January 1986). For comparison, according to Radach and Meier Raimer (1975), 1 mg l − 1 chlorophyll is equal to 0.5 mM nitrogen.

1973). Walsh and Dugdale (1971) have pioneered the construction of spatial models in the Peruvian upwelling ecosystem. For the Arabian Sea, Mc Creary et al. (1996) studied the processes that determines the annual cycle of biological activity using a coupled, physical – biological 212 layer model. At the Cabo Frio upwelling, Valentin and Coutinho (1990) made a preliminary approach of modelling maximum chlorophyll at the sea surface using the Stella program starting from assumptions suggested by data and incorporating as much as possible of what was known about the

processes. Valentin (1992) presented a time dependent one-dimensional model for the vertical distribution of primary biomass of the Cabo Frio upwelling in which the light and grazing pressure were forcing functions on phytoplankton. The results confirm that vertical distribution of phytoplankton maximum depends on depth of the thermocline. The intention of this paper is to contribute in the understanding of the space-time response of coastal areas in relation to the primary production associated with upwelling events. Therefore, a deterministic model based on the conservation

138

C.A.A. Carbonel H., J.L. Valentin / Ecological Modelling 116 (1999) 135–148

equations of physical and biological variables is proposed to study the phenomenon of phytoplankton bloom in upwelling regions, applied to the coast of Cabo Frio (Brazil). We use a coupled, physical – biological 112 layer model (Fig. 2), in which the conservation equations are solved numerically in space and time for a dynamic upper layer. The physical model is based on the 112 layer model described by Carbonel (1998). Among other things, this model provides a satisfactory representation of the rapid changes of SST and current system in the Cabo Frio Coast. Additionally, the implemented weakly reflective open boundary conditions allow to treat adequately open and limited coastal areas.

The biological model consists of a set of transport equations for the upper layer determining the three biological concentrations: nutrients, phytoplankton and zooplankton. These equations were adapted for the 112 layer model formulation and are based on system equations used by Walsh and Dugdale (1971) and Walsh (1975) for simulating nitrogen flow in the Peruvian upwelling system, by O’Brien and Wroblewski (1973) for studying trophic levels off West Florida, Valentin and Coutinho (1990) and Valentin (1992) in the Cabo Frio (Brazil) upwelling. Steele and Henderson (1981), Andersen (1985) and Sekine et al. (1991) were also consulted for building our equations. The coupled physical–biological model is relatively simpler because only one dynamic layer is considered, but its simplicity allow us to evaluate the influence and role of different factors in a simple way. The influence of physical and biological processes is analyzed for study cases forced by transient upwelling—favorable NE winds. Particularly, the influence on phytoplankton bloom intensity of transient wind duration and some biological parameter, is discussed.

2. The coupled physical–biological model We consider two layers of incompressible fluid flow in a coastal ocean domain V with boundary G(G =G0 @G1; G0: land type boundary; G1: ocean type boundary). The planar Cartesian coordinates xi (i= 1, 2) and the summation convention with repeated indices are used in this paper.

2.1. Physical model The physical model is based on the hydrodynamic and thermodynamic model described by Carbonel (1998). The conservation of momentum, mass and heat is described by the vertical integrated equations for the dynamic upper layer in an open and limited coastal upwelling region,

112

Fig. 2. The coupled physical–biological layer model: (a) The physical structure in the space coordinates; and (b) The biological processes which take place in the upper layer of the model.

!

"

(Ui (uj Ui (h hu (T + + oij Uj + gh s + + nUi (t (xj (xi 2m¯ (xi −

ti =0 ru

(1)

C.A.A. Carbonel H., J.L. Valentin / Ecological Modelling 116 (1999) 135–148

(h (Ui + − we = 0 (t (xi

(2)

(T 1 (T + ui + (q− Q)= 0 (xi h (t

(3)

where oij =



0 f

−f 0



s =(r l −r u)/r l,

m¯ = m/(m −s) and ui Ui = uih we v h ti cw Wi 6 f r u, r l rair q T Tl T* u Q H ts

terface between the two layers from surfacing adopting the smooth function (Mc Creary and Kundu, 1988), we = (He − h)2/teHe, we = 0,

m =r u/r l,

are the velocity components denotes the corresponding flux in the upper layer is the entrainment velocity is a vertical upwelling velocity defined equal to (Ui /(xi represents the upper layer thickness the wind stress components defined equal to cw rairWi W is the wind drag coefficient is the wind velocity component is the Rayleigh friction coefficient is the Coriolis parameter are ocean water densities in the upper and lower layer respectively. is the air density source and sink of SST equal to v T* when T]T l and equal zero otherwise is the SST is the constant temperature in the lower layer is T value at the interface between the 2 layers is the coefficient of thermal expansion is the surface heat flux defined as (T u–T)H/ts, (T u initial SST) initial thickness of upper layer is the surface heat flux time scale

In the momentum Eq. (1) the mixing terms were included into the Raleight friction term (6Ui ) parameterizing in this way all the dissipative losses. In the continuity Eq. (2), in the third term, the entrainment velocity we only prevent the in-

139

for hBHe

for h]He

(4) (5)

where He is the entrainment thickness and, te is the entrainment time scale. In the equation of heat conservation (3), the last term included the source and sink of SST q and the surface heat flux Q. The definition of the surface heat flux Q follows the parameterization proposed by Mc Creary and Kundu (1988), which is similar to the one proposed by Haney (1971).

2.2. Biological model The biological model for the upper layer of the coastal ocean is a vertically integrated version in space and time based on the primary production local model described by Valentin (1992). The time dependent conservation equations for the vertical integrated variables in the upper layer such as nutrient N, phytoplankton P and zooplankton Z are described by: (N (N SN +ui − (eZ− aP)− =0 (xi h (t

(6)

(P (P SP + ui − (aP− gzpZP−wsrP)− = 0 h (t (xi

(7)

(Z (Z SZ +ui − (gzpZP−dZ− eZ)− = 0 (t (xi h

(8)

where e a gzp wsr d S N, S P, S Z

is the zooplankton excretion rate, the phytoplankton assimilation rate, the grazing rate from zooplankton on phytoplankton phytoplankton sinking rate (natural death) zooplankton mortality sources and sinks of nutrients, phytoplankton and zooplankton, defined equal to vN*, vP*, vZ*, respectively

C.A.A. Carbonel H., J.L. Valentin / Ecological Modelling 116 (1999) 135–148

140

N*, P*, Z* nutrients, phytoplankton and zooplankton at the interface between the two layers

are defined at rest, At t= 0: Ui = ui = 0; h= H, T= T u,

The assimilation rate a is light (I) and nutrient (N) dependent, being calculated from the Michaelis–Menten equation: a=

6mNI (kn +N)(ki + I)

(9)

where 6m is the constant of optimum photosynthesis rate; kn the half-saturation constant for N uptake; ki the half-saturation constant for light. The light intensity depends on its diurnal variation and is described by: I(t) = − Imax cos for

I(t) =0

    2pt tday

cos

for

cos

 

2pt .B 0 tday

2pt .\ 0 tday

(10)

where tday = 24 h. The grazing by zooplankton on phytoplankton is a process which removes phytoplankton at the rate gzp P and is represented by the quantity gzp which is defined as the fraction of P per unit mass of phytoplankton during each unit time. The quantity gzp depends on algal density and is described by the following relations: gzp = 0 gzp =

when P BPmin

gmax(P− Pmin) (Popt − Pmin)

gzp = gmax

when Pmin 5P 5Popt

when P\ Popt

(11)

where Pmin is the concentration of phytoplankton at which feeding begin and Popt is the concentration of phytoplankton at which grazing rate is maximum (gzp =gmax). The quantities S N, S P, S Z parameterize the transfer of biological concentrations between the layers in function of the upwelling velocity v.

2.3. Initial and boundary conditions To solve the resulting Eqs. (1) – (3) and Eqs. (6) – (8), we need to establish initial and boundary conditions. The initial conditions of the problem

N= 0, P= 0, Z=0

(12)

On land type boundaries G0, non-slip conditions are prescribed for the velocities and homogeneous conditions for the other variables are assumed, Ui = ui = 0

(13)

(T (N (P (Z = = = =0 (xn (xn (xn (xn

(14)

where the subscript n indicates the component normal to the open boundary. On ocean type boundaries G1, the following weakly reflective conditions are prescribed in an axis normal to the boundary (xn ), in the following form: ((Un 9 ch) ((Un 9 ch) +c + G=0 (t (xn

(15)

where, c= gsh,

and

G= F9cb

gh u (T tn ((unUn ) − u + 6Un + + o*, 2m¯ (xn r (xn 2

F=

b= −we The factor o* is the Coriolis term in the normal direction. On the open boundary, the in-going characteristic of the presented equations is used. Additionally, for the temperature and biological variables, homogeneous conditions are prescribed, (T (N (P (Z = = = =0 (xn (xn (xn (xn

(16)

2.4. Numerical model The numerical model follows the scheme for opened regions presented by Carbonel (1998) but with the additional biological conservation equations. The numerical solutions are found on a uniform grid system, with all the variables located at the nodes of the grid. The derivatives in the conservation equations are approximated using forward differences in time, with a dissipative interface and centered differences in space.

C.A.A. Carbonel H., J.L. Valentin / Ecological Modelling 116 (1999) 135–148 Table 1 Physical parameters used in the standard run Parameter

Notation

Value

Initial upper layer thickness The entrainment thickness Upper layer density Lower layer density Coriolis parameter Wind draw coefficient Rayleigh friction coefficient Initial upper layer temperature Lower layer temperature Thermal expansion coefficient The entrainment time scale The surface heat flux time scale Grid step Time step

H He ru rl f cw 6 Tu

30 m 30 m 1024 kg m−3 1025 kg m−3 −5.68×10−5 s−1 2.25×10−3 1.8×10−6 s−1 23°C

Tl u

15°C 2.6×10−4 (°C)−1

te ts

1 4 day 40 days

Dx Dt

1/60° 600 s

On the open boundaries, the weakly-reflective conditions are implemented for the evaluation of the upper layer thickness. The velocities are evaluated using laterally finite differences instead of centered differences in space.

Table 2 Biological parameters used in the standard run Parameter

Notation

Value

Lower layer nutrients Lower layer phytoplankton Lower layer zooplankton Maximal surface light intensity Optimum photosynthesis rate Phytoplankton for initial grazing Phytoplankton for max. grazing Half-saturation for NO3 uptake Half -saturation for light Maximum grazing rate Phytoplankton sinking rate Zooplankton death rate Zooplankton excretion rate

Nl Pl Zl lmax

8.00 mM 0.24 mM 0.10 mM 90 ly h−1

6m

0.22 h−1

pmin

0.7 mM

popt

0.96 mM

kn

1.0 mM

ki gmax wsr d e

2.0 ly h−1 0.04 mM−1 h−1 0.02 h−1 0.01 h−1 0.0026 h−1

141

2.5. Parameter choices and forcing functions Unless specified otherwise, the physical and biological parameter values used in the model are those listed in Tables 1 and 2. The indicated values are based on available data and values reported by several authors. The effects of varying some of the parameters are discussed in Section 3. In Table 1, the values of r u, r l, T u, T l, H are representative of observations (Ikeda et al., 1974; Valentin, 1984; Valentin et al., 1987a; IEAPM, 1988). T u is a typical value of SST in the offshore side of Cabo Frio and T l is a typical temperature just below the thermocline. H is consistent with observed depths of colder waters (T B18°C). The value of He is defined equal to H assuming that entrainment is active when the upper layer thickness hB He. The values of the Rayleigh friction coefficient 6 and the wind draw coefficient cw are adequate calibration parameters to describe upwelling events in the coastal area of Cabo Frio (Carbonel 1998). The time scales te and ts were taken from Mc Creary et al. (1991) and Carbonel (1998). In Table 2, the values of the biological parameters were taken from the literature. The values of N l, P l, Z l are representative of observed conditions in the region of Cabo Frio (Valentin, 1984; Valentin et al., 1987b; IEAPM, 1988). The optimum photosynthesis rate was set at 0.22 h − 1 for a first run and at 0.15 h − 1 for a second run, which are equivalent to the values measured by Gonzalez-Rodriguez (1994) at Cabo Frio (15 mgC mgChlora − 1 h − 1, assuming 50–100 carbon/ chlorophyll ratio) and close to the values used by Jamart et al. (1977) in Pacific ocean. Half-saturation for light (ki = 2.0) and nitrate uptake (kn = 1.0) were taken from Lazzaro (1978) and agree with Walsh (1975) at the Peru upwelling. Photoinhibition was not considered, as Walsh (1975) did for ki = 1.8–2.9 ly h − 1. Phytoplankton output from the system was included in the sinking rate (wsr = 0.02 h − 1, Steele and Henderson, 1981) and in the grazing pressure by zooplankton. Maximum grazing rate (gmax = 0.04 mM − 1 h − 1) and phytoplankton abundance for initial (Pmin = 0.7 mM) and optimal (Popt = 0.96 mM) grazing was taken from Wroblewski (1977). He applied these

C.A.A. Carbonel H., J.L. Valentin / Ecological Modelling 116 (1999) 135–148

142

t= 0, the variables of N, P and Z in the upper layer are set equal zero, whereas in the lower layer N l, P l and Z l are prescribed. During the time integration, for wind-favorable upwelling, the injection of deeper biological concentrations into the upper layer takes place, due the terms S N, S P, S Z.

3. Results

3.1. The standard run. basic north-east wind solution Fig. 3. The non-dimensional time dependent function x(t).

values as constant grazing by the herbivorous copepod Pseudocalanus in Oregon upwelling. At Cabo Frio, zooplankton is also composed mainly of larvae and adults of herbivorous copepods (Valentin et al., 1987a). Zooplankton death (d = 0.01 h − 1) and excretion (0.0026 h − 1) rates are from Andersen (1985). The model is forced by a wind field uniform in space and varying in time, such that the wind components are described by W1 =W01x(t),

W2 =W02x(t),

(17)

where W01 and W02 are the amplitude of the wind velocity components and fixed equal to −8 m s − 1. The function x(t) represent the non-dimensional time dependent structure of the wind, defined by x(t) =





2pt 1 1− cos , PW 2

x(t) = 0,

for t BPW for t ]PW

(18)

where PW is the time duration of x(t). Fig. 3 presents the functions x(t) used in the experiments of Section 3. The diurnal light variability was included in Eq. (10). With this choice, the light intensity I is zero from 18:00–06:00, it reaches a maximum value Imax at noon. The coupled model is integrated from a state of rest (no motion) for a period of 10 days. At time

In this part, we describe the standard run in detail. This solution uses a transient wind with PW = 5 days and the parameter values listed in Tables 1 and 2. The time response of the model is evaluated during 10 days. Fig. 4 illustrates the response of the coupled non-linear model to the transient uniform wind forcing, showing the current system and the SST distribution at day 2, 4 and 6. A band of upwelled cool water is generated along the coast of Cabo Frio, reaching colder temperatures along the zonal coastline. Warmest temperatures occur well offshore approaching the initial value. During the spin-up stage of upwelling the band increases in width, decreasing slowly during the spin-down stage (after day 4). The velocity field shows a current flowing along the coast. At the offshore side of the zonal coastline, an eastward flow component is noted, only during the first days. Fig. 5 presents the spatial configurations of the nutrients, phytoplankton and zooplankton fields on different days. At day 2 the nutrient field is concentrated in a band along the zonal coastline. At day 4 the nutrients fields is more intense and reaches the meridional coastline. The phytoplankton field is perceptible at day 4 with a similar pattern to the nutrient field. At day 6 the nutrient concentrations are minimal when the phytoplankton bloom is practically finished (the maximum intensity of the phytoplankton occurred at day 5) and the maximal concentrations are located offshore, whereas the zooplankton field has an important intensity along the zonal coastline. The results presented for the standard run indicate

C.A.A. Carbonel H., J.L. Valentin / Ecological Modelling 116 (1999) 135–148

143

Fig. 4. Upper layer currents and SST fields for the standard run at day 2 and 4. The SST field (in °C) shows an upwelling band along the coast with colder waters at the zonal coastline. The velocity field shows a coastal jet well developed along the coast at day 4. The arrows are normalized to the maximum velocity in each case (0.19 m s − 1 at day 2 and 0.11 m s − 1 at day 4). Contour interval for the SST is 0.5°C.

that biological concentrations show a spatial distribution along the zonal coastline similar to the SST field presented in Fig. 4, during the spin-up stage of the biological variables. Additionally, biological concentration centers are noted in front of Saquarema, an influence of the upwelling dynamic in this coastal area. These centers could be stronger when wind patch structures are observed in this region, as occurred with the plumes of SST (Carbonel, 1998), which intensify at about 30 nautical miles in southwest direction (Ikeda et al., 1974).

In Fig. 6 the time-dependent variation of the biological components in two control points of the model located at the coast (Saquarema and Cabo Frio) is shown. The nutrient concentrations in the upper layer increases from zero reaching a maximum at the beginning of day 4 and then decreasing to zero at day 5, indicating the strong dependence of nutrient source (function of the upwelling velocities) resulting from wind forcing. The phytoplankton concentrations are perceptible only on day 3 reaching a maximum value of 4.7 mM at day 5 at the Saquarema point and drop-

144

C.A.A. Carbonel H., J.L. Valentin / Ecological Modelling 116 (1999) 135–148

Fig. 5. Biological fields (in mM) for the standard run at different days. At day 2 only the nutrients have perceptible concentrations. At day 4, the nutrient and phytoplankton fields are well developed. At day 6, the phytoplankton is weak and shows an interesting maximum offshore, whereas the zooplankton field has a marked presence along the zonal coastline and the nutrient concentrations are not perceptible. The contour interval for the concentrations is 0.5 mM except for the phytoplankton at day 6 (0.2 mM).

ping to values in the range of 0.5 – 0.8 mM on day 6. At the Cabo Frio point, the phytoplankton reaches the maximum value of 2.5 mM on day 5

and is equivalent to 5 mg l − 1chlorophyll (1 mg l − 1 chlorophyll is equal to 0.5 mM nitrogen according to Radach and Meier Raimer, 1975), similar to

C.A.A. Carbonel H., J.L. Valentin / Ecological Modelling 116 (1999) 135–148

145

the maximum value observed in January 1976 very near of Cabo Frio (Fig. 1b). The maximum of phytoplankton (phytoplankton bloom) coincides with the rapidly drops of nutrients and the beginning of zooplankton increase which reaches a maximum only on day 6, indicating that grazing by zooplankton play an important role in limiting the persistence of the phytoplankton bloom. The time-dependent evolution of the phytoplankton is marked by the presence of diurnal oscillations, due the influence of diurnal light variation in the phytoplankton assimilation rate a. The diurnal oscillations are stronger during the spin-up of the bloom and have amplitudes of around 1 mM. The model forced by NE winds shows results qualitatively compatible with observed colder waters and important concentrations of biological variables along the zonal coastline of Cabo Frio (Ikeda et al., 1974; Valentin, 1984).

3.2. The influence of wind forcing duration (PW) To investigate how the time duration of wind events influences the biological response, we carried out two calculations with a NE wind with PW equal to 4 and 6 days. Increasing PW to 6 days (see Fig. 7, lower panel), the explosion of nutrient concentrations in the upper layer has a greater Fig. 7. The influence of the time dependence of wind forcing in the biological variables (in mM) at control points, when the duration of the wind forcing PW is changed. Upper panel: PW =4 days. Lower panel PW =6 days.

Fig. 6. Time variability of the biological variables (in mM) at control points of the model (Cabo Frio and Saquarema) for the standard run. The explosion of biological variables at different times is the main characteristic of the solution.

amplitude at the control points maintaining a stock (persistence) during a longer time (drops to zero at day 7). The peaks of phytoplankton are approximately 6% larger as in the standard run. In this experiment the zooplankton concentrations are more important reaching amplitudes similar to those of phytoplankton. When PW is decreased to 4 days (see Fig. 7, upper panel), the time response of nutrient and phytoplankton concentrations in the upper layer had reduced amplitude. The phytoplankton peak reaches a maximum of 2.7 mM but the presence of the bloom is maintained for a longer time due to the lower concentration of zooplankton. The decay of

146

C.A.A. Carbonel H., J.L. Valentin / Ecological Modelling 116 (1999) 135–148

the phytoplankton bloom then results from the sinking rate (natural death). It is interesting to note that, whereas the instant of maximum of nutrients concentration changes due variation of PW, the peak of phytoplankton and zooplankton remains unchanged, indicating that it depends on the beginning of favorable winds. Another feature demonstrated in the results is the persistence (duration) of the nutrient explosion, which is associated with the increase of PW. Additionally, it is noted that the diurnal oscillation of phytoplankton is dominant when PW decreases.

3.3. Sensibility of biological parameters The influence and sensitivity of the biological parameters in the time dependent solutions is studied by performing some runs changing the values of optimum photosynthesis rate 6m and maximum grazing rate gmax by zooplankton in the standard run.

3.3.1. Optimum photosynthesis rate 6m Changes of the optimum photosynthesis rate coefficient 6m altered the solution obtained in the standard run in expected ways. Fig. 8, shows the response of the biological variables when 6m decreases to 0.15. In this case, the phytoplankton

Fig. 8. Time evolution of biological variables (in mM) at control points, when the optimum photosynthesis rate vm is reduced to 0.15 h − 1.

Fig. 9. Time evolution of biological variables (in mM) at control points, when the maximum grazing rate gmax is reduced to 0.02 mM − 1 h − 1.

assimilation rate a is reduced and the evolution time of the nutrient concentration changes with the consumption of nutrients dropping slowly to zero only at day 7, whereas the maximum amplitude of nutrients remains practically unchanged (at day 4). The phytoplankton bloom decreases its amplitude ( 40%), reaching the maximum value at day 7 (2 days later than occurred in the standard run). In a local model of primary production the value of 6m = 0.15 is enough to generate a peak of phytoplankton (Valentin, 1992), but here an important dispersion of concentrations is produced due the presence of a current system with a stronger flow divergence near the zonal coastline.

3.3.2. Maximum grazing rate gmax Here the maximum grazing rate gmax = 0.04 mM − 1 h − 1, used in the standard run, is reduced by 50% to investigate its influence. The time response of the biological variables is shown in Fig. 9. In this case, the time evolution of nutrients keeps practically unchanged, but the phytoplankton bloom is a little more intense with a more longer persistence due the weakness of the zooplankton concentration (there is less grazing), which in this case is minimal. This time response of the biological variables is very similar to that occurred when PW = 4 days (Fig. 7, upper panel),

C.A.A. Carbonel H., J.L. Valentin / Ecological Modelling 116 (1999) 135–148

where the phytoplankton bloom decay was due to natural death of phytoplankton and not to grazing of zooplankton. Zooplankton grazing is an important parameter in the biological dynamics. These results verified the fact that maximum grazing rate gmax influences directly the strength and duration of the phytoplankton bloom.

4. Summary and conclusions In this paper a mathematical model based on conservation equations of physical and biological variables is proposed to study the current system and primary biomass of upwelling regions. The model is applied to the coast of Cabo Frio (Brazil) with the purpose of investigating the phenomenon of phytoplankton blooms observed in this upwelling region. The physical–biological model is a gravity reduced 112 layer model, composed by a set of conservation equations of a dynamic upper layer, for momentum, mass, heat, nutrients, phytoplankton and zooplankton. These equations, are solved numerically in space and time by the finite difference method. The physical model is forced by upwelling favorable winds acting at the sea surface to generate areas of divergence at the coast (upwelling centers), accompanied by reduction of upper layer thickness and cooling of the SST. The biological model which has three-components (nutrients, phytoplankton and zooplankton) is forced by the primary biomass concentration upwelled into the upper layer by upwelling favorable winds and by the light producing the photosynthetic activity of phytoplankton. Numerical solutions forced by transient winds are arranged initially for a standard run developing an offshore Ekman drift along the coast and causing coastal upwelling. A westward coastal jet and eastward velocity components in the offshore side are observed. The spatial distribution of phytoplankton has a similar configuration as those of nutrient, zooplankton and temperature along the zonal coastline of

147

Cabo Frio. During the spin-up stage of winds, the primary biomass concentrations are injected into the upper layer due to upwelling. A maximum of nutrient concentrations takes place (generated by the nutrients injection in the upper layer), decayed by the increasing assimilation by phytoplankton and by reduction of nutrient supply (spin-down stage of winds). One day after, the phytoplankton bloom occurs with amplitudes similar to those observed in Cabo Frio (IEAPM 1988; Valentin and Coutinho, 1990). This bloom is very short (1 day), due the increase of grazing by zooplankton and compare qualitatively well with the existing time-dependent observations of phytoplankton bloom in the region of Cabo Frio, suggesting that the physical–biological model presented here contains much of this phenomenon. The peak amplitude of nutrient explosions depends mainly on the wind driven upwelling. The change of transient wind duration affects primarily the nutrient stock in the upper layer, changing the instant of the maximum of nutrient concentration, but remaining unchanged at the instant of maximum of phytoplankton. The intensity of the phytoplankton blooms depends on the factors influencing phytoplankton growth and grazing by zooplankton. In the experiments, the reduction of the constant of optimum photosynthesis rate affects the bloom of phytoplankton reducing its amplitude and delaying in time. The maximum grazing rate is an essential factor in limiting the persistence and amplitude of the phytoplankton bloom. When it is reduced, the amplitude and duration of the phytoplankton bloom increase, whereas the zooplankton stock decreases strongly. The instant of maximum phytoplankton depends directly on the beginning of the upwelling-favorable winds, when the duration of the wind forcing changes. This study has given us insights into the space and time mechanisms controlling the phytoplankton bloom in the coastal region of Cabo Frio. With this knowledge we can now explore the model testing others effects and processes which could be important in particular cases.

148

C.A.A. Carbonel H., J.L. Valentin / Ecological Modelling 116 (1999) 135–148

References Andersen, V., 1985. Mode´lisation d’ eosyste`mes pe´lagiques. ´ tude de processus. The`se de Doctorat d’ Etat, Universite´ E P&M. Curie, Paris VI, 135 pp. Anderson, D., 1984. An advective mixed layer model with application to the diurnal cycle of the low-level East African Jet. Tellus 36, 278–291. Anderson, D., Mc Creary, J.P., 1985. Slowly propagating disturbances in a coupled ocean-atmosphere model. J. Atmos. Sci. 42, 615 – 629. Carbonel, C., 1998. Modelling of upwelling in the coastal Area de Cabo Frio (Rio de Janeiro, Brazil). Braz. J. Oceanogr. 46 (1), 1 – 17. Gonzalez-Rodriguez, E., 1994. Yearly variation in primary productivity of marine phytoplankton from Cabo Frio (RJ, Brazil) region. Hydrobiologia 294, 145–156. Haney, R.L., 1971. Surface thermal boundary condition for ocean circulation models. J. Phys. Oceanogr. 1 (4), 241– 248. IEAPM, 1988. Comissa˜o Cabo Frio VIII. Doc. Tec. 208. Instituto de Estudos do Mar Alte. Paulo Moreira. Arraial do Cabo, Rio de Janeiro, 21 pp. Ikeda, Y., Miranda, L.B., Rock, N.J., 1974. Observations on stages of upwelling in the Region of Cabo Frio (Brazil). Bol. Inst. Oceanogr. Sao Paulo 23, 33–46. Jamart, B.M., Winter, D.F., Banse, K., Anderson, G.C., Lam, R.K., 1977. A theoretical study of phytoplankton growth and nutrient distribution in the Pacific Ocean off the northeastern US coast. Deep-Sea Res. 24 (8), 753–773. Jones, R., Henderson, E.W., 1987. The dynamics of energy transfer in marine food chains. In: Paine, A., Gulland, J.A., Brink, K.H. (Eds.), Benguela and Comparable Ecosystems, S. Afr. J. Mar. Sci. 5, 447–465. Lavoie, R.L., 1972. A mesoscale numerical model of lake-effect storms. J. Atmos. Sci. 29, 1025–1040. Lazzaro X., 1978. Mode´lisation sur le the`me production primaire en mer Ligure au printemps. DEA Oceanogr. Biol., Univ. P and M Curie, Paris VI. Mascarenhas, A.S., Miranda, L.B., Rock, N.J., 1971. A study of oceanographic conditions in the region of Cabo Frio. In: Costlow Jr. (Ed.), Fertility of the Sea. Gordon and Breach, New York, pp. 285–308. Mc Creary, J. P., Kundu, P., 1988. A numerical investigation of the Somali Current during the southwest monsoon. J. Mar. Res. 46, 25 – 58. Mc Creary, J.P., Fukamachi, Y., Kundu, P., 1991. A numerical investigation of jets and eddies near an eastern ocean boundary. J. Gheophys. Res. 96, 2515–2534. Mc Creary, J. P., Kohler, K., Hood, R.R., Olson, D.B., 1996. A four-component ecosystem model of biological activity in the Arabian Sea. Prog. Oceanogr. 37, 193–240.

.

O’Brien, J.J., Wroblewski, J.S., 1973. A simulation of the mesoscale distribution of the lower marine trophic levels off West Florida. Invest. Pesq. 37, 193 – 244. Patten, B.C., 1968. Mathematical models of plankton production. Int. Ver. Ges. Hydrobiol. 53, 357 – 408. Radach, G., Meier Raimer, E., 1975. The vertical structure of phytoplankton growth dynamics. A mathematical model. Mem. Soc. Roy. Sci. Liege 6, 113 – 146. Riley, G. A., 1965. A mathematical model of regional variations in plankton. Limnol. Ocean. 10, 202 – 215. Ripa, P., 1993. Conservation laws for primitive equations models with inhomogeneous layers. Geophys. Astrophys. Fluid Dyn. 70, 85 – 111. Sekine, M., Nakanishi, H., Ukita, M., 1991. A shallow-sea ecological model using an object-oriented programming language. Ecol. Model. 57, 221 – 236. Schop, P., Cane, M., 1983. On equatorial dynamics, mixed layer physics and a sea surface temperature. J. Phys. Oceanogr. 13, 917 – 935. Steele, J.H., Henderson, E.W., 1981. A simple plankton model. Am. Nat. 117, 676 – 691. Valentin, J.L., 1984. Analyse des parame`tres hydrbiologiques dans la remonte´e de Cabo Frio (Bre´sil). Mar. Biol. 82, 259 – 276. Valentin, J.L., Lins da Silva, N.M., Bastos, C.T.B.N., 1985. Les diatome´es dans l’upwelling de Cabo Frio (Bre´sil): liste d’espe`cies et e´tude e´cologique. J. Plankton Res. 7, 313 – 337. Valentin, J.L., Andre, D.L., Jacobs, S.A., 1987a. Hydrobiology in the Cabo Frio (Brazil) upwelling: two-dimensional structure and variability during a wind cycle. Contin. Shelf Res. 7 (1), 77 – 88. Valentin, J.L., Monteiro-Ribas, W.M., Mureb, M.A., Pessotti, E., 1987b. Sur quelques zooplanctontes abondants dans l’upwelling de Cabo Frio (Bre´sil). J. Plank. Res. 9 (6), 1195 – 1216. Valentin, J.L., Coutinho, R., 1990. Modelling maximum chlorophyll in the Cabo Frio (Brazil) upwelling: a preliminary approach. Ecol. Model. 52, 103 – 113. Valentin, J.L., 1992. Modeling of the vertical distribution of marine primary biomass in the Cabo Frio upwelling region. Cienc. Cult. 44 (2/3), 178 – 183. Walsh, J.J., 1975. A spatial simulation model of the upwelling ecosystem. Deep-Sea Res. Oceanogr. Abstr. 22, 201 – 236. Walsh, J.J., Dugdale, R.C., 1971. A simulation model of the nitrogen flow in the Peruvian upwelling system. Inv. Pesq. 35, 309 – 330. Wroblewski, J.S., 1977. A model of phytoplankton plume formation during variable Oregon upwelling. J. Mar. Res. 35, 357 – 394. Wroblewski, J.S., 1980. A simulation of the distribution of Acartia clausi during Oregon upwelling, August 1973. J. Plank. Res. 2, 43 – 68.