Delf3D is a three-dimensional simulation software package developed by the Deltares Institute in the Netherlands. It is based on the Navier-Stokes equations (N-S equations) that describe hydrodynamics. It uses the finite difference method to solve the N-S equations, combined with the material balance equations that describe the state of sediments, to achieve sediment transport and sedimentary landform evolution. In the horizontal direction, Delf3D provides three types of coordinates: Cartesian coordinates (intuitively represent the distribution of objects); orthogonal coordinates (reduce the discrete error of some curved edge parts); while in the vertical direction, Delf3D mainly provides two differential vertical grid systems: α coordinate system and Cartesian coordinate system, both of which are used to characterize undulating terrain and moving free surfaces.

Horizontal coordinate system.

1) Cartesian coordinate system (rectangular coordinate system)

(x _i , y _i ) are the coordinates of any point. The distance formula d between two points, such as A (x ₁ , y ₁ ) and B (x ₂ , y ₂ ) is as described above.

2) Polar coordinate system (a two-dimensional orthogonal coordinate system)

Any point P in the plane is represented by the polar radius r and the polar angle θ , that is P (r, θ ), the conversion relationship with rectangular coordinates is shown in 2) above.

Vertical coordinate system: σ = z − ζ d + ζ = z − ζ H   ;   H = d + ζ

z is the vertical coordinate in physical space; ζ is the elevation of the free surface relative to the reference plane (z=0); d is the water depth under the reference plane (m); H is the total water depth (m),

The key components of Delft3D include: Hydrodynamics: The model calculates the flow field through simplified Navier-Stokes equations, assuming that vertical acceleration is neglected and the pressure distribution is hydrostatic. The model can handle water flows caused by tides, currents and wind. Sediment transport: Delft3D simulates the transport process of cohesive sediments (such as clay) and non-cohesive sediments (such as sand). It takes into account bed load and suspended load transport, and uses Van Rijn and other formulas to describe the erosion and deposition process of sediment. Landform evolution: The model calculates the change of riverbed elevation based on sediment transport and simulates erosion and deposition rates. These changes affect the hydrodynamics, forming a feedback loop between flow, sediment transport and geomorphic changes. Meshing: Delft3D uses structured curved or rectilinear grids to represent the physical domain. The accuracy of the results depends on the resolution of the grid. The finer the grid, the more accurate the simulation results, but the calculation time will also increase. Numerical solver: The model uses iterative solvers such as Newton’s method to solve nonlinear equations. The convergence and accuracy of the solution depend on the time step, grid resolution and solver settings. Error and stability: Numerical dispersion and truncation errors are important considerations that affect the reliability of the simulation. The stability of the simulation is ensured by appropriate time steps and CFL condition checks.

This study mainly uses the Flow module in Delf3D. The Flow module is mainly used to solve two-dimensional (depth average) and three-dimensional unsteady shallow water equations. Taking into account the effects of water flow and meteorology, the software adopts a curve grid discrete format with good boundary fitting and reasonable structure. The resolution of a single grid is set by the simulator. The curve grid can be used to calculate unsteady flow and material transport phenomena, and is particularly suitable for simulating the sediment deposition process of rivers. The simulation process follows the law of conservation of fluid momentum and the law of conservation of mass. The completed simulation mainly includes four steps: (a) Solve the Navier-Stokes equation (N-Sequations) to obtain the flow state of the water body in the simulation grid; (b) Calculate sediment transport, deposition and erosion; (c) Obtain relevant parameters such as sedimentary landforms based on the results of steps (a) and (b); (d) Save the results and update the simulation data. The Flow module realizes continuous dynamic simulation of the sediment deposition process through iterative calculations of steps (a) to (d). In order to accurately simulate the transport, deposition and erosion of muddy and sandy sediments in rivers and lakes, the Engelund-Hansen formula is used for sediment transport calculation. The main calculation formulas 1–8 is as follows:

1) Flow field calculation:

Wherein, x is the coordinate in the downstream direction; y is the coordinate in the direction across the water flow; z is the vertical coordinate; Z _w is the water level; u is the flow velocity in the x direction, m/s; v is the flow velocity in the y direction, m/s; w is the flow velocity in the z direction, m/s; h is the water depth; C is the Chezy roughness, m^1/2/s; g is the acceleration of gravity, m/s²; V _h is the horizontal eddy viscosity coefficient, m²/s; V _v is the vertical eddy viscosity coefficient, m²/s.

2) Suspended sediment transport and deposition calculations:

Where: c ^(l) is the concentration of sediment component (l), kg/m³; u, v, w are the flow velocities in the x, y, z directions, m/s; Ꜫs,x(l), Ꜫs,y(l), Ꜫs,z(l)) are the eddy diffusion rates of sediment component (l); w _s (l) is the settling velocity of the obstructed sediment component (l), m/s.

3) Non-cohesive material deposition and erosion equation:

Where: C l is the sediment concentration (l) (kg/m³); C ∂ l is the reference sediment concentration (l) (kg/m³); ∂ is the Van Rijn reference height; h is the water depth; z is the height from the bottom bed; A l is the Rouse number

4) Sediment settling velocity equation:

Where: C s o i t is the reference density (input); ω s , 0 l is the basic sedimentation ratio at a specific sedimentation velocity; C s t o t is the total sedimentation mass of different sediment components.

5) Depth average kinetic energy calculation:

Wherein, ζ is the height of the water surface; h is the water depth; u and v are the depth-averaged flow velocities in the downstream x direction and along the coast y direction, m/s; g is the gravitational acceleration, m/s²; t is the flow time, s; μ is the eddy viscosity; C is the Chèzy viscosity coefficient, m^1/2/s.

Poyang Lake, situated in Jiangxi Province (28°24′–29°46′N, 115°49′–116°46′E), extends 110 km from north to south and 50–70 km from east to west, narrowing to just 5–8 km in the northern region. Covering a total area of 3,210 km², it is characterized by relatively flat terrain and is recognized as the largest freshwater lake in China. The lake features well-developed delta plains, delta fronts, and prodeltas. The upper and lower delta plains are expansive and widely distributed, whereas the delta front is relatively narrow, with a clear distinction between the upper and lower delta plains (Figure 1) (Li et al., 1990; Jin et al., 2014).

The hydrological system of Poyang Lake is extensively developed and primarily fed by five major rivers: the Ganjiang, Xiuhe, Xinhe, Raohe, and the Jiangxi River. Among these, the Ganjiang River has the most complex and well-established water system. Hydrological data from the Ganjiang River, recorded by a station located in its upper reaches since the establishment of the People’s Republic of China, indicate an annual average flow ranging from 1,000 to 1,800 m³/s, a sediment concentration of approximately 0.15 kg/m³, and a flow velocity between one and 2.5 m/s (Wang, 2002). Additionally, Poyang Lake is hydraulically connected to the Yangtze River in its northern region. During the flood season, water from Poyang Lake flows into the Yangtze River, while in the dry season, there is a possibility of backflow from the Yangtze into the lake. Consequently, Poyang Lake is classified as a typical shallow lake with river connectivity (Jin et al., 2014).

The annual water level of Poyang Lake exhibits significant variability. Data from the Duchang hydrological station, which recorded the average monthly water levels from 1953 to 1992, reveal a seasonal difference of up to 7 m between the wet and dry seasons. During the summer wet season (May–October), the lake’s water level is notably higher, with specific values illustrated in Figure 2. Conversely, in the autumn and winter dry season (December–March), the water level decreases substantially, reducing the lake’s surface area to less than 1,000 km². During this period, the delta sedimentary area advances northward and retreats southward by approximately 30–50 km annually (Liu et al., 2016; Zhang et al., 2016; Gao et al., 2016).

The boundary conditions of the simulation are crucial to ensure the accuracy of the calculation and the reliability of the results. In this study, the setting of the experimental work area scale is mainly based on the modern Poyang Lake sedimentary area in China. There are many typical shallow water delta deposits in this area. The study takes a delta in the middle reaches of the Ganjiang River as an example (Blue frame line), and its sedimentary area is about 20 km * 20 km by measuring (Figure 1 blue area). Based on this, the size of the simulated work area in this study is set to 21 km long and 17.5 km wide (Figure 3). This configure provides sufficient space for the dynamic evolution of the delta under different lake levels. At the same time, the initial sedimentary bed slope was set to 0.07° (Table 1). In order to accurately capture the geometry of the delta, the grid size used in the simulation was defined as 70 m × 70 m, which minimized the error of the delta geometry and ensured that each diversion channel occupied at least one to three grid cells. Previous studies have shown that grid sizes ranging from 25 m to 100 m produce the best simulation results.

Regarding boundary conditions, the left and right boundaries of the simulation domain were designated as closed boundaries, while the upper boundary was set as an open boundary. Additionally, a water supply port was included at the bottom of the domain to maintain a constant total flow.

For the sediment dynamics parameters, settings were based on sediment dynamics data from the Ganjiang Delta of Poyang Lake and similar simulation studies (Edmonds and Slingerland, 2010; Caldwell and Edmonds, 2014; Xu et al., 2017). The simulation included both sandy (non-cohesive) sediments and muddy (cohesive) sediments. The sandy sediments had median particle sizes of 300 mm, 150 mm, and 80 mm. The sand-to-mud ratio was set at 3:2, with a river supply rate of 1,500 m³/s and a flow velocity range of 0.8–2.0 m/s, consistent with the hydrological and sediment data of the Ganjiang Delta (Table 2).

To compare and analyze the scale characteristics and sedimentation processes of shallow deltas under varying lake levels, this study conducted three simulation experiments, designated as A1, A2, and A3. These simulations were based on real lake level data from Poyang Lake, with lake levels set at minimum, average, and maximum values. Specifically, the lowest lake level of A1 was set to −1 m to simulate delta sedimentation under low lake conditions. The highest lake level of A3 was set to 6 m to simulate delta evolution under high lake conditions. Experiment A2 represented a transitional lake level, reflecting average conditions to examine the growth and development of shallow deltas under intermediate lake levels. In all three experiments, only the lake level values were varied, while the remaining simulation parameters were kept constant.

To ensure the accuracy of the simulation experiments, the time step was set to 0.4 min. To expedite the simulation process, previous research has demonstrated that incorporating a morphological scale factor can be effective. Specifically, when the morphological scale factor is less than 200, it does not adversely impact the simulation results (Caldwell and Edmonds, 2014). Consequently, this study set the morphological scale factor to 150, resulting in a 150-fold acceleration of the simulation. Additional parameters are detailed in Table 2.

Delta Length (L): This parameter represents the maximum distance the delta extends along the coastline, as described by Caldwell and Edmonds (2014) and Xu et al. (2021). It indicates the delta’s capacity to prograde into the basin. A specific schematic diagram illustrating delta length is shown in Figure 4A.

Delta Width (W): Defined as the maximum span of the delta perpendicular to its length, this metric reflects the delta’s spreading capability (Caldwell and Edmonds, 2014). Figure 4A provides a schematic representation of delta width.

Delta Area (A): This parameter is determined by accumulating and summing the two-dimensional plane data of sediments over the simulation period. By binarizing the simulation results based on the two-dimensional water depth data and corresponding lake surface values, the area of the shallow water delta can be accurately measured. A plane schematic diagram of the delta area is presented in Figure 4A.

Delta Shoreline (SH): The shoreline length, representing the outer boundary of the delta front, serves as a key parameter for characterizing delta shoreline morphology. A longer shoreline suggests a more complex delta morphology with increased branching and frequent flow velocity changes, whereas a shorter shoreline indicates a simpler and more stable delta form. The specific schematic diagram is shown in Figure 4A.

Delta Thickness (SE): The sediment thickness distribution reflects the topographic variations of the delta. This study extracts two-dimensional sediment thickness data from different time series and analyzes it using MATLAB software. By calculating the standard deviation of sediment thickness through custom code, the variability in terrain thickness is quantified. A larger standard deviation indicates greater topographic variation. The plane and 3D schematic diagrams of delta thickness are shown in Figure 4C.

Depth-averaged velocities: in Delft3D numerical simulation, the two-dimensional depth-averaged velocities of each time series on a certain grid is saved. This study counted the depth-averaged velocities of each grid in the work area under all simulation time series conditions and averaged them to obtain a two-dimensional plane distribution map of the depth-averaged velocities in the simulation area.

Delta Land Area: In this study, the delta sedimentary area is defined as the land portion of the delta where the water depth is less than 0.5 m. Figure 4A presents a schematic diagram of the delta area, from which the land area of the delta is measured.

Distributary Channel: The classification of distributary channels is based on flow velocity. Two-dimensional plane diagrams of flow velocity for various time series in the simulated sections are analyzed. Channels with an average flow rate greater than 1 m/s are classified as active distributary channels. These channels are characterized by higher flow velocities, indicating a longer survival time, a greater sediment and water transport capacity, and a significant influence on delta growth and development. Channels with flow velocities ranging from 0 m/s to 1 m/s are classified as breach distributary channels, indicating shorter survival times and quick abandonment after formation. Active distributary channels are depicted in red, while breach distributary channels are shown in blue in Figure 4B.

As illustrated in Figure 5, the scale of the delta exhibits considerable variation across different lake levels. Within the same simulation period, the length, width, and area of the delta are inversely related to the lake level. At the maximum lake level of 6 m and a simulation duration of 100 h, the delta attains its smallest dimensions, with a length of 2.09 km, a width of 1.79 km, and an area of 1.28 km², assuming a flower-like shape. At an intermediate lake level of 2.5 m, the delta’s dimensions increase slightly, reaching a length of 2.26 km, a width of 4.41 km, and an area of 6.99 km². Conversely, at the minimum lake level of −1 m, the delta reaches its largest scale, with a length of 4.43 km, a width of 12.56 km, and an area of 30.87 km².

As simulation time progresses, the scale of deltas under varying lake levels generally increases, although deltas at higher lake levels remain the smallest. In comparison to other lake levels, deltas under lower lake levels exhibit a pronounced growth. After 500 h of simulation, the rate of increase in length, width, and area of deltas under high lake levels is slower compared to those under medium and low lake levels. As illustrated in Figure 6A, the growth rate of delta length, width, and area at high lake levels is the least pronounced, with dimensions of 4.12 km in length, 7.58 km in width, and an area of 22.82 km². Conversely, deltas at medium and low lake levels grow at a faster rate. The disparity in scale between deltas at medium and low lake levels and those at high lake levels becomes more pronounced over time. At the 500-h mark, deltas under medium lake levels measure 5.46 km in length, 11.76 km in width, and cover an area of 45.12 km², while deltas under low lake levels reach 8.49 km in length, 18.39 km in width, and encompass an area of 100.97 km².

As the delta scale continues to expand, the growth rate of the delta’s dimensions gradually diminishes, leading to the maximum disparity in scale between deltas at different lake levels. After 800 h of simulation, the delta at a low lake level achieved dimensions of 10.64 km in length, 21.03 km in width, and 139.79 km² in area. These measurements are approximately twice the length and width and four times the area of the delta observed at high lake levels, highlighting a substantial scale difference.

Further analysis of the delta shoreline values provides additional insights into the morphological changes under varying lake levels over time. The study reveals significant differences in delta morphology across lake levels. The shoreline length of the delta exhibits a gradual increase over time, with the growth rate being slower at high lake levels compared to low lake levels. Shoreline measurements taken at 100, 500, and 800 h indicate that at low lake levels, the delta shoreline extends to 17.22 km, 28.75 km, and 36.02 km, respectively. In contrast, at high lake levels, the shoreline values are substantially lower, measuring 5.80 km, 16.33 km, and 18.69 km. This data underscores that lower lake levels significantly promote the expansion of the delta shoreline.

The MATLAB code developed for this study measures the thickness of the delta at various lake levels. As illustrated in Figure 7, there is considerable variation in delta thickness across different lake levels. Specifically, the average thickness of the delta increases with rising lake levels. Correspondingly, the standard deviation of delta thickness also increases, indicating that higher lake levels lead to more pronounced topographic variations and a more rugged surface terrain.

At a low lake level (lake level = -1m), the distribution of delta thickness is relatively uniform, with sediment thickness predominantly ranging between 10 and 15 m. The most common sediment thickness is approximately 13 m, which constitutes the largest proportion at 5.7%. As shown in Table3, at this low lake level, the average delta thickness is the smallest, approximately 12.66 m, and the maximum thickness is also relatively modest, around 16.73 m. The delta exhibits minimal topographic variation at this level, as evidenced by the smallest standard deviation of thickness, which is approximately 1.22 m.

At a high lake level (lake level = 6m), the delta thickness predominantly ranges between 15 and 20 m. The most common sediment thickness is approximately 18.5 m, accounting for 25% of the measurements. As indicated in Table3, at this high lake level, the average delta thickness is the greatest, around 18.03 m, with a maximum thickness reaching up to 21.03 m. However, the distribution of delta thickness is uneven, as evidenced by a relatively high standard deviation of approximately 2.5 m.

At a medium lake level (lake level = 2.5 m), the delta thickness primarily falls between 14 and 17 m, with a maximum thickness of 18.03 m and an average thickness of 14.84 m. This average thickness is intermediate between the values observed at low and high lake levels. The standard deviation of delta thickness at this level is 2.03 m, which is slightly lower than the standard deviation observed at high lake levels.

As illustrated in Figure 8, the characteristics of the delta distributary channels vary significantly across different lake levels. Specifically, the number of distributary channels is inversely proportional to the lake level height. At low lake levels, the number of distributary channels is the highest, followed by medium lake levels, with the fewest distributary channels observed at high lake levels. Conversely, the number of active distributary channels remains relatively constant regardless of the lake level height. This observation suggests that variations in lake level have a more pronounced effect on breach distributary channels than on active distributary channels.

Similarly, the distribution range of distributary channels is inversely proportional to lake level height. This study uses the distribution angle of active distributary channels as a parameter to quantify the distribution range. At low lake levels, the development angle of active distributary channels is approximately 178.8°, covering nearly the entire lake basin. At medium lake levels, the development angle decreases to 144.31°, with reduced development observed on the right side of the lake basin. At high lake levels, the distribution angle of active distributary channels reaches a minimum of 125.81°, with active distributary channels primarily concentrated near the river mouth. Consequently, higher lake levels correspond to a narrower distribution range of distributary channels and a more restricted riverine function.

This study also assessed the length of active distributary channels across different lake levels. As depicted in Figure 9, the length of active distributary channels is inversely related to lake level height. At low lake levels, the average length of active distributary channels is approximately 10.33 km, with channels remaining distinct and separated by significant distances. At medium lake levels, the average length decreases to 7.02 km, with multiple active distributary channels becoming more compact. At high lake levels, the average length of active distributary channels is reduced to 4.52 km, with channels closely distributed in the direction of the river source.

Under different lake levels, the growth and development process of the delta determines the differences in its scale. Numerical simulation reproduces the complete development process of the river-controlled delta from early to middle to late. This time, we will take low and high lake levels as examples to compare and analyze the differences in its evolution process.

Previous researchers believed that the supply of river materials was the key factor affecting the scale of river-controlled deltas. Previous studies were based on bird-foot deltas, such as the modern Mississippi River Delta and the Yellow River Delta, whose sediment ratios were 11.2 and 2.8, respectively; while typical tongue-shaped deltas, such as the Wax Lake Delta and the Mossy Delta in the United States, had sediment ratios less than 1 (Edmonds and Slingerland, 2010). This confirmed that the supply of fine-grained, highly viscous sediments can promote the accumulation of natural levees, enhance the cohesion of natural levees, and be conducive to the formation of bird-foot deltas. Such deltas have a smaller geometric scale and a larger sediment thickness; while coarse-grained, low-viscosity sediments lead to the rapid deposition of estuary bars and the diversion of distributary channels, promoting the formation of tongue-shaped deltas with distributary sand bars as the main source. Mainly, this type of delta has a larger geometric scale and a relatively smaller sediment thickness than that of a bird-foot delta (Edmonds and Slingerland, 2007; Caldwell and Edmonds, 2014; Burpee et al., 2015; Orton and Reading, 1993).

In addition to supplying river materials, river discharge is also an important factor affecting the formation of delta scale. Based on numerical simulation of sedimentation, Xu et al. (2021) found that under high river discharge, river-controlled deltas are prone to develop diversion and breach and form large and thin tongue-shaped deltas; under low discharge, the diversion channel in the river-controlled delta is relatively stable, and it is easy to form a small and thick bird-foot-shaped delta, and pointed out that the discharge value of 1,000 m³/s is a key threshold.

However, based on the numerical simulation method, this paper found that the height of the lake level also has a great influence on the scale of the delta by comparing three different lake level values: high, medium and low. When the lake level is high, the delta is bird-foot-shaped, with a small plane scale but a large sediment thickness (Figures 5G, H, L; Figures 6A–D); when the water level is low, the delta is tongue-shaped, with a large plane scale but a thin sediment thickness; when the water level is between the two, the scale of the delta is also between the two.

The size of the delta reflects the strength of the river’s energy supply. The larger the delta, the stronger the river’s energy supply (Galloway, 1975; Olariu, 2014). When the lake level is low, the delta is large, indicating that the lake basin has low energy, and the corresponding supply river has relatively high energy. When the lake level is low, the limited accommodation space and the smaller A/S ratio lead to the supply river being mainly erosive, and the sediments can be transported farther. As the supply river’s ability to carry sediments gradually weakens with the increase of the extension distance, the coarser sand bodies are deposited first because they are less affected by the support of the top of the lake surface; the finer sand bodies are more affected by the support of the top of the lake surface, and expand outward to form a larger delta plane. On the contrary, when the lake level is high, the accommodation space is larger, resulting in a larger A/S ratio, and the downcutting ability of the supply river is weakened, mainly vertical accumulation, resulting in a thicker delta and a smaller scale at high lake levels. The study characterized the water velocity data corresponding to the delta plane at different lake level values ​​under a simulation time of 506h (Figures 14A–C). As can be seen from Figure 14, at high lake level water levels, the water flow rate is slow and the flow range is small. Most sediments are deposited quickly and it is difficult to migrate further. At low water levels, there is flow to varying degrees in all parts of the delta front (i.e., the underwater delta part), which leads to the continuous transportation of sediments and their rapid diffusion to the surrounding areas, resulting in a large and rapid increase in the size of the delta.

Consistent with delta size, the number of distributary channels is also affected by multiple factors, including the flow and sediment characteristics of the source river and the depth of the lake basin. Edmonds et al. (2009) demonstrated through numerical simulation of sediment that changes in river flow affect the stability of shallow deltas. Specifically, a decrease in river flow leads to a proportion of existing distributary channels being abandoned, resulting in a decrease in their number. Conversely, a 60% increase in river flow may lead to the formation of new distributary channels. In addition, sediments with high viscosity, low sand-mud ratio and fine grain texture generally lead to fewer distributary channels, while sediments with low viscosity, high sand-mud ratio and coarse grains are associated with more distributary channels (Burpee et al., 2015).

This study also found that different water level values ​​also have a great influence on the number of distributary channels in the delta, but the impact of water level on distributary channels is mainly concentrated in breach distributary channels, while it has little effect on active distributary channels (Figure 8). At low lake levels, the supply capacity of the river is relatively strong, resulting in the formation and dispersion of distributary channels. This condition increases the likelihood of breach distributary channels. In contrast, during high lake levels, the lake supports the feeding rivers, concentrating the flow at the front of the feeding channels. As a result, these channels tend to migrate laterally, which reduces breaching and the likelihood of breach distributary channels.