Today, due to the increased computational power of advanced computers, interest in numerical techniques increased dramatically. The solution of the equations of fluid mechanics on computers has become so important that it now occupies the attention of perhaps a third of all researchers in fluid mechanics and the proportion is still increasing. This field is known as computational fluid dynamics (CFD). Contained within it are many subspecialties. We shall discuss only a small subset of methods for solving the equations describing fluid flow and related phenomena.

Flows and related phenomena can be described by partial differential equations, which cannot be solved analytically except in special cases. To obtain an approximate solution numerically, we have to use a discretization method which approximates the differential equations by a system of algebraic equations, which can then be solved on a computer.

The approximations are applied to small domains in space and/or time so the numerical solution provides results at discrete locations in space and time. Much as the accuracy of experimental data depends on the quality of the tools used, the accuracy of numerical solutions is dependent on the quality of the discretizations used. Contained within the broad field of computational fluid dynamics are activities that cover the range from the automation of well-established engineering design methods to the use of detailed solutions of the Navier-Stokes equations as substitutes for experimental research into the nature of complex flows.

• In order to solve a given PDE by numerical methods, the partial differentials of the dependent variable in PDEs must be approximated by finite difference relations (algebraic equations).
• The solution of a steady PDE in a two-dimensional rectangular domain with initial and boundary conditions.
• Division of domain in uniform/non-uniform mesh.

• Formulation of algebraic equations for each grid point/control volume (matrix system Ax = b).
• Methods for finite difference approximations.

1. Taylor series expansions.
2. Finite Difference by Polynomials.

• Divide the continuous domain into a number of discrete subdomains (control volumes) by a grid. The grid defines the boundaries of a control volume, while the computational node lies at the center of each control volume.

• For each sub-domain, derive governing algebraic equations from the governing differential equations.
• Obtain a system of algebraic equations from above.
• Solve the above system of algebraic equations to obtain values of the dependent variables at identified discrete points (computational nodes).

The model usually used to describe the free surface flows is based on the two-dimensional Saint-Venant equations. These quite known models in the literature equations are obtained by integrating the vertical dimensional incompressible Navier-Stokes equations under the assumptions of hydrostatic pressure and averaged velocities in the vertical level. Terms from the turbulence, viscosity, and the Coriolis forces are not considered in this study. The system can be set as conservative form:

where u and v are the depth-averaged water velocities in x and y direction, h the water depth, g the gravitational acceleration,  Sox and Soy are respectively slopes following the direction x and y. They are defined by Sox=dZ/dx and Soy=dZ/dy.

We consider here a diffusive flux F(x,t) of form F(x,t)=−Λ(x,t)∇u(x,t), where Λ(x,t) is a linear operator from Rd to Rd, and u is one of the continuous unknowns (uj)j=1,…, N of the problem. Such an expression of the flux is encountered when using e.g. Fourier's law (heat conduction), Darcy's law (flow in porous media), or Fick's law (chemical diffusion).

In some cases, a simple, conservative, and consistent discrete expression for F(n)K,σ may be obtained. This is the case when Λ(x,t) is the identity operator and when the grid satisfies the following "orthogonal property": in each control volume KM , there exists a particular point xK, called the "center" (or centroid) of the control volume, such that, for any adjacent control volume L to K with common face σ, the straight line (xK, xL) is orthogonal to σ.

Figure 2: Two control volumes with a common face, satisfying an orthogonality property

This orthogonality property is satisfied by several types of meshes, such as triangles, rectangles, Voronoï boxes. Then a simple finite difference expression for the approximation of  

Consider the partial differential equation

∂tu(x,t)+∇⋅F(u(x,t),x,t)=0.

This is the conservation equation (1) with A(x,t)=u(x,t) , F(x,t)=F(u(x,t),x,t) and S(x,t)=0 . The linear transport equation is a particular case, with F(u,x,t)=uV(x,t) , and V is defined from Ω×[0,T] to Rd . The discrete unknowns are the values u(n)K , which are expected to be approximations of u in the control volumes K∈M at time t(n) . Since the convection flux is of order 0 (no derivative involved) its approximation does not require any assumption on the mesh, contrary to the diffusion flux. Nevertheless, the fluxes must be carefully approximated in order to ensure the stability (and convergence) of the scheme. Consider the flux ϕ(n)K,σ through an interface σ at time t(n) , seen as a function of the solution u