Articles | Volume 12, issue 4
Research article
16 Nov 2021
Research article |  | 16 Nov 2021

Accounting for surface waves improves gas flux estimation at high wind speed in a large lake

Pascal Perolo, Bieito Fernández Castro, Nicolas Escoffier, Thibault Lambert, Damien Bouffard, and Marie-Elodie Perga

The gas transfer velocity (k) is a major source of uncertainty when assessing the magnitude of lake gas exchange with the atmosphere. For the diversity of existing empirical and process-based k models, the transfer velocity increases with the level of turbulence near the air–water interface. However, predictions for k can vary by a factor of 2 among different models. Near-surface turbulence results from the action of wind shear, surface waves, and buoyancy-driven convection. Wind shear has long been identified as a key driver, but recent lake studies have shifted the focus towards the role of convection, particularly in small lakes. In large lakes, wind fetch can, however, be long enough to generate surface waves and contribute to enhance gas transfer, as widely recognised in oceanographic studies. Here, field values for gas transfer velocity were computed in a large hard-water lake, Lake Geneva, from CO2 fluxes measured with an automated (forced diffusion) flux chamber and CO2 partial pressure measured with high-frequency sensors. k estimates were compared to a set of reference limnological and oceanic k models. Our analysis reveals that accounting for surface waves generated during windy events significantly improves the accuracy of k estimates in this large lake. The improved k model is then used to compute k over a 1-year time period. Results show that episodic extreme events with surface waves (6 % occurrence, significant wave height > 0.4 m) can generate more than 20 % of annual cumulative k and more than 25 % of annual net CO2 fluxes in Lake Geneva. We conclude that for lakes whose fetch can exceed 15 km, k models need to integrate the effect of surface waves.

1 Introduction

Lakes are universally regarded as significant sources of CO2 to the atmosphere; however, the accurate quantification of the magnitude of such emissions currently remains challenging (Cole et al., 2007; Tranvik et al., 2009; Raymond et al., 2013). CO2 fluxes can be directly measured with floating chamber or eddy covariance systems (Vachon et al., 2010; Vesala et al., 2006). However, both approaches have their own constraints. The former suffers from limited time and space integration (from minutes to hours and from centimetres to metres respectively; Klaus and Vachon, 2020), whereas the latter remains technically difficult and can be influenced by non-local processes (entrainment from the shore or advection; Vachon et al., 2010; Esters et al., 2021). Thus, long-term direct flux measurements are mostly restricted to small lakes (Huotari et al., 2011), and fluxes remain mostly estimated with models. CO2 fluxes at the surface of lakes operate through a net diffusive transport, obeying the first Fickian law:

(1) F = k α Δ p CO 2 ,

where F (molm-2s-1 but often expressed as µmolcm-2h-1) is the CO2 gas flux, α is the CO2 solubility coefficient (µmolcm-3µatm-1), ΔpCO2 is the gradient of partial pressure of CO2 (pCO2) between the water and the atmosphere corrected for altitude (µatm), and k is the gas transfer velocity (cm h−1).

Therefore, lake carbon emissions are primarily driven by the gradient of partial pressure of CO2 between the surface lake water and the atmosphere, but the gas transfer velocity controls the rate of CO2 exchange across the lake–atmosphere interface. Assessing the amount of lake CO2 emissions to the atmosphere has been a major issue, starting with the Cole and Caraco (1998) seminal paper, with debates regarding both the representativeness of the measurements and the optimal conceptual model for air–water gas transfer (e.g. MacIntyre et al., 2001; Borges et al., 2004). As recent developments in sensor technologies allow continuous and accurate measurements of aqueous CO2 concentrations, the gas transfer velocity currently remains the main source of uncertainties, which hinders attempts to achieve full carbon budgets (Dugan et al., 2016) or to quantify greenhouse gas emissions by lakes, at local, regional, or worldwide scales (Maberly et al., 2013; Raymond et al., 2013; Engel et al., 2018).

Figure 1Conceptual scheme of the four main processes driving gas transfer velocity (k) in a large lake induced by wind and cooling events. These four processes are split into two types of k: k-bubble for the bubble formation (kB=kb) and k-no bubble for the convective mixing, wind shear, and wave action term which are added (kNB=kc+ku+kw). Below this scheme, a non-exhaustive review of the conceptual approaches of k models used in first Fickian law is given. From left to right, the increase in the complexity level of k models and their study site (limnological to oceanic case) is visible. All of these variables are described in Sect. 2.4 and Table 1.


Table 1Summary of the characteristics of kSc models for predicting the air–water gas transfer velocity based on wind speed (CC98 and CW03) and lake size (VP13), surface renewal model (T14, R12, and S07), COARSE approach (DM18), and both adapted models, namely SD21 and SD21-fit, from a combination of S07 and DM18.

 Adaptation of εWave action for a large lake: CC98 – Cole and Caraco (1998), CW03 – Crusius and Wanninkhof (2003), VP13 – Vachon and Prairie (2013), T14 – Tedford et al. (2014), R12 – Read et al. (2012), S07 – Soloviev et al. (2007), and DM18 – Deike and Melville (2018).

Download Print Version | Download XLSX

k is inherently tied to turbulent mixing within the surface boundary layer, which enhances the diffusive gas exchange by renewing the surface mass content (Zappa et al., 2007). At the lake–atmosphere interface, turbulent mixing is the product of wind shear (ku), buoyancy flux (kc), and wind-driven surface waves, whose effect can be split into wave action (kw) and wave breaking (kb or kB), with the latter producing air bubble and water spray (Fig. 1; Wüest and Lorke, 2003, Soloviev et al., 2007). Regarding the prominent role of wind action on surface turbulence, first quantitative models have empirically scaled k to wind speed (referenced at a 10 m height; U10), as a proxy for the level of wind-driven turbulence (Fig. 1; Cole and Caraco, 1998; Crusius and Wanninkhof, 2003). The parameterisations of the k–wind relationships vary between authors (e.g. Klaus and Vachon 2020), as a likely consequence of the local characteristics of the lakes used in the calibration datasets (Table 1). However, all studies suggested a polynomial relationship between U10 and k with an exponent larger than 1. Further development of empirical models integrated the lake surface area as a second parameter in the k–wind relationships, to account for the role of the fetch length for wind action (Vachon and Prairie, 2013). Generally, empirical wind-based models tend to underestimate fluxes, especially at low wind speed (e.g. Schubert et al., 2012; Heiskanen et al., 2014; Mammarella et al., 2015) where turbulent mixing through buoyancy flux is expected to take over wind shear. Moreover, these empirical models require a proper calibration each time they are applied in a new system with different characteristics, i.e. a new set of lakes and/or meteorological conditions (Klaus and Vachon, 2020), thereby limiting their universal applicability.

In parallel to empirical wind-based models, process-based models attempt to link k directly to near-surface turbulence. The surface renewal model (SRM) is one of the first, and still most widely used, theories (Danckwerts 1951; Lamont and Scott, 1970) with k depending on the product of the turbulent kinetic energy dissipation rate (ε) and the kinematic viscosity of water (υ), both to a power of one-quarter as follows:

(2) k = a 1 ( ε υ ) 1 / 4 Sc - 1 / 2 ,

where a1 is a calibration constant parameter, and Sc is the Schmidt number. Recently, Lorke and Peeters (2006) and Katul and Liu (2017) demonstrated that this relationship, to which different approaches converge, can be seen as a universal scaling. As opposed to the practical empirical models presented above, process-based models have the potential to predict k using the turbulent dissipation rate over a wide range of environmental conditions extending beyond those encountered in the calibration dataset (Zappa et al., 2007). As for lakes, SRM k models have so far considered the friction velocity at the water side (u*, wat) and the turbulence created by thermal convection using the buoyancy flux at the surface (B0) (Fig. 1; Eugster et al., 2003; MacIntyre et al., 2010; Read et al., 2012; Tedford et al., 2014; Heiskanen et al., 2014). It is noteworthy that the SRM approach leads to k being related to u*, wat (or U10) to the first order (Wanninkhof, 1992; Lorke and Peeters 2006; see Sect. 2), whereas the empirical models described above predict a higher-order polynomial relationship. This inconsistency is tentatively solved in oceanography by adding another source of gas exchange associated with wind-waves' whitecaps. Early gas flux parameterisation already accounted for wind and buoyancy-driven turbulence as well as surface waves (Fig. 1; Woolf et al., 1997; Soloviev et al., 2007; Fairall et al., 2011). However, the buoyancy-driven contribution can often be neglected in oceanography, and recent efforts have been dedicated to a better parameterisation of the bubble enhancement term (Fig. 1; Deike and Melville, 2018). In lakes, wind fetch can be long enough to generate surface waves (Wanninkhof, 1992; Frost and Upstill-Goddard, 2002; Borges et al., 2004; Guérin et al., 2007), implying that surface waves could be a significant driver of k and subsequent CO2 fluxes (Schilder et al., 2013; Vachon and Prairie, 2013). Thus, the role of surface waves has been essentially empirically accounted for in lake k models, through the polynomial scaling to U10 in wind-based models, and most often neglected in studies using process-based parameterisations (mainly SRM) (e.g. Read et al., 2012). While this approximation may be appropriate for small-shielded lakes, it is likely to be insufficient in larger, long-fetched lakes.

Herein, we aim to identify the most adequate k model for Lake Geneva, a large, clear, hard-water lake in the Swiss Alps, to assess k values over a full annual cycle. We compare the performances of different models of gas transfer velocity, in their original or slightly modified published formulations from the limnological and oceanic literature. This set of models includes different levels of complexity, ranging from empirical models integrating wind speed and lake size to process-based models including wind shear, convection, and surface waves. Continuous ΔpCO2 measurements by in situ automated sensors and CO2 fluxes, obtained from a new generation of automated (forced diffusion) flux chamber, were collected during specific periods of intensive field survey covering a wide range of natural conditions. Empirical k values computed from chamber data are then compared to outputs from the different k models. Owing to the size of Lake Geneva, we anticipate that models accounting, implicitly or explicitly, for the four key exchange drivers (i.e. wind shear, convective mixing, wave action, and bubble formation) will show the highest accuracy and precision in their estimation of k and that a precise integration of surface wave effects in such a large system should enhance model predictions. Thereafter, the relative distribution of these components is computed over a full year and analysed in the scope of the temporal variability of the gas transfer velocity. Finally, we expect that extreme wind and associated wave events should contribute disproportionately to accumulated k values over the year. In such a case, episodic weather events could generate large CO2 fluxes over very short timescales that should be accounted for when computing annual CO2 emission budgets.

Figure 2Location and map of Lake Geneva with the two prevailing winds (left) also depicted by the wind rose (top right). The wave rose highlights the highest wave field generated at the sampling location by the southwest wind with a larger fetch (bottom right). Both the wind and wave roses are computed with annual data from 13 June 2019 to 12 June 2020 at LéXPLORE.

2 Material and methods

2.1 Study site

Lake Geneva is a peri-Alpine lake defining part of the Swiss–French border, at 372 ma.s.l. (metres above sea level) (4626 N, 633 E). Its surface area (582 km2) and its maximum depth (309 m) make it the largest freshwater body in western Europe, with a volume of 89 km3 (Fig. 2). Lake Geneva is monomictic. The two prevailing winds are almost diametrically opposed and come from the southwest and northeast respectively (Fig. 2). The lake water has been surveyed monthly or fortnightly since the late 1950s (OLA-IS, AnaEE-France, INRAE of Thonon-les-Bains, and CIPEL; Rimet et al., 2020). The surface CO2 concentrations, as computed from the routine temperature, alkalinity, and pH measurements (Stumm and Morgan, 1981), show a typical seasonal cycle with high, supersaturated values during winter mixing and values below saturation in summer (Perga et al., 2016).

2.2 Field data at LéXPLORE

All field data were collected from the LéXPLORE platform, a 10 m by 10 m pontoon equipped with high-tech instrumentation and installed on Lake Geneva in 2019 (Wüest et al., 2021). LéXPLORE is moored at a 110 m depth, 570 m off the northern lake shore (Fig. 2).

On LéXPLORE, local weather conditions (air temperature, wind speed and direction, relative humidity, short-wave radiation, and atmospheric pressure) were continuously recorded (at 10 min intervals) by a Campbell Scientific automatic weather station. Lake surface temperature was measured every minute at 50 cm depth using a Minilog II-T (VEMCO, resolution 0.01 C). The partial pressure of water surface CO2 (pCO2) was also measured at 50 cm depth during specific surveys (see Sect. 2.3) using a mini CO2 sensor (Pro-Oceanus Systems Inc.) with an accuracy of ± 30 ppm. Values of pCO2 in parts per million (ppm) were converted into micro-atmospheres (µatm) following the basic equation correcting for altitude (Russell and Denn, 1972). Therefore, we assume that the concentration and the temperature are homogeneous over the first 50 cm.

Fetch distance (m) from LéXPLORE to the lake shores considering wind direction was computed using data from the Federal Office of Topography online portal (Swisstopo geoportal:, last access: 14 October 2020). The position of LéXPLORE is particularly relevant for this study as the fetch ranges from  0.5 to  30 km for the two prevailing winds. Significant wave height, Hs (in m), was computed after Hasselmann et al. (1973) according to the following equation:

(3) H s = 1.6 × 10 - 3 U 10 ( Fetch / g ) 1 / 2 ,

where g is the gravitational constant. This variable Hs is defined as the average height of the highest one-third of the waves (crest to trough) corresponding to the thickness over which the wind can push laterally (Wüest and Lorke, 2003). This equation is equivalent to the formulation by Carter (1982) that is more widely used in the oceanic literature. Simon (1997) tested the model for significant wave heights in Lake Neuchâtel (a lake close to Lake Geneva) with a fetch distance of 9 km. These results showed that the significant wave height in this lake was consistent with this oceanic formulation. However, Simon (1997) highlighted that the Joint North Sea Wave Project (JONSWAP) wave breaking parameterisation did not hold for winds greater than 5 m s−1, producing faster wave breaking and with a higher probability in the case of surface waves that were not fully developed. Such lake waves are characterised by steeper slopes that favour their wave breaking and wave action (Wüest and Lorke, 2003).

The net CO2 flux at the lake–atmosphere interface, F, was directly measured with an automated (forced diffusion) floating CO2 flux chamber (eosFD, Eosense: environmental gas monitoring; Fig. A1; Risk et al., 2011), which was originally developed for soil flux studies. The flux chamber had a detection limit close to 0.05 µmolcm-2h-1 and measured F every 15 min in summer and 30 min in winter for battery-saving purposes. The standard floating chambers require quiet surface conditions (e.g. Cole et al., 2010; Vachon et al., 2010; Bastviken et al., 2015), thereby limiting studies to low to moderate wind speed conditions. One typical problem with floating chambers arises from the possible atmospheric leakage under rough surface conditions (Fig. A2a). To work around this problem, Vachon et al. (2010) recommend the creation of 10 cm long-edges entering the water (Fig. A2b), and this design also reduces artificial turbulence generated by the chamber's walls at surface. A second typical issue with this method is potential flux enhancement by artificial (chamber-generated) turbulence. This was also studied by Vachon et al. (2010), who demonstrated that the overestimations due to this effect can be as high as 1000 % at low wind but less than 50 % when the wind speed exceeds 4 m s−1 in large lakes. At even higher wind speed, this overestimation should decrease further because the surface water turbulence becomes much greater than that produced by the floating chamber. Thus, our flux chamber was specifically conceived to increase stability under calm and windy conditions and to limit artificial turbulence, but we do not exclude a bias under low and moderate wind conditions (Figs. A1, A2c).

Regarding the operation of the eosFD, it has two independent cavities: one for the chamber and one for the atmosphere (Fig. A1). These are connected to the same CO2 sensor by a pump which sends either the chamber gas or the air gas to the sensor at regular intervals (about 20 s) and then completely flushes the chamber cavity according to the programmed measurement time step (15 or 30 min). Therefore, the advantage of this new instrument is to have a constant monitoring of the chamber's variation but also of the atmosphere. In addition, the use of the same CO2 sensor for the two measurements limits the need for intercalibration between CO2 sensors. We tested the performance of the floating chamber by comparing the standard deviation of the CO2 concentrations of the atmosphere and in the chamber estimated from two separated cavities (Fig. A1; Risk et al., 2011). We did not observe any difference in the standard deviation between high and low wind conditions (Fig. A3), suggesting that the measured fluxes remained reliable at high wind speed without leakage of the chamber.

We assessed the performances of our flux chamber during five specific periods over the annual cycle: 13–14 June 2019, 27–28 August 2019, 1–5 October 2019, 18–20 December 2019, and 20–26 February 2020. To select the most robust dataset for comparison with k estimates derived from models, we discarded flux data that were below the detection limit as well as CO2 gradients that ranged within the uncertainty of sensors (i.e. ± 20 ppm for air and ± 30 ppm for water, resulting in a ± 50 ppm uncertainty) (Fig. A2). Accordingly, we were able to retain the most robust data points during the following deployment periods: 18–20 December 2019 and 20–26 February 2020. Finally, all these field data were standardised at a 1 h time step.

2.3 Computed k values from field data

k values (cm h−1) from field observations (kobs) were computed from the gas transfer velocity equation:

(4) k obs = F / ( α Δ p CO 2 ) ,

where F is the measured CO2 flux (µmolcm-2h-1); α is the gas solubility coefficient (µmolcm-3µatm-1), which depends on the measured water temperature (Wanninkhof, 1992); and ΔpCO2 is the differential of pCO2 measured at 0.5 m below the surface (pCO2-water) and pCO2 at saturation (pCO2-sat; in ppm) measured from the flux chamber corrected by altitude (µatm). It is noteworthy that the chemical enhancement factor (Wanninkhof and Knox, 1996) was not considered in this equation, as the fluxes retained corresponded to conditions of moderate pH (i.e. < 8). kobs was then standardised in k600 using the dimensionless Schmidt number (Sc) of CO2: k600-obs=kobs(600/Sc)-1/2 (600 for freshwater standardised at 20 C).

2.4 Models for air–water gas transfer velocity

After years of debate, a consensus has begun to emerge on the relationship linking k, the intensity of turbulence, and Sc (Eq. 2), even when starting from different physical assumptions (see Katul et al., 2018). In this study, we selected six parameterisations widely used in limnology and oceanography, combining specific calibration characteristics (Table 1). We first show that they can all be expressed following Eq. (2) for wind shear and convection, despite their different formulations. We then develop the effects of surface waves from oceanic models and adapt the wave action for a large lake. The final lake model integrating the wave effect is ultimately calibrated using our field data (Table 1).

2.4.1 Wind shear stress

We start with the case where near-surface dynamics are driven by a weak to moderate wind, in the absence of heat exchange. In this case, the contribution of surface waves can be neglected and the wind stress (τ0=ρairC10U102, where ρair is air density and C10 is the drag coefficient at 10 m) is equal to the tangential shear stress (τt=ρwatu,wat2, where ρwat is water density). The relationship between ε and the sheared velocity on the water side, u,wat, is then derived from a law-of-the-wall scaling for the velocity profile: ε=u,wat3/κz(0), where κ is the von Kármán constant (= 0.41), and z(0) is the thickness of the diffusive boundary layer. This relationship leads to

(5) k NB = a 1 ν u , wat 3 / κ z ( 0 ) 1 / 4 Sc - 1 / 2 .

The challenge is then to define z(0). Tedford et al. (2014) followed an ad hoc observational approach and chose z(0)= 0.15 m as the shallower depth where ε was measured. In contrast, theoretical studies have linked z(0) to the thickness of the diffusive or viscous sublayer ( 0.1–1 cm). In line with theory, we scale this layer as z(0)=cν/u,wat (Wüest and Lorke, 2003; Lorke and Peeters, 2006) with c as a constant value. Taking c= 114 (Soloviev et al., 2007), the thickness of this layer typically ranges from 0.04 to 0.14 m under a wind regime of 10 to 1 m s−1. Therefore, we modify Eq. (5) to compute the interfacial (no bubble, NB) exchange coefficient:

(6a) k NB = a 1 u , wat ( 1 / κ c ) 1 / 4 Sc - 1 / 2 ,


(6b) k NB = a 1 ( ρ air / ρ wat ) C 10 U 10 ( 1 / κ c ) 1 / 4 Sc - 1 / 2 .

These equations show that the SRM formulation (Table 1, Fig. B1; Soloviev et al., 2007; Read et al., 2012) is analogous to the Coupled Ocean–Atmosphere Response Experiment Gas transfer algorithm (COAREG) flux (Fairall et al., 2011) and to the formulation used in Deike and Melville (2018) with the sheared velocity on the atmosphere side, u,atm (Table 1: DM18). Indeed, when equating the expression by Deike and Melville (2018),

(7) k NB = A NB u , atm ( Sc / 600 ) - 1 / 2 = A NB u , wat ( ρ wat / ρ atm ) 1 / 2 ( Sc / 600 ) - 1 / 2 .

With Eq. (6a), we find that the coefficient a1= 0.29 in Soloviev et al. (2007) and Read et al. (2012) is essentially equivalent to the coefficient ANB=a1(κc)-1/4(1/600)1/2(ρwat/ρair)-1/21.5×10-4 in Deike and Melville (2018) (Fig. B1), which, in turn, was found to be equal to the coefficient of A= 1.5 in Fairall et al. (2011). These results agree with Lorke and Peeters (2006), who derived a unified relation for interfacial fluxes (air–water and water–sediment) through a linear relationship of u*, wat with k, especially at the bottom interface where shear is the only relevant process. Furthermore, Eq. (6) has a similar (i.e. quasi-linear) k–wind relationship to the data-driven parameterisation from VP13 but cannot explain the higher-order polynomial relationship reported in CW03 and CC98.

2.4.2 Convection

A second source of dissipation at the surface is the convection (εc) resulting from surface cooling. In the SRM formulation, only the negative buoyancy flux is considered when this term directly enters into the turbulent kinetic equation as a production term. The combination of wind shear and free convection near a boundary is described by the Monin–Obukhov similarity theory (MOST) with a general form derived from a turbulent kinetic energy balance (Lombardo and Gregg, 1989; Tedford et al. 2014):

(8) ε ( z ) = ε u ( z ) c u + c c z L MO ,

where LMO is the Monin–Obukhov length scale defined as LMO=u*, wat3/κB0, including ε(z)=cuεu+ccεc in Eq. (2). The latter expression can be rearranged as follows:

(9) k NB = a 1 ε u ( c u + c c B 0 / ε u ) 1 / 4 Sc - 1 / 2 ,

where a1 ranges in the literature from 0.2 to 1.2 (Soloviev et al., 2007; MacIntyre et al., 2010; Tedford et al., 2014; Heiskanen et al. 2014; Winslow et al., 2016), cu ranges from 0.84 to 1 (Winslow et al., 2016), and cc ranges from 0.37 to 2.5 (Wyngaard and Coté, 1971; Tedford et al., 2014). Hereafter, we use the following set of values: cu= 1 and cc= 1. Fairall et al. (2011) used an essentially equivalent approach but formulated in terms of a Richardson number to describe the partitioning between dissipation from convection and wind shear, expressing the wind shear in terms of the air-side friction velocity: Rf=B0υ/u*, atm4, which can be integrated into Eq. (9) as

(10) k NB = a 1 ρ atm ρ wat 1 / 2 u *, atm c u κ c 1 + R f R f , c 1 / 4 Sc 600 - 1 / 2 ,

where Rf,c=cuρatm2ccρwat2κc. The details of this demonstration can be seen in Soloviev and Schlüssel (1994).

2.4.3 Wave action

The effect of surface waves is commonly implemented in oceanography but barely considered in limnology. All process-based models rely on the same parameterisation of energy dissipation by wind shear and convection. However, they differ in how they parameterise energy dissipation by wave action and wave breaking.

The contribution of the wave action (εw) is, accordingly, added as a third source of turbulence (Fig. 1). In the presence of surface waves, the balance between τt and τ0 no longer holds. Therefore, Soloviev and Schlüssel (1994) added a corrective factor, φ, using the Keulegan number (Ke=u*, wat3/(gν)), in order to decrease the component τt=τ0φ (where φ=1/(1+Ke/Kec)), with the critical Keulegan number (Kec) defined in Soloviev and Lukas (2006). As a result, the equation for shear-driven dissipation εu(z) is as follows:

(11) ε u ( z ) = u *, wat 4 κ c υ φ 2 .

Following this step, the turbulent kinetic energy dissipation rate from wave action (εw) is added and defined with the Keulegan number by Soloviev et al. (2007) as follows:

(12) ε w = α W 3 B S q 1 / 2 ( Ke / Ke c ) 3 / 2 ( 1 + Ke / Ke c ) 3 / 2 u *, wat g 0.062 κ C T ( 2 π A w ) 3 / 2 ρ atm ρ wat ,

where CT=(z0/Hs). z0 is the surface roughness scale from the water side, and the CT value is set as a constant at 0.6 (more details are given in Soloviev et al., 2007). This definition does not hold for closed basins because, in the case of incompletely developed waves, the dissipation of energy from wind shear transmitted to the waves is not fully redistributed in the water body (Simon, 1997). Hence, for the application in Lake Geneva, we followed Terray et al. (1996), who defined a varying CT:

(13) C T = 1.38 × 10 - 4 U 10 C p 2.66 ,

where Cp is the peak speed of the wave spectrum defined in Deike and Melville (2018) according to Toba (1972, 1978). This leads to CT 1. This allows one to increase the effect of εw (inversely proportional to CT in Eq. 12) on k. Here, we used this formulation to adapt the S07 ocean model for a large lake (closed basin). Henceforth we refer to the adapted uncalibrated and calibrated models as SD21 and SD21-fit respectively (Table 1). Finally, these three terms of ε (εu, εc, and εw) can be added before computing the SRM (Eq. 2) for determining k-no bubble (kNB).

2.4.4 Bubble enhancement

Additional deviations from the linear relationship to U10 are explained by the gas transfer resulting from bubbles and sprays during wave breaking. This mechanism is accounted for by adding a k-bubble (kB) term to the previously mentioned kNB. Soloviev et al. (2007) used the empirical k-bubble parameterisation from Woolf et al. (1997):

(14) k B = W 2450 O s 1 + 1 14 O s Sc - 0.5 1 / 1.2 1 / 2 ,

where W is the fractional whitecap coverage only expressed as a function of wind (3.84×10-6U103.41), and Os is the Ostwald gas solubility. This formulation does not take wave height into account (Fig. B1). Nevertheless, a recent study (Deike and Melville, 2018) performed a new numerical process-based parameterisation for gas transfer velocity from bubble enhancement considering Hs through the following equation:

(15) k B = A B O s u , atm 5 / 3 ( g H s ) 2 / 3 Sc 600 - 1 / 2 ,

where AB is an empirical factor with dimension (= 10−5m−2 s2), and Os is defined by the ideal gas constant (R), the surface water temperature (T0), and the CO2 solubility coefficient in freshwater (α) (Reichl and Deike, 2020). The gas transfer velocity is the expressed as a sum of the no bubble kNB and bubble kB components (Table 1: S07, DM18, SD21, and SD21-fit) following Keeling (1993) and Woolf et al. (1997, 2005). Our adapted model for lake includes a refined parameterisation of the wave action term εw from S07 along with the bubble term from DM18. For these reasons, the model will be called SD21 in the rest of the paper. In addition, for the appellation SD21-fit, the a1 parameter of Eq. (2) and the AB parameter of Eq. (15) were fitted to the k600 observations (a1= 0.33 and AB= 3 × 10−5m−2 s2).

With this review of existing and adapting parameterisations, we show that (i) there is a discrepancy between an SRM-based model with shear stress as the only energy source and empirical parameterisations with a polynomial (order > 1) wind-based relationship. Such a discrepancy is tentatively resolved by adding the effect of convection and surface waves. (ii) We further highlight that most fitting parameters from the different SRM-based models are in good agreement. (iii) We finally recall that it is possible to provide a unifying parameterisation of k with an SRM model including wind shear, wind-induced waves, and convection with only a few input parameters such as U10, B0, and Fetch.

Figure 3The annual distribution of the three main components used to compute k600 models: wind speed at 10 m (orange), buoyancy flux at the surface during cooling (blue), and significant wave height (turquoise). These survey data observed during CO2 flux measurements after quality control (+) are also shown.


3 Results

3.1 Observed and predicted k

After quality check, our dataset contains 94 discrete CO2 flux observations. We first assess the representativeness of our sampling by comparing the survey-specific and annual distributions of the three main inputs for k models: U10 (all models), B0 during convective periods (T14, S07, SD21, and SD21-fit), and Hs (S07, DM18, SD21, and SD21-fit) (Fig. 3; the temporal evolution of these three terms is shown in Fig. C1). From 13 June 2019 to 12 June 2020, the average wind speed over Lake Geneva is 2.9 m s−1 with a mode at 2.5 m s−1; very low wind speeds (< 1 m s−1) are encountered 12 % of the year, whereas high (> 5 m s−1) to very high (> 10 m s−1) wind events represent 15 % and 2 % of the year respectively. The sampling surveys covered the full annual range of U10. Average and modal values of B0 over the year are close to 0.25 × 10−7m2 s−3. However, the sampling covered only the lowest 50 % of the annual distribution and undersampled conditions of potentially strong convection. Considering that the dissipation by buoyancy flux, as parameterised in the process-based models, is already well known in the literature and that it is not the central point of our study, we posit that the undersampling of B0 is not expected to significantly affect our analysis. The predicted modal Hs value is 0.15 m over the year. Events of high Hs (> 0.4 m) represent 6 % of the year, with a maximum Hs of 1.1 m. As for U10, the surveys covered the full range of annual Hs.

Figure 4(a) k600 observed as a function of U10 and coloured according to Hs (colour bar), showing the error bars produced by the uncertainty of ΔCO2 (± 50 ppm) as well as the u-k600 linear regression (solid line; see also Fig. B1); (b) k600 wind-based models (CC98, CW03, and VP13); (c) k600 process-based models (T14, S07, DM18, SD21, and SD21-fit) computed with annual data. Observed k600 derived from CO2 flux chamber measurements is shown using the “+” symbol.


Table 2RMSE of k600 models for all wind speed (U10), U10< 5 m s−1 (i.e. LW), and U10 5 m s−1 (i.e. SW).

Download Print Version | Download XLSX

k600 values based on observations are shown in Fig. 4a with their error bars corresponding to the uncertainties of the pCO2 in air and in water (± 50 ppm). We notice that all of the measurements with a wave height > 0.4 m were observed for wind speeds > 5 m s−1, and the corresponding k600 values are located above the linear function (i.e. from a linear regression against wind shear velocity; Fig. B1) scaling k600 to u (i.e. first-order relationship). We then compare the k600 observed during the specific surveys to the values computed with all k600 models throughout the annual cycle, in relation to U10 (Fig. 4b, c). Table 2 provides the root-mean-square errors (RMSEs) for all model estimates compared to kobs during the flux surveys, (i) for the full dataset (All Wind), and split (ii) for low wind (< 5 m s−1, LW) and (iii) strong wind conditions ( 5 m s−1, SW). The three empirical wind-based models only depend on wind (Fig. 4b). Both CC98 and CW03 were originally calibrated for small lakes, using a mass balance calibration method (Table 1). However, they lead to divergent gas transfer velocities, particularly above 5 m s−1, illustrated by a RMSE for SW as high as 22.8 cm h−1 for CC98, whereas CW03 performs better (RMSE SW = 12.8 cm h−1). Furthermore, both models underestimate k600 at low wind (Fig. 4a), with a higher deviation for CW03 (Table 2). The k values predicted by VP13 are closer to those of the process-based models that explicitly integrate wave actions (S07 and DM18) (Fig. 4b, c), demonstrating that lake size integration in the empirical model captures at least part of the wave action on k. Performances of VP13 at strong winds (RMSE SW = 12.7 cm h−1) were better than those of the ocean-derived models integrating surface waves (RMSE SW = 13–15.9 cm h−1). However, VP13 shows a positive offset during calm periods, along with the highest RMSE of the set of models at low wind speed (Fig. 4b, Table 2).

The process-based models (Fig. 4c) provide different k600 values for a given wind speed, owing to the integration of additional environmental components (i.e. the varying drag coefficient, the convective mixing in R12 and T14, and the effect of waves in S07, DM18, SD21, and SD21-fit). All process-based models are similar at low winds, as they share a common physical basis for parameterisation of wind shear and convection. Therefore, they lead to similar RMSEs (2.9–3.5 cm h−1) under such conditions, where surface waves are negligible (Table 2). Divergences occur at higher wind speeds. T14, initially developed for small lakes with limited wind exposure, performed the worst (RMSE = 19.8 cm h−1). This increased k underestimation at high winds can be attributed to (i) dissipation by wave action and bubble formation not being considered (in R12 and T14) and (ii) to the use of a constant z(0)= 0.15 m in the T14 model (Eq. 5). This approximation of the diffusive layer is consistent with low wind speed but is almost 1 order of magnitude too large under strong wind speed. Other process-based models, designed for greater wind range (> 10 m s−1), integrate surface waves and, as a result, lead to better estimates than R12 and T14 (RMSE = 10.4–15.9 cm h−1). However, the ocean wave model of DM18 shows lower performances at strong winds than CW03 and VP13 (Fig. 4b, c). Finally, the specific fit parameterisation of the SD21-fit model improves the performance at high wind speeds by  30 % (RMSE = 10.5 cm h−1), outperforming all the other methods.

Figure 5The U10 vs. k600 relationship modelled and coloured according to Hs (colour bar) in panels (a)(e) as well as coloured according to fetch distance (colour bar) in panel (f): (a) R12 integrating wind shear and convection; (b) S07 integrating wind shear, convection, and wave action for fully developed waves; (c) S21 integrating wind shear, convection, and wave action for waves that are not fully developed; (d) SD21 is similar to S21 but the k-bubble term of DM18 is added; (e, f) SD21-fit is similar to SD21 but with a1 and AB fitted to observed k.


3.2 Surface wave integration

Herein, we scrutinise how those varying parameterisations ultimately alter the shape of relationship between k600 and U10. In R12 (Fig. 5a), wind is only included through wind shear, resulting in a linear relationship between k600 and U10, as already anticipated. Adding the wave action (no bubble) through the S07 parameterisation (Fig. 5b) does not lead to any significant departure from the minimal R12 model. Adapting wave action by decreasing CT (for the no bubble term) leads to a departure from the wind shear linear relationship for Hs> 0.4 m (Fig. 5c).

Figure 6(a) Cumulative k600 modelled over an annual cycle; (b) cumulative k600 for wind > 5 m s−1; (c) cumulative k600 for wind  5 m s−1.


Adding the k-bubble term related to wave breaking of DM18 further increases this deviation from the linear k600U10 relationship but also scatters k600 estimates for a given U10 (Fig. 5d). Finally, the fitting with observationally based k600 improves the estimation for strong wind (Fig. 5e, Table 2). Given the range of wind fetch from  0.5 to  30 km, the contribution of waves varies for a given wind speed depending on the fetch, as evidenced by the scattering of the parameterised k600 for a given U10 (Fig. 5f). A significant modification of k600 by wave action and wave breaking occurs for a fetch length > 15 km and U10> 5 m s−1 (Fig. 5f), in the case of Lake Geneva, generating wave of Hs> 0.4 m (Fig. 5e).

Compared with k600 estimated by R12 (Fig. 5a), the SD21 and SD21-fit models provide k estimates that are 20 %–50 % higher for U10= 10 m s−1 respectively and 40 %–70 % higher for U10= 15 m s−1 respectively. Therefore, adapting the surface waves, through the change in the wave action for incompletely developed waves and the fitting to observed data encountered in local lake conditions, leads to better performances of the SD21 models. SD21-fit reached the lowest RMSE at all wind speeds and was thereafter used as a reference for the modelling of the annual gas transfer velocities.

3.3 Annual cumulative gas transfer velocity and the effect of extreme conditions

We show above that models accounting for surface waves better represent the non-linear increase in k600 at high winds. Because high-wind events remain rare, we test whether a better representation of k600 during rare, high-wind events affects the local estimates of k600 over a full year. To this end, cumulative sums of hourly k600 were computed over a full annual cycle (13 June 2019–12 June 2020) for all k models (Fig. 6). The annual dynamics, such as the annually averaged k600, were compared using SD21-fit as a new reference model.

Cumulated k600 computed for SD21-fit shows some episodic steep increases between December and March, due to the winterly prevalence of high-wind events (the winter average wind speed, 3.25 m s−1, was greater than the summer mean, 2.55 m s−1, by 25 %) and greater significant wave height (the winter average wave height, 0.15 m, was greater than the summer value, 0.10 m, by 50 %) (Fig. C1). The average hourly k600 by the SD21-fit model is 7.3 ± 7.4 cm h−1 (mean ± SE, Fig. 6a). Periods of high winds, although accounting for only 15 % of data points, contribute 44 % of annually cumulative k600 in the SD21 and SD21-fit models (Fig. 6b), whereas periods of high waves (Hs 0.4 m), accounting for only 6 % of data points, contribute to more than 20 % of annually cumulative k600. The wind-based models are those for which cumulative k600 diverges the most from the SD21-fit reference model, with the lowest annual averaged k600 for CC98 and CW03 (3.9 ± 2.7 and 4.8 ± 9.3 respectively) and the highest for VP13 (9.9 ± 6.1). These divergences arise from the low performances of these models at low wind regimes (Figs. 4a, 6b; Table 2), which represent 85 % of annual data points. All of the other process-based models have relative dynamics of cumulative k600 similar to that of the SD21-fit model and end up with annually averaged k600 that are 15 % lower than for the SD21-fit. The representation of k600 at low wind speeds is similar for all process-based models, and the divergence arises from the representation of the rarer high-wind-speed episodes, which contribute to 43 %–46 % of annual cumulative k600 (Fig. 6b).

4 Discussion

The history of k models, simulating the gas transfer velocity for surface waters, dates back from the early 1990s. k models have been developed and tested in small lakes sheltered from winds (e.g. Crusius and Wanninkhof, 2003; Tedford et al., 2014), large lakes under low to moderate wind speed (Vachon and Prairie, 2013), and oceans (e.g. Soloviev et al., 2007; Fairall et al., 2011; Esters et al., 2017; Deike and Melville, 2018). While the effects of surface waves on k can be neglected in small lakes, we question whether this assumption holds for large lakes such as Lake Geneva, in which surface waves are frequently observed (Figs. 2, C1). We evaluated the performance of different experimental-based and process-based models to estimate k600 in the large Lake Geneva. We show that integrating the effect of wave formation at high wind speeds and long fetch better represents the sharp increase in the k600 values during such episodic windy events.

4.1 Choice of k models

Wind-based models have long been known to misestimate k600 at low wind speeds (Eugster et al., 2003; MacIntyre et al., 2010; Erkkilä et al., 2018). Consistently, wind-based models showed the lowest performances for Lake Geneva, especially at low wind speeds (CW03 and VP13), which resulted in large discrepancies in annually averaged and cumulative k600 over the full year. They are, however, easy to compute, require few inputs (only U10), and remain by far the most used to estimate lakes' CO2 emissions worldwide (e.g. Raymond et al. 2013). Another possibility is to broadly adopt process-based models. The presented process-based models require input data that are currently more easily accessible and are routinely acquired at high frequency in many lakes: wind speed, heat flux, and wind fetch (i.e. distance from the shore). The development of R packages, such as “LakeMetabolizer” (Winslow et al., 2016), in which the calculations of process-based models are implemented also alleviates their computational difficulty. Both increased data availability and computational tools should foster the use of process-based k models, which hold great potential to obtain more accurate global k600 estimates.

The analysis of the models adapted from the existing literature to account for the effect of lake surface waves, SD21 and SD21-fit (Fig. 5d, e, f), showed that the wave contribution to k becomes significant for Hs> 0.4 m, corresponding, for Lake Geneva, to winds blowing at 5 m s−1 from the southwest where the fetch length is maximal (> 15 km) with respect to the measurement site. A significant contribution to the gas transfer velocity by surface waves is expected in lakes where Hs> 0.4 m is not infrequent. Wave heights beyond this threshold value of Hs are frequently encountered in lakes that are larger than or a similar size to Lake Geneva (6 % of annual time at the LéXPLORE platform). In the Great Lakes of North America, Hubertz et al. (1991) showed that the mean wave heights of all of these lakes were > 0.4 m in summer and close to 1 m in winter with a maximum of up to 5 m. Hs> 0.4 m can also form over elongated lakes of smaller size, such as smaller Swiss lakes (e.g. Lake Neuchâtel and Lake Bienne) (Amini et al., 2017). As SD21-fit is a process-based model integrating the four main processes in a mathematically coherent way, we would expect that it can be applied to lakes experiencing Hs> 0.4 m and can improve the accuracy of k estimates. Because waves can physically damage inshore and offshore infrastructure, many large lakes benefit from wave forecasts. Hs data from those forecasting systems (e.g. National Data Buoy Centre – NOAA and Wave Atlas from; Amini et al., 2017) could allow one to test whether the SD21-fit models can be applied to those lakes and whether kNB and kB through a1 and AB need to be recalibrated or fitted to the local context if flux measurement data are available, as for this study. Energy dissipation during high-wave events increases the gas transfer velocity well beyond the linear relationship derived for wind shear alone. Therefore, we expect that computed gas fluxes at the air–water interface should be significantly improved by the integration of surface waves into the k models.

4.2 Implication of four components on the annual k estimation and the annual CO2 fluxes

4.2.1 Seasonal and hourly distribution of kCO2

Converting k600 to kCO2 using the Schmidt number (Wanninkhof, 1992) highlights the importance of water temperature in gas exchange dynamics. Indeed, the seasonal distribution of the cumulative k600 is  20 % and  30 % for the warm (spring and summer) and cold (fall and winter) seasons respectively. Once the temperature effect is accounted for, this distribution increases to 26.1 % for summer and decreases to 24.9 % for winter but remains unchanged for spring and fall. While R12 only uses wind shear and convective terms, the selected process-based model (SD21-fit) allows a decomposition into the four main drivers of the gas transfer velocity, thereby paving the way to a better understanding of the implication of these processes throughout an annual cycle.

Figure 7(a) Distribution of kCO2 generated by two main processes (ku and kc) in R12 and four main processes (ku, kc, kw, and kb) in SD21-fit for each season: spring (April–May–June), summer (July–August–September), fall (October–November–December), and winter (January–February–March). The height of the bar represents the cumulative kCO2 by season for both models (R12 and SD21-fit). (b) Distribution of four k generated by wind shear, convection, wave action, and bubble enhancement (ku, kc, kw, and kb respectively) along the annual cycle. For the SD21-fit model, ku=SRM(εu), kc=SRM(εu+εc+εw)-SRM(εu+εw), kw=SRM(εu+εc+εw)-SRM(εu+εc), and kb.


Wind shear remains the dominant component of the gas exchange velocity over the different seasons (Fig. 7a). The annual contribution of surface waves (wave action and bubble formation) is limited to 9 %–10 % of the cumulative k in fall and winter. The contribution of the buoyancy flux at the surface to k is even smaller for both models (R12 and SD21-fit) at this seasonal scale. However, both the buoyancy flux and the surface waves can significantly increase k during episodic events, during which they can contribute disproportionately to k at hourly (up to 80 % for convection) and daily (up to 25 % for surface waves) timescales (Fig. 7b). Several studies have emphasised the disproportionate contribution of episodic mixing events to the annual flux, bringing CO2 back to lake surfaces such as after ice break in dimictic lakes (Karlsson et al., 2013; Finlay et al., 2019) or during fall mixing in a eutrophic deep lake (Reed et al., 2018). Process-based k models integrating both the buoyancy flux and the wind-induced waves offer the opportunity to mechanistically investigate how much these episodic events contribute to annual emissions through short-term modifications of the gas exchange velocity.

Table 3Seasonal to annual CO2 flux estimation (mmolCm-2d-1) from k models (with, denoted using “-CE” (italic font), and without (roman font) chemical enhancement), the monthly ΔCO2 average (µatm), and the models' deviation from SD21-fit. CE was only considered for seasons when pH > 8.4.

Download Print Version | Download XLSX

4.2.2 Consequences of the choice of k model on the seasonal to annual CO2 flux estimation

We produced coarse estimates of monthly CO2 fluxes with the objective of scaling the effects of wave integration at seasonal and annual scales. Monthly fluxes were computed based on k estimates at LéXPLORE from the different models at an hourly time step as well as the monthly average of water temperature and recorded pCO2 at the lake surface (OLA-IS, AnaEE-France, INRAE of Thonon-les-Bains, CIPEL; Rimet et al., 2020; Perga et al., 2016) as well as a constant pCO2 in the atmosphere (400 µatm). For months for which surface pH > 8.4,k values were computed with and without considering the chemical enhancement (CE; Wanninkhof and Knox, 1996) (Table 3). The dependency of the chemical enhancement factor on k (Wanninkhof and Knox, 1996) might generate a further uncertainty in estimated CO2 fluxes (e.g. with a greater k value being related to a lower chemical enhancement factor).

As predicted by Fick's law, the highest outgassing fluxes occur in fall and winter, when water mixing brings CO2 up to the lake surface, whereas low gas uptake fluxes occur in spring and summer, when primary production depletes surface CO2 below saturation. However, annual estimates of net CO2 outgassing vary from 14.7 to 37.1 mmolCm-2d-1 (Table 3) depending on the k model used for computation. Consistently, differences between model estimates are relatively low in summer, as both the ΔpCO2 gradient (100–200 µatm) and wave occurrence are limited. Adding CE during the spring and summer months causes an increase in influx of the order of 5 %–128 % depending on the k model used. For example, using CC98 in summer leads to a 93 % increase in the estimated influx (compared with fluxes without CE), whereas the increase is only 39 % for SD21-fit. This chemical enhancement factor deserves more attention in future studies, especially for high pH lakes and summer seasons. However, the variability introduced by the CE at the annual scale remains low (10 %) compared with that introduced by the choice of the k model. Estimated fluxes are also strongly dependent on the chosen k model in winter when both ΔpCO2 (475 µatm) and surface wave occurrence are higher (Fig.7, Table 3). Therefore, while high-wave events represent only 6 % of the total surface wave occurrence (Hs> 0.4 m), an incomplete consideration and description of their contribution may lead to an annual flux underestimation of about 20 %–25 %.

The weak contribution of convection is at odds with observations in small lakes, although not unexpected, as large lakes are exposed to stronger winds, such that wind-shear-driven εu often outpaces convectively driven εc (Read et al., 2012). However, the limited impact of the buoyancy flux on k does not rule out its contribution to CO2 exchange. Indeed, convective mixing plays a central role in the deepening of the mixed layer, allowing the export of the CO2 stored in the hypolimnion towards the surface during the cold period and, thus, controlling the pCO2 gradient (Zimmermann et al., 2019) and the observed wintertime outgassing. Altogether, both surface oversaturated CO2 concentrations (as a result of convective mixing) and wind-induced waves are more relevant in fall and wintertime for the monomictic Lake Geneva, leading to most of the annual outgassing during this season (Table 3). As for many monomictic lakes, these seasons drive most of the annual CO2 budget of Lake Geneva (Perga et al., 2016), whereas they usually correspond to those where direct measurements are the scarcest. An improved quantification of k values through SRM models including wind-induced waves should contribute to refining the overall estimation of large lakes' contribution to regional CO2 emissions.

Figure 8(a, b) Wind fields from COSMO-1 for two episodes of northeast and southwest wind directions; (c, d) wave fields on Lake Geneva considering the two prevailing winds (northeast and southwest); (e, f) gas transfer velocity from SD21-fit; (g) boxplots of the spatial variability, at the lake scale, of k values computed from CC98, R12, and SD21-fit under both meteorological conditions. Diamonds represents the spatial mean, and the cross (+) represents the k value computed from the averaged fetch distance (NE: 9.5 km; SW: 9.3 km).


4.3 Wind and wave field on Lake Geneva and their impact on the spatial integration of k600

SD21 and SD21-fit were built on the basis of a single measurement point on the lake, just as for most of the existing k models. Therefore, the question of the extrapolation of the model to the whole lake remains essential. Herein, we showcase two snapshot situations of high wind (i) to illustrate how process-based models could enable spatially resolved estimates of k values and (ii) to exemplify how much k can vary as a result of the spatially variable wave and wind fields during a single episode. Two events of high and similar wind speed (11 m s−1) but different directions (NE on 30 March 2020 at 08:00 UTC and SW on 10 February 2020 at 06:00 UTC; Fig. 8a, b) were extracted from the 0.01 hourly resolved numerical weather model of the Swiss Federal Office of Meteorology and Climatology (COSMO-1, MeteoSwiss). The two fetch distances from both prevailing wind directions (NE, SW) were then measured at each pixel of the grid (n=583), from which the wave height field was mapped (Fig. 8c, d) considering Eq. (3) and the two wind grids. The maps were qualitatively consistent with previous studies on wind waves for Lake Geneva using the spectral wave model (SWAN) for wave height (Amini et al., 2016). Spatially resolved k values were computed from the fetch and wind grids, using the wind-based model CC98, the process-based model without lake wave implementation R12, and the SD21-fit model containing the lake wave parameterisation.

Taking Hs> 0.4 m as a threshold for the significant effect of waves on k, the two prevailing winds show opposite wave responses, with long fetch and higher wave heights affecting either the northern or the southern shores for the SW and NE winds respectively. Under both conditions, more than 60 % of the full lake surface area experiences Hs values > 0.4 m, leading to k values as high as 68 cm h−1 as computed by SD21-fit (Fig. 8g). The eastern part of lake experiences the lowest wind speeds and wave heights in both situations, as a consequence of the orographic effect of the Alps surrounding the Grand Lake.

Both the range and the mean of estimated k values increase with increasing model complexity. Accounting only for wind speed, through the wind-based model CC98, leads to a spatially integrated k value of 15.7 cm h−1 (range of 1–22 cm h−1). k values computed from R12, accounting for both the wind shear and buoyancy flux, are on average 55 % greater (24.3 cm h−1) with a moderate effect on the range of spatial variability. Finally, including surface waves results in a spatially averaged k value that is more than double the k value computed from wind speed only (35.5 cm h−1), with a variability that is almost 3 times greater (range of 0–68). It is noteworthy that the spatial average of the k values computed by the wind grids (Fig. 8g, diamonds) is equivalent to the average of k values computed from the spatial mean fetch (9.5 km for NE and 9.3 km for SW; Fig. 8g, crosses) under these specific weather conditions. Thus, the application of an average fetch would be relevant to estimate a spatially averaged k value.

For all models, the integration of spatially resolved wind fields may improve the accuracy of k at the lake scale, but accounting for wind only would underestimate both the average gas piston velocity and its spatial variability. Therefore, a better understanding of wave behaviour in large lakes, using different approaches such as field and laboratory measurements, new physical models, and technical development, would improve the accuracy of gas exchange estimates at the lake–water interface, at both temporal and spatial scales. Further, estimates of lake-scale CO2 fluxes would also require one to account for the spatial variability of CO2. The question of the spatial variability of the ΔCO2 is still open and remains difficult to analyse at high frequency in large lakes.

5 Conclusions

Investigations of the four main processes generating the gas transfer velocity in the large Lake Geneva demonstrated the importance of considering surface waves during episodic windy events responsible for more than 44 % of annual cumulative k600. The in-depth study of the behaviour of the process-based models has enabled us to underscore their consistent predictions under low and strong wind conditions, especially considering the new combination and adaptation model, SD21-fit. This last model significantly improves the estimation of the CO2 flux when these three thresholds appear in the field, U10> 5 m s−1, Fetch > 15 km, and Hs> 0.4 m, making it applicable to a wide range of lake sizes. Furthermore, SD21-fit is assembled on solid theoretical bases coming from limnological and oceanic literature and allows one to analyse the distribution of these four main terms (ku, kc, kw, and kb) across a variety of timescales depending on the kind of study.

To conclude, the development of high-resolution gridded atmospheric models such as COSMO-1 is an asset for future estimates of the gas transfer velocity and the spatial heterogeneity of lake biogeochemical cycles. Moreover, this study sheds light on the complexity of large lakes located at the interface between small, sheltered lakes and the open oceans, experiencing a combination of processes relevant for both small and large systems. The possibility of using process-based models in a fairly simple way with few inputs to improve the precision of the gas transfer velocity and, therefore, the gas flux should be supported in future research. In addition, this approach is very promising with respect to uncovering long-term trends of CO2 emissions from lakes as well as for finer estimation of fluxes during more intense episodic events.

Appendix A

Figure A1Schematics of eosFD operation (, last access: 18 August 2020), its mini-platform construction, and its positioning for measurements in the field (Lake Geneva at LéXPLORE platform). The raft design also complies with recommendations to minimise artificial turbulence induced by the chamber's walls, with 10 cm long-edges entering the water (Vachon et al., 2010).

Figure A2(a) Classical floating chamber; (b) floating chamber with 10 cm long-edges; (c) platform design used in this study: 10 cm long-edges, rounded-edges, and flat and long water wings.


Figure A3Visualisation of 304 observed k600 values during the five periods of flux measurements: 13–14 June 2019, 27–28 August 2019, 1–5 October 2019, 18–20 December 2019, and 20–26 February 2020.


Figure A4(a) Raw outputs of the eosFD during one period of CO2 flux measurements, showing the ΔCO2 between both cavities of measures (atmosphere cavity and chamber cavity) (blue line); the standard deviation of each cavity between two automated flushing events (30 min of interval) – chamber cavity (red line) and atmosphere cavity (yellow line); and the CO2 flux (black dash line). (b) Temporal evolution of U10 and Hs during the same period as the CO2 flux measurements. The increase in flux on 25 February corresponds to the increase in wind speed and waves.


Appendix B

Figure B1Panel (a) shows the comparison of Soloviev et al. (2007) and Deike and Melville (2018) for the first-order function of friction velocity at the water side (u,wat) (blue points) and at the atmosphere side (u,atm) (green points) with their linear regression (black line), the linear function of Vachon and Prairie (2013) for a lake size of 582 km2 (yellow points), and the linear regression from u*,wat or u*,atm. Panel (b) is a visualisation of S07 with empirical parameterisation of the bubble term (Woolf, 1997) regardless of wave height as a function of wind speed at 10 m. Panel (c) is a visualisation of DM18 as a function of wind speed, only using the effect of the bubble term from 10 m s−1.


Appendix C

Figure C1Annual evolution of three main inputs of k models: wind speed at 10 m (U10), buoyancy flux at surface (B0), and significant wave height (Hs).


Data availability

Water temperature, buoyancy flux at the surface, and meteorological data are available from the datalakes open research data publishing platform (, Bouffard et al., 2021). Water and air CO2 concentration and CO2 flux measurements are available on Zenodo (; Perolo, 2021).

Author contributions

PP, MEP, and DB designed the study. PP, BFC, NE, and TL collected field data and carried out data pre-processing. PP developed the model code and performed the simulation with contribution from DB. PP, MEP, and DB drafted the paper and all co-authors contributed to the final submitted article.

Competing interests

The contact author has declared that neither they nor their co-authors have any competing interests.


Publisher's note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Special issue statement

This article is part of the special issue “Modelling inland waters in a changing climate (GMD/ESD/TC inter-journal SI)”. It is a result of the 6th Workshop on Parameterization of lakes in Numerical Weather Prediction and Climate Modelling, Toulouse, France, 22–24 October 2019.


We would like to thank the entire team of the LéXPLORE platform for their administrative and technical support and for the LéXPLORE core dataset. We also acknowledge the five partner institutions involved with LéXPLORE: Eawag, EPFL, the University of Geneva, the University of Lausanne, and CARRTEL (INRAE-USMB). This study was supported by the CARBOGEN project (SNF 200021_175530), which is linked to the LéXPLORE project (SNF R'Equip, P157779). The authors thank Sébastien Lavanchy, the chief technical officer (APHYS-EPFL), and Aurélien Ballu, a member of the technical pool (IDYST-UNIL) of LéXPLORE platform, for their technical and field support. Bieito Fernández Castro was supported by the European Union's Horizon 2020 Research and Innovation programme under a Marie Skłodowska-Curie action (grant agreement no. 834330; SO-CUP).

Financial support

This research has been supported by the Schweizerischer Nationalfonds zur Förderung der Wissenschaftlichen Forschung (grant no. SNF 200021_175530; CARBOGEN).

Review statement

This paper was edited by Tom Shatwell and reviewed by two anonymous referees.


Amini, A., Dhont, B., and Heller, P.: Wave atlas for Swiss lakes: modeling design waves in mountainous lakes, Journal of Applied Water Engineering and Research, 5, 103–113,, 2016. 

Bastviken, D., Sundgren, I., Natchimuthu, S., Reyier, H., and Gålfalk, M.: Technical Note: Cost-efficient approaches to measure carbon dioxide (CO2) fluxes and concentrations in terrestrial and aquatic environments using mini loggers, Biogeosciences, 12, 3849–3859,, 2015. 

Borges, A. V., Vanderborght, J.-P., Schiettecatte, L.-S., Gazeau, F., Ferron-Smith, S., Delille, B., and Frankignoulle, M.: Variability of the gas transfer velocity of CO2 in a macrotidal estuary (the Scheldt), Estuaries, 27, 593–603,, 2004. 

Bouffard, D., Sukys, J., Runnals, J., and Odermatt, D.: Datalakes: Heterogeneous data platform for operational modelling and forecasting of Swiss lakes, available at:, last access: 23 February 2021. 

Carter, D. J. T.: Prediction of wave height and period for a constant wind velocity using the JONSWAP results, Ocean Eng., 9, 17–33,, 1982. 

Cole, J. J. and Caraco, N. F.: Atmospheric exchange of carbon dioxide in a low-wind oligotrophic lake measured by the addition of SF6, Limnol. Oceanogr., 43, 647–656,, 1998. 

Cole, J. J., Prairie, Y. T., Caraco, N. F., McDowell, W. H., Tranvik, L. J., Striegl, R. G., Duarte, C. M., Kortelainen, P., Downing, J. A., Middelburg, J. J., and Melack, J.: Plumbing the Global carbon cycle: Integrating inland waters into the terrestrial carbon budget, Ecosystems, 10, 172–185,, 2007. 

Cole, J. J., Bade, D. L., Bastviken, D., Pace, M. L., and de Bogert, M. V.: Multiple approaches to estimating air–water gas exchange in small lakes, Limnol. Oceanogr.-Meth., 8, 285—293,, 2010. 

Crusius, J. and Wanninkhof, R.: Gas transfer velocities measured at low wind speed over a lake, Limnol. Oceanogr., 48, 1010–1017,, 2003. 

Danckwerts, P. V.: Significance of liquid-film coefficient in gas absorption, Ind. Eng. Chem., 43, 1460–1467,, 1951. 

Deike, L. and Melville, W. K.: Gas transfer by breaking waves, Geophys. Res. Lett., 45, 482–492,, 2018. 

Dugan, H. A., Woolway, R. I., Santoso, A. B., Corman, J. R., Jaimes, A., Nodine, E. R., Patil, V. P., Zwart, J. A., Brentrup, J. A., Hetherington, A. L., Oliver, S. K., Read, J. S., Winters, K. M., Hanson, P. C., Read, E. K., Winslow, L. A., and Weathers, K. C.: Consequences of gas flux model choice on the interpretation of metabolic balance across 15 lakes, Inland Waters, 6, 581–592,, 2016. 

Engel, F., Farrell, K. J., McCullough, I. M., Scordo, F., Denfeld, B. A., Dugan, H. A., de Eyto, E., Hanson, P. C., McClure, R. P., Nõges, P., Nõges, T., Ryder, E., Weathers, K. C., and Weyhenmeyer, G. A.: A lake classification concept for a more accurate global estimate of the dissolved inorganic carbon export from terrestrial ecosystems to inland waters, Sci. Nat.-Heidelberg, 105, 25,, 2018. 

Erkkilä, K.-M., Ojala, A., Bastviken, D., Biermann, T., Heiskanen, J. J., Lindroth, A., Peltola, O., Rantakari, M., Vesala, T., and Mammarella, I.: Methane and carbon dioxide fluxes over a lake: comparison between eddy covariance, floating chambers and boundary layer method, Biogeosciences, 15, 429–445,, 2018. 

Esters, L., Landwehr, S., Sutherland, G., Bell, T. G., Christensen, K. H., Saltzman, E. S., Miller, S. D., and Ward, B.: Parameterizing air-sea gas transfer velocity with dissipation: Dissipation-based k-parametrization, J. Geophys. Res.-Oceans, 122, 3041–3056., 2017. 

Esters, L., Rutgersson, A., Nilsson, E., and Sahlée, E.: Non-local impacts on eddy-covariance air-lake CO2 fluxes, Bound.-Lay. Meteorol., 178, 283–300,, 2021. 

Eugster, W., Kling, G. W., Jonas, T., McFadden, J. P., Wüest, A., MacIntyre, S., and Chapin, F. S.: CO2 exchange between air and water in an Artic Alaskan and midlatitude Swiss lake: Importance of convective mixing, J. Geophys. Res., 108, D12,, 2003. 

Fairall, C. W., Yang, M., Bariteau, L., Edson, J. B., Helmig, D., McGillis, W., Pezoa, S., Hare, J. E., Huebert, B., and Blomquist, B.: Implementation of the Coupled Ocean-Atmosphere Response Experiment flux algorithm with CO2, dimethyl sulfide, and O3, J. Geophys. Res., 116, C00F09,, 2011. 

Finlay, K., Vogt, R. J., Simpson, G. L., and Leavitt, P. R.: Seasonality of pCO2 in a hard-water lake of the northern Great Plains: The legacy effects of climate and limnological conditions over 36 years, Limnol. Oceanogr., 64, 118–129,, 2019. 

Frost, T. and Upstill-Goddard, R. C.: Meteorological controls of gas exchange at a small English lake, Limnol. Oceanogr., 47, 1165–1174,, 2002. 

Guérin, F., Abril, G., Serça, D., Delon, C., Richard, S., Delmas, R., Tremblay, A., and Varfalvy, L.: Gas transfer velocities of CO2 and CH4 in a tropical reservoir and its river downstream, J. Marine Syst., 66, 161–172,, 2007. 

Hasselmann K., Barnett T. P., Bouws, E., Carlson, H, Cartwright D. E., Enke, K., Ewing, J. A., Gienapp, H., Hasselmann, D. E., Kruseman, P., Meerburg, A., Müller, P., Olbers, D. J., Richter, K., Sell, W., and Walden, H.: Measurements of wind-wave growth and swell decay during the Joint North Sea Wave Project (JON- SWAP), Dtsch. Hydrog. Z. Suppl. A, 8, 1–95, 1973. 

Heiskanen, J. J., Mammarella, I., Haapanala, S., Pumpanen, J., Vesala, T., MacIntyre, S., and Ojala, A.: Effects of cooling and internal wave motions on gas transfer coefficients in a boreal lake, Tellus B, 66, 22827,, 2014. 

Hubertz, J. M., Driver, D. B., and Reinhard, R. D.: Wind waves on the Great Lakes: A 32 year hindcast, J. Coastal Res., 7, 945–967, 1991. 

Huotari, J., Ojala, A., Peltomaa, E., Nordbo, A., Launiainen, S., Pumpanen, J., Rasilo, T., Hari, P., and Vesala, T.: Long-term direct CO2 flux measurements over a boreal lake: Five years of eddy covariance data, Geophys. Res. Lett., 38, L18401,, 2011. 

Karlsson, J., Giesler, R., Persson, J., and Lundin, E.: High emission of carbon dioxide and methane during ice thaw in high latitude lakes, Geophys. Res. Lett., 40, 1123–1127,, 2013. 

Katul, G. and Liu, H.: Multiple mechanisms generate a universal scaling with dissipation for the air–water gas transfer velocity, Geophys. Res. Lett., 44, 1892–1898,, 2017. 

Katul, G., Mammarella, I., Grönholm, T., and Vesala, T.: A structure function model recovers the many formulations for air–water gas transfer velocity, Water Resour. Res., 54, 5905–5920,, 2018. 

Keeling, R. F., Najjar, R. P., Bender, M. L., and Tans, P. P.: What atmospheric oxygen measurements can tell us about the global carbon cycle, Global Biogeochem. Cy., 7, 37–67,, 1993. 

Klaus, M. and Vachon, D.: Challenges of predicting gas transfer velocity from wind measurements over global lakes, Aquat. Sci., 82, 53,, 2020. 

Lamont, J. C. and Scott, D. S.: An eddy cell model of mass transfer into the surface of a turbulent liquid, AIChE J., 16, 513–519,, 1970. 

Lombardo, C. P. and Gregg, M. C.: Similarity scaling of viscous and thermal dissipation in a convecting surface boundary layer. J. Geophys. Res., 94, 6273–6284,, 1989. 

Lorke, A. and Peeters, F.: Toward a unified scaling relation for interfacial fluxes, J. Phys. Oceanogr., 36, 955–961,, 2006. 

Maberly, S. C., Barker, P. A., Stott, A. W., and De Ville, M. M.: Catchment productivity controls CO2 emissions from lakes, Nat. Clim. Change, 3, 391–394,, 2013. 

MacIntyre, S., Eugster, W., and Kling, G. W.: The critical importance of buoyancy flux for gas flux across the air–water interface, in: Geophysical Monograph Series, edited by: Donelan, M. A., Drennan, W. M., Saltzman, E. S., and Wanninkhof, R., American Geophysical Union, Washington, DC,, pp. 135–139, 2001. 

MacIntyre, S., Jonsson, A., Jansson, M., Aberg, J., Turney, D. E., and Miller, S. D.: Buoyancy flux, turbulence, and the gas transfer coefficient in a stratified lake: Turbulence and gas evasion in lakes, Geophys. Res. Lett., 37, L24604,, 2010. 

Mammarella, I., Nordbo, A., Rannik, Ü., Haapanala, S., Levula, J., Laakso, H., Ojala, A., Peltola, O., Heiskanen, J., Pumpanen, J., and Vesala, T.: Carbon dioxide and energy fluxes over a small boreal lake in Southern Finland: CO2 and Energy Fluxes Over Lake, J. Geophys. Res.-Biogeo., 120, 1296–1314,, 2015. 

Perga, M.-E., Maberly, S. C., Jenny, J.-P., Alric, B., Pignol, C., and Naffrechoux, E.: A century of human-driven changes in the carbon dioxide concentration of lakes: 150 years of human impacts on lakes CO2, Global Biogeochem. Cy., 30, 93–104,, 2016. 

Perolo, P.: CO2 flux measurements in Lake Geneva, Zenodo [dataset],, 2021. 

Raymond, P. A., Hartmann, J., Lauerwald, R., Sobek, S., McDonald, C., Hoover, M., Butman, D., Striegl, R., Mayorga, E., Humborg, C., Kortelainen, P., Dürr, H., Meybeck, M., Ciais, P., and Guth, P.: Global carbon dioxide emissions from inland waters, Nature, 503, 355–359,, 2013. 

Read, J. S., Hamilton, D. P., Desai, A. R., Rose, K. C., MacIntyre, S., Lenters, J. D., Smyth, R. L., Hanson, P. C., Cole, J. J., Staehr, P. A., Rusak, J. A., Pierson, D. C., Brookes, J. D., Laas, A., and Wu, C. H.: Lake-size dependency of wind shear and convection as controls on gas exchange: Lake-size dependency of u and w, Geophys. Res. Lett., 39, L09405,, 2012. 

Reed, D. E., Dugan, H. A., Flannery, A. L., and Desai, A. R.: Carbon sink and source dynamics of a eutrophic deep lake using multiple flux observations over multiple years: Carbon sink and source dynamics, Limnol. Oceanogr. Lett., 3, 285–292,, 2018. 

Reichl, B. G. and Deike, L.: Contribution of sea-state dependent bubbles to air-sea carbon dioxide fluxes, Geophys. Res. Lett., 47, L087267,, 2020. 

Rimet, F., Anneville, O., Barbet, D., Chardon, C., Crépin, L., Domaizon, I., and Monet, G.: The Observatory on LAkes (OLA) database: Sixty years of environmental data accessible to the public, J. Limnol., 78, 164–178., 2020. 

Risk, D., Nickerson, N., Creelman, C., McArthur, G., and Owens, J.: Forced Diffusion soil flux: A new technique for continuous monitoring of soil gas efflux, Agr. Forest Meteorol., 151, 1622–1631,, 2011. 

Russell, T. W. F. and Denn, M. M.: Introduction to chemical engineering analysis, Wiley, New York, USA, 1972. 

Schilder, J., Bastviken, D., van Hardenbroek, M., Kankaala, P., Rinta, P., Stötter, T., and Heiri, O.: Spatial heterogeneity and lake morphology affect diffusive greenhouse gas emission estimates of lakes: Spatial heterogeneity of diffusive flux, Geophys. Res. Lett., 40, 5752–5756,, 2013. 

Schubert, M., Paschke, A., Lieberman, E., and Burnett, W. C.: Air–Water Partitioning of 222Rn and its Dependence on Water Temperature and Salinity, Environ. Sci. Technol., 46, 3905–3911,, 2012. 

Simon A.: Turbulent mixing in the surface boundary layer of lakes, PhD thesis no. 12,272, Swiss Fed. Inst. Technol. (ETH), Zurich, 1997. 

Soloviev, A. and Lukas, R.: The Near-Surface Layer of the Ocean: Structure, Dynamics, and Applications, Springer, Dordrecht, NL, 572 pp.,, 2006. 

Soloviev, A. and Schlüssel P.: Parametrization of the cool skin of the ocean and of the air-ocean gas transfer on the basis of modeling surface renewal, J. Phys. Oceanogr., 24, 1339–1346, 1994. 

Soloviev, A., Donelan, M., Graber, H., Haus, B., and Schlüssel, P.: An approach to estimation of near-surface turbulence and CO2 transfer velocity from remote sensing data, J. Marine Syst., 66, 182–194,, 2007. 

Stumm, W. and Morgan, J. J.: Aquatic Chemistry: An introduction emphasizing chemical equilibria in natural waters, 2nd edn., John Wiley and Sons Ltd, New York, USA, 1981. 

Tedford, E. W., MacIntyre, S., Miller, S. D., and Czikowsky, M. J.: Similarity scaling of turbulence in a temperate lake during fall cooling, J. Geophys. Res.-Oceans, 119, 4689–4713,, 2014. 

Terray, E. A., Donelan, M. A., Agrawal, Y. C., Drennan, W. M., Kahma, K. K., Williams III, A. J., Hwang, P. A., and Kitaigorodkii, S. A.: Estimates of kinetic energy dissipation under breaking waves, J. Phys. Oceanogr., 26, 792–807, 1996. 

Toba, Y.: Local balance in the air-sea boundary processes, J. Oceanogr., 28, 109–120,, 1972. 

Toba, Y.: Stochastic form of the growth of wind waves in a single-parameter representation with physical implications, J. Phys. Oceanogr., 8, 494–507,<0494:SFOTGO>2.0.CO;2, 1978. 

Tranvik, L. J., Downing, J. A., Cotner, J. B., Loiselle, S. A., Striegl, R. G., Ballatore, T. J., Dillon, P., Finlay, K., Fortino, K., Knoll, L. B., Kortelainen, P. L., Kutser, T., Larsen, Soren., Laurion, I., Leech, D. M., McCallister, S. L., McKnight, D. M., Melack, J. M., Overholt, E., Porter, J. A., Prairie, Y., Renwick, W. H., Roland, F., Sherman, B. S., Schindler, D. W., Sobek, S., Tremblay, A., Vanni, M. J., Verschoor, A. M., von Wachenfeldt, E., and Weyhenmeyer, G. A.: Lakes and reservoirs as regulators of carbon cycling and climate, Limnol. Oceanogr., 54, 2298–2314,, 2009. 

Vachon, D. and Prairie, Y. T.: The ecosystem size and shape dependence of gas transfer velocity versus wind speed relationships in lakes, Can. J. Fish. Aquat. Sci., 70, 1757–1764,, 2013. 

Vachon, D., Prairie, Y. T., and Cole, J. J.: The relationship between near-surface turbulence and gas transfer velocity in freshwater systems and its implications for floating chamber measurements of gas exchange, Limnol. Oceanogr., 55, 1723–1732,, 2010.  

Vesala, T., Huotari, J., Rannik, Ü., Suni, T., Smolander, S., Sogachev, A., Launiainen, S., and Ojala, A.: Eddy covariance measurements of carbon exchange and latent and sensible heat fluxes over a boreal lake for a full open-water period, J. Geophys. Res., 111, D11101,, 2006. 

Wanninkhof, R.: Relationship between wind speed and gas exchange over the ocean, J. Geophys. Res., 97, 7373–7382,, 1992. 

Wanninkhof, R. and Knox, M.: Chemical enhancement of CO2 exchange in natural waters, Limnol. Oceanogr., 41, 689–697,, 1996. 

Winslow, L. A., Zwart, J. A., Batt, R. D., Dugan, H. A., Woolway, R. I., Corman, J. R., Hanson, P. C., and Read, J. S.: LakeMetabolizer: An R package for estimating lake metabolism from free-water oxygen using diverse statistical models, Inland Waters, 6, 622–636,, 2016. 

Woolf, D. K.: Bubbles and their role in gas exchange, in: The Sea Surface and Global Change, edited by: Liss, P. S. and Duce, R. A., Cambridge University Press, Cambridge, UK,, pp. 173–206, 1997. 

Woolf, D. K.: Parametrization of gas transfer velocities and sea-state-dependent wave breaking, Tellus B, 57, 87–94,, 2005. 

Wüest, A. and Lorke, A.: Small-scale hydrodynamics in lakes, Annu. Rev. Fluid Mech., 35, 373–412,, 2003. 

Wüest, A., Bouffard, D., Guillard, J., Ibelings, B. W., Lavanchy, S., Perga, M.-E., and Pasche, N.: LéXPLORE: A floating laboratory on Lake Geneva offering unique research opportunities, Wires Water, 8, e1544,, 2021. 

Wyngaard, J. C. and Coté O. R.: The budgets of turbulent kinetic energy and temperature variance in the atmospheric surface layer, J. Atmos. Sci., 28, 190–201, 1971. 

Zappa, C. J., McGillis, W. R., Raymond, P. A., Edson, J. B., Hintsa, E. J., Zemmelink, H. J., Dacey, J. W. H., and Ho, D. T.: Environmental turbulent mixing controls on air–water gas exchange in marine and aquatic systems, Geophys. Res. Lett., 34, L10601,, 2007. 

Zimmermann, M., Mayr, M. J., Bouffard, D., Eugster, W., Steinsberger, T., Wehrli, B., Brand, A., and Bürgmann, H.: Lake overturn as a key driver for methan oxidation, CSH. Lab. bioRXiv,, 2019. 

Short summary
Wind blowing over the ocean creates waves that, by increasing the level of turbulence, promote gas exchange at the air–water interface. In this study, for the first time, we measured enhanced gas exchanges by wind-induced waves at the surface of a large lake. We adapted an ocean-based model to account for the effect of surface waves on gas exchange in lakes. We finally show that intense wind events with surface waves contribute disproportionately to the annual CO2 gas flux in a large lake.
Final-revised paper