## Abstract

Measurements of deep convection from fully Lagrangian floats deployed in the Labrador Sea during February and March 1997 are compared with results from model drifters embedded in a large eddy simulation (LES) of the rapidly deepening mixed layer. The deep Lagrangian floats (DLFs) have a large vertical drag, and are designed to nearly match the density and compressibility of seawater. The high-resolution numerical simulation of deep convective turbulence uses initial conditions and surface forcing obtained from in situ oceanic and atmospheric observations made by the R/V *Knorr.* The response of model floats to the resolved large eddy fields of buoyancy and velocity is simulated for floats that are 5 g too buoyant, as well as for floats that are correctly ballasted. Mean profiles of potential temperature, Lagrangian rates of heating and acceleration, vertical turbulent kinetic energy (TKE), vertical heat flux, potential temperature variance, and float probability distribution functions (PDFs) are compared for actual and model floats.

Horizontally homogeneous convection, as represented by the LES model, accounts for most of the first and second order statistics from float observations, except that observed temperature variance is several times larger than model variance. There are no correspondingly large differences in vertical TKE, heat flux, or mixed layer depth. The augmented temperature variance may be due to mixing across large-scale temperature and salinity gradients that are largely compensated in buoyancy. The rest of the DLF statistics agree well with the response of correctly ballasted model floats in the lowest 75% of the mixed layer, and are less consistent with results from buoyantly ballasted model floats.

Other differences between observation and simulation in the mean profiles of heat flux, vertical TKE, and Lagrangian heating and vertical acceleration rates are confined to the upper quarter of the mixed layer. These differences are small contributions to layer-averaged quantities, but represent statistically significant profile features. Larger observed values of heat flux and vertical TKE in the upper quarter of the mixed layer are more consistent with model floats ballasted light. Float buoyancy, however, cannot fully account for the observed PDFs, temperature profiles, and Lagrangian rates of heating and acceleration. A test of Lagrangian self-consistency comparing vertical TKE and Lagrangian acceleration also shows that DLF measurements are not significantly affected by excess float buoyancy. These upper mixed layer features may instead be due to the interaction of wind-driven currents and baroclinicity.

## 1. Introduction

The Labrador Sea has in recent decades yielded manifold observations of open-ocean deep convection (Lazier 1973; Gascard and Clarke 1983; Lilly et al. 1999; Lab Sea Group 1998). These field experiments, and hydrographic records from Ocean Weather Station Bravo (Lazier 1980) have documented the annual formation of Labrador Sea Water during wintertime storm forcing in a central gyre 500–600 km in diameter, reaching depths greater than 2 km.

In February 1997, 11 Lagrangian drifters were deployed in the Labrador Sea during the cruise of the R/V *Knorr* (Steffen and D'Asaro 2002), in the first winter of a two-year study of the Labrador Sea under the Deep Convection Accelerated Research Initiative of the Office of Naval Research (Lab Sea Group 1998). In situ observations of the ocean and atmosphere during this cruise were subsequently used to initialize and force a large eddy simulation (LES) of the deepening mixed layer. This paper compares the Lagrangian drifter measurements to an ensemble of virtual measurements from model drifters embedded in the LES. This work is motivated by the need to understand the dynamics of rotating turbulent convection in the field, and to lay a foundation for parameterizations of deep water formation in large-scale oceanic circulation models.

Comparing LES time-averaged Eulerian statistics with the ensemble-averaged model drifter measurements can elucidate the response of drifters to convection. But comparing real drifters to model drifters embedded in the LES simultaneously compares both the response of real versus modeled drifters and the dynamics of real versus modeled convection. Consequently, statistically significant differences between observed and modeled quantities derived from drifter time series may arise from several distinct sources: 1) the suitability of the computational problem posed in the LES model as a representation of convection in the Labrador Sea, 2) the suitability of the drifter model as a representation of the float response to ambient conditions, 3) the numerical accuracy of the LES and embedded model drifters in simulating the given physical problem, and 4) uncertainties pertaining to the field experiment due to incomplete knowledge of environmental conditions and equipment characteristics. By combining the comparisons of model floats to model convection and of model floats to real floats, for several distinct statistical measures of convection, we seek to build a more unequivocal conceptual model of deep convection and to provide a critical assessment of the performance of the real drifters, the LES model, and the drifter models.

## 2. Methods

### a. Lagrangian drifter observations

Measurements of Lagrangian trajectories were made using deep Lagrangian floats (DLFs), a variety of the Lagrangian float described by D'Asaro et al. (1996), but designed to work to a depth of 2 km. The DLF aims to accurately follow the three-dimensional motion of water parcels by matching the density and compressibility of seawater and by having a large vertical drag. Compressibility is controlled by a specially designed hull and by the extrusion of a small piston. The drag is increased by a cloth drogue of approximately 1 m^{2} area. The float's vertical motion is measured by its pressure, and its horizontal motion is measured by low-frequency acoustic (RAFOS) tracking (Rossby et al. 1986). The temperature of the water at the top of the float is measured to an accuracy of 10^{−3} °C. Temperature and pressure are recorded every five minutes.

Lagrangian floats follow the water motion only imperfectly. Most importantly, the float's density differs slightly from that of the water. The compressibility of the 25-kg float is somewhat less than that of seawater. If the float density is adjusted to exactly match that of seawater at 500 db, the float will be 0.5 g heavy at the surface and 0.5 g light at 1000 db. In reality, it is difficult to control the float's density to this accuracy. The uncertainty in the float's absolute mass and volume is equivalent to about ±0.01% uncertainty in the buoyancy. Furthermore, the float's mass increases by several grams per month due to chemical reaction between the hull and seawater. The float attempts to correct for these deficiencies by varying its volume. These corrections are imperfect. For these reasons, the 1997 floats may be up to 5 g too light, equivalent to a steady-state upward velocity of 7 mm s^{−1} relative to the water. The buoyancy offset probably increases with time.

Other effects are probably less important. The floats are of finite size and are therefore insensitive to the motion of water on scales smaller than about a meter. The floats are asymmetrical and are thus potentially subject to lift forces. The float's temperature sensor is offset from its center and thus does not follow a Lagrangian trajectory, even if the float does. Finally, the float's buoyancy increases very near the surface by the weight of 3–8 g due to air pockets trapped within its crevices.

### b. Large eddy simulation

The nonhydrostatic LES model employed here was developed by Deardorff (1970, 1973), in studies of atmospheric boundary layers, and modified by Moeng (1984) to take advantage of fast Fourier transforms. Garwood et al. (1994) adapted it for oceanic convection, introducing a nonlinear equation of state that includes thermobaric effects through a depth-dependent thermal expansion coefficient, *α*(*z*) = *α*_{0} + *α*_{1}*z.* For the Labrador Sea, the LES uses *α*_{0} = 8.7848 × 10^{−5}° C^{−1} and *α*_{1} = 2.4964 × 10^{−8} (°C m)^{−1}; *β* = 7.7733 × 10^{−4} psu^{−1} is expansion coefficient for salinity. Similar LES models have been used by Denbo and Skyllingstad (1996), Harcourt et al. (1998), and Harcourt (1999) to study rotating oceanic deep convection with uniform surface forcing.

Some changes have been made to the LES of Moeng (1984) and Garwood et al. (1994). The spectral cutoff filter now removes TKE from horizontal wavenumbers *k*_{ρ} = *k*^{2}_{x} + *k*^{2}_{y} > *k*_{c}, where *k*_{c} is two-thirds the Nyquist wavenumber, rather than only filtering |*k*_{x}| > *k*_{c} and |*k*_{y}| > *k*_{c}. Since this isotropic filter removes TKE from 65% of wavenumbers versus 56% for the square spectral cutoff, closure constants are adjusted: In the notation of Moeng and Wyngaard (1988), subgrid eddy viscosity *K*_{M} = *c*_{K}Δ*e*^{1/2} and dissipation *ε* = *c*_{e}*ε*^{3/2}/Δ are now related to unresolved TKE *e* of length scale Δ = 1.33(Δ*x*Δ*y*Δ*z*)^{1/3} by *c*_{e} = 2.29 and *c*_{K} = 0.128. Scalars are not filtered, the turbulent Prandtl number is 1, and vertical scalar advection uses an upwind correction near sharp interfaces. Subsequent to these and other minor changes explained in Harcourt (1999), TKE and scalar variances each resolve part of their *k*^{−5/3}_{ρ} inertial range spectra in the turbulent mixed layer for both LES cases described below.

The LES dynamically models unresolved subgrid TKE with a scalar equation for ∂*e*/∂*t.* A critical test of an LES is that the sum of explicitly resolved large eddy (LE) TKE and parameterized subgrid (SG) TKE should be insensitive to changes in numerical resolution. The ability of an LES model to satisfy this test rests upon the validity of parameterizations used to effect turbulence closure in the subgrid TKE equation. Absent limiting stratification, turbulence closure is achieved by assuming that unresolved energy lies in the turbulent inertial subrange where Kolmogorov similarity theory can be invoked. The requirement that model resolution be sufficient to resolve a portion of the inertial subrange combines with practical limitations of computational resources to constrain model domain size.

To demonstrate the insensitivity of *net* (LE + SG) turbulence statistics to numerical resolution, simulations were carried out with two different grid spacings. The lower-resolution domain is 6400 m × 6400 m wide and 2373 m deep, with a 128 × 128 horizontal grid at 50 levels. The higher-resolution domain (Fig. 1) is 3200 m × 3200 m wide and 1804 m deep, with a 128 × 128 horizontal grid at 76 levels. Abbreviations referring to low-resolution results with 50-m horizontal grids are in lowercase italics (e.g., *le, **sg,* etc.), and uppercase abbreviations refer to (25 m) high-resolution results. The boundary condition for momentum is slip at the upper surface, with a prescribed surface wind stress. A semi-slip boundary condition at the bottom gives no mean horizontal stress, but is no-slip for velocity fluctuations with respect to the horizontal mean. Vertical velocity is zero at both upper and lower boundaries. The numerical time step Δ*t* is 2 min for the low-resolution case, and 40 s at high resolution.

The model domains are small compared to the Labrador Sea basin scale, but large compared to convective plumes. While the simulated volume is considered to be advected by the barotropic circulation, the model assumes a horizontally homogeneous local environment at scales larger than 6 km. This approach is reasonable for regimes with fairly uniform atmospheric forcing and no significant oceanic frontal features. Where interactions with mesoscale or basin-scale features are important, this homogeneous case serves as an important reference point.

The scalar fields are initialized with profiles of temperature and salinity matching those obtained from the central Labrador Sea (near 57°N, 54°W) on 12 February 1997 during the hydrographic survey of the R/V *Knorr* (Lab Sea Group 1998). Time-dependent surface wind stress and net heat flux are specified as a horizontally uniform surface boundary condition. Hourly meteorological observations from the R/V *Knorr* were used to construct the time series of surface wind stress and net heat flux. The LES surface forcing is shown in Fig. 2, with profiles from two subsequent CTD casts near the same position on 25 February and 10 March. The LES surface wind stress (Fig. 2) was deduced from observed 15-m winds using a constant drag coefficient *c*_{D} = 1.4 × 10^{−3}. Further details of the meteorological observations can be found in Lab Sea Group (1998) and in Renfrew and Moore (1999).

Wind stress and net heat flux derived from in situ observations were interpolated to the model time step to specify the LES surface boundary conditions for the 28-day period simulated, 12 February–12 March. The change in heat content of the water column between the first CTD on 12 February and the final repeat station on 10 March implies an average surface heat flux of 427 W m^{−2}, assuming no significant contribution from horizontal advection. By comparison, the average surface heat flux used for the LES was 383 W m^{−2} over the 26-day period between the two CTDs, 367 W m^{−2} for the entire 4 weeks simulated, and 370 W m^{−2} for the comparison period while DLFs were actually sampling the mixed layer. Mean surface wind stress and 15-m wind speed were 0.240 N m^{−2} and 11.27 m s^{−1} for the 4 weeks simulated, and 0.238 N m^{−2} and 11.19 m s^{−1} for the comparison period. The initial velocity field was derived from an earlier simulation of steady-state convection with comparable surface heat loss in a nonentraining turbulent mixed layer of similar depth (Harcourt et al. 1998; Harcourt 1999).

### c. Lagrangian drifter simulations

The model drifters respond only to resolved motion, not to unresolved subgrid TKE. Resolved velocity is a local volume average defined by LES resolution, and a model drifter thus moves with the mean velocity of point drifter trajectories through this finite volume. A model drifter is in this sense effectively the size of the numerical resolution, although it does not disturb the ambient flow. If the effective dimensions of either real or model drifters are small compared to the integral scale of motion, then most of the kinetic energy is resolved, and any effects due to finite drifter size should be correspondingly minor. For model Lagrangian drifters embedded in the LES, the unresolved TKE and vertical fluxes along drifter trajectories can be combined with the TKE and fluxes from the time series of depth and temperature to calculate what would be observed by an arbitrarily small float, and to test the sensitivity of the net (LE + SG) model float results to LES resolution.

Harcourt (1999) reports that a pure Lagrangian drifter model that is unaffected by the small amount of buoyancy specified in the DLF design is a poor representation of the actual float behavior. A perfectly Lagrangian model float will occasionally become embedded in the stable pycnocline and may for some time mix downward in the lowest part of the entrainment zone ahead of the advancing mixed layer. This behavior was not observed for floats in the Labrador Sea in 1997. Subsequently, model drifter results discussed here incorporate at least the small amount of float buoyancy inherent in the DLF design that causes floats to gravitate toward an equilibrium level within the mixed layer.

Large eddy fields are spectrally interpolated in two dimensions to horizontal drifter positions, then vertically spline interpolated to match the fields and their first-order finite-difference vertical derivatives at adjacent grid levels. The drifters are initially distributed in regular horizontal arrays spaced uniformly in the vertical within the initial mixed layer. Drifter trajectories are computed in tandem with the numerical simulation of the large eddy fields, using the same Adams–Bashforth time stepping for the drifter position as was used in the LES numerical scheme, and the same iteration time step Δ*t.* Using the interpolated horizontal velocity *u,* the trajectory iteration in the *x* direction for time step *i* is *x*_{i+1} = *x*_{i} + (3*u*(*x*_{i}, *t*_{i}) − *u*(*x*_{i−1}, *t*_{i−1}))Δ*t*/2.

Horizontal float velocity (*ẋ, **ẏ*) is the local resolved velocity (*u, **υ*), but *ż* for buoyant model floats differs from the LES *w.* In the low-resolution case, small amounts of float buoyancy inherent in the DLF design are accounted for kinematically in the centrally ballasted (cbl) model float. The kinematic assumption, that the vertical velocity of the fluid relative to the float *w*^{r} = *w*′ − *ż* balances drag force *F*^{D} against float buoyancy *F*^{B}, is replaced in the high-resolution LES by a dynamic centrally ballasted (CBL) model. In the dynamic model, the net force acting on the drifter and entrained fluid determines the relative deceleration of the surrounding fluid in its Lagrangian reference frame (see, e.g., Isaacson and Sarpkaya 1981),

The drag plate of the *m*_{0} = 25 kg drifter has area *A* = *πr*^{2} = 1 m^{2} and entrains fluid with density *ρ* and added mass *m*_{a} = (8/3)*ρr*^{3} = 19*m*_{0}. At time step *i,* the dynamic model solves

for *w*^{r}_{i}, with *F*^{D}_{i} = −*ρAC*_{D}|*w*^{r}_{i}|*w*^{r}_{i}, and drag coefficient *C*_{D} = 1. The buoyancy force

is the sum of the intended return force centered about *z*_{0} = 500 m with *m*_{c} = 0.5 g km^{−1}, the unintended buoyancy force *F*^{B}_{0} = −Δ*m*_{0}*g* due excess mass Δ*m*_{0}, and the buoyancy force *F*^{B}_{Θ,S} due to changes in ambient *θ* and *S.* Assuming *F*^{B}_{θ,S} = 0 at temperature *θ*_{0}, salinity *S*_{0}, and depth *z*_{0} when a float is initially deployed,

accounts for different changes in the buoyancy of the float and that of the volume *V* = 0.15 m^{3} of displaced fluid, relative to the initial condition. Here the thermal expansion coefficient for the drifter was taken to be that of aluminum (6061), *α*_{float} = 6.9 × 10^{−5}°C^{−1}, while *α*(*z*) and *β* are the expansion coefficients for potential temperature and salinity in the LES equation of state.

To test the possible effects of excess buoyancy in the 1997 DLFs, simulations were performed using both centrally ballasted (CBL and cbl) and buoyantly ballasted (BBL, bbl) float models. The CBL model has *F*^{B}_{0} = −Δ*m*_{0}*g* = 0, while the BBL model has excess buoyancy from Δ*m*_{0} = 5 g, making it initially 4.5 g too light near the surface. This value of Δ*m*_{0} is at the high end of the estimated range for the possible difference between float and fluid buoyancy. The BBL model float cannot be trapped at the rigid-lid surface boundary because it is prevented from rising above *z* = −Δ*z*/2, as though it were physically the size of a grid volume.

## 3. Comparison of real and model Lagrangian floats

This comparison section is organized as follows: First, section 3a describes observed and simulated time series and data processing, and discusses differences between Eulerian and Lagrangian statistics. Sections 3b–g compare profiles of turbulence statistics (Figs. 5–10) from the Eulerian LES fields, the real (DLF) floats and correctly ballasted CBL model, for both high- and low-resolution LES cases. Profiles from incorrectly ballasted BBL float models also appear in Figs. 5–10, but are not presented or discussed until section 3h. Section 4 summarizes bulk layer-integrated quantities in Table 3. Table 1 gives brief descriptions for common abbreviations appearing in tables and figures.

### a. Data and processing

To illustrate the Lagrangian time series from observations and LES-embedded drifter simulations, Fig. 3 shows one sample DLF vertical trajectory, and one each from CBL and BBL models. The domain sampled by each drifter ensemble is indicated in juxtaposition to the mixed layer depth *H*_{LES} and in the Eulerian time series of horizontally averaged *θ* from the LES. Mixed layer depth *H*_{LES} is determined every 50 time steps from Eulerian profiles as the depth of maximum vertical gradient in a spline-interpolated profile of *θ*. During the course of the simulation, the mixed layer deepens from 546 to 1122 m and reaches 684 m at the beginning of the indicated 3-week comparison period. For more details about the DLF ensemble elements, the autoballasting period preceding the comparison, or the restratification that follows, see Steffen and D'Asaro (2002).

Each real and simulated time series of drifter *θ* is filtered identically to decompose potential temperature into gradually changing background mean 〈*θ*〉_{LGR} and higher-frequency fluctuations *θ*′. Each time series is first unfolded at the edges, reflecting the data about the endpoints to accommodate the filter, then detrended by subtracting a linear best fit *θ*_{fit} before filtering. A high-pass finite impulse response filter *F*_{4d}(*t* − *t*′) removes fluctuations with periods longer than 4 days,

and 〈*θ*〉_{LGR} = *θ* − *θ*′. The vertical velocity record is taken to consist only of fluctuations, *w* = 〈*w*〉 + *w*′ ≃ *w*′, since for drifters in the deepening mixed layer 〈*w*〉 ≃ 0.5*dH*_{LES}/*dt* ≪ *w*_{rms}.

Every 50 time steps, Eulerian fluctuations are obtained by subtracting the horizontal average from a resolved LES field, for example, *θ*′ = *θ* − *θ*. This option is generally not available to the interpretation of Eulerian data from individual moorings, which requires the separation of the slowly changing background mean from turbulent fluctuations by a temporal method of decomposition such as the one used here for drifters. The bulk Eulerian mean,

uniformly averages together all depths within the mixed layer. In general, the overbar *θ* denotes an Eulerian average over a statistically homogeneous ensemble of turbulent fluctuations, while angle brackets 〈*θ*〉_{EUL} indicate that the ensemble may be inhomogeneous. The bulk Eulerian average is similar in this sense to the Lagrangian drifter average 〈*θ*〉_{LGR}, but their long-time means will generally differ when drifter sampling is not uniform within the mixed layer.

In Figs. 4–11, the vertical coordinate *z* of all DLF and model float time series is scaled continuously by the time-varying mixed layer depth *H*_{LES}. These data points are then binned in the scaled vertical space *ζ* = *z*/*H*_{LES} for statistical analysis. Mean profiles of Lagrangian statistics are constructed by averaging in vertical bins *ζ*_{i} ± Δ*ζ*/2 of width Δ*ζ* = 0.05, or *H*_{LES}/20 in *z.* Mean profiles of Eulerian statistics are first vertically interpolated to 1-m resolution, then scaled on the mixed layer depth and averaged over time into the same bins. The width of the bins *H*_{LES}/20 is 34 m and the beginning of the comparison period, increasing to 56 m by the end.

Turbulence statistics based upon Lagrangian and Eulerian fluctuations have certain fundamental differences. Potential temperature variance from Eulerian fluctuations *θ*′*θ*′ at a given depth will generally not agree with variance from Lagrangian fluctuations *θ*′*θ*′ averaged together at that depth since the latter includes variance due to vertical transits of the inhomogeneous profile. The problem can sometimes be overcome by constructing an Eulerian profile of Lagrangian fluctuations *θ* − 〈*θ*〉_{LGR} from the float time series by averaging measurements across depths in a fixed reference frame. If the error in estimating this profile of fluctuations is small, the Eulerian variance *θ*′*θ*′ obtained directly from the LES fields may then be compared with the mean profile of Lagrangian variance from which the inhomogeneous contribution (*θ* − 〈*θ*〉_{LGR})^{2} has been removed.

With unsteady surface forcing and entrainment, comparisons between second-order Eulerian and Lagrangian statistics are further hampered by limitations of finite time series analysis. Only part of the Eulerian temperature fluctuations are obtained in the Lagrangian decomposition because the high-pass filter removes slower fluctuations. Furthermore, the finite impulse response filter does not effect a sharp cutoff in the temperature spectrum and some higher-frequency fluctuations are lost as well. Consequently, variances and covariances involving temperature fluctuations are only expected to be comparable between modeled and observed drifters.

### b. Float probability distributions

In Fig. 4 the probability density functions (PDFs) of simulated CBL and BBL floats, and of DLF observations, are shown as a function of scaled depth *z*/*H*_{LES} within the mixed layer. Using bins of width Δ*ζ* = 1/20, the probability density function is

where *H* is the Heaviside unit step function and *N* is the number of data points.

The shaded regions ±*σ*^{CBL}_{DLF} that bracket CBL profiles in Figs. 4–9 estimate the uncertainty of mean profiles from a CBL ensemble restricted to the sample size of the DLF data. The DLF dataset begins with 11 functioning drifters at the beginning of the comparison period. A steady rate of subsequent failures leaves only three floats after three weeks. The resulting sample is equivalent to *N*_{DLF} = 4.37 drifters operating continuously for the full period, compared with *N*_{CBL} = 100 floats in each model ensemble. Uncertainty estimates *σ*^{CBL}_{DLF} also account for the difference between *N*_{CBL} and the number of model floats *n*_{CBL}(*z*) that sample a given bin. The expected standard deviation *from* the CBL mean in a subensemble of DLF sample size is

where *σ*^{CBL}_{1} is the standard deviation from the mean for elements in a set of *n*_{CBL}(*z*) profiles generated by each drifter sampling the depth *z.* For bootstrap error estimates derived from field results, see Steffen and D'Asaro (2002).

The meaning of *σ*^{CBL}_{DLF} is that ensemble mean descriptions from model drifters are of limited use in predicting the exact observations of DLFs in the Labrador Sea. Even if oceanic convection is dynamically well represented by LES, the comparison must contend with the fundamentally stochastic nature of turbulent convection. Furthermore, *σ*^{CBL}_{DLF} is only a lower bound on the expected scatter of DLF observations from the simulated CBL ensemble mean. Absent systematic differences between the LES and the drifters' environment, systematic experimental error, or significant additional sources of random error, it is most probable that only (1 − *e*^{−1}), or 63%, of the profile elements will lie within *σ*^{CBL}_{DLF} of the mean CBL profile.

Uncertainty in mean CBL profiles themselves is the standard deviation *of* the mean, approximately *σ*^{CBL}_{DLF}/5. No uncertainty estimates for DLF measurements are given for levels below *z* = −1.05*H*_{LES} not sampled by the CBL ensemble. Moreover, since DLF sampling of *z* < −1.05*H*_{LES} is very low, measurements from these depths cannot be considered statistically meaningful on their own. They are omitted from some plots by the choice of axes, but they contribute to entrainment and bulk averages (Tables 2 and 3).

Given the anticipated dispersion in the observed drifter PDF (Fig. 4), CBL model results are consistent with the DLF profile. The distribution of drifters is relatively uniform in the central and upper mixed layer, tapering off as *z* approaches −*H*_{LES}. Fluid in the lowest 25%–40% of the mixed layer is a graded mixture of “old” mixed layer water tracked by the floats and recently entrained parcels not seeded with drifters. The agreement through these depths supports the working hypothesis that entrainment processes are well represented by the LES.

Summing up elements of the normalized DLF distribution, 88% of the *N* = 2.78 × 10^{4} measurements were in the upper 3*H*_{LES}/4, 11% in the lowest 1/4 of the LES mixed layer, and 1% below *z* = −*H*_{LES}. For the CBL model floats, *N* = 4.59 × 10^{6}, and the corresponding breakdown of the distribution is 87%, 13%, and 0.05%. Most of the low-resolution kinematic cbl values lie very close to the CBL profile, demonstrating that increased resolution and greater float model complexity do not significantly affect results.

### c. Mean temperature profiles

Figure 5 shows profiles of *θ* − 〈*θ*〉,the mean deviation from the bulk mixed layer temperature. The average Eulerian difference between *θ* and the layer average 〈*θ*〉_{EUL} differs significantly from the mean profile of Lagrangian fluctuations *θ*′. There are three reasons for this. First, the two deviations are defined with respect to different means. Mean 〈*θ*〉_{LGR} is to first approximation an average of *θ* weighted by the PDF (Fig. 4):

while 〈*θ*〉_{EUL}weights depths above *H*_{LES} uniformly. Second, some fluctuating Lagrangian signal remains in 〈*θ*〉_{LGR} after filtering. As a result, drifter profiles have smaller vertical gradients than the LES Eulerian one. Third, drifters that have been in the mixed layer for some time preferentially sample older mixed layer parcels at depths near the pycnocline. These drifters are unlikely to sample either recently entrained fluid or the warm and salty dense parcels temporarily displaced above the mean Eulerian entrainment zone at *H*_{LES} by internal waves in the pycnocline.

The observed DLF profile is consistent with those of the correctly ballasted CBL floats, and there is no significant effect from LES resolution or float model indicated in the comparison of CBL and cbl profiles. The level of disagreement between model and observation is also consistent with expectations for the scatter of data, as 60% of the DLF points fall within *σ*^{CBL}_{DLF} of the CBL mean, and the balance within 2*σ*^{CBL}_{DLF}.This comparison indicates that both the surface heat flux boundary condition and the vertical transport of heat away from the surface in the LES model (i.e., the *net* vertical eddy diffusivity) is consistent with the experimental result.

### d. Vertical TKE

Both correctly ballasted (CBL, cbl) model floats are able to reproduce the corresponding large eddy (LE, le) profile of resolved vertical TKE *w*^{2} with only minor differences. Small differences are caused by the centrally directed ballasting, and by a small bias in *w* of [u[u1012]*dH*_{LES}/*dt* relative to Eulerian statistics (Harcourt 1999). Since some TKE is unresolved, profiles from real floats much smaller than LES resolution are best compared with profiles that combine large eddy and subgrid components.

In Fig. 6, a contribution to *w*^{2} from unresolved inertial subrange turbulence is combined with resolved Eulerian *ww*. Given the LES closure assumption of isotropy in the subgrid TKE, (*uu* + *υυ* + *ww*)/2 = 3*ww*/2, the subgrid contribution to *ww* is taken to be two-thirds of the horizontally averaged subgrid TKE in the LES model. This accounts for the difference between high-resolution Eulerian profiles of resolved (LE) and net (LE + SG) *w*^{2}. The profiles of *w*^{2} from the two model drifter profiles (CBL + SG, BBL + SG) in Fig. 6 also include this contribution, determined from Eulerian subgrid TKE values along the Lagrangian drifter trajectories.

Eulerian profiles (LE + SG, le + sg) of net TKE from the high- and low-resolution cases are very close, demonstrating the insensitivity of the net *ww* profiles to model resolution. The differences are small compared with the remaining unresolved (SG) component, and are comparable to the sampling error of up to 5% in the mean Eulerian profiles (instantaneous LE + SG is only sampled every 50 time steps).

The DLF profile of observed *w*^{2} in Fig. 6 is larger at shallow depths than the net CBL profile from the model (CBL + SG), but generally agrees for 0.3*H*_{LES} ≤ −*z* ≤ −0.9*H*_{LES}. The observed peak in *w*^{2} occurs much closer to the surface (between 0.1*H*_{LES} and 0.15*H*_{LES}) in the observations than it does for the net *w*^{2} in the model profiles (LE + SG, CBL + SG). The most significant differences between modeled and observed profiles of *w*^{2} appear near the surface, where DLF values exceed the net vertical TKE of both the CBL and the BBL models. Some differences near the bottom of the mixed layer may also be significant, despite the larger uncertainty due to reduced PDFs.

### e. Vertical heat flux

Eulerian profiles of net vertical heat flux (le + sg, LE + SG) are obtained by integrating the successive changes in the potential temperature profiles *θ*(*z*_{i}) upward from the no-flux bottom boundary. The average Eulerian profile of net heat flux between two mean profiles of *θ* separated by time interval *δt* is

where the no-flux boundary condition applies to the bottom of grid volumes at *i* = *N*_{z}. This method eliminates sampling error and other small analysis errors associated with different vertical interpolations and the LES numerical scheme. The large eddy (LE) profile is obtained by subtracting the subgrid component *τ*_{wθ} from the net flux profile. The difference between the net Eulerian heat flux profiles from high- and low-resolution LES averages 2 W m^{−2} in the mixed layer interior. The entrainment zone below *z* = −*H*_{LES} is sharper in the high-resolution case, in part because the horizontal domain is smaller. The LES fields explicitly resolve 96% (94%) of the layer-averaged net heat flux of 273 W m^{−2} (271 W m^{−2}) in the high (low) resolution case, with most of the unresolved subgrid component *τ*_{wθ} located adjacent to the surface. Given the surface heat flux of 370 W m^{−2} over the period of drifter comparison, the bulk Eulerian heat flux implies an average entrainment heat flux of 176 W m^{−2} (172 W m^{−2} for the lower-resolution case).

The covariance of *w* with 4-day high-passed *θ*′ determines the vertical heat flux profiles for drifters in Fig. 7. The CBL heat flux profile averages only 84% of resolved (LE) Eulerian heat flux. This reduction stems from the high-pass filter's 4-day cutoff period, and the spectral width of the cutoff needed to filter finite DLF time series with a 5-min sampling interval. Increasing the cutoff period or simply detrending to obtain *θ*′ includes more of the cospectrum of *w* and *θ*′, but increases uncertainty in heat flux values. Removing only the mean temperature may generate large random errors. Where longer time series or steady cooling permit good resolution of the slowest modes, the ratio of bulk LES covariance to bulk CBL covariance converges exponentially to unity with the square of cutoff frequency (Harcourt 1999).

The Eulerian profile of subgrid vertical heat flux *τ*_{wθ} (the difference between the LE and LE + SG plots) is added to resolved covariance flux for high-resolution CBL and BBL profiles in Fig. 7, and correspondingly for low-resolution cbl + sg. For CBL floats, this is a good approximation, since the filter primarily reduces longer timescale covariance and because CBL float sampling of turbulent fluctuations appears reasonably unbiased. High- and low-resolution (CBL + SG, cbl + sg) profiles are also close, demonstrating that the effect of different float models on covariance heat flux is small. Observations and model results all agree adjacent to the surface, confirming that LES surface heat flux is close to the actual value.

The DLF observations are consistent with CBL float results for depths below *z* = −0.2*H*_{LES}, and adjacent to the surface. For 0.2 ≥ −*z*/*H*_{LES} > 0.05, DLF results differ significantly from CBL values. In this range, the values are up to twice as large as expected, and the vertical gradient of *wθ*′ is significantly different from the trend in the central mixed layer.

### f. Temperature variance

Comparison of Lagrangian and Eulerian temperature variance requires that variance along drifter trajectories due to vertical inhomogeneity in the mixed layer be removed from the Lagrangian time series. In Fig. 8, this is done for each drifter by removing the squared mean profile of fluctuations (Fig. 5). Resulting CBL profiles of temperature variance are reduced from Eulerian values by filtering, but would approach them if the filter cutoff frequency were reduced (at the expense of accuracy) to include slower modes of fluctuation. In the upper half of the mixed layer, CBL *θ*′*θ*′ averages 3/4 of the Eulerian value, but this fraction drops below 1/3 in the bottom 20%. Higher-resolution CBL values exceed lower-resolution cbl results near the surface because more temperature variance is resolved and because no subgrid contribution is included here. Nonetheless, all model Eulerian and Lagrangian *θ*′*θ*′ profiles are similar in shape and magnitude.

In contrast, the DLF temperature variance profile is very different from any model profile. At the uppermost level, observed temperature variance is four times the CBL model drifter value. These much larger values of *θ*′*θ*′ are not directly associated with either larger heat flux or increased TKE production (Figs. 6, 7). The variance observed in the central half of the layer is two to four times the model's, without any significant enhancement to *w*^{2} or *wθ*′ at the same depths. Clearly, additional sources of temperature variance production are at hand in the Labrador Sea.

### g. Lagrangian heating and acceleration rates

Lagrangian derivatives along drifter trajectories provide important additional insights. Figures 9 and 10 show profiles of the (unfiltered) Lagrangian rates of heating *Dθ*/*Dt* and vertical acceleration *Dw*/*Dt*, respec-tively. In the layer interior, contributions to *Dθ*/*Dt* and *Dw*/*Dt* from fluid parcels that decelerate or lose heat cancel those from parcels that accelerate or gain heat.

Adjacent to the surface, heat is carried upward, out of drifter reference frame, by turbulence and diffusion at length scales smaller than the float. The effect of this surface cooling is not balanced in *Dθ*/*Dt* by equivalent subdrifter-scale flux from below. Near-surface CBL *Dθ*/*Dt* closely balances the gradient of subgrid vertical heat flux (Fig. 9), since for horizontally homogeneous convection

and because heat flux **τ**_{θ} = (*τ*_{uθ}, *τ*_{υθ}, *τ*_{wθ}) at scales effectively smaller than model floats is carried by subgrid diffusion.

Resolution effects are small but apparent in the comparison of near-surface CBL and cbl *Dθ*/*Dt*. Lagrangian cooling is more concentrated in the uppermost *H*_{LES}/20 bin for the CBL profile than for the lower-resolution (cbl) case. Since DLFs are much smaller than *H*_{LES}/20, the mean divergence of heat flux at subdrifter length scales should only contribute to observed *Dθ*/*Dt* in the uppermost bin. Given model surface heat flux 370 W m^{−2}, DLF cooling is expected to be 20*wθ*′_{z=0}/*H*_{LES} = 1.9 × 10^{−6}°C s^{−1} in the top *H*_{LES}/20, compared with 1.8 and 1.7(× 10^{−6}°C s^{−1}) in the CBL and cbl profiles, respectively.

No corresponding signal of the entrainment flux appears in the profile of −∇·*τ*_{θ}. Ultimately this is because mixed layer deepening occurs by entrainment, not diffusion. Parcels of fluid initially below the layer are lifted into it, producing temperature and salinity variance. Evidently this variance production is balanced by turbulent dissipation processes that are on average nearly isotropic in the vertical, and hence contributions to −∇·*τ*_{θ} largely cancel.

Observed *Dθ*/*Dt* (Fig. 9) is similar to model drifter profiles, but with substantially greater scatter than expected. In the central *H*_{LES}/2, only 2 out of 10 DLF bin averages fall within ±*σ*^{CBL}_{DLF} of the CBL mean, without suggesting a different trend. But in the deepest quarter of the layer, most DLF observations lie within *σ*^{CBL}_{DLF} of CBL results. Comparing this pattern of higher stochastic variability above *z* = −3*H*_{LES}/4 to observed *θ*′ variance in Fig. 8, it is apparent that larger DLF *θ*′*θ*′ values relative to the CBL profile correspond to greater uncertainty in observations of *Dθ*/*Dt* relative to *σ*^{CBL}_{DLF}. The presence of about twice the model's *θ*_{rms} probably results in an equivalent increase in random error for the DLF measurement of *Dθ*/*Dt*.

While this may explain why most DLF values lie between one and two *σ*^{CBL}_{DLF} away from the CBL profile, it does not account for larger differences in the uppermost *H*_{LES}/4. In this region, observed *θ*_{rms} is at most 6 of CBL values, while DLF and CBL values for *Dθ*/*Dt* differ by more than 4*σ*^{CBL}_{DLF}. This upper mixed layer feature of *Dθ*/*Dt* clearly corresponds to differences found between observed and modeled heat flux (Fig. 7) in the upper *H*_{LES}/4.

One of the most useful features of the profile of *Dθ*/*Dt* is that the mean Lagrangian heating rate vertically decomposes the surface and entrainment heat flux signals, allowing them to be assessed separately (Fig. 9). The arbitrary aspect of using this feature to estimate surface and entrainment^{1} heat fluxes, *τ*_{wθ}|_{z=0} and *w*′*θ*′|_{z=HLES}, is specifying which fraction of the mixed layer profile of *Dθ*/*Dt* to attribute to either source. Table 2 lists the values obtained when the depth of separation is arbitrarily taken to be the center of the mixed layer,

where *P*(*ζ*_{i}) is the fraction of data from each drifter ensemble falling into the bin at depth *ζ*_{i} = *z*_{i}/*H*_{LES}. The CBL and cbl values obtained from *Dθ*/*Dt* for the surface and entrainment heat fluxes are 32 − 55 W m^{−2} more than the corresponding net Eulerian values (Table 2). This method appears to overestimate fluxes from the model drifters, but improves with increasing LES resolution, converging toward the DLF results. The same estimates from *Dθ*/*Dt* for the observations are both lower than CBL and cbl values, and imply a bulk heat flux (*τ*_{wθ}|_{z=0} + *w*′*θ*′|_{z=−H})/2 = 281 W m^{−2} that is only about 10 W m^{−2} larger than the net LES Eulerian bulk flux. Differences between model and observed heat flux values are comparable to the uncertainty estimate *σ*^{CBL}_{DLF}. Furthermore, this uncertainty should scale with *θ*_{rms}, so *σ*^{CBL}_{DLF} probably underestimates DLF error.

The observed Lagrangian mean rate of vertical acceleration *Dw*/*Dt* in Fig. 10 is consistent with the CBL model profiles over the lowest 75% of the mixed layer. The larger scatter in the data seen in the Lagrangian heating rate is not evident here because, unlike *θ*′*θ*′, the observed variance in *w* is consistent with the LES over these depths. In the upper *H*_{LES}/4 of the layer, the DLF profile of *Dw*/*Dt* shows another manifestation of the near-surface feature that also affects *ww*, *wθ*′, and *Dθ*/*Dt*. Near-surface DLF acceleration is concentrated in the uppermost 5% of the layer, with *Dw*/*Dt* more than twice the CBL value adjacent to the surface. Elsewhere in the top *H*_{LES}/4, observed downward acceleration is significantly less than the CBL model prediction.

### h. Excess float buoyancy

Initial DLF deployment field reports and comparisons with LES-embedded model floats raised considerable concern that DLF observations were affected by unintended excess float buoyancy. Observations of vertical TKE and heat flux (Figs. 6, 7) in the upper *H*_{LES}/4 are more consistent with the incorrectly ballasted BBL model floats than the correctly ballasted CBL model. Physically, values of *w*^{2} and *wθ*′ are larger for buoyant drifters because only the colder, more energetic plumes generated in strong surface convergence zones are able to overcome the buoyancy bias and convect these floats well below the surface. Once at depth, the extra buoyancy will lift them out of more quiescent zones toward the surface, inhibiting the sampling of these low TKE regions with a correct Lagrangian probability.

However, better agreement between DLF and correctly ballasted CBL results for *ww* and *wθ* in the central mixed layer indicates that any excess DLF buoyancy only affected sampling near the surface. Such near-surface buoyancy might arise from bubbles adhering to the float or from compression of flexible seals in the drifter hull. The observed profile of (*θ* − 〈*θ*〉) (Fig. 5) is not inconsistent with this possibility, since the overall shift in BBL values is due to sampling biases in the lower mixed layer, where ∂*θ*/∂*z* is large. Observed *θ*′*θ*′ (Fig. 8) is equally inconsistent with both BBL and CBL model results.

The contrast between the PDFs (Fig. 4) for CBL and BBL model floats shows that excess DLF buoyancy from Δ*m*_{0} = 5 g should cause a large upward shift in the observed float distribution. The BBL profile in Fig. 4 shows oversampling by a factor of 2 in the upper 10% of the mixed layer, and a PDF reduction spread over the lower 2*H*_{LES}/3. The DLF distribution is consistent with that of the correctly ballasted floats, and inconsistent with the excessively buoyant BBL model results. The PDF (not shown) for the low-resolution bbl model was actually less skewed toward the surface than the BBL profile, even though less TKE is resolved. The absence of such an effect in observed PDFs strongly refutes the possibility that the DLF sampling is biased by incorrect drifter ballasting at any depth.

The Lagrangian heating and vertical acceleration profiles (Figs. 9, 10) from BBL model floats are consistent with the scenario of excess near-surface buoyancy at some depths, but inconsistent in the uppermost bin, where the large observed DLF cooling and downward acceleration is not replicated in BBL model results. No trend toward such behavior is found for increasing LES resolution when comparing bbl (not shown) and BBL results.

The mean Lagrangian rate of vertical acceleration *Dw*/*Dt* (Fig. 10) can then be integrated downward from the surface to provide a test of Lagrangian self-consistency. Assuming that Eulerian ∂*w*/∂*t* = 0, and that horizontal stresses are horizontally homogeneous, ∂*uw*/∂*x* + ∂*υw*/∂*y* = 0, then

where *ww*(0) = 0 (surface waves should not directly contribute to DLF pressure fluctuations). If the averaging of *ww* and *Dw*/*Dt* is into bins of scaled depth *ζ* = *z*/*H*_{LES} over the comparison period, then

Vertical TKE from integrated *Dw*/*Dt* may disagree with directly computed variance in *w* (Fig. 11) if float trajectories are not truly Lagrangian, but they may also disagree if either shear stress convergence in the horizontal plane or mean Eulerian vertical acceleration are not negligible. For the horizontally homogeneous LES simulations, ∂*w*/∂*t* = ∂*uw*/∂*x* = ∂*υw*/∂*y* = 0 for averages over the Eulerian field. Even so, the two profiles are inconsistent for CBL floats in the lower 40% of the mixed layer. This is because the mean float PDF (Fig. 4) is spreading downward with entrainment. Because of this downward flux of drifters, the *three-dimensional* PDF is correlated with *w*′ at these depths, contributing to terms neglected in Eq. (14). But in the upper 60% of the layer where the mean PDF profile is more uniform, integrated *Dw*/*Dt* agrees to ±10% with the variance of CBL velocity. Small differences may be due to random sampling error or to the intended DLF buoyancy characteristics incorporated in the CBL model.

The two BBL profiles in Fig. 11 show that trajectories of excessively buoyant floats fail the test of Lagrangian self-consistency at all levels. While the BBL variance of *w* is 24% larger than the corresponding Eulerian profile (Fig. 6), the integrated Lagrangian acceleration is reduced from the Lagrangian variance by 30% to over 100% in the upper 0.6*H*_{LES}. The effects of excess float buoyancy are cumulative in the BBL integral of *Dw*/*Dt*. Differences between the two DLF profiles in the upper 0.6*H*_{LES} are more random in sign and smaller, indicating that DLF trajectories are nearly Lagrangian in the upper mixed layer. Good agreement at middepths indicates no net excess buoyancy in the upper layer, and the sign of the difference for DLFs in the uppermost bin is opposite that expected for excessively buoyant floats. It seems rather unlikely that shear stress convergence in the horizontal plane or mean Eulerian vertical acceleration are significant, but most unlikely that they would compensate almost exactly for the effect of any near-surface excess DLF buoyancy in the integral of *Dw*/*Dt* at middepths.

## 4. Bulk layer averages

Table 3 summarizes bulk layer averages of vertical TKE, heat flux, and potential temperature variance from the LES model, the embedded drifter models, and the 1997 Labrador Sea Lagrangian drifter observations. Bulk mixed layer averages are given both as averages over float ensembles (〈 〉_{LGR}), and as averages of mean profiles over *z* > −*H*_{LES} (〈 〉_{EUL}). Eulerian layer-averaged values are given for both resolved (LE) and net (LE + SG) turbulence statistics. Anticipated errors in bulk DLF quantities are computed as *σ*^{CBL}_{DLF} = *σ*^{CBL}_{1}/*N*_{DLF}, using the CBL ensemble standard deviation from the mean *σ*^{CBL}_{1}.

Corresponding profiles (Figs. 6–8) should be considered when interpreting bulk averages in Table 3. For example, proximity of bulk DLF potential temperature variance 〈*θ*′^{2}〉_{EUL} (Table 3) and the Eulerian (LE) value is misleading, since real differences for *z* > −3*H*_{LES}/4 cancel sampling differences from the lowest *H*_{LES}/4 in the bulk average.

Observed DLF 〈*w*^{2}〉 and 〈*wθ*′〉 are marginally higher than CBL results, and less than BBL model values. Although the differences between CBL and DLF values are statistically significant, the fractional differences in the *rms* turbulent velocity scale are small. The observed vertical velocity scale *w*_{rms} = 〈*w*^{2}〉^{1/2}_{LGR} of 2.31 cm s^{−1} is only 13% above the CBL model value of 2.03 cm s^{−1}.

The net bulk heat flux from the high (low) resolution LES can be combined with the surface heat flux of 366 W m^{−2} in the model to estimate the entrainment heat flux of 180 W m^{−2} (176 W m^{−2}). This value is generally consistent with estimates of the entrainment flux obtained from all three drifter profiles of covariance (Fig. 7) by extrapolating the trend across the well-sampled central mixed layer to *z* = −*H*_{LES}. This estimate of *wθ*′|_{z=−H} is lower but consistent with the entrainment heat flux derived from the CBL (cbl) model Lagrangian heating rate, and only 8 W m^{−2} (10 W m^{−2}) more than the entrainment heat flux estimated by the same method from DLF observations.

## 5. Discussion

The comparisons point to two distinct patterns of difference between DLF observations and CBL model results. One of these patterns appears only in the greatly augmented level of temperature variance observed throughout the upper 3/4 of the mixed layer, while the other pattern appears in a set of interrelated features for *Dθ*/*Dt*, *Dw*/*Dt*, *wθ*′, and *w*^{2} that are confined to the upper quarter of the layer.

Given the incongruity at middepths between the comparisons of model and experiment for *w*^{2} or *wθ*′ (Figs. 6, 7) on the one hand, and for *θ*′*θ*′ (Fig. 8) on the other, it appears that there are two separate processes that generate these two patterns of differences. Although both of these processes may arise from large-scale horizontal gradients in temperature and salinity, the dynamics behind these two features exhibit distinctly different vertical length scales. The near-surface feature appearing in DLF profiles of *w*^{2} and *wθ*′ exhibits coherence on vertical scales less than 250 m, while the process responsible for the large temperature variance may be simply limited by the depth of the mixed layer.

Except for observed temperature variance, DLF profiles are in good agreement with CBL model floats over the lowest 75% of the layer. Small yet significant differences in the upper 25% of *ww*, *wθ*′, *Dθ*/*Dt*, and *Dw*/*Dt* profiles are clearly related to a process that is either hydrodynamic and not included in the physical problem modeled by the LES, or a dynamic aspect of the DLF response to convection that is not included in the embedded drifter model.

Taken together, the comparisons between BBL model floats and experimental results (Figs. 4–11) discount the possibility that any unintended excess float buoyancy significantly affected DLF observations. Although elevated levels of heat flux and vertical TKE (Figs. 6, 7) in the upper mixed layer could have been caused by unintended float buoyancy, the effect of 5 g in excess model float buoyancy on the mean Lagrangian rate profiles adjacent to the surface (Figs. 9, 10), on the mean profile of *θ*′, and on the PDF of the floats (Fig. 4) is inconsistent with the observed DLF response. Furthermore, a test of Lagrangian self-consistency (Fig. 11) indicates that DLF trajectories were not significantly affected by excess buoyancy. In the uppermost 0.15*H*_{LES} some inconsistency between *w*^{2} and the integrated profile of *Dw*/*Dt* might be generated by lift forces, but this effect is not much larger than random error in DLF profiles. Both DLF profiles in Fig. 11 differ significantly from CBL *w*^{2}, indicating that non-Lagrangian float behavior is not the source of this disagreement.

The much larger level of temperature variance and the absence of a corresponding increase in vertical heat flux show that most temperature fluctuations are uncorrelated with vertical velocity. If not associated with vertical heat flux, the large observed temperature variance may be produced by horizontal heat flux (Legg and Williams 2000). This horizontal heat flux could be generated in regions with large-scale horizontal gradients of temperature and salinity when baroclinic currents are destabilized by vertical mixing across the baroclinic shear (Stone 1999). Such large variances may also be generated by convection in regions with mean horizontal temperature and salinity gradients that are compensating in buoyancy, or when horizontal gradients of buoyancy and spiciness (the state variable orthogonal to buoyancy) are not aligned.

Various mechanisms of baroclinic instability have been found for convection with localized surface heat loss simulated numerically by Visbeck et al. (1996) and in laboratory models by Maxworthy and Narimousa (1994) and by Coates et al. (1995), and also found by Legg et al. (1998) for convection localized by a mesoscale eddy. The horizontal heat fluxes carried by baroclinic eddies in these models may be involved here in production of additional temperature variance, but the absence of a drifter-observed feature of *wθ*′ or *w*^{2} in the central mixed layer unambiguously related to these eddies makes their presence difficult to demonstrate. Buoyancy-compensated vertical temperature and salinity gradients in the pycnocline result in a well-observed entrainment flux signal (Fig. 7), which clearly demonstrates that a cessation of entrainment deepening has not been brought about by a horizontal flux of heat borne by baroclinic eddies.

A small reduction of entrainment deepening is still a possibility; comparison of CBL and DLF entrainment flux from mean Lagrangian heating rates (Table 2) suggests that the observed ratio 170/391 = 0.43 between entrainment and surface heat fluxes is lower than the simulation value of 215/402 = 0.53 for the 25-m-resolution LES, and 227/410 = 0.55 for the 50-m-resolution case. However, these differences are comparable to *σ*^{CBL}_{DLF}, which probably underestimates the measurement error. Moreover, agreement between the simulated CBL *w*′*θ*′ (Fig. 7) profile and the observations in the central and lower mixed layer shows that if there is a process of horizontal advection responsible for reducing the entrainment flux by 40–50 W m^{−2}, it is confined vertically to the entrainment zone. One element of the profile of vertical TKE (Fig. 6) is significantly larger than anticipated in the entrainment zone, supporting the hypothesis that some additional process not realized in the LES may be affecting dynamics in the pycnocline.

The DLF measurement of vertical heat flux adjacent to the surface in Fig. 7 is consistent with in situ meteorological observations that determine *wθ*′|_{z=0} in the LES. If the DLF response to convection is sufficiently Lagrangian, then the extra cooling found adjacent to the surface in *Dθ*/*Dt* (Fig. 9) can only occur if fluid tracked by the float mixes horizontally with fluid that is colder. Mean Lagrangian heating due horizontal mixing is

where unresolved fluxes *τ*_{iθ} in the DLF context refer to subdrifter-scale diffusion. Excess DLF heat loss relative to CBL is balanced by warming from horizontal mixing at depths between 0.05*H*_{LES} and 0.25*H*_{LES}. This feature of the Lagrangian heating rate is clearly responsible for increased values of both *wθ*′ and its vertical gradient over the same depths. The DLF profile alone is inconsistent with observed water column changes and conservation of heat in the horizontally homogeneous LES. It is therefore not possible that near-surface differences found between simulated CBL and observed DLF profiles are primarily due to numerical grid size, lack of a free surface, or treatment of subgrid TKE in either the LES or drifter models.

The connection between the manifestation of this near-surface feature in profiles of *Dθ*/*Dt* and *wθ*′ and its homologous appearance in profiles of *Dw*/*Dt* and *ww* suggests that the horizontal advection of heat in the drifter reference frame described above is also accompanied by an extra buoyancy loss at the surface that is regained in the depth interval 0.05 ≥ −*z*/*H*_{LES} > 0.25. This extra gain and loss of buoyancy in the upper *H*_{LES}/4 could account for the additional downward acceleration found in the DLF profile of *Dw*/*Dt* (Fig. 10) adjacent to the surface beyond that found for the CBL model result. This would also account for the significant increase in vertical TKE over the CBL values (Fig. 6) in the upper quarter of the layer. But, it is not clear in this scenario whether the fluid parcels tracked by DLFs are losing additional buoyancy at the surface because of the extra cooling, or simply that they are surrounded by increasingly buoyant ambient fluid that is colder and fresher.

If the near-surface difference between CBL and DLF profiles is due to a dynamic process in the upper ocean, then it represents the effects of horizontal currents or advection in the presence of large-scale temperature and buoyancy gradients not represented in the LES model. These horizontal processes are confined to the upper mixed layer (and possibly near the pycnocline as well) and do not significantly affect *w*^{2}, *wθ*′ or the vertical gradient of temperature in the central mixed layer. Horizontal advection of heat into and out of the drifters' reference frame could be generated by Ekman circulation in the vertical plane or by other ageostrophic currents horizontally transporting fluid over large-scale gradients of heat and buoyancy.

In a horizontally homogeneous environment, the mean Lagrangian heating rate in a layer *dz* adjacent to the surface is *Dθ*/*Dt*_{SG} = −*τ*_{wθ}|_{z=0}/*dz,* and zero at depths below *z* = −*dz,* if *dz* is much larger than the smallest scale of resolved TKE. Using surface heat flux *τ*_{wθ}|_{z=0} = 391 W m^{−2} estimated from the DLF heating rate (Table 2), and averaging into the LES grid using bins of width *dz* = 23.74 m, Fig. 12a shows the difference Δ(*Dθ*/*Dt*)_{SG} = (*Dθ*/*Dt*)_{DLF} − (*Dθ*/*Dt*)_{SG} between the DLF-observed heating rate and the profile expected in the homogeneous environment. The uncertainty estimate *σ*_{DLF} in this plot is *σ*^{CBL}_{DLF} rescaled by the ratio between observed and modeled CBL *θ*^{′}_{rms} (see Figs. 8, 9, and related text). The profile Δ(*Dθ*/*Dt*)_{CBL} = (*Dθ*/*Dt*)_{DLF} − (*Dθ*/*Dt*)_{CBL} of this Lagrangian heating anomaly computed directly from DLF and CBL profiles is also shown, and differs somewhat from Δ(*Dθ*/*Dt*)_{SG} because of model resolution.

The variation of *Dθ*/*Dt*_{SG} with depth is similar to the profile of wind-driven current in certain directions. The mean LES surface wind stress for this period was 0.238 N m^{−2} in the direction *ϕ* = 113°, measured clockwise from the north. Figure 12b shows the correlation coefficient over the 10 uppermost LES grid levels, *γ* = 〈*Dθ*/*Dt*_{SG}*U*(*ϕ*)〉/(〈(*Dθ*/*Dt*_{SG})^{2}〉〈*U*(*ϕ*)^{2}〉)^{1/2}, between the anomalous Lagrangian heating and the mean horizontal velocity *U*(*ϕ*) as a function of direction *ϕ.* The anomalous Lagrangian heating and *U*(*ϕ*) are most correlated at *ϕ* = 272° and most anticorrelated for *ϕ* = 92°. The high correlation indicates that a linear relationship between *Dθ*/*Dt*_{SG} and *U*(*ϕ*) in these directions can potentially account for 91% of the anomalous Lagrangian heating.

A Lagrangian parcel with velocity **U** moving in horizontal temperature gradient ∇*θ*_{G} becomes warmer or colder than ambient fluid at a rate −**U** · ∇*θ*_{G}. In steady state, subdrifter-scale mixing balances this change in relative *θ.* If the average effect on *Dθ*/*Dt* from wind-driven transport across a large-scale temperature gradient *θ*_{G} is simply **U** · ∇*θ*_{G}, the slope of the best fit to the horizontal current in heading *ϕ* estimates the temperature gradient. At *ϕ* = 272° (or *ϕ* = 92°), the gradient indicated is ±1°C in 13 km. The current profile in the direction *ϕ* = 150° (or *ϕ* = 330°) of mean LES surface transport at 3.9 cm s^{−1} is somewhat less correlated (*γ* = 0.80), but implies a smaller temperature gradient of ∓1°C in 23 km. In the region of DLF deployment, the largest surface temperature gradients observed in the vicinity during float deployment are about 1°C in 20 km (Pickart et al. 1998).

These estimates suggest that wind-driven or other ageostrophic mechanisms may be extracting potential energy from a baroclinic density field. Further simulations of deep convection with substantial wind forcing in the presence of baroclinic density gradients will be necessary to properly investigate these possibilities. In particular, it is probably incorrect to assume that baroclinicity will not alter that part of the mean current profile due to surface wind stress. Surface-intensified ageostrophic circulations induced by the interaction of free convection and baroclinicity can make similar contributions of the Lagrangian heating rate. Such secondary circulations may carry significant buoyancy across large-scale density gradients with horizontal velocities averaging 0.1 − 1 cm s^{−1} (Legg et al. 1998; Stone 1999).

It is also possible that near-surface anomalies found in *Dθ*/*Dt* and *w*′*θ*′ are unrelated to differences between DLF and CBL profiles of *w*′*w*′ and *Dw*/*Dt* in the upper mixed layer. If buoyant TKE production from elevated levels of *w*′*θ*′ is compensated by a larger salinity flux, homology between *Dθ*/*Dt* and *Dw*/*Dt* anomalies would be coincidental. Interaction of wind-driven currents and compensated temperature gradients would result in the observed nonlinear heat flux profile, while differences in TKE profiles may be due to shortcomings of the LES.

Net TKE in the LES model might be low near the surface because assumptions for the LES subgrid parameterization are less valid there. However, the absence of a trend toward better agreement in *w*^{2} with increased LES resolution indicates no such shortcoming in the subgrid parameterization. Larger observed *w*^{2} may simply be due to incorrect LES surface forcing, as higher surface wind stress or surface heat flux would also produce more vertical TKE. The 11% larger surface heat flux deduced from the heat budget of the water column between CTDs on 12 February and 10 March, and the 6% larger value derived from the Lagrangian heating rate, can account for at least some difference in TKE production.

The lack of surface gravity waves in the LES might also reduce TKE production, as suggested by Skyllingstad and Denbo (1995). Also, the effect of breaking waves on pressure measurements in a Lagrangian reference frame may be substantial at depths on the order of 1–10 m below the surface. Another possibility is that large-scale vorticity may either augment or counteract TKE reduction due to planetary rotation. This effect might depend on depth, altering both the magnitude and shape of the vertical TKE profile.

Such separate explanations of why observed TKE is larger than in the LES must still contend with the potential effects of baroclinicity. Even if large-scale gradients of temperature and salinity are compensated at one level in a deep mixed layer, the effect of thermobaricity produces horizontal density gradients at other depths. If DLF heat flux measurements imply the horizontal advection of a substantial temperature gradient, the associated baroclinicity will also affect TKE.

## 6. Conclusions

Simulations of drifting sensors embedded in a large eddy simulation (LES) of a small domain in the 1997 wintertime Labrador Sea have demonstrated the effectiveness of deep Lagrangian floats (DLFs) in measuring the quasi-steady-state dynamical statistics of a deeply convecting mixed layer. By comparing results from buoyantly ballasted model floats to 1997 DLF measurements of the rapidly deepening mixed layer in the central Labrador Sea, it was possible to deduce that these observations were unaffected by any unintended excess buoyancy, despite concerns to that effect.

Disagreements between mean profiles of statistical measures of turbulence from DLF observations and from correctly ballasted model floats embedded in the LES fall into two distinct patterns. The first is that observed temperature variance is much greater in the upper 3/4 of the layer than it is in the model. This difference appears to have no direct impact on the TKE budget of the layer and may be due to mixing across large-scale buoyancy-compensating gradients of temperature and salinity that are not represented in the LES.

The second pattern of disagreement is confined to the upper ocean in the mean profiles of vertical turbulent kinetic energy (TKE), heat flux, and the Lagrangian derivatives of heat and vertical momentum. The systematic interrelationship between these upper ocean features and their coincidence with a component of the wind-driven Ekman spiral suggest that they may be due to the interaction of wind-driven horizontal currents with a large-scale mean gradient of temperature that is uncompensated in buoyancy by a coincident gradient in salinity.

The otherwise good agreement between simulated and actual drifters demonstrates that the numerical technique of LES with embedded drifter models can be successfully used to simulate the high–Reynolds number turbulence of deep convection and the response of drifting sensors. The areas of disagreement encountered in the comparison appear to be small perturbations to the underlying mixed layer budgets of TKE and buoyancy variance. To understand the effects of large or mesoscale features on the small-scale vertical mixing, the horizontally homogeneous scenario modeled and presented here serves as an indispensable reference point.

## Acknowledgments

This research was supported by grants from the Office of Naval Research (Codes N001400WR20250 and N00014-94-1-0025) and by a grant from the National Science Foundation (OPP-9530530). LES computations were accomplished with a grant of HPC resources from the DOD Office of High Performance Computing.

## References

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

## Footnotes

*Corresponding author address:* Dr. Ramsey R. Harcourt, Dept. of Oceanography/Code OC/Ha, Naval Postgraduate School, 833 Dyer Road, Bldg. 232, Rm. 328, Monterey, CA 93943. Email: harcourt@oc.nps.navy.mil

^{1}

Entrainment flux reference level is figurative, invoking an equivalent layer model that is homogeneous above *z* = −*H*_{LES}.