Interaction Between Hydraulic and Numerical Models for the Design of Hydraulic Structures

The constant evolution of the calculation capacity of the modern computers implies in a permanent effort to adjust the existing numerical codes, or to create new codes following new points of view, aiming to adequately simulate fluid flows and the related transport of physical properties. Additionally, the continuous improving of laboratory devices and equipment, which allow to record and measure fluid flows with a higher degree of details, induces to elaborate specific experiments, in order to shed light in unsolved aspects of the phenomena related to these flows. This volume presents conclusions about different aspects of calculated and observed flows, discussing the tools used in the analyses. It contains eighteen chapters, organized in four sections: 1) Smoothed Spheres, 2) Models and Codes in Fluid Dynamics, 3) Complex Hydraulic Engineering Applications, 4) Hydrodynamics and Heat/Mass Transfer. The chapters present results directed to the optimization of the methods and tools of Hydrodynamics.

equalization between chambers, for which a one-dimensional numerical model was used; (b) the calculation of the amplitude of free surface oscillations in the lock chambers (which leads to increments in the hawser forces) due to close-to-resonance conditions under interaction with an oscillation in a flow partition component of the filling/emptying system (triggered by large turbulent eddies), for which a full three-dimensional numerical modeli.e., a CFD approach -was applied.

Equalization times
The equalization time, named as Filling/Emptying (F/E) time during the study, is a key parameter in establishing the system performance of a Lock Complex, as it has a direct impact on the vessel throughput, measured as the number of vessels passing through the system per day. Contractual requirements existed for the design of the Third Set of Locks of the Panama Canal, imposing maximum allowable F/E times for different scenarios. Minimization of these F/E times, through the reduction of local head losses, was the main strategy used during the design optimization process. Consequently, scale effects affecting these F/E times were carefully studied.

Description of the physical system
The Third Set of Locks of the Panama Canal, presently under construction, comprises twin lock complexes located near each ocean. Each complex has three lock chambers in series ( Fig. 2.1). These lock complexes allow vessels to be transported up or down between Gatun Lake and the oceans, spanning a 27 m water level difference. Each lock chamber has three side pools, called Water Saving Basins (WSB). This WSBs store part of the water used during the equalization operations, that otherwise would be flushed downstream towards the ocean. This stored water is then utilized to refill part of the lock, allowing a reduction of freshwater consumption during dry hydrological seasons. The lock chambers are connected through two longitudinal main culverts located within the lateral walls, running from the lake to the ocean. Four sets of valves are located along these culverts, which can isolate each chamber ( Fig. 2.1). The operation of these valves allows the successive equalization of the water level between chambers by gravity flow. In order to allow an even filling or emptying of the chambers, thus minimizing longitudinal water surface slopes and consequently hawser forces, water enters or exits each lock chamber through 20 ports located on each lateral wall. They are connected to secondary culverts which, in turn, connect to the main culvert at the midpoint of each chamber, through a carefully designed hydraulic component called Central Connection (CC) (Fig. 2.2). The CC was designed so as the flow coming from the main culvert is evenly split between both secondary culverts, before filling the chamber through the ports. Additionally, during the opposite operation, i.e., chamber emptying, equal discharges should flow through both secondary culverts, before its confluence towards the main culvert. The symmetry of the flow with respect to the chamber midpoint, for both the filling and emptying operations, is what ultimately guarantees low longitudinal water slopes and, consequently, low longitudinal hawser forces.

Fig. 2.2. Detail of the F/E system and Central Connection
The two branches of the CC are called 'U' and 'S' branches, in reference to the trajectory followed by the incoming flow ( Fig. 2.2). Note that the role as a U or a S branch depends on the hydraulic operation (filling or emptying). A conduit arising from the midpoint of each secondary culvert, and running below the chambers, connects them with the WSBs, after traversing a 'trifurcation' and a set of valves ( Fig. 2.2).

Description of numerical model
The hydrodynamic model of the F/E system was built using Flowmaster V7 (http://www.flowmaster.com/), a commercial code which solves one-dimensional transient flow over a network of conduits. Incompressible fluid and rigid pipe hypothesis were made, without compromising accuracy. For convenience, the whole F/E system was divided into sub models, one for each type of operation, such as equalization of the Upper Chamber with the Lake, equalization of two contiguous chambers, equalization of chamber with its Water Saving Basins, and so on.
The F/E system was represented as a network of interconnected component elements, namely:  Reservoirs, representing the lock chambers, WSBs, the lake, and the oceans. Level-Area relations were specified for each one of them. The lake and oceans were considered as infinite area constant level reservoirs.  Rigid rectangular pipes, representing primary and secondary culverts, WSB conduits, etc. The calculation of friction losses was made using Darcy-Weisbach and Colebrook-White equations, as a function of the flow Reynolds number and the effective roughness height of the conduit walls.  Local energy losses parameterized with a cross-section area and a head loss coefficient, representing most of the special hydraulic components, such as bends, bifurcations, transitions, etc.  Local energy losses expressed as laws for cross-section area and head loss coefficient in terms of a control parameter, representing valves for which the control parameter is the aperture.

Numerical modeling of physical model
The Third Set of Locks has been subject to physical modeling, both during the development of the conceptual design, and later during the design for the final project. Both physical models where commissioned to the Compagnie Nationale du Rhône (CNR), Lyon, France. The physical models were built at a 1/30 scale, comprising 2 chambers and one set of three WSBs. Extensive tests were made for various normal and special operations, measuring water levels, discharges, pressures and water slopes in the chambers. Some tests included the presence of a design vessel model, measuring hydraulic longitudinal and transversal forces over its hull. Based on these tests, a correlation between forces on the ship, and water surface slopes in the chamber in the absence of the ship (easier to measure and allegedly more repeatable), was established. This correlation was used to impose maximum values to the longitudinal and lateral water surface slopes, as contractual requirements. The flow in the hydraulic model was numerically simulated. Real physical dimensions of the physical model components (culverts, conduits, chambers) were used. Local head loss coefficients for the special hydraulic components were obtained through steady CFD modeling (see Section 3 for more details on CFD modeling), by calculating the difference between upstream and downstream mechanical energy, and subtracting energy losses due to wall friction. Most parts of the physical model were made out of acrylic (with a 0.025 mm roughness height), which behaves as a hydraulically smooth surface, for which the roughness height is completely submerged within the viscous sublayer (White, 1974). Some of the special hydraulic components, though, were built with Styrofoam (enclosed inside of acrylic boxes), as the initial expectations were that many alternative geometries would have to be tested, so this system would allow swapping with relative ease (very few alternatives were finally tested, due the great success of the optimization process carried out with CFD models). As it was later demonstrated that Styrofoam behaves as hydraulically rough at the physical model scale, most of it had to be coated with a low roughness layer of paint in order to avoid a spurious response (a scale effect in itself). The results obtained with the numerical model (water levels, discharges, pressures) showed a very good agreement with physical model measurements, for different operations and conditions. As an illustration, Figs 2.3 and 2.4 show comparisons for a typical Lock to Lock operation, with maximum initial head difference. All comparisons are presented with results scaled up to prototype dimensions.

Numerical modeling of the prototype
Practical knowledge exists about the discrepancies between F/E times as measured in a physical model and those effectively occurring at the prototype. For instance, USACE manual on hydraulic design of navigation locks (2006) states: ..."A prototype lock filling-and-emptying system is normally more efficient than predicted by its model"..."The difference in efficiency is acceptable as far as most of the modeled quantities are concerned (hawser forces, for example) and can be accommodated empirically for others (filling time and over travel, specifically)."... In the specific commentaries about F/E times, it suggests quantitative corrections: ..."General guidance is that the operation time with rapid valving should be reduced from the model values by about 10 percent for small locks (600 ft or less) with short culverts; about 15 percent for small locks with longer, more complex culvert systems; and about 20 percent for small locks (Lower Granite, for example) or large locks having extremely long culvert systems."... The alternative, rigorous strategy proposed in the present paper is to numerically simulate the flow in the prototype. This means using the physical dimensions of the prototype, the corresponding local head loss coefficients for the special hydraulic components, and the roughness height for concrete. Though the concrete wall also behaves as hydraulically smooth, the friction coefficient for smooth pipes is a function of the flow Reynolds number, as indicated by the "smooth pipe" curve in the Moody chart ( Fig. 2.5). For example, the Reynolds number in the primary culvert (in which most of the friction losses are produced) changes in time following the flow hydrograph, from zero to the peak discharge, and back to zero again. The peak discharge for 21 m initial head difference in a Lock to Lock operation is around 425 m 3 /s (the corresponding flow velocity is 7.87 m/s). This leads to a Reynolds number of around 6.5 10 7 for the prototype. When scaled to the physical model, the Reynolds number is only 3.9 10 5 , i.e., a drop of more than two orders of magnitude. The associated friction coefficients are then below 0.008 for the prototype, and about 0.014 for the physical model. The consequently higher friction losses produced in the physical model, exclusively due to scale effects, reduce the flow velocities, then increasing the F/E times. The numerical model contemplates the variation of frictional losses with the Reynolds number. Hence, it allows to be used in order to extrapolate the physical model results to those expected for the prototype, overcoming the distortion introduced by scale effects in the physical model results. For the Panama Canal Third Set of Lock, the validated 1D model was scaled up to prototype dimensions. Variations in local head loss coefficients, indicated by 3D models, were also introduced. Relatively little effects were observed in the simulations because of the change in local head loss coefficients. On the contrary, friction losses decreased significantly, as already explained. Consequently, for a typical Lock to Lock operation with maximum initial head difference, F/E times showed a 10% decrease (61 seconds Additionally, a 5% increase in the peak discharge of the main culverts was also observed ( Fig. 2.7). This has an effect over the pressures on the vena contracta, downstream of the main culvert valves (Fig. 2.8), which had to be contemplated during the design stage, as air intrusion had to be avoided (for contractual reasons), and because piezometric levels downstream of the valves were close to the roof level of the culvert for various special operating conditions. So avoiding scale effects was also significant to correctly deal with these two limitations.

Free surface oscillations
Free surface oscillations in the lock chambers leads to forces in the hawsers. Based on results from the physical model constructed during the development of the conceptual design, a correlation was found between these forces and the free surface slope in the absence of the vessel, as already mentioned in Section 2. Hence, the free surface slope was used as an indicator for the hawser forces. As a design restriction, a maximum value of 0.14 ‰ was contractually established for the longitudinal water surface slope. A 2D (vertically averaged) hydrodynamic model, based on code HIDROBID II developed at INA (Menéndez, 1990), was used to simulate the surface waves. It was driven by the inflow from the ports, specified as boundary conditions through time series for each one of them, that were obtained with the 1D model described in the Section 2. Fig. 3.2 shows the comparison between the calculated longitudinal free surface slope (using the dimensions of the physical model) and the recorded one at the physical model, for a case www.intechopen.com with a relatively low initial head difference (9 m in prototype units) between the Lower Chamber and the Ocean. The agreement is considered as very good, taking into account that the numerical model does not include the resolution of turbulent scales (which introduce a smaller-amplitude, higher-frequency oscillation riding on the basic oscillation). However, the 2D model completely fails to correctly predict the longitudinal free surface slope for higher initial head differences, as observed in Fig. 3.3 for a Lock to Lock operation with an initial head difference of 21 m. More specifically, the recorded oscillation indicates a quite more irregular response, with a much higher amplitude than the one calculated with the www.intechopen.com numerical model. This indicates that turbulence scales are exerting a significant influence, so a more elaborated theoretical approach is needed. Hence, 3D modeling of the combination Central Connection + Secondary Culvert + Ports + Lock Chamber (actually, only half of the chamber, assuming that the flow is symmetrical with respect to the longitudinal axis) was undertaken using a Large-Eddy Simulation (LES) approach (Sagaut, 2001).

Improved theoretical approach
As sub-grid scale (SGS) model for the LES approach, a sub-grid kinetic energy equation eddy viscosity model was used (Sagaut, 2001). Deardorff's method was selected to define the filter cutoff length (Sagaut, 2001). A wall model was considered to treat the boundary conditions at solid borders; Spalding law-of-the-wall -which encompasses the logarithmic law (overlap region), but it holds deeper into the inner layer -was selected for the velocity (White, 1974), while a zero normal gradient condition was taken for the remaining variables. At the inflow boundary, in addition to the ensemble-averaged velocity (which arises from the 1D model), the amplitude of the stochastic components were provided (Sagaut, 2001): 4% for the longitudinal component, and 1.3% for the transversal one, values associated to a fully developed flow, very appropriate for the present problem; additionally, a weighted average of the previous and present generated stochastic components was imposed in order to add some temporal correlation; for the turbulent kinetic energy, a zero normal gradient was taken. For the free surface at the Chamber, the rigid-lid approximation was used, where uniform pressure was imposed, together with zero normal gradient conditions for the remaining quantities. The model was implemented using OpenFOAM (Open Field Operation And Manipulation), an open source toolbox for the development of customizable numerical solvers and utilities for the solution of continuum mechanics problems (Weller et al., 1998). The model solves the integral form of the conservation equations using a finite volume, cell centered approach in the spirit of Rhie and Chow (1983). PISO (Pressure Implicit with Splitting of Operators) algorithm is used for time marching (Ferziger & Peric, 2001). Fig. 3.4 presents a view of the model domain. The computational mesh was composed by 1.5 million elements. Special considerations were made for the mesh near the wall, as the center of the first cell has to lie within a distance range to the wall -30  y+  300 -to rigurously apply the logaritmic velocity profile as boundary condition (Sagaut, 2001). Typical computing times for stabilization with a steady discharge, in a Core i7 PC running 8 parallel processes, were 3 to 8 days. When complete hydroghaphs were simulated (of approximately 550 secs), 15 to 30 days of computing time were required. By parallelizing the simulation using more than one PC, computing times were reduced, though non-linearly.

www.intechopen.com
Note that the rigid-lid approximation implies that the free surface oscillations are not solved by the 3D model; this was done in order to avoid extremely high computing times. Instead, the 3D model provided the time series of the flow discharge for each port, which were used to drive the 2D model of the chamber. Alternatively (and less costly in post-processing), the time series of the discharges at the U and S branches of the Central Connection, provided by the 3D model, were used to feed the 1D model, from which the discharge distribution among ports was obtained, and used to feed the 2D model. Fig. 3.5 shows the longitudinal water surface slope obtained with the two approaches (using the dimensions of the physical model), and their comparison with the results from the physical model, for the high initial head difference case. It is observed that both numerical simulations are now able to capture the high amplitude oscillations, indicating that large eddies must be responsible for this amplification phenomenon. Note that the numerical results with input straight from the 3D model show oscillations, associated to large eddies, which are not present in the ones with input through the 1D model (which filters out those oscillations), but they are quite compatible between them. The differences between the numerical results and the measurements at the physical model are due essentially to the variability of the system reponse (variations in amplitude and phase of the oscillations), under the same driving conditions, due to the stochastic nature of turbulence. This was verified both experimentally ( Fig. 3.6a) and numerically (Fig. 3.6b) by repeating the same test (in the case of the numerical model, using the 'through 1D model' approach, and different initializations for the stochastic number generator). This behavior puts a limit to the degree of agreement that can be attained between the results from the numerical and physical models. In any case, the maximum amplitudes for any of the experimental or numerical realizations are relatively consistent among them. Before proceeding to simulate prototype conditions, it is relevant to analyze the response provided by the numerical model, in order to be confident about using this tool to make such a prediction. Specifically, the physical mechanisms involved in the present problem should be fully understood. This is performed in the following. Fig. 3.7a shows the time series of the discharges through the U and S branches of the Central Connection (in prototype units), according to the 3D numerical model. It is observed that, for the higher discharges, they present oscillations, which seem coherently out-of-phase. The difference between those discharges is shown in Fig. 3.7b (together with the total discharge, i.e., the one through the Main Culvert). It is effectively observed that this difference oscillates, and that during the time window of higher discharges (above about 250 m 3 /s) there is a dominant period of oscillation which spans from 40 to 80 seconds, approximately.

www.intechopen.com
Now, these periods are close to, and include, the period of free surface oscillations in the Chamber (around 70 seconds), indicating that conditions close to resonance are achieved, thus resulting in an amplification of the free surface oscillation, which is the observed effect on the water surface slope. As in the numerical simulation the free surface was represented like a rigid lid, the oscillation in the discharge difference between the two branches of the Central Connection is not influenced at all by free surface oscillations themselves, i.e., the dominant period arises from the flow properties in the Central Connection. This dominant period must then be associated to the largest, energy-containing eddies (the ones resolved with the LES approach). Before pursuing with the analysis, it is worth to confirm that the close-to-resonance conditions are responsible for the amplification of the free surface oscillations. Hence synthetic hydrographs for the U and S branches of the Central Connection were built, introducing a purely sinusoidal oscillation to their difference during the higher-discharges time window, as indicated in Fig. 3.8a for 2 m 3 /s amplitude of oscillation, and two different www.intechopen.com periods: 60 and 120 seconds. Fig. 3.8b presents the results from the 2D model for the two different periods. It is clearly observed that amplitude amplification occurs for the 60 seconds period (during the time window of forced discharge oscillation), which is under close-to-resonance conditions. On the contrary, the amplitude attenuates for the 120 seconds period, which is far from the resonant period. In Fig. 3.9 the results of the 2D model with the synthetic hydrographs, for the 60 seconds case, are compared with the physical model measurements, indicating a quite reasonable agreement, providing an extra validation to the physical explanation of the observed phenomenon. The 2D model results are much 'cleaner' than the measurements because the triggering signal (discharge difference) has a single frequency, in lieu of the set of frequencies associated to the turbulent eddies. Now, the relation between the discharge oscillation and the larger, energy-containing eddies generated at the wake zones, in the U and S branches (Fig. 3.10), is analyzed. The characteristics of those large eddies are quantified based on an analysis of scales (Tennekes & Lumley, 1980). The size of these eddies, the so called 'integral scale' of turbulence in the wake region, is limited by the physical dimensions of the Secondary Culvert. Hence, it is of the order of the conduits widths (4.5 m for the U branch, and 3.1 m for the S branch). On the other hand, the relation between the velocity-scale of the largest eddies and the ensemblemean of the incoming velocity is of the order 10 -2 . It is assumed that this relation is 1% if the section-averaged velocity at the Secondary Culvert (which changes with the total discharge) is taken as a reference. The relation between the integral scale and the velocity scale provides a scale for the period of the largest eddies. Fig. 3.11 shows the variation of the period-scale of the largest eddies, for the two branches of the Central Connection (which differ between them due to the different incoming velocities), with the total discharge (i.e., the one through the Main Culvert). It is claimed that the interaction of the largest eddies of the U branch with those of the S branch is responsible for the generation of the coherent out-of-phase oscillations in the discharges through each branch (as explained below). When this oscillation has a period close to the Chamber free surface oscillation period, also represented in Fig. 3.11, amplification occurs, as already explained. From Fig. 3.11, it is observed that close-to-resonance conditions should be expected for total discharges higher than about 200 m 3 /s, and up to at least 500 m 3 /s. This is completely consistent with the numerical and physical model results obtained for high initial head difference.  In order to complete the analysis, an explanation for the mechanism of interaction between the largest eddies of the U and S branches, leading to the coherent out-of-phase oscillations in the discharges through each branch, is undertaken, inspired in the one for a von Karman vortex street (Sumer & Fredsoe, 1999). Vortices (largest eddies) are shed from the separation points. Subject to small disturbances, one of those vortices, for example the one on the U branch, grows larger, increasing the blockage effect in that branch; as a result, the discharge through the U branch decreases, leading to an increase of the discharge through the S branch (in order to maintain the total discharge). Now, the next vortex shed in the S branch is of higher intensity, due to the increased incoming flow velocity in this branch; but this has the effect of increasing the blockage of the S branch, then producing a decrease in the discharge through that branch, and a consequent increase of the discharge through the U www.intechopen.com branch. Then, the phenomenon described for the S branch now occurs in the U branch, leading to a cyclic behavior, as observed.

Numerical modeling of the prototype
Having understood and numerically modeled, with a reasonable degree of satisfaction, the oscillatory phenomenon which develops in the Central Connection for the physical model, the flow in the prototype was simulated in order to determine the behaviour at that scale. The calculation with the 3D model was undertaken using the same (rescaled) mesh as for the physical model. Though the condition on the location of the first node, in order to correctly represent the wall shear stress, is not fulfilled, it is considered that this should not significantly affect the results, based on the fact that tests performed in the physical model including triggering devices indicated that the appearance of the oscillatory phenomenon is not conditioned by the location of the separation point. Fig. 3.12a shows the evolution of the longitudinal water surface slope arising from the results of the 3D model. It is compared with the numerical results for the physical model; the ones arising from the 1D modeling approach (no 3D LES model) are also represented, as a reference. Note that the prototype response is significantly more noisy than the physical model response, as it includes a higher range of turbulent frequencies. It is observed that, though the amplification effect manifest in the prototype (the amplitude of oscillation is higher than the one predicted by the 1D model), its amplitude is definitely smaller than the one for the physical model. In fact, the oscillation in the discharge difference, presented in Fig. 3.12b, is sensitively less significant for the prototype than for the physical model (compare with Fig. 3.6b). It is especulated that this should be due to differences in the energy spectrum: the larger eddies of the prototype would contain less energy than the corresponding ones in the physical model. It is concluded that, for this problem, scale effects tend to increase the amplification effects.

Conclusions
The proposed strategy for the design of hydraulic structures, consisting in a first stage where the flow in the physical model is numerically simulated, in order to validate the numerical model, and in a second stage where the flow in the prototype is numerically simulated, in order to extrapolate the results to this scale, has been shown to be effective in correcting for the scale effects present in the physical model. This has been illustrated for the particular case of the design of the Third Set of Locks of the Panama Canal, for two problems with quite different levels of complexity. The first problem was the determination of the time for water level equalization between chambers, using a one-dimensional numerical model. Friction losses are shown to be overrepresented in the physical model, leading to larger equalization times. Differences of the order of 10% are calculated for a case with maximum initial head difference. The second problem was the calculation of the amplitude of free surface oscillations in the lock chambers, due to close-to-resonance conditions, under interaction with an oscillation in a flow partition component of the filling/emptying system, using a full three-dimensional numerical model with a LES approach. Differences in the energy spectrum lead to a significant amplification of the amplitude of oscillation in the physical model. The paper indirectly stresses, through an in-depth analysis of the involved physical mechanisms for the case studies, the necessity of thoroughly understanding the responses provided by the numerical model, in order to be confident in using the tool to make predictions at the prototype scale.

Acknowledgments
In addition to the present authors, Emilio Lecertúa, Martín Sabarots Gerbec, Fernando Re and Mariano Re were part of the numerical modeling team for the Panamá Canal Project. The team worked under the coherent supervision of Nicolás Badano (MWH), responsible for the hydraulic studies, with the help of Mercedes Buzzela. The smooth interaction with the responsible for the physical model, Sébastien Roux, from CNR, was fundamental in order to achieve the goals of the study. The constant evolution of the calculation capacity of the modern computers implies in a permanent effort to adjust the existing numerical codes, or to create new codes following new points of view, aiming to adequately simulate fluid flows and the related transport of physical properties. Additionally, the continuous improving of laboratory devices and equipment, which allow to record and measure fluid flows with a higher degree of details, induces to elaborate specific experiments, in order to shed light in unsolved aspects of the phenomena related to these flows. This volume presents conclusions about different aspects of calculated and observed flows, discussing the tools used in the analyses. It contains eighteen chapters, organized in four sections: 1) Smoothed Spheres, 2) Models and Codes in Fluid Dynamics, 3) Complex Hydraulic Engineering Applications, 4) Hydrodynamics and Heat/Mass Transfer. The chapters present results directed to the optimization of the methods and tools of Hydrodynamics.