Serviços Personalizados
Journal
Artigo
Indicadores
- Citado por SciELO
- Acessos
Links relacionados
- Citado por Google
- Similares em SciELO
- Similares em Google
Compartilhar
Journal of the Brazilian Society of Mechanical Sciences and Engineering
versão impressa ISSN 1678-5878
J. Braz. Soc. Mech. Sci. & Eng. vol.32 no.1 Rio de Janeiro jan./mar. 2010
http://dx.doi.org/10.1590/S1678-58782010000100007
TECHNICAL PAPERS
A finite element model for the ocean circulation driven by wind and atmospheric heat flux
Carlos A. A. Carbonel H^{I}; Augusto C. N. R. Galeão^{II}
^{I}carbonel@lncc.br
^{II}acng@lncc.br, Lab. Nacional de Computação Científica - LNCC 25651-075 Petropolis, RJ, Brazil
ABSTRACT
A finite element model for solving ocean circulation forced by winds and atmospheric heat fluxes is presented. It is vertically integrated 1 1/2 layer model which solves the motion and continuity hydrodynamic equations coupled with the advection-diffusion transport equation for temperature. A space-time Petrov-Galerkin formulation is used to minimize the undesirable numerical spurious oscillation effects of unresolved boundary layers solution of classical Galerkin method. The model is employed to simulate the circulation of the eastern South Pacific Ocean represented by a mesh covering the area between the equator and 30ºS and between the 70ºW to 100ºW. Monthly climatological data are used to determine the wind and heat fluxes forcing functions of the model. The model simulates the main features of observed sea surface temperature (SST) pattern during the summer months showing the upwelling along the coastal boundary with currents oriented northwestward and the presence of southward flows and warm water intrusion in offshore side. The calculated SST fields are compared with the mean observed SST showing that the coastal processes and the interaction with the equatorial band are physically better resolved.
Keywords: ocean circulation model, stabilized finite element formulation, atmospheric forcing
Introduction
The finite element method using unstructured computational grids of triangular elements of variable size allows high resolution of local features in irregular domains. Consistent and stable weak formulation method lead to reliable accurate numerical schemes to solve the governing equations of the problem, allowing to rigorous integration of local and far-field interactions of the processes included in the mathematical formulation of the real problem.
The application of finite element methods in oceanography requires, as well, non-uniform refined mesh to take into account the irregular coastline shape. The spatial variability in grid spacing gives rise to unphysical wave scattering, which (under certain circumstances) could be of minor importance for steady state engineering problems, and to hydrodynamic problems at relatively short integration time; but might become of crucial importance as meso to large-scale ocean circulation modeling is focused. Besides, in ocean modeling, the numerical treatment of advection is still a very important issue, because the approximation of advection-dominated problems present numerical difficulties related to the lack of stability of the discretized operator.
There are several problems using finite element models for transient and convection dominated problems. When advection dominates diffusion, the classical Galerkin finite element methods generate unstable approximations, which usually exhibit spurious oscillations. Additionally, the models based on classical Galerkin formulations have shown some misbehavior in time-dependent propagation problems. Only for small values of the time integration step Δt it is possible to obtain numerical schemes that retain the high accuracy of the element-based spatial discretizations. For quite modest values of Δt, the accuracy and phase-propagation properties of the Galerkin formulation when applied to parabolic problems are lost, and the stability range is reduced in comparison with finite difference schemes (Donea and Quartapelle, 1992). In the case of hydrodynamic problems, it is more representative the features of the solution related to the accuracy/stability dependence on the dimensionless Courant number, instead of the size of the time step Δt. These disadvantages have motivated, in the last years, the development of alternatively finite element formulations such as Taylor-Galerkin, least-squares Galerkin and Petrov-Galerkin methods to improve the time-accuracy approximations of time evolution problems, particularly for those for which advection is the dominant behavior.
Finite element formulations based on the Streamline Upwind Petrov-Galerkin (SUPG) formulation and similar operators can solve the problem of unphysical wave scattering and unresolved boundary layer. These formulations were initially proposed for advection-diffusion and compressible flow problems (e.g., Brooks and Hughes, 1982; Shakib, 1988; Almeida et al., 1996). Formulations using a symmetric form of the shallow water wave equations were reported in last years (e.g., Bova and Carey, 1995; Carbonel et al., 2000). Also, for hydrodynamic baroclinic circulation, an approach was also presented (Carbonel and Galeão, 2006) focusing the wind-driven dynamic response in open and limited areas.
The ocean models range in dynamical complexity from simple linear surface-layer model to sophisticated, non-linear, continuously stratified systems. The physical inhomogeneous layers models were introduced by Lavoie (1972) and later important contributions using the same idea were published (e.g., Schop and Cane, 1983; Anderson, 1984; Anderson and McCreary, 1985; McCreary and Kundu, 1988; Ripa, 1993). The 1 ½-layer model presented here is the simplest possibility to represent physical processes in a real stratified ocean. The model has one active layer overlaying a quiescent, deep ocean, where pressure gradients are assumed to vanish. The water is allowed to move across the interface between the layers. The system is thermodynamically active so that the temperature in the upper layer changes its response according to the surface heat flux, the interface cooling flux, and the horizontal advection.
In this paper, a Petrov-Galerkin formulation is developed for the coupled system of equations of motion, continuity and transport of heat energy in the 1 ½-layer domain. These equations apply only to the active oceanic surface layer, while the deep layer is assumed to be motionless, arbitrarily deep, and separated from the upper layer by a density discontinuity. The formulation uses a stabilizing operator in space and time to improve the classical Galerkin approach. This methodology allows accurate solutions, precluding the use of artificial filters, mass diffusion and dampers, commonly used in several popular models. The treatment of the open boundary conditions is a very critical issue in a limited-area domain. For the model developed here, the treatment of the open boundaries is based on the formulation of weakly reflective open boundary conditions derived from the characteristic equations (Verboom, 1982; Verboom et al., 1983). The effectiveness of those kinds of boundary condition types for coastal circulation problems has been previously demonstrated (e.g., Carbonel, 2003; Carbonel and Galeão, 2006).
In the eastern ocean boundaries, the upwelling-favorable winds force offshore drift at the ocean surface, leading to upwelling of cold subsurface water at the coast, the formation and the offshore advection of coastal fronts, the generation of alongshore currents. Other commonly observed phenomena are the appearance of an equatorward surface jet, and a poleward undercurrent. In the present paper, the ability of the model to predict the above mentioned behavior is tested for the eastern South Pacific Ocean, focusing the simulation of the upper layer circulation and SST dispersion during the summer months. During this period the coastal upwelling remains weakly, but warming waters from the north invade the region, mainly in the offshore side. The simultaneously occurrence of colder waters near the coast and a southward warmer intrusion offshore configures a scenario with strong SST gradients, which is a challenger modeling ocean process in a limited area domain. The model is forced by monthly mean wind and solar radiation fields, but it develops its own fluxes of sensible and latent heat.
Nomenclature
A_{H}, A_{V} | = | horizontal and vertical eddy viscosity respectively, m^{2}/s |
b | = | b=, m/s |
c | = | wave celerity, m/s |
c_{0} | = | referential wave celerity, m/s |
C_{S}, CL | = | Stanton and Dalton parameters for the sensible heat and latent heat respectively |
c_{pa} | = | specific heat of the air, J/(kg^{ o}C) |
c_{p} | = | specific heat of the surface water, J/(kg^{ o}C) |
d | = | d=2c m/s |
f | = | Coriolis parameter, sec^{-1} |
H | = | initial upper layer thickness, m |
H_{e} | = | depth scale of the entrainment process, m |
h | = | upper layer thickness, m |
L | = | latent heat of evaporation, J/kg |
N | = | number of finite element sub-domains |
P^{k} | = | set of piecewise polynomials of degree < k |
Q | = | net heat flux across the surface, Watt/m^{2} |
Q_{I} | = | heat flux across the interface, Watt/m^{2} |
Q_{s} | = | net short wave radiation, Watt/m^{2} |
Q_{lw} | = | outgoing long wave radiation, Watt/m^{2} |
Q_{h} | = | sensible heat flux, Watt/m^{2} |
Q_{e} | = | latent heat flux, Watt/m^{2} |
q | = | saturation specific humidity |
R^{2} | = | the space of (x,y) points |
R^{±} | = | Riemann invariant |
S_{n} | = | space-time finite element "slab" |
T | = | sea surface temperature (SST), ^{o}C |
T^{u} | = | referential initial upper layer temperature, ^{o}C |
T^{l} | = | lower layer temperature, ^{o}C |
t_{e} | = | time scales of the entrainment process, sec |
u, v | = | velocity components in the upper layer, m/s |
u_{n} | = | velocity component normal to the open boundary, m/s |
= | set of admissible trial functions | |
= | finite element subspace of weighting functions | |
V^{h} | = | trial function |
= | weighting function | |
w_{ed} | = | entrainment-detrainment velocity, m/s |
ω | = | divergence velocity, m/s |
x_{n} | = | axis normal to the open boundary |
x, y | = | Cartesian coordinates |
Greek Symbols | ||
α | = | coefficient of thermal expansion, ^{o}C^{-1} |
Γ_{L} | = | land type boundary |
Γ_{O} | = | open ocean boundary. |
π^{h , }^{Δ}^{t} | = | pace time partition |
Ω_{e} | = | finite element sub-domain |
Ω | = | spatial domain |
Γ_{e} | = | finite element boundary sub-domain |
= | matrix of intrinsic time scale | |
ϑ_{T} | = | thermal diffusivity, m^{2}/s |
Φ | = | the empty set |
ρ_{air}, ρ^{u}, ρ^{l } | = | densities of air, upper layer and lower layer respectively (kg/m^{3}) |
τ_{x },τ_{y} | =_{} | wind stress components, N/m^{2} |
The Ocean Model
The Governing Equations
We consider an incompressible fluid, which satisfies the Boussinesq approximation in a coastal ocean. An apparent temperature (Foffonof, 1962) is utilized and constant eddy coefficients are used to represent the turbulent exchange processes. Under these assumptions and defining a set of Cartesian coordinates (x, y, z) oriented positively to the East, North and upwards directions respectively, the governing equations describing the hydro-thermodynamical processes can be written as:
In the preceding equations (u, v, w) are the velocity components in (x, y, z) directions respectively, t is the time and p is the pressure. ρ is the density and T is the apparent temperature with constant reference values ρ^{ref} and T^{ref} respectively. The Coriolis parameter that depends on the geographical coordinates is represented by f and g is the gravity acceleration. A_{H}, A_{V} and K_{H}, K_{V} are the constant horizontal and vertical eddy coefficients for momentum and temperature respectively and α is the thermal expansion coefficient.
The turbulent eddy viscosities and diffusivity constants are difficult to calculate accurately for oceanic flows, but some referential values could be considered for applications purposes (Stewart, 2004).
In the present paper we are mainly interested in the dynamical description of the upper layer ocean of density ρ^{u}(x,y,t), thickness h (x,y,t) and temperature T(x,y,t), overlying a deep inert layer of constant density ρ^{l} and temperature T^{l} where the horizontal pressure gradients are set to zero. The lower layer is arbitrarily deep and its thickness H^{l}(x,y,t) is measured from a deep reference level up to the interface. The interface between these two layers represents the thermocline (see Figure 1).
In a layered medium, boundary and interface conditions should be specified to close the system of equations governing the model behavior. For our particular model, let us introduce an interface level function η, such that η = H^{l} + h identifies the free surface, and η = H^{l} defines the common interface between the active upper layer and the inert lower layer. Under those assumptions the following set of conditions must be fulfilled:
where ω is a vertical velocity due to mass transfer that should be zero at the free surface and equal to the entrainment-detrainment velocity w_{ed} at the interface between the upper and lower layers. At the free surface the heat source Θ = Q (the heat flux between the ocean and the atmosphere), and in the internal interface Θ = Q_{I} (the gain or loss of heat across this surface). The shear stress Υ_{x} , Υ_{y} are equal to the wind stresses (τ_{x} , τ_{y}) at the surface; and are assumed to be zero in the lower interface according to the proposed 1 1/2 layer model definition.
Integrating Eqs. (1), (2), (4), and (5), from the interface to the free surface, using the Leibnitz rule, and considering the boundary and kinematics conditions, we obtain
Equations (9) to (12) describe the upper layer dynamics of the proposed 1 ½ layer model in a space time domain Ω x (0,T), where the spatial domain Ω R^{2} with land type boundary Γ_{L}_{ } and an ocean type boundary Γ_{O}_{ }such that (empty set) and , the boundary of Ω.
For simplicity, the same previous notation u, v, T is now being used to represent the mean horizontal velocity components and the mean temperature in this layer. For a generic variable φ(x,y,z,t), its mean value is:
The vertical velocity w_{ed} which appears in Eq. (11), represents the entrainment-detrainment velocity given by
where H_{ed} and t_{ed}_{ }are the depth and time scale of the entrainment-detrainment process. This definition is similar to the entrainment equation proposed by McCreary and Kundu (1988). In the present case, the velocity w_{ed} avoids the interface from surfacing when (H_{ed} -h) > 0 and also avoids the interface from deepens excessively when (H_{ed} -h) < 0. It is also possible to define a similar expression using different depths and time scales for the entrainment and detrainment (e.g., Roed and Shi, 1999).
According to Eq. (3), the pressure distributions in the upper layer p^{u} and lower layer p^{l} are given by
Using the anzatz
corresponding to the specification of ρ^{ref} = ρ^{l} and T^{ref} = T^{l} in Eq. (6), and eliminating the thickness H^{l} we arrived at
where , . Here, the influence of the salinity is not included, which does not mean a loss of generality, because an apparent temperature could be estimated modifying the constant α (Hantel, 1971; Fofonoff, 1962). The presented model is classified as an Inhomogeneous Layer Primitive Equations Model (ILPEN) valid for low-frequency dynamics (Ripa, 1993).
Forcing Functions, Boundary Conditions and Initial State
To completely define the proposed ocean model we need to characterize the applied external sources of energy, the time dependent boundary conditions and the initial state from which the relevant hydrothermodynamical response will be obtained.
The hydrodynamic behavior of the model is forced by the wind stresses evaluated from the wind velocity fields using the quadratic relations
where ρ_{air} is the air density and is the drag coefficient (Enriquez and Friehe, 1997). W_{x}, W_{y}_{ }are the wind velocity components and |W| is the wind module.
The surface heat flux Q is the heat input at the ocean surface defined as
where Q_{sw} is the net short wave radiation; Q_{lw} the resultant outgoing longwave radiation; Q_{S} the sensible heat and Q_{L}_{ }the latent heat. The sensible heat is calculated from the bulk formula
where c_{pa} is the specific heat of the air, C_{S} an empirical coefficient, and T and T_{air} are the sea surface temperature and air temperature respectively. The latent heat is estimated from
where L is the latent heat of vaporization, C_{L }an empirical coefficient, and q, q_{air} are the saturation specific humidity of the sea surface and air respectively. For a constant lower layer temperature T^{l} the heat flux across the interface is given by:
showing that the gain or loss of heat across this interface, depends on the dynamical convergence or divergence of the flows ω in the upper layer, as expressed by the parameter (Carbonel, 2003). Concerning the boundary conditions it is assumed that land type boundaries give rise to the classical non-slip conditions:
and for the upper thickness h and for the temperature T homogeneous conditions are prescribed.
On the sea side boundaries, we prescribe a weakly reflective boundary condition, based on the characteristic method. This is done assuming that on a normal direction (x_{n}) to the open boundary
where the Riemann invariant R^{±}= u_{n} ± 2c is valid along the characteristics dx_{n} /dt = u_{n} ± c, normal to the open boundary. u_{n} represents the velocity component normal to this boundary and . The weakly reflective conditions on the open boundary are defined by the in-going characteristic of the presented equations. For the temperature at the open boundaries, homogeneous condition is assumed. Finally, the initial state requires the knowledge of a set of: initial velocity components
( u^{o}(x,y), v^{o}(x,y) ); initial upper layer thickness h^{o}(x,y) and an initial temperature field T^{o}(x,y), such that for
A Quasi-Symmetric Form of the Governing Equations
Symmetric forms are always desirable because they possess better stability properties, an important feature to be used in the design of approximate solutions of advection dominated problems. For the particular case of fluid flow problems these aspects are pointed out in Shakib (1988), Almeida et al., (1996), Bova and Carey, 1995, Ribeiro et al. ( 1996). For the hydrothemodynamical problem treated in this paper a complete symmetrization is not achieved. Nevertheless, the hydrodynamic part of the equations governing the ocean circulation could be symmetrized. Furthermore, the whole set of unknown variables retains the same physical dimensions, which means that the transformed equations are adequately scaled.
This quasi-symmetrization is done introducing a new set of kinematic variables defined as:
which have the dimensions of velocity, and are function of the upper layer thickness h and the upper layer temperature T respectively. The parameter is a referential celerity function of the referential initial upper layer thickness H and T^{u} is the referential initial temperature in the upper layer. After this change of variables, the set of equations for the unknown vector V = (u,v,d,b) could be written as:
where
and . The resulting equation system is an extension of the symmetric form of the hydrodynamic part of the model.
Space Time Petrov-Galerkin Method (STPG)
In this section we present the space-time Petrov-Galerkin finite element approximation for the continuous hydrothermodynamical ocean model described by Eq. (29).
We will assume that: the spatial region of interest is mathematically represented by an open bounded set Ω R^{2}, with boundary Γ; and the time period of observation is the interval R^{+}. We now consider that for n = 0,1,2...M an ordered partition of time levels
is constructed. Denoting by I_{n} = (t_{n+1} - t_{n}) the nth interval, we will say that, for each n, the space-time domain of interest is the "slab" , with boundary . Then for n = 1,2...N we will define a finite element triangulation of Ω such that
In accordance with the previous definitions, we will assume that the finite element subspace is the set of continuous piecewise polynomials in , which may be discontinuous across the slab interface at t_{n} , i.e.:
where is the set of polynomials of degree less than or equal to k. For a prescribed boundary condition on , a general trial function is then an element of , where:
For the present numerical ocean model, continuous interpolation in space and time will be adopted. Under this assumption, we will say that the approximate solution of Eq. (29) is the particular element V^{h} , if it satisfies
where
represents the residual associated to (29)
is a space-time operator, and at the initial time
where V_{0} (Ω) is the initial condition in Ω.
The first integral appearing in Eq. (34) reproduces the classical Galerkin formulation whereas the integral under the summation index comes from the space-time Petrov-Galerkin formulation. is the matrix of intrinsic time scales, containing free parameters, normally known as upwind functions (Hughes and Mallet, 1986). This matrix is symmetric and positive semidefinite and could be interpreted as an array of "scaling" parameters adequately chosen to "adjust" the weighting function effect. An optimal choice for is an open problem since convergence analysis and error estimates can provide some design conditions for , but are insufficient to determine an unique definition. Some procedures commonly used to estimate the intrinsic time scale matrix (e.g., Shakib, 1988; Bora and Carey, 1995; Ribeiro et al., 1996) do not give the best parameters for general applications in hydrodynamics and water waves. For these cases, some specific design procedures must be used to improve the physical and numerical mechanisms involved in practical problems.
For complex system of equations, the definition is complicated due to the presence of multiple wave components. For the particular problem focused in this paper, we propose
where is the representative length of the triangular spatial element , and are the local coordinates. Since is symmetric, it may be computed in terms of its eigenvalues and eigenvectors, from the expansion
where e_{k} are the eigenvectors and
are the eigenvalues. The intrinsic time scales components of the matrix defined in Eq.(38) will always be greater than Δt. The resulting system of algebraic equations is solved using the direct method of Gauss.
Numerical Results
A numerical experiment was conducted to evaluate the performance of the proposed discrete model to predict the main hydrothermodynamical features in a limited open region of the eastern South Pacific Ocean. The most familiar characteristic of eastern boundary ocean dynamics is the importance of equatorward surface currents and coastal upwelling resulting from a predominantly equatorward wind stress. The coexistence of countercurrents and undercurrents contributes to an even more complex dynamical behavior (Yoshida, 1967).
The region of interest (Fig. 2a) extends in the north-south direction from equator to the 30ºS and in the east-west direction from the 70ºW to 100ºW. This region was discretized by an unstructured mesh of 5684 surface triangular elements and 2994 nodes, as illustrated in Fig. 2b. Near the coastline, a finer mesh was employed which gradually becomes coarser and coarser into the offshore region. The smaller element sides along the coastline are about 0.2º, whereas in the offshore left outer boundary they reach 0.9º.
The model was forced at the surface by the monthly mean climatological wind velocity, net solar radiation (Q_{n} = Q_{sw} - Q_{lw }), the air temperature T_{air}, and the air specific humidity q_{air}. These forcing data were extracted from the Atlas of Da Silva et al. (1994). It is necessary to remark that for the evaluation of the surface heat flux Q (Eq. 21), the contributions Q_{S} and Q_{L} were obtained using the water temperature T calculated by the model, instead of the observed monthly SST (Fig. 3). The discrete fields Q_{sw}, Q_{lw},_{ }T_{air}, q_{air} were obtained interpolating the original data at each node of the finite element mesh.
The parameters of the model are mostly physically realistic values. The used densities of the upper layer, lower layer and air, ρ^{u }= 1022.5 kg m^{-3}, ρ^{l }= 1025.5 kgm^{-3} and ρ_{air }= 1.1 kg m^{-3}, respectively, are representative of the observed conditions in the region. The coefficient A_{H} is 800m^{2}s^{-1 }and K_{H} is 200m^{2}s^{-1}. The adopted entrainment-detrainment time scale t_{ed }= 1/4 day ensures that the upper layer thickness does not vanish in coastal upwelling areas, and also avoids the upper layer from becoming excessively thick in convergence regions.
In the literature, the parameters C_{L} , C_{S }used for the evaluation of the heat fluxes Q_{L} and Q_{S} range from 0.0008 to 0.002 (Bunker, 1976; Mc Creary et al., 1989; Kraus and Bussinger, 1994; Rosenfeld et al., 1994). Here, they were taken equal to 0.001. The specific heat of the air is c_{pa} = 1004 J/kg. The saturation specific humidity q was evaluated in function of the saturation vapor pressure e (see Rosenfeld et al., 1994). The latent heat of evaporation was fixed equal to L = 2.44 x 10^{6 }J/kg. The value of the coefficient of thermal expansion is α = 2.4 x 10^{3 o}C^{-1}. The time step is defined as Δt = 3 hrs.
The model response was started at rest in the month of January. The initial upper layer thickness was H = 50 m and the entrainment-detrainment depth H_{ed} = 50m. The initial upper layer and lower layer temperatures T^{u} and T^{l} were taken as 28ºC and 15ºC respectively. The model was integrated forward for a period of 1 year using the monthly forcing functions from January to December, and the results shown in Fig. 4 to 6 correspond to the summer period of the second year.
Figure 4 shows the SST fields calculated by the model for the summer months (corresponding to end January, end February and end March). The resulting patterns of the temperature contours are similar to the observed ones.
The numerical solution shows colder waters (lower than 20ºC), along the coastline and reaching also the southern offshore side of the region. This cooling process occurring in that region (from 20ºS to 30ºS) is not connected to a coastal upwelling phenomenon, but the cooling of the waters offshore is a consequence of the heat exchange fluxes between the atmosphere and the ocean dominated by the evaporation.
The calculated SST contour of 22ºC clearly describes the evolution of the warm water intrusion observed in the data. The sequence of snapshots in the summer months shows the gradual advancing of the isotherm of 22ºC to Southeast (January to February). During March, this 22ºC isotherm returns to the north. The temperature contours show that the warm water intrusion reaches almost all the coast around the 20ºS, diminishing the colder pattern in the sector comprising the 19ºS to 25ºS.
Figure 5 shows the calculated current fields during the summer months. A noticeable feature is the presence of an alongshore northward coastal jet associated with the coastal upwelling.
In the northern side (from equator to 5ºS) a westward current is developed during January, reaching maximal values of 0.36 ms^{-1}. This current weakens during February and March due to the weak winds and the reduction of evaporative cooling, leading to higher temperatures and owing to an unstable current pattern with the presence of southward flows. Note that a coastal countercurrent develops in the northeast. The southward flows are more intense in the west side of the study region.
The current fields also reveal a sector of relative calm located in the southeast side of the domain, more precisely in the sector where the meridional coastline direction changes to the northwest direction.
Finally, to highlight the performance of the proposed model, we compare the calculated SST fields showed in Fig. 4 with the observed mean monthly SST fields presented in Fig. 3. Figure 6 shows the error fields (%). Note that in the region along the eastern boundary and the equatorial band, the relative errors are less than 5%.
The largest errors (greater to 15%) occur in small sectors of the model domain located in the northeast side (January) and in the southern side (February and March).
In general, the errors during the total period of integration were positive, indicating excessive cooling produced by the model.
The mesh used in the present section is adequate for the evaluation of the hydrothermodynamic study. A comparison of the results using three meshes were conducted and the results are presented in the appendix A.
Conclusions
A 1 ½ finite element layer model for open and limited area domain is designed to describe the coupled hydrodynamic and thermodynamic processes of the surface layer of the ocean.
The numerical results obtained with this model for the hydrothermodynamics of the eastern South Pacific Ocean, during the summer months, are in a good agreement with the observed data. In particular, the occurrence of strong temperature gradients due to the intrusion of warm water to the south, and the presence of cool wind-driven coastal waters are well represented. Also, it is confirmed that the hydrothermodynamics of the upper layer in the studied region depends mainly on the interaction with the atmosphere. The relative errors are mostly less than 5% along the eastern boundary and the largest errors (~15%) occur in limited sectors.
These results indicate that the proposed model can be a useful tool to describe the most important transient circulation phenomena of the surface layer in the eastern South Pacific Ocean considering a limited area domain.
A numerical investigation of the annual cycle of the SST in this region is being prepared. Further implementation of the model with advective-diffusive equations of biological tracers will be of interest for marine biology.
Acknowledgments
This work was supported by normal fund of LNCC/MCT, FAPERJ/PRONEX. E26/171.199/2003 and by FAPERJ grant E-26/150.93/2007.
References
Anderson, D., 1984, "An advective mixed layer model with application to the diurnal cycle of the low-level East African jet (Indian Ocean)", Tellus - Series A, Vol. 36A, No. 3, pp. 278-291. [ Links ]
Anderson, D., McCreary, J.P., 1985, "Slowly propagating disturbances in a coupled ocean-atmosphere model", Journal of the Atmospheric Sciences, Vol. 42, No. 6, pp. 615-629. [ Links ]
Bova, S.W., Carey G.F., 1995, "An entropy variable formulation and Petrov-Galerkin methods for the shallow water equations". In: G.F. Carey (ed.), "Finite Element Modeling of Environmental Problems-Surface and Subsurface Flow and Transport", John Wiley: London, pp. 85-114. [ Links ]
Brooks A.N., Hughes T.J., 1982, "Streamline upwind/Petrov-Galerkin formulations for convection-dominated flows with particular emphasis on the incompressible Navier-Stokes equations", Computer Methods in Applied Mechanics and Engineering, Vol. 32, No. 1-3, pp. 199-259. [ Links ]
Bunker, A. F., 1976, "Computations of surface energy flux and annual air-sea interaction cycles of the North Atlantic Ocean", Monthly Weather Review, Vol. 104, No. 9, pp. 1122-1140. [ Links ]
Carbonel, C.A.A.H., Galeão, A.C.N., 2006, "A stabilized finite element model for the hydrothermodynamical simulation of the Rio de Janeiro coastal ocean", Communications in Numerical Methods in Engineering, Vol. 23, No. 6, pp. 521-534. [ Links ]
Carbonel, C.A.A.H., 2003, "Modelling of upwelling-downwelling cycles caused by variable wind in a very sensitive coastal system", Continental Shelf Research, Vol. 23, No. 16, pp. 1559-1578. [ Links ]
Carbonel, C., Galeão, A.C., Loula, A., 2000, "Characteristic response of Petrov-Galerkin formulations for the shallow water wave equations", Journal of the Brazilian Society of Mechanical Sciences, Vol. 22, No. 2, pp. 231-247. [ Links ]
Da Silva , A.M.M., Young, C.C., Levitus, S., 1994, "NOAA Atlas NESDIS 6", U.S. Department of Commerce, NOAA, NESDIS. [ Links ]
Donea, J., Quartapelle, L., 1992, "An introduction to finite element methods for transient advection problems", Computer Methods in Applied Mechanics and Engineering, Vol. 95, No. 2, pp. 169-203. [ Links ]
Enriquez, A.G., Friehe, C., 1997, "Bulk parameterization of momentum, heat, and moisture fluxes over a coastal upwelling area", Journal of Geophysical Research, Vol. 102, No. C3, pp. 5781-5798. [ Links ]
Fofonoff, N.P., 1962. "Dynamics of ocean currents". In: M.N. Hill (ed), "The Sea, Vol. 1", Interscience Publishers: New York, pp. 323-395. [ Links ]
Hantel, M., 1971, "Zum Einfluss des Entrainmentprozesses auf the Dynamik der Oberflachenschicht in einem tropisch-subtropischen Ozean", Deutsche Hydrographische Zeitschrift, Vol. 24, No. 3, pp. 120-137. [ Links ]
Hughes, T.J.R., Mallet, M., 1986, "A new finite element formulation for computational fluid dynamics: III. The generalized streamline operator for multidimensional advective-diffusive systems", Computer Methods in Applied Mechanics and Engineering, Vol. 58, No. 3, pp. 305-328. [ Links ]
Kraus, E.B., Businger, J.A., 1994, "Atmosphere-Ocean Interaction", 2^{nd} ed. Clarendon Press: Oxford, 362 p. [ Links ]
Lavoie, R.L., 1972, "A mesoscale numerical model of lake-effect storm", Journal of Atmospheric Sciences, Vol. 29, No. 6, pp. 1025-1040. [ Links ]
McCreary, J.P., Kundu, P.K., 1988, "A numerical investigation of the Somali Current during the Southwest Monsoon", Journal of Marine Research, Vol. 46, No. 1, pp. 25-58. [ Links ]
McCreary, J., Kundu, P., 1989, "A numerical investigation of sea surface temperature variability in the Arabian Sea", Journal of Geophysical Research, Vol. 94, No. C11, pp. 16097-16114. [ Links ]
Schop, P., Cane, M., 1983, "On equatorial dynamics, mixed layer physics and a sea surface temperature", J. Phys. Oceanogr., Vol. 13, pp. 917-935. [ Links ]
Ribeiro, F.L., Galeão, A.C., Landau, L., 1996, "A space-time finite element formulation for the shallow water equations". In: Proceedings of the International Conference on Development and Application of Computer Technics to Environmental Studies, pp 403-414. [ Links ]
Ripa, P., 1993, "Conservation laws for primitive equations models with inhomogeneous layers", Geophysical and Astrophysical Fluid Dynamics, Vol. 70, No. 1-4, pp. 85-111. [ Links ]
Rosenfeld, L.K., Schwing, F.B., Garfield, N., Tracy, D.E., 1994, "Bifurcated flow from na upwelling center: a cold water source for Monterey Bay", Continental Shelf Research, Vol. 14, No. 9, pp. 931-939, 941-964. [ Links ]
Roed, L.P., Shi, X.B., 1999, "A numerical study of the dynamics and the energetics of cool filaments, jets, and eddies off the Iberian Peninsula", Journal of Geophysical Research, Vol. 104, No. C12, pp. 29817-29841. [ Links ]
Shakib, F., 1988, "Finite Element Analysis of the Compressible Euler and Navier Stokes Equations", Phd. Thesis, Stanford University. [ Links ]
Stewart, R.H., 2004, "Introduction To Physical Oceanography Department of Oceanography", Texas A & M University. [ Links ]
Verboom, G.K., 1982, "Weakly-reflective Boundary Conditions for the Shallow Water Equations", Delft Laboratory, Publication No. 26. [ Links ]
Verboom, G.K., Stelling, G.S., Officier, M.J., 1983, "Boundary conditions for the shallow water equations". In: M.B. Abbott and J.A. Cunge (eds.), "Engineering Applications of Computational Hydraulics, Vol. I", Pitman Publ.: London. [ Links ]
Yoshida, K., 1967, "Circulation in the eastern tropical ocean with special references to upwelling and undercurrents", Japan Journal of Geophysics, Vol. 4, No. 2, pp. 1-75. [ Links ]
Paper accepted September, 2009.
Technical Editor: Francisco Ricardo Cunha
Appendix A
The influence of the mesh in the Hydrothermodynamic response is studied for three different discretizations. In Figure A1, it is presented the SST results for these three meshes. It is observed that for increasing number of nodes the SST get a similar pattern. For meshes with more than 2600 nodes the results show very small differences.