Modeling the Suspended Sediment Transport in a Very Wide, Shallow, and Microtidal Estuary, the Río de la Plata, Argentina

The impact of the diverse processes responsible for the distribution of suspended sediments in the wide, shallow, and microtidal Río de la Plata estuary and the adjacent shelf is studied by means of a set of process‐oriented numerical simulations. Model results show that the large width and the geometry of the estuary play a major role in the sedimentation processes. The widening and deepening of the estuary drives a significant reduction in current speeds at (i) the confluence of the tributaries and (ii) downstream Barra del Indio Shoal. Thus, sediment deposition enhances downstream those areas. Even though tides are of small amplitude in the study area, they have a significant impact on lateral mixing and resuspension of the bottom sediments. Resuspension augments the concentration of fine sediments in the layers close to the bottom, but tidal energy is not enough to lift them to the surface. Winds (which can be quite strong over this area) enhance horizontal mixing, smoothing the pattern produced by tides. Wind waves increase the concentration of sediment by vertical mixing, and their effect is most evident along the southern coast where wind waves lift the sediments resuspended by tides to the surface. The estuarine circulations associated to the bottom salinity front acts retaining sediments upstream the Barra del Indio Shoal, where estuarine currents and flocculation play an important role in sediment deposition.


Introduction
The Río de la Plata Estuary (RdP Estuary; Figure 1) is one of the most turbid estuaries in the world, with extreme sediment concentrations above 400 mg/L (Framiñan & Brown, 1996;Moreira et al., 2013). The input of suspended sediments to the world's oceans is between 15,000 and 20,000 Mt/year (Walling & Webb, 1996), of which more than 1% (between 80 and 160 Mt/year) is contributed by the sediment load of the RdP Estuary (Menéndez & Sarubbi, 2007;Urien, 1972). These sediments reach the RdP Estuary as suspended sediment depending on the sediment size and the flow characteristics (van Rijn, 1984a). About 90% of the sediments that flow into the RdP Estuary as suspended sediments are silts and clays; the remaining 10% corresponds to very fine sands (Sarubbi, 2007). Bottom deposits are composed of sands in the upper estuary, silt in the intermediate estuary and clayey silts in the exterior zone .
The RdP Estuary has large social, ecological, and economic importance for its coastal countries, Argentina and Uruguay. The capitals of both countries, Buenos Aires and Montevideo, respectively, and a number of harbors, resorts, and industrial centers are located on its margins and influence zone. The estuary is the main source of drinking water and an important amusement area for the millions of inhabitants in the region. Because it is rich in nutrients, the RdP Estuary has abundant and diverse fauna and it is home to large fisheries. Moreover, it has the unusual feature of being both a spawning and nursery area for several coastal species (Acha & Macchi, 2000;Berasategui et al., 2004Berasategui et al., , 2006. Samborombón Bay is one of the most important wetlands in Argentina and is home to a number of species of fish, turtles, crabs and migratory birds (Canevari et al., 1998;Lasta, 1995;Volpedo et al., 2005).
The fresh water plume of the RdP Estuary flows more than 500 km to the north (Campos et al., 1999) carrying sediments, carbon, and nutrients, which have obvious impact on the adjacent shelf and ocean processes. In addition to these land-ocean exchanges, many environmental characteristics and processes in the RdP Estuary and the adjacent shelf are linked to the distribution, variability, and major drivers of the dynamics of suspended sediments (Moreira et al., 2013). The most significant issues related to the load of sediments in the RdP Estuary include the need for optimizing dredging operations (Cardini et al., 2002), designing effective flood alerts (Moreira et al., 2014), understanding geomorphological changes (Codignotto et al., 2011), pollution (Avigliano et al., 2015;Colombo et al., 2005Colombo et al., , 2007, benthic ecology (Gómez-Erache, 1999), primary productivity (Gómez-Erache et al., 2004;Huret et al., 2005), and fisheries (Acha et al., 2008;Jaureguizar et al., 2003Jaureguizar et al., , 2003Jaureguizar et al., , 2008. Despite their importance, sedimentological processes in this area still need to be fully understood. Based on observations collected along the navigation channels of the RdP Estuary, Ottmann and Urien (1966) considered for the first time the dynamic phenomena that might occur in the RdP. They were, in their view, the fluvial discharge, the astronomic tide, and the coastal drift. Parker et al. (1986b), based on observations of the bottom sediments, speculated that at the upper and intermediate RdP Estuary the dispersion of sediments would be determined by the loss of energy of the flow of the tributaries, the waves effect, and the tide. This would originate the Paraná River Delta, which extends to the southeast in the form of lobes formed by medium to fine sands and silts, without clay elements. They also postulated that between La Plata and El Codillo, erosion increases due to high tides. For the first time they noted the difference in the characteristics of the sediments transported by the diverse tributaries to the RdP Estuary and related them to the bottom sediments distribution.
One of the main difficulties regarding the study of suspended sediments dynamics in the RdP Estuary has been the historical lack of consistency of the available observations. The samples had been collected at a few places (mainly along the navigation channels) and different moments (in occasions with decades of difference), using diverse techniques, what makes impossible their comparison in such a dynamical system. This situation was in part overcome when a new set of hydro-sedimentological observations was collected during six oceanographic cruises in the frame of the FREPLATA/FFEM Project (Environmental Protection of the Río de la Plata and its Maritime Front/French Fund for the Global Environment Experiment) between November 2009 and December 2010 (Simionato, Moreira, Piedra Cueva, et al., 2011;Simionato, Moreira, Re, et al., 2011). Moreira et al. (2013) analyzed those data, together with 10 years of daily intermediate resolution (1 km) Moderate Resolution Imaging Spectroradiometer (MODIS)-Aqua observations processed for surface suspended matter. Those authors showed that the mean surface concentration of suspended matter and turbidity maximizes at the upper estuary and along the southern coast of the upper and intermediate RdP Estuary. At the outer RdP Estuary the concentration and the turbidity decay seaward. Maximum concentrations are observed close to Punta Piedras and Punta Rasa (Moreira et al., 2013). This way, two significant sediments concentration gradients are observed at the RdP Estuary: (i) one oriented across the estuary axis at the upper and intermediate RdP Estuary and (ii) another one along the estuary axis at its the outer part.
Several hypotheses were raised by Moreira et al. (2013) to explain the observed mean distribution pattern; in particular: 1. The sediments concentration at upper and intermediate estuary would be mainly controlled by the solid discharge, the morphology, and the tides, with a lesser effect of the winds and wind waves. 2. The relative maximum observed along the southern coast of the upper estuary would be related to the combined action of the tributaries to the RdP Estuary, which discharge 73% of the solid load through the Paraná River (at the south), and the stronger tidal currents  along the southern coast. 3. The maximum observed close to Punta Piedras and Punta Rasa is located over the one of the areas of maximum tidal currents and tidal dissipation by bottom friction (Simionato, Dragani, Meccia, & Nuñez, 2004). Nevertheless, Moreira et al. (2013) suggest that tides act resuspending the sediments near the bottom, but wind waves mix the water column vertically lifting the sediments to the surface. 4. At the outer estuary, in the region of the bottom salinity front, decantation and flocculation would become more important; here tides and waves would act resuspending sediments and mixing them, but playing a secondary role because depth is larger. 5. The secondary maximum of the mean suspended matter concentration at the outer estuary, close to Punta Rasa, would be connected to another area of maximum tidal currents and tidal dissipation by bottom friction. In this area wave action would become more important, as depth decreases.
Even though the above mentioned hypotheses can hardly be "proved" in a huge and variable system as the RdP Estuary, numerical models are excellent tools to test them. Models constitute a "laboratory" that permit "including" and/or "removing" the different forcings, or modifying the natural environmental conditions (bathymetry, geometry, etc.), in order to analyze their individual or combined effects. In this sense, the aim of this paper is, by means of process oriented numerical simulations, to analyze the role of the diverse forcings and the environmental conditions of the estuary, which determine the distribution of suspended sediment in the RdP Estuary, in order either verify or reject the conjectures of the previous works through a more robust tool. This work also aims to provide some insight on several processes, which effect is almost unknown, as the estuarine circulation, the effect of the large estuary widening at the imaginary line between Punta Piedras and Montevideo, and the presence of the Barra del Indio Shoal.
The hydro sedimentological model used for the simulations is a regional implementation of the Model for Applications at Regional Scale (MARS; Lazure & Dumas, 2008) originally built during the FREPLATA/ FFEM Project (Simionato, Moreira, Re, & Fossati, 2011) and further developed, calibrated, and validated for the RdP Estuary by Fossati (2013), Fossati et al. (2014), and Moreira et al. (2016). In the version applied in this paper, the concentration of three different kinds of textures (fine sands, silts, and clays) is considered. The model is run under realistic bathymetry and geometry but idealized forcing conditions, adding the forcings (discharge, tides, winds, and wind waves) one at a time. Then two additional simulations ignoring the upper and intermediate estuary bottom topography and the effect of the estuarine circulation are run.
This paper is organized as follows: section 2 provides a description of the most relevant features of the study area; section 3 describes the characteristics and implementation of the model and the different datasets used; the results of the simulations are presented in section 4. Finally, in section 5 we summarize and discuss the results and present a conceptual model of the dominant sedimentological processes that occur in the different parts of the estuary in relation to action of the different forcing.

The Study Area
The RdP Estuary is one of the largest estuaries in the world (Shilkomanov, 1998), located on the eastern coast of southern South America at approximately 35°S ( Figure 1). It has a northwest to southeast oriented funnel shape of approximately 300 km long, which narrows from 220 km at its mouth to 40 km at its upper end (Balay, 1961). The estuarine area is 35,000 km 2 and the fluvial drainage area is 3.1×10 6 km 2 (Framiñan et al., 1999). The RdP Estuary displays a complex geometry and bathymetry. A complete description of its morphology can be found in Depetris and Griffin (1968), López Laborde (1987a, 1987b, Ottman and Urien (1966), Parker et al. (1986aParker et al. ( , 1986b, Parker and López Laborde (1989), and Urien (1966Urien ( , 1972. The estuary is naturally divided into two regions by the Barra del Indio Shoal, a shallow area that runs perpendicular to the estuary axis between Punta Piedras and Montevideo ( Figure 1). The upper region is mainly occupied by fresh water. Seaward the shoal is the Maritime Channel, a wide depression with depths ranging from 12 to 14 m in the north to 20 m in the south. It separates Samborombón Bay to the west from a region of banks known as Alto Marítimo (with depths ranging from 6 to 8 m) and the Rouen Bank (with depths between 10 and 12 m). The Oriental Channel, the deepest zone in the estuary with depths of up to 25 m, extends along the Uruguayan coast. Samborombón Bay is a very shallow and extensive area with depths ranging from 2 to 10 m extending south of Punta Piedras and up to Punta Rasa.
The RdP Estuary is the second largest basin in South America, after the Amazon. It exhibits a mean discharge of around 22,500 m 3 /s, and extreme discharges as high as 90,000 m 3 /s and as low as 8,000 m 3 /s, related to the El Niño-Southern Oscillation cycles event (Jaime & Menéndez, 2002). The main tributaries to the RdP Estuary are (i) the Paraná River, which reaches the estuary after forming a large delta with two main branches, the Paraná Guazú-Bravo through the north, and the Paraná de las Palmas through the south, and (ii) the Uruguay River with lower liquid and solid discharges, which converges to the estuary to the north of the Paraná Guazú-Bravo River (Menéndez & Re, 2009). The Uruguayan coast is mostly affected by waters of the Uruguay River and the Argentinean coast, by waters of the Paraná de las Palmas, whereas the waters of the Paraná Guazú-Bravo occupy the central part (Simionato et al., 2009).
Water stratification is controlled by the confluence of highly buoyant continental discharge advecting offshore, lying on top of denser shelf water that intrudes into the estuary as a topographically controlled salt wedge. This salt wedge is typically between 100 and 250 km long (Guerrero et al., 1997) and defines a bottom salinity front, over the Barra del Indio Shoal following approximately the 10 m isobath (Guerrero et al., 1997).
The RdP Estuary is a microtidal system. Tidal waves associated with the South Atlantic amphidromes reach the Continental Shelf while propagating northward (Glorioso & Flather, 1995O'Connor, 1991;Simionato, Dragani, Meccia, & Nuñez, 2004). Over the shelf, the geographic setting modifies their direction and they enter the estuary mainly from the southeast (Simionato, Dragani, Meccia, & Nuñez, 2004). Tidal amplitudes are not amplified towards the upper part. The estuary is long and converges only at its innermost part, where it is extremely shallow and bottom friction plays a fundamental role in controlling wave amplitude (Simionato, Dragani, Meccia, & Nuñez, 2004). Tidal amplitudes are larger and currents are stronger along the southern coast of the estuary than along the northern one, with maximum tidal current speeds and energy dissipation by bottom friction at the northernmost and southernmost limits of Samborombón Bay, Punta Piedras, and Punta Rasa respectively (Simionato, Dragani, Meccia, & Nuñez, 2004).
Subtidal transport in the RdP Estuary is highly influenced by the geometry and bathymetry, the rotation of the Earth, and particularly, by the wind . At the upper estuary, after discharge, the flow concentrates along the deep North and Intermediate channels. As the fresh water plume reaches the central part of the estuary, it is deflected to the north by the rotation of the Earth (Coriolis effect). Even though the Arquimedes and English banks divide the flow into two branches in the outer estuary, in the absence of winds, they meet again after flowing through this region. The transport increases (reduces) for higher (lower) runoff conditions, but the described pattern is preserved (Simionato, Dragani, Nuñez, & Engel, 2004).
A schematic representation of the path of the plumes of the main tributaries to the Río de la Plata, which in turn carry the sediments, was presented by Simionato et al. (2009). Whereas the waters of the Paraná de las Palmas flow to the south, along the southern coast, the waters of the Paraná Guazú-Bravo and Uruguay partially mix after the confluence and flow along the northern coast. Part of the waters of the Paraná Guazú-Bravo flow between the other two masses, particularly when the discharge of Uruguay is high.
The hydrodynamic module of MARS code solves the primitive equations of an incompressible fluid, with the assumption of hydrostatic balance and the Boussinesq approximation. It is based on traditional finite differences on an Arakawa C grid (vertical and horizontal). The vertical coordinate is of generalized sigma-type. The originality of the MARS code relies in the treatment of the barotropic mode, which is evaluated semiimplicitly, thus eliminating time splitting and allowing direct coupling between the barotropic and baroclinic modes.
In addition, open boundary conditions are expressed in the center of the grid (forcing in level) because MARS was originally applied to model tide driven currents.
The sediment dynamics module integrated into MARS derives from the SIAM model (Le Hir et al., 2011). It is based on a simple transport equation that estimates sediment advection and dispersion based on the concentration and falling speed.
The model includes the processes of erosion, transport in suspension, and simultaneous deposition of different classes of sedimentological particles: gravel, sand, silts, and clays. The bed is described as thin layers characterized by a distribution of those classes. Sediment composition and concentration change with the action of erosion, deposition, and consolidation processes. Two unique features characterize the code: (i) first, the treatment of mixtures of different size classes to better simulate natural sediments, usually heterometric, and (ii) secondly, the inclusion of an algorithm of simultaneous compaction of the different sediment classes, which allows for instance to reproduce deposits that alternate layers of sand or muds dominated sediments. The sedimentological module also incorporates a variant to the classic treatment of sand-mud mixtures. Finally, the sediment settling velocity is expressed as a function of sediment concentration as a way of taking flocculation processes into consideration (Le Hir et al., 2001). The reader is referred to the above mentioned papers and the on line documentation (http://www.ifremer.fr/docmars/html/documentation.html) for further details on model formulation and equations.
The implementation for the RdP Estuary was made using two nested grids ( Figure 2). The father grid "Rank 0" covers the southwestern South Atlantic Ocean, from 69.35 to 45.50°W and from 25.50 to 54.80°S, with resolutions of 0.10°(~10 km) and 0.15°(~12 km) in latitude and longitude, respectively, and one vertical layer. This barotropic and bidimensional implementation permits the free propagation of tidal waves and provides open boundary conditions for free surface height to the child domain ("Rank 1"). The latter covers the RdP Estuary and a large part of the adjacent continental shelf, spanning from 32.900 to 38.139°S in latitude and 60.670 to 50.620°W in longitude, with a regular grid of 0.027°(~3 km) resolution. In the vertical, 10 not equidistant sigma levels were used, corresponding to 0.05, 0.2, 0.35, 0.5, 0.65, 0.75, 0.85, 0.95, and 0.97 of the total depth at every grid point. The child model is three dimensional and baroclinic and was rotated in the direction of the continental shelf to make the application of boundary conditions easier and simulations cheaper for the computational point of view by avoiding large depths outside the shelf break.
The bathymetry of the Río de la Plata and the continental shelf (Figures 1 and 2) was obtained from data provided by the Hydrographic Service of Argentina, which come from digitalization of nautical charts for the continental shelf and the RdP Estuary, and from the National Oceanic and Atmospheric Administration (NOAA)'s National Geophysical Data Center (http://www.emrl.byu.edu/gsda/data_bathymetry_obtain. html) ETOPO5 for the deep ocean. The high-resolution coastline was obtained from NOAA (http://www. ngdc.noaa.gov/mgg/gdas) at 1:250.000 scale.
To provide realistic weather conditions, surface winds and sea level pressure data from the 6-hourly Reanalyses of the National Center for Environmental Prediction/National Center for Atmospheric Research (NCEP/NCAR; Kalnay et al., 1996) were used. The spatial resolution of this data set is 2.5°× 2.5°. A constant wind drag coefficient of Cd = 0.0016 was used.

Journal of Advances in Modeling Earth Systems
Liquid and solid discharges of the three main tributaries, Uruguay, Paraná Guazú-Bravo, and Paraná de las Palmas (Figure 1), were kept constant in time in all the simulations. Data were obtained from the "Base de Datos Hidrológica Integrada" of the "Subsecretaría de Recursos Hídricos de La Nación" of Argentina (http:// bdhi.hidricosargentina.gov.ar). They were collected at the city of Concepción for the Uruguay River, at Brazo Largo for the Paraná Guazú-Bravo, and at Zárate for the Paraná de las Palmas, at nonregular intervals during the period ranging from February 1993 to February 2016. Mean values of liquid and solid discharge used in the simulations, and their standard deviation are shown in Table 1.

Journal of Advances in Modeling Earth Systems
Three types of sediment were considered in the simulations: (i) fine sands, with settling velocities between 0.001 and 0.002 m/s; (ii) silts, with settling velocities between 0.000015 and 0.000075 m/s; and (iii) clays, with settling velocities between 0.000010 and 0.000015 m/s. The shear stress of erosion and deposition was fixed at 0.1 N/m 2 for fine sand and at 0.07 N/m 2 for silt and clay, following previous studies for this and other areas of the world (e.g., Gibbs, Mathews, & Link, 1971;van Rijn, 1984;Fossati, 2013). When considered, the bottom erosion parameter E 0 was set to 0.0001 Kg/m 2 s The mean liquid and solid discharges were distributed into the textures considered in the simulations following the relative percentages suggested by Sarubbi (2007); results are shown in Table 1. According to these observations, the mean runoff, was 26,209 m 3 /s, 16% above the long term (more than a century) mean discharge of 22,500 m 3 /s reported by Jaime & Menéndez (2002). The mean solid discharge (75 Mt/year) is similar to the mean value reported by Urien (1972), but almost half of the estimates by Menéndez & Sarubbi (2007).
To initialize the MARS model, it is necessary to specify the initial percentage of every texture that composes the bottom layer. Here also three different types of sediments were considered (clay, silt, and sand), according to the characteristic average diameter suggested by Wentworth (1922). The gravel was not taken into account in the simulations since its concentration is not significant in the RdP Estuary and it responds to different transport processes. The initial percentage distribution was obtained from the data collected during the FREPLATA/FFEM Project analyzed by Moreira et al. (2016) and complemented with the historical distribution reported by López Laborde (1987b).
MARS in this configuration was calibrated and validated for the hydrodynamic module in both domains by Simionato, Moreira, Re, andFossati (2011), Fossati (2013), and Moreira (2016). The sedimentological module was also applied, calibrated, and validated by Fossati et al. (2014) and Moreira (2016). The simulations of Fossati (2013) and Fossati et al. (2014) show that the model results well compared to time series collected during the FREPLATA/FFEM Experiment at three fixed points of the RdP Estuary and with observations gathered at Montevideo Harbor.
Wind waves had to be necessarily included our simulations because of their crucial role in forcing sedimentological processes. Given the lack of long time series of wave data in the RdP Estuary, we derived an empirical relationship between wind speed and direction data recorded by an Oceanographic Buoy and (nondirectional) wave data recorded at Torre Oyarvide (red dot in Figure 1). These data were collected in the frame of the FREPLATA-FFEM Experiment during 2012 (Simionato, Moreira, Re, & Fossati, 2011), and the recorded period spanned from 26 October to 17 December 2012, with a total of 412 observations. To obtain the empirical equation linking wind speed to wave height, wind speed collected by the Oceanographic Buoy was divided into 15 intervals of 1 m/s, between 0 to 15 m/s, and, the average wave height observed in Torre Oyarvide was calculated for every interval. Then, a polynomial was adjusted by minimum squares. The obtained relation was where H is the wave height in meters, and Wndsp is the wind speed in m/s. The correlation coefficient was 0.922, which indicates good agreement. This equation was added to the MARS code in order to have wave heights for every time step and grid point of the domain. Because the wave data available are not directional, the direction of wave propagation was assumed to be equal to the direction of wind. This condition is characteristic within the RdP Estuary, where waves are mainly local (Dragani & Romero, 2004). Since the wave gauge was located at the intermediate region of the RdP Estuary, Torre Oyarvide, it is possible that the empirical relationship overestimates wave height in the upper estuary and underestimates it in the exterior estuary. As the aim of this paper is not to realistically reproduce any particular situation, but to understand the influence of the diverse forcings and, in particular, of wind waves on the suspended sediment dynamics, this approach can be regarded as valid.
Incorporating wind waves in the simulations also requires the parameterization of the maximum orbital velocities, that is, the speeds of the wave generated currents produced just above the bottom boundary layer. The equation for the bottom orbital velocity was derived from the theory proposed by Soulsby and Smallman (1986) following the exponential approximation (equation (2)) for shallow waters proposed by Soulsby (1987). Accordingly, the orbital velocity in the bottom (u b ) is where H is the wave height in meters, h is the depth in meters, T is the wave period in seconds, and g is the acceleration of gravity. The period T was set to 4.56 s, the average period of the measured wave periods.
The simulations discussed in this paper were started from rest, and a ramp of 1 day was used to introduce the forcings, to prevent the excessive development of inertia-gravity waves at the beginning of the simulation and, thereby, to reduce the spin-up time. Rank 0 was run for an equivalent of 24 months; after about 10 days, the simulation reached a steady state, but the first 4 months were discarded. The model was run twice, one run considering atmospheric forcing and tides, and the other one with tides only. Hourly results for free surface elevation and barotropic currents at every grid points were stored. To ensure maximum efficiency throughout the simulation, MARS adjusts time steps automatically ensuring stability continuously, according to the CFL (Courant et al., 1928;Kowalik & Murty, 1993) criteria. Time steps ranged from 200 to 600 s (3 to 10 m). Rank 1 was also initialized from rest, and a ramp was applied; however, boundary conditions came from bilinear interpolation of the solutions of Rank 0. The model was run for 20 months in baroclinic mode (ranging from 1 April 2009 to 31 December 2010), considering an initial temperature of 10°C and salinity of 34 PSU; freshwater river discharge was introduced with the same temperature. Sediment distribution reached a steady state after approximately 4 months of integration. Nevertheless, only the hourly solutions of the last simulation year were used for the analysis.

Results
A set of process oriented numerical simulations was performed to study the effect of different forcings (runoff, tides, winds, and wind waves) on the sedimentological process in the RdP Estuary. To favor interpretation, some of the simulations were repeated with and without bottom erosion. Consolidation of the sediments deposited at the bottom was allowed. To simplify comparison, the hourly numerical solutions for the last year were averaged. Continental discharge, tides, winds, and waves were incorporated one at a time in the simulations. This way the different cases are Case 1 (continental discharge), Case 1 (continental discharge + tide), Case 3 (continental discharge + tide + wind), and Case 4 (continental discharge + tide + wind + waves). They are summarized in Table 2.

Comparison of the Simulated Mean Field of Suspended Particulate Matter to Observations
Simulation results permit, for every case, the computation of the suspended particulate matter (SPM) concentration as the addition of the concentrations (in mg/L) of sands, silts, and clays. This variable is convenient because it can be directly compared to SPM observations, more abundant than direct estimations of the concentration of the diverse textures in the RdP Estuary. As a way of validating the model configuration applied in this paper and its capability to reproduce the observed spatial patterns of the distribution of suspended sediments (that is the interest of this work), such a comparison is displayed in Figure 3.
Upper panel of Figure 3 shows (i) the average concentration of SPM derived from direct observations collected during the six oceanographic cruises of the FREPLATA/FFEM Experiment (from Moreira et al., 2013, upper left panel); the red dots indicate the sites where the samples were gathered; (ii) the mean concentration of SPM from 10 years of MODIS satellite images processed by IFREMER (from Moreira et al., 2013, upper central panel); and (iii) the average solution of our simulations with all the forcings included (Case 4, upper right panel). Lower panel of Figure 3 shows the surface SPM concentration along a transect spanning for the estuary head to the estuary mouth along the center of the estuary channel, for every of the cases; the depth (in m) and the observed SPM concentration (from MODIS observations) have been superimposed in the figure. When interpreting the figures, it should be taken into account that differences in the absolute concentrations between the three panels of the upper panel can derive, in a large extent, of the different sources of the data and their limitations. In situ values correspond to an average of six measurements along one year (November 2009 to December 2010), satellite estimations are an average of daily observations (when they were available) over a period of ten years (2002)(2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012), and the model results are an average over one year of idealized simulation with a constant liquid and solid discharges. The largest

10.1029/2018MS001605
Journal of Advances in Modeling Earth Systems differences between the model results and the observed fields occur close to Punta Piedras, where the model provides concentrations that exceed 240 mg/L, what seems to be large when compared to the observed fields. Nevertheless, almost no direct observations are available close to that very shallow area (see the red dots in the upper left panel of Figure 3) and satellite observations of SPM that were obtained from the channel at 550 nm for clear waters and from the channel at 670 nm for turbid waters. A strict flagging is applied to the data through the NASA flags. As the saturation effect can also limit the validity of SPM at 670 nm in the highly turbid waters of the RdP Estuary, the values superior to 98 g/m 3 have been flagged in the procedure (Moreira et al., 2013). It must be noted that there is a possible negative bias in the most turbid areas as the highest concentrations of SPM tend to be flagged (Dogliotti et al., 2011). In spite of the above mentioned limitations, the good qualitative correspondence of the spatial patterns is evident in both panels of the  between the head of the estuary (strongly influenced by the tributaries) and up to about 58°W (imaginary line between Buenos Aires and Colonia); here, the concentration of sediment decreases with the distance from the source; (ii) in the intermediate RdP (until just before the Barra del Indio Shoal, about 57°W) the concentration of SPM increases; (iii) finally, over and downstream the Barra del Indio Shoal, the maximum concentration gradient is observed, and the concentration of SPM decreases to less than 10 mg/L; this region defines the beginning of the exterior RdP, where the concentration of SPM eventually decreases to 0. Note that the Case 4 reproduces the main described features. In general, simulated concentrations are higher than the observed ones over the areas of maximum erosion, what could be explained by the limitations of the data sets used as well as by the idealized nature of the simulations. Keeping that in mind, it can be concluded that the simulations well reproduce the observations and, therefore, are useful to analyze the dominant sedimentological processes in the RdP Estuary. for the bottom layers for Case 2. In this case, forcings were constant continental solid and liquid discharges and tides. The comparison of these results with those of Case 1 (without tide, Figure 4, upper panel) reveals that the intensification of the currents due to tides plays a very important role in the distribution of suspended sediment in the RdP Estuary. In a period of integration equal to that considered in the case without tides, both clays and silts are transported over a much greater distance and reach the outer estuary. The reason for this is the combined effect of enhanced currents, which inhibit sediment deposition and re-suspend them, and the tide associated lateral mixing (or stirring). The highest concentrations of silts at the surface layer (mid upper left panel in Figure 4) appear close to the mouth of the tributaries. Concentrations decrease toward the outer estuary, and a strong gradient is observed over the Barra del Indio region. A relative maximum is observed in the proximity of Punta Rasa. At the bottom layer, a similar distribution pattern is observed, although concentration in the upper estuary is higher. A relative maximum is observed to the east of Punta Piedras. Even with the mixing effect introduced by tides, the paths of the plumes of the tributaries remain separate, although in this case maximum concentrations occur at the mouth of the southernmost branch of the Paraná River, the Paraná de las Palmas. This suggests that tide is very effective in the northern part of the upper estuary, where it mixes the waters of the Paraná Guazú-Bravo and Uruguay Rivers, and inhibits the deposition of fine sediments. Clays at the surface layer (mid upper right panel of Figure 4) are distributed throughout the upper and intermediate RdP Estuary, and concentrations decrease offshore and from the southern to the northern coast.

Sensitivity Studies to the Forcings
A marked gradient parallel to the estuary axis occurs over the Barra del Indio Shoal. A small portion of clays remains in suspension downstream the shoal. Clay concentration is low in the outer RdP Estuary. Similarly to silt, maximum clay concentrations are observed at the mouth of the Paraná de las Palmas River, which points to very effective mixing between the Paraná Guazú-Bravo and the Uruguay rivers, and to the effect of tidal currents keeping sediment in suspension. A relative maximum of clay concentration is also observed at the surface layer in a wide region between Punta Rasa and Punta Médanos, which has no clear connection with other areas. Such lack of connection indicates that this maximum is due to resuspension of sediments rather than to advection. The distribution of the mean suspended clay concentration at the bottom layer (mid upper right panel in Figure 5) shows a similar pattern, although values are higher. This is particularly evident at the Barra del Indio Shoal region, with a relative maximum to the east of Punta Piedras. Seaward Barra del

Journal of Advances in Modeling Earth Systems
To identify the regions where bottom tidal erosion processes are more significant, a new run (Case 2b) was performed in which erosion was neglected (E 0 = 0). Figures 6 and 7 show the difference between the cases with erosion (Case 2a) and without erosion (Case IIb) for the surface (upper panel in Figure 6) and bottom (lower panel in Figure 7) layers. As in the previous figures, results are an average over the last year of the simulation. The difference on the concentration of suspended maximizes in two regions: (i) close to Punta Piedras and extending to Barra del Indio Shoal and (ii) from Punta Médanos to Punta Rasa. At both areas, erosion causes an increment of the suspended sediment concentrations at both the surface and the bottom. In those two regions tidal currents are the largest in the estuary and dissipation by bottom friction maximizes (Simionato, Dragani, Meccia, & Nuñez, 2004).

Case III Simulations Forced With Continental Discharge, Tides, and Winds
Figures 4 and 5 show the average distribution over 1 year of simulation, of the concentration of suspended silts and clays at the surface (mid lower left and right panels in Figure 4) and at the bottom (mid lower left and right panels in Figure 5) for Case 3, in which constant continental discharge, tides, and winds were included as forcings, and erosion of bottom sediments was allowed. Note that results do not correspond to any particular condition of wind direction and speed, because they are an average of one year of simulation; this way, results represent the integrated effect of the forcing over a relatively long period compared to the typical scales of synoptic wind variability (4 days in the study region; Vera et al., 2002). Qualitatively, the distribution of suspended silt concentration at the surface layer (mid lower left panel in Figure 4) is

10.1029/2018MS001605
Journal of Advances in Modeling Earth Systems similar to that of the case without winds (mid upper left panel in Figure 4). Maximum concentrations are observed at the mouth of the tributary rivers. A gradient parallel to the estuary axis is observed over the Barra del Indio Shoal, where concentration decreases from 20 to 0 mg/L over a very short distance. The concentration of suspended silts is relatively larger along the estuary axis and decreases towards the coasts, being lower on the northern coast of intermediate RdP Estuary. The concentration of suspended clays shows a similar pattern, even though concentrations are larger in general. This texture also presents larger concentration along the southern coast and reaches Punta Piedras and Punta Médanos, as in Case 2.
Another way of comparing the cases with (Case 3) and without wind (Case 2a) is by computing the differences between them (middle panel in Figures 6 and 7 for surface and bottom layer, respectively). Winds cause the increment of the concentration of suspended sediments (particularly silts) at the surface over almost the entire estuary as a result of enhanced vertical mixing (mid left panel in Figure 6). In addition, stronger currents inhibit sediment deposition and increase erosion. In particular, the concentration of sediments at the surface is much higher over the Barra del Indio Shoal. Close to the bottom, the concentration of silt is smaller at the southern part of the Barra del Indio and greater at the northern part of Samborombón Bay (mid left panel in Figure 7). These changes result from the action of wind on the fresh water plume of the RdP Estuary that produces significant lateral motion (see, e.g., Simionato, Meccia, Dragani, Guerrero, & Nuñez, 2006;Simionato, Meccia, Guerrero, Dragani, & Nuñez, 2007), which in turn advects sediment and contributes to lateral and vertical mixing.
In summary, the effects of wind are diverse in the different parts of the RdP Estuary. In the upper estuary it enhances vertical mixing and helps maintaining fine sediments in suspension. Consequently, in most of the  Figure 4). At the upper and intermediate RdP Estuary a marked gradient perpendicular to the estuary axis is observed, with higher concentrations along the southern coast. This gradient is stronger and much more defined than in the previous simulation. Over the Barra del Indio Shoal, a marked gradient parallel to the estuary axis occurs, with concentrations decreasing offshore. Sediment distribution at the layer close to the bottom (lower panels in Figure 5) is similar to that of the surface layer at the upper and intermediate regions, but concentrations have higher values. The concentration peak close to Punta Piedras is stronger and extends over the Barra del Indio Shoal for both textures.
To better analyze the effect of waves, lower panels of Figures 6 and 7 show the differences between Cases 4 (waves) and III (no waves). When waves are included in the simulation, higher concentration of suspended sediment are observed throughout the RdP Estuary. This provides evidence of the importance of the effect of wind waves in keeping sediment suspended and/or lifting sediments to the surface. However, the increments in the concentration of both suspended silt and clay are especially important along the entire southern coast of the RdP Estuary (maximizing close to Punta Piedras and the northern part of Samborombón Bay) and along the northern coast, approximately offshore Juan Lacaze. Suspended sediment concentration is also higher over the Barra del Indio Shoal. At the northern coast of the upper and intermediate RdP Estuary, the difference in concentration between Cases 4 and 3 decreases considerably, showing that this region is substantially less influenced by waves. The main characteristics of the suspended sediments distribution reported by Moreira et al. (2013) and obtained from in situ and remote observations were similar to the results obtained by the Case 4 with waves ( Figure 3). In fact, it can be concluded that after tides, wind waves are the main driver of the dynamics of suspended sediment in the RdP Estuary.

Idealized Simulations With a Flat Bottom and Without Salinity
To complement the study and to contribute to a better understanding of the physical processes involved in the suspended sediment dynamics within the RdP Estuary, two additional idealized simulations were made:  , and a maximum between Punta Piedras and the northern portion of Samborombón Bay, are also observed in the FB simulation. This is particularly significant as it suggests that the zone of sediment deposition between Montevideo and Punta Piedras is independent from bathymetry and is mainly linked to the Journal of Advances in Modeling Earth Systems geometry of the RdP Estuary. In this region, the estuary widens significantly, leading to a reduction of the barotropic transport associated to continental discharge. In this simulation the bottom salinity front (not shown) locates following the isobath that connect the 8 meter flat bottom with the natural bathymetry.

Case "Without Salinity"
The second idealized simulation (Case WS) consisted in neglecting seawater salinity and, therefore, removing the salt wedge and the associated estuarine circulation. For this case, the initial salinity condition was set to 0 UPS throughout the domain. Bathymetry was realistic. Figure 8 shows the distribution of silt (left) and clay (right) concentrations hourly averaged over the last year of simulation for the surface (upper panel) and bottom (lower panel) model layers. It is interesting to compare this figure with Case 4 (Figure 4 lower panel) in which the only difference is the presence of seawater on the Continental Shelf and outer RdP Estuary. When ocean salinity is neglected the most important features of the suspended sediment distribution at the intermediate and upper parts of the RdP Estuary remain unchanged, as could be expected. At the outer RdP Estuary, however, in the absence of salinity both silt and clay are advected offshore Barra del Indio, particularly over its central part. The effect of gravitational or estuarine circulation, even when is weak (Simionato et al., 2007), is to maintain the particles of sediments upstream the shoal, In the Case WS, the sediments exported out the Barra del Indio Shoal deposit more slowly. The gradient along the estuary axis still occurs, although it is much weaker.
This simulation shows that the occurrence of the salt wedge and the associated estuarine currents play an important role on sediment dynamics in the RdP Estuary. The sediments that reach the Barra del Indio Shoal deposit there as a result of the combination of a number of factors: (i) the weakening of current speed due to the broadening and deepening of the estuary, (ii) the physicochemical flocculation processes occurring at the region of the bottom salinity front, and (iii) the estuarine circulation.

Conclusions
In this paper, we analyze the effect of diverse drivers (continental discharge, tides, winds, and winds waves), of the morphology and of the estuarine circulation associated to the adjacent ocean on the distribution of fine suspended sediments of the RdP Estuary by means of a set of process oriented numerical simulations. The applied model is the IFREMER's hydro-sedimentological MARS model, implemented for three types of sediments: fine sand, silt, and clay. This model was tuned and validated for this estuary by Fossati (2013), Fossati et al. (2014), and Moreira et al. (2016), and in the particular implementation made in this paper, it qualitatively well reproduces the main features of the distribution pattern of fine sediments in the RdP Estuary, so as the magnitude of the concentrations and gradients observed there. The SPM concentration obtained from the simulations presents the two gradients observed on in situ and remote data (i) one with northeast to southwest direction in the intermediate RdP Estuary and (ii) a second one oriented from southeast to northwest in the outer estuary. It also reproduces the observed relative maximum concentrations around Punta Piedras and Punta Rasa. In this sense, although there are natural differences between the observations and the idealized numerical simulations run in this paper, the model can be considered to perform satisfactorily and the process oriented simulations made in this work, oriented to (qualitatively more than quantitatively) understand the hydro-sedimentological processes that occur in this estuary, can be regarded as valid.
The results of the simulations of this work complement previous works based on observations (i.e., Moreira et al., 2013;Parker et al., 1986b;Parker & López Laborde, 1989) and contribute to building a more robust conceptual model of the processes that determine the distribution of suspended sediments in the RdP Estuary. In what follows, we summarize and discuss the conclusions drawn from the simulations.
In general, and as expected, in all the simulations the fine sediments (fine sand, silt, and clay) transported to the RdP Estuary by the main tributary rivers (Paraná and Uruguay) deposit from the head of the estuary to its mouth, according to their diameter and settling velocity. However, the dynamic processes occurring in this wide and shallow estuary significantly affect deposition. Thicker and with higher settling velocity, fine sand deposits immediately at the mouth of the tributary leading to the progress of the Delta of the Paraná River, in agreement with other authors (Parker et al., 1987). Also, in agreement with what postulated by several authors (Moreira et al., 2013;Ottman & Urien, 1966;Parker et al., 1986b), silt deposits mostly along the upper and intermediate RdP Estuary and presents low concentrations over the Barra del Indio Shoal. Clays, however, can remain in suspension over a longer distance, because they are smaller than the others textures, and have lower settling velocity. Therefore, clay is found throughout the estuary, even downstream the Barra del Indio and in Samborombón Bay. In summary, the hydrodynamic conditions in every area of the estuary affect the concentration of textures in the water column.
The results obtained from simulations incorporating one forcing at a time led to the following conclusions:

Solid and Liquid Discharge
The discharge of the three main tributaries (Uruguay, Paraná de las Palmas, and Paraná Guazú-Bravo) has its greatest impact at the upper and intermediate estuary. The differences in both the liquid discharge and the concentration of the different types of sediments between the tributaries set up the gradient oriented across the estuary axis in the upper and intermediate RdP Estuary. It results from the fact that the Paraná (and particularly the Paraná de las Palmas) River carries more sediment than the Uruguay River. This fact was postulated by Parker et al. (1986b) and Parker and López Laborde (1989).

Tides
The simulations show that tides play a key role enhancing lateral mixing and increasing kinetic energy over the entire RdP Estuary. Consequently, sediment deposition is inhibited, and sediments can be more rapidly transported to the area of the Barra del Indio. The erosion associated with tidal currents causes a general increment of the sediments concentration over all the RdP Estuary that, of course, maximizes in the areas where the dissipation of tidal energy by bottom friction is highest: around Punta Piedras and close to Punta Rasa. Results also suggest that the northeast-to-southwest oriented gradient of sediment concentration caused by the differences in the discharges of the tributaries, strengthens by the tidal mixing. This hypothesis was postulated for the first time by Moreira et al. (2013).

Winds
Winds help to maintain sediments in suspension and increase lateral mixing throughout the estuary. This forcing has its greatest impact on the Barra del Indio Shoal and to the south of Punta Piedras.

Waves
The effect of wind on the concentration of suspended sediments is amplified indirectly through wind waves. Bottom friction and vertical mixing due to waves increase the concentration of suspended sediments throughout the water column over a significant area of the estuary. Waves enhance lateral mixing, particularly in the shallower regions of the RdP Estuary, homogenizing the distribution of the sediments discharged by the tributaries. It also increases the concentration gradient perpendicular to the estuary axis on the uppermost and central parts of the estuary. When waves are included in the simulation, the largest increment in suspended sediment concentration occurs at the areas where tidal currents maximize and erode bottom sediments (particularly, along the southern coast of the intermediate and outer RdP Estuary). Those sediments are moved to the surface by wave-driven vertical mixing. This hypothesis was posed from observations by Moreira et al. (2013). Thus, both observations and simulations suggest that whereas the sediment distribution pattern at the RdP Estuary is geographically determined by tides, tidal energy is enough to re-suspend sediments only in the lower layers but waves raise them to the surface. Offshore Barra del Indio, both the wind and wind waves are still important, but given that there is little suspended sediment remaining, their greatest impact occurs on the coasts and shallow areas, especially at the southern coast and in Samborombón Bay.

Geometry and Estuarine Circulation
When the model is forced with runoff and tides it succeeds in reproducing an abrupt reduction of the concentration of suspended sediments in the region of Barra del Indio. Simulations results of course improve when the winds and waves are included, but the area where the longitudinal gradient occurs is always present after the imaginary line between Punta Piedras and Montevideo. There, several factors converge: (i) an abrupt increment of depth, (ii) an abrupt increment of the width of the estuary, and (iii) the occurrence of the bottom salinity front associated to the salt wedge. The first two factors contribute increasing the estuarine transversal section and, thus, reducing the strength of the mean transport. This, in turn, increases the sedimentation capacity. The third factor works containing the fresh water and sediment discharges of the tributaries at the area of the Barra del Indio Shoal and, presumably, being the reason of its location. In the region of the bottom salinity front, flocculation processes also occur due to the physicochemical associated processes and the large increment of the concentration.
The results of the simulations do not only aid on the construction of a conceptual model but also the fact of identifying the areas where the diverse processes become relatively more important could be used in the future for a better design of new field experiments. The so obtained observations, in turn, should be used to improve the model performance by fine-tuning its parameters. The capability of the model of reproducing the main features observed in the RdP Estuary even on simple and idealized simulations like those of this paper augurs that it will be possible to improve implementations with more and better information on the forcings and their variability. It is envisaged to continue this research focusing in that direction. In this sense, an important limitation in our simulations is the simple empirical wave model that we applied. Although numerical models could be applied to estimate waves characteristics from wind data, observations are still indispensable for model calibration and validation. This poses a serious limitation, as the scarcity of wave observations in the RdP Estuary has limited the description of their features particularly at its upper and intermediate parts. Therefore, efforts should be made to improve this aspect first.