## Abstract

Coupling models of different dimensions is one of the most important yet under-represented challenges. This paper introduces a new modeling strategy to streamline a more flexible and effective integrated one-dimensional (1D)/two-dimensional (2D) model for floodplains along lowland rivers. The 1D model, utilizing the finite volume method, solves the Saint–Venant equations, while the 2D mesh employs unstructured quadrilateral elements. The two strategies couple the 1D/2D models: direct 1D/2D connection by the law of mass conservation at supernode, and lateral 1D/2D model connection by spillways at riverbank. The coupling strategy in F28 guarantees the water balance and the conservation of momentum at the integrated 1D/2D nodes. The model was applied to the Mekong Delta to address the capacity of hydrodynamic simulations integrating various water infrastructures. Results showed that the developed model has a strong potential to be applied to other lowland rivers worldwide with complex infrastructures.

## HIGHLIGHTS

Innovative 1D–2D model coupling enhances floodplain flexibility.

Integrated nodes ensure water balance and momentum conservation.

Unstructured 2D mesh improves river discharge simulation.

Model employs explicit finite volume, semi-analytical coupling, and parallel approach for stability and precision.

## INTRODUCTION

Numerical models have long been utilized to investigate various water-related processes, including hydrology and hydrodynamics, to test hypotheses, and improve decision-making on water resources, water quality, and related hazard management of governance under global climate change and economic development (Fang *et al.* 2016; Shah *et al.* 2018). For over a century until the middle of the 1960s, mathematical modeling has been primarily employed for the development of concepts, theories, and models of individual components of the hydrologic cycle, such as overland flow, channel flow, infiltration, depression storage, evaporation, interception, subsurface flow, and base flow (Green & Ampt 1911; Horton 1935; Lighthill & Whitham 1955; Singh & Woolhiser 2002). However, due to the lack of data on physical basin characteristics and the limitations of computing power, conceptual and often lumped-parameter models mainly focus on the simulation of aggregated output (total streamflow) at the watershed outlet. Hence, these models are unable to provide sufficient details of dynamic hydrology processes and accurate estimation on a wide range of temporal and spatially heterogeneous domains (Clarke 1973; Chen *et al.* 2007; Guo *et al.* 2021). With the computer revolution in the 1960s, numerical hydrology has dramatically improved simulation duration and presented a more visually dynamic hydrology process, by integrating different component models of the entire hydrologic cycle in various scenarios (Crawford & Linsley 1966; Feldman 1981; Arnold *et al.* 1998; DHI Water Environment 2003). Therefore, hydrology modeling has been widely used as a tool to help define problems and identify some possible solutions for stakeholders (Anees *et al.* 2016; Singh *et al.* 2021).

Various numerical models for flood investigations have been developed based on input variables, accuracy of water depth in the specific temporal and spatial scales, and efficiency computation duration (Asselman *et al.* 2009; Teng *et al.* 2017). For example, 0D models have adapted simple rules to simulate water flow over floodplains in short run times based on limited data (Pender 2006; Kauffeldt *et al.* 2016). As for one-dimensional (1D) cell-based flows, they have been highly efficient tools to simulate water depth and average velocity at each cross-section in sewer systems or channels (Krysanova *et al.* 1999; Yoshida & Dittrich 2002; Papanicolaou *et al.* 2004; Hicks & Peacock 2005; Martínez-Aranda *et al.* 2019). The 1D channel flow and simulation of inundation are integral to the solution of the St Venant equations. Nevertheless, the technique is disadvantageous in that it assumes floodplain flow to be in one direction parallel to the main channel, while completely neglecting the momentum transfer (Pender 2006; Kauffeldt *et al.* 2016; Morales-Hernández *et al.* 2016). As a result, 1D model tends to produce inaccurate results in meandering rivers (Sellin 1964; Bousmar *et al.* 2004; Betsholtz & Nordlöf 2017). In contrast, the floodplain modeling is carried out via two-dimensional (2D) models, which apply shallow water equations to calculate water depths and velocities across the finer grid or mesh defined from topographic information (Néelz 2009; Shustikova *et al.* 2019). However, 2D models can be time-consuming, and large variations in cell sizes may occur because of the complex basin underlying surface or local refinement (Fan *et al.* 2017; Martínez-Aranda *et al.* 2019; Echeverribar *et al.* 2020). This leads to lower efficiency because the time step is limited by numerical stability and is determined mainly by the smallest cells in the mesh (Mignot *et al.* 2006; Fan *et al.* 2017). Besides, flood risk assessments in floodplain areas rely on the accuracy of supercritical flow representation that can be offered by a promising numerical model that simulates fluid dynamics (3D model). Still, flood studies using 3D models are limited by the need for experimental laboratory data and high computational cost (Li *et al.* 2006; Teng *et al.* 2017).

Under rapid urbanization and economic growth, developing countries have been changing the configuration of waterwork infrastructures along rivers (Playán *et al.* 2018). However, incompatible infrastructure has led to the degradation of stream ecosystems by triggering changes in the long-term floodplain system (Sholtes *et al.* 2018; Tran *et al.* 2021). Although 2D models are robust in simulating floodplains, complex water infrastructure necessitates intensive input data, which impacts computational time, and produces inaccurate water levels at meandering rivers. Consequently, researchers have proposed combining a 1D numerical schematization for the main channel, and floodplain 2D mesh for overbank flow (Werner 2004; Finaud-Guyot *et al.* 2011b; Bladé *et al.* 2012; Timbadiya *et al.* 2015). In 1975, the first integrated 1D/2D models were developed for the Mekong River Delta by solving the Saint–Venant equations with the Preissmann Scheme for a 1D model of looped channel flow and integrating, with a storage cell algorithm using the mass conservation equation to link domains (Cunge 1975). Subsequently, the method was adopted in the first versions of Mike-11 and the present version of Hec-Ras (DHI 2009; Bladé *et al.* 2012; USACE-HEC 2013). Nonetheless, this integrated scheme is unable to include the model's front wave advance and recession over the floodplain (Tayefi *et al.* 2007; Bladé *et al.* 2012). Several other approaches are linked to the model's separate 1D model and 2D domain via water level compatibility. However, the conveyance's water level in the main river channel in dynamic time cannot be accounted for during the model runs. Thus, the predicted overbank flows onto the floodplain will be erroneous (Kharat 2009; Deltares 2017). Besides, the most challenging integrated 1D/2D is the momentum transfer from the highly meandering main channel section to the floodplains (Bousmar *et al.* 2004; Aricò *et al.* 2016; Betsholtz & Nordlöf 2017). Therefore, the overbank flow of the main channel into the floodplain could overestimate the discharge capacity of the main channel (Kharat 2009; Finaud-Guyot *et al.* 2011a).

This paper introduces a new modeling approach to address the more flexible and effective integrated 1D/2D model for floodplains. The 1D model solves both Saint–Venant equations and shallow water equations using the finite volume method. At the same time, the 2D mesh is unstructured with quadrilateral elements. The 1D/2D models are coupled by two strategies: the balancing of the water level at supernodes; and the lateral exchange of fluxes between 2D nodes to 1D sections via the spillway. The coupling strategy guaranteed the integrated notes' water balance and conservation of momentum. The significant new model approach solves the interaction boundary 1D/2D fully transmitted, provides more reliable dynamic flux exchange between 1D and 2D models, and reduces the coupling time process. We successfully demonstrate the developed F28 in the heavily engineered floodplains of the Vietnamese Mekong Delta (VMD). Since flood discharge in the lowland is influenced by the overflow between roads, embankments, and mainstream, the coupling 1D and 2D of F28 could provide reduced computation time and more accurate results. The advantage of the new coupling model will boost a better early warning system, particularly for complex infrastructure floodplains.

## GOVERNING EQUATIONS

### Governing equations for 1D F28 hydraulic model

*η*is the water level (m);

*Q*,

*A*, and

*K*are the flow rate (m

^{3}/s), the cross-section area (m

^{2}), and the flow rate module for 1D flow (m

^{3}/s);

*s*is the space coordinate along the channel axis;

*q*

_{l}and

*u*

_{l}– the lateral inflow along the river (m

^{3}/s/m) and its axial component of velocity (m/s/m). At the tributary nodes (Figure 1), the volume conservation equation is expressed as follows:

where *W _{J}* is the volume of the

*J*th river node (m

^{3});

*Q*is the flow out from the tributary node into the

_{i}*i*th control cross-section (m

^{3}/s);

*L*

_{i}is the length of the river segment (m), divided from the tributary node to the

*i*th cross-section.

### Governing equations for 2D F28 hydraulic model

*D*is the water depth (m); ∇ is the differential operator; is the flux vector of flow rate; and is the vector of external forces. The flux vector and external forces vector have the following forms:where

*f*is the Coriolis parameter; (

*τ*,

_{wx}*τ*) is the shear stress on the water surface due to the wind (N/m

_{wy}^{2}); (

*τ*,

_{bx}*τ*) is the bottom shear stress (N/m

_{by}^{2});

*A*

_{H}is the eddy viscosity (m

^{2}/s);

*q*and

_{v}*u*,

_{a}*v*are the lateral inflow per unit surface area (m

_{a}^{3}/s/m

^{2}) and its velocity components.

*n*is the Manning's roughness coefficient; is the bottom shear stress velocity (m/s);

*g*is the gravity acceleration (m

^{2}/s).

At boundaries connecting with other flow models, appropriate connection conditions will be applied.

## NUMERICAL SOLUTION PROCEDURE

### Computational mesh

### 1D unsteady flow

*Δs*

_{j}, and the segment between the river's cross-sectional points is represented by

*δs*in Figure 3.

_{j}*A*=

*f*(

*η*). From Equation (3), the volume of the junction node is also calculated at timestep

*n*+ 1 (Figure 1). Here, we can also determine the water level (m) at the node from the function

*W*=

*f*(

*η*). The water level at the junction node is assigned for the nearby cross-sections:

### 2D unsteady flow

*n*+ 1/2 is computed. Equation (7) is integrated over the control area surrounding the edge, for instance, edge ict (Figure 4(b)). Upon applying the Green's integral transformation formula, we obtain:where

*S*and

*L*are control area (m

^{2}) and control perimeter (m);

*q*

_{n}is the component in the normal direction of the unit flow rate on the control perimeter.

*u*

_{s}and

*u*

_{b}are the two components along the

*s*and

*b*directions, respectively, of the inflow velocity. The convective term is interpolated using the Upwind scheme, while the friction term is estimated using the difference formula.

*n*+ 1, Equation (6) is integrated over the control area surrounding the node, such as node C in Figure 4(a). After transforming the integral using the Green's formula, we obtain:

*n*+ 1/2, we get the expression to calculate the node water level at the timestep

*n*+ 1:where

*L*is side length

_{j}*j*of the control perimeter;

*q*is the component in the normal direction of the unit flow rate on

_{nj}*l*, which is calculated as follows:

_{j}### Coupling 1D/2D strategy in F28 hydraulic model

The integration of the 1D/2D coupling establishes a continuous flow linkage between the two submodels, 1D and 2D, guided by the principles of conservation. These principles encompass the preservation of mass and momentum, articulated, respectively, through the continuity and momentum equations. When amalgamating submodels into a cohesive framework, the inter-submodel flow must adhere to these conservation principles. This entails a seamless transition of water and momentum between the submodels, ensuring that water leaving one submodel enters the other without loss, and momentum from one submodel is entirely transferred to the other (if neglecting energy losses). Notably, the inclusion of momentum transfer has been validated to enhance simulation precision, as corroborated by the findings of Finaud-Guyot *et al.* (2011a). Within the F28 model, two coupling types are offered: through the utilization of supernodes and spillways.

CalQ1: This subroutine is responsible for computing the 1D discharge.

CalQ2: This dedicated subroutine is designed to calculate the 2D flow rate.

CalQlk: This subroutine focuses on determining the discharge at riverbanks.

CalE1: This specialized subroutine carries out two primary tasks. Firstly, it calculates the volume of 1D nodes using Equation (23). Secondly, it computes the area of 1D cross-sections as defined by Equation (22). This procedure aids in deriving water levels at nodes and cross-sections within the 1D model, excluding those associated with the supernode coupling mechanism.

CalE2: This unique subroutine is responsible for estimating water levels for 2D nodes based on Equation (19), while excluding nodes engaged in the supernode coupling mechanism.

CalElkn: This subroutine is specifically designed to calculate the volume of water levels at linkage supernodes.

CalQlkn: This subroutine is tailored to interpolate the discharge at the link interface (Figure 5).

### Direct 1D/2D connection by the law of mass conservation at supernode

*W*is the volume of supernode

_{J}*J*;

*B*is the width of 2D interface;

_{2}*q*is normal to the 2D interface component of flow rate per unit width;

_{n}*Q*

_{1D}is the discharge at 1D interface; and

*Q*

_{BC}is the boundary discharge at node

*J*.

*Q*

_{int}, is estimated through interpolation using values from both the 1D and 2D interfaces. The apportionment of the 1D discharge to the 2D flow rate per unit width at the link interface is executed with the underlying assumption of a consistent hydraulic gradient

*J*across all points along the coupling sections:where

*B*represents the width of the link interface (m), and

*D*indicates the depth.

*Q*

_{int}and are used to evaluate the term in Equation (2) and in Equation (7), both pertaining to the 1D section and 2D edges located at the interfaces of the super node. This approach facilitates the transfer of momentum from one submodel to the other.

### Lateral 1D/2D model connection by spillways at riverbank

*q*

_{t}) is calculated using the spillway formula, and its unit is m

^{3}/s/m. When the water level at the 2D node is higher than the one at the 1D section and by neglecting the velocity head, the discharge is:where:

*η*

_{2D}and

*η*

_{1D}are the water levels at 2D node and 1D section;

*z*

_{ct}is the crest elevation;

*m*and

*φ*are the discharge and velocity coefficients. From calculated

_{n}*q*, mass transfer from 2D note to 1D section is evaluated:and momentum transfer from 2D edge to 1D section:where MAE is the mass transfer from 2D note to 1D section (m

_{t}^{3}/s) and MOE is momentum transfer from 2D edge to 1D section (m

^{4}/s

^{2});

*L*is the length of link interface between 2D edge and 1D section (m);

*L*′ is the length of link interface between 2D node and 1D cross-section (m);

*q*′

*is the discharge per unit width evaluated at the middle of link interface (m*

_{t}^{3}/s/m);

*V*

_{1}stands for the 1D velocity (m/s); and

*u*

_{s}is the tangent component of the 2D velocity at link edge. Considering these exchanges, flowrate at 1D section and discharge per unit width at 2D edge adjacent to the link interface now are calculated instead of Equations (17) and (31):and the areas at 1D cross-section and at 2D node located on link interface are calculated, instead of Equations (22) and (34):

### Model results analysis and comparison

*n*is the total number of observations,

*y*is the observed

_{i}*i*th value, is the respective simulated

*i*th value, and is the mean value of the

*n*observations.

## THE MODEL APPLICATION

### Research domain

The Mekong Delta is one of the world's largest and most intensively used estuaries, and the Mekong Delta provides ecological and food security for its inhabitants (Hung *et al.* 2012; Li *et al.* 2017). The most significant part of the Mekong Delta is located in southern Vietnam, where the Mekong River drains into the South China Sea. The region encompasses an area of 39,000 km^{2} and is the homeland of 17 million people. The VMD region provides 40% of the nation's agricultural production and accounts for more than 50% of agro-food exports (Bolay *et al.* 2019). The climate is influenced by two monsoon systems, the southwest Indian monsoon, and the northwest Pacific monsoon, causing two distinct seasons; the dry season from December to the end of April and the flood season from May to November (Holmes *et al.* 2009). The flow is confined to the rivers and channels in the dry season. However, in flood season, a large delta area is inundated, and an essential part of the flow in the floodplain deals with complex irrigation structures (Hung *et al.* 2012). The Mekong river flow benefits millions of habitants but causes some inconveniences, especially in the delta.

Flow study in this delta is needed and has attracted scientists' attention over the years. The most popular way to calculate the flow in the Mekong delta is to use 1D models with lumped storage and artificial branches (Nguyen, 1986, 2004; To *et al*. 2008; Van *et al*. 2012; Dung *et al*. 2013). Although well-calibrated 1D models can provide many valuable outputs, they are inadequate for flood inundation warnings. Since they cannot represent actual flow features on floodplains and hardly simulate the effects of various infrastructures such as roads, dikes, and ditches, some authors have used 2D models instead (Zanobetti & Lorgeré 1968; Hai *et al*. 2008). However, using a 2D model for the whole Mekong delta is impossible as the channel system in Vietnam is too complex to be reproduced by a 2D mesh. An integrated 1D/2D model is the best choice in this case. In this model, flow in the river is considered 1D, flow in the floodplain is 2D, and they are connected at river banks. The model of Dutta *et al*. (2007) is one of this kind. In our opinion, the limitation of expanding this model to Vietnam is using 2D rectangular mesh.

### Model setup

### Numerical results

^{3}/s in calibration and verification, which means the error is less than 5% of the observed river discharge. Besides, the NSE at all stations is relatively good with a value of 0.8. Although the downstream Mekong river in Vietnam is characterized by complex irrigation systems that influence the discharge of the Mekong river (Park

*et al.*2020), F28 could simulate flow in the floodplain stabilized from upstream to downstream Vietnamese Mekong river. A comparison of the results in Can Tho station (Figure 11) shows that the integrated two submodels' simulated flood flow and floodplain water depth agree with the observation. Figure 12 indicates that the coupled model presented reliable flood discharge in the complex irrigation basin. The simulated flows are slightly higher than the observations during the dry period. The groundwater flux exchange has not been accounted for in the 1D model, which may cause the discrepancy in the dry season discharge of the 2D mesh. However, the model results are considered certain overall.

*R*

^{2}values at all observation stations are notably high, approaching unity. A more detailed comparison at the Can Tho station (Figure 13) underscores the alignment between the simulated floodplain water depths generated by the integrated submodels and the observed data. This outcome underscores the model's reliability in forecasting flood water levels within the intricate irrigation basin.

## DISCUSSION

Seasonal flooding in VMD not only brings enormous benefits to local inhabitants and the regional economy but also constitutes a potential safety risk (Hoa *et al.* 2008; MONRE 2013; Dang *et al.* 2021). Over recent decades, the VMD has witnessed extensive development of human-made water control infrastructure, especially high-dike systems in the upper part of the delta, to fully protect agricultural productions from floods under climate change and enable intensive agriculture production (three rice crops per year) (Hung *et al.* 2012; Park *et al.* 2020, 2022). Besides, the riverbed incision rate is likely to increase due to the accelerating excavation activities to meet the growing demands for construction sands from both domestic and international markets (Park *et al.* 2020; Thi Kim *et al.* 2020; Van Binh *et al.* 2020). Consequently, the flow regime of the VMD has been altered and has caused more flood hazards in the region (Dang *et al.* 2016; Duc Tran *et al.* 2018; Xuan *et al.* 2022). Understanding delta flood dynamics is an essential yet highly challenging task to ensure the region's safety and sustainable flood management (Hoang *et al.* 2018; Triet *et al.* 2020). To fill this gap, the delta's flow regime needs to be accessed by a reliable flood model that integrates local water infrastructures. The F28 proposes a new approach to address the more flexible and effective integrated 1D/2D model for lowland basins, which can improve the accuracy of flood results and reduce simulation time.

The F28 model introduces a novel approach to enhance flood dynamics and reduce computation time in local water infrastructures. It couples its 1D/2D submodels using two strategies: water level balancing at supernodes and lateral flux exchange between 2D nodes and 1D sections through spillways. This coupling strategy ensures mass balance and momentum between the 1D grid and 2D mesh without the need for a looped algorithm at the integration segment. The comparison of daily discharge in VMD indicates that F28 excels in peak flow performance without compromising total flow volume. Additionally, simulated discharges maintain consistency and reliability throughout the watercourse. Our study demonstrates superior NSE for discharge compared to a comparable 1-D-quasi2-D river network model by Duc Tran *et al.* (2018) using the MIKE Model. While Duc Tran *et al.* (2018) reported NSE values of 0.88, 0.92, 0.9, and 0.92 for Tan Chau, Chau Doc, Vam Nao, and Can Tho, respectively, our study consistently achieved NSE values approaching 1, indicating a stronger alignment between observed and modeled discharge and water level. Moreover, the spatial distribution of inundated depth in October 2011 closely resembles the maximum flood extent from Song & Bakker (2011). The results show that the F28 can properly simulate flux exchange between 1D grid and 2D mesh in a complex basin such as VMD. In addition, The F28 model adeptly simulates flux exchange between the 1D grid and 2D mesh in complex basins such as VMD. Furthermore, F28 leverages a finite volume solution algorithm, a semi-analytical coupling method, and a parallel 1D/2D approach, significantly enhancing numerical stability and precision. This enables both 1D and 2D submodels to run efficiently with larger time steps of up to 45 s. Consequently, the computational time is impressively reduced to just 3 h for a single run spanning 375 days, using a 2D mesh with 38,791 quadrilateral elements averaging 1 km and 50,654 nodes, on a PC equipped with an i9 core clocked at 5.2 GHz. Although a direct performance comparison was not conducted in this study, F28's exceptional performance is evident. Unlike other coupling models that require smaller grid sizes and computation time steps to prevent overreliance on the smallest mesh cells, F28 stands out (Morales-Hernández *et al.* 2016; Fan *et al.* 2017). Other models of similar scale demand grid sizes from 0.05 m to a few meters, with computational times spanning 2–120 h for convergence (Cardoso *et al.* 2020). Hence, given its efficacy and robustness, F28 stands as a suitable choice for simulating floods in densely developed infrastructure basins.

We observe that F28 performs relatively better in flood conditions than in drought. This could have been due to the model's attention to simulating flood dynamics. In contrast, groundwater interaction with surface water dominates during the dry season. We are currently developing a module for groundwater interaction with surface water that can adopt the F28 model for better water resources management.

## CONCLUSIONS

The study developed a new coupling strategy for 1D/2D hydraulic model to simulate the hydrodynamics of a large river and floodplain more flexibly, effectively, with reasonable accuracy and reduced time. With two robust coupling strategies: supernodes and the lateral exchange of fluxes, the new model approach archived solving the water balance and conservation of momentum at interaction boundary 1D/2D, providing more reliable dynamic floodplains with less computational time. Furthermore, by employing an unstructured 2D mesh, the hydrology model demonstrates remarkable versatility in representing various types of infrastructure within large rivers, including sluice gates, dikes, and roads, without the constraints of limited mesh grid size. Thus, the model could be a new helpful tool for further flood study in the lowland, unique with complex irrigation basins. Besides, the F28 should account for the groundwater flux exchange to the main river for simulating discharge in dry periods.

## ACKNOWLEDGEMENTS

This research is funded by Ministry of Science and Technology under grant number KC-4.0-09/19-25.

## DATA AVAILABILITY STATEMENT

All relevant data are included in the paper or its Supplementary Information.

## CONFLICT OF INTEREST

The authors declare there is no conflict.

## REFERENCES

*Mike 11 reference manual*

*.*