Application of Optical Flow Analysis to Shadowgraph Images of Impinging Jet

Masato Hijikuro, Masayuki Anyoji

Interdisciplinary Graduate School of Engineering Sciences, Kyushu University, Fukuoka, Japan.

DOI: 10.4236/jfcmv.2020.84011                

Abstract

In this study, we apply the optical flow method to the time-series shadowgraph images of impinging jets using a high-speed video camera with high spatial and temporal resolution. This image analysis provides quantitative velocity vector fields in the object space without tracer particles. The analysis results clearly capture the details of the coherent vortex structure and its advection from the shear layer of the free jet. Although the results still leave challenges for the quantitative validation, the results show that this analysis method is effective for understanding the details of the physical phenomenon based on the quantitative values extracted from the shadowgraph images.

Keywords

Impinging Jet, Optical Flow, Shadowgraph, Image Analysis, Velocity Vector

1. Introduction

Particle Image Velocimetry (PIV) [1] [2] and Molecular Tagging Velocimetry (MTV) [3] have been widely employed for optical, quantitative, and non-contact surface or spatial fluid analysis. In particular, the development of three-dimensional velocity distribution measurement using stereo PIV has been remarkable with the rapid evolution of cameras and PC performance [4]. However, in the fluid analysis method using tracer particles, fluid analysis methods such as PIV using tracer particles have some problems, i.e. traceability of the tracer particles at high speeds accompanied by shock waves. As conventional qualitative optical visualization methods, the shadowgraph and schlieren optical visualization methods without tracer particles have been applied to understand high-speed flow phenomena with large changes in density. In contrast, Liu et al. [5] devised an algorithm for quantitative image analysis using the optical flow method [6] [7] [8] [9] [10]. The image analysis technique that calculates the moving velocity of an object on the image plane, so-called optical flow, based on the change in the brightness between consecutive images, is defined as the optical flow method or optical flow analysis. Besides, they derived that the constraint formula of optical flow and the theoretical equations of shadowgraph or schlieren methods were similar and revealed the physical meaning of applying this algorithm to these qualitative visualization methods. The feature of this algorithm is that this method enables pixel-level analysis within a captured image and allows to calculate the object moving velocity of objects on the image at the same spatial resolution as the acquired image resolution. Indeed, Liu et al. [11] applied this algorithm to particle images of impinging jets and compared the accuracy with the cross- correlation method used in PIV. Their results indicate that the flow field closer to the theoretical solution can be extracted in the region where the complex interaction occurs between the vortex and the boundary layer than the PIV result by the cross-correlation method. However, their test only applied to low-speed flows with a flow velocity of 4.5 m/s. Furthermore, as far as we have searched, there were no examples of applying the optical flow method in the high-speed flow visualized by a shadowgraph method without tracer particles. When applying this algorithm to high-speed flow phenomena with compressive effects, the camera frame rate and image resolution relative to the velocity of the flow phenomenon are important factors that determine the effectiveness and accuracy of the analysis.

In this study, we performed shadowgraph visualizations of impinging jet in transonic flows by applying the algorithm proposed by Liu et al. Numerous studies have revealed the typical flow field of the impinging jet as shown in Figure 1 [12]. For a two-dimensional nozzle, the two-dimensionality of the phenomenon and the flow symmetry in the time-averaged flow field are assumed. In general, the flow field is roughly divided into three areas: i) a free jet region [13], ii) an impingement or stagnation region [14], and iii) a wall jet region [15]. In the free jet region, primary vortices with counterclockwise rotation are formed by the Kelvin–Helmholtz [KH] instability [16] of the shear layer and advect along the shear layer. In the wall jet region, secondary vortices with clockwise rotation are formed near the wall. This study aims to extract a quantitative velocity vector field by the image analysis for these well-known characteristic vortices captured by the shadowgraph method. The image analysis of the high-speed impinging jet requires both high spatial and temporal resolution for shadowgraph images. We carried out the image analysis of the unstable and complicated high-speed impinging jet using FASTCAM SA-Z, which meets both requirements.

2. Theory

2.1. Projected-Motion Equation for Shadowgraph Image and Optical Flow Equation

The projected-motion equations and physics-based optical flow equation were derived by Liu and Shen [5]. Figure 2 illustrates the perspective projection from

Figure 1. Typical flow field of impinging jet.

Figure 2. Projection from fluid flow on image plane [5].

a fluid medium on an image plane. The perspective center is associated with a camera/lens system. The x₁-axis and x₂-axis are defined as the image coordinate system. The object coordinates X = (X₁, X₂, X₃) is denoted as the projections of the object-space position vector from the perspective center in this frame.

The projected-motion equation for shadowgraph visualization method is expressed by the following equation as a function of the second-order derivative of the fluid density ρ [17]:

$\frac{I-I_T}{I_T}=C{\displaystyle {\int }_{{\Gamma }_1}^{{\Gamma }_2}{\nabla }_{12}^2\rho \text{d}{X}_3}.$ (1)

where, I is the image intensity, I_T is the initial image intensity, C is a coefficient related to the setting of the shadowgraph system, and ${\nabla }_{12}^2$ is the Laplace operator expressed by Equation (2):

${\nabla }_{12}^2=\frac{{\partial }^2}{\partial X_1^2}+\frac{{\partial }^2}{\partial X_2^2}.$ (2)

The visualized domain is confirmed by two control surfaces, $X_3={\Gamma }_1=const.$, $X_3={\Gamma }_2=const.$, corresponding to optical windows. Taking partial differentiation with respect to time in Equation (1) and using the continuity equation expressed by Equation (3), we obtain Equation (4):

$\frac{\partial \rho }{\partial t}+\nabla \cdot \left(\rho U\right)=0.$ (3)

$-\frac{1}{C}\frac{\partial \left(I/I_T-1\right)}{\partial t}={\nabla }_{12}^2\left[{\nabla }_{12}\cdot {\displaystyle {\int }_{{\Gamma }_1}^{{\Gamma }_2}\rho U_{12}\text{d}{X}_3}\right].$ (4)

where, U₁₂ = (U₁, U₂) are the projected two-dimensional velocity components on the coordinate plane (X₁, X₂) of the fluid.

From Equation (1), the Poisson equation for the integration of the fluid density is obtained:

${\nabla }_{12}^2{\displaystyle {\int }_{{\Gamma }_1}^{{\Gamma }_2}\rho \text{d}{X}_3}=\frac{1}{C}\left(I/I_T-1\right).$ (5)

By putting $g=C{\displaystyle {\int }_{{\Gamma }_1}^{{\Gamma }_2}\rho \text{d}{X}_3}$, the density integral can be obtained by solving the Poisson equation ${\nabla }_{12}^{2}g=I/I_T-1$ with appropriate boundary conditions. From Equation (4), we have

$\frac{\partial g}{\partial t}+C\cdot {\nabla }_{12}\cdot {\displaystyle {\int }_{{\Gamma }_1}^{{\Gamma }_2}\rho U_{12}\text{d}{X}_3}=0.$ (6)

The path-averaged velocity weighted with the fluid density is defined as

${\langle U_{12}\rangle }_{\rho }=\frac{{\displaystyle {\int }_{{\Gamma }_1}^{{\Gamma }_2}\rho U_{12}\text{d}{X}_3}}{{\displaystyle {\int }_{{\Gamma }_1}^{{\Gamma }_2}\rho \text{d}{X}_3}}.$ (7)

Substituting Equation (7) into Equation (6), we have

$\frac{\partial g}{\partial t}+{\nabla }_{12}\cdot \left(g{\langle U_{12}\rangle }_{\rho }\right)=0.$ (8)

We note that Horn and Schunck [18] originally suggested the brightness constraint equations $\partial g/\partial t+u\cdot \nabla g=0$ determine the visual motion from a time sequence of images. Comparing Equation (8) with the brightness constraint equations, only for $\nabla \cdot u=0$, Equation (8) has the same form as the Horn-Schunck optical flow equation. However, in general case, $\nabla \cdot u\ne 0$.

2.2. Variational Formulation with Smoothness Constraint

To solve the above optical flow problem, a variational formulation with a smoothness constraint has been proposed by Horn and Schunck [18], which is the first-order form of Tikhonov’s formulation for so-called ill-posed problems (Tikhonov and Arsenin [19]). Corpetti et al. [10] proposed a regularization functional based on the gradients of divergence and curl to preserve the spatial characteristics of divergence and vorticity in flows that is conceptually more plausible to preserve fine structures in turbulent flows. However, these forms of regulation are not derived based on the principles of fluid mechanics. Liu and Shen [5] suggested a physics-based constraint for optical flow computation by projecting the Navier-Stokes equation on the image plane. Their method is to minimize a functional J(u) defined by the following equation:

$J\left(u\right)={\displaystyle {\int }_A{\left[\frac{\partial g}{\partial t}+\nabla \cdot \left(g{\langle U_{12}\rangle }_{\rho }\right)\right]}^2\text{d}{x}_1\text{d}{x}_2}+\alpha {\displaystyle {\int }_A\left({\left|\nabla u_1\right|}^2+{\left|\nabla u_2\right|}^2\right)\text{d}{x}_1\text{d}{x}_2}.$ (9)

where α is the Lagrange multiplier. To minimize J(u), by introducing an arbitrary smooth function $v=\left(v_1,v_2\right)$, computing ${\text{d}J\left({\langle U_{12}\rangle }_{\rho }+pv\right)/\text{d}p|}_{p=0}$, and using the Green’s theorem where the Neumann condition $\partial {\langle U_{12}\rangle }_{\rho }/\partial n=0$ is imposed on the domain boundary $\partial A$, the Euler-Lagrange equation is given as follows:

$g\nabla \left[\frac{\partial g}{\partial t}+\nabla \cdot \left(g{\langle U_{12}\rangle }_{\rho }\right)\right]+\alpha {\nabla }^2{\langle U_{12}\rangle }_{\rho }=0.$ (10)

To solve Equation (10), Jacobi’s block-wise iteration is used. The solution of Horn & Schunck’s equation is used as an initial approximation for Equation (10) for faster convergence.

3. Experimental Setup and Condition

3.1. Impinging Jet Device

Figure 3 illustrates an overall view of the experimental equipment. The working fluid is air, which is filled into a high-pressure gas tank with a volume of 5.0 m³. The compressed air in the tank is released to the atmosphere from the nozzle through a high-pressure pipe, a manual valve, and a plenum chamber. The stagnation pressure in the plenum chamber is measured by a pressure transducer (JTEKT PMS-5M-21M). Atmospheric pressure is measured by a high-precision barometer (GE Druck, DPI740). Figure 4 shows the schematic image of the impinging jet system. The diameter of the nozzle outlet denoted by D is 8 mm. The nozzlehasa two-dimensional shape with a span length of 80 mm and an aspect ratio to the height of the nozzle outlet of 10, and side walls with optical windows are installed on both sides. Despite slight wall interference, this arrangement allows us to carry out the two-dimensional test. The nozzle is made of stereolithography resin, and the surface has a smooth finish. The collision wall was composed of a 300 mm × 80 mm × 15 mm iron plate and a 5 mm thick aluminum plate to improve the structural strength. The collision wall surface was also polished to reduce the effect of the wall surface roughness on the flow field. The distance from the nozzle outlet to the aluminum plate surface was set to 34 mm [11]. In order to evaluate the quantitative validity of the image analysis, an unsteady wall-pressure measurement was performed using a semiconductor pressure sensor (Kulite, XCQ-062-255G) at (x, y) = (15, 6) where the characteristic vortex shedding was observed. Conclusively, the correlation between the frequency of the wall pressure and the one of the vortex shedding based on the brightness change on the image was evaluated.

Figure 3. Experimental equipment.

Figure 4. Impinging jet system.

3.2. Optical Setup

Figure 5 shows the optical setup for the shadowgraph method. The light source was a metal halide lamp (KYOWA, MID-25FC) with a maximum output of 250 W. Two concave mirrors with a diameter of 150 mm and a focal length of 1000 mm were used. Time series images were captured by a high-speed video camera (Photron, FASTCAM SA-Z) with a zoom lens (Nikon, EDAF NIKKOR 80-200).

3.3. Experimental Condition

The nozzle pressure ratio (NPR) is defined by Equation (13) as the airflow condition:

$\text{NPR}=\frac{P_0}{P_b}.$ (13)

where P₀ is the stagnation pressure and P_b is the atmospheric pressure.

The NPR was fixed at 1.5, corresponding to the Mach number at the nozzle exit of 0.78 and the flow velocity of 255 m/s. The frame rate was 300 kfps with an exposure time of 1.0 μsec and a bit count of 8 bits. The image resolution was 256 × 128 pixels, equivalent to the actual spatial resolution of 0.14 mm/pixel.

3.4. Image Analysis Method

Considering the vertical symmetry of the phenomenon with respect to the mainstream axis, the image measurement area was limited to the upper half area including the nozzle outlet and a part of the wall surface (area enclosed by the dash line in Figure 4) to measure with higher spatial resolution. The image analysis area was further reduced by masking a part. The masking was performed on the following two areas: 1) a region in the jet where the main flow velocity is much more dominant than the complicated impinging jet and 2) an area where the change in brightness is very small on the shadow graph image due to the low density change in the flow field. In this image analysis, one snapshot solution is obtained from two image pairs with a certain time interval (d_t) between the two images. The d_t was set to the minimum time interval of 3.33 μsec.

Figure 5. Shadowgraph optical system.

4. Experimental and Analysis Results

4.1. Shadowgraph Image

Figure 6 shows the time-series shadowgraph images at NPR = 1.5. The vertical and horizontal axes are normalized by the nozzle diameter (D).The initial time, t = 0 μs, was arbitrarily defined.

As shown in Figure 6, the change in the brightness corresponding to the change in density is particularly noticeable in the free jet region and the wall jet region. We note that the algorithm is based on a change in the image intensity associated with the flow field on the image. This indicates that the optical flow cannot be calculated in a region where there is no change in brightness due to partial saturation in the schlieren image (i.e. vortex shedding). Hence, shadowgraph images were used for analysis.

Figure 6(a) shows that, in the shear layer of the free jet, a vortex due to the KH instability [16] is generated of the shear layer at around x/D = 2.0. This vortex develops into a coherent structure as it approaches the wall (Figure 6(b) and Figure 6(c)), and then its traveling direction is changed along the wall in the impinging jet region (Figure 6(d)). This study focuses on the advection of the coherent KH vortex generated from the shear layer, and we extract quantitative velocity fields by applying the image analysis using optical flow to time-series shadowgraph images.

4.2. Image Analysis Results

Figures 7-10 show velocity vectors, streamlines, velocity contours, and vorticity contours, respectively, obtained by applying the image analysis using optical flow to the time-series shadowgraph images. Some vectors and streamlines are intentionally decimated for viewability. The velocity contour shown in Figure 9 is normalized by the nozzle exit velocity (U₀ = 255 m/s).

Figure 7 and Figure 8 very clearly capture the spiral structure of the KH vortex and the counterclockwise rotation of the vortex. Also, the speed of rotation gradually increases with the development of the vortex as the center of the vortex approaches the wall. In Figure 7(c) and Figure 7(d), it is observed that the flow in the shear layer of the free jet is entrained and accelerated by the rotation of the KH vortex. Figure 9 shows the accelerated flow field more significantly. The velocity in the vicinity of the center of the KH vortex is relatively higher than that around it. We note that the velocity on the right side of the vortex is higher than that of the spiral center compared to Figure 7. In Figure 9(d), the accelerated flow due to the confluence of the KH vortex and the entrained flow becomes more prominent. The red and blue areas in Figure 10 indicate clockwise and counterclockwise vortex, respectively. Alternating clockwise and counterclockwise vortices appear near the shear layer of the free jet, corresponding small vortices due to the KH instability of the shear layer. We note that no clear secondary vortices are observed. We consider that the secondary vortices can be captured by expanding the visualization area. The velocity distribution around the KH vortex and the small vortices associated with the KH instability of the shear layer are reasonable compared to past studies [7] [20].

(a) t = 0 μsec(b) t = 33 μsec(c) t = 66 μsec(d) t = 99 μsec

Figure 6. Time-series shadowgraph images.

(a) t = 0 μsec(b) t = 33 μsec(c) t = 66 μsec(d) t = 99 μsec

Figure 7. Velocity vector fields.

(a) t = 0 μsec(b) t = 33 μsec(c) t = 66 μsec(d) t = 99 μsec

Figure 8. Streamlines.

(a) t = 0 μsec(b) t = 33 μsec(c) t = 66 μsec(d) t = 99 μsec

Figure 9. Velocity contour.

(a) t = 0 μsec(b) t = 33 μsec(c) t = 66 μsec(d) t = 99 μsec

Figure 10. Vorticity contour.

Figure 11 compares the analysis results of the spectral of the wall pressure fluctuation of the jet and the brightness change on the shadowgraph image intensity at the same location. The horizontal axis of Figure 11 is the frequency [Hz], the left vertical axis shows the wall pressure PSD (power spectrum density) [kPa²/Hz], and the right vertical axis indicates the PSD [count²/Hz] of the count value of the camera corresponding to the brightness on the image. The peak frequencies of the vortex shedding based on the wall pressure and image analysis are 7162 [Hz] and 6800 [Hz], respectively. There is a good correlation between both, indicating the quantitative validity of the image analysis.

5. Conclusions

In order to extract quantitative velocity vectors by image analysis without using tracer particles, we applied the optical flow analysis based on the physical relationship between the density change and the brightness change to the time-series shadowgraph images of the impinging jet. We carried out flow visualizations of transonic impinging jet using the high-speed video camera with high spatial and temporal resolution. Since the impinging jet system has a two-dimensional nozzle with a high aspect ratio and side walls, two-dimensionality of the phenomenon is assumed.

We succeeded in capturing the shear layer fluctuations due to KH instability, the detailed spiral structure of the KH vortices, and complicated interference flows near the collision wall by extracting the quantitative velocity vectors with high spatial resolution. The vector fields and streamlines allow us to clearly define the vortex center that is difficult to understand from qualitative shadowgraph images. Besides, we observed that the flow in the shear layer of the free jet is entrained and accelerated by the rotation of the KH vortex. In addition, there was a correlation between the spectral analysis results of both the wall pressure fluctuation and the brightness value fluctuation of the shadowgraph image.

Although there is a still problem in quantitative validation by comparison with PIV measurements and numerical calculations, the analysis results are reasonable compared with previous studies. Our results indicate that this analysis

Figure 11. Spectrum.

method is effective for quantitatively extracting velocity fields and understanding flow physics from shadowgraph images.

Acknowledgements

We would like to sincerely thank T. Liu of Western Michigan University for his many suggestions on this experiment and the analytical method.

Nomenclature

ρ: Fluid density

I: Image intensity

I_T: Initial image intensity

C: Shadowgraph system coefficient

α: Lagrange multiplier

u: Optical flow velocity

U: Local instantaneous velocity

U₀: Nozzle outlet velocity

D: Nozzle outlet width

Ω: Vorticity

P₀: Stagnation pressure

P_b: Atmospheric pressure

X: Object coordinate

x: Horizontal coordinate

y: Vertical coordinate

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1]	Nayler, J.L. and Frazer, B.A. (1917) Preliminary Report upon an Experimental Method of Investigating, by the Aid of Kinematographic Photography, the History of Eddying Flow Past a Model Immersed in Water. Technical Report of the Advisory Committee for Aeronautics for 1917-18.
[2]	Adrian, R.J. (1991) Particle-Imaging Techniques for Experimental Fluid Mechanics. Annual Review of Fluid Mechanics, 23, 261-304. https://doi.org/10.1146/annurev.fl.23.010191.001401
[3]	Handa, T., Mizuta, S. and Imamura, K. (2012) Velocity Measurements in Gaseous Flows Using MTV: Velocity Measurements in Supersonic Microjets. Journal of the Visualization Society of Japan, 32, 76-81. https://doi.org/10.3154/jvs.32.26
[4]	Prasad, A.K. and Adrian, R.J. (1993) Stereoscopic Particle Image Velocimetry Applied to Liquid Flows. Experiments in Fluids, 15, 49-60. https://doi.org/10.1007/BF00195595
[5]	Liu, T. and Shen, L. (2008) Fluid Flow and Optical Flow. Journal of Fluid Mechanics, 614, 253-291. https://doi.org/10.1017/S0022112008003273
[6]	Quenot, G.M., Pakleza, J. and Kowalewski, T.A. (1998) Particle Image Velocimetry with Optical Flow. Experiments in Fluids, 25, 177-189. https://doi.org/10.1007/s003480050222
[7]	Ruhnau, P., et al. (2005) Variational Optical Flow Estimation for Particle Image Velocimetry. Experiments in Fluids, 38, 21-32. https://doi.org/10.1007/s00348-004-0880-5
[8]	Yuan, J., Schnorr, C. and Memin, E. (2007) Discrete Orthogonal Decomposition and Variational Fluid Flow Estimation. Journal of Mathematical Imaging and Vision, 28, 67-80. https://doi.org/10.1007/s10851-007-0014-9
[9]	Corpetti, T., Memin, E. and Perez, P. (2002) Dense Estimation of Fluid Flows. IEEE Transactions on Pattern Analysis and Machine Intelligence, 24, 365-380. https://doi.org/10.1109/34.990137
[10]	Corpetti, T., et al. (2006) Fluid Experimental Flow Estimation Based on an Optical Flow Scheme. Experiments in Fluids, 40, 80-97. https://doi.org/10.1007/s00348-005-0048-y
[11]	Liu, T. and Merat, A. (2015) Comparison between Optical Flow and Cross-Correlation Methods for Extraction of Velocity Fields from Particle Images. Experiments in Fluids, 56, 166-189. https://doi.org/10.1007/s00348-015-2036-1
[12]	Jambunathan, K., et al. (1992) A Review of Heat Transfer Data for Single Circular Jet Impingement. International Journal of Heat and Fluid Flow, 13, 106-115. https://doi.org/10.1016/0142-727X(92)90017-4
[13]	Yule, A.J. (1978) Large-Scale Structure in the Mixing Layer of a Round Jet. Journal of Fluid Mechanics, 89, 413-432. https://doi.org/10.1017/S0022112078002670
[14]	Tummers, M.J., Jacobse, J. and Voorbrood, S.G. (2011) Turbulent Flow in the Near Field of a Round Impinging Jet. International Journal of Heat and Mass Transfer, 54, 4939-4948. https://doi.org/10.1016/j.ijheatmasstransfer.2011.07.007
[15]	Van Hout, R., Rinsly, V. and Grobman, Y.G. (2018) Experimental Study of a Round Jet Impinging on a Flat Surface: Flow Field and Vortex Characteristics in the Wall Jet. International Journal of Heat and Fluid Flow, 70, 41-58. https://doi.org/10.1016/j.ijheatfluidflow.2018.01.010
[16]	Greco, C.S., et al. (2018) Effects of the Stroke Length and Nozzle-to-Plate Distance on Synthetic Jet Impingement Heat Transfer. International Journal of Heat and Mass Transfer, 117, 1019-1031. https://doi.org/10.1016/j.ijheatmasstransfer.2017.09.118
[17]	Goldstein, R.J. and Kuehn, T.H. (1996) Optical Systems for Flow Measurement: Shadowgraph, Schlieren, and Interferometric Techniques. In: Fluid Mechanics Measurement, Second Edition, Taylor & Francis, Abingdon-on-Thames, Chapter 7.
[18]	Horn, B.K. and Schunck, B.G. (1981) Determining Optical Flow. Artificial Intelligence, 17, 185-204. https://doi.org/10.1016/0004-3702(81)90024-2
[19]	Tikhonov, A.N. and Arsenin, V.Y. (1977) Solutions of Ill-Posed Problems. Wiley, Hoboken, Chapter 2.
[20]	Arthurs, D. and Ziada, S. (2012) Self-Excited Oscillations of a High-Speed Impinging Planar Jet. Journal of Fluids and Structures, 34, 236. https://doi.org/10.1016/j.jfluidstructs.2012.06.002