- Email: [email protected]
Modelling bund overtopping using Shallow Water Theory
Journal of Loss Prevention in the Process Industries 23 (2010) 662e667
Contents lists available at ScienceDirect
Journal of Loss Prevention in the Process Industries journal homepage: www.elsevier.com/locate/jlp
Modelling bund overtopping using Shallow Water Theory D.M. Webber, M.J. Ivings* Health and Safety Laboratory, Harpur Hill, Buxton SK17 9JN, UK
a r t i c l e i n f o
a b s t r a c t
Article history: Received 4 February 2010 Received in revised form 7 July 2010 Accepted 7 July 2010
A model is presented of liquid flow under gravity, over horizontal and sloping surfaces, based on the Shallow Water Equations. The aim is to facilitate assessment of the efficacy of bunds surrounding liquid storage tanks, in cases where the liquid front encountering the bund carries significant momentum. The model, called SPLOT, incorporates terms representing the effects of turbulent and laminar friction, and a sub-model for the effect of vertical walls. Having produced a robust numerical solution method, and verified its performance in cases where analytic solutions are known, it is the objective of this work to test the model’s validity against experimental data on flows overtopping bunds. SPLOT is found to be very satisfactory in a wide range of circumstances involving both vertical and sloping bunds. Physical reasons are identified for the overly conservative predictions found in the remaining cases, and a simple parameter is derived which indicates whether the model is being used outside its range of optimum validity. Crown Copyright Ó 2010 Published by Elsevier Ltd. All rights reserved.
Keywords: Shallow water equations CFD modelling Bund overtopping Validation
1. Introduction Assessment of the safety of liquid storage tanks requires an estimation of how much liquid may overtop a surrounding bund in the event of a sudden and significant failure of a tank (Thyer, Hirst, & Jagger, 2002). In such cases it is not sufficient that the bund may simply have the capacity to hold the entire contents of the tank, as a rapidly advancing liquid front may not be entirely held back by the bund wall. In general therefore, it is necessary to understand the dynamics of an advancing liquid front, and what happens when it encounters a wall. This paper summarises a mathematical model of this process, and presents a comparison with some well defined laboratory-scale experiments using water, which measure the quantity of liquid escaping from the bund in a number of different geometries and at different scales. In general, a bunded area may be of arbitrary shape and surround one or more tanks. A starting point for understanding the flow is therefore provided by the two-dimensional Shallow Water Equations for the flow of shallow pools in two horizontal dimensions e see for example Landau and Lifshitz (1959), Webber and Brighton (1986) and Toro (2001). However, two specific limitations of the standard Shallow Water Equations are that they describe inviscid flow, and they are only strictly valid when vertical accelerations can be neglected. These are addressed as follows. A model of the effect of friction (in laminar and turbulent flow) is introduced at the outset and the importance of vertical * Corresponding author. Tel.: þ44 (0) 1298 218133; fax: þ44 (0) 1298 218840. E-mail address: [email protected] (M.J. Ivings).
acceleration is discussed in some detail after the model has been compared with the experimental data, where it will be seen to provide limits on where the model can be expected to perform optimally. The model is defined by the equations:
vh v huj ¼ 0 þ vxj vt
(1.1)
vðhui Þ v hui uj vh vb þ gh ¼ fi gh þ vxj vt vxi vxi
(1.2)
with
cn fi ¼ ui Cjujui h
(1.3)
where h(x, t) is the depth of liquid, u(x, t) is the 2-component horizontal liquid velocity, b(x) is the ground level (but note that only its gradient is important), f(x, t) is the friction term, n is the kinematic viscosity of the liquid (w1.0 106 m2/s for water), c is a laminar friction coefficient, C is a turbulent skin friction coefficient and (t, x1, x2) are time and the two horizontal space coordinates Note that where friction is negligible, the model has no free parameters. The friction model introduces two free parameters c and C which have been discussed in detail by Webber and Jones (1987) who point out that the parameter c, governing friction in laminar flow should take the value 3.0 in order to reproduce the lubrication theory results of Huppert (1982) who studied thin
0950-4230/$ e see front matter Crown Copyright Ó 2010 Published by Elsevier Ltd. All rights reserved. doi:10.1016/j.jlp.2010.07.002
D.M. Webber, M.J. Ivings / Journal of Loss Prevention in the Process Industries 23 (2010) 662e667
liquid films between parallel plates1. For turbulent flow, parameter C is the more important one and is2 of the order of 103. This model, known as “SPLOT”, has been discussed by Ivings and Webber (2002), Ivings, Toro, and Webber (2003) and Ivings and Webber (2007) highlighting some of the problems which must be overcome to produce a numerical solution algorithm which is robust enough to be used routinely in the analysis of industrial hazards, and which provides accurate numerical solutions in cases where analytic properties of the solution are known. The final step, necessary to enable confident use of the model, is to test its validity against experimental data for dynamic bund overtopping. This is the main objective of this work but, in comparing the model with experimental data, we have also investigated its sensitivity both to grid resolution and to the size of the friction coefficients. We have also performed a verification of the inviscid model (c ¼ C ¼ 0) against the expectation from dimensional arguments that the fraction of liquid overtopping the bund should not change if all the length dimensions of a scenario are scaled by the same factor. The experiments considered here all involved complete failure of a cylindrical tank with a concentric cylindrical bund. The model SPLOT is not constrained to this geometry, but it does produce a very useful test case which, as we shall show, allows analysis of the reasons underlying the way in which the model performs. It should also be emphasised that the current work is aimed at validation of the spreading model and the model of the interaction between the spreading liquid and the bund. Vaporising pools are outside the scope of this paper. 2. Experimental data Two series of experiments have been used here to validate SPLOT: those reported by Greenspan and Johansson (1981) and those commissioned by the UK Health and Safety Executive (HSE) (see Appendix) at Liverpool John Moores University (LJMU) and reported by Atherton (2005). Of these the latter are at a somewhat larger scale and use vertical bund walls, and so we shall address them first. The Greenspan and Johansson (1981) experiments used bunds sloping away from the containment area and require a different modelling approach by SPLOT. A detailed (and of necessity somewhat lengthy) tabulated summary of the experimental parameter space (tank and bund geometries), for both series of experiments, is given by Webber, Ivings, Maru, and Thyer (2009). 2.1. The LJMU experiments Of the LJMU experiments (Atherton, 2005), 84 trials were carried out with instantaneous releases of water from a cylindrical tank with a concentric cylindrical bund at some radial distance from it. The scenario addressed is shown in Fig. 1. In order to maximise the scale of the experiments in the laboratory space available, advantage was taken of the expected cylindrical symmetry of the flow, and the whole experiment was conducted between walls meeting at a right angle at the centre of the tank; thus the whole experiment was effectively a quarter of the modelled cylinder.
663
H hb
R
rb
Fig. 1. The geometric parameters in experiments involving a cylindrical tank surrounded by a concentric circular bund.
The salient parameters are the radius of the tank, R, the depth of the water in the tank, H, the radial position of the bund, rb, and the height of the bund, hb. The volume of the liquid released (for the full cylinder), V, and the capacity of the bund, vb, are then given by
V ¼ pR2 H and vb ¼ prb2 hb
(2.1)
In each experiment it was arranged that vb > V so that the bund would be sufficient to retain all of the liquid if it were released slowly. All measured overtopping is therefore a dynamic effect owing to the rapid release of liquid from the tank. The primary point of comparison of the model with the data is the fraction of the released liquid which overtops the bund, Q. The experimental geometry gives three independent dimensionless ratios of the parameters R, H, rb and hb. As long as the scale of the experiments is large enough, it is reasonable to expect the (dimensionless) overtopping fraction Q to be independent of the scale, and therefore populate a three-dimensional domain of results defined, for example, by the values of rb/R, hb/H, H/R. Other possible dimensionless combinations may have more physical significance. For example (when vb < V) 1 vb/V provides an obvious lower bound on the value of Q. We discuss other physically significant dimensionless combinations below. The experiments were performed with R ¼ 0.3 m using instantaneous high, medium and low aspect releases with H ¼ 0.6 m, 0.3 m and 0.12 m respectively. Bund heights hb were in the range 0.03 m e 0.72 m (the last representing a high collar bund) and at radial distances rb in the range 0.315 m e 1.897 m but always such that the bund could in principle contain the entire volume of the tank e i.e. vb > V. 2.2. The Greenspan and Johansson (1981) experiments Greenspan and Johansson (1981) also conducted instantaneous releases of water from a cylindrical tank into a concentric cylindrical bund, but with the bund consisting of an upward slope away from the tank at angles of 30 and 60 to the horizontal. In these experiments the linear scale was about a factor three smaller than those done by Atherton (2005) with tank radius R ¼ 0.0953 m and the fill level H (extracted indirectly from the data by Hirst, 2006) in the range 0.07 me0.26 m. A particular point of interest in these experiments is that the sloping bund may be expected to have different retention properties from those of a vertical wall. 3. Comparison of model and data
1
Huppert was concerned with a thin liquid film between parallel plates, and c ¼ 3 is what is required to match Huppert’s NaviereStokes solution using the above laminar friction term added to the inviscid Shallow Water Equations. It is not entirely obvious that the same value of c will be optimal for a front advancing on to dry ground, but it is a useful starting point. 2 This is typical for a skin friction coefficient. Note that in the simple friction model presented, the turbulent term dominates the laminar when the Reynolds number Re ¼ uh/n >> c/C.
3.1. The base case The base case for comparison with experimental data comprises a grid of 100 100 cells covering the experimental domain to just outside the bund (exceeding it by around 5%), with friction switched off in the model by setting c ¼ C ¼ 0. Sensitivity of the
664
D.M. Webber, M.J. Ivings / Journal of Loss Prevention in the Process Industries 23 (2010) 662e667
predictions to variations away from this base case is examined in conjunction with the comparisons below. As a “standard” friction model we set c ¼ 3.0 and C ¼ 0.001. These values are used as scales by which to vary friction away from the inviscid base case. 3.2. The LJMU releases into cylindrical bunds 3.2.1. Comparison with the inviscid flow model Fig. 2 shows the inviscid predictions of SPLOT with the data on the overtopped liquid fraction. It shows results from the base case 100 100 cell grid and other resolutions ranging from 50 50 to 400 400 cells. In each case the model was run until it was clear that no further bund overtopping would occur. This was generally after about two to three reflected waves had interacted with the bund. Following the experimental set-up, the model was run for one quarter of the cylindrical geometry, applying symmetry boundary conditions along the bounding radii. The marginally better fit with lowest resolution, 50 50 cells, is fortuitous and indicates that care must be taken in ensuring convergence. The results indicate that 100 100 grid is of sufficient resolution to obtain converged results. The three peaks in the overtopping fraction Q plotted on the vertical axis correspond from left to right to high medium and low aspect ratio releases. Trials 1 to 28 are high aspect ratio releases (H/R ¼ 2.0); trials 29 to 56 are medium (H/R ¼ 1.0); trials 57 to 84 are low (H/R ¼ 0.4). For each tank aspect ratio the trials are ordered in increasing hb/H and neither the experimental data nor SPLOT results show the bund overtopping fraction, Q, monotonically decreasing as hb/H increases, though SPLOT comes closer with increasing grid resolution. Looking at the converged results (for grids 100 100 and finer), broad agreement is found for larger values of hb/H for all tank aspect ratios. The trend is well represented, and the slight under-prediction of the overtopping fraction is not significant: the maximum error in Q is approximately 0.1. It is worth noting that SPLOT’s bund model is simple and, whilst one can think of a number of ways in which more liquid might be retained by the bund, it is difficult to imagine how more could get out. We shall return to this in more detail below. For smaller values of hb/H, SPLOT over-predicts the liquid fraction which escapes the bund and most significantly so for the low aspect ratio tanks. A low bund appears to be much better at retaining liquid than SPLOT predicts. As this is an over-prediction, one might ignore it as simply being conservative, but the effect is so clear that it requires further investigation, especially because something as simple (though completely unjustified) as smoothing
out the depth and velocity fields (as is the effect of the coarser resolution) seems to improve the fit so dramatically. 3.2.2. Friction Clearly a possible explanation for the over-prediction of the overtopping fraction seen above is due to the lack of a friction model. In fact it seems a plausible conjecture that friction effects (loss of momentum from the flow due to turbulence or laminar viscosity or both) might slow the front down enough as it encounters the bund to decrease SPLOT’s prediction of the escaped liquid fraction. As discussed earlier, the friction model is an extension of the Shallow Water Equations. Its introduction leads to a degree of freedom in choosing the exact size of the friction terms as long as they remain within physically meaningful limits. The friction term f (equation (1.3)) is modelled as a sum of laminar and turbulent terms, of which in many circumstances (and points in the flow) one or other will dominate. Fig. 3 shows the SPLOT predictions using a range of friction coefficients with all results computed on a grid of 100 100 cells. In the inviscid case the friction coefficients c and C are both set to zero as before. The base case friction coefficient parameters are taken as (c, C) ¼ (3, 0.001) and the results are denoted simply as ‘Friction’ in Fig. 3. The results in Fig. 3 that are labelled as ‘Friction a/b’ multiply the laminar and turbulent friction terms by a and b respectively, i.e. (c, C) ¼ (3.0 a, 0.001 b). For the results computed using the base case friction coefficient values, there is a slight improvement over the inviscid predictions for the low aspect releases (small H/R) where rb/H is large, but no improvement for the high aspect (large H/R) where rb/H is small. A better fit is achieved if the laminar term is multiplied by 10 and the turbulent term by 2, for example. However, there is no justification for such large values, and in any case it does essentially nothing at all to improve the results from the high aspect ratio releases. The inviscid model predicts that if hb/H and rb/R are constant, then the overtopping is also approximately constant. However, the data show that the overtopping fraction is significantly affected by rb/H. Introducing friction terms into the model qualitatively reproduces this physical phenomenon, although in cases where hb/ H is small then SPLOT continues to overpredict Q, but to a lesser extent. It can therefore be concluded that friction is not the reason for the discrepancy in the predictions where hb/H is small. It is possible that the friction coefficients may be better prescribed by fitting less energetic releases if appropriate data become available. Importantly, in cases where the inviscid model under-predicts the overtopping fraction, the introduction of friction does not make
1 0.9
1
0.8
0.9
0.7
0.8
Experiment
0.6 Experiment
Q
0.6
SPLOT 400x400
0.5
Q
0.7
SPLOT inviscid
0.5
Friction
0.4
Friction 10/1
SPLOT 200x200 SPLOT 100x100
0.4
SPLOT 50x50
friction 10/2
0.3
Friction 1/3
0.2
0.3
0.1 0.2
0 0.1
0 80
70
60
50
40
30
20
10
0
0 trial
Fig. 2. SPLOT predictions of the overtopping fraction Q with different grid resolutions.
20
40
60
80
trial
Fig. 3. SPLOT predictions of bund overtopping fraction using different friction coefficients compared against experimental data where ‘Splot inviscid’ uses (c,C) ¼ (0,0), ‘Friction’ uses (c,C) ¼ (3, 0.001) and ‘Friction a/b’ uses (c,C) ¼ (3.0 a, 0.001 b).
D.M. Webber, M.J. Ivings / Journal of Loss Prevention in the Process Industries 23 (2010) 662e667
the predictions significantly worse, even where the over-large friction coefficients are used. The main preliminary conclusion is that SPLOT gives reasonable predictions where hb/H is large (with a discrepancy in the value of Q generally no larger than approximately 0.1), but over-predicts the overtopping (is conservative) where hb/H is small. A more detailed quantitative assessment is made of the SPLOT predictions below, after considering other data sets. 3.3. The Greenspan and Johansson (1981) releases SPLOT was used to simulate the experiments of Greenspan and Johansson (1981) using a computational domain fixed at 0.3 m by 0.3 m. Three computational grids were used with 50 50, 100 100 and 200 200 cells respectively. Again the last two of these grids gave very similar results, indicating a converged solution. The sloping bund was modelled simply as a variation in the bed height, b, leading to non-zero source terms in the momentum equations (1.2) dependent on the bed gradient. The results of the inviscid model and using the base case friction coefficient values are shown in Fig. 4. The data are ordered as follows: the first two peaks are for bunds angled at 30 , the next three for bunds at 60 . Within that, the data are ordered in decreasing bund radius, each of the five “bumps” in Fig. 4 thus corresponding to a given angle and radius of the bund. Within each bump, the data are ordered in decreasing initial liquid level. The agreement is reasonable, with Q generally being over-predicted by between 0.0 and 0.2. 3.4. Scaling To test the scaling properties of the model the LJMU cylindrical runs were scaled down by a factor of approximately 3 to match the radius of the tank in the experiments of Greenspan and Johansson (1981) and also scaled up by a factor of 100 to represent a typical petrochemical storage tank. In the inviscid case, the values of the overtopping fraction obtained barely differed from those obtained in the original runs, confirming expectations based on dimensional arguments, though of course the time taken for the liquid to escape did change. The same tests with the standard friction values of C and c, also yielded results only very marginally different from the originals. This is another confirmation that friction (especially the laminar term) is not important in these cases. Finally simulations of the Greenspan and Johansson (1981) experiments were carried out scaled up to the dimensions of the LJMU experiments. Again the overtopping fraction Q did not differ noticeably from the originals.
1.0
30° bunds
0.9
60° bunds
0.8 0.7
Q
0.6 0.5
Experiment
0.4
SPLOT inviscid SPLOT friction
0.3 0.2 0.1
665
3.5. The interaction of the liquid front with the bund For the LJMU experiments SPLOT over-predicts the overtopping fraction for bunds that are low with respect to the tank fill height. In other words, low bunds seem to be able to retain liquid in reality better than the model predicts. To provide a better understanding of the discrepancies between the model predictions and the experimental data it is useful to look at the assumptions underlying the Shallow Water Equations e in particular the assumption that vertical accelerations are small. The Shallow Water Equations include a depth which depends on time, but paradoxically they ignore all vertical liquid motion. It is instructive to look at the energy balance in the Shallow Water equations. The derived kinetic, T, potential, P, and total energy, E, (ignoring the constant density) are
T ¼
1 hðu$uÞ 2
P ¼ ghðb þ h=2Þ
E ¼ T þP
(3.1)
And from the Shallow Water Equations (1.1, 1.2) it follows that
vT v Tuj vðh þ bÞ þ ghuj ¼ uj fj þ vxj vt vxj 2 vP v Puj vðh þ bÞ 1 v uj h ghuj ¼ þ vxj vxj vt vxj 2
(3.2)
(3.3)
and hence
2 vE v Euj 1 v uj h ¼ uj f j þ vt 2 vxj vxj
(3.4)
In equation (3.4) the fact that f is always in the opposite direction to u guarantees that the friction terms dissipate energy. The final term on the right is a pure divergence, which therefore redistributes energy within the pool but doesn’t lose or gain any. The g-terms in (3.2) and (3.3) convert potential energy to kinetic where the top surface (at height h þ b) is not level, but again do not lose any energy. The above describes a consistent energy balance model, but does not incorporate kinetic energy associated with vertical motion. Now consider what happens when the front arrives at a bund wall. What we would like to see happen, so that the assumptions made in the derivation of the Shallow Water Equations remain valid, is something gradual, along the lines of Fig. 5a. For a sufficiently slowly moving front, this may not be far from the truth. However, the reality for a rapidly moving front more closely resembles Fig. 5b, where the front is redirected violently upwards. The important aspect of the flow here is not the large value of the depth gradient jVhj, but the large value of vh=vt.3 In some of the model calculations the depth profile looks very much like Fig. 5b, capturing the rapid increase in the liquid depth, despite the fact that it has no vertical component of velocity. Of course Shallow Water Theory is only valid in the limit of negligible vertical acceleration and so, as we have seen in the energy balance, it cannot be expected to handle situations in which there is a lot of kinetic energy associated with the vertical motion. It is therefore necessary to understand the limits imposed by the adoption of Shallow Water Theory, and whether (as we did by incorporating friction terms) it is possible to extend them. If the front arrives at the bund with velocity u and the bund height is hb,
0.0 0
5
10
15
20
25
30
Trial
Fig. 4. SPLOT simulations of the Greenspan and Johansson experiments with and without friction.
3 Note that the initial conditions for the dam break problem (i.e. the instantaneous release of a column of liquid) involve a large value of jVhj, but the Shallow Water Equations have been shown to approximate this flow very well.
666
D.M. Webber, M.J. Ivings / Journal of Loss Prevention in the Process Industries 23 (2010) 662e667
Fig. 5. (A) A schematic sequence showing the front slowly approaching a bund wall in keeping with the shallow water approximation. The liquid is reflected from the wall but builds up slowly until it flows over the top. The horizontal arrows emphasise that the flow is always considered horizontal in Shallow Water Theory. (B) A schematic sequence showing what may happen in practice when a fast liquid front hits a bund. If the front is moving rapidly when it encounters the bund, then a strongly vertical flow will ensue.
then simple energy considerations tell us that if the front is directed sharply upwards and u2/2 > ghb then the front may be deflected to a height higher than that of the bund. We therefore define a dimensionless “kinetic factor”
k ¼ u2 =ð2ghb Þ
(3.5)
and adopt the hypothesis that Shallow Water Theory will work well for k < O(1) but start to go wrong if k >> 1. In the case of a concentric cylindrical bund, k can be estimated based on the difference between the initial potential energy of a cylinder of liquid with volume V and radius R and the kinetic energy when it encounters the bund with radius rb. Thus
gV V gV V V ¼ þ u2 2 pR2 2 pr 2 4 b
(3.6)
where the ¼ in the kinetic term on the right arises from assuming a linear variation of velocity from the centre of the pool to the edge. Hence
2 ! R u ¼ 2gH 1 rb 2
(3.7)
and 1 0.9 0.8 0.7
Q
0.6 Experiment
0.5
SPLOT inviscid
0.4
k/20
0.3 0.2 0.1 0 0
20
40
60
80
trial
Fig. 6. SPLOT predictions (neglecting friction) and LJMU measurements of overtopping fraction compared with the kinetic factor (divided arbitrarily by 20 to fit it on the same vertical axis).
2 ! H R 1 kh hb rb
(3.8)
In Fig. 6 we show this value of k for the various trials alongside the data and inviscid predictions using SPLOT. The agreement between the SPLOT predictions and the experimental data is reasonable in all cases with k < 10 (k/20 < 0.5) but less so when k > 10. It follows that to improve the fit to this data it would be desirable to have a better bund model, in which a low bund can hold back a fast moving liquid front by a ‘splash-back’ mechanism. Such a sub-model is suggested by the fact that it is precisely in this region where the assumptions underlying Shallow Water Theory are most strongly violated. It is now also clearer why too coarse a grid resolution may give an improved fit. The coarse grid smears out velocity and depth gradients, taking the edge off the violence with which the model front encounters the bund, partly, at least, mitigating the effect of large k. However these improvements must be assumed to apply only in very limited circumstances, and selecting a coarse grid cannot be taken to be a substitute for a better bund model. 4. Conclusions The model SPLOT provides robust numerical solutions of the Shallow Water Equations, enhanced by friction terms, and allows the introduction of vertical bund walls. It has been compared here with a variety of bund overtopping experiments. In every experimental trial the bund was large enough to retain the entire volume of liquid released, had that liquid been released slowly and gradually. In most cases, however, the releases were sufficiently sudden that some nonzero fraction Q of the released volume overtopped the bund. Various sensitivity tests were undertaken examining the key modelling assumptions including the grid resolution and friction model. For the experiments considered a grid domain of 100 100 cells provided good evidence of convergence. This allowed simulations to be carried out in a minute or so on a modest personal computer, whereas 400 400 cells takes noticeably longer. For many releases the inviscid predictions of the overtopping fraction were very good, especially in view of the fact that the inviscid model contains no free parameters. Introducing friction made little overall difference to the predictions, although significant improvements were made to a small number of cases where SPLOT had been significantly over-predicting the overtopping fraction. Encouragingly, the introduction of friction
D.M. Webber, M.J. Ivings / Journal of Loss Prevention in the Process Industries 23 (2010) 662e667
does not make the best inviscid predictions any worse. Setting friction coefficients much larger than can be justified did improve the fit slightly, but not in all cases. Thus while friction does introduce two parameters which we have some freedom to adjust, the results are largely insensitive to them, and it appears that reasonable sized friction terms are neither going to make the best predictions worse, nor the worst ones much better. Therefore friction does not appear to be a particularly important factor in these releases. This is confirmed by the scaling properties of the model, where not only the inviscid results but also those with friction gave the same values of the overtopping fraction, independent of quite a range of scaling factors applied to all the linear dimensions of the problem. The “kinetic factor”, defined in equation (3.5), has been shown to distinguish between the cases where SPLOT is in reasonable agreement with the experimental data and where it over-predicts the overtopping fraction considerably. For k < 10 a reasonably good fit was found whereas for k > 10 the fit was poorer. We therefore conclude that the reasons for very conservative predictions in some cases are associated with the model of the flow as a strong front hits the bund wall. Furthermore, it is in exactly in these cases (releases with large values of k), that standard Shallow Water Theory is inappropriate in the region close to the bund wall. It may therefore be possible to produce a better bund model which will account for some ‘splash-back’ in cases of large k, where SPLOT currently overestimates the overtopping, though that is outside the scope of this work. We conclude with the recommendation that SPLOT be used for estimating hazards of bund overtopping, but that its results should be accepted as conservative in cases where the kinetic factor k is estimated to be greater than around 10. Acknowledgements The Authors would like to thank the UK Health and Safety Executive (HSE) for supporting this work.
HSE’s motivation in procuring the LJMU experiments was twofold. In the short term the test results were seen as providing a check on the validity of a simple correlation, based on the small-scale Greenspan and Young (1978) tests that HSE was then using to predict bund overtopping at large scale. In the longer term they were seen as providing a sound base from which the development of more advanced predictive methods, such as SPLOT, could be guided. In the event, HSE learned that their simple correlation performed only moderately well at the scale of the LJMU tests, and they derived a new correlation giving a better fit. The new correlation is
Q ¼ 1:0255 0:1886* rHb 2:9951* hHb þ 0:3842* HR 2 2 2 þ 0:0140* rHb þ 2:7535* hHb 0:0637* HR 3 3 0:0005* rHb 0:8595* hHb :
The new correlation, known as OVERTOP, which HSE now uses routinely for its predictive work, reproduces the overtopping fraction measured in each of the eighty-four LJMU tests with a rootmean-square error of 0.036 and a maximum error of 0.10 in its prediction of the overtopping fraction Q.
References Atherton, W. (2005). An experimental investigation of bund wall overtopping and dynamic pressures on the bund wall following catastrophic failure of a storage vessel. Health and Safety Executive Research Report 333. HMSO, ISBN 0 7176 2988 0. Greenspan, H. P., & Johansson, A. V. (1981). An experimental study of flow over an impounding dike. Studies in Applied Mathematics, 64, 211e223. Greenspan, H. P., & Young, R. E. (1978). Flow over a containment dyke. Journal of Fluid Mechanics, 87(1), 179e192. Hirst, I. (2006). Selected reconstructed data from the experiments of Greenspan & Johansson (1981). Private Communication. Huppert, H. E. (1982). The propagation of two-dimensional and axisymmetric viscous gravity currents over a rigid horizontal surface. Journal of Fluid Mechanics,121, 43e58. Ivings, M. J., Toro, E. F., & Webber, D. M. (2003). Numerical schemes for the 2D shallow layer equations including dry fronts. Computational Fluid Dynamics Journal, 12(1), 41e52. Ivings, M. J. & Webber, D. M. (2002). Development of SPLOT: a numerical model for the two dimensional shallow layer equations. Health and Safety Laboratory Report CM/02/08. Ivings, M. J., & Webber, D. M. (2007). Modelling bund overtopping using a shallow water CFD model. Journal of Loss Prevention in the Process Industries, 20(1), 38e44. Landau, L. D., & Lifshitz, E. M. (1959). Fluid mechanics. Pergamon, ISBN 0 08 009104 0. Thyer, A. M., Hirst, I. L., & Jagger, S. F. (2002). Bund overtopping d the consequence of catastrophic tank failure. Journal of Loss Prevention in the Process Industries, 15 (5), 357e363. Toro, E. F. (2001). Shock-capturing methods for free-surface shallow flows. Wiley. Webber, D. M., Ivings, M. J., Maru, W-A., & Thyer, A. M. (2009). Validation of the Shallow Water model “SPLOT” against experimental data on bund overtopping. Health and Safety Executive Research Report RR755. Webber, D. M., & Jones, S. J. (1987). A model of spreading vaporising pools. In J. Woodward (Ed.), International conference on Vapor Cloud modeling, Boston Massachusetts. USA (November 1987) proceedings. AIChE. Webber, D. M., & Brighton, P. W. M. (1986). Inviscid similarity solutions for slumping from a cylindrical tank. Journal of Fluids Engineering, 108, 238e240.
Nomenclature
A. Appendix
ðA:1Þ
667
b: bed height c, C: laminar and turbulent friction coefficients E: total energy fi: friction components f: friction vector g: gravitational constant (9.81 ms2) h: liquid depth hb: bund height H: tank fill height k: ‘kinetic factor’ P: potential energy Q: bund overtopping fraction rb: bund radius R: tank radius t: time T: kinetic energy ui: velocity components u: velocity vector V: tank volume vb: bund volume xi: Cartesian co-ordinate directions m: kinematic viscosity
Contents lists available at ScienceDirect
Journal of Loss Prevention in the Process Industries journal homepage: www.elsevier.com/locate/jlp
Modelling bund overtopping using Shallow Water Theory D.M. Webber, M.J. Ivings* Health and Safety Laboratory, Harpur Hill, Buxton SK17 9JN, UK
a r t i c l e i n f o
a b s t r a c t
Article history: Received 4 February 2010 Received in revised form 7 July 2010 Accepted 7 July 2010
A model is presented of liquid flow under gravity, over horizontal and sloping surfaces, based on the Shallow Water Equations. The aim is to facilitate assessment of the efficacy of bunds surrounding liquid storage tanks, in cases where the liquid front encountering the bund carries significant momentum. The model, called SPLOT, incorporates terms representing the effects of turbulent and laminar friction, and a sub-model for the effect of vertical walls. Having produced a robust numerical solution method, and verified its performance in cases where analytic solutions are known, it is the objective of this work to test the model’s validity against experimental data on flows overtopping bunds. SPLOT is found to be very satisfactory in a wide range of circumstances involving both vertical and sloping bunds. Physical reasons are identified for the overly conservative predictions found in the remaining cases, and a simple parameter is derived which indicates whether the model is being used outside its range of optimum validity. Crown Copyright Ó 2010 Published by Elsevier Ltd. All rights reserved.
Keywords: Shallow water equations CFD modelling Bund overtopping Validation
1. Introduction Assessment of the safety of liquid storage tanks requires an estimation of how much liquid may overtop a surrounding bund in the event of a sudden and significant failure of a tank (Thyer, Hirst, & Jagger, 2002). In such cases it is not sufficient that the bund may simply have the capacity to hold the entire contents of the tank, as a rapidly advancing liquid front may not be entirely held back by the bund wall. In general therefore, it is necessary to understand the dynamics of an advancing liquid front, and what happens when it encounters a wall. This paper summarises a mathematical model of this process, and presents a comparison with some well defined laboratory-scale experiments using water, which measure the quantity of liquid escaping from the bund in a number of different geometries and at different scales. In general, a bunded area may be of arbitrary shape and surround one or more tanks. A starting point for understanding the flow is therefore provided by the two-dimensional Shallow Water Equations for the flow of shallow pools in two horizontal dimensions e see for example Landau and Lifshitz (1959), Webber and Brighton (1986) and Toro (2001). However, two specific limitations of the standard Shallow Water Equations are that they describe inviscid flow, and they are only strictly valid when vertical accelerations can be neglected. These are addressed as follows. A model of the effect of friction (in laminar and turbulent flow) is introduced at the outset and the importance of vertical * Corresponding author. Tel.: þ44 (0) 1298 218133; fax: þ44 (0) 1298 218840. E-mail address: [email protected] (M.J. Ivings).
acceleration is discussed in some detail after the model has been compared with the experimental data, where it will be seen to provide limits on where the model can be expected to perform optimally. The model is defined by the equations:
vh v huj ¼ 0 þ vxj vt
(1.1)
vðhui Þ v hui uj vh vb þ gh ¼ fi gh þ vxj vt vxi vxi
(1.2)
with
cn fi ¼ ui Cjujui h
(1.3)
where h(x, t) is the depth of liquid, u(x, t) is the 2-component horizontal liquid velocity, b(x) is the ground level (but note that only its gradient is important), f(x, t) is the friction term, n is the kinematic viscosity of the liquid (w1.0 106 m2/s for water), c is a laminar friction coefficient, C is a turbulent skin friction coefficient and (t, x1, x2) are time and the two horizontal space coordinates Note that where friction is negligible, the model has no free parameters. The friction model introduces two free parameters c and C which have been discussed in detail by Webber and Jones (1987) who point out that the parameter c, governing friction in laminar flow should take the value 3.0 in order to reproduce the lubrication theory results of Huppert (1982) who studied thin
0950-4230/$ e see front matter Crown Copyright Ó 2010 Published by Elsevier Ltd. All rights reserved. doi:10.1016/j.jlp.2010.07.002
D.M. Webber, M.J. Ivings / Journal of Loss Prevention in the Process Industries 23 (2010) 662e667
liquid films between parallel plates1. For turbulent flow, parameter C is the more important one and is2 of the order of 103. This model, known as “SPLOT”, has been discussed by Ivings and Webber (2002), Ivings, Toro, and Webber (2003) and Ivings and Webber (2007) highlighting some of the problems which must be overcome to produce a numerical solution algorithm which is robust enough to be used routinely in the analysis of industrial hazards, and which provides accurate numerical solutions in cases where analytic properties of the solution are known. The final step, necessary to enable confident use of the model, is to test its validity against experimental data for dynamic bund overtopping. This is the main objective of this work but, in comparing the model with experimental data, we have also investigated its sensitivity both to grid resolution and to the size of the friction coefficients. We have also performed a verification of the inviscid model (c ¼ C ¼ 0) against the expectation from dimensional arguments that the fraction of liquid overtopping the bund should not change if all the length dimensions of a scenario are scaled by the same factor. The experiments considered here all involved complete failure of a cylindrical tank with a concentric cylindrical bund. The model SPLOT is not constrained to this geometry, but it does produce a very useful test case which, as we shall show, allows analysis of the reasons underlying the way in which the model performs. It should also be emphasised that the current work is aimed at validation of the spreading model and the model of the interaction between the spreading liquid and the bund. Vaporising pools are outside the scope of this paper. 2. Experimental data Two series of experiments have been used here to validate SPLOT: those reported by Greenspan and Johansson (1981) and those commissioned by the UK Health and Safety Executive (HSE) (see Appendix) at Liverpool John Moores University (LJMU) and reported by Atherton (2005). Of these the latter are at a somewhat larger scale and use vertical bund walls, and so we shall address them first. The Greenspan and Johansson (1981) experiments used bunds sloping away from the containment area and require a different modelling approach by SPLOT. A detailed (and of necessity somewhat lengthy) tabulated summary of the experimental parameter space (tank and bund geometries), for both series of experiments, is given by Webber, Ivings, Maru, and Thyer (2009). 2.1. The LJMU experiments Of the LJMU experiments (Atherton, 2005), 84 trials were carried out with instantaneous releases of water from a cylindrical tank with a concentric cylindrical bund at some radial distance from it. The scenario addressed is shown in Fig. 1. In order to maximise the scale of the experiments in the laboratory space available, advantage was taken of the expected cylindrical symmetry of the flow, and the whole experiment was conducted between walls meeting at a right angle at the centre of the tank; thus the whole experiment was effectively a quarter of the modelled cylinder.
663
H hb
R
rb
Fig. 1. The geometric parameters in experiments involving a cylindrical tank surrounded by a concentric circular bund.
The salient parameters are the radius of the tank, R, the depth of the water in the tank, H, the radial position of the bund, rb, and the height of the bund, hb. The volume of the liquid released (for the full cylinder), V, and the capacity of the bund, vb, are then given by
V ¼ pR2 H and vb ¼ prb2 hb
(2.1)
In each experiment it was arranged that vb > V so that the bund would be sufficient to retain all of the liquid if it were released slowly. All measured overtopping is therefore a dynamic effect owing to the rapid release of liquid from the tank. The primary point of comparison of the model with the data is the fraction of the released liquid which overtops the bund, Q. The experimental geometry gives three independent dimensionless ratios of the parameters R, H, rb and hb. As long as the scale of the experiments is large enough, it is reasonable to expect the (dimensionless) overtopping fraction Q to be independent of the scale, and therefore populate a three-dimensional domain of results defined, for example, by the values of rb/R, hb/H, H/R. Other possible dimensionless combinations may have more physical significance. For example (when vb < V) 1 vb/V provides an obvious lower bound on the value of Q. We discuss other physically significant dimensionless combinations below. The experiments were performed with R ¼ 0.3 m using instantaneous high, medium and low aspect releases with H ¼ 0.6 m, 0.3 m and 0.12 m respectively. Bund heights hb were in the range 0.03 m e 0.72 m (the last representing a high collar bund) and at radial distances rb in the range 0.315 m e 1.897 m but always such that the bund could in principle contain the entire volume of the tank e i.e. vb > V. 2.2. The Greenspan and Johansson (1981) experiments Greenspan and Johansson (1981) also conducted instantaneous releases of water from a cylindrical tank into a concentric cylindrical bund, but with the bund consisting of an upward slope away from the tank at angles of 30 and 60 to the horizontal. In these experiments the linear scale was about a factor three smaller than those done by Atherton (2005) with tank radius R ¼ 0.0953 m and the fill level H (extracted indirectly from the data by Hirst, 2006) in the range 0.07 me0.26 m. A particular point of interest in these experiments is that the sloping bund may be expected to have different retention properties from those of a vertical wall. 3. Comparison of model and data
1
Huppert was concerned with a thin liquid film between parallel plates, and c ¼ 3 is what is required to match Huppert’s NaviereStokes solution using the above laminar friction term added to the inviscid Shallow Water Equations. It is not entirely obvious that the same value of c will be optimal for a front advancing on to dry ground, but it is a useful starting point. 2 This is typical for a skin friction coefficient. Note that in the simple friction model presented, the turbulent term dominates the laminar when the Reynolds number Re ¼ uh/n >> c/C.
3.1. The base case The base case for comparison with experimental data comprises a grid of 100 100 cells covering the experimental domain to just outside the bund (exceeding it by around 5%), with friction switched off in the model by setting c ¼ C ¼ 0. Sensitivity of the
664
D.M. Webber, M.J. Ivings / Journal of Loss Prevention in the Process Industries 23 (2010) 662e667
predictions to variations away from this base case is examined in conjunction with the comparisons below. As a “standard” friction model we set c ¼ 3.0 and C ¼ 0.001. These values are used as scales by which to vary friction away from the inviscid base case. 3.2. The LJMU releases into cylindrical bunds 3.2.1. Comparison with the inviscid flow model Fig. 2 shows the inviscid predictions of SPLOT with the data on the overtopped liquid fraction. It shows results from the base case 100 100 cell grid and other resolutions ranging from 50 50 to 400 400 cells. In each case the model was run until it was clear that no further bund overtopping would occur. This was generally after about two to three reflected waves had interacted with the bund. Following the experimental set-up, the model was run for one quarter of the cylindrical geometry, applying symmetry boundary conditions along the bounding radii. The marginally better fit with lowest resolution, 50 50 cells, is fortuitous and indicates that care must be taken in ensuring convergence. The results indicate that 100 100 grid is of sufficient resolution to obtain converged results. The three peaks in the overtopping fraction Q plotted on the vertical axis correspond from left to right to high medium and low aspect ratio releases. Trials 1 to 28 are high aspect ratio releases (H/R ¼ 2.0); trials 29 to 56 are medium (H/R ¼ 1.0); trials 57 to 84 are low (H/R ¼ 0.4). For each tank aspect ratio the trials are ordered in increasing hb/H and neither the experimental data nor SPLOT results show the bund overtopping fraction, Q, monotonically decreasing as hb/H increases, though SPLOT comes closer with increasing grid resolution. Looking at the converged results (for grids 100 100 and finer), broad agreement is found for larger values of hb/H for all tank aspect ratios. The trend is well represented, and the slight under-prediction of the overtopping fraction is not significant: the maximum error in Q is approximately 0.1. It is worth noting that SPLOT’s bund model is simple and, whilst one can think of a number of ways in which more liquid might be retained by the bund, it is difficult to imagine how more could get out. We shall return to this in more detail below. For smaller values of hb/H, SPLOT over-predicts the liquid fraction which escapes the bund and most significantly so for the low aspect ratio tanks. A low bund appears to be much better at retaining liquid than SPLOT predicts. As this is an over-prediction, one might ignore it as simply being conservative, but the effect is so clear that it requires further investigation, especially because something as simple (though completely unjustified) as smoothing
out the depth and velocity fields (as is the effect of the coarser resolution) seems to improve the fit so dramatically. 3.2.2. Friction Clearly a possible explanation for the over-prediction of the overtopping fraction seen above is due to the lack of a friction model. In fact it seems a plausible conjecture that friction effects (loss of momentum from the flow due to turbulence or laminar viscosity or both) might slow the front down enough as it encounters the bund to decrease SPLOT’s prediction of the escaped liquid fraction. As discussed earlier, the friction model is an extension of the Shallow Water Equations. Its introduction leads to a degree of freedom in choosing the exact size of the friction terms as long as they remain within physically meaningful limits. The friction term f (equation (1.3)) is modelled as a sum of laminar and turbulent terms, of which in many circumstances (and points in the flow) one or other will dominate. Fig. 3 shows the SPLOT predictions using a range of friction coefficients with all results computed on a grid of 100 100 cells. In the inviscid case the friction coefficients c and C are both set to zero as before. The base case friction coefficient parameters are taken as (c, C) ¼ (3, 0.001) and the results are denoted simply as ‘Friction’ in Fig. 3. The results in Fig. 3 that are labelled as ‘Friction a/b’ multiply the laminar and turbulent friction terms by a and b respectively, i.e. (c, C) ¼ (3.0 a, 0.001 b). For the results computed using the base case friction coefficient values, there is a slight improvement over the inviscid predictions for the low aspect releases (small H/R) where rb/H is large, but no improvement for the high aspect (large H/R) where rb/H is small. A better fit is achieved if the laminar term is multiplied by 10 and the turbulent term by 2, for example. However, there is no justification for such large values, and in any case it does essentially nothing at all to improve the results from the high aspect ratio releases. The inviscid model predicts that if hb/H and rb/R are constant, then the overtopping is also approximately constant. However, the data show that the overtopping fraction is significantly affected by rb/H. Introducing friction terms into the model qualitatively reproduces this physical phenomenon, although in cases where hb/ H is small then SPLOT continues to overpredict Q, but to a lesser extent. It can therefore be concluded that friction is not the reason for the discrepancy in the predictions where hb/H is small. It is possible that the friction coefficients may be better prescribed by fitting less energetic releases if appropriate data become available. Importantly, in cases where the inviscid model under-predicts the overtopping fraction, the introduction of friction does not make
1 0.9
1
0.8
0.9
0.7
0.8
Experiment
0.6 Experiment
Q
0.6
SPLOT 400x400
0.5
Q
0.7
SPLOT inviscid
0.5
Friction
0.4
Friction 10/1
SPLOT 200x200 SPLOT 100x100
0.4
SPLOT 50x50
friction 10/2
0.3
Friction 1/3
0.2
0.3
0.1 0.2
0 0.1
0 80
70
60
50
40
30
20
10
0
0 trial
Fig. 2. SPLOT predictions of the overtopping fraction Q with different grid resolutions.
20
40
60
80
trial
Fig. 3. SPLOT predictions of bund overtopping fraction using different friction coefficients compared against experimental data where ‘Splot inviscid’ uses (c,C) ¼ (0,0), ‘Friction’ uses (c,C) ¼ (3, 0.001) and ‘Friction a/b’ uses (c,C) ¼ (3.0 a, 0.001 b).
D.M. Webber, M.J. Ivings / Journal of Loss Prevention in the Process Industries 23 (2010) 662e667
the predictions significantly worse, even where the over-large friction coefficients are used. The main preliminary conclusion is that SPLOT gives reasonable predictions where hb/H is large (with a discrepancy in the value of Q generally no larger than approximately 0.1), but over-predicts the overtopping (is conservative) where hb/H is small. A more detailed quantitative assessment is made of the SPLOT predictions below, after considering other data sets. 3.3. The Greenspan and Johansson (1981) releases SPLOT was used to simulate the experiments of Greenspan and Johansson (1981) using a computational domain fixed at 0.3 m by 0.3 m. Three computational grids were used with 50 50, 100 100 and 200 200 cells respectively. Again the last two of these grids gave very similar results, indicating a converged solution. The sloping bund was modelled simply as a variation in the bed height, b, leading to non-zero source terms in the momentum equations (1.2) dependent on the bed gradient. The results of the inviscid model and using the base case friction coefficient values are shown in Fig. 4. The data are ordered as follows: the first two peaks are for bunds angled at 30 , the next three for bunds at 60 . Within that, the data are ordered in decreasing bund radius, each of the five “bumps” in Fig. 4 thus corresponding to a given angle and radius of the bund. Within each bump, the data are ordered in decreasing initial liquid level. The agreement is reasonable, with Q generally being over-predicted by between 0.0 and 0.2. 3.4. Scaling To test the scaling properties of the model the LJMU cylindrical runs were scaled down by a factor of approximately 3 to match the radius of the tank in the experiments of Greenspan and Johansson (1981) and also scaled up by a factor of 100 to represent a typical petrochemical storage tank. In the inviscid case, the values of the overtopping fraction obtained barely differed from those obtained in the original runs, confirming expectations based on dimensional arguments, though of course the time taken for the liquid to escape did change. The same tests with the standard friction values of C and c, also yielded results only very marginally different from the originals. This is another confirmation that friction (especially the laminar term) is not important in these cases. Finally simulations of the Greenspan and Johansson (1981) experiments were carried out scaled up to the dimensions of the LJMU experiments. Again the overtopping fraction Q did not differ noticeably from the originals.
1.0
30° bunds
0.9
60° bunds
0.8 0.7
Q
0.6 0.5
Experiment
0.4
SPLOT inviscid SPLOT friction
0.3 0.2 0.1
665
3.5. The interaction of the liquid front with the bund For the LJMU experiments SPLOT over-predicts the overtopping fraction for bunds that are low with respect to the tank fill height. In other words, low bunds seem to be able to retain liquid in reality better than the model predicts. To provide a better understanding of the discrepancies between the model predictions and the experimental data it is useful to look at the assumptions underlying the Shallow Water Equations e in particular the assumption that vertical accelerations are small. The Shallow Water Equations include a depth which depends on time, but paradoxically they ignore all vertical liquid motion. It is instructive to look at the energy balance in the Shallow Water equations. The derived kinetic, T, potential, P, and total energy, E, (ignoring the constant density) are
T ¼
1 hðu$uÞ 2
P ¼ ghðb þ h=2Þ
E ¼ T þP
(3.1)
And from the Shallow Water Equations (1.1, 1.2) it follows that
vT v Tuj vðh þ bÞ þ ghuj ¼ uj fj þ vxj vt vxj 2 vP v Puj vðh þ bÞ 1 v uj h ghuj ¼ þ vxj vxj vt vxj 2
(3.2)
(3.3)
and hence
2 vE v Euj 1 v uj h ¼ uj f j þ vt 2 vxj vxj
(3.4)
In equation (3.4) the fact that f is always in the opposite direction to u guarantees that the friction terms dissipate energy. The final term on the right is a pure divergence, which therefore redistributes energy within the pool but doesn’t lose or gain any. The g-terms in (3.2) and (3.3) convert potential energy to kinetic where the top surface (at height h þ b) is not level, but again do not lose any energy. The above describes a consistent energy balance model, but does not incorporate kinetic energy associated with vertical motion. Now consider what happens when the front arrives at a bund wall. What we would like to see happen, so that the assumptions made in the derivation of the Shallow Water Equations remain valid, is something gradual, along the lines of Fig. 5a. For a sufficiently slowly moving front, this may not be far from the truth. However, the reality for a rapidly moving front more closely resembles Fig. 5b, where the front is redirected violently upwards. The important aspect of the flow here is not the large value of the depth gradient jVhj, but the large value of vh=vt.3 In some of the model calculations the depth profile looks very much like Fig. 5b, capturing the rapid increase in the liquid depth, despite the fact that it has no vertical component of velocity. Of course Shallow Water Theory is only valid in the limit of negligible vertical acceleration and so, as we have seen in the energy balance, it cannot be expected to handle situations in which there is a lot of kinetic energy associated with the vertical motion. It is therefore necessary to understand the limits imposed by the adoption of Shallow Water Theory, and whether (as we did by incorporating friction terms) it is possible to extend them. If the front arrives at the bund with velocity u and the bund height is hb,
0.0 0
5
10
15
20
25
30
Trial
Fig. 4. SPLOT simulations of the Greenspan and Johansson experiments with and without friction.
3 Note that the initial conditions for the dam break problem (i.e. the instantaneous release of a column of liquid) involve a large value of jVhj, but the Shallow Water Equations have been shown to approximate this flow very well.
666
D.M. Webber, M.J. Ivings / Journal of Loss Prevention in the Process Industries 23 (2010) 662e667
Fig. 5. (A) A schematic sequence showing the front slowly approaching a bund wall in keeping with the shallow water approximation. The liquid is reflected from the wall but builds up slowly until it flows over the top. The horizontal arrows emphasise that the flow is always considered horizontal in Shallow Water Theory. (B) A schematic sequence showing what may happen in practice when a fast liquid front hits a bund. If the front is moving rapidly when it encounters the bund, then a strongly vertical flow will ensue.
then simple energy considerations tell us that if the front is directed sharply upwards and u2/2 > ghb then the front may be deflected to a height higher than that of the bund. We therefore define a dimensionless “kinetic factor”
k ¼ u2 =ð2ghb Þ
(3.5)
and adopt the hypothesis that Shallow Water Theory will work well for k < O(1) but start to go wrong if k >> 1. In the case of a concentric cylindrical bund, k can be estimated based on the difference between the initial potential energy of a cylinder of liquid with volume V and radius R and the kinetic energy when it encounters the bund with radius rb. Thus
gV V gV V V ¼ þ u2 2 pR2 2 pr 2 4 b
(3.6)
where the ¼ in the kinetic term on the right arises from assuming a linear variation of velocity from the centre of the pool to the edge. Hence
2 ! R u ¼ 2gH 1 rb 2
(3.7)
and 1 0.9 0.8 0.7
Q
0.6 Experiment
0.5
SPLOT inviscid
0.4
k/20
0.3 0.2 0.1 0 0
20
40
60
80
trial
Fig. 6. SPLOT predictions (neglecting friction) and LJMU measurements of overtopping fraction compared with the kinetic factor (divided arbitrarily by 20 to fit it on the same vertical axis).
2 ! H R 1 kh hb rb
(3.8)
In Fig. 6 we show this value of k for the various trials alongside the data and inviscid predictions using SPLOT. The agreement between the SPLOT predictions and the experimental data is reasonable in all cases with k < 10 (k/20 < 0.5) but less so when k > 10. It follows that to improve the fit to this data it would be desirable to have a better bund model, in which a low bund can hold back a fast moving liquid front by a ‘splash-back’ mechanism. Such a sub-model is suggested by the fact that it is precisely in this region where the assumptions underlying Shallow Water Theory are most strongly violated. It is now also clearer why too coarse a grid resolution may give an improved fit. The coarse grid smears out velocity and depth gradients, taking the edge off the violence with which the model front encounters the bund, partly, at least, mitigating the effect of large k. However these improvements must be assumed to apply only in very limited circumstances, and selecting a coarse grid cannot be taken to be a substitute for a better bund model. 4. Conclusions The model SPLOT provides robust numerical solutions of the Shallow Water Equations, enhanced by friction terms, and allows the introduction of vertical bund walls. It has been compared here with a variety of bund overtopping experiments. In every experimental trial the bund was large enough to retain the entire volume of liquid released, had that liquid been released slowly and gradually. In most cases, however, the releases were sufficiently sudden that some nonzero fraction Q of the released volume overtopped the bund. Various sensitivity tests were undertaken examining the key modelling assumptions including the grid resolution and friction model. For the experiments considered a grid domain of 100 100 cells provided good evidence of convergence. This allowed simulations to be carried out in a minute or so on a modest personal computer, whereas 400 400 cells takes noticeably longer. For many releases the inviscid predictions of the overtopping fraction were very good, especially in view of the fact that the inviscid model contains no free parameters. Introducing friction made little overall difference to the predictions, although significant improvements were made to a small number of cases where SPLOT had been significantly over-predicting the overtopping fraction. Encouragingly, the introduction of friction
D.M. Webber, M.J. Ivings / Journal of Loss Prevention in the Process Industries 23 (2010) 662e667
does not make the best inviscid predictions any worse. Setting friction coefficients much larger than can be justified did improve the fit slightly, but not in all cases. Thus while friction does introduce two parameters which we have some freedom to adjust, the results are largely insensitive to them, and it appears that reasonable sized friction terms are neither going to make the best predictions worse, nor the worst ones much better. Therefore friction does not appear to be a particularly important factor in these releases. This is confirmed by the scaling properties of the model, where not only the inviscid results but also those with friction gave the same values of the overtopping fraction, independent of quite a range of scaling factors applied to all the linear dimensions of the problem. The “kinetic factor”, defined in equation (3.5), has been shown to distinguish between the cases where SPLOT is in reasonable agreement with the experimental data and where it over-predicts the overtopping fraction considerably. For k < 10 a reasonably good fit was found whereas for k > 10 the fit was poorer. We therefore conclude that the reasons for very conservative predictions in some cases are associated with the model of the flow as a strong front hits the bund wall. Furthermore, it is in exactly in these cases (releases with large values of k), that standard Shallow Water Theory is inappropriate in the region close to the bund wall. It may therefore be possible to produce a better bund model which will account for some ‘splash-back’ in cases of large k, where SPLOT currently overestimates the overtopping, though that is outside the scope of this work. We conclude with the recommendation that SPLOT be used for estimating hazards of bund overtopping, but that its results should be accepted as conservative in cases where the kinetic factor k is estimated to be greater than around 10. Acknowledgements The Authors would like to thank the UK Health and Safety Executive (HSE) for supporting this work.
HSE’s motivation in procuring the LJMU experiments was twofold. In the short term the test results were seen as providing a check on the validity of a simple correlation, based on the small-scale Greenspan and Young (1978) tests that HSE was then using to predict bund overtopping at large scale. In the longer term they were seen as providing a sound base from which the development of more advanced predictive methods, such as SPLOT, could be guided. In the event, HSE learned that their simple correlation performed only moderately well at the scale of the LJMU tests, and they derived a new correlation giving a better fit. The new correlation is
Q ¼ 1:0255 0:1886* rHb 2:9951* hHb þ 0:3842* HR 2 2 2 þ 0:0140* rHb þ 2:7535* hHb 0:0637* HR 3 3 0:0005* rHb 0:8595* hHb :
The new correlation, known as OVERTOP, which HSE now uses routinely for its predictive work, reproduces the overtopping fraction measured in each of the eighty-four LJMU tests with a rootmean-square error of 0.036 and a maximum error of 0.10 in its prediction of the overtopping fraction Q.
References Atherton, W. (2005). An experimental investigation of bund wall overtopping and dynamic pressures on the bund wall following catastrophic failure of a storage vessel. Health and Safety Executive Research Report 333. HMSO, ISBN 0 7176 2988 0. Greenspan, H. P., & Johansson, A. V. (1981). An experimental study of flow over an impounding dike. Studies in Applied Mathematics, 64, 211e223. Greenspan, H. P., & Young, R. E. (1978). Flow over a containment dyke. Journal of Fluid Mechanics, 87(1), 179e192. Hirst, I. (2006). Selected reconstructed data from the experiments of Greenspan & Johansson (1981). Private Communication. Huppert, H. E. (1982). The propagation of two-dimensional and axisymmetric viscous gravity currents over a rigid horizontal surface. Journal of Fluid Mechanics,121, 43e58. Ivings, M. J., Toro, E. F., & Webber, D. M. (2003). Numerical schemes for the 2D shallow layer equations including dry fronts. Computational Fluid Dynamics Journal, 12(1), 41e52. Ivings, M. J. & Webber, D. M. (2002). Development of SPLOT: a numerical model for the two dimensional shallow layer equations. Health and Safety Laboratory Report CM/02/08. Ivings, M. J., & Webber, D. M. (2007). Modelling bund overtopping using a shallow water CFD model. Journal of Loss Prevention in the Process Industries, 20(1), 38e44. Landau, L. D., & Lifshitz, E. M. (1959). Fluid mechanics. Pergamon, ISBN 0 08 009104 0. Thyer, A. M., Hirst, I. L., & Jagger, S. F. (2002). Bund overtopping d the consequence of catastrophic tank failure. Journal of Loss Prevention in the Process Industries, 15 (5), 357e363. Toro, E. F. (2001). Shock-capturing methods for free-surface shallow flows. Wiley. Webber, D. M., Ivings, M. J., Maru, W-A., & Thyer, A. M. (2009). Validation of the Shallow Water model “SPLOT” against experimental data on bund overtopping. Health and Safety Executive Research Report RR755. Webber, D. M., & Jones, S. J. (1987). A model of spreading vaporising pools. In J. Woodward (Ed.), International conference on Vapor Cloud modeling, Boston Massachusetts. USA (November 1987) proceedings. AIChE. Webber, D. M., & Brighton, P. W. M. (1986). Inviscid similarity solutions for slumping from a cylindrical tank. Journal of Fluids Engineering, 108, 238e240.
Nomenclature
A. Appendix
ðA:1Þ
667
b: bed height c, C: laminar and turbulent friction coefficients E: total energy fi: friction components f: friction vector g: gravitational constant (9.81 ms2) h: liquid depth hb: bund height H: tank fill height k: ‘kinetic factor’ P: potential energy Q: bund overtopping fraction rb: bund radius R: tank radius t: time T: kinetic energy ui: velocity components u: velocity vector V: tank volume vb: bund volume xi: Cartesian co-ordinate directions m: kinematic viscosity