Application of Finite Difference Schemes to 1D St. Venant for Simulating Weir Overflow

. Depth averaged equations are commonly used for modelling hydraulics problems. Nevertheless, the model may not be able to accurately assess the flow in the case of different flow regimes, such as hydraulic jump. The model requires appropriate numerical method or other numerical treatments in order to simulate the case accurately. A finite volume scheme with shock capturing may provide a good result, but it is time consuming as compared to the commonly used finite difference schemes. In this study, 1D St. Venant equation is solved using Artificial Viscosity Lax-Wendroff and Mac-Cormack with TVD filter schemes to simulate an experiment case of weir overflow. The case is chosen to test each scheme ability in simulating flow under different flow regimes. The simulation results are benchmarked to the observed experimental data from previous study. Additionally, to observe the scheme efficiency, the simulation time between the models are compared. Therefore, the most accurate and efficient scheme can be determined.


Introduction
Open channel flow is a channel that drains water with free surface. There are many problems that occur in the flow in the open channel, e.g. is the failure of the structure on the weir. Therefore, early prevention is needed at the planning stage. On the planning stage, structure performance can be simulated by physical model or a numerical model, therefore structure failure can be prevented.
The physical model is created when the physical phenomenon of the problems in the prototype can be made on a smaller scale with adequate congruence. Numerical model is created when the problems to be studied can be formulated with the formulation / mathematical equations in detail.
However, physical modelling is usually more expensive than numerical models. Therefore, making a physical model is not economical if the existing problem can be solved using a numerical model and as stated in Toombes and Chanson [7], real-life situations are frequently too complex to solve without the aid of numerical models. Numerical models also offer ease of workmanship that can save time in very wide study that can be produced.
Geometry parameters, calibration process, and validation of simulation results in this study are benchmarked to the observed experimental data from previous study by Binilang [1]. On that previous study, flow is observed in multiple discharge and downstream water level.
This research develops numerical model using one dimensional shallow water equation with TVD filter applied MacCormack and Artificial Viscosity Lax Wendroff discretization scheme. Numerical model results then compared to previous experimental study for validation. Additionally, those schemes are compared to each other to determine the most accurate and efficient scheme.

Governing Equation
The system of partial differential equations which describe the 1-D incompressible flow in an open channel of arbitrary cross section used in this model is continuity and momentum of Saint-Venant equations. It derived from mass and momentum conservation laws shown in equations [1] and [2].

Continuity Equation
(1) where A is cross-sectional area, Q is discharge, t and x are time and space respectively. where g is gravitational acceleration, h is water depth, Sf and S0 is energy slope and bed slope respectively.

MacCormack Discretization Scheme
The application of MacCormack method to the hyperbolic equation proceeds in two steps; a predictor step which is followed by a corrector step. In predictor step, spatial and temporal derivative were made by forward difference scheme. In corrector step, forward scheme is used for temporal derivation and backward scheme for spatial derivation. Final step for (t+1) time discretization is computed as average of predictor and corrector results. In mathematical form, MacCormack discretization scheme described in equations [3] to [7]. For predictor step: Where U is velocity, t and x are time and space respectively. Farid et al. [8] has studied discretization of continuity and momentum equation in [1] and [2] using MacCormack scheme as shown in equations [8] to [15].

Predictor Step Discretization
For continuity equation: ) ∆t (9) where A is cross-sectional area, Q is discharge, t and x are time and space respectively.
where A is cross-sectional area, Q is discharge, t is time, x is space respectively, g is gravitational acceleration, h is water depth, n is manning roughness coefficient, R is hydraulic radius, Sf and S0 is energy slope and bed slope respectively.

Corrector Step Discretization
For continuity equation:

Total Variation Diminishing (TVD) Filter and Averaged Results
As stated in Tadmor [4], every discretization method is associated with an appropriate finite-dimensional projection. It is known that linear projections which are monotone, are the most first order accurate. The lack of monotonicity for higher order projections is reflected by spurious oscillations near jump discontinuities. These are evident with the second -order (and higher) centered difference, whose dispersive nature is responsible to the formation of binary oscillations. TVD schemes avoid spurious oscillations to this end they use the necessarily nonlinear projection. TVD filters, instead, suppress spurious oscillations.
Lal and ASCE [2] stated that the MacCormack method shows oscillatory behaviour near steep gradients and shocks which can be suppressed using a TVD method.

Lax-Wendroff Discretization Scheme
Lax-Wendroff Scheme is modification of FTCS with artificial diffusion. This scheme consists of Lax and Wendroff steps as shown in Figure 1. Both of them are using central difference to get space derivative and forward difference for time derivative. It is needed to apply, first, one step "Lax-step", and later half a time step.

Lax Step
For continuity equation: where A is cross-sectional area, Q is discharge, t is time, x is space respectively, g is gravitational acceleration, h is water depth, n is manning roughness coefficient, R is hydraulic radius, Sf and S0 is energy slope and bed slope respectively.

Wendroff Step
For continuity equation:

Lax Wendroff Artificial Viscosity
The equations of hydrodynamics can be modified to include additional terms as stated in Neumann and Richtmeyer [6] which greatly simplify the procedures needed for stepwise numerical solution of the equations in problems involving shocks. A set of difference equations is given below. Consider: Replace flux by: Giving: The right term in equation (35) simulates the effect of the physical viscosity on the grid scale, concentrated around discontinuity and negligible elsewhere.
Approximation of artificial viscosity with Lax Wendroff discretization gives:

Numerical Parameters
This model developed with cell size of 1 centimetres of 3.61 meters domain, thus resulting 362 number of cell nodes. Courant-type stability criterion is used in this model to calculate the maximum allowed time-step size.
As stated in study conducted by Haun et al. [11], the Courant Number implies how fast the fluid passes through a cell. If the Courant Number is greater than 1, the velocity of a particle is so high that it passes through a cell in less than one-time step. This leads to numerical instabilities. Stability criteria and model trial lead to time step 0.001 seconds avoids instabilities in the computation and makes satisfied model result.

Initial and Boundary Conditions
Boundary condition of this model is static discharge on the upstream and static water depth on the downstream, using value from previous experimental research.

Convergence Criteria
The solution in this model is satisfy if the discharge has converged into upstream boundary value resulting water depth no longer change in time. Discharge is considered to be converged if the sum of discharge deviation of all cell is less than 0.01 m3/s. After the convergence criteria is met, the model should stop running and value of water depth and velocity are recorded.

Model Calibration
Calibrated parameter in this model is manning roughness coefficient. Three main criteria of this calibration process are hydraulic jump position, water depth in supercritical region, and upstream water depth. Calibration is taken using long section profile with discharge 0.00255 m3/s and downstream water level 0.0715 m from previous experimental research. After calibration process is taken, it is obtained different value of manning roughness coefficient on weir body and floor. On the weir body, manning roughness coefficient varies gradually between 0.0175 on the upstream and 0.00835 on the downstream, and constant 0.00835 on the weir floor.

Previous Experimental Research
Previous experimental research which is benchmarked to this research consists of two set of discharge. Each discharge tested with 8 downstream different water depths. L1, L2, and L3 are conditions of non-submerged hydraulic jump; P1 and P2 are conditions of transitional hydraulic jump; T1, T2, and T3 are conditions of submerged hydraulic jump. Value of discharge and downstream water depth are shown in Table 1.

Result and Discussion
The weir experiment consists of a 3.61 m long, 0.075 m wide flat channel with a streamline weir located on the upstream. The model is observed for two discharge values and for each discharge, there are eight downstream water depth variables. Computed water profile compared and examined to the result from previous experimental research as shown in Figure 2 and Figure  The close agreement between the computed results and the experiment for hydraulic jump location, supercritical region and upstream water depth shows that the proposed method is comparatively accurate, which value of error shown in Table 2. However, both models fail to compute length of hydraulic jump since 1D SWE gives a linear height profile before the jump as stated in Singha et al. [10]. In addition, models are not valid for submerged hydraulic jump condition and there is significantly error on the upstream as roller effect influence flow significantly, which is 2D phenomenon. Comparative result between MacCormack and Lax Wendroff schemes for L2, P2, and T2 in Q3 case as shown in Figure 6 shows that Lax Wendroff scheme provides more oscillation result relative to MacCormack. Overall for both schemes, hydraulic jump location, supercritical and upstream water depth match each other with only slight difference due to oscillation on Lax Wendroff's result.

Fig. 6. MacCormack and Lax Wendroff schemes comparison (Q3 Case)
Numerical result is verified using Belanger formula as stated in Liu et al. [12] for the case of L2 (dark red line in Figure 2 to 5). Verification result satisfies with maximum error 0.003 meter, as shown in Table 3. Simulation times for TVD MacCormack scheme to obtain all results in 5 discharges and 8 downstream water depth cases are 56 minutes and 54 minutes for Lax Wendroff scheme.

Conclusion
Numerical model has been developed for 1D shallow water equations (Saint-Venant equations) to compute weir overflow including hydraulic jump phenomenon at downstream. The 1D shallow water equations are solved with the MacCormack and Lax Wendroff scheme discretization with addition of a symmetric total variation diminishing (TVD) term to the corrector step of MacCormack, and artificial viscosity term to Wendroff step of Lax Wendroff. Both schemes show good result to estimate hydraulic jump location, upstream and supercritical water depth. In Lax Wendroff's result, there is small oscillation around bed slope transition and hydraulic jump. Simulation time between both schemes shows almost same compared each other. Therefore, it can be concluded that TVD MacCormack scheme more accurate and efficient to model weir overflow and hydraulic jump compared to Lax Wendroff with Artificial Viscosity, it is because TVD is more effective to filter shocks.
Numerical simulation for the weir overflow and hydraulic jump indicates that this model is fairly effective and accurate for non-submerged condition and can be used to estimate hydraulic jump location for a planning stage of weir construction. Additionally, this model is able to simulate weir overflow and hydraulic jump with dynamic discharge and water depth. High-resolution mathematical model can be further modified and extended to multi-dimensional hyperbolic conservation laws to effectively simulate 2D and 3D weir overflow, in condition of submerged hydraulic jump occur at downstream.