From the source to the reservoir and beyond—tracking sediment particles with modeling tools under climate change predictions (Carpathian Mts.)

The study tracks spatial and temporal distribution of sediment particles from their source to the deposition area in a dammed reservoir. This is particularly important due to the predicted future climate changes, which will increase the severity of problems with sediment transport, especially in catchments prone to erosion. Analyses were performed with a monthly step for two mineral and one mineral/organic sediment fractions delivered from the Carpathian Mts. catchment (Raba River) to the drinking water reservoir (Dobczyce) by combining SWAT (Soil and Water Assessment Tool), and AdH/PTM (Adaptive Hydraulics Model/Particle Tracking Model) modules on the digital platform—Macromodel DNS (Discharge Nutrient Sea). To take into account future changes in this catchment, a variant scenario analysis including RCP (representative concentration pathways) 4.5 and 8.5, and land use change forecasts, was performed. The differences between the two analyzed hydrological units (catchment and reservoir) have been highlighted and showed a large variability of the sediment load between months. The predicted climate changes will cause a significant increase of mineral fraction loads (silt and clay) during months with high flows. Due to the location and natural arrangement of the reservoir, silt particles will mainly affect faster loss of the first two reservoir zones capacities. The increased mobility of finer particles (clay) in the reservoir may be more problematic in the future, mainly due to their binding pollutant properties, and the possible negative impact on drinking water abstraction from the last reservoir zone. Moreover, the study shows that the monthly approach to forecasting the impact of climate change on sediment loads in the reservoir is recommended, instead of a seasonal one.


Introduction
Natural processes, such as erosion, sediment transport, and its deposition, can shape the land surface and influence the structure and function of river catchment ecosystems worldwide (Guillén Ludeña et al. 2017;Vercruysse et al. 2017;Wu et al. 2018;Stähly et al. 2020;Hu et al. 2021;Feng and Shen 2021). However, these processes are very sensitive to both climate change and land use, especially in mountainous areas (Hohmann et al. 2018;Dibike et al. 2018;Zhao et al. 2018;Aksoy et al. 2019;Zhang et al. 2020). Already, many mountain catchments are struggling with an increase in soil loss due to intensification of rainfall, and consequently, surface runoff (Halecki et al. 2018a, b;Borrelli et al. 2018;Mostowik et al. 2019;Berteni and Grossi 2020;Ciampalini et al. 2020). Moreover, the forecasted changes in precipitation and temperature suggest that the intensity of this phenomenon may increase significantly (Borrelli et al. 2020;Gianinetto et al. 2020;Szalińska et al. 2020). It also turns out that in such areas, even favorable changes in land use (LU) associated with a significant increase in the area of forests, which naturally stabilize the soil, may turn out to be insufficient to stop the negative effects of the climate changes (Orlińska-Woźniak et al. 2020a).
The consequences will be particularly severe in catchments with dammed reservoirs, trapping most of the sediment particles. Sediment entrapment has long been a major factor in reducing the capacity and the deterioration of water quality in many reservoirs around the world, endangering their durability, human health, and safety (Zarfl and Lucía 2018;Bilali et al. 2020;Huang et al. 2021). Multiple studies have shown today, as a result of the increase in sediment loads transported by rivers into reservoirs, that each year about 1% of their total capacity in the world is lost (Jain 2005;Rahmani et al. 2018). In order to effectively stop, or at least mitigate this negative trend, it is necessary to precisely predict the effects of climate and land use changes in river catchment systems.
Environmental models can play a key role in this, and research with their use has been conducted for many years, allowing, i.e., to estimate sediment loads introduced into a reservoir, the effectiveness of the reservoir as a sediment trap, loss of reservoir volume, and finally to identify areas of sediment deposition (Croley et al. 1978;Banasik et al. 1993;Garg and Jothiprakash 2010;Wisser et al. 2013;Charafi 2019;Bladé Castellet 2019). Currently, a wide range of sediment models are available that differ in complexity, accuracy, and spatial and temporal scales . Nevertheless, there are still significant limitations in these models, especially when a dam reservoir is localized within a catchment. Then, important elements like sorting of grain fractions, role of a reservoir backwater, and movement of sediments downstream from a reservoir are not taken into account. Although a river and a reservoir together form an integral dynamic system with deep interactions, at the same time, they are two separate entities in terms of hydrology, with a complex and stochastic nature of the processes of transport and deposition of sediment particles (Idrees et al. 2021). Until now, both of these elements of the catchment area have usually been studied using tools operating at different spatial and temporal scales, often treating them as separate units or largely generalizing the relationship between them. The main novelty of our work is the combination of two macroscale models on a digital platform-Macromodel DNS, one for the watershed-toriver-basin sediment transport (SWAT) (Arnold et al. 2012), and the other for sediment transport and deposition into the reservoir (AdH/PTM) (Thornton et al. 1996;Green et al. 2015). In this way, we created a unique tool to precisely track sediment particles of individual fractions from their source, through transport, and ultimately to settlement in a designated reservoir zone (Kimmel et al. 1990; Thornton et al. 1996;Green et al. 2015). Moreover, the conducted analyses included four detailed variant scenarios which took into account both monthly changes of climate (precipitation and temperature) in short-and long-time horizons, as well as changes in land use (forest and urban areas). This enabled the precise indication of critical periods for the sediment load increase in the reservoir, which is of great interest, since storage capacity of dammed reservoirs is being lost worldwide due to sedimentation (Kondolf et al. 2014). The implemented approach can thus be a valuable tool for the assessment of reservoir sedimentation problems with regard to climate change effects and implementation of appropriate countermeasures. The developed tool was used for the mountain catchment of the upper Raba River, located in the Polish part of the Carpathian Mountains, and flowing into the Dobczyce Reservoir. This area is distinguished by intensity of meteorological phenomena and water erosion (Mikuś et al. 2019;Orlińska-Woźniak et al. 2020a;Szalińska et al. 2020). Moreover, land use forecasts indicate a gradual increase in forest cover in this area Price et al. 2017). The reservoir itself is the main source of drinking water for one of the largest agglomerations in Poland, and also serves as flood protection. As the Raba estuary to the reservoir is an effective trap of large particles, like sand and gravel, only suspended sediment fractions were analyzed in this study. Ultimately, two mineral SILT and CLAY, and one organic-mineral SMAG fractions have been tracked, taking into account the different nature of the catchment area and reservoir, and their different sensitivities to future changes. As a consequence, we were able to show the effects of the implementation in the variant scenario model, both for the catchment area, and subsequent zones of the dammed reservoir in a monthly time step.

Study area
The study area consists of two parts, the upper Raba River catchment (RR) and the Dobczyce Reservoir (DR) (Fig. 1). The upper Raba River flows 60 km from its source, located in the Carpathian Gorce Mts. (780 m a.s.l.), to the Dobczyce dammed reservoir (265 m a.s.l.). The RR catchment covers an area of 768 km 2 (Mazurkiewicz-Boroń 2016; Operacz 2017; Mikuś et al. 2019) and has a typically mountainous character with almost 43% of the area covered by slopes over 25% and overgrown by forest (Fig. 2). The mountainous nature of this catchment also manifests it with a fast reaction to precipitation. The river features very large amplitudes of water levels and flows, with maximum values recorded in the snow melt period, and early-summer long-lasting rainfall events. The studied catchment is covered generally by clays, and its lower slopes are used mostly for agricultural activities (over 40%). The entire area is well known as being particularly vulnerable to water erosion (Partyka 2002;Kijowska-Strugała et al. 2017;Halecki et al. 2018a, b).
The Dobczyce Reservoir was created as a result of dam construction undertaken in 1986 to provide drinking water for approx. half a million people, and to regulate the river flows (Mazurkiewicz-Boron 2016; Hachaj 2019). It is a multipurpose dam reservoir supplying people with drinking water, protecting against floods and droughts, Fig. 1 The research area with division into the upper Raba River catchment (RR) and the Dobczyce Reservoir (DR) Fig. 2 The upper Raba River catchment (RR) with soil, land use type, and terrain slopes producing electricity, and being a place for fish farming. The created reservoir has an area of approx. 10.7 km 2 (length of 8 km and width of 1.6 km), and an average depth of 12 m (maximum of 35 m) Zemełka et al. 2019;Wilk-Woźniak et al. 2021). The reservoir is divided into four zones (Fig. 3) including a riverine part (A), backwater (B), Myślenice basin (C), and Dobczyce basin (D) closed by the dam cross section (OUT). The entire area of the upper Raba River catchment belongs to the Carpathian climatic region of Poland (Kędra and Szczepanek 2019;Hachaj and Kołodziejczyk 2020) where the topography has a great influence on climatic conditions (Wypych et al. , 2018. Average annual air temperature is around 7 °C; however, the annual temperature amplitude reaches 21 °C. The average rainfall ranges from about 700 mm to about 1600 mm at the higher altitudes (Gorczyca et al. 2018). The growing season lasts about 200-210 days (April to mid-October).

Modeling tool
To achieve the goal of the study, i.e., to track sediment particles from their source to the deposition area, the DNS (Discharge Nutrient Sea) digital platform (Macromodel DNS) has been used as a modeling tool. The Macromodel DNS, developed at the Institute of Meteorology and Water Management-National Research Institute (Instytut Meteorologii i Gospodarki Wodnej-Państwowy Instytut Badawczy, IMGW-PIB) (Wilk et al. 2018b;Wilk and Orlińska-Woźniak 2019), provides an interactive platform allowing for integration of the SWAT (Soil and Water Assessment Tool) module (version 2012) (Arnold et al. 2012;Abbaspour et al. 2015) with other modeling tools (modules) to track different processes of the sediment/contaminant transport in a catchment (Wilk et al. 2018a;Szalińska et al. 2020Szalińska et al. , 2021Orlińska-Woźniak et al. 2020a). The RR catchment module has been created in the SWAT module ( Fig. 4) with the use of the following data: Subsequently, the created RR module has been used to simulate sediment yields from the catchment (land phase) and sediment loads (river phase) delivered to the reservoir (DR). Moreover, the SWAT module has been also used to estimate the sediment loads discharged into the DR in future time horizons under climate and land use changes. To simulate sediment transport and deposition processes into the DR, two additional modules have been incorporated into the DNS digital platform. Hydraulic conditions in the DR for sediment load scenarios have been modeled with use of the AdH (Adaptive Hydraulics Model) module, while the particle transport within the DR has been simulated with the use of the PTM (Particle Tracking Model) module ( Fig. 4)

SWAT module
The RR module applied in the current catchment study has been adopted from the model built for the entire Raba River catchment area, calibrated and validated for flow and sediment for the upper (upstream from DR) and lower (downstream from DR) parts. Calibration, verification, and validation of the model were performed using the SWAT-CUP (Abbaspour 2013), and the SUFI-2 algorithm (Khalid et al. 2016). The Latin Hypercube (LH-OAT) was used to identify the most influential model parameters (Jaiswal et al. 2017). Model calibration and validation were performed for the calculation profiles closing both parts of the catchment (  (Pool et al. 2018), whose respective value ranges are shown in Table S1. For the flow calibration statistical measures, R 2 , NSE, and KGE classified the model performance as good and very good, and PBIAS obtained satisfactory values. Sediment calibration for the upper Raba River catchment indicated a satisfactory model performance according to R 2 and NSE, and very good and good according to PBIAS and KGE. Better results were obtained for the lower part of the catchment, where the efficiency of the model can be considered very good (R 2 ) and good (NSE, PBIAS, and KGE). The validation of the flow model (according to R 2 , PBIAS, and KGE) showed good performance of the model, while only NSE had obtained a satisfactory value. For sediments, according to R 2 , NSE, and KGE, the efficiency of the model can be considered satisfactory and good according to PBIAS (Table A2).
The upper part of the catchment ( Fig. 1) has been separated using the GIS tool and the catchment boundaries. The calculation profile closing this area is Myślenice located directly above the DR. The land phase of the model (sediment yields) has been simulated for all RR 17 sub-catchments, delineated as hydrological response units (HRUs) combining unique land use, soil type, and slope features. The sediment yield estimations (SYLD) have been based on the modified universal soil loss equation (MUSLE), embedded in the SWAT module, according to Eq. (1) (Williams 1975;Vigiak et al. 2015): where sed is the sediment yield on a given day (metric tons), Q surf is the surface runoff volume (mmH 2 O ha −1 ), q peak is the peak runoff rate (m 3 s −1 ), area HRU is the area of the HRU (ha), K USLE is the USLE soil erodibility factor, C USLE is the USLE cover and management factor, P USLE is the USLE support practice factor, LS USLE is the USLE topographic factor, and CFRG is the coarse fragment factor.
The SWAT model offers four methods of tracking sediments in the riverbed phase (Bagnold, Kodoatie, Molinas, and Wu, Yang) of which Kodoatie is considered to be particularly effective for suspended and small sediment particles (Simons et al. 2004;Neitsch et al. 2011;Yen et al. 2017). Since a negligible amount of large particles, such as sand and gravel, reaches the DR, as the RR mouth is an effective trap for them, the current study focuses only on particles classified as suspended sediments (up to 0.062 mm). Therefore, the Kodoatie method (Kodoatie 2000) was used to simulate the channel phase of sediment transport, according to Eq. (2): where conc sed,ch,mx is the maximum sediment concentration (t m 3 ); v b ch is mean flow velocity (m s −1 ); y is mean flow depth (m); S is energy slope (m m −1 ); Q in is water entering the reach (m 3 ); a, b, c, and d are regression coefficients for different bed materials; W is channel width at the water level (m), and W btm is bottom width of the channel (m).
The climate change predictions were based on data from Euro-CORDEX, RCM models (Rummukainen 2016;Dosio 2016), and Global Climate Models-GCM (Yang et al. 2019). They included RCP4.5 and RCP8.5 emission scenarios (stabilization and raising pathways) for the short-term (H1, average for 2026-2035), and long-term (H2, average for 2046-2055) time horizons. The scenarios employed to estimate sediment load discharged into the DR were adopted from the Polish Development of Urban Adaptation Plans (UAP) (MPA 2021a, b), which have been described in detail in previous studies (Orlińska-Woźniak et al. 2020b;Szalińska et al. 2021). Since they have been prepared with a monthly time step, the impact of meteorological parameter changes could be tracked in detail, which was not possible when seasonal forecasts were given previously (Orlińska-Woźniak et al. 2020a). Briefly, besides the baseline scenario, four variant scenarios (VS1-VS4) have been created covering predicted temperature and precipitation changes under the short-term (H1) and long-term (H2) perspectives (Fig. 5a)  Along with the climate change predictions, the LU future changes have also been taken into consideration, and were included in each VS since they occur concurrently. The estimates projected on LU change impacts on sediment loads in the RR were based on results of the FORECOM project (http:// www. gis. geo. uj. edu. pl/ FOREC OM/ index. html), and transposed into the analyzed area with use of the DYNA-Clue model for the future time horizon of 2060 Price et al. 2017). LU scenarios included future changes in the forest cover and growth of urban areas. Depending on the scenario, one of the following forecasts was assigned to each RR sub-catchment: trend (growth of forest and urban areas by 23% and 10%, respectively), or liberal (growth of forest and urban areas, respectively by 30% and 15%) (Fig. 5b, c) (Orlińska-Woźniak et al. 2020a).

AdH/PTM module
To execute the sediment transport and deposition simulations in the DR, the AdH/PTM model Hachaj 2018Hachaj , 2019 was adapted in the current study. AdH is a 2D depth-averaged, finite element, modeling tool for simulating hydrodynamic conditions in river systems, but also in bay, reservoir, and lacustrine systems (Tate et al. 2014;Herrera-Díaz et al. 2017;Liu et al. 2021).
It uses the finite element method to solve two-dimensional momentum conservation equations for water in the Eulerian frame for both x and y directions (Berger 2010) according to Eq. (3):  works (2001-2012) commissioned by the Waterworks of the City of Kraków and reported by Nachlik et al. (2009) and Bojarski et al. (2012). This comparison resulted in a fine agreement between model simulations and reservoir observations. Based on the planar velocity field calculated by the AdH module, the particle's behavior over time (entrainment, advection, diffusion, settling, deposition, burying, etc.), the DR was simulated with the PTM module. Missing in the AdH module a vertical velocity component was calculated in the PTM module, based on the bed and water surfaces elevation profiles and the 2D velocity field provided by AdH. At each time step, PTM performs calculations to determine local characteristics of the environment (Euler frame, meshbased), and the behavior of each tracked representative particle (Lagrange frame, particle-based). The potential transport rate of sediment particles was calculated by the Soulsby-van Rijn method based on Eq. (5), while the other equations used by PTM, including the calculation of the particle trajectory is presented in (Neil et al. 2006).
where q t -total transport rate; A s -coefficient dependent on grain size; U -depth-averaged planar velocity; U w -average wave orbital velocity; C D -wave drag coefficient; and U cr -critical (threshold) velocity for motion/suspension regimes U.
The calculation profile of Myślenice, previously used as the closure of the RR part in the Eulerian SWAT module, was selected as the source of the sediment for the Lagrangian PTM module simulation in the DR part. Since the Eulerian approach is based on tracking a mass load while the Lagrangian one uses a concept of tracking individual particles, thus, at the connection point, a given mass load has to be transformed into a given number of representative particles. This issue has been solved by generating a number of representative particles for each suspended sediment fraction to be proportional to the appropriate mass load calculated by the SWAT module for a given month. To ensure that all the model-generated representative particles will start their movement in the water, the source area was set as an ellipsoid with the center point at 1.9 m above the bed and with the vertical and horizontal radii of 1 and 25 m, respectively. These dimensions were determined by the shape of the riverbed in the place where the sediment where v x , v y -velocity components in x and y directions, In the current approach, the boundary conditions of the DR AdH model, i.e., inflow rate and water surface elevation level (WSE), were adapted to include monthly average high flow (AHQ) and normal water surface level 269.9 m a.s.l. The simulation was performed with the default bottom roughness parameter and water density, as well as the neglected impact of wind, wave, and Coriolis forces. The AHQ values (Wrzesiński and Sobkowiak 2020;Luong et al. 2021) were calculated using flow data simulated by the SWAT module for each VS according to Eq. (4): where AHQ-average high flow for a given month, Y-number of years within the data series, and max (Qd)-the highest daily flow observed within a calendar month ( m index of AHQ and the max operator) of a given year ( y index of the max operator).
To verify the AdH model performance, the sensitivity analyses were first carried out for main variables accountable for the model response to channel roughness and wind forcing. As for the generalized Manning coefficient, the differences among the modeled flow velocities along the main current line were negligible for the coefficient range of 0.0025-0.05, and only noticeable when the parameter value was set at an unreasonably high level of 0.25 (Fig.  S1). Therefore, the AdH model default value of 0.025 was adopted in the current study. As for the model response to the wind forcing, testing was performed to verify the wind conditions impact on the resulting velocity field. The impact of the wind appeared to be important for low inflows (at the range of the modal inflow, approx. 2 m 3 s − 1 ), but negligible for values higher than the mean value for the DR (approx. 10 m 3 s − 1 ) (Hachaj 2019). The AdH model simulation results for the DR were validated through comparison with previously observed reservoir features, like appearance of current patterns and stagnant areas, as well as typical transport times/velocities in the DR selected basins. The hydrodynamics of the DR was investigated during the field Y source is located and allowed for the maximum filling of the selected profile with the shape of an ellipsoid. Simulations were performed using sediment density of 2650 kg m −3 and 1240 kg m −3 for mineral particles and mineral/organic aggregates respectively (Czuba et al. 2015). As for the particle diameter the following values have been adopted: CLAY d = 0.002 mm (dev = 1.25), SILT d = 0.05 mm (dev = 1.2), SMAG d = 0.03 mm (no deviation).
The PTM model validation was performed through comparison between results of sediment deposition simulation and grain distribution analyses for the DR bottom sediments and localization of progressive shallowing areas observed in satellite images for over 20 years (1997-2017) (Bojarski et al. 2012;Landsat8). The granulometry data contained samples (0-5 cm) collected in the eight reservoir cross sections (Szarek-Gwiazda and Sadowska 2010), as well the authors own granulometry analyses performed in the FEPE CUT laboratory (Szlapa 2019) on samples collected from the additional 5 cross sections in the DR backwater area (0-10 cm). The comparison results showed a good reflection of sedimentation processes (localization of progressive shallows, and grain segregation along the reservoir in various hydrodynamic conditions). Comparison of granulometric analysis results with the PTM simulations for SILT, included in the Supplementary Information (SI), (Fig. S2), revealed PBIAS value of − 21%, which signifies very good model performance (Table S2).

Sediment loads
The average sediment monthly loads (tons per month, t m −1 ) for the Myślenice calculation profile produced with use of the SWAT module have been presented in Fig. 6. For the baseline scenario, high load variability in all three analyzed sediment fractions (SMAG, SILT, and CLAY) is clearly visible with the coefficient of variation (CV) ranging from 85 to 191%. Following monthly distribution, two periods can be distinguished, one from October until April, and the second from May to September. During the first, average loads reached about 6 t m −1 for the SILT fraction, and 220 t m −1 for CLAY, while the loads for SMAG were negligible for this period. During the second period, a distinct increase of the loads could be noticed, up to 2.6 t m −1 for SMAG, 223 t m −1 for silt, and 785 t m −1 for CLAY on average. Generally, the CLAY fraction dominates the sediment loads introduced into the DR, reaching over 66% of the total load in June, when the maximum values during the year are observed.
Under the variant scenarios (VS1-VS4), the temporal pattern of sediment loads has been generally maintained, with elevated values during the May-October period, and their decrease during the remaining months of the year (Fig. 6). However, differences in the scenarios' impact on particular sediment fractions should be noticed. For the SMAG fraction, a decrease of the loads under all discussed scenarios is visible when compared to the baseline scenario. For the remaining fractions, deviations from this pattern are apparent, especially for CLAY. An increase of the clay fraction loads is evident during the winter and spring months (except for May), reaching extreme values in April, elevated by 197-265% when compared to the baseline scenario. During the remaining months, an increase in September is also noticeable, by 49.5-131.6% of the baseline load. A similar pattern of the load changes under the variant scenarios has been also detected for SILT, with the April loads for this fraction growing from 5.3 t m −1 , up to over 110 t m −1 under VS3 and VS4.

The average high discharges
The monthly AHQ values for the Myślenice calculation profile, subsequently used for the AdH/PTM module simulations, were obtained from the SWAT for the baseline and variant scenarios (Fig. 7). In the baseline scenario, the increase of the flow values can be observed from May to September, reaching up 147 m 3 s −1 in June. During the remaining months of the year (October-April), the AHQ flows vary from 16 m 3 s −1 for December to 58 m 3 s −1 in October. Under the variant scenarios, extension of the period with elevated AHQ values can be observed due to the increase of the flows from October to April, even by 45 m 3 s −1 for VS3 and VS4. During the remaining months, the decrease of the AHQ can be observed, even by 49 m 3 s −1 in August for VS2). Since the obtained range of the AHQ values of the in-flow to the studied reservoir is similar to the conditions in the riverine system, the use of the AdH/PTM model was assumed as appropriate.

The Dobczyce Reservoir simulations
As a result of the AdH/PTM simulations, the percentages of each SMAG, SILT, and CLAY fractions that were deposited in individual zones of the reservoir, under the baseline and variant scenarios, have been determined (Fig. 8). Generally, on average over 86% of all fractions flowing into the DR have been deposited within the reservoir in the baseline scenario; however, this process displayed a strong seasonal pattern, and also varied from fraction to fraction. For SMAG, the presence of particles belonging to this fraction has been detected generally between May and September, with 26-47% of these particles deposited only in zone B of the DR, and 52-70% in C. The SMAG portion introduced into zone D was small, at a level of 0.2-3.8%, while the outflow of these particles from the DR (OUT) has been observed only in June, reaching 0.7% of the total SMAG particles. For the SILT fraction, the  A and B); however, the share of particles in each zone displayed a seasonal pattern. From October to April, SILT particles were mostly retained in the first zone of the DR (A; 42-71%), while their share in the next zone was relatively smaller (B; 27-48%). In this period, the share of this fraction reaching zone C was at the level of 2-9%, while their further transport into zones D and OUT was almost negligible. In the remaining months of the year, May-September, the SILT particles have been mainly deposited into zone B (67-72%), and partially into C (16-23%). The transport and deposition of this fraction into zones D and OUT was at a level of 1-4%. For CLAY, a similar seasonal pattern could be observed; however, the localization of these particles' final deposition zones differs from the SILT fraction. Generally, between October and April, these particles are deposited into zones B and C (23-47% and 30-40%, respectively), while the share reaching D is estimated at 4-13%. Also, the elevated share of CLAY particles transported downstream from the reservoir (OUT) should be noticed in October, November, and March (17-35%). For the May-September period, this process of particle transport downstream from the reservoir dominates, with 46-60% of the total number of CLAY particles estimated in the OUT zone. The remaining share of these particles is distributed between zones B, C, and D, with higher values for C (19-23%). The impact of the adopted variant scenarios on the SMAG is expressed by the extension of the period when these particles are transported and deposited into the DR. For the baseline scenario, these processes were confined to May-September, while under the VS1-VS3 scenarios, also April and October are marked with SMAG deposition into the DR zones B and C. While under VS4, March is also included in this process. As for other zones of the DR, the share of these particles was either 0 (A), or at the level similar to the baseline scenario (0.1-3.3%). For the SILT fraction, similar patterns of temporal changes could be expected under VS1-VS4. Although the majority of these particles will still be deposited into the first two DR zones, the period when particles are retained in zone A will be shortened from October-April to November-March and favor a higher percentage of particles settling into zone B (on average by 21%). For the CLAY fraction, the variant scenario impact is expressed mainly by extension of the period in which these particles flow downstream of the DR (OUT). In the baseline scenario, this process lasted from March to November, while in the VS1-VS3 scenarios, December and January are also marked as meaningful for the CLAY particle transport to the lower part of the Raba River. In turn, from May to September, the amount of the CLAY fraction reaching OUT will decrease by an average of 9%, depositing these particles mainly into zones B, C, and D.

Discussion
Dammed reservoirs pose a huge impact on entire aquatic ecosystems disrupting many natural processes, including sediment transport within the river through its effective entrapment (e.g., Sundborg 1992;Vörösmarty et al. 2003;Sedláček et al. 2017;Dong et al. 2019). Excessive sediment deposition not only threatens reservoir capacity and functionality, but also creates an additional threat to the reservoir water quality since sediment particles are considered effective carriers of contaminants (e.g., Szalińska et al. 2013;Zhu et al. 2019;Jalowska and Yuan 2019;Aradpour et al. 2020;Babek et al. 2020). For holistic assessment of sediment inflow impact on a reservoir, exhaustive information on its catchment as a sediment source is required. In the current study, such information has been collected through previous research performed on the Raba River catchment for the land and river phases of the sediment transport (Orlińska-Woźniak et al. 2020a, b;Szalińska et al. 2020). The mountainous part of this catchment, located upstream from the Dobczyce Reservoir, is very prone to erosion; therefore, multiple debris flow barriers were constructed in the 1970s along the main river and its tributaries (Wyżga et al. 2020;Orlińska-Woźniak et al. 2020a). However, currently, the majority of these structures are overfilled and require renovation (Korpak et al. 2008;Lenar-Matyas et al. 2015). The mountainous character of this area also demonstrates the hydrometeorological conditions typical for the entire area of the Polish Carpathian Mts., with average rainfall values reaching 1600 mm accompanied by frequent torrential rains, especially during the summer. These conditions combined with a high share of slopes above 25%, and many areas still used for agriculture (41%), cause intensive loss of soil fractions from the land phase of the catchment. Moreover, they are the reasons behind large fluctuations in flows and frequent floods affecting sediment transport in the riverbed (Kijowska-Strugała and Kiszka 2018). These features had a decisive impact on sediment loads introduced into the DR simulated for the Myślenice calculation profile with the first module of the DNS digital platform-SWAT.
The obtained results showed high loads of SMAG, SILT, and CLAY in the baseline scenario amounting to an average of 184 t m −1 and showing elevated temporal variability, with an average coefficient of variation (CV) of 153%. During the year, two periods can be distinguished: May-September when over 76% of the total annual sediment load flows through the Myślenice profile, and October-April when the rest of this load is delivered into the DR. As for the share of the individual fractions, mineral fine particles constitute over 99% of the total load, with CLAY and SILT contributing 82% and 17%, respectively. Such a dominant share of these fractions results from the natural conditions in the studied catchment (Szarek-Gwiazda and Sadowska 2010), where clay soils cover more than 95% of the area (Fig. 2). Moreover, frost lasting more than 100 days per year in this part of the Carpathian Mts. also promotes erosion of soils, especially of those finegrained. The high content of fine fractions facilitates groundwater capillary forces, causing exfoliation of the surface layers during periods of freezing temperatures (Augustowski and Kukulak 2017). Apart from the two mineral fractions (CLAY and SILT), a small amount of mineral-organic agglomerates (SMAG) has also been observed (0.2%), mostly during the May-August period. The period of their appearance coincides not only with the period of high flows, but also favorable conditions for the biomass production in the river (e.g., temperature, solar radiation, and agricultural activities) (Lee et al. 2019;Hoffman et al. 2020).
The subsequent fate of the sediment particles downstream from the Myślenice calculation profile has been tracked with the second module of the DNS digital platform, AdH/PTM. Due to the RR width increase in the outlet section and the water flow decrease, zone A favors deposition of larger particles. Although due to the choice of the Kodoatie method, fractions above 0.062 mm have not been considered in our study, their presence in this system is visible in the form of vast sandbanks in the mouth of the RR (Fig. 3). This phenomenon is caused by a combination of multiple natural features and processes, including presence of less fine-grained soils in the last sub-catchment of the RR (Fig. 2), and urban use in this part of the area, supplying sand particles due to construction work and winter maintenance of the local roads. Moreover, other mineral particles lose their driving force of advection and begin to settle in this zone, especially during the October-April period, when the average flow (avg. 35 m 3 s −1 ) favors the deposition of 60% of all sediment particles (56% of SILT and 4% of CLAY) (Fig. 9). With the increasing flow, up to 114 m 3 s −1 on average during the May-September period, sediment deposition is limited to only 6% of their total load.
In zone B, the DR water level fluctuations result in alternating cycles of sedimentation and resuspension in this part of the DR (backwater) . As a consequence, deposition of particles not settled in zone A is favored here, especially of the SILT larger fractions (Fig. 9), which are retained up to 70% during the May-September period. For the finer particles of the CLAY fraction, this zone is also the major deposition zone (34%), especially during the lower flow period (October-April). Furthermore, over 40% of the SMAG fraction is deposited during the May to September period in this zone. Due to their lower density, SMAG particles are larger than mineral particles of the same mass which results in their slower settling velocity and higher drag force (Hoffmann et al. 2020). As of which, SMAG particles can be transported much further than mineral particles of the same mass, which explains the lack of their deposition in zone A.
In zone C, with depth increasing to over 9 m and flow velocity continuing to decline, the sediment-carrying potential becomes even weaker. This zone is the last which receives large amounts of SILT fractions (up to 18% on average in May-September), and accumulates over 57% of the SMAG particles. Moreover, an apparent pattern of "bedforms," created by the particle sedimentation, is noticeable in this zone (Fig. 9). Localization and distribution of these bedforms correspond to the variations in reservoir-floor relief visible in a high-resolution bathymetric map ( Fig. S3 and animated GIF file); however, confirmation of their presence requires further research.
The last zone of the DR (D) is characterized by the greatest depth (over 16 m) and the presence of the drinking water intake. During high flows (May-September), even more than 60% of the CLAY fraction load flows into this zone, which requires special attention due to the adsorption of contaminants on fine particles (Szalińska et al. 2013;Palma et al. 2015;Rügner et al. 2019). Since only 10% of this load settles into this part of the DR, it is not really prone to the loss of capacity. However, the remaining particles are transported downstream of the DR (OUT) to the lower part of the RR catchment.
The impact of the future climate changes, for the land and river phases of the sediment transport in the RR catchment, has underlined the role of the spring and winter season, mainly due to the increase of precipitation (Orlińska-Woźniak et al. 2020a, b;Szalińska et al. 2020). Moreover, the detailed simulations (Hachaj et al. 2021) have shown that such increases of sediment yields and loads could not be attenuated by the forecasted land use changes in this area, e.g., 30% increase of forest area at the expense of agricultural land use. In the current study, where the climate change predictions with a monthly step were applied, a more detailed temporal pattern of sediment transport has been elucidated. In particular, a considerable increase of precipitation observed especially in April (up to 81% under VS4),but also in December,February,and October (up to 76% under VS1) in all the discussed scenarios (Fig. 5), will ultimately affect the mineral and organic fraction of sediments delivered into the Myślenice calculation profile. Currently (baseline scenario), during these 4 months, a total of just over 18 t m −1 of SILT and 934 t m −1 of CLAY flow through this profile, whereas the future climate changes will lead in these months to an average SILT load increase by 154 t m −1 (VS3), and by more than 2000 t m −1 (VS4) for CLAY. Since the pattern of precipitation changes will be reversed during May-September (Fig. 5), a decrease of sediment loads in the Myślenice profile can be expected in these months, when compared to the baseline scenario. Indeed, the load of transported fractions will be reduced in these months by 430 and 1160 t m −1 for SILT and CLAY, respectively. These monthly changes require further attention; however, it has been already recognized that this area is likely to experience substantial climate exposure leading to further alterations in the ecosystem (Hlasny et al. 2016).
The climate change scenarios applied to the reservoir simulations emphasize a fundamental difference between both parts of the discussed system. While even a small change of parameters (such as precipitation, temperature, or land cover) causes a rapid hydrological reaction of the mountainous river catchment and consequently alters the intensity of erosion, and sediment transport and fate, the response of the reservoir is very disparate. The flow velocity field, which is crucial for the range of the sediment particles transport, remains sensitive to flow rate changes only in zones A and B of the DR (river and backwater). Therefore, these two zones will remain the main depositional sections of the DR, especially for larger particles, even under flood conditions (Szlapa 2019). Directly behind the backwater zone, the sudden extension of the reservoir's cross-sectional area and a slowdown of the water flow radically change the sediment transport possibilities in further parts of the DR. Therefore, when the SILT fraction behavior in the DR is discussed, only changes in particle numbers and a slight shift in temporal pattern in zones A and B can be noticed. Regardless of the selected scenario, a pronounced increase, even by 800%, of the SILT fraction can be expected in these two zones (Fig. 10). As for the monthly distribution, since the high flow period has been extended from May-September in the baseline scenario, to April-October in the V1-V4 scenarios, therefore, a similar extension in retention effectiveness of SILT particles in zone B is visible (Fig. 10). This leads towards increasing the rate of the reservoir volume loss in time, and consequently causes a transition of the backwater zone (B) towards the dam at the cost of the intermediate zone (C). The lacustrine zone (D) stays unaffected by this SILT rate increase. Moreover, a similar pattern is also noticeable for the SMAG particles. Although their spatial pattern remains the same under the adapted scenarios, i.e., their appearance is restricted only to zones B and C, their temporal distribution is noticeably different, as it is extended to include even March under V4. It should be also noted that an extended period of allochthonous organic material delivery to the reservoir can promote further development of the autochthonous organic matter due to the favorable conditions in the DR . As for the CLAY particles, the reservoir response to the applied scenarios will completely differ. The altered high flow period, April-October, will increase the mobility of these particles, pushing them to flow to the last reservoir zone (D), and even downstream from the reservoir (OUT). Generally, the share of this fraction flowing to zone D will be increased by even 30% (December) when compared to the baseline scenario, but only an average of 10% of them will settle there. Since this zone is used as a source of drinking water, the extended presence in this part of the reservoir should be further investigated due to their role as a contaminant carrier, as previously discussed .
Although the adopted modeling approach has been proved as very suitable in tracking sediment particles from their source to the deposition place under a variety of scenarios, one should remember about possible limitations. As for the SWAT model, the land phase of simulations has been rarely contested by the scientific community, whereas a representation of the stream processes in its bed phase is considered as oversimplified (Sarkar et al. 2019).
It should be also noted that in our simplified approach, the nonlinear phenomena affecting suspended sediment transport, such as density currents, has not been taken into consideration. However, their presence and impact will be assessed through further research. Moreover, due to the semi-empirical character of this model an impact of adopted calculation equations and corresponding coefficients (Yen et al. 2017) could be a possible source of bias. Another source of uncertainties for sediment tracking, under the adopted scenarios, are climate change predictions. The currently available forecasts are usually characterized by a large spatial scale and a low time resolution (Fatichi et al. 2016). Therefore, their application for catchments with limited areas where extreme weather events occur should be only used as a general indication of predicted changes.
The main limitations of the applied modeling approach for the reservoir part are related to the characteristics of both modules. As for the hydrodynamic AdH model, it generates a planar, two-dimensional velocity field, impacting the sediment transport and altering grain deposition zones. However, in our approach, the vertical velocity component has been calculated in the PTM module, therefore, the 3D velocity field could be used to track representative particle transport. As for the PTM module, its biggest limitation is an inability to simulate mass transport. With this tool, one can observe the transport of characteristic sediment particles and indicate the places of their deposition, but no direct information about the transported or deposited mass could be provided. In order to switch from a qualitative to a quantitative approach, it is necessary to characterize the sediment inflow carried from the river to the reservoir, which will be the goal of further research.

Conclusions
In the current study, we have tracked the selected fractions of sediment particles (SMAG, SILT, and CLAY) from their source of origin (Raba River catchment) to the deposition area in a dammed reservoir (Dobczyce Reservoir). This research was possible due to the combined performance of two models (SWAT and AdH/PTM) under the umbrella of the Macromodel DNS digital platform. Such an approach could be easily adopted in similar catchments, featuring dam reservoirs. Especially, since the proposed solution allows to overcome the insufficient performance of the reservoir module included in the SWAT model itself. The obtained results underlined the fundamental difference between both hydrological units, i.e., river catchment and dammed reservoir. Although the studied river catchment is extremely prone to erosion, and under the forecasted climate change will respond in increasing sediment loads, this response will eventually be attenuated by the reservoir. Due to the very fortunate location and natural setup of the studied reservoir, the two first zones will maintain their trapping role for larger particles (SILT), even during periods of highly increased sediment delivery periods under the adopted climate and land use scenarios, although special attention must be paid towards acceleration of the capacity loss for this part of the reservoir. As for the finer particles (CLAY), their increased mobility into the reservoir is clearly visible both under shortand long-term scenarios which raises concerns due to their contaminant binding affinity, and possible impact on drinking water. Moreover, a monthly resolution of our scenarios provided detailed analysis of the studied grain size fractions, which are necessary when appropriate measures to evaluate sediment fate in a catchment need to be discussed.