ADVANCES IN LAND SURFACE HYDROLOGY REPRESENTATION IN INM RAS EARTH - - PowerPoint PPT Presentation

ADVANCES IN LAND SURFACE HYDROLOGY REPRESENTATION IN INM RAS EARTH SYSTEM MODEL

Victor Stepanenko1

1Lomonosov Moscow State University, Research Computing Center

Contributors: V.Yu.Bogomolov, V.N.Lykossov, E.M.Volodin, I.Mammarella, H.Miettinen, A.Ojala, T.Vesala, S.P.Guseva, A.Medvedev

International Young Scientists School and Conference, Zvenigorod, 5 September 2017

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 1 / 54

Пример: снегопады над Великими Американскими озерами (lake-effect snow)

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 2 / 54

При холодных вторжениях континентального воздуха интенсивное испарение и конвекция приводят к образованию облачности и осадков. "Озерные снегопады" парализуют дорожную ситуацию, закрываются школы, отменяются полеты и т.д. В течение XX в. наблюдается тренд увеличение суммы снежных осадков в данном районе, +1.9 см год−1.

Example: convection over Great African Lakes

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 3 / 54

Nocturnal convection over Victoria accounts for annual fishers death toll ∼ 5000.

Thiery et al. 2015, J. of Climate, DOI: 10.1175/JCLI-D-14-00565.1

Example: cloudiness over the Ladoga Lake

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 4 / 54

Ice-free lake evaporates, and resulting stratiform clouds are advected to Finland. Cloudiness increases the surface net radiation, and 2m-temperature rises by 15-20◦C

Eerola et al. Tellus A 2014, 66, 23929, http://dx.doi.org/10.3402/tellusa.v66.23929

Freshwaters in global carbon cycle

Total freshwater methane emission is 104 Tg yr−1, i.e. 50% of global wetland emission (177-284 Tg yr−1, IPCC, 2013) greenhouse warming potentials from freshwater-originating CO2 and CH4 are roughly equal

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 5 / 54

(Tranvik et al. 2009) (Bastviken et al. 2011)

CO2 emissions by lakes and rivers

Raymond et al., 2013, Nature

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 6 / 54

Водоемы Водотоки global emission of CO2 by freshwaters is 2.1 Pg C yr−1 lake emission is 0.3 Pg C yr−1, river emissions is 1.8 Pg C yr−1 significant contribution of Volga hydropower reservoirs

Эмиссия парниковых газов из водохранилищ

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 7 / 54

Затопленные экосистемы подвергаются длительному разложению в преимущественно анаэробных условиях В отличие от естественных водоемов, имеется дополнительный путь для эмиссии метана в атмосферу – через турбины

Global warming of lakes

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 8 / 54

The majority of lakes are warming at a rate higher than T2m.

O’Reilly et al., 2015, GRL, doi:10.1002/ 2015GL066235

1D lake model framework

1D equations result from boundary-layer approximation

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 9 / 54

1D heat and momentum equations k − ǫ turbulence closure Monin-Obukhov similarity for surface fluxes Beer-Lambert law for shortwave radiation attenuation Momentum flux partitioning between wave development and currents (Stepanenko et al., 2014) Soil heat and moisture transfer including phase transitions Multilayer snow and ice models 1D concept does not suffice the greenhouse gas modeling task, as it does not take into account differences between CH4 & CO2 emissions at deep and shallow sediments

k − ǫ turbulence closure

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 10 / 54

1D+ framework

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 11 / 54

Traditional 1D model concept 1D+ model concept 1D+ model includes friction, heat and mass exchange at the lateral boundaries Heat, moisture and gas transfer are solved for each soil column independently

In 1D+ model horizontally averaged quantity f obeys the equation: ∂f ∂t = 1 A ∂ ∂z A kf ∂f ∂z + F(z, t, f, A ) + Hf 1 A dA dz .

Coupling heat transport in water and soil

Lake body

zs0 zs0

Soil column 5

Ts1 Fs1

Soil column 1

zs1 Ts2 Fs2

Soil column 2

zs2Ts3 Fs3

Soil column 3

zs3Ts4 Fs4

Soil column 4

zs4 Ts1 Fs1

Soil column 1

zs1 Ts2 Fs2

Soil column 2

zs2 Ts3 Fs3

Soil column 3

zs3 Ts4 Fs4

Soil column 4

zs4

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 12 / 54

Boundary conditions: at soil-water interface

Continuity of temperature (gas) Continuity of flux

1D equations for enclosed basins

Horizontally-averaged 3D equations for basic prognostic quantities:

cwρw ∂T ∂t = · · · 1 A ∂ ∂z

• A (λm + cwρwνT ) ∂T

∂z

• −

− 1 A ∂AS ∂z + 1 A dA dz [Sb + FT,b(z)], – heat conservation equation (1) ∂u ∂t = · · · −

• 1

ρw ∂p ∂x

• + 1

A ∂ ∂z

• A(ν + νm)∂u

∂z

• +

+ 1 A dA dz Fu,b(z) + fv, – momentum equation for x-speed component (2) ∂v ∂t = · · · −

• 1

ρw ∂p ∂y

• + 1

A ∂ ∂z

• A(ν + νm)∂v

∂z

• +

+ 1 A dA dz Fv,b(z) − fu – momentum equation for y-speed component (3)

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 13 / 54

Barotropic pressure gradient and seiches

Lake surface

Lx,0 Ly,0 v u y x

Vertical cross-section, y = 0

Wind hx1 hx2 u z x

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 14 / 54

Mass conservation

dhN

dt A0(t) = − dhS dt A0(t) = 2 1 0 vLW −Ehdξ, dhE dt A0(t) = − dhW dt A0(t) = 2 1 0 uLS−Nhdξ,

Barotropic pressure gradient force

• g ∂hs

∂x ≈ gπ2 4 hE−hW LW −E,0 ,

g ∂hs

∂y ≈ gπ2 4 hN −hS LS−N,0 .

Barotropic (surface) seiches are lake surface and related velocity oscillations after strong wind events.

Surface oscillations in the model Turbulent kinetic energy profile (modeled), June 2013, Kuivajarvi Lake, seiches produce TKE near bottom

Biogeochemical processes in the model

O2 CO2 CH4 Biochemical

• xygen

demand (BOD) Sedimentary

• xygen

demand (SOD) Photosynthesis Respiration Methane

• xidation

Methane production Turbulent diffusion Bubble transport

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 15 / 54

Photosynthesis, respiration and BOD are empirical functions of temperature and Chl-a (Stefan and Fang, 1994) Oxygen uptake by sediments (SOD) is controlled by O2 concentration and temperature (Walker and Snodrgass, 1986) Methane production ∝ P0qT −T0

10

, P0 is calibrated (Stepanenko et al., 2011) Methane oxidation follows Michaelis-Menthen equation

Sinks and sources of gases in a lake

Model validation for Seida Lake

Guseva et al., Geogr. Env. Sust., 2016

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 16 / 54

Seida lake location Bubble flux (starting from 01.07.2007)

Kuivaj¨ arvi Lake (Finland)

Point of measurements

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 17 / 54

Mesotrophic, dimictic lake Area 0.62 km2 (length 2.6 km, modal fetch 410 m) Altitude 142 m a.s.l. Maximal depth 13.2 m, average depth 6.4 m, depth at the point of measurements 12.5 m Catchment area 9.4 km2 Secchi depth 1.2 – 1.5 m

Measurements

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 18 / 54

Measurement raft Footprint of the raft measurements

Conducted since 2009 by University of Helsinki Ultrasonic anemometer USA-1, Metek GmbH Enclosed-path infrared gas analyzers, LI-7200, LI-COR Inc. Four-way net radiometer (CNR-1) relative humidity at the height of 1.5 m (MP102H-530300, Rotronic AG) thermistor string of 16 Pt100 resistance thermometers (depths 0.2, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, 5.0, 6.0, 7.0, 8.0, 10.0 and 12.0 m) Turbulent fluxes were calculated from 10 Hz raw data by EddyUH software

Water temperature

Mixed layer depth and surface temperature (RMSE=1.54 ◦C) are well reproduced Stratification strength in the thermocline is overestimated Model results lack frequent temperature oscillations in the thermocline

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 19 / 54

Measurements Model

Oxygen

Stepanenko et al., Geosci. Mod. Dev., 2016 Seasonal pattern is well captured: oxygen is produced in the mixed layer and consumed below Oxygen concentration in the mixed layer is underestimated by 1-1.5 mg/l, and more significantly during autumn overturn

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 20 / 54

Measurements Model

Carbon dioxide concentration

Stepanenko et al., Geosci. Mod. Dev., 2016 Seasonal pattern is simulated realistically: carbon dioxide is consumed by photosynthesis in the mixed layer and produced in the thermocline and hypolimnion by aerobic organics decomposition Sudden CO2 increase prior to autumn overturn is absent in the model

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 21 / 54

Measurements Model

Methane

Stepanenko et al., Geosci. Mod. Dev., 2016 Methane starts to accumulate near bottom in the late summer when oxygen concentration drops to low values Surface methane concentration is very small leading to negligible diffusive flux to the atmosphere, consistent with measurements

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 22 / 54

Measurements Model

Further development: dissolved and particulate carbon

Adopting approach from Hanson et al., 2004

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 23 / 54

The Hanson et al. model is reformulated to explicitly reproduce vertical distribution of DOC, POCL, POCD (instead of using mixed-layer and hypolimnion pools, as in original paper) The horizontal influx from catchment is to be included

Extended biogeochemical model

∂CCH4 ∂t = DifA(CCH4 ) + BCH4 − OCH4 , (1) ∂CO2 ∂t = DifA(CO2 ) + BO2 + PO2 − RO2 − DO2 − SO2 − OO2 , (2) ∂CDIC ∂t = DifA(CDIC ) + BCO2 − PCO2 + RCO2 + DCO2 + SCO2 + OCO2 , (3) ∂ρDOC ∂t = Dif(ρDOC ) + EP OCL − DDOC , (4) ∂ρP OCL ∂t = Dif(ρP OCL) + PP OCL − RP OCL− EP OCL − Dh,P OCL , (5) ∂ρP OCD ∂t = Dif(ρP OCD) − wg h ∂ρP OCD ∂ξ − DP OCD + Dh,P OCL . (6)

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 24 / 54

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 25 / 54

Motivation for inclusion of rivers in ESMs

• river runoff affects thermohaline circulation
• river runoff is the most precisely measured component of the land water balance
• rivers are considered as an substantial player in land carbon cycle
• the level and ice regimes of rivers can become the one of the most in-demand output of

ESMs

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 26 / 54

All values are in Pg C yr−1 (Battin et al., 2009)

River runoff in INMCM model

54 major basins surface and subsurface runoff are integrated over basins and instantaneously "added" to oceans in salinity equation no river tile in the surface energy balance calculations

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 27 / 54

River routing for Earth System Models

Exemplified by (Yamazaki et al., 2009)

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 28 / 54

Stream upscaling

External parameters for river model: flow direction riverbed slope parameters of cross-section geometry

Riverflow dynamic equations

Saint-Venant system: ∂S ∂t + ∂SU ∂x = Er, ∂SU ∂t + ∂SU 2 ∂x = −g ∂(hb + hr) ∂x − gU 2 C2(R)R + ∂ ∂xνr ∂U ∂x , hr = f(S).

Highlighted are inertia terms that can be omitted if Fr . =

U2 g(∆hb+∆hr) ≪ 1

| ∂hr

∂x | ≪ | ∂hb ∂x | at Fr < 0.1 (Dingman, 1984)

Longitudinal viscosity effects are also considered small Using ∆hb = ∂hb

∂x ∆x = s∆x, Froude number criterium becomes

∆x > 10U 2 gs ∼ 100 m for plain rivers. Under these conditions comes Manning’s equation: U = 1 nR2/3s1/2 .

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 29 / 54

Outlook

• > introducing of POC and DOC dynamics in lakes to improve CO2

simulations

• > introducing nutrient dynamics in lakes
• > simulations of future climate with lakes embedded in INMCM ESM (with

both thermodynamic and biogeochemical coupling)

• > testing river dynamics module in INMCM land surface scheme

The work is supported by grants RSF 17-17-01210 and RFBR 17-05-01165

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 30 / 54

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 31 / 54

Soil columns in the model

Horizontal projection

Soil columns are geometric figures of the same vertical dimension and with horizontal sections confined by sequential isobaths: z = 0 z = h z = z1 z = z2 z = z3 z = z4 y x

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 32 / 54

Freshwaters in global carbon cycle

Total freshwater methane emission is 104 Tg yr−1, i.e. 50% of global wetland emission (177-284 Tg yr−1, IPCC, 2013) greenhouse warming potentials from freshwater-originating CO2 and CH4 are roughly equal

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 33 / 54

(Tranvik et al. 2009) (Bastviken et al. 2011)

Thermocline thickness

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 34 / 54

Thermocline thickness is defined as a depth difference between 8 ◦C and 14 ◦C isotherms

Methane budget in the surface mixed layer

Mixed layer

The diffusive flux through thermocline is negligible compared to other terms

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 35 / 54

TKE profiles

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 36 / 54

TKE TKE balance terms

Significance of Coriolis force for Kuivaj¨ arvi Lake

Rossby deformation number, Ro = NH

f

≈ √

gρ−1 ∆ρ hML f

The lake’s length

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 37 / 54

Rotational effects are comparable with those of stratification.

The effect of barotropic seiches on methane

Neglecting barotropic seiches leads to TKE ≈ 0 below thermocline, less

• xygen flux from above and earlier accumulation of methane near bottom

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 38 / 54

Control simulation Seiches excluded

Effect of weak turbulent mixing in the thermocline

Oxygen diffuses downwards, oxidizing methane

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 39 / 54

Control simulation Increased minimal diffusion coefficient (10 ∗ λw0)

Figure: Эволюция глубины перемешанного слоя в эксперименте К.-Ф., дополненном учетом силы Кориолиса и параметризацией бароклинных сейш (результаты моделирования), при горизонтальных размерах озера: 300 км×300 м и 300 км×300 км

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 40 / 54

Observations

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 41 / 54

Measurement raft Footprint of the raft measurements

Conducted since 2009 by University of Helsinki Ultrasonic anemometer USA-1, Metek GmbH Enclosed-path infrared gas analyzers, LI-7200, LI-COR Inc. Four-way net radiometer (CNR-1) relative humidity at the height of 1.5 m (MP102H-530300, Rotronic AG) thermistor string of 16 Pt100 resistance thermometers (depths 0.2, 0.5, 1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, 5.0, 6.0, 7.0, 8.0, 10.0 and 12.0 m) Turbulent fluxes were calculated from 10 Hz raw data by EddyUH software

Internal seiche mixing parameterization in k − ǫ model

Goudsmit et al. 2002

Shear production is generalized to include seiches P = νtM 2 + Ps ; TKE production by seiche-induced shear at lake’s margins Ps = −

1−Cdiss√ Cd,bot ρw0cAb

γ 1

A dA dz N 2E3/2 s

, Es - seiche energy; Seiche energy is derived from wind forcing:

dEs dt = αA0ρaCd(u2 + v2)3/2 − γE3/2 s

Stationary Richardson number (Burchard, 2002) may be derived for this case as Rist =

P r∆cǫ21 ∆cǫ23−ν−1 P rCs∆cǫ21(u2+v2)3/2 ≈ 0.30 for typical wind speed

Ri ≫ 1

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 42 / 54

Kato-Phillips experiment

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 43 / 54

no heat and radiation flux at the top and bottom boundaries constant surface wind stress 0.01 N/m2 linear stable initial temperature profile, 2 K/m no morphometry no rotation depth 7 m, 60 vertical computational layers 10 days of the model integration

Kato-Phillips experiment results: Standard k-ǫ model

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 44 / 54

Time-depth distribution of velocity The deepening of the mixed layer follows well the known formula HML = 1.05u∗wN−1/2 t1/2 (Price, 1979) Temperature and eddy conductivity profiles, 9-th day After complete mixing of temperature the flow is classical Cuette flow

Kato-Phillips experiment: k-ǫ model + barotropic seiches

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 45 / 54

Time-depth distribution of velocity Temperature and eddy conductivity profiles, 9-th day

The bottom return flow creates another mixing layer, having very thin extremely stratified interface with the upper one

Kato-Phillips experiment: k-ǫ model + baroclinic seiches

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 46 / 54

Time-depth distribution of velocity Temperature and eddy conductivity profiles, 9-th day

The response of velocity to wind stress is waves, with dominating 1-st vertical mode, ∼ 1 day period. The thermocline is preserved, with both surface and bottom mixed layer present

Mixed-layer depth

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 47 / 54

K-P – Kato-Phillips experiment, K-P+kor – Kato-Phillips experiment with Coriolis force, K-P+bts – Kato-Phillips experiment with barotropic seiches, K-P+bts – Kato-Phillips experiment with baroclinic seiches

Rotation and seiching impose similar suppressing effect on vertical mixing Barotropic seiche parameterization is not enough to produce "correct" mixing

Mixed-layer depth

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 48 / 54

Kato-Phillips experiment with Coriolis force and baroclinic seiches at different lake sizes: 300 m×300 m, LR × LR and 300 km×300 km (LR ≈2.77 km) Coriolis force playes significant role in mixing compared to seiching only for the lake size L ≫ LR The effect of Coriolis force for very large lakes is similar in magnitude to that of seiching for small lakes

Why the mixing is suppressed by rotation and seiching?

In classical Kato-Phillips setup, the friction is zero at the base of mixed layer, leading to continuous increase of total momentum in mixed layer (under constant momentum flux from the atmosphere), the shear production of TKE and mixed-layer deepening until complete mixing of temperature and achieving stationary Cuette flow (where the momentum flux at the top is compensated by friction at the bottom) In both cases of rotation and seiching quasi-stationary oscillatory velocity patterns are established where Coriolis and pressure gradient terms (respectively) "consume" the constant momentum flux from the atmosphere

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 49 / 54

Пограничные слои в водоемах (данные наблюдений)

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 50 / 54

верхний перемешанный слой – эпилимнион: ТКЕ генерируется в основном за счет сдвига скорости (∼ напряжение трения) средний слой – металимнион (термоклин): очень устойчиво стратифицированный нижий слой – гиполимнион: генерация ТКЕ за счет сдвига внутренних циркуляций (сейши, волн Кельвина)

W¨ uest and Lorke 2003, Annu.Rev.Fluid Mech.

Критерий справедливости одномерного приближения

Показано на примере оз.Куйваярви (Финляндия)

W = g∆ρh2

1

ρ0u2

∗L

Running means

Wcr = 1

2

Wedderburn number

Shintani et al., 2010

LN = 2(zm−zv)V ρ0gh1

zvτA0L

LN,cr = 1

Lake number

Imerito, 2015

Thermocline displacement is negligible compared to mixed-layer depth

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 51 / 54

Схемы маршрутизации водотоков (на примере TRIP)

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 52 / 54

Вычислительно простые схемы, достаточные для воспроизведения средних расходов Диагностические формулы для расхода рек -> не воспроизводят экстремальные явления Нет расчета термодинамики и льдообразования Не учитываются биогеохимические процессы

Модели водоема в климатических моделях и системах прогноза погоды

V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 53 / 54

Модель пузырька

Ci Mi, Pi vb rb

dissolution, exsolution V.M.Stepanenko (MSU) Advances in land surface hydrology 5 September 2017 54 / 54

For shallow lakes (several meters), bubbles reach water surface not affected, for deeper lakes bubble dissolution has to be taken into account. Five gases are considered in a bubble: CH4, CO2, O2, N2, Ar Bubbles are composed of CH4 and N2 when they are emitted from sediments The velocity of bubble, vb, is determined by balance between buoyancy and friction The molar quantity of i-th gas in a bubble, Mi, changes according to gas exchange equation (McGinnis et al., 2006): dMi dt = vb ∂Mi ∂z = −4πr2

bKi(Hi(T)Pi − Ci).

Gas exchange with solution is included in conservation equation for i-th gas : ∂Ci ∂t = 1 A ∂ ∂z Ak ∂Ci ∂z + 1 A ∂ABCi ∂z + F(z, t, Ci, A) + (HCi − BCi,b ) 1 A dA dz .

Methane ebullition from different soil columns