close

Вход

Забыли?

вход по аккаунту

?

307

код для вставкиСкачать
INTERNATIONAL JOURNAL FOR NUMERICAL METHODS IN ENGINEERING
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
Development of physically based meshes for two-dimensional
models of meandering channel #ow
M. S. Horritt*,s
School of Geographical Sciences, ;niversity of Bristol, Bristol BS8 1SS, ;.K.
SUMMARY
The choice of mesh generation and numerical solution strategies for two-dimensional "nite element models
of #uvial #ow have previously been based chie#y on experience and rule of thumb. This paper develops
a rationale for the "nite element modelling of #ow in river channels, based on a study of #ow around an
annular reach. Analytical solutions of the two-dimensional Shallow Water (St. Venant) equations are
developed in plane polar co-ordinates, and a comparison with results obtained from the TELEMAC-2-D
"nite element model indicates that of the two numerical schemes for the advection terms tested, a #ux
conservative transport scheme gives better results than a streamline upwind Petrov}Galerkin technique. In
terms of mesh discretization, the element angular deviation is found to be the most signi"cant control on the
accuracy of the "nite element solutions. A structured channel mesh generator is therefore developed which
takes local channel curvature into account in the meshing process. Results indicate that simulations using
curvature-dependent meshes o!er similar levels of accuracy to "ner meshes made up of elements of uniform
length, with the added advantage of improved model mass conservation in regions of high channel
curvature. Copyright ( 2000 John Wiley & Sons, Ltd.
KEY WORDS:
mesh generation; shallow water #ow; model error
INTRODUCTION
Finite element models solving the two-dimensional shallow water equations are currently the
most sophisticated tools capable of predicting free-surface #ows in compound channels at the
reach scale (typically up to 30 km long). While these models are capable of making distributed
predictions of free-surface height and depth-averaged velocities, the modelling process is fraught
with uncertainties. The governing equations are incapable of capturing three-dimensional
processes, and therefore do not describe the #ow physics completely. The high Reynolds numbers
associated with such #ows mean that some practical technique of dealing with turbulence is also
* Correspondence to: M. S. Horritt, School of Geographical Sciences, University of Bristol, University Road, Bristol
B58 1SS, U.K.
s E-mail: Matt.Horritt@bris.ac.uk
Contract/grant sponsor: U.K. Natural Environment Research Council; contract/grant number GR3 CO/030
Copyright ( 2000 John Wiley & Sons, Ltd.
Received 5 February 1999
Revised 18 June 1999
2020
M. S. HORRITT
needed. The "nite element technique only delivers an approximate solution to these equations,
and the accuracy of this solution may be a function of "nite element mesh discretization and the
exact formulation of the numerical solver used. There is also uncertainty present in the data input
to these models [1] (boundary friction parameterization, boundary conditions, bed topography
[2]) which will a!ect model output.
Fluvial process studies [3] can shed light on the "rst of these issues, by identifying the more
signi"cant #ow mechanisms and structures and assessing their importance. Recent advances in the
provision of distributed input data (especially through remote sensing techniques [4]) can also
reduce uncertainty in the modelling process, and may encourage a move away from calibration
studies, where these parameters are adjusted to give a best "t between model predictions and
real-world observations, to more physically based parameterizations. Fluvial #ow models have
been calibrated and validated through comparison with environmental data [5], and this will tend
to result in "nal model uncertainties which stem from all three sources listed above. It is thus
di$cult to pinpoint exactly where a model is going wrong, and a more appropriate approach would
be to separate these three sources of uncertainty and hence quantify each of them.
The aim of this paper is thus an assessment of the accuracy of "nite element solutions, i.e. the
ability of the "nite element model to give accurate solutions of the shallow water equations.
Fundamental to this problem is that of the "nite element mesh discretization. While much
literature concerned with the generation of unstructured triangular meshes [6}8] exists, and
techniques for adapting meshes to the spatial properties of the "nite element solutions (for
example [9]) are available, mesh generation for #uvial applications has been governed by
practical considerations. The need for at least three elements across a channel to capture, even in
a crude fashion, variations in velocity and the channel bathymetry, generally leads to meshes of
a prohibitively high number of elements if such triangular meshing strategies are adopted for
reaches which can be up &10 km in length, but with channels only &10 m wide. A typical
method of reducing the computational burden is to discretize the channel with elements elongated in the #ow direction, thus allowing a reasonable representation of channel cross-section
and cross-channel variations in predicted variables without an excessive number of elements.
A region of mesh of this type is shown in Figure 1, where a channel made up of elements with an
aspect ratio of 5 : 1 is seen to be embedded in a #oodplain of more equilateral elements [5]. The
choice of this aspect ratio is essentially arbitrary, and its e!ect on the quality of model solutions
for #ow in the channel has yet to be assessed.
Figure 1. Part of a "nite element mesh for #ood #ow modelling. The channel is described by
elongated elements (with a length-to-breadth ratio of 5) embedded in the approximately
equilateral elements of the #oodplain.
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
DEVELOPMENT OF PHYSICALLY BASED MESHES
2021
The mesh size a!ects the model's ability to represent second-order process gradients, which can
a!ect bulk #ow properties, as seen in [10]. It is unclear, however, whether the use of a "ner mesh
merely a!ects the detail at which the solutions are represented or can make a more fundamental
systematic di!erence to model predictions. The choice of mesh resolution is often based on
computational considerations, opting for the "nest mesh feasible within the constraints of the
hardware available and reasonable simulation times. The choice of numerical solver used has also
been seen to a!ect bulk #ow properties [10], but previous studies have tended to concentrate
upon the solver's computational e$ciency (again a pragmatic consideration) and intercomparisons between the results from di!erent solvers or "nite element models [11], rather than using
objective tests of model accuracy. This has mainly been due to a lack of detailed validation data
against which model results can be compared.
To overcome this lack of validation data, we here develop analytical solutions of the model
governing equations to aid an assessment of the e!ects of mesh discretization and solver
formulation when compared with the "nite element solution. The development of analytical
solutions of the shallow water equations raises a dichotomy with respect to the choice of test
strategy. A modelling situation is required which is on the one hand complex enough to be
non-trivial (and interesting) but on the other simple enough to be analytically soluble for the set of
simultaneous non-linear partial di!erential equations. The analytical approach does, however,
have the advantage of being able to invoke a number of symmetry arguments which will not be
apparent in the model formulation. For example, a transformation of the governing equations
into plane polar co-ordinates can reveal rotational symmetry in the #ow, a symmetry that will not
be apparent to the "nite element model working in an essentially cartesian system. Thus, the
problem of #ow around a channel which is rectangular in cross-section but annular in planform
may succumb to analytical methods and yet provide enough variations in velocity magnitude and
direction and free surface height to make an interesting "nite element problem.
The result will be a rationale for the "nite element discretization of the channel region for
#uvial #ow models, based on the accuracy of model solutions rather than rule of thumb. The
technique also isolates and quanti"es one of the sources of uncertainty in the modelling process,
an advantage over comparisons with "eld or laboratory data, where numerical e!ects may be
masked by poor process representation and model parameterization.
MODELLING AND ANALYSIS OF AN ANNULAR REACH
The TELEMAC-2-D model [10,12] is used throughout this study, having undergone extensive
development in its application to #uvial #ow problems and validation/calibration with respect to
point hydrometric data [5]. The technique will, however, be applicable to any model solving the
two-dimensional shallow water equations. The model aims to solve the two-dimensional shallow
water equations, which can be written as a vector equation for the conservation of momentum
and a scalar equation for the conservation of mass:
Lv
l
gDvDv
#(v ) +)v#g+(z #h)! t +(h+ ) v)#
"0
0
Lt
h
C2h
(1)
Lh
#+(vh)"0
Lt
(2)
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
2022
M. S. HORRITT
where v describes the two-dimensional depth-averaged velocity, h is the water depth, z the bed
0
elevation, g the acceleration due to gravity, l the turbulent kinematic eddy viscosity, and C is
5
Chezy's coe$cient of bed friction. In this study, in order to obtain a comparison with an
analytical solution, the equations are simpli"ed by neglecting the turbulent momentum di!usion
term. This may seem to be an oversimpli"cation, but the e!ect of l on reach length studies has
5
been found to be small [10], and the value of l poorly de"ned [5]. A further simpli"cation is that
5
this study is limited to steady-state #ows. In this case, the time dependence of Equation (1) is used
as a surrogate for an iterative procedure to "nd steady-state solutions, and Equation (2) reduces
to a statement of the non-divergence of the mass #ux.
To achieve an analytic solution, these equations (simpli"ed by assuming steady state and
l "0) can be written in polar form in terms of the depth h and v , v , the radial and tangential
5
r h
velocity components. The assumption of radial symmetry (for depth, velocity and bed elevation)
means all derivatives with respect to the polar angle, h, are eliminated, and the continuity
equation becomes
1 d
(rv h)"0
r dr r
(3)
With all variables now simply functions of the radius, r, the partial derivatives are replaced with
normal derivatives. The derivation of the polar versions of the spatial derivatives found in
Equations (1) and (2) can be found in [13]. If the side walls of the channel are impermeable, v will
r
be zero at the walls, and from Equation (3) zero across the whole cross-section (assuming
non-zero depth). Thus, all radial velocities can be eliminated, and the equations for radial and
tangential momentum become
dh v2
! h "0
dr
r
(4)
g dz
gv2
0# h "0
C2h
r dh
(5)
g
These two equations can be interpreted physically: the "rst equates the free surface slope in a radial
direction with the centripetal acceleration required to maintain the circular motion of the #uid, the
second balances down reach slope against bed friction. If we assume constant friction and down
reach slope, v2 is easily eliminated and the equations integrated to give the analytical solutions:
h
C
C2 dz
0
h"D exp !
gr dh
C
D
dz DC2
C2 dz
0
exp !
v2" 0
h
dh r
gr dh
(6)
D
(7)
where D is a constant of integration. It will be useful to derive the volumetric #ow rate, Q, down
the reach from these expressions
r1
Pr hvh dr
Q"
(8)
0
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
DEVELOPMENT OF PHYSICALLY BASED MESHES
2023
where the limits r and r are the inner and outer radii of the channel. The expression de"es
0
1
analytic integration, but can be integrated numerically (using the trapezium rule) or can be
written as a Taylor expansion in powers of the dimensionless channel width, u"(r !r )/r . The
1
0 0
power series converges to within 1 per cent with terms up to and including O(u3). This enables the
constant of integration D to be derived from a given #ow rate Q.
The experimental methodology runs as follows. Eleven "nite element meshes were generated of
channels 20 m in width and of inner radii 50 and 80 m (Figure 2). These parameters were taken as
being representative of a small U.K. meandering river, such as the upper Thames. Meshes both of
unstructured equilateral elements and structured elongated elements are included. The down
reach slope was set at 0.001, giving values of dz /dh of !0.06 m for the smaller meshes and
0
!0.09 m for the larger. The Chezy friction coe$cient was set at 40 m1@2 s~1 and the #ow rate at
40 m3 s~1, a typical "gure for bankful #ow in such a channel. The boundary conditions for the
TELEMAC-2-D model were the imposed #ow rate at the upstream end (#ow is in an anticlockwise direction) and periodic boundary conditions were applied at the downstream end, the water
depths at this boundary being "xed as equal to those at the equivalent nodes on the upstream
boundary. These boundary conditions leave the model free to "x both velocities and water
depths, but still gives a well-posed numerical problem.
The TELEMAC-2-D model incorporates a number of options for the numerical solution of the
shallow water equations [12], and two of these are investigated here. TELEMAC-2-D generally
uses a fractional step technique: advection terms are solved initially, with propagation, di!usion
and source terms solved in a second step, with these steps being repeated to achieve iterative
convergence. The "rst scheme uses a conservative streamline-upwind Petrov}Galerkin (SUPG)
technique for the advection of depth in the continuity equation, and the method of characteristics
for the advection of velocity. This scheme solves the continuity equation exactly (no fractional
step method is used), and should, in theory, give exact mass conservation. The second scheme uses
a #ux conservative transport (FCT) method for the advection both of depth and velocity, and
although not guaranteed to be totally mass conservative, mass error can be reduced by using
a number of sub-iterations of the advection solving step. The mass conservative nature of the
SUPG scheme without the need for sub-iterations makes it a popular choice for reach length
simulations [5]. The combination of the 11 computational meshes and the 2 advection schemes
gives a total of 22 simulations.
Figures 3 and 4 show the modelling results using the SUPG and FCT schemes in comparison
to the analytical solutions for mesh 3, the "nest mesh used in the study. The FCT scheme gives
markedly better results, both in terms of the absolute errors in depth and velocity and with
respect to the cross channel velocity gradient. The poorest results were those from mesh 6
(Figures 5 and 6), where #ow depth has been over-predicted by '50 per cent for the SUPG
scheme. In terms of mesh discretization, it was found that the ratio of mesh size to channel
planform radius (or the elemental angular deviation) was the most signi"cant control on model
prediction errors, having a far greater in#uence than the number of elements across the channel,
the choice of structured or unstructured mesh or the channel radius. Figures 7 and 8 summarize
this result, demonstrating that both schemes converge with decreasing angular deviation, but
with FCT giving the lower overall errors. The dependence of model accuracy on element angular
deviation could be a result of two e!ects: the angular deviation of the streamline over an element
may a!ect the "nite element model, and the representation of the channel itself is di!erent. Each
mesh represents the boundary of the circular channel as two concentric polygons, and for the
"ner meshes, these polygons coincide more closely with the true circular boundaries.
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
2024
M. S. HORRITT
Figure 2. Finite element meshes for modelling #ow around an annular reach.
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
DEVELOPMENT OF PHYSICALLY BASED MESHES
2025
Figure 2. Continued.
Two obvious conclusions can be drawn from these results. Firstly, the FCT scheme gives more
accurate solutions than the SUPG scheme, despite the better mass conservation properties
theoretically achievable with the SUPG scheme. This may in part be a property of this relatively
simple test case, and when dynamic e!ects and turbulence come to be taken into account, the
signi"cance of the advection terms may be reduced, and the two schemes may give more similar
results. Given the lack of validation data for more complex #ows, however, this is mere
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
Copyright ( 2000 John Wiley & Sons, Ltd.
Figure 4. (a) Comparison between model predicted depth and true depth for the FCT scheme operating on mesh
3. (b) Comparison between model predicted velocity and true velocity for the FCT scheme operating on mesh 3.
Figure 3. (a) Comparison between model predicted depth and true depth for the SUPG scheme operating on mesh
3. (b) Comparison between model predicted velocity and true velocity for the SUPG scheme operating on mesh 3.
2026
M. S. HORRITT
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
Copyright ( 2000 John Wiley & Sons, Ltd.
Figure 6. (a) Comparison between model predicted depth and true depth for the FCT scheme operating on mesh
6. (b) Comparison between model predicted velocity and true velocity for the FCT scheme operating on mesh 6.
Figure 5. (a) Comparison between model predicted depth and true depth for the SUPG scheme operating on mesh
6. (b) Comparison between model predicted velocity and true velocity for the SUPG scheme operating on mesh 6.
DEVELOPMENT OF PHYSICALLY BASED MESHES
2027
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
2028
M. S. HORRITT
Figure 7. Height prediction error as a function of element angular deviation
for the SUPG and FCT schemes.
Figure 8. Velocity prediction error as a function of element angular deviation
for the SUPG and FCT schemes.
speculation, and the poor treatment of the advection term by the SUPG scheme would indicate
that it should be abandoned in favour of the FCT scheme. Secondly, the strong dependence of
solution accuracy upon the elemental angular deviation would imply that local channel curvature
should be taken into account in mesh generation strategy, rather than simply using meshes with
uniform element lengths. This is explored further in the next section.
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
DEVELOPMENT OF PHYSICALLY BASED MESHES
2029
A CURVATURE-DEPENDENT MESH GENERATION STRATEGY
The results given above can be used to formulate a channel meshing strategy given the level of
accuracy required of the "nite element solution. Mesh generation with the constraint of uniform
angular deviation would seem to be a better use of computational resources than elements of
uniform length, concentrating computational resources in bend regions where the solutions are
likely to be less accurate due to the high local streamline curvature. To this end, a structured "nite
element channel mesh generator is developed, which calculates element lengths according to local
channel curvature.
The local curvature, i of an arc can be de"ned as
dh
i"
ds
(9)
where h is the angle between the tangent to the curve and some arbitrary direction, and s is the arc
length parameter along the curve. i can thus be thought of as factor for converting between
lengths and angular deviations, and a mesh generator should aim to produce elements not of
equal length, l, but with equal DiDl. Since information on the channel planform is likely to be in the
form of a series of nodes de"ning the channel banks, a more useful "nite di!erence approximation
of curvature at a node is given by
(v1]v2)
2
i"
'
Dv1 D Dv2 D (Dv1 D#Dv2 D)
(10)
where v and v are vectors linking the previous and successive nodes to the one in question. The
1
2
"rst term gives the angular deviation dh (using the small angle approximation sin h&h), the
second approximates the change in arc length ds.
The channel form is likely to be subject to a certain degree of noise (uncertainty in the position
of the individual nodes) due to the digitization process and errors in the surveying of the channel,
and this may generate high values of DiD, obscuring the underlying channel curvature. Coherent
smoothing of i (observing the sign of i) along the channel can eliminate this noise, but this can
have the drawback of producing zero curvature in cross-over regions, where the contributions
from convex and concave sections before and after cancel. Further incoherent smoothing (DiD is
smoothed rather than i) is used to avoid this problem, allowing DiD to di!use into these cross-over
regions and promoting the use of shorter elements therein. This smoothing procedure is illustrated in Figure 9.
Having established an appropriate measure of curvature along each of the banks of
the channel, the meshing algorithm proceeds as follows. A cross-section moves along the reach,
in the direction of a vector determined from the mean of the direction of the two banks. The
distance from the start of the reach, l, is measured, and if the quantity DiDl exceeds a threshold, the
end of the "rst set of elements is inserted. The process is then repeated, inserting elements
whenever the length from the last set, l is such that DiDl exceeds the threshold. A second threshold
value for l itself is also set, to prevent overly long elements forming in straight sections of the
reach.
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
2030
M. S. HORRITT
Figure 9. Schematic showing action of coherent and incoherent smoothing on the curvature measure
for a meandering channel.
This mesh generation strategy is tested on a &1 km reach of the upper Thames, in
southern England, the results being shown in Figure 10, along with the results of applying
two other mesh generation strategies commonly adopted in #uvial #ow studies [5]. Meshes
1 and 2 were formed using the older strategy of equal size elements (e!ectively "xing i"0
along the reach) with lengths 10 and 30 m. The elemental angular deviations are up to 903 for
mesh 2. Meshes 3 and 4 use the new curvature-dependent strategy with angular deviations
between elements of 0.4 and 0.1 rad. The technique ensures that the meander regions of high
curvature are represented by smaller elements, and that the angular deviation between sets of
elements remains approximately constant (Figure 11). Finally, two meshes (5 and 6) were also
generated using unstructured near-equilateral elements of size 2 and 10 m. These parameters,
along with the number of nodes and elements for each mesh, are summarised in Table I. The six
meshes thus represent coarse and "ne examples of the results of each of the three generation
strategies.
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
DEVELOPMENT OF PHYSICALLY BASED MESHES
2031
Figure 10. Finite element meshes of the Thames reach. Mesh 1 is shown in its entirety, for the others a detail
of a meander region is shown.
PERFORMANCE OF CURVATURE-DEPENDENT MESHES
Steady state simulations were performed using the six meshes of Figure 10 at approximately
bankful discharge (40 m3 s~1) with both the SUPG and FCT schemes. Tests were also performed
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
2032
M. S. HORRITT
Figure 11. Detail of mesh 4, demonstrating the constant element angular deviation property.
Sets of edges between elements from regions of the channel with di!erent curvatures subtend
similar angles at points a and b.
Table I. Parameters for the six Thames meshes.
Mesh
1
2
3
4
5
6
Element size
(m)
Angular
deviation (rad)
10
30
0.4
0.1
2
10
Nodes
Elements
378
130
146
470
7685
522
562
190
214
700
14370
770
on meshes of intermediate resolution, but the results were predictable from those from the coarse
and "ne meshes, and have therefore been omitted for the sake of clarity. The results, in the form of
water depths as functions of down reach distance, are shown in Figures 12 and 13. Cross-reach
variations in water level were eliminated by taking the mean depth across the channel. The
predicted water levels are more consistent between meshes for the FCT scheme, where some level
of convergence can be seen (the similarity between prediction from meshes 4 and 5 for example).
The "gures demonstrate the advantages of the new curvature dependent meshing scheme: the
results from mesh 4 are similar to those from mesh 5 but with fewer elements, and the advantage
of shorter computation times. This advantage will become greater when such channel meshes
are embedded in larger meshes for the simulation of #oodplain #ows. Because process
gradients on the #oodplain are lower than those in the channel, generally larger elements can be
used, but there is often an upper limit imposed on their size by the need to couple with the mesh in
the channel. The curvature-dependent scheme will impose small #oodplain elements on the inside
of meander bends, but should make for larger elements surrounding straighter sections of the
reach.
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
DEVELOPMENT OF PHYSICALLY BASED MESHES
2033
Figure 12. Water depths (averaged across the channel) predicted using the SUPG scheme as a function of
along reach distance, for the six "nite element meshes.
Figure 13. Water depths (averaged across the channel) predicted using the FCT scheme as a function of
along reach distance, for the six "nite element meshes (note the change in scale of the vertical axis).
A measurement of mass conservation error is often used as a measure of con"dence in model
predictions in situations such as these where no validation data is available. Mass conservation
may be de"ned on a global basis, in terms of mass #uxes at the boundaries and the rate of change
of volume of #uid within the domain (which will be zero at steady-state), but it may be more
instructive to look at mass balance errors on a local basis. For a steady-state simulation, the
divergence of the mass #ux vector (vh) should be zero, and its actual value will be a measure of
local mass balance error. This can be calculated for each element by adding the #uxes across its
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
2034
M. S. HORRITT
Table II. Mass balance errors for the Thames reach simulations.
Mesh/advection
scheme
1
1
2
2
3
3
4
4
5
5
6
6
Final d</dt
(m3 s~1)
: +(vh) dA
(m3 s~1)
[1 :(+(vh))2 dA]1@2
A
(cm s~1)
0.16
0.90
!0.65
0.011
!0.65
0.45
0.22
!0.15
(0.001
(0.001
(0.001
(0.001
1.36]10~3
!5.68]10~1
!2.76]10~2
4.71]10~1
!5.82]10~2
!7.72]10~1
!3.13]10~3
!3.50]10~1
5.52]10~4
!8.85]10~1
!8.27]10~3
!1.23]10~1
2.38
2.70
4.24
4.02
4.11
3.80
2.52
3.04
0.96
1.03
1.27
1.28
SUPG
FCT
SUPG
FCT
SUPG
FCT
SUPG
FCT
SUPG
FCT
SUPG
FCT
three boundaries and dividing by the element area (according to the divergence theorem). This is
essentially a "nite volume measure of mass balance. The integral of the divergence over the area of
the model domain is presented in Table II. Some of this error may be assigned to the fact that the
simulations are not true steady state, but are really dynamic simulations with steady-state
boundary conditions, run for long enough for transients to die down. Rates of change of volume
in the domain are also given in Table II, measured over the last two time steps of each simulation.
For a simulation that conserved mass, d</dt and the integral of the divergence would be equal in
magnitude and opposite in sign (a result obtained by integrating the continuity equation over the
domain), but there appears to be no such correlation in Table II. This means that some of the
mass balance error is genuine, and not due to transients still present in the simulation. A second
measure of the mass balance error can be found by calculating the root-mean-square value of the
divergence over the reach, and is also shown in Table II. The "gure can be interpreted as
the volume error per unit area over the reach, or the depth of water lost or created averaged over
the mesh. All the errors are of the same order, but the errors for the triangular meshes are lower
than those for the regular meshes.
As a further demonstration of potential problems with mass conservation, Figures 14 and 15
show down reach #ow (which should be uniformly 40 m3 s~1) for meshes 1 (10 m elements) and
4 (angular deviation of 0.1 rad). The #ow for mesh 1 varies by up to 5 per cent (2 m3 s~1), mainly
in regions of high curvature. Mesh 4 performs better in the main section of the reach, but there are
large #uctuations in #ow towards the upstream boundary, perhaps associated with the large
elements in this part of the reach. Elsewhere, the discrepancy is generally less than 0.5 m3 s~1.
Both numerical advection schemes (SUPG and FCT) give similar mass balance errors. The better
mass conservation properties of the curvature dependent mesh are a further advantage of its use,
but some special treatment at the upstream end of the reach may be required if adequate mass
conservation is required in this region.
The #uctuations in along reach #ow are mainly due to velocity variations, the depths predicted
by TELEMAC-2-D being relatively much smoother and free from node-to-node oscillations (the
SUPG and FCT schemes are very good at damping oscillations in depth). Both numerical
schemes seem to encounter problems dealing with the non-linear velocity advection term,
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
DEVELOPMENT OF PHYSICALLY BASED MESHES
2035
Figure 14. Flow rates as functions of along reach distance for the two numerical schemes
operating on mesh 1.
Figure 15. Flow rates as functions of along reach distance for the two numerical schemes
operating on mesh 4.
a problem not encountered with the study of the annular reach. The problem only arises in the
Thames reach probably because irregularities in channel planform excite these modes of oscillation, especially in regions of high curvature and at the upstream boundary of the domain. The
modes of node to node oscillations will tend to be aligned with the channel direction for meshes
1}4 and at approximately 303 to the channel for the triangular meshes (Figure 16). This may
explain the poorer mass conservation properties (at least according to the RMS measure) of the
regular meshes, as the oscillations are able to propagate freely along the channel.
The poor performance of even the conservative SUPG scheme is partly due to the "nite volume
method used to assess mass balance in this study. The conservative scheme ensures mass
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
2036
M. S. HORRITT
Figure 16. (a) a possible mode of node to node oscillation for a regular mesh (e.g. meshes 1}4) of elongated
elements. (b) the mode of oscillation for approximately equilateral elements is no longer aligned with
the channel direction.
conservation integrated over the "nite element test functions, allowing velocity variations to
cancel and conserve mass on a global basis (which is the usual measure used). The results do,
however, demonstrate that even though the conservative SUPG scheme gives much better results
in terms of global mass balance, both schemes perform similarly in terms of local mass balance.
This in some way negates one of the advantages of the conservative SUPG scheme, in that its
globally conservative property is not re#ected in the local mass balance, and even simulations
that appear to conserve mass at the boundaries of the domain may have local mass balance
problems, caused by unphysical spatial #uctuations in the predicted velocity "eld.
CONCLUSIONS
Uniform #ow around an annular reach has proved to be a good test of the TELEMAC-2-D
model when compared with analytical solutions of the shallow water equations, discriminating
between di!erent numerical schemes and "nite element meshing strategies. The #ux conservative
transport scheme for the advection of velocity and depth gives more accurate solutions than the
streamline upwind Petrov}Galerkin technique when compared to the analytical solutions. The
solution accuracy is a strong function of the elemental angular deviation, and this measure should
be used to drive any mesh generation strategy, rather than simply prescribing element length. The
analysis of this simple scenario has thus given a rationale for the choice of numerical scheme and
"nite element mesh. With this rationale in place, a mesh generator for meandering #uvial
channels has been developed, with element sizes dependent on a measure of local channel
curvature. Meshes generated using this technique appear to give better results than those of
uniform element size, when compared to irregular meshes of small triangular elements. Despite
the lack of validation data for this reach, the FCT scheme still seems to o!er an advantage over
the SUPG scheme, in terms of the apparent convergence of the solutions as progressively "ner
meshes are used. Both numerical schemes su!er from mass conservation errors due to unphysical
node-to-node oscillations in the predicted velocity "eld, but these are smaller in magnitude for the
curvature-dependent meshes, o!ering a further incentive for their use.
It must be noted that this modelling study makes certain simpli"cations to facilitate the
analytical solution of the shallow water equations. Turbulence has been neglected, but its e!ect
on #oodplain #ow studies is as yet only poorly understood, and this o!ers some justi"cation for
its omission. The assumption of a constant eddy viscosity is in itself a gross oversimpli"cation,
and the use of a more sophisticated model (such as the k}e model) may have a bearing on the
mesh generation strategy. Model accuracy may change when #ow through more complex channel
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
DEVELOPMENT OF PHYSICALLY BASED MESHES
2037
geometries (changes in curvature, variable widths, complex bathymetry) are simulated, and the
modelling of overbank #ows (with the implication of complex interactions between channel and
#oodplain #uid) may introduce further inaccuracies. The simple relationship between model
accuracy and streamline curvature found here may then no longer hold. It is also unclear from
this study whether model accuracy is solely dependent on streamline curvature, or whether the
representation of the model boundary is as important, but either way, the curvature-dependent
mesh generation strategy is still valid. Dynamic simulations may pose more problems, as the "nite
element mesh may have an e!ect on the subtle changes that will need to be represented to model
the passing of a #oodwave that may last hours or days.
These results will have some practical implications for future #uvial modelling studies. Given
the error as a function of angular deviation (as in Figures 7 and 8), and the accuracy required of
model predictions, an appropriate mesh resolution can be selected. Conversely, using a given
elemental angular deviation, the uncertainty in the "nite element solutions due to the approximate way in which the governing equations are solved can be estimated. In studies comparing
modelled results to laboratory or "eld data, this numerical error can be discounted, and the
remaining discrepancies assigned to poor process representation and model parameterization.
This will shed light on how well the governing equations describe the physics of the problem, and
aid in the calibration and validation process.
ACKNOWLEDGEMENTS
This research has been supported by the UK Natural Environment Research Council grant GR3 CO/030.
Thanks go to Paul Bates (University of Bristol) and Jean-Michel Hervouet (Laboratoire National
d'Hydraulique, ElectriciteH de France) for comments on the text.
REFERENCES
1. Anderson MG, Bates PD. Evaluating data constraints on two dimensional "nite element models of #oodplain #ow.
Catena 1994; 22:1}15.
2. Bates PD, Horritt MS, Hervouet J-M. Investigating two-dimensional, "nite element predictions of #oodplain
inundation using fractal generated topography. Hydrological Processes 1998; 12:1257}1277.
3. Bathurst JC. Flow processes and data provision for channel #ow models. In Modelling Geomorphological Systems,
Anderson MG (ed.). Wiley: Chichester, 1988; 127}150.
4. Bates PD, Horritt MS, Smith CN, Mason DC. Integrating remote sensing observations of #ood hydrology and
hydraulic modelling. Hydrological Processes 1997; 11:1777}1795.
5. Bates PD, Stewart MD, Siggers GB, Smith CN, Hervouet J-M, Sellin RJH. Internal and external validation of
a two-dimensional "nite element code for river #ood simulations. Proceedings of the Institution of Civil Engineers,
=ater, Maritime and Energy 1998; 130:127}141.
6. Frey WH. Selective re"nement: a new strategy for automatic node placement in graded triangular meshes.
International Journal for Numerical Methods in Engineering 1987; 24:2183}2200.
7. Rebay S. E$cient unstructured mesh generation by means of Delaunay triangulation and Bowyer}Watson algorithm.
Journal of Computational Physics 1994; 110:23}38.
8. Frey WH, Field DA. Mesh relaxation: a new technique for improving triangulations. International Journal for
Numerical Methods in Engineering 1991; 31:1121}1133.
9. Knupp P. Mesh generation using vector "elds. Journal of Computational Physics 1995; 119:142}148.
10. Bates PD, Anderson MG, Price DA, Hardy RJ, Smith CN. Analysis and development of hydraulic models for
#oodplain #ow. In Floodplain Processes, Anderson MG, Walling DE, Bates PD (eds). Wiley: Chichester, 1996;
215}254.
11. Bates PD, Anderson MG, Hervouet J-M. An initial comparison of two 2-dimensional "nite element codes for river
#ood simulation. Proceedings of the Institution of Civil Engineers, =ater, Maritime and Energy, 1995; 112:238}248.
12. Hervouet J-M, Van Haren L. Recent advances in numerical methods for #uid #ow. In Floodplain Processes, Anderson
MG, Walling DE, Bates PD (eds). Wiley: Chichester, 1996; 183}214.
13. Arfken G. Mathematical Methods for Physicists (3rd edn). Academic Press: London, 1985.
Copyright ( 2000 John Wiley & Sons, Ltd.
Int. J. Numer. Meth. Engng. 2000; 47:2019}2037
Документ
Категория
Без категории
Просмотров
4
Размер файла
315 Кб
Теги
307
1/--страниц
Пожаловаться на содержимое документа