Global cold and wet ecosystems, such as permafrost, northern peatlands, Arctic tundra, boreal wetlands, and alpine swamp meadows, store large amounts of soil organic carbon (SOC) and are typically water-rich. While it is well recognized that these ecosystems are highly vulnerable to climate warming as it accelerates SOC decomposition, how soil water levels regulate SOC decomposition and CO2 emissions specifically by constraining oxygen (O2) availability during the growing season remains poorly understood at large spatial scales. Here, we integrate field observations, global data sets, and process-based models to quantify how soil water dynamics influence CO2 emissions by regulating O2 diffusion across these ecosystems. Results from 107 field sites consistently reveal a protective effect of high soil water levels against SOC decomposition and CO2 release during the peak growing season (representing one-third of the year), suggesting that, beyond warming, the loss of this protection due to soil moisture reduction under near-saturated conditions is a critical driver of increased CO2 emissions. However, global data-driven data sets and process-based model simulations show divergent correlations between soil CO2 emissions and soil water levels. Many models, which rely on simplified soil moisture scalars or empirical functions to represent soil water level effects, inadequately reproduce the effects of soil water levels on SOC decomposition and CO2 emissions in cold and wet ecosystems. Our study underscores the urgent need to incorporate more mechanistic representations of soil water-mediated SOC decomposition into models and to improve spatially explicit hydraulic parameters, thereby enhancing projections of soil carbon-climate feedback.

This study demonstrates that the magnitude and spatial organization of pore-scale heterogeneity jointly govern the emergence of asymptotic hydrodynamic dispersion in porous media. Two-dimensional digitally constructed pore networks representing graded, laminated, and random sedimentary architectures were generated with systematically varied heterogeneity index σg. Flow and solute transport simulations were performed over longitudinal domains exceeding 20 cm across multiple Péclet numbers to quantify the evolution from pre-asymptotic to asymptotic transport behavior. Results show that spatial architectural organization—rather than the magnitude of heterogeneity alone—exerts a first-order control on emergent continuum-transport behavior. Bulk velocity variance remained comparatively insensitive to σg, yet the asymptotic hydrodynamic dispersion coefficient and accompanying longitudinal dispersivity varied by orders of magnitude following the hierarchy laminated > graded > random, reflecting how spatially coherent velocity segregation imposed by sedimentary architecture governs emergent continuum-transport behavior. In the graded and laminated architectures, both transport metrics increased exponentially with σg, while the laminated architecture additionally exhibited power-law relationships with the permeability contrast between coarse and fine sediment layers, defining architecture-specific constitutive relationships absent in random architectures. The transition length to asymptotic Fickian transport depended on both sedimentary architecture and σg, exhibiting only a weak dependence on the Péclet number in the laminated and graded domains. In contrast, the relative increase from early-time to asymptotic dispersion showed no systematic Péclet number dependence, indicating primary structural control. Once asymptotic conditions were reached, longitudinal dispersivity became independent of the Péclet number, consistent with classical hydrodynamic dispersion theory. Asymptotic dispersivity values reproduced both the magnitude and vertical scatter of compiled field dispersivity–length datasets, supporting the interpretation that much of the observed field-scale variability at a given length likely reflects the impact of geological architectural heterogeneity rather than measurement uncertainty.

Abstract Understanding baseflow dynamics is essential for sustainable ecosystem and water resource management. While most studies have documented climate‐baseflow relationships at regional scales, observational evidence of how groundwater influences baseflow dynamics across different spatiotemporal scales remains limited, substantially constraining mechanistic understanding. To address this gap, we conducted a systematic continental‐scale analysis of groundwater‐baseflow linkages across 18 major US watersheds from 1980 to 2020. We analyzed daily discharge from 7,016 stations and estimated baseflow using a recursive digital filter. We also identified 470 paired shallow groundwater well‐gauge sites within 5‐km proximity to evaluate hydraulic connectivity across multiple time scales. Key findings were as follows: (a) Annual baseflow trends revealed pronounced east‐west spatial clustering, with predominantly increasing trends in eastern watersheds and decreasing trends in western watersheds. (b) Seasonal analysis indicated strong synchronicity between groundwater level fluctuations and baseflow dynamics, with matched trend directions from 68.75% (fall) to 92.86% (summer) of paired sites. Spearman correlations further demonstrated robust positive seasonal relationships (ρ = 0.686–0.791, all p < 0.005). (c) At the monthly scale, linear regression and Spearman correlation analyses showed weak associations between baseflow and climate variables (R2 = 0.02–0.31), but substantially stronger associations with upstream groundwater levels (R2 = 0.61–0.95). (d) At the daily scale, strong associations between baseflow dynamics and groundwater levels (R2 = 0.63–0.95) further support direct hydraulic connectivity between groundwater and river systems. Collectively, these findings indicate that groundwater is the dominant control on baseflow temporal variability across multiple timescales, whereas climate influences may emerge at annual or longer time scales due to groundwater travel times and long system memory. Our results underscore the importance of groundwater monitoring for improving baseflow prediction and water resources management.

Understanding how cooling efficiency (CE) trends of tree covers respond to background climate change and local urbanization is essential for optimizing urban greening strategies. However, previous studies have largely focused on city-wide CE trends, often with a limited number of cities or with limited geographic scope. Consequently, the spatial heterogeneity in CE trends across urban-rural gradients, and their underlying drivers remain poorly understood, especially across cities worldwide. Here we quantified summer daytime CE trends from 2003 to 2020 across six urban-rural gradients (i.e., urban cores, new towns, urban fringes, suburbs, rural fringes, and rural backgrounds) in more than 5,000 global cities, utilizing MODIS-derived land surface temperature and tree cover data. We observed an inverse-V shape in CE trends along urban-rural gradient, with peak values in urban fringes (0.10 f 0.01 degrees C/%/century, mean f one standard error), followed by new towns, cores, and three rural gradients. The trends in CE were generally positive across most climate zones, yet arid regions exhibited a decline (-0.06 f 0.06 degrees C/%/century). CE trends strengthened with increasing city size in urban fringes, yet they decreased in cores and new towns. Using a LightGBM-SHAP algorithm, we found that the macro-scale background climate dominated the CE trends in urban areas (43%), whereas micro-scale local surface properties emerged as the primary contributors in rural areas (48%). Our findings provide critical insights into the spatial heterogeneity of CE trends of urban tree covers at a global scale.

This study investigated the Representative Elementary Volume (REV) for solute transport using pore-scale computational simulations and laboratory column experiments in homogeneous media up to 1.83 m long. Non-Fickian tailing was observed at transport distances shorter than the solute transport REV, becoming more pronounced with increasing Peclet number (Pe). Non-Fickian 'tails' transitioned to Fickian 'Gaussian' characteristics between 1 and 1.83 m. The apparent dispersion coefficient converged towards a steady hydrodynamic dispersion coefficient between 1 and 1.5 m, both signifying the emergence of pore-scale solute transport physics into a continuum and defining the solute transport REV, rather than the 'scale effect'. This challenges the assumption that Darcy's law, with its much smaller REV (2.8 cm), is adequate for defining the continuum for solute transport. When the transport length was shorter than the solute transport REV, there was an overestimation of the Pe-dependent dispersion coefficient due to solute stretching from non-Fickian tailing, leading to nonlinear relationships. As the transport length approached the solute transport REV, the nonlinear Pe-dependence converged towards a linear relationship, conforming with hydrodynamic dispersion theory in the mechanical transport regime. Similarly, the longitudinal dispersivity coefficient showed Pe dependence until the transport reached the solute transport REV, beyond which it became Pe-independent, confirming theoretical expectations once again. These findings suggest a need to reconsider previously proposed exponents that did not account for the REV specific to solute transport. The differences observed between computational simulations and laboratory experiments, such as the persistent memory of subtle tails in simulations, suggest areas for future research. However, this study sets a benchmark for the methodology of estimating dispersion and dispersivity coefficients, including their Pe dependencies, in accordance with theoretical expectations through REV analysis specific to solute transport.

Bedform-induced variations in pressure along streambeds can stimulate vertical hyporheic exchange in the hyporheic zone (HZ). The role of HZ in supporting riverine ecosystems depends on the spatial extent of HZ and the magnitude of mass and energy exchange rates. However, the quantification of the spatial extent of bedforminduced HZ remains less constrained, especially for the heterogeneous subsurface settings. This study therefore develops a streamline-based method to map the HZ. The identification process is fulfilled by analyzing the subsurface flow field based on simulating both of surface and subsurface fluid flow. Given that there is no flux exchange across the interface of HZ and underflow, the interface can be interpreted as the global minimum positive value of vertically-integrated flux for the entire domain. Like traditional particle tracking method, our method is tested to be effective to delineate the HZ for homogeneous cases, where the volume, mean residence time, and maximum depth of HZ can be directly estimated. Notably, our approach outperforms traditional methods for the heterogeneous porous and fractured streambeds, since the new method is viable to cover the entire stagnation zone within the HZ, which is critical for anaerobic reactions.

Land subsidence is widely present across the globe and brings catastrophic hazards. The well-acknowledged mechanism of subsidence is groundwater pumping, which leads to the reduction of hydraulic head and hence increases the effective stress, resulting in the vertical compaction of unconsolidated sediment. Here, we propose a hypothesis that subsidence in the coastal areas might be caused by osmotic effects, given the presence of seawater intrusion. The hypothesis is corroborated by simulating fluid flow, solute transport, and elastic deformation of multi-layered aquifer-aquitard systems. The simulations potentially cover a variety of natural environments by varying concentration, hydraulic head, thickness of aquitard, and hydraulic conductivity. We find that osmotic effects due to seawater intrusion play a non-negligible role in controlling subsidence in our studied cases, suggesting that future work on subsidence in areas impacted by seawater intrusion should fully incorporate osmotic effects to improve our understanding and prediction of subsidence.

Ecosystem water use efficiency (eWUE = gross primary production (GPP)/evapotranspiration (ET)) is widely used to characterize the coupling of ecosystem water and carbon processes. To investigate how eWUE responds to environmental changes, we compared environmental controls on annual, monthly, and daily GPP, ET, and eWUE from diverse ecosystems and climate regimes and quantified their daily relative impacts with machine learning techniques. Similar to GPP and ET, eWUE was strongly related to environmental variables at daily and monthly scales than at annual scales, indicating the tighter interplays of ecosystem processes with surroundings at shorter timescales. More critically, daily GPP and ET variations were primarily driven by net radiation (Rn) at most sites; whereas, vapor pressure deficit (VPD) dominated daily eWUE variations from humid to semi-arid sites, leaf area index (LAI) controlled eWUE variations at arid sites. It was largely attributed to the asynchronous responses of daily GPP and ET to environmental variables: the positive responses (though with different degrees) of daily GPP and ET to Rn weakened the Rn impact on eWUE; whereas, the opposite responses of daily GPP (negative) and ET (positive) to VPD enhanced the VPD impact on eWUE. The greater LAI impact on daily eWUE at arid sites was due to the dominant control of LAI on GPP variations under arid conditions. Unlike early eWUE models that incorporate VPD, our data showed that Rn could significantly improve eWUE models. This work provides valuable insights into understanding the controlling mechanisms of eWUE and ameliorating the representation of GPP and ET coupling.

Coastal wetlands are crucial ecosystems with high fragility to environmental alterations. Despite early efforts, interactions of water and salt with vegetation in coastal wetlands have rarely been investigated especially with increasing anthropogenic disturbances. Here, field data and numerical simulations were utilized to understand the influences of vegetation, evapotranspiration (ETa), tide, and inundation events on porewater and salt transport in highly regulated tidal riparian wetlands. At the site, the presence of Phragmites australis significantly enhanced soil hydraulic conductivity (KS) and facilitated hydraulic connections between groundwater and surface soil, resulting in the rapid response of soil salinity to groundwater level fluctuations. Simulation results further demonstrated the combined effects of tides and ETa on the salinity distributions, revealing the formation of a high salinity zone in the middle region due to the inflow of water from the tides and ETa-driven upward movements of water and salt; however, as salinity levels increased within the high salinity zone, the growth of Phragmites australis was constrained, leading to reduced transpiration (Ta). Tidal parameters such as amplitudes and mean water levels could affect salt accumulation processes in soils and further the high-salinity zone. Finally, inundation events from extreme upstream discharge had a flushing and diluting effect on accumulated salts, alleviating salinity stress and promoting vegetation Ta and growth; however, further increases in ETa led to salt re-accumulation, potentially reinstating the salinity distribution pattern after inundation. These findings enhance our understanding of ecohydrological processes in coastal wetlands and have important implications for wetland management and restoration.

Evaporative fraction (EF) is a useful measure for quantifying land surface energy partitioning processes and determining evaporative regimes; however, its influencing factors remain highly uncertain. Here, global data sets were compiled to disentangle the effects of environmental variables on EF variations along climate and land surface gradients. We found that (a) at annual timescales, ecosystem-level EF could be expressed as a power law function of aridity index. The relationships of mean annual soil water content (SWC) and leaf area index (LAI) with mean annual EF resembled the traditional evaporative regime theory; (b) at daily timescales, the boosted regression tree method quantitatively revealed that the impacts of environmental variables (including meteorological variables) on EF showed equal importance, especially at humid sites, primarily due to the different response direction and magnitude of latent heat (LE) and sensible heat (H) fluxes to environmental changes. Particularly, the contrasting responses of LE (positive) and H (negative) to SWC, LAI, and relative humidity enhanced the positive effects of those influencing variables on EF; whereas, the correlations between EF and energy-related factors (i.e., net radiation-Rn and air temperature-Ta) deteriorated as both LE and H showed positive response patterns to those variables; (c) meteorological factors were also found to have nonlinear effects on daily EF, further modified by climatic conditions. Rn near 150 W/m2 and Ta near 15 degrees C appeared to be important energy-partitioning thresholds at drier and humid sites, respectively. Moreover, changing interactions among environmental variables with climates were demonstrated to be important for better explaining EF variations. Evaporative fraction (EF; defined as latent heat flux-LE divided by the sum of LE and sensible heat flux-H) can vary under different climatic and land surface conditions, but its influencing factors remain poorly understood. To this end, we explored the effects of environmental variables on annual and daily EF variations at different sites around the globe. We found that mean annual EF decreased with increasing aridity index and increased with mean annual soil water content and leaf area index. By comparison, the boosted regression tree method quantitatively showed that environmental variables exerted equally important roles in regulating daily EF, especially at humid sites. The complex interplays of daily EF with environmental variables were mainly due to the different responses of LE and H to surrounding environments and the strong nonlinear and interactive effects of environmental variables. These results are important for understanding the driving mechanisms of EF and land surface energy partitioning processes along climate and land surface gradients. Environmental variables exerted equally important roles in regulating daily evaporative fraction, especially at humid sites Synchronous responses of latent and sensible heat to surroundings determined how environmental variables affected evaporative fraction Environmental factors had strong nonlinear and interactive effects on daily evaporative fraction, which was further modified by climates

Controlling groundwater table decline could mitigate land subsidence and induced environmental hazards in over-explored areas. Nevertheless, this becomes a challenge in the multi-layered porous system as (in)elastic deformation simultaneously occurs due to vast spatiotemporal variability in the groundwater table. In this study, SBAS-InSAR was used to estimate annual land deformation during 2017–2022 in a specific region of North China Plain, in which aquifers are composed of many layers of fine-grained compressible sediments and the groundwater table has experienced a prolonged decline. The random forest (RF) was applied to establish the nonlinear relationship between accumulated deformation and its potential driving factors, including the depth to the groundwater table (GWD) and its change rate, and the compressible sediment thickness. Results show that the marked subsidence and uplift co-exist in the region even though the groundwater table has risen widely since the South–North Water Diversion Project. The land subsidence is attributed to inelastic compaction of the thick compressible deposits in depression cone centers, where the GWD is over 40 m and 90 m in the shallow and deep aquifers, respectively. In contrast, the marked uplift is primarily attributed to fast rising of the groundwater table (e.g., −2.44 m/a). The RF predictions suggest that, to control the subsidence, the GWD should be less than 20 and 70 m in the shallow and deep aquifers, respectively, and the rising rate of the GWD should increase to 2–5 times of current rates in the depression cones. To mitigate the marked uplift, the rising rate of the GWD should reduce to 1/2–1/5 of the current rates in the shallow aquifers. The uneven deformations of sediments in the depression cone centers and uplift in their boundaries may exacerbate geohazards. Therefore, it is vital to implement appropriate governance of groundwater recovery in the multi-layered porous system.

Understanding and quantifying the solute tailing process in fractured rocks is critical for various subsurface applications but remains challenging. Particularly, the effects of flux exchange between fracture and matrix on the tailing process are considered to be important but have not been explicitly investigated. To fill this knowledge gap, we estimate the power-law scaling of late-time tracer tailing, n, based on three dimensional direct numerical simulations that consider both solute diffusion and advection across fracture-matrix interface. The dependence of solute decay rate n on fracture-matrix interaction is elucidated by simulating cases over a wide range of matrix permeability (km). In the case with matrix diffusion but without advection between fracture and matrix, we find that matrix diffusion significantly enhances tracer tailing. Then, the consideration of advection does not show significant effect on n up to the matrix permeability km of 10-12 m2. When k is greater than 10-12 m2, n declines rapidly down to 1.0 at k = 10-9 m2. We finally demonstrate that the upscaled model of continuous time random work can capture the late-time tailing behavior. This study provides a mechanistic explanation for the heavy tailing process commonly observed in fractured rocks.

Hydraulic heterogeneity leads to non‐Fickian transport characteristics, which cannot be entirely accounted for by the continuum‐scale advection‐dispersion equation. In this pore‐scale computational study, we investigate the combined effects of flow rate (i.e., Peclet number, Pe) and first‐order hydraulic heterogeneity, that is, resulting from intrapore geometry exclusively, on the transition from non‐Fickian to Fickian dispersion. A set of intrapore geometries is designed and quantified by a dimensionless pore geometry factor ( β ), which accounts for a broad range of pore shapes likely found in nature. Navier‐Stokes and Advection‐Diffusion equations are solved numerically to study the transport phenomenon using velocity variance, residence time distribution, and coefficients of hydrodynamic dispersion and dispersivity. We determine the length scale (i.e., the linear distance in flow direction) for each pore shape and Pe when non‐Fickian features transition to the Fickian transport regime by incrementally extending the length, that is, the linear array of pores. We show how velocity distribution and variance ( σ 2 ) depend on β , and directly control the transition to Fickian dispersion. Pores with a larger β , that is, complex pore shapes with constricted pore‐body or with “slit‐type” attributes, result in a substantial non‐Fickian characteristics. The magnitude of non‐Fickian characteristics gets amplified with an increase in Pe requiring a significantly longer length scale, that is, up to 1 m or a linear array of 500 pores to transition to the Fickian transport regime. We find the hydrodynamic dispersion coefficient ( D h ) exponentially depends on the pore shape factor β , with its exponent dependent on flow rate or Pe. We determine constitutive relations to quantify how σ 2 , β , and Pe, contribute to the degree of non‐Fickian characteristics, the length scale needed for the transition to Fickian transport regime, asymptotic D h , and the length‐scale dependence of longitudinal dispersivity.

In arid area, the liquid water and water vapor states in soil profiles and fluxes at the upper and bottom interfaces are extremely complex due to heterogeneity of soil textures and the driving forces of heat and matrix potential. In this study, we used Hydrus-1D to simultaneously simulate liquid water, water vapor, and heat transports based on the observed data of atmosphere, soil and groundwater at three soil profiles in an arid area of northwest China. Comparison and contrast of the observed and simulated results at the three soil profiles show that there are diurnal vapor entry and outlet fluxes at the dry surface layer (DSL) of 30 cm in the summer season. The vapor entry and re-evaporation account for about 14% of annual precipitation for the heterogeneity soil profile with a mean groundwater depth of 210 cm. Because of limited soil moisture in this arid area, vapor induced re-evaporation occurs shortly in the early daytime. Moreover, the extent of vapor entry, condensation and re-evaporation strong depends on soil properties and water table depths. The lower water table produces the drier soil surface, allowing more vapor entry, condensation and re-evaporation. Whereas the finer grained soil layers benefits the vapor fixation to produce zero fluxes that substantially inhibit the upward liquid water and vapor fluxes, and thereby reduces soil evaporation. The reduced soil evaporation correspondingly decreases the capillary effect on phreatic evaporation, proven by that soil evaporation decreases slowly with decline of water table and the large extinct depth of phreatic evaporation for the finer grained soil profiles. The estimated extinct depth is 180 cm and 200 cm for the soil profiles consisting of silt loam and loamy sand, respectively, much larger than 100 cm of the sandy soil profile. Additionally, as water table is higher and lower than the extinct depth, the models neglecting the vapor - heat function could respectively overestimate and underestimate soil evaporation.

Estimation of actual evapotranspiration (ET a ) is challenging due to its complex interactions with surrounding environments; thus, understanding the complexity and predictability of ET a along with its influencing factors is essential to improve ET a estimation. Based on the FLUXNET and ChinaFLUX datasets, we first examined whether daily ET a exhibited chaotic behaviors, and then investigated how daily ET a complexity and predictability varied across climate regimes and ecosystems using the chaos theory. The results of recurrence plot and correlation dimension (CD) analysis suggested the existence of chaotic behaviors in daily ET a , implying the system controlling ET a processes was deterministic with a limited number of controlling variables. The recurrence quantification analysis further revealed a varying degree of complexity for ET a processes (e.g., determinism‐DET with a range from 0.01% to 40.7%) across the sites. Specifically, daily ET a featured stronger deterministic properties (e.g., higher DET and lower CD values) and thus higher predictability in humid and arid regions; whereas, higher degrees of daily ET a complexity emerged in sub‐humid and semi‐arid regions due to more complex interactions with environmental factors. Moreover, primary environmental controls on daily ET a complexity varied with climatic and land surface conditions. Results showed that energy and water supplies along with average daily temperature and relative humidity were the primary controls on the ET a complexity across the sites, while leaf area index also played an important role at drier sites. As a first attempt, this study provides additional avenues to understand the complexity of ET a processes from a global perspective using the chaos theory.

Solute transport process in fractured rocks is central to many engineering applications. However, how matrix diffusion affects non-Fickian transport behavior, which is typically characterized by the late time tailing in breakthrough curves and residence time distributions, remains underexplored, especially considering the coupling effects of the finite matrix thickness (h) and transport regimes quantified by Peclet number (Pe). In this work, the fluid flow and solute transport were simulated by directly solving the Navier-Stokes equations and advection-diffusion equation in a two-dimensional fracture-matrix system considering matrix diffusion. In total, eighty cases considering variations of h and Pe were conducted to investigate the impacts of h and Pe on the late time solute tailing characterized by the exponent (n) in the power law function. Numerical results demonstrate that n decreases with an increase of h and Pe, where n eventually reduces to 1.5 that is a theoretical value for an infinite matrix domain. Moreover, the relationship between n and (h, Pe) is established via a polynomial function, which is further validated by the three-dimensional heterogeneous fracture-matrix system. The newly-established function can be used to accurately assess the power law tailing process given the knowledge of measurable h and Pe.

EDITORIAL article Front. Environ. Sci., 26 May 2023Sec. Freshwater Science Volume 11 - 2023 | https://doi.org/10.3389/fenvs.2023.1221837