Complex bog landscape model (COMBOLA) as an integrated tool for modeling of biotic turnover and peat deposit processes

Biotic cycling in ecosystems consists of live organic matter production and dead organic matter destruction. The latter is accompanied by the emission of greenhouse gases into the atmosphere. In peatland landscapes, additional conditions are imposed due to the presence of a water table depth (WTD), under which the destruction is anaerobic with methane generation, while above the WTD it is aerobic, and part of the diffusing methane is consumed by the methanotrophic bacteria. Hence, due to the complexity of heat and water transfer processes in the peat deposit and the nonlinearity of biological turnover, it is necessary to make a combination of their models. The COmplex Model of BOg LAndscapes (COMBOLA) is a set of dynamic models of carbon and nitrogen turnover, net ecosystem exchange, water balance, heat and water transport, generation and transfer of CO2 and CH4 in a peat deposit on annual, seasonal, and daily time scales. The main component includes a series of biotic turnover models – from a mass-balance equation on an annual time scale to a NEE dynamics model on a daily one. Biotic turnover can be represented by a single carbon cycle, a single nitrogen cycle or both. Another important component of the COMBOLA system is a one-dimensional vertical model of heat-water-gas exchange in a peat deposit. Thus, a number of interconnected modules constitute an integrated mathematical model of peatland landscapes adapted to any given initial information.


Introduction: models of biological turnover, heat and mass transfer in peatland landscapes
One of extremely large carbon storages on the Earth is contained in peat deposits of wetland ecosystems: by various estimates, about 120-455 GtC [1,2]. The wetlands in Western Siberia occupy almost 50% of the area containing 36% of the whole soil carbon pool in Russia [3]. Carbon pool assessment of Western Siberian peatlands varies from 55 GtC [4] to 70 GtC [5]. Natural and anthropogenic transformations of peatlands are reflected in the functioning of their biogeochemical cycles, impact on the regional climatic processes, the greenhouse gas emission and the feedback between the climate and biosphere. Having a unique ability to conserve atmospheric carbon in peat, peatlands simultaneously are among the key emitters of methane whose effect is more than 20 times stronger than the greenhouse effect in the atmosphere compared to carbon dioxide. Thus, the reaction of peatlands to climate change and anthropogenic perturbations is determined by both the biological turnover including microbe destruction of dead organic matter and heterotrophic respiration and the processes of generation, transport, and absorption of methane in peat layers that depend on the water table depth regulating the aerobic/anaerobic layers' relation for a peat deposit. The resulting 2 greenhouse gas emission from the peatland surface becomes the product of complex biochemical interactions and depends on a landscape type. The main carbon reservoirs in peatland ecosystems are concentrated in the phytomass and dead organic matter in a peat deposit. The carbon balance for peatland ecosystems is formed by photosynthetic fixation of atmospheric CO 2 in the gross primary production (GPP), the autotrophic and heterotrophic respiration, the emission of carbon dioxide and methane from the peat deposit, and the organic matter output by run-off. Since a variety of ecological factors can influence the carbon exchange, the actual climate change can perform them with transforming both the peatland landscapes and their carbon flows and balance [6]. The contemporary climate change seems extremely hard in the northern hemisphere and, in particular, in the boreal zone of Northern Eurasia where west-Siberian peatlands are located. Under them the peatland ecosystems can become an essential source of carbon to the atmosphere, since warming stimulates not only the NPP increase, but also dead organic matter decay [7]. A number of investigations show complicated reactions of biotic turnover in peatlands to the air temperature change [1,8], as well as of carbon dioxide emission to the water table depth oscillations [7,9]. Also, the net primary productivity (NPP) is definitely sensitive to the hydrological and temperature conditions initiating changes in the productive branch of the carbon cycle under the climate change [10]. All of these reasons make necessary an integrated approach in mathematical modeling of the biological turnover and heat and mass transfer in a peat deposit.
Due to complexity of the structure and functioning and the lack of experimental data, mathematical models of the biological turnover in peatlands are essentially weaker developed than similar models for the forest dynamics. Most of them -Wetland-DNDC [11], MWM [12], CaMP [13] -follow the simulation modeling principles and operate on the daily or less time scale providing for the need in large amount of information on the species-specific (for vegetation) or highly dispersed parameters. In the MILLENNIA model [14], an age-cohort approach with a variable water table depth is used allowing validation with the peat age data and testing of the peat soil organic carbon stock response to climate change scenarios without biotic turnover dynamics. In this paper, we describe a model system COMBOLA (Complex Model of Bog Landscapes) for combined modeling of the biological turnover (C-and N-cycles both or separately) and the processes of water, heat, and carbon gases transfer in a peat deposit on various time scales.
2. General design scheme for the model system COMBOLA: layer-module structure In Figure 1, an elementary unit of peatland landscape with vertical structure is shown. The biotic cycling processes in ecosystems consist of the live organic matter producing branch and the dead organic matter destruction one. The latter is accompanied by the emission of greenhouse gases into the atmosphere. Since the destruction of the dead organic matter takes place over the whole depth profile, the overall emission is closely connected with the generation, transport, and consumption of gases accompanied by the heat and water transport. In peatlands, additional conditions are imposed due to the presence of water table depth (WTD) under which the destruction is anaerobic with methane generation, while above the WTD it is aerobic, and a part of the diffusing methane is consumed by methanotrophic bacteria. Hence, due to complexity of the heat and water transfer processes in a peat deposit as well as nonlinearity of the biotic turnover functioning, it is necessary to develop combined models of them on various time scales in order to assess possible evolution of peatland landscapes under the climate change and economic exploitation.
The COmplex Model of BOg LAndscapes (COMBOLA) is a set of dynamic models of carbon and nitrogen turnover, net ecosystem exchange, water balance, heat and water transport, generation and transfer of CO 2 and CH 4 in a peat deposit on annual, seasonal, and daily time scales. The first main component is a number of biotic turnover modelsfrom the mass-balance equation on the annual scale to the NEE dynamics model on the daily time scale. Another important component of the COMBOLA system is the one-dimensional model of heat, water, and gas exchange in a peat deposit. On the whole, the COMBOLA can be characterized by a layer-module structure that provides the selected problem to be solved using a set of programming modules from different layers depending on 3 the modeling goal, selected time scale (annual, seasonal, daily), and the presence of initial data for calibration. A layer is called a set of identical phenomena on all time and spatial scales. A module is a part of the layer distinguished for the description of concrete process on the selected time scale. Thus, a number of interconnected modules constitute an integrated mathematical model of the peatland landscape adopted for given initial information. For example, if the user has annual data on the carbon storages and fluxes in a peatland ecosystem together with information on the WTD and components of the water balance, there is an opportunity to make a dynamic model of carbon turnover with water cycle on the annual time scale. Provided other information on the biophysical processes in the peat deposit, the user can include the corresponding modules of heat and/or water and/or CO 2 transfer into it for more adequate calculation of the overall annual carbon emission.

Multiple time scale modeling of biogeochemical cycles in peatland landscapes
A number of ecosystem models and approaches for analyzing the dynamics of biological turnover form a basis for the COMBOLA system. On the annual time scale, the dynamic compartment approach is applied, where for the set of carbon and/or nitrogen pools exchanging by matter fluxes between each other and the environment, dynamic mass-balance equations are designed with the flow functions depending on the storages in the interacting pools in a simply linear (donor or recipient type) or nonlinear (Lotka-Volterra or rational types) way [15]. The coefficients of the flow functions depend on the climate induced parameters like the annual mean air temperature and the total precipitation. Since peatlands are ecosystems with a special role of water, the standard compartment approach is added by elements or the full water cycle with the coefficients also depending on the mean air temperature. The dynamic compartment equations can be calibrated by really measured data on the storages and fluxes in peatland landscapes that serve as equilibria for the compartment model, since they represent observable states stably existing for a rather long time. In Figure 2 one can see a set of 1 Water where C i are the carbon storages in pools (gC/m 2 ), w is the water content of the unsaturated peat layer (%), NPP = C 1 φ 1 (C 1 ) φ 2 (w) is the net primary production (gC/m 2 /year), P is the total annual precipitation, E 0 is the potential evaporation, d 2 c is the dead organic matter decay rate, d 2 p is the peat formation intensity, q i are the input flows, and α 12 is the litterfall rate. Calibration of the parameters in the water balance equation is carried out by WTD data and the peat deposit depth (Hilbert et al. 2000 [19]).
If sufficient data are available for a multi-component static turnover scheme at one or more time moments or for various steady states, similar dynamic models can be easily obtained (e.g., [20]). Bifurcation analysis of this type of dynamic equations can give parametric domains of stability for steady states reflecting the observed landscapes and allows one to calculate critical parameter values which initiate change of the landscape or peatland type [15].
Modification of this approach is possible if there is a combined carbon-nitrogen cycle in a biological turnover with corresponding data available. The interaction between carbon and nitrogen serves as a basis for the aggregated dynamic model of a combined carbon-nitrogen cycle in various types of peatland ecosystems. This interaction includes two main biochemical principles realized in bog models for the first time in [18,20]: 1) the litterfall intensity is proportional to the C/N ratio in the living phytomass that reflects nitrogen starvation of plants; the decay rate for dead organic matter decreases with increase of the C/N ratio of dead organic matter. Figure 3 shows three types of bog landscapes from Western Siberian middle taiga with a twocompartment structure of the storages and fluxes. The pools include living and dead organic matter. Dynamic models for a combined carbon-nitrogen turnover with similar principles were developed for forests in [21,22], where oscillatory processes initiated by different types of felling are considered. Including the nitrogen cycle module allows the researcher to study how the deficit of mineral nitrogen affects the evolution of bog landscapes under the climate change. On the monthly and daily time scales, another approach is used as a keystone for modeling the biological turnover. For inclusion into the COMBOLA system, a new dynamic model of Net Ecosystem Exchange (NEE) is developed and calibrated on data measured at the low forested bog and sedge-sphagnum bog landscape types. The net ecosystem exchange models are used to study the relative importance of different environmental factors, fill gaps in the time series to calculate the daily and annual carbon budget [25]. The measured NEE was partitioned into Gross Primary Production (GPP) and Ecosystem Respiration (ER): NEE = ER -GPP.
The total ecosystem respiration was modeled using the exponential equation widely used for explanation of ER variations [26,27]: where T a is the air temperature ( o C), E H and E R are the base respiration levels at 0 o C, k H and k A are the temperature sensitivity factors. The air temperature was used as an explanatory factor. We tested the soil temperature at a depth of 10 cm and the air temperature as explanatory variables in bare soil experiments. In each case the air temperature explained the observed variations in the CO 2 fluxes better than the soil temperature due to smoothed oscillations in the soil temperatures. Autotrophic respiration is absent at the beginning of the vegetation season, it started to increase proportionally to the LAI with rise of the green vegetation. It was shown in many experiments that the photosynthetic response to low light intensities has a linear character and photosynthetic saturation is observed at high light intensities. A rectangular hyperbolic function was used for the light response of the NEE in the daytime [25][26][27]: GPP = α×PAR×G m / (α×PAR + G m ), where α is the initial slope of the light response curve at low light (photosynthetic efficiency), and G m is the theoretical maximum rate of photosynthesis at infinite PAR (photosynthetic capacity). For its calibration, the proposed model requires measurements of the CO 2 emission, air temperature, photosynthetically active radiation at least in a part of the vegetation period of a chosen year. Then the model allows one to calculate the carbon dioxide fluxes for the particular peatland landscape. Figure 4 shows the calculated daily values of the main carbon balance components. The total ecosystem respiration increases from 0.64 gCm -2 day -1 at the beginning of April to 2.37 gCm -2 day -1 in the middle of June. The maximal value of the GPP in 2016 (5.02 gCm -2 day -1 ) is obtained for June 23 and the minimal one (0.39 gCm -2 day -1 ), at the beginning of April. Due to essential variations of the control parameters, the calculated NEE value is oscillating in the interval from 0.89 to -2.86 gCm -2 day -1 with an average value of -0,81 gCm -2 day -1 . The carbon assimilation in vegetation usually is greater than the total ecosystem respiration. The total annual accumulation of carbon in the sedge-sphagnum fen in 2016 constituted 148,6 gCm -2 and the net primary productivity (NPP = GPP -AR), 300 gCm -2 . The carbon balance calculations by monitoring of the meteorological parameters result for the studied sedge-sphagnum fen to be a stable sink of atmospheric carbon accumulating from 112 (2011) to 148 (2016) gCm -2 .

Modeling heat, carbon dioxide, and methane transfer processes in a peat deposit
One of the key components of the COMBOLA system is a one-dimensional model of heat, water, and gas transfer in a heterogeneous peat deposit consisting of a pre-assigned number of peat layers with particular bulk density, moisture, temperature profile, and WTD on the annual and monthly time scales through the vegetation period. The simplest case of this model is represented by a dynamic model of heat and carbon dioxide transfer in a two-layer peat deposit. The model is based on an approach from [28] with changes allowing one to combine it with the mass-balance model of biological turnover.
The one-dimensional dynamic equation of CO 2 diffusion in a two-layer peat deposit separated by WTD has the following form: In this system, C CO2 is the total concentration of carbon dioxide at depth z, T(z) is the peat temperature at depth z, D CO2 (w) is the combined diffusion coefficient, C a , T a are the air concentration of CO 2 and the surface air temperature, F d (T, w) is the factor of temperature and water content impact on the dead organic matter C s decay rate, λ is the heat conductivity coefficient, and c eff is the effective depends on the diffusion of CO 2 in the air and in the water, the peat porosity P or , the solubility of CO 2 in the water α CO2 (T) and the tortuosity τ. This problem with the initial and boundary conditions is solved by different numerical methods together with the standard initial value problem for differential equations of biological turnover. For the calibration and verification, data from measurements at some typical landscapes of the Bakchar bog in the Western Siberian south taiga is used.
There is rather large number of methane emission models for peatlands ( , where all mechanisms of methane vertical motions are accounted for: Q b is the bubble transport, P CH4 (z), O CH4 (z) are the production and oxidation at depth z, and B CH4 (z) is the consumption by the plant roots. Similar to CO 2 the onedimensional PDE problem with initial and boundary conditions is now tested on data obtained in 2016-2017 in stationary experimental areas in the middle and southern taiga of Western Siberia. Combining this model with the carbon dioxide dynamics and the biological turnover, we will make a feasible modeling tool for the researcher not only to estimate local emissions of greenhouse gases from the peatland surface, but also to make predictions on variability of the peatland landscapes under the impact of climate and human activities.

Conclusions
Although not all parts and mechanisms of the model system COMBOLA have so far been realized in a combined programming code, a general approach has been developed, and the essential components work together, allowing researchers to make specific separate and combined models of peatland landscapes. Depending on the available data, the researcher can use annual, monthly, or daily time scale models consisting of software modules describing the biological turnover and heat-gas transfer processes in a peat deposit on each of these time scales. On the annual time scale, a bifurcation parametric analysis is possible for estimating changes in the stability boundaries of the biological turnover, steady states under climate change and, in perspective, economic perturbations. The modular structure provides a sufficient degree of freedom for the researcher to design a set of models describing the complexity and data availability of particular peatland cases, as well as time scalability with a possibility of expert oriented determination of the dependencies in biological turnover. The first numerical experiments with COMBOLA components on data measured in southern taiga peatland landscapes of Western Siberia [32] have shown good agreement between a model of biological turnover (carbon cycle) combined with a CO 2 emission model and the observed values of the emission and other carbon cycle components.