Paper received: 18.02.2009 Paper accepted: 03.04.2009 Quantification of Flow Kinematics Using Computer-Aided Visualization Tom Bajcar *- Brane Širok - Matjaž Eberlinc University of Ljubljana, Faculty of Mechanical Engineering, Slovenia The paper presents a study of the ability of the visualization methods to measure velocity field in the fluid flow from the concentration field of a passive tracer that is introduced to the fluid. The method is based on the advection diffusion equation. A sequence of successive flow images was processed to obtain the required unknown terms. Concentration of the passive tracer was acquired through greyscale intensity on images at chosen measurement points. Standard numerical procedures were applied to obtain spatial- and time derivatives of the passive tracer concentration. The method was tested in an experiment with an air jet flowing out of a pipe, where the velocity field in a jet was measured with a hotwire anemometer. Comparison of the results of both methods shows good agreement, albeit with certain deviations that are a consequence of simplified approach of the visualization method. The results prove the suitability of the visualization method for quantitative measurements of flow kinematics. Further studies and improvements of the method could provide a simple and relatively cheap yet reliable tool not only for preliminary assessment of the flow velocity field, but also for velocity measurements in complex flow structures. © 2009 Journal of Mechanical Engineering. All rights reserved. Keywords: computer-aided visualization, flow kinematics, velocity measurment 0 INTRODUCTION Flow visualization, common or computer-aided, is a valuable tool in fluid dynamics for observation and monitoring of the fluid flow structures. Structures in the flow become visible by introducing a pollutant - passive tracer into the flow, which is usually in a form of dye or smoke. The pollutant should be properly chosen for each fluid by means of density and should exactly follow the fluid flow. Visualization then focuses on passive tracer flow and its spatial- and time changes of concentration, which helps to identify the qualitative properties of a particular fluid flow. The images of the fluid/pollutant flow acquired by a camera can then be further processed for different purposes (e.g. vortex separation analysis [1], spatial determination (location) of recirculation zones [2], monitoring of cavitation structures [3], influence of the fluid flow on flow deposits [4], etc.). The question that arises at this stage is whether it is possible to acquire the quantitative properties of the flow by means of flow visualization, i.e. extracting the kinematic quantities of the flow such as velocity field from a set of flow images with a passive tracer. Up to now, computer-aided visualization (CAV) could be applied to measure the velocity field in the flow only by introducing solid particles in the flow, a common example being the Particle Image Velocimetry (PIV) method. In order to research qualitative as well as quantitative properties of the flow by means of computer-aided visualization, at least two methods/systems should be applied simultaneously [5] and [6]; one with the pollutant in the flow for qualitative measurements, such as Planar Laser Induced Fluorescence (PLIF), and the other for qualitative measurements, such as PIV method. Such combined systems are often complex since two systems require two cameras with appropriate filters and two types of pollutants. PIV method itself requires pulsing laser light sheet, while the number of added particles in the flow is usually limited by the ability of the image processing software to recognize a particular particle in the crowd on two successive images. Other measurement methods, such as Hot-Wire Anemometry (HWA) or Laser-Doppler Anemometry (LDA) are capable of measuring velocity field in one measurement point at a time. They are either intrusive methods that cause physical disturbance in the flow (HWA) or quite complex and require tedious preparations (LDA). Apart from that, the majority of common measurement methods are unable to measure velocities in non-stationary or *Corr. Author's Address: University of Ljubljana, Faculty of Mechanical Engineering, Aškerčeva 6, 1000 Ljubljana, Slovenia, tom.bajcar@fs.uni-lj.si periodical flow structures, such as von Karman vortices. The aim of this paper is to research examine in practice the possibility of measuring the flow velocity field only by applying a visualization method on a fluid flow with a passive tracer. A method should be based on known physical relations to combine the field of pollutant concentrations with the field of flow/pollutant velocities. The obtained values should be the result of processing a sequence of successive flow images. Such a method would have several advantages over the existing measurement methods mentioned above for it could measure the velocity field in arbitrary many measurement points at virtually the same time since the same sequence of flow images can be used for all the measurement points. Other advantages would be simple and readily available measurement equipment (camera, PC) and ease of handling. 1 METHODOLOGY The most important step of the described methodology is to link a vector field of velocities with a scalar field of concentrations, the latter represented by greyscale images. By introducing a pollutant with the concentration N (i.e. the number of pollutant molecules per unit volume of fluid) to the fluid, the gradient of N in the direction x can be maintained by supplying the pollutant molecules at a rate [7]: 0 = - D-, dx (1) where 0 is a current in the direction of x, and D denotes the molecular diffusivity. Eq. 1 represents Fick's law of diffusion. If unsteady transport processes are considered where the rate of change of concentration N in an elementary volume is represented by the difference between the inflow and the outflow of the pollutant through the boundaries of such a volume, then the Eq. 1 can be rewritten as [7]: dN = DV2N . dt (2) Where fluid movement (i.e. convection) is present, Eq. 2 can be further extended by replacing the partial time derivative with the total derivative: dN d( Nvi ) 2 -—^ = DV 2N, dt dxi where the second term on the left-hand-side of Eq.3 involves the effects of fluid convection. Eq. 3 is known as advection diffusion equation and represents the basic connection between tracer / pollutant concentration and the flow kinematics. By knowing the values of concentration N and, consequently, the values of its derivatives along with the molecular diffusivity D of the pollutant in the fluid, the only unknown that remains in Eq. 3 is velocity v, i.e. its components vi in directions x. While the molecular diffusivity D is a property of matter and is defined for different pairs of substances in accordance with pressure and temperature, the concentration N of the pollutant changes constantly and should be instantaneously measured; in this case, its value should be extracted from greyscale images of the flow. A greyscale image of the fluid flow with the introduced tracer forms a matrix of pixels with different greyscale levels. Normally, integer greyscale levels extend from 0 (black) to 255 (white); they could be normalized on the interval between 0 and 1. When the pollutant introduced to the fluid is illuminated by external lighting, its greyscale values are usually higher than those of the fluid. The higher the pollutant concentration in the fluid the higher the greyscale values of the pollutant in the greyscale images. Therefore, A x N, (4) where A denotes an average greyscale value in the selected window of pixels on a greyscale image. It was shown by other authors [8] and [9], that the proportionality between greyscale values and pollutant concentration in Eq. 4 is linear. Average greyscale level in the window of pixels can be obtained by a simple relation [2] : l m l m ' i=i 1=1 (5) where E(i,j) presents the grey intensity of the observed single pixel at a position (ij) inside the window of l x m pixels. The numerical calculation of the time derivative of concentration (dN/dt) is obtained from the average greyscale values A of a window whose position is fixed on successive images (Fig. 1). By knowing the time between two successive images At (easily obtained from the image acquisition frequency of the camera), the term dN/dt can be approximated by the expression: dN ^AA dt ~ At ' (6) where AA denotes the difference in the average greyscale level in a fixed window between two successive images that were taken At apart. The spatial derivative of the pollutant concentration (BN/ck) can be calculated from a single image of the flow. It is now the window that moves in the selected direction xi in order to get the spatial difference of the average greyscale level A, as depicted in Fig. 2. A calculation of spatial derivatives can then be performed using standard numerical techniques, such as a method of central differences [10] as shown in Eqs. 7a and 7b (neglecting the truncation error): dN _ A(xt + Axt) - A(xt - Axt) Bx- 2 ■ Ax, (7a) d 2N dxf A(xi + Ax;) - 2 • A(x,) + A(xt - Axt) (Axi) 2 (7b) In a two dimensional system, Eq. 3 represents a set of two linear differential equations with two unknowns, i.e. two velocity components. Solving such a set of equations requires the initial conditions to be specified, which could often be quite tedious and difficult due to the fact that kinematic properties of the flow are generally not known unless measured by some other method. Instead, the set of differential equations (Eq. 3) can be converted to a set of ordinary linear equations with velocity components as well as their spatial derivatives as unknowns. Such sets are generally easy to solve numerically and work well even when overdetermined (i.e. the number of equations is larger than the number of unknowns). The number of unknowns therefore doubles in comparison with a set of differential equations, which implies that for a two dimensional system at least four ordinary equations are required. Each equation requires two successive images (because of time derivative calculation) of the flow. In other words, at least five successive images of the flow are needed to get the information on two velocity components of the flow. The above described methodology is based on known physical relations and reflects basic principles of determining the flow kinematics by means of computer-aided visualization. Its validity can be tested by comparison with other more conventional measurement methods. For this purpose, a simple experiment was conducted using hot-wire anemometry as an alternative measurement method. Fig. 1. Fixed square window in three successive images of the fluid flow for calculation of the time derivative of concentration N Fig. 2. Movement of the window on a single image in order to get the spatial difference in concentration N 2 EXPERIMENTAL SET-UP Fig. 3 shows the experimental set-up for measurements of a velocity field along an air jet coming out of a vertical pipe. Reynolds number of the flow was determined in the pipe at its outlet and amounted Re = 1300. The flow in the pipe was therefore laminar. Measurements were performed with HWA and CAV. Velocity distribution was measured with the two-component-hot-wire anemometer in the selected measuring points. The two-component hot-wire anemometer was used for measurements of instantaneous velocity fluctuations in both streamwise (y) and lateral (x) direction. Measurements were performed in accordance with the procedure by Bruun [11] and Jorgenson [12]. At selected operating point, 153 measurement points were selected on a mesh 9 x 17 measuring points (Fig. 3, point 2), where the distance between each measuring point was set to 5 mm with the positioning table (Fig. 3, point 8). For measurements with hot-wire anemometry, a DANTEC 55P62 constant temperature hot-wire anemometer sensor (CTA) system was used. The sensor was made out of Pt-plated tungsten wire and had the length and diameter of 1.25 mm and o5 ^m, respectively. The anemometer amplifier cut-off frequency was set to 10 kHz. The hot-wire anemometer working temperature was 250 °C. The positioning of the two-component-hot-wire anemometer sensor on the wind tunnel was such that it allowed the measurement of instantaneous velocities in two directions. The positioning of the anemometer sensor was performed using a PC controlled precision positioning device. For the acquisition of the hot-wire anemometer signal, a 16-bit data acquisition board from National Instruments was used. Fig. 3. Experimental set-up of the test section. (1. illumination; 2. measuring points; 3. smoke generator; 4. pipe 500 mm*100 mm; 5. high-speed camera; 6. inlet of the fluid volume flow; 7. hot-wire anemometer; 8. positioning table; 9. computer) The acquisition time was 6 s and the acquisition frequency was 50 kHz. The output of the anemometer was filtered prior to the A/D conversion with SCXI module using a low pass 4th order Bessel filter with frequency of 10 kHz. Lab View software was used for data acquisition and storage. Calibration of the hot-wire anemometer was performed on the measuring station for calibrating anemometers. In accordance with King's law [11], the constants for the first hot-wire A1 = 1.7325, B1 = 0.8552, and n1 = 0.4517 and for the second wire A2 = 2.1419, B2 = 0.6819, and n2 = 0.4928 were determined. Temperature correction was applied during calibration and during the measurements. Temperature was measured using a Pt-100 class A resistance thermometer with four wire connections and Agilent 34970A instrument. Actual streamwise velocity vy and lateral velocity vx were calculated from the hot-wire anemometer output according to King's law equation [11] and equation procedure by Jergenson [12]. Measurements were performed with the same measurement equipment, which was calibrated prior to measurements [13]. The main sources of total measurement uncertainty [14] result from the uncertainty of selecting the operating point and the fluctuations, from selecting the starting point of positioning, rotational speed of the fan, calibration, positioning of the hot-wire anemometer sensor with a positioning table, linearization, A/D resolution, data acquisition, temperature correction and humidity and hot-wire anemometer limited frequency response. A separating transformer was used to separate fan power supply from the measurements' power supply. Overall measuring uncertainty of instantaneous velocity was estimated to 2.8% of the measured value. For CAV measurements, passive tracer smoke consisting of vaporized paraffin oil was introduced into fully developed flow in order to visualize the flow [15]. CAV measurements were performed in the same measurement points as previously with hot-wire anemometry. A passive tracer was set inside of the pipe, Fig. 3, point (4), Camera Dragonfly Express IEEE-13 94b with Edmund Optics 75 mm Double Gauss lens, (Fig. 3, point 5) lens was used for image acquisition, with image acquisition frequency of 400 Hz. The recorded images had 8-bit grey level depth, resolution of 600*250 pixels, and were captured with shutter speed of 0.6 ms. The camera was placed perpendicular to surface in the distance of 1 m from the test section of the wind tunnel. For illumination, a Vega Velum DC 150 W (Fig. 3, point 1) continuous flicker-free light source with the line light guide was used. Software package LabView was used for camera set-up, acquisition and storage of digital images. For each measurement point, 1000 successive images were taken. The number of images was limited by the time of uniform passive tracer smoke generation. Overall measuring uncertainty of instantaneous passive tracer concentration was estimated to 3.8% of the measured value [15]. 3 RESULTS The sequence of images acquired during the experiment was processed as described in Section 2 using Matlab program package. Both velocity components were determined by computer-aided visualization in the same measurement points in a x-y plane (Fig. 3) as in the case of experimental (hot-wire) measurements. The value for molecular diffusivity D of paraffin was taken from [16] for decane vapour (C10H22) in air (T = 23 oC, p = 1 bar): D = 7.0110-12 m2/s. A sequence of 30 successive images was used to define the flow velocity field in x-y plane; this number of images in the sequence assured that the set of Eq. 3 was well defined in all measurement points, even if the passive tracer was scarce (outer edges of the measurement plane x-y). On the other hand, averaging the velocity components over larger sequences of images require substantially more computing time and may as well result in certain loss of information at a given acquisition frequency. The value of a step (Axi) in Eq. 7 was chosen the smallest possible in order to minimize the truncation error in numerical computation of derivatives. The step therefore amounted 1 pixel, its size being bounded by the image resolution. Fig. 4 shows qualitatively the flow velocity field measured by both, the HWA method (Fig. 4a) and CAV method (Fig. 4b). It can be seen from Fig. 4 that the agreement between the results of CAV method and those of HWA method is good in the most measurement points. The results from both methods show that the velocity of the jet decreases away of the orifice and that the flow at the outer edges of the jet tends to move radially outward (the jet spreads in an axial as well as in radial direction). Nevertheless, there are some differences in length and orientation of velocity vectors between both methods, which are mostly encountered at the outer edges of the jet. At these locations, the paraffin vapour is scarce, which adversely affects the possibilities of CAV method, since it needs sufficient passive tracer concentration on every image in the sequence in order to work properly. Fig. 4. Qualitative comparison of velocity measurement results between HWA method (a) and CAVmethod (b) Quantitative comparison between the results of both methods is presented for measured velocity components in streamwise (y) and normal (x) direction. Normalized values for x-and y-direction are used for the presentation of results, where x*, y* e [0,1]. x* = 0 and x* = 1 correspond to the leftmost and rightmost measuring points, respectively, on the measurement plane (Fig. 3). Similarly, y* = 0 and y* = 1 corresponds to the closest and the farthest measuring points from the pipe orifice, respectively. Fig. 5 shows measured values of velocity component vy at normalized streamwise stations y* = 0.125, 0.375, 0.625 and 0.875 for both measurement methods. A comparison of results presented in Fig. 5 confirms that velocity profiles in y-direction have similar form regardless of the applied measurement method. However, the results of the CAV method show more intense fluctuations of vy along the lateral direction x in comparison with the results of the HWA method, which are somehow more pronounced far away from the pipe orifice (Fig. 5d). This phenomenon may well be the consequence of the vorticity of the flow together with the data acquisition frequency of a particular measurement method. At selected Re number of the flow from the pipe, a formation of intense eddies is observed at streamwise stations y* > 0.3. 0 0.2 0.4 0.6 0.! 0 0.2 0.4 0.6 0.8 1 Fig. 5. Velocity component in y-direction at different streamwise stations y*: a) y * = 0.125, b) y * = 0.375; c) y* = 0.625; d) y* = 0.875 Data acquisition frequency of HWA method was about 50 kHz and the measurement in a single point with this method took about 2 s. The results obtained by the HWA method were thus averaged along this period of time, which was long enough to ensure that the rotation of large vortices did not significantly affect the mean velocity of the flow. For the same time period, the CAV method would need a sequence of about 700 images, which would affect the processing / computing time considerably. On the other hand, the acquisition frequency of the CAV method in this case is more than 140-fold lower than the acquisition frequency of the HWA method. Increasing the sequence of images in order to get more accurate result would, therefore, be applicable only if significant deviations from the mean velocity (fluctuations) in the flow last long enough to be sensed by the CAV method. In other words, a higher acquisition frequency of the camera is needed for the CAV method to be more closely comparable with the HWA method. Table 1 presents the relative difference in mean values of vy (calculated along the lateral distance x) between both methods for each of the streamwise stations presented in Fig. 5. Table 1. Relative difference (according to HWA method) in mean values of vy between both methods at four streamwise stations y* Streamwise station Relative difference y* in mean vy values 0.125 0.095 0.375 0.092 0.625 0.099 0.875 0.108 The results in Table 1 show the average relative difference of about 10% at each of the four selected streamwise stations. It can, therefore, be concluded that the CAV method can of predicting the values of vy fairly well in comparison with the HWA method even if a short sequence of images is used. The values of velocity component in x-direction are shown in Fig. 6. Similar conclusions as in the case of velocity component values in y-direction (Fig. 5) can be drawn for the velocity components in x-direction (Fig. 6). The maximum absolute values of vx do not change much along the streamwise y-direction, the exception being the stations close to the pipe orifice (Fig. 6a). X* X* Fig. 6. Velocity component in x-direction at different streamwise stations y*: a) y* = 0.125, b) y* = 0.375; c) y* = 0.625; d) y* = 0.875 The jet obviously spreads more pronouncedly in the positive x-direction (positive vx values), except at streamwise station y* = 0.375 (Fig. 6b), where spreading in both, +x and -x directions is approximately equal; this is most probably the result of formation of large vortices mentioned earlier. Nevertheless, the values of vx are much lower than the values of vy (Fig. 5) and this was sensed by both methods. The difference in mean values of vx of the CAV method relative to those of the HWA method is shown in Table 2. Table 2. Relative difference (according to HWA method) in mean values of vx between both methods at four streamwise stations y* Streamwise station y* Relative difference in mean vx values 0.125 0.143 0.375 0.090 0.625 0.139 0.875 0.128 The relative differences in mean vx in Table 2 are somewhat larger than those of vy (Table 1), which is mainly the consequence of much smaller values (close to 0 and up to 0.04 m/s) compared to values of vy (Fig. 5). Nevertheless, the relative difference does not exceed 15%, which could be considered as a good result for the basic CAV method. 4 CONCLUSIONS A basic principle for quantification of kinematic properties of the flow by means of computer-aided visualization was presented in the paper. A two-dimensional velocity field of a flow is measured from a sequence of successive images of the pollutant (passive tracer) in the flow. The methodology itself is based on advection diffusion equation, which couples the velocity field with the field of concentrations. The calculation of needed parameters was carried out by image processing and common well known numerical procedures. The method itself was tested with a simple experiment where a conventional measurement method - hot wire anemometry - was applied in parallel with the computer-aided visualization. The comparison of results from both methods showed that the visualization method, even in its basic form, is capable of calculating both, the direction of the flow (i.e. of the passive tracer in the flow) as well as the magnitude of the flow velocity with fair accuracy. Apart from that, the velocity component profiles at different streamwise stations are comparable between both methods. The relative differences in mean values are within the range of 10 and 15% for vy and vx, respectively. The difference between both results is low enough to conclude that the principle of measuring the flow velocity field by means of computer-aided visualization is feasible, which was the main purpose of the study. On the other hand, the visualization method at this stage is still inferior to other common measurement methods, such as HWA, PIV or LDA, in terms of accuracy. Future research should therefore involve several refinements not only in the field of numerical calculations, where the influence of the step size and the application of different methods for derivative computation on the end result would be interesting, but also in the field of measurement preparation, where particular studies could determine the influence of illumination, three-dimensional flow structures, non-stationarity, number of images in the sequence, etc. on the uncertainty of the results. The method according to the presented algorithms and equipment requirements (appropriate image resolution and acquisition frequency) should be applicable to turbulent flows with higher Re numbers as well. Nevertheless, further studies in this area should prove the qualities and limitations of the presented method. 5 REFERENCES [1] Širok, B., Potočar, E., Novak, M. (2000) Analysis of the flow kinematics behind a pulsating adaptive airfoil using computer-aided visualisation, Strojniški vestnik -Journal of Mechanical Engineering 46(6), p. 330-341. [2] Širok, B., Bajcar, T., Dular, M. (2002) Reverse flow phenomenon in a rotating diffuser, Journal of Flow Visualization and Image Processing 9, p. 193-210. [3] Širok, B., Blagojevic, B., Bajcar, T., Trenc, F. (2003) Simultaneous study of pressure pulsation and structural fluctuations of a cavitated vortex core in the draft tube of a Francis turbine, Journal of Hydraulic Research 41, (4), p. 541-548. [4] Bajcar, T., Blagojevic, B., Širok, B., Dular, M. (2007) Influence of flow properties on a structure of a mineral wool primary layer, Experimental thermal and fluid science 32(2), p. 440-449. [5] Hu, H., Saga, T., Kobayashi, T., Taniguchi, N. (2004) Analysis of a Turbulent Jet Mixing Flow by Using a PIV-PLIF Combined System, Journal of Visualization 7(1), p.33-42. [6] Reungoat, D., Rivière, N., Fauré, J.P. (2007) 3C PIV and PLIF Measurement in Turbulent Mixing, Journal of Visualization 10(1), p. 99-110. [7] McComb W.D. (1996) The Physics of Fluid Turbulence, Oxford Engineering Science, Clarendon Press, Oxford. [8] Aider, J.L., Westfreid, J.E. (1995) Visualizations and PDF of the fluctuations of a passive scalar in a turbulent Goertler flow, Experimental and Numerical Flow Visualization 218, p. 123-130. [9] Simoens, S., Ayrault, M. (1994) Concentration flux measurements of a scalar quantity in turbulent flows, Experiments in Fluids 16, p. 273-281. [10] Mathews, J.H., Kurtis, K.F. (2004) Numerical Methods Using Matlab, 4th Edition, Prentice-Hall, Upper Saddle River. [11] Bruun H.H. (1995) Hot-wire Anemometry -Principles and Signal analysis, Oxford University Press, New York. [12] Jorgensen F.E. (2005) How to measure turbulence with hot-wire anemometers, A practical guide, Dantec Dynamics. [13] Coleman, H.W., Steele, W.G. (1999) Experimentation and Uncertainty Analysis for Engineers, John Wiley & Sons Inc., New York. [14] Eberlinc, M., Širok, B., Hočevar, M., Dular, M. (2009) Numerical and experimental investigation of axial fan with trailing edge self-induced blowing, J. Eng. Research 73 (in print). [15] Eberlinc, M., Dular, M., Širok, B., Lapajne, B. (2008) Influence of Blade Deformation on Integral Characteristic of Axial Flow Fan, J. of Mech. Eng. 54, p. 159-169. [16] Won, D., Shaw, C.Y. (2003) Determining coefficients for mass-transfer models for volatile organic compound emissions from architectural coatings, Proceedings of the 7th International Conference of Healthy Buildings, 7-11 Dec. 2003, p. 1-6.