**Three Dimensional Simulation of Gas-Radiation Interactions in Gas Lasers**

Timothy J. Madden

*US Air Force Research Laboratory, Directed Energy Directorate Kirtland Air Force Base, New Mexico USA* 

#### **1. Introduction**

174 Fluid Dynamics, Computational Modeling and Applications

[8] Beckhaus, K. and Larisch, M., Alehossein, H., Northey, S., Ney, P., Lucas, G.,

[9] Beckhaus, K. and Larisch, M., Alehossein, H., Northey, S., Ney, P., Lucas, G.,

[10] Alehossein, H. 2009, "A triangular model of caving subsidence", *J. Applied Earth Science* 

[11] Shen, B. and Alehossein, H. 2009. *ACARP Project* C16023, Australia – *Subsidence Control* 

[15] Fox, R.W., Pritchard, P.J. and McDonald, A.T. 2009. *Introduction to fluid mechanics*, Don

[16] Alehossein, H. and Hood, M. 1999, "Application of linearized dimensional analysis to

rock cutting", *International Journal of Rock Mechanics and Mining Sciences*, V. 36, 1999,

on "Subsidence control using overburden grout injection technology". [12] Alehossein, H., Poulsen, B. A. 2010, "Stress analysis of longwall top coal caving". *International Journal of Rock Mechanics and Mining Sciences*, 47(1), 30-1. [13] Alehossein, H. 2010, "Mechanics of slurries for rock factures & voids", WCHFA-2010,

foundations", ICAGE conference, Perth 2011.

*(Trans. Inst. Min. Metall. B)*, Vol 118, No 1, pp. 1-4.

[14] R.F. Craig, 1997. *Soil Mechanics*. E & FN SPON, Routledge. New York.

2011, Leipzig, Germany

August 2010, Xian, China

pp. 701–709.

Fowley (Wiley & Sons), New York.

Vanderstaay L. 2011, "New approach for tremie concrete used for deep

Vanderstaay L. 2011, "Performance-based specifications for concrete, 14 – 15 June

*Using Coal Washery Waste*. An extension research project to ACARP C12019 in 2003

The spectroscopically measured lineshape of an atomic transition provides a wealth of useful information relative to diagnosing the state of a gas. The center of the lineshape is specific to a particular transition of a specific atom. The width of the lineshape indicates the amount of broadening of the transition, due to the effects of both collisions with other particles in the gas and Doppler shift due to the movement of the atom. Since the Doppler shift is proportional to velocity, the width of the transition can be used to estimate the degree of random molecular motion in the gas, expressed macroscopically as temperature. A Doppler shift to the frequencies in the transition can also occur through the bulk motion of the gas, and this can be used to examine the velocity field of the gas. The astronomy community was the first to recognize the utility of these concepts in practical applications, stemming to the early 1930s.1,2 In 1934 Stuve and Elvey3 showed that by including a bulk gas velocity Doppler broadening term in the Voigt equation for the transition lineshape in addition to the random thermal motion term, it was possible to estimate whether a stellar atmosphere was 'turbulent' or not based upon the fit of the Voigt equation to the measured transition lineshapes for that atmosphere. Taking the theory a step further, using estimated optical paths or length scales for the stellar atmospheres that they were measuring, they were able to estimate median gas velocities and correlate increasing velocity magnitude with increasing temperature of the stellar atmosphere.

More recently, continuously tunable diode lasers have been applied to lineshape measurement of transitions within species in the chemical oxygen-iodine laser (COIL) flowfield as a means to determine number density and laser gain on the I 2P1/2→2P3/2 transition. Davis and Allen et al4,5,6,7 applied lineshape measurement and Voigt fitting methods to various COIL species as a mechanism to determine concentration and translational temperature in the COIL flowfield as an experiment diagnostic. While the effect of the bulk gas velocity was not taken into account in these investigations, the diagnostics developed by these investigations provided the means to do so. Nikolaev et al,8 applied the same type of continuously tunable diode laser to the I 2P1/2→2P3/2 transition to investigate the influence of mean gas velocity through a COIL device upon the laser gain. They showed that by varying the angle at which the diode laser beam passes through the COIL flowfield, the line center gain can be varied through bulk velocity induced Doppler broadening. Using this mechanism, they were able to determine the mean flow velocity for the COIL and convert the laser gain

Three Dimensional Simulation of Gas-Radiation Interactions in Gas Lasers 177

<sup>1</sup> 3 \* *Oa I OX I* 2 22 2

 1\* 3 2 <sup>2</sup> 2 2 3/2 *Oa I OX I P* 2

12 32 *Oa I P OX I P* <sup>2</sup> 3/2 2 1/2

The final product of the COIL chemical kinetics processes is the 2P1/2 electronically excited state of atomic iodine, a second population inversion in atomic iodine that results from the first in O2. The flow unsteadiness associated with the jet-in-crossflow interaction of the secondary He/I2 flow with the primary O2(a1) carrying flow leads to the generation of vortices that are transported downstream into the supersonic region where the laser cavity is located. These vortices through their rotation have associated velocity components that are parallel to the transit of the optical field. The motion of the bulk flow in the direction of the optical field leads to Doppler broadening of the atomic iodine (I) 2P1/2→2P3/2 transition lineshape. This effect is important in two ways: the first is the decrease in the line center gain of the 2P1/2→2P3/2 transition and the second is in the breadth of the lineshape. The first effect is important relative to the performance of the laser, but the second effect is most important relative to the application of spectroscopic lineshape measurement. Understanding the influence of the vortices upon the breadth of the I 2P1/2→2P3/2 transition lineshape is therefore important to understanding the accuracy of the temperature diagnostic in COIL experiments. Additionally, the magnitude of the flow fluctuations, their temporal variation, and their scale all are related to the question of whether the flow is fully turbulent or not. Relative to the modeling of the gas flows through COIL devices, these flows can best be described as the translation of particles of different chemical composition with collisional interactions occurring between the particles and between the particles and the photons within the radiation field. Mathematically, this flow of particles is treated as a continuum and is approximated by Navier-Stokes continuity equations for mass, momentum, and

( ) ˆ

The production rates of the species due to chemical reactions are determined for each species continuity equation using the Arrhenius rate law for the forward rates, with the backward rates determined from the equilibrium and forward rates. The computational fluid dynamics (CFD) code GASP from AeroSoft, Inc. is employed to solve the Navier Stokes equations for the simulations performed here. The GASP COIL model utilizes 10 to 13 species mass conservation equations for the chemically reacting components of the COIL flow, depending upon the COIL model used, in addition to the base conservation equations for momentum and energy. An effective binary diffusion model is used to describe concentration and pressure contributions to mass diffusion, an important process in the low density COIL flowfield. 13 species, 52 reactions and 10 species, 20 reactions finite-rate chemistry mechanisms11 are used to model the gas phase chemical kinetic processes that generate the population inversion in atomic iodine necessary for laser oscillation in

*Q dvol F F ndA S dvol*

 *<sup>v</sup> vol A vol*

 2 2 1/2 2 3/2 2 *I P HO I P HO*

energy. In integral form, these equations are given by:

*t*

diagnostic into a flow diagnostic. Similarly, they recognized the propensity for 'turbulence' to contribute to the broadening of the lineshape. Invoking an assumption of isotropic turbulence to characterize the bulk vorticity in the flows that they were measuring, they examined the lineshape data using basically the same theory developed by Stuve and Elvey and estimated a 10 K effect for a mean turbulent fluctuation velocity of 30 m/s. Given experimental evidence that temperatures determined from experimentally measured lineshapes for the I 2P1/2→2P3/2 transition in COILs may be as much as 50 K higher than what is expected based upon thermodynamic predictions, the implication is that either the magnitude of the flow fluctuations is much larger than 30 m/s or some other process such as H2O condensation is liberating heat and increasing the temperature.

The work presented here examines the issue of the interplay between an optical field in a COIL and the bulk motion of the gas through simulation of the physics. Previous work by Madden and Miller9 presented data from COIL simulations indicating that the flowfield should be marked by substantial flow unsteadiness, as documented in the experiments of Fric and Roshko10 for fluid dynamically similar types of flows. These simulations provide the capability to examine the optical field to gas interaction in the lineshape directly, including the full physics of the fluid dynamics, chemistry, heat release, and flow instability without resorting to statistical turbulence theory. Direct prediction of the fluid dynamic physics allows for examination of the magnitude, and spatial and temporal characteristics of the flow fluctuations, an important consideration given the presence of both strong compressibility which stabilizes and heat release which may either stabilize or de-stabilize the flow. From this perspective, it is not only possible to characterize the optical field to gas interaction through the I 2P1/2→2P3/2 transition lineshape, but also to assess whether the flow is actually turbulent or not.

#### **2. Problem and methodology**

Gas lasers operate through a series of processes, often exothermic, that lead to a population inversion, or non-Boltzmann energy distribution, within the energy states of the species in the gas. In COILs, two separate but coupled non-Boltzmann energy distributions are created. The first is through the liquid phase reaction of chlorine with the hydroxyl ion generating a population inversion within oxygen with the electronically excited state of singlet delta oxygen:

$$\text{Cl}\_2 + 2\text{HO}\_2^- \rightarrow \text{O}\_2\text{(}a^1\Delta\text{)} + \text{H}\_2\text{O}\_2 + 2\text{Cl}^-$$

The O2(a1) subsequently diffuses out of the liquid phase into a carrier gas, He, and is transported toward the mixing and reaction region upstream of the resonator. During transport, energy pooling processes erode the O2(a1) population through pooling and deactivation. Further downstream, a mixture of molecular iodine, I2, and a diluent gas, He, is injected into the O2(a1) carrying flow and reacts with the singlet oxygen. The combination of these pooling and reaction processes can be briefly summarized as:

$$O\_2\left(a^1\Delta\right) + O\_2\left(a^1\Delta\right) \to O\_2\left(b^1\Sigma, \upsilon = 1, 2\right) + O\_2\left(X^3\Sigma\right)$$

$$O\_2\left(b^1\Sigma, \upsilon = 1, 2\right) + H\_2O\left(a^1\Delta, \upsilon\right) \to H\_2O$$

diagnostic into a flow diagnostic. Similarly, they recognized the propensity for 'turbulence' to contribute to the broadening of the lineshape. Invoking an assumption of isotropic turbulence to characterize the bulk vorticity in the flows that they were measuring, they examined the lineshape data using basically the same theory developed by Stuve and Elvey and estimated a 10 K effect for a mean turbulent fluctuation velocity of 30 m/s. Given experimental evidence that temperatures determined from experimentally measured lineshapes for the I 2P1/2→2P3/2 transition in COILs may be as much as 50 K higher than what is expected based upon thermodynamic predictions, the implication is that either the magnitude of the flow fluctuations is much larger than 30 m/s or some other process such as H2O condensation is

The work presented here examines the issue of the interplay between an optical field in a COIL and the bulk motion of the gas through simulation of the physics. Previous work by Madden and Miller9 presented data from COIL simulations indicating that the flowfield should be marked by substantial flow unsteadiness, as documented in the experiments of Fric and Roshko10 for fluid dynamically similar types of flows. These simulations provide the capability to examine the optical field to gas interaction in the lineshape directly, including the full physics of the fluid dynamics, chemistry, heat release, and flow instability without resorting to statistical turbulence theory. Direct prediction of the fluid dynamic physics allows for examination of the magnitude, and spatial and temporal characteristics of the flow fluctuations, an important consideration given the presence of both strong compressibility which stabilizes and heat release which may either stabilize or de-stabilize the flow. From this perspective, it is not only possible to characterize the optical field to gas interaction through the I 2P1/2→2P3/2 transition lineshape, but also to assess whether the flow

Gas lasers operate through a series of processes, often exothermic, that lead to a population inversion, or non-Boltzmann energy distribution, within the energy states of the species in the gas. In COILs, two separate but coupled non-Boltzmann energy distributions are created. The first is through the liquid phase reaction of chlorine with the hydroxyl ion generating a population inversion within oxygen with the electronically excited state of

> <sup>1</sup> 2 22 2 2 *Cl HO O a H O Cl* 2 2

 11 1 <sup>3</sup> 22 2 2 *Oa Oa Ob v OX* , 1,2

> 1 1 2 22 *O b v HOa v HO* , 1,2 ,

The O2(a1) subsequently diffuses out of the liquid phase into a carrier gas, He, and is transported toward the mixing and reaction region upstream of the resonator. During transport, energy pooling processes erode the O2(a1) population through pooling and deactivation. Further downstream, a mixture of molecular iodine, I2, and a diluent gas, He, is injected into the O2(a1) carrying flow and reacts with the singlet oxygen. The

combination of these pooling and reaction processes can be briefly summarized as:

liberating heat and increasing the temperature.

is actually turbulent or not.

singlet delta oxygen:

**2. Problem and methodology** 

$$O\_2\left(a^1\Delta\right) + I\_2 \to O\_2\left(X^3\Sigma\right) + I\_2$$

$$O\_2\left(a^1\Delta\right) + I\_2\stackrel{\*}{\to} O\_2\left(X^3\Sigma\right) + 2I\left(^2P\_{3/2}\right)$$

$$O\_2\left(a^1\Delta\right) + I\left(^2P\_{3/2}\right) \to O\_2\left(X^3\Sigma\right) + I\left(^2P\_{1/2}\right)$$

$$I\left(^2P\_{1/2}\right) + H\_2O \to I\left(^2P\_{3/2}\right) + H\_2O$$

The final product of the COIL chemical kinetics processes is the 2P1/2 electronically excited state of atomic iodine, a second population inversion in atomic iodine that results from the first in O2. The flow unsteadiness associated with the jet-in-crossflow interaction of the secondary He/I2 flow with the primary O2(a1) carrying flow leads to the generation of vortices that are transported downstream into the supersonic region where the laser cavity is located. These vortices through their rotation have associated velocity components that are parallel to the transit of the optical field. The motion of the bulk flow in the direction of the optical field leads to Doppler broadening of the atomic iodine (I) 2P1/2→2P3/2 transition lineshape. This effect is important in two ways: the first is the decrease in the line center gain of the 2P1/2→2P3/2 transition and the second is in the breadth of the lineshape. The first effect is important relative to the performance of the laser, but the second effect is most important relative to the application of spectroscopic lineshape measurement. Understanding the influence of the vortices upon the breadth of the I 2P1/2→2P3/2 transition lineshape is therefore important to understanding the accuracy of the temperature diagnostic in COIL experiments. Additionally, the magnitude of the flow fluctuations, their temporal variation, and their scale all are related to the question of whether the flow is fully turbulent or not.

Relative to the modeling of the gas flows through COIL devices, these flows can best be described as the translation of particles of different chemical composition with collisional interactions occurring between the particles and between the particles and the photons within the radiation field. Mathematically, this flow of particles is treated as a continuum and is approximated by Navier-Stokes continuity equations for mass, momentum, and energy. In integral form, these equations are given by:

$$\frac{\partial}{\partial t} \iiint\_{\text{vol}} \{\mathcal{Q}\} \, dvol + \iint\_{A} (\vec{F} - \vec{F}\_{\text{v}}) \cdot \hat{n} dA = \iiint\_{\text{vol}} \{S\} \, dvol$$

The production rates of the species due to chemical reactions are determined for each species continuity equation using the Arrhenius rate law for the forward rates, with the backward rates determined from the equilibrium and forward rates. The computational fluid dynamics (CFD) code GASP from AeroSoft, Inc. is employed to solve the Navier Stokes equations for the simulations performed here. The GASP COIL model utilizes 10 to 13 species mass conservation equations for the chemically reacting components of the COIL flow, depending upon the COIL model used, in addition to the base conservation equations for momentum and energy. An effective binary diffusion model is used to describe concentration and pressure contributions to mass diffusion, an important process in the low density COIL flowfield. 13 species, 52 reactions and 10 species, 20 reactions finite-rate chemistry mechanisms11 are used to model the gas phase chemical kinetic processes that generate the population inversion in atomic iodine necessary for laser oscillation in

Three Dimensional Simulation of Gas-Radiation Interactions in Gas Lasers 179

section with one large and two small injector orifices that issue reactants into the primary flow through this nozzle. It was previously shown15 that this particular configuration of the 'unitdomain' as opposed to one half as wide splitting the large injector orifice is necessary to properly capture unsteady fluctuations in the lateral or Z direction relative to Figure 1. Fig. 1 illustrates the unit-domain computational grid within the context of a surface rendering of the experiment's mixing nozzle. The orifices inject a sonic mixture of He and I2 into a subsonic primary flow composed of He, O2(a1Δ), O2(X3Σ), H2O, and Cl2 inducing the complex 3-D flow structure as the jet issuing from the orifice interacts with the crossflow. It is the combination of the interfacial area created by the complex 3-D flow structure and molecular diffusion that mixes the two flows. Boundary conditions accomplish the unit-domain approximation through the enforcement of planar symmetry at the nozzle centerline in the vertical direction and periodicity at the side boundaries in the lateral direction. No-slip constant temperature boundary conditions are used at the wetted surfaces of the nozzle and orifices, with the temperature fixed at 400 K at the orifice region walls and 300 K at the nozzle walls. The nozzle

Fig. 1. Computational grid used in the 3-D GASP simulations of COIL hardware, shown in relation to a rendering of the wetted surfaces in the COIL laser hardware. Note that the

resolution show here is 1/64th that used in the execution of the simulation.

COILs. These models utilize an empirical construct for I2 dissociation in the presence of O2(a1Δ) developed by Heidner12 that has been utilized within the COIL community for 25 years. Recent work by Azyazov et al13 indicates that the I2 dissociation process is somewhat more complex than the simplified construct developed by Heidner. However, work by Waichman14 indicates that the difference in power predictions from CFD simulations using the new mechanism from Azyazov and the original Heidner mechanism is of order 10% for what is termed 'high pressure operation,' i.e. essentially optimal operation with regard to the I2/O2 ratio, and similar to the conditions considered in this investigation. Form this standpoint, the 10% difference between the I2 dissociation mechanisms for the 'high pressure,' optimal COIL operational conditions simulated here is not expected to substantially alter interpretation of the results given here.

The laser gain and the associated lineshape for the I 2P1/2→2P3/2 transition is determined from the flowfield calculated variables using the equations:

$$\begin{aligned} \label{E} \alpha \left(\nu\_{0}\right) &= \frac{7}{12} \left(\frac{\Delta l\_{0}^{-2}}{8\pi}\right) \phi \left(\nu\_{0}\right) \left(N\_{\mathcal{H}^{-2}\_{\text{V}\_{\text{V}\_{\text{L}}}}} - \frac{1}{2} N\_{\mathcal{H}^{-2}\_{\text{V}\_{\text{V}\_{\text{L}}}}}\right) \\\\ \phi \left(\alpha \nu(\nu)\right) &= \left[ \left(\frac{a^{2} 4 \ln 2}{\pi^{3}}\right)^{\frac{1}{2}} \frac{1}{\Delta \nu\_{D}} \right] \int\_{-\infty}^{\infty} \frac{e^{-\gamma^{2}}}{\left(\alpha - \gamma\right)^{2} + a^{2}} d\nu \\\\ \alpha &= \frac{2 \left(\nu - \nu\_{0} \left(1 + \frac{\mathcal{W}}{c}\right)\right)}{\Delta \nu\_{D}} \sqrt{\ln 2} \\\\ a &= \frac{\Delta \nu\_{L}}{\Delta \nu\_{D}} \sqrt{\ln 2} \\\\ \chi &= \frac{2 \nu\_{\circ} \nu\_{0}}{c \Delta \nu\_{D}} \sqrt{\ln 2} \\\\ \Delta \nu\_{D} &= \frac{2}{\lambda} \sqrt{\frac{2RT \ell \nu \mathcal{Q}}{m\_{I}}} \\\\ \Delta \nu\_{L} &= \frac{T\_{\nu \mathcal{I}}}{T} P \sum\_{i=1}^{N} a\_{i} \mathcal{X}\_{i} \end{aligned}$$

A difference between the approach to including the effect of the bulk flow velocity components on the lineshape used here and that used by Stuve and Elvey and later by Nikolaev is that they are not assumed to be described by a Gaussian probability distribution. Here, the velocity distribution is directly determined through prediction of the velocity components including their scale, magnitude, and temporal variation.

The computational grid used for the COIL simulations consists of 29 blocks and 8.1 million grid cells. The computational domain which this grid discretizes represents the smallest geometrically similar element within the flowfield in the COIL experiment hardware, denoted a 'unit-domain.' The unit-domain consists of a supersonic (M~2.2) converging-diverging nozzle

COILs. These models utilize an empirical construct for I2 dissociation in the presence of O2(a1Δ) developed by Heidner12 that has been utilized within the COIL community for 25 years. Recent work by Azyazov et al13 indicates that the I2 dissociation process is somewhat more complex than the simplified construct developed by Heidner. However, work by Waichman14 indicates that the difference in power predictions from CFD simulations using the new mechanism from Azyazov and the original Heidner mechanism is of order 10% for what is termed 'high pressure operation,' i.e. essentially optimal operation with regard to the I2/O2 ratio, and similar to the conditions considered in this investigation. Form this standpoint, the 10% difference between the I2 dissociation mechanisms for the 'high pressure,' optimal COIL operational conditions simulated here

The laser gain and the associated lineshape for the I 2P1/2→2P3/2 transition is determined

2 2 1 2 3 2

*<sup>A</sup><sup>λ</sup> α ν <sup>ν</sup> N N <sup>π</sup>*

 

<sup>0</sup> 2 1

 

*<sup>D</sup> m*

*ref <sup>L</sup> P T T*

A difference between the approach to including the effect of the bulk flow velocity components on the lineshape used here and that used by Stuve and Elvey and later by Nikolaev is that they are not assumed to be described by a Gaussian probability distribution. Here, the velocity distribution is directly determined through prediction of the

The computational grid used for the COIL simulations consists of 29 blocks and 8.1 million grid cells. The computational domain which this grid discretizes represents the smallest geometrically similar element within the flowfield in the COIL experiment hardware, denoted a 'unit-domain.' The unit-domain consists of a supersonic (M~2.2) converging-diverging nozzle

*<sup>D</sup>*

 

*a*

velocity components including their scale, magnitude, and temporal variation.

*/ / I( P ) I( P )*

<sup>2</sup> <sup>1</sup> 2 2

*W c*

ln 2 *<sup>L</sup> D*

<sup>0</sup> <sup>2</sup> ln 2 *<sup>z</sup> D v <sup>c</sup> y*

2 2*RTn*2

1

*i*

*N*

*I*

 

*i i*

3 2 2 4ln 2 1 *<sup>y</sup> D a e dy*

ln 2

*y a*

7 1 12 8 2 

is not expected to substantially alter interpretation of the results given here.

2 0 0 0

from the flowfield calculated variables using the equations:

 section with one large and two small injector orifices that issue reactants into the primary flow through this nozzle. It was previously shown15 that this particular configuration of the 'unitdomain' as opposed to one half as wide splitting the large injector orifice is necessary to properly capture unsteady fluctuations in the lateral or Z direction relative to Figure 1. Fig. 1 illustrates the unit-domain computational grid within the context of a surface rendering of the experiment's mixing nozzle. The orifices inject a sonic mixture of He and I2 into a subsonic primary flow composed of He, O2(a1Δ), O2(X3Σ), H2O, and Cl2 inducing the complex 3-D flow structure as the jet issuing from the orifice interacts with the crossflow. It is the combination of the interfacial area created by the complex 3-D flow structure and molecular diffusion that mixes the two flows. Boundary conditions accomplish the unit-domain approximation through the enforcement of planar symmetry at the nozzle centerline in the vertical direction and periodicity at the side boundaries in the lateral direction. No-slip constant temperature boundary conditions are used at the wetted surfaces of the nozzle and orifices, with the temperature fixed at 400 K at the orifice region walls and 300 K at the nozzle walls. The nozzle

Fig. 1. Computational grid used in the 3-D GASP simulations of COIL hardware, shown in relation to a rendering of the wetted surfaces in the COIL laser hardware. Note that the resolution show here is 1/64th that used in the execution of the simulation.

Three Dimensional Simulation of Gas-Radiation Interactions in Gas Lasers 181

the more complex structure associated with the break-up. Since molecular diffusion and chemical reactions are strongly correlated with spatial gradients of reactant concentration, the presence of these unsteady flow structures has considerable impact upon the generation and variation of I(2P3/2) and I(2P1/2). Additionally, the gain is directly proportional to the number densities of I(2P3/2) and I(2P1/2) and the gain will correlate with the flow structure as

Fig. 3. Fixed time snap shot of the three-dimensional gain variation within the simulation

The underlying sources of the flow unsteadiness are the 'wake vortices' in the flow behind the He/I2 jets. These vortices have their origin in the primary flow boundary layer fluid that wraps around the orifice flow and projects away from the wall with rotation, subsequently interacting with the jet fluid. Fig. 4 illustrates this with plots of streamlines originating in the boundary layer fluid upstream of the jet shown in conjunction with vortex cores. The cores were extracted using the vortex core identification functions of Sujudi and Haimes18 implemented in the Tecplot visualization software from Amtec Engineering. The red vortex core traces can be found at the center of the swirls in the streamlines, identifying the location of the vortices. Two wake vortices are seen projecting away from the wall immediately behind the large orifice, and a third projects downward from the wall in the center of the domain between the large and small injector orifices. Similarly, additional vortices are seen projecting away from the wall behind the small injector orifices. Temporal variation of the flow fluctuations and the Fourier analysis of these fluctuations for the two simulations is shown in Figs. 5, 6. Fourier analysis of the fluctuations from the 10 species, 22 reaction simulation yielded frequencies of 56, 94, 131, and 150 kHz, whereas the 13 species, 52 reaction simulation yielded frequencies of 12.1, 24.9, 42.3, 60.4, 96.6, 133.9, 151, and 292 kHz. The peak fluctuation amplitudes for the lateral velocities are 175 m/sec in both simulations. The presence of the unsteady fluctuations within the jet structure is consistent with the experimental literature on jet-in-crossflow interactions. In particular, the capture of the wake vortices as well as their location is entirely consistent with the experimental investigations,

illustrated in Fig. 3.

volume.

exemplified in the work of Frik and Rhosko.

and injector subsonic inflow boundary conditions fix the total pressure, total temperature, and the species fractions at constant values, while the derivative of the static pressure is set to 0. The nozzle outflow boundary condition sets the second derivative of the dependent variables to 0 as is appropriate for supersonic flows. For execution of this model, References providing additional details beyond the scope of this article may be found in Madden et al.16,17

#### **3. Results**

The 3-D GASP model for the COIL flowfield was executed in time accurate mode utilizing 1st order accurate time integration and 3rd order spatial accuracy. The computation was advanced with a physical time step of 4.0x10-9 sec over 159,500 time steps to a physical time of 0.00064 sec. Analysis of this simulation indicated the need to increase the resolution of the grid in the streamwise direction within the nozzle expansion, and a second simulation was performed that quadrupled the number of grid cells in this region bringing the total number of grid cells to 14.7 million. This simulation was continued for an additional 124,000 time steps with a time step of 2. 0x10-9 sec., and added 0.000247 sec of physical time to the simulation.

The time accurate execution of the GASP COIL model generated a prediction of flow unsteadiness that did not decay over the time intervals that the computations were advanced. The unsteadiness was found to extend from the jet/primary interaction region at the point of the He/I2 transverse injection and continue downstream undiminished. Fig. 2 illustrates the resulting impact of the flow unsteadiness upon the structure of the flow within the He/I2 jet. Here a volumetric rendering of the I2 mole fraction is plotted within 3- D space. The rendering demonstrates the presence of complex, periodic structures associated with the unsteady vortex generation, and indicates that there is a strong interaction between the fluid from the large and small orifices. The break-up of the jets occurs as the fluid from the small orifices begins to interact with that from the large, illustrated by the change in the vortical structure from smooth, periodic fluctuations to

Fig. 2. Fixed time snap shot of the unsteady fluid dynamic structures manifested in the He/I2 jet from the 10 species, 22 reactions 3-D GASP COIL simulation. Note that the view is reflected about the symmetry plane at the nozzle centerline to provide an enhanced perspective.

and injector subsonic inflow boundary conditions fix the total pressure, total temperature, and the species fractions at constant values, while the derivative of the static pressure is set to 0. The nozzle outflow boundary condition sets the second derivative of the dependent variables to 0 as is appropriate for supersonic flows. For execution of this model, References providing

The 3-D GASP model for the COIL flowfield was executed in time accurate mode utilizing 1st order accurate time integration and 3rd order spatial accuracy. The computation was advanced with a physical time step of 4.0x10-9 sec over 159,500 time steps to a physical time of 0.00064 sec. Analysis of this simulation indicated the need to increase the resolution of the grid in the streamwise direction within the nozzle expansion, and a second simulation was performed that quadrupled the number of grid cells in this region bringing the total number of grid cells to 14.7 million. This simulation was continued for an additional 124,000 time steps with a time

The time accurate execution of the GASP COIL model generated a prediction of flow unsteadiness that did not decay over the time intervals that the computations were advanced. The unsteadiness was found to extend from the jet/primary interaction region at the point of the He/I2 transverse injection and continue downstream undiminished. Fig. 2 illustrates the resulting impact of the flow unsteadiness upon the structure of the flow within the He/I2 jet. Here a volumetric rendering of the I2 mole fraction is plotted within 3- D space. The rendering demonstrates the presence of complex, periodic structures associated with the unsteady vortex generation, and indicates that there is a strong interaction between the fluid from the large and small orifices. The break-up of the jets occurs as the fluid from the small orifices begins to interact with that from the large, illustrated by the change in the vortical structure from smooth, periodic fluctuations to

Fig. 2. Fixed time snap shot of the unsteady fluid dynamic structures manifested in the He/I2 jet from the 10 species, 22 reactions 3-D GASP COIL simulation. Note that the view is reflected about the symmetry plane at the nozzle centerline to provide an enhanced perspective.

additional details beyond the scope of this article may be found in Madden et al.16,17

step of 2. 0x10-9 sec., and added 0.000247 sec of physical time to the simulation.

**3. Results** 

the more complex structure associated with the break-up. Since molecular diffusion and chemical reactions are strongly correlated with spatial gradients of reactant concentration, the presence of these unsteady flow structures has considerable impact upon the generation and variation of I(2P3/2) and I(2P1/2). Additionally, the gain is directly proportional to the number densities of I(2P3/2) and I(2P1/2) and the gain will correlate with the flow structure as illustrated in Fig. 3.

Fig. 3. Fixed time snap shot of the three-dimensional gain variation within the simulation volume.

The underlying sources of the flow unsteadiness are the 'wake vortices' in the flow behind the He/I2 jets. These vortices have their origin in the primary flow boundary layer fluid that wraps around the orifice flow and projects away from the wall with rotation, subsequently interacting with the jet fluid. Fig. 4 illustrates this with plots of streamlines originating in the boundary layer fluid upstream of the jet shown in conjunction with vortex cores. The cores were extracted using the vortex core identification functions of Sujudi and Haimes18 implemented in the Tecplot visualization software from Amtec Engineering. The red vortex core traces can be found at the center of the swirls in the streamlines, identifying the location of the vortices. Two wake vortices are seen projecting away from the wall immediately behind the large orifice, and a third projects downward from the wall in the center of the domain between the large and small injector orifices. Similarly, additional vortices are seen projecting away from the wall behind the small injector orifices. Temporal variation of the flow fluctuations and the Fourier analysis of these fluctuations for the two simulations is shown in Figs. 5, 6. Fourier analysis of the fluctuations from the 10 species, 22 reaction simulation yielded frequencies of 56, 94, 131, and 150 kHz, whereas the 13 species, 52 reaction simulation yielded frequencies of 12.1, 24.9, 42.3, 60.4, 96.6, 133.9, 151, and 292 kHz. The peak fluctuation amplitudes for the lateral velocities are 175 m/sec in both simulations. The presence of the unsteady fluctuations within the jet structure is consistent with the experimental literature on jet-in-crossflow interactions. In particular, the capture of the wake vortices as well as their location is entirely consistent with the experimental investigations, exemplified in the work of Frik and Rhosko.

Three Dimensional Simulation of Gas-Radiation Interactions in Gas Lasers 183

Fig. 6. Fourier analysis of temporal variation of the streamwise or Z direction velocity (*W*)

While the predicted presence of the unsteady fluctuations within the jet structure is entirely consistent with the experimental literature on jet-in-crossflow interactions, the question of validation relative to COIL flow conditions is still relevant. This issue was addressed by conducting a non-reacting flow experiment in which a He/I2 mixture was injected into a He/O2 primary flow with the nozzle geometry and flow conditions matching the hardware simulated here. A single large and two small orifices were used on the upper surface of the nozzle, thus differing from the hardware simulated here which used an array of orifices on both the upper and lower surfaces. Planar laser induced fluorescence (PLIF) was used to image the jet fluid, providing the ability to visualize the jet structure. Fig. 7 shows the PLIF image

jet from the large orifice and causing fluorescence within the I2 molecules, illuminating the jet structure; Fig. 8 shows a 2-D cut from the simulation at the large injector centerline as a point of comparison. The combination of the primary flow velocity, the spacing between vortices in the jet fluid, and the size of the vortices gives information to compare directly with the simulations. The size of the first 3 vortices in the PLIF image is 0.8, 1.9, and 2.3 mm which give an indication of the spatial growth rate of the vortices as they are shed. The size of the first 3 predicted vortices in the simulation are 0.8, 1.1, and 1.3 mm, a somewhat smaller growth rate than that in the PLIF image. Observing the difference in the jet trajectory between the experiment and the simulation, attributable to the presence of an opposing set of identical jets across the vertical centerline in the simulation, it is possible that the smaller growth rate in the simulation is to due to interference from the opposite jets. However, the frequency determined

m thick sheet of 565 nm laser light passing through the centerline of the

component in Fig. 5.

resulting from a 600

Fig. 4. 3-D streamtraces (in black) initiated upstream of the large injector orifice 2% of the channel height from the wall in conjunction with vortex cores (in red) extracted from the 10 species, 22 reaction 3-D GASP COIL simulation.

Fig. 5. Time variance of the lateral (Z) velocity component (*W*) at various monitoring points in the jet nearfield from the COIL simulation.

Fig. 4. 3-D streamtraces (in black) initiated upstream of the large injector orifice 2% of the channel height from the wall in conjunction with vortex cores (in red) extracted from the 10

Fig. 5. Time variance of the lateral (Z) velocity component (*W*) at various monitoring points

species, 22 reaction 3-D GASP COIL simulation.

in the jet nearfield from the COIL simulation.

Fig. 6. Fourier analysis of temporal variation of the streamwise or Z direction velocity (*W*) component in Fig. 5.

While the predicted presence of the unsteady fluctuations within the jet structure is entirely consistent with the experimental literature on jet-in-crossflow interactions, the question of validation relative to COIL flow conditions is still relevant. This issue was addressed by conducting a non-reacting flow experiment in which a He/I2 mixture was injected into a He/O2 primary flow with the nozzle geometry and flow conditions matching the hardware simulated here. A single large and two small orifices were used on the upper surface of the nozzle, thus differing from the hardware simulated here which used an array of orifices on both the upper and lower surfaces. Planar laser induced fluorescence (PLIF) was used to image the jet fluid, providing the ability to visualize the jet structure. Fig. 7 shows the PLIF image resulting from a 600 m thick sheet of 565 nm laser light passing through the centerline of the jet from the large orifice and causing fluorescence within the I2 molecules, illuminating the jet structure; Fig. 8 shows a 2-D cut from the simulation at the large injector centerline as a point of comparison. The combination of the primary flow velocity, the spacing between vortices in the jet fluid, and the size of the vortices gives information to compare directly with the simulations. The size of the first 3 vortices in the PLIF image is 0.8, 1.9, and 2.3 mm which give an indication of the spatial growth rate of the vortices as they are shed. The size of the first 3 predicted vortices in the simulation are 0.8, 1.1, and 1.3 mm, a somewhat smaller growth rate than that in the PLIF image. Observing the difference in the jet trajectory between the experiment and the simulation, attributable to the presence of an opposing set of identical jets across the vertical centerline in the simulation, it is possible that the smaller growth rate in the simulation is to due to interference from the opposite jets. However, the frequency determined

Three Dimensional Simulation of Gas-Radiation Interactions in Gas Lasers 185

where *ttot* is the total time interval of the integration, *ztot* is the total distance of integration in

simulations performed here, the exit plane of the nozzle as shown in Fig. 1 was exported every 1.0e-07 sec of physical time, and the Z integration was performed along lines of constant Y. Fig. 10 shows a composite gain lineshapes versus normalized frequency with the frequency axis normalized to the line center of the I 2P1/2→2P3/2 transition at 1.31527

Here *ttot* is 6.410-04 sec and *ztot* is the full dimension in the z direction. This plot shows that the lineshape including the influence of the gas velocity is broadened somewhat and the line center gain is lower, indicating lower amplification for the regions of the flowfield where the

*,y,z,t)* is determined using the equations listed above. For the

m.

the Z axis direction, and

*(*

broadening effects of the velocity field are significant.

Fig. 9. Time, space averaged gain lineshape from the COIL simulation.

equation, it is seen that it consists of a shape function

amplitude given by:

The transition lineshapes predicted with the above methodology can next be fit using the gain equations from above to characterize temperature extraction when strong flow fluctuations impact the lineshape. The spatially, temporally averaged lineshape is made using the assumption that all of the Doppler broadening be due to random, thermal motion, allowing the lineshape to be treated as synthetic representation of experimentally measured lineshapes. The absence and the presence of the bulk gas velocity influence in the lineshape provides an independent parameter for quantifying the effect of bulk gas velocity in terms of temperature. To fit the lineshape equation to the spatially, temporally averaged lineshapes reduced from the unsteady simulations, the MINPACK19 implementation of the Levenberg-Marquardt20,21 nonlinear fit algorithm is used22 to determine the best fit of to the data as a function of a given set of fit parameters. Observing the functional form of the lineshape

> *(*

*)* given by the Voigt function, and an

from the spacing between the vortices in the PLIF image and the primary flow velocity is 132 kHz, differing by less than 1% from the 131 and 133.9 kHz frequencies extracted from the simulation associated with the vortex shedding.

Fig. 7. PLIF experiment image of a jet-in-crossflow at flow conditions and with a configuration similar to that in the simulations performed here. This image was used to extract frequency and scale for validation of the models unsteady prediction characteristics. Note that the nozzle throat blocks the PLIF laser light, darkening the upper 1/3 of the image.

Fig. 8. Planar cut at the centerline of the large injector orifice from a single time step in the COIL simulation.

The end question to be answered here is 'what effect does the unsteady flow structure have upon the predicted lineshape?' This question is answered by accumulating a dataset containing the spatial and temporal variation of the quantities necessary to calculate the gain. This dataset is first spatially averaged in the Z axis direction and then temporally averaged, consistent with the behavior of a tunable diode small signal gain probe or lasing action over a period of time. Mathematically, this is expressed as:

$$\alpha\_{\text{av}}\left(\nu,\nu\right) = \frac{1}{t\_{\text{tot}}z\_{\text{tot}}} \iint \alpha\left(\nu,\nu,z,t\right)dzdt$$

from the spacing between the vortices in the PLIF image and the primary flow velocity is 132 kHz, differing by less than 1% from the 131 and 133.9 kHz frequencies extracted from the

Fig. 7. PLIF experiment image of a jet-in-crossflow at flow conditions and with a configuration similar to that in the simulations performed here. This image was used to extract frequency and scale for validation of the models unsteady prediction characteristics. Note that the nozzle

Fig. 8. Planar cut at the centerline of the large injector orifice from a single time step in the

The end question to be answered here is 'what effect does the unsteady flow structure have upon the predicted lineshape?' This question is answered by accumulating a dataset containing the spatial and temporal variation of the quantities necessary to calculate the gain. This dataset is first spatially averaged in the Z axis direction and then temporally averaged, consistent with the behavior of a tunable diode small signal gain probe or lasing

> <sup>1</sup> , ,,, *ave tot tot*

*t z*

*y y z t dzdt*

 

throat blocks the PLIF laser light, darkening the upper 1/3 of the image.

action over a period of time. Mathematically, this is expressed as:

 

simulation associated with the vortex shedding.

COIL simulation.

where *ttot* is the total time interval of the integration, *ztot* is the total distance of integration in the Z axis direction, and *(,y,z,t)* is determined using the equations listed above. For the simulations performed here, the exit plane of the nozzle as shown in Fig. 1 was exported every 1.0e-07 sec of physical time, and the Z integration was performed along lines of constant Y. Fig. 10 shows a composite gain lineshapes versus normalized frequency with the frequency axis normalized to the line center of the I 2P1/2→2P3/2 transition at 1.31527m. Here *ttot* is 6.410-04 sec and *ztot* is the full dimension in the z direction. This plot shows that the lineshape including the influence of the gas velocity is broadened somewhat and the line center gain is lower, indicating lower amplification for the regions of the flowfield where the broadening effects of the velocity field are significant.

Fig. 9. Time, space averaged gain lineshape from the COIL simulation.

The transition lineshapes predicted with the above methodology can next be fit using the gain equations from above to characterize temperature extraction when strong flow fluctuations impact the lineshape. The spatially, temporally averaged lineshape is made using the assumption that all of the Doppler broadening be due to random, thermal motion, allowing the lineshape to be treated as synthetic representation of experimentally measured lineshapes. The absence and the presence of the bulk gas velocity influence in the lineshape provides an independent parameter for quantifying the effect of bulk gas velocity in terms of temperature. To fit the lineshape equation to the spatially, temporally averaged lineshapes reduced from the unsteady simulations, the MINPACK19 implementation of the Levenberg-Marquardt20,21 nonlinear fit algorithm is used22 to determine the best fit of to the data as a function of a given set of fit parameters. Observing the functional form of the lineshape equation, it is seen that it consists of a shape function *()* given by the Voigt function, and an amplitude given by:

Three Dimensional Simulation of Gas-Radiation Interactions in Gas Lasers 187

lineshape introduced errors because the integral bounds are limited to the frequency domain of ±1.75 GHz about the line center used to construct the lineshape numerically and are and not infinite as required by the Voigt function. The effect of this error on the fit is

> *(*

K error in the temperature determination, and was subsequently determined to be avoidable

, , (, , )

temperature can be found using the gain equations from above in conjunction with an independent determination of the pressure and the mole fractions. Here, the simulation data were spatially and temporally averaged using the same procedure as for the gain determine average pressures and species mole fractions. When performing this analysis with experimentally measured lineshapes, wall pressures and mean mole fractions based upon the measured molar flow rates may be the only values available and may lead to errors in

The degree of agreement between the Doppler and Lorentzian determined temperatures and the average thermodynamic temperature from the simulation provides a means to assess the lineshape fit and the interpretation of the fit. To perform this comparison, the thermodynamic temperatures were averaged from the simulation data using the same integration procedure as used for the gain. As a test of the procedure, a fit was performed to a lineshape that was calculated without the functional dependence upon the bulk velocity. The comparison between the Doppler, Lorentzian, and thermodynamic temperatures in this test is shown in Fig. 10 is very close, to within 1%. Next, the temperatures are determined for the lineshape in which the bulk velocity was included in the functional dependence of the lineshape equation. As shown in Fig. 11, the Doppler temperature diverges from the Lorentzian and thermodynamic temperatures in the flow regions closer to the nozzle centerline, with the maximum difference being 14 K. This trend is consistent with the observation that the unsteady eddies increase in size toward the center of the nozzle as illustrated in Fig. 2. Thus the Doppler temperature is registering the clear influence of the bulk gas velocity, whereas the Lorentzian temperature remains consistent with the thermodynamic temperature. Given the results illustrated here, a suggested method of fitting experiment lineshapes Where poor signal-to-noise can interfere with the nonlinear fitting would be to fix the Lorentzian width based upon an independent temperature measurement such as a rotational temperature, and then determine the Doppler width. Then the degree to which the Doppler temperature differs from the rotational or other temperature measurement would be an indicator of the magnitude of flow rotation. It should be noted that fitting procedures in which the Doppler and Lorentzian component temperatures are forced to be equal, as is commonly done, could very well be masking information regarding the magnitude of eddy strength and the degree of unsteadiness in the experiments. Another test would be to re-visit existing lineshape temperature determinations for experiments in which flow unsteadiness could potentially be present, and perform the nonlinear fit without fixing either the Doppler and Lorentzian component. The difference in temperatures could then be correlated against fluid dynamic parameters

 

*<sup>D</sup> <sup>L</sup> D L ν*

*)* via is 0.99, not 1.0. This translated into a 5

*<sup>L</sup>* determined, separate values for the

*)* were treated as a separate fitting parameter. In practical

 

*<sup>L</sup>* are the fitting parameters used in the nonlinear fit. With the values for

such that evaluation of the infinite integral of

*(*

 

*D* and

such as Reynolds number that are relevant to eddy generation to identify trends.

if the infinite integral over

the extracted temperatures.

where , *D*, and

terms, the fitting function becomes:

the Doppler and Lorentzian widths

$$\left(\frac{7}{12}\left(\frac{A\lambda\_0}{8\pi}\right)\left(N\_{I\_{I\_\circ}^{\circ\_2}R\_{\mathbb{Z}\_2}}-\frac{1}{2}N\_{I\_{I\_\circ}^{\circ\_2}R\_{\mathbb{Z}\_2}}\right)\right)$$

Given that the Voigt represents the probability versus frequency, the infinite integral is 1.

$$\int\_{-\infty}^{\infty} \rho(\nu)d\nu = \mathbb{I}$$

Integrating allows for determination of and separation of the Voigt from the lineshape:

$$\begin{split} \int\_{-\infty}^{+\infty} a(\nu) \, d\nu &= \int\_{-\infty}^{+\infty} \frac{7}{12} \left( \frac{A \dot{\lambda}\_0^{-2}}{8\pi} \right) \phi(\nu) \left( N\_{I\_{\parallel}^{-2} P\_{\parallel \perp 2}} - \frac{1}{2} N\_{I\_{\parallel}^{-2} P\_{\parallel \perp 2}} \right) d\nu \\ &= \frac{7}{12} \left( \frac{A \dot{\lambda}\_0^{-2}}{8\pi} \right) \left( N\_{I\_{\parallel}^{-2} P\_{\parallel \perp 2}} - \frac{1}{2} N\_{I\_{\parallel}^{-2} P\_{\parallel \perp 2}} \right) \prod\_{-\infty}^{+\infty} \phi(\nu) d\nu \end{split}$$

and

$$\alpha\left(\nu\right) = \phi(\nu) \int\_{-\infty}^{+\infty} a\left(\nu\right) d\nu$$

Fig. 10. Comparison of temperature values from the 13 species, 52 reaction simulation; Doppler and Lorentzian values were determined based upon the lineshape in Fig. 9 without the velocity dependence and the thermodynamic values averaged from the simulation data.

Thus the individual values of the gain *()* at each can be divided by the infinite integral of *()* to recover the Voigt. However, in practice it was found that fitting this normalized

12 8 */ /* 2 *I( P ) I( P ) <sup>A</sup><sup>λ</sup> N N*

Given that the Voigt represents the probability versus frequency, the infinite integral is 1.

 () 1

*<sup>d</sup>*

7 1 12 8 2

*α ν <sup>d</sup> <sup>ν</sup> N Nd <sup>π</sup>*

 *ν αν <sup>d</sup>*

<sup>2</sup> 7 1 <sup>0</sup>

Integrating allows for determination of and separation of the Voigt from the lineshape:

2 0

7 1 12 8 2

> ( )

Fig. 10. Comparison of temperature values from the 13 species, 52 reaction simulation; Doppler and Lorentzian values were determined based upon the lineshape in Fig. 9 without the velocity dependence and the thermodynamic values averaged from the simulation data.

*)* at each

*)* to recover the Voigt. However, in practice it was found that fitting this normalized

*(*

Thus the individual values of the gain

2 0

*Aλ*

*π*

and

of *(* 2 2 1 2 3 2

 

can be divided by the infinite integral

2 2 1 2 3 2

*/ /*

*I( P ) I( P )*

2 2 1 2 3 2

*/ /*

 

*I( P ) I( P )*

*<sup>A</sup><sup>λ</sup> N N <sup>ν</sup> <sup>d</sup> <sup>π</sup>*

lineshape introduced errors because the integral bounds are limited to the frequency domain of ±1.75 GHz about the line center used to construct the lineshape numerically and are and not infinite as required by the Voigt function. The effect of this error on the fit is such that evaluation of the infinite integral of *()* via is 0.99, not 1.0. This translated into a 5 K error in the temperature determination, and was subsequently determined to be avoidable if the infinite integral over *()* were treated as a separate fitting parameter. In practical terms, the fitting function becomes:

$$\alpha\left(\nu, \Delta\nu\_D, \Delta\nu\_L\right) = \beta\rho(\nu, \Delta\nu\_D, \Delta\nu\_L)$$

where , *D*, and *<sup>L</sup>* are the fitting parameters used in the nonlinear fit. With the values for the Doppler and Lorentzian widths *D* and *<sup>L</sup>* determined, separate values for the temperature can be found using the gain equations from above in conjunction with an independent determination of the pressure and the mole fractions. Here, the simulation data were spatially and temporally averaged using the same procedure as for the gain determine average pressures and species mole fractions. When performing this analysis with experimentally measured lineshapes, wall pressures and mean mole fractions based upon the measured molar flow rates may be the only values available and may lead to errors in the extracted temperatures.

The degree of agreement between the Doppler and Lorentzian determined temperatures and the average thermodynamic temperature from the simulation provides a means to assess the lineshape fit and the interpretation of the fit. To perform this comparison, the thermodynamic temperatures were averaged from the simulation data using the same integration procedure as used for the gain. As a test of the procedure, a fit was performed to a lineshape that was calculated without the functional dependence upon the bulk velocity. The comparison between the Doppler, Lorentzian, and thermodynamic temperatures in this test is shown in Fig. 10 is very close, to within 1%. Next, the temperatures are determined for the lineshape in which the bulk velocity was included in the functional dependence of the lineshape equation. As shown in Fig. 11, the Doppler temperature diverges from the Lorentzian and thermodynamic temperatures in the flow regions closer to the nozzle centerline, with the maximum difference being 14 K. This trend is consistent with the observation that the unsteady eddies increase in size toward the center of the nozzle as illustrated in Fig. 2. Thus the Doppler temperature is registering the clear influence of the bulk gas velocity, whereas the Lorentzian temperature remains consistent with the thermodynamic temperature. Given the results illustrated here, a suggested method of fitting experiment lineshapes Where poor signal-to-noise can interfere with the nonlinear fitting would be to fix the Lorentzian width based upon an independent temperature measurement such as a rotational temperature, and then determine the Doppler width. Then the degree to which the Doppler temperature differs from the rotational or other temperature measurement would be an indicator of the magnitude of flow rotation. It should be noted that fitting procedures in which the Doppler and Lorentzian component temperatures are forced to be equal, as is commonly done, could very well be masking information regarding the magnitude of eddy strength and the degree of unsteadiness in the experiments. Another test would be to re-visit existing lineshape temperature determinations for experiments in which flow unsteadiness could potentially be present, and perform the nonlinear fit without fixing either the Doppler and Lorentzian component. The difference in temperatures could then be correlated against fluid dynamic parameters such as Reynolds number that are relevant to eddy generation to identify trends.

Three Dimensional Simulation of Gas-Radiation Interactions in Gas Lasers 189

The first point of examination of the effect of increased grid resolution is the vorticity. The isosurface of constant total vorticity at a value of 2.0·105 s-1 is generated for the new grid and plotted in Fig. 13. The first and most obvious feature distinguishing the two plots is that the vorticity loss within the downstream end of the nozzle expansion found in the simulation using the original grid does not appear in the simulation using the higher cell density grid, confirming that the grid is impacting the transport of the flow fluctuations in this region.

Fig. 12. An isosurface of constant total vorticity with a value of 2.0·105 (1/s). This simulation

Fig. 13. An isosurface of constant total vorticity with a value of 2.0·105 (1/s). This simulation used 512 points in the streamwise direction in the grid within the nozzle expansion. Note

While vorticity is a measure of flow rotation and is associated with the velocity components in the lateral direction that influence the lineshape, it is also correlated with mixing and the grid resolution effect shown here would be expected to extend beyond the broadening to the

the continued presence of the vorticity through the nozzle exit.

used 128 points in the streamwise direction in the grid within the nozzle expansion.

The temperature determination from the averaged lineshapes illustrates the sensitivity of the temperature determination process to the magnitude of the flow velocity components. Observing the decrease in grid density within the nozzle expansion in Fig. 1 and the dissipation of the vortices within Fig. 2, it is reasonable to question whether the increase in grid cell size with streamwise distance in the expansion is filtering the smaller eddies as well as decreasing the strength of the larger eddies through numerical dissipation. Fig. 12 illustrates this point through the plotting of a constant vorticity isosurface with a value of 2.0·105 s-1. The figure shows that this particular vorticity magnitude, while persistent within the freestream before the nozzle throat, disappears rapidly as the flow passes through the throat. This decrease is consistent not just with the decrease in grid density but also with the physical effect of the flow expansion and the commensurate stretching of the vortices. To separate the effect of flow dilation upon vorticity from that of grid dissipation, it is necessary to decrease grid dissipation by increasing the grid resolution.

Fig. 11. Comparison of temperature values from the simulation; Doppler and Lorentzian values were determined based upon the lineshape in Fig. 9 with the velocity dependence and the thermodynamic values averaged from the simulation data.

The effect of grid dissipation upon the transport of vorticity in the nozzle expansion and the subsequent effect upon the lineshape is examined by quadrupling the number of grid cells in the streamwise direction within the nozzle expansion while keeping the number of cells in the transverse directions constant. This increases the number of cells in the streamwise direction from 128 to 512. This new grid is used to continue the simulation and accumulate sufficient temporal information regarding the influence of the unsteadiness upon the lineshape. Data necessary for the lineshape analysis were collected within the same plane transverse to the flow as in the first simulation.

The temperature determination from the averaged lineshapes illustrates the sensitivity of the temperature determination process to the magnitude of the flow velocity components. Observing the decrease in grid density within the nozzle expansion in Fig. 1 and the dissipation of the vortices within Fig. 2, it is reasonable to question whether the increase in grid cell size with streamwise distance in the expansion is filtering the smaller eddies as well as decreasing the strength of the larger eddies through numerical dissipation. Fig. 12 illustrates this point through the plotting of a constant vorticity isosurface with a value of 2.0·105 s-1. The figure shows that this particular vorticity magnitude, while persistent within the freestream before the nozzle throat, disappears rapidly as the flow passes through the throat. This decrease is consistent not just with the decrease in grid density but also with the physical effect of the flow expansion and the commensurate stretching of the vortices. To separate the effect of flow dilation upon vorticity from that of grid dissipation, it is

necessary to decrease grid dissipation by increasing the grid resolution.

Fig. 11. Comparison of temperature values from the simulation; Doppler and Lorentzian values were determined based upon the lineshape in Fig. 9 with the velocity dependence

The effect of grid dissipation upon the transport of vorticity in the nozzle expansion and the subsequent effect upon the lineshape is examined by quadrupling the number of grid cells in the streamwise direction within the nozzle expansion while keeping the number of cells in the transverse directions constant. This increases the number of cells in the streamwise direction from 128 to 512. This new grid is used to continue the simulation and accumulate sufficient temporal information regarding the influence of the unsteadiness upon the lineshape. Data necessary for the lineshape analysis were collected within the same plane

and the thermodynamic values averaged from the simulation data.

transverse to the flow as in the first simulation.

The first point of examination of the effect of increased grid resolution is the vorticity. The isosurface of constant total vorticity at a value of 2.0·105 s-1 is generated for the new grid and plotted in Fig. 13. The first and most obvious feature distinguishing the two plots is that the vorticity loss within the downstream end of the nozzle expansion found in the simulation using the original grid does not appear in the simulation using the higher cell density grid, confirming that the grid is impacting the transport of the flow fluctuations in this region.

Fig. 12. An isosurface of constant total vorticity with a value of 2.0·105 (1/s). This simulation used 128 points in the streamwise direction in the grid within the nozzle expansion.

Fig. 13. An isosurface of constant total vorticity with a value of 2.0·105 (1/s). This simulation used 512 points in the streamwise direction in the grid within the nozzle expansion. Note the continued presence of the vorticity through the nozzle exit.

While vorticity is a measure of flow rotation and is associated with the velocity components in the lateral direction that influence the lineshape, it is also correlated with mixing and the grid resolution effect shown here would be expected to extend beyond the broadening to the

Three Dimensional Simulation of Gas-Radiation Interactions in Gas Lasers 191

is fit with a Gaussian distribution function, shown with the solid blue line in Fig. 15, the Gaussian distribution shows a somewhat lower probability in the region from -40 to -100 m/sec and 40 to 100 m/sec. The reason for this difference can be explained in terms of the turbulent flow assumption underlying the use of the Gaussian velocity distribution. The assumption of turbulence implies a character to the flow fluctuations that is consistent with a Gaussian distribution, i.e. random and having an exponential range of energies and scales. Fig. 15 indicates that the predicted velocity probability distribution does not have this type of character, implying that the random, stochastic character of isotropic turbulence that allows an exponential range of energies and scales to develop is not present in this

Fig. 15. Velocity probability as determined from statistics collected in the simulation

The question of turbulence can be examined in a preliminary sense by examining the temporal variation of the lateral velocity component at a point within the 2D plane. The point chosen is within a region in this plane marked by high fluctuation magnitudes, is associated with the largest differences between the Doppler and Lorentzian temperatures, and is more likely to be turbulent. Fig. 16 shows the variation of the *W* velocity component with time at this point. As can be seen, the velocity variation is regular and sinusoidal in character, with some random fluctuation in the amplitude within the first 7x10-5 sec and near 2.4x10-4 sec. Comparing with Fig. 5 which plots the *W* velocity versus time in the region surrounding the injectors where the flow unsteadiness originates, it is seen that the sinusoidal character of the fluctuations is preserved from the injection region, through the nozzle throat, and downstream to the end of the isentropic expansion region where the data for Fig. 16 are extracted. If the flow were to transition to fully isotropic turbulence, the sinusoidal character of the fluctuations would be lost between Figs. 5 and 16 as the larger

simulation. Simply said, the flow is not turbulent.

compared with a Gaussian fit to this data.

general prediction of the COIL flowfield. An examination of the Doppler and Lorentzian temperature components based upon fitting of the lineshapes for the original and higher resolution grid simulations provides the final point of comparison. Fig. 14 shows the Doppler, Lorentzian, and thermodynamic temperatures for the new grid simulation compared with the Doppler temperature from the original grid simulation. Comparing the Doppler temperatures, the values from the higher density grid simulation are somewhat higher over a greater section of the nozzle than the original grid, consistent with the analysis performed above. Additionally, the maximum difference between the Doppler and Lorentzian temperatures for the new grid simulation is 25 K, an increase from 14 K for the original grid simulation. This confirms the hypothesis that the grid resolution impacts the lineshape analysis and strengthens the argument that flow unsteadiness is a factor in temperature extraction from lineshapes in these flows.

Fig. 14. Comparison of temperature values from both the original and high density grid simulations. Note that the Dopper temperature is somewhat higher with the high density grid simulation.

The result of a 25 K effect of the unsteady flow structures on the temperature extraction from the Doppler and Lorentzian components naturally leads to questions about the velocity distribution predicted here and the Gaussian distribution assumed by Stuve and Elvey and later by Nikolaev. To analyze this question, the time dependent *W* velocity component data from the 2D planes used above to generate lineshapes was separated into 1 m/sec bins ranging from -200 to 200 m/sec, and temporal and spatial frequency of occurrence was counted. The probability of occurrence for a particular 1 m/sec bin was calculated as the frequency of occurrence divided by the total frequency of occurrence for all velocity bins within the freestream portion of the 2D plane across the entire time of integration. The plot of this data is shown in Fig. 15. The distribution shows a symmetric and finite probability of occurrence of the *W* velocity from -100 to 100 m/sec. When this data

general prediction of the COIL flowfield. An examination of the Doppler and Lorentzian temperature components based upon fitting of the lineshapes for the original and higher resolution grid simulations provides the final point of comparison. Fig. 14 shows the Doppler, Lorentzian, and thermodynamic temperatures for the new grid simulation compared with the Doppler temperature from the original grid simulation. Comparing the Doppler temperatures, the values from the higher density grid simulation are somewhat higher over a greater section of the nozzle than the original grid, consistent with the analysis performed above. Additionally, the maximum difference between the Doppler and Lorentzian temperatures for the new grid simulation is 25 K, an increase from 14 K for the original grid simulation. This confirms the hypothesis that the grid resolution impacts the lineshape analysis and strengthens the argument that flow unsteadiness is a factor in

Fig. 14. Comparison of temperature values from both the original and high density grid simulations. Note that the Dopper temperature is somewhat higher with the high density

The result of a 25 K effect of the unsteady flow structures on the temperature extraction from the Doppler and Lorentzian components naturally leads to questions about the velocity distribution predicted here and the Gaussian distribution assumed by Stuve and Elvey and later by Nikolaev. To analyze this question, the time dependent *W* velocity component data from the 2D planes used above to generate lineshapes was separated into 1 m/sec bins ranging from -200 to 200 m/sec, and temporal and spatial frequency of occurrence was counted. The probability of occurrence for a particular 1 m/sec bin was calculated as the frequency of occurrence divided by the total frequency of occurrence for all velocity bins within the freestream portion of the 2D plane across the entire time of integration. The plot of this data is shown in Fig. 15. The distribution shows a symmetric and finite probability of occurrence of the *W* velocity from -100 to 100 m/sec. When this data

temperature extraction from lineshapes in these flows.

grid simulation.

is fit with a Gaussian distribution function, shown with the solid blue line in Fig. 15, the Gaussian distribution shows a somewhat lower probability in the region from -40 to -100 m/sec and 40 to 100 m/sec. The reason for this difference can be explained in terms of the turbulent flow assumption underlying the use of the Gaussian velocity distribution. The assumption of turbulence implies a character to the flow fluctuations that is consistent with a Gaussian distribution, i.e. random and having an exponential range of energies and scales. Fig. 15 indicates that the predicted velocity probability distribution does not have this type of character, implying that the random, stochastic character of isotropic turbulence that allows an exponential range of energies and scales to develop is not present in this simulation. Simply said, the flow is not turbulent.

Fig. 15. Velocity probability as determined from statistics collected in the simulation compared with a Gaussian fit to this data.

The question of turbulence can be examined in a preliminary sense by examining the temporal variation of the lateral velocity component at a point within the 2D plane. The point chosen is within a region in this plane marked by high fluctuation magnitudes, is associated with the largest differences between the Doppler and Lorentzian temperatures, and is more likely to be turbulent. Fig. 16 shows the variation of the *W* velocity component with time at this point. As can be seen, the velocity variation is regular and sinusoidal in character, with some random fluctuation in the amplitude within the first 7x10-5 sec and near 2.4x10-4 sec. Comparing with Fig. 5 which plots the *W* velocity versus time in the region surrounding the injectors where the flow unsteadiness originates, it is seen that the sinusoidal character of the fluctuations is preserved from the injection region, through the nozzle throat, and downstream to the end of the isentropic expansion region where the data for Fig. 16 are extracted. If the flow were to transition to fully isotropic turbulence, the sinusoidal character of the fluctuations would be lost between Figs. 5 and 16 as the larger

Three Dimensional Simulation of Gas-Radiation Interactions in Gas Lasers 193

variation based upon location in the flow. Based upon the observed decrease in vorticity within the nozzle expansion in conjunction with decreased grid density, an additional simulation was performed that increased the grid in the streamwise direction in this region by a factor of 4. Comparing the results from the simulations from different grids, the higher density grid dramatically decreased dissipation of vorticity, captured higher magnitude velocity amplitudes, and showed a 11 K increase in the difference between the Doppler and Lorentzian temperatures over the original grid simulation, illustrating the impact of numerical considerations upon the predictions of the model. Based upon the observed trends in temperature prediction, it is suggested that the difference between the temperatures extracted from the fit values for the Doppler and Lorentzian widths could be correlated against changes in fluid dynamic parameters such as Reynolds number to

The velocity probability distribution for the flow at the exit of the nozzle expansion was examined, in the location where the lineshape examination was performed. Analysis of this velocity distribution showed it to be non-Gaussian, with higher probabilities for higher magnitude velocity fluctuations. Further analysis illustrates that the fluctuations are sinusoidal in nature as is characteristic of laminar unsteady flows, as opposed to the random, stochastic fluctuations associated with turbulent flows. This indicates that a Gaussian probability distribution function, as is commonly used in the astronomical community for Voigt lineshape analysis of spectral lineshapes of atoms in stellar atmospheres that are assumed to have turbulent bulk flow velocities, is not strictly

The author would like to thank Drs. Kevin Hewett and John McCord for fruitful discussions regarding the methods used in experimental lineshape measurement and temperature extraction for COIL devices. Mr. Jeremy Stanford generated the computational grids used in the simulations discussed here. Computational resources were provided by the DoD High Performance Computing Modernization Challenge Allocation at the Aeronautical Systems Center, Army Research Laboratory, Engineering Research and Development Center, and

[1] Struve, O. and Morgan, W. W., "On the Intensities of Stellar Absorption Lines,"

[2] Elvey, C. T., "The Intensities of Some Multiplets of fe II and TI II in Stellar Spectra,"

[3] Struve, O. and Elvey, C. T., "The Intensities of Stellar Absorption Lines," *Astrophysical* 

[4] Allen, M.G., Carleton, K.L., Davis, S.J., Kessler, W. J., and McManus, K. R., "Diode laser-

*Proceedings of the National Academy of Sciences of the United States of America*, 18, No.

based measurements of water vapor and ground state oxygen in chemical oxygen iodine lasers," AIAA-1994-2433, *25th AIAA Plasmadynamics and Lasers Conference*,

determine trends.

appropriate for these flows.

**5. Acknowledgements** 

**6. References** 

Naval Oceanographic Major Shared Resource Centers.

*Astrophysical Journal*, 79, pp. 263-269, 1934.

9, pp. 590-594, 1932.

*Journal*, 79, pp. 409-440, 1934.

Colorado Springs, CO, 1994.

scale vortices interact to generate smaller vortices, leading to the cascade of temporal and spatial scales that characterizes turbulence.

Returning to the issue of lineshape fitting and the use of a Gaussian function to describe the bulk velocity probability distribution, the results of this work indicate that for the low Reynolds number flows simulated here that are typical of COIL operation, a Gaussian velocity distribution function does not describe the bulk flow fluctuations because the flow is not turbulent. However, this statement has broader implications as statistical turbulence descriptions such as the *k*, and *k*, models are now frequently employed in COIL simulation. Since the COIL chemistry, particularly I2 dissociation, is sensitive to the local concentrations of reactants and the concentrations are directly determined by mixing and chemistry, the predictions given by COIL models that use turbulence models to predict the mixing should be viewed cautiously. This statement holds particularly for efforts to develop and prove new COIL chemistry mechanisms, where accurate prediction of the reactant mixing characteristics provides a basis for accurate examination of the predictions from the chemistry model.

Fig. 16. Example of *W* velocity variation with time within the nozzle exit plane for the higher density grid over the period of the simulation.

#### **4. Summary and conclusions**

In conclusion, the results provided here predict that the bulk gas flow velocity components associated with vortex rotation will broaden the lineshape for quantum transitions of species within the gas flow. The studies performed here show that if the bulk gas velocity influence is ignored when fitting the lineshape and the vortices are present, the effect will be a divergence between the reduced Doppler and Lorentzian values of the temperature. The results presented here indicate that this difference will be as high as 25 K, with some variation based upon location in the flow. Based upon the observed decrease in vorticity within the nozzle expansion in conjunction with decreased grid density, an additional simulation was performed that increased the grid in the streamwise direction in this region by a factor of 4. Comparing the results from the simulations from different grids, the higher density grid dramatically decreased dissipation of vorticity, captured higher magnitude velocity amplitudes, and showed a 11 K increase in the difference between the Doppler and Lorentzian temperatures over the original grid simulation, illustrating the impact of numerical considerations upon the predictions of the model. Based upon the observed trends in temperature prediction, it is suggested that the difference between the temperatures extracted from the fit values for the Doppler and Lorentzian widths could be correlated against changes in fluid dynamic parameters such as Reynolds number to determine trends.

The velocity probability distribution for the flow at the exit of the nozzle expansion was examined, in the location where the lineshape examination was performed. Analysis of this velocity distribution showed it to be non-Gaussian, with higher probabilities for higher magnitude velocity fluctuations. Further analysis illustrates that the fluctuations are sinusoidal in nature as is characteristic of laminar unsteady flows, as opposed to the random, stochastic fluctuations associated with turbulent flows. This indicates that a Gaussian probability distribution function, as is commonly used in the astronomical community for Voigt lineshape analysis of spectral lineshapes of atoms in stellar atmospheres that are assumed to have turbulent bulk flow velocities, is not strictly appropriate for these flows.

#### **5. Acknowledgements**

192 Fluid Dynamics, Computational Modeling and Applications

scale vortices interact to generate smaller vortices, leading to the cascade of temporal and

Returning to the issue of lineshape fitting and the use of a Gaussian function to describe the bulk velocity probability distribution, the results of this work indicate that for the low Reynolds number flows simulated here that are typical of COIL operation, a Gaussian velocity distribution function does not describe the bulk flow fluctuations because the flow is not turbulent. However, this statement has broader implications as statistical

COIL simulation. Since the COIL chemistry, particularly I2 dissociation, is sensitive to the local concentrations of reactants and the concentrations are directly determined by mixing and chemistry, the predictions given by COIL models that use turbulence models to predict the mixing should be viewed cautiously. This statement holds particularly for efforts to develop and prove new COIL chemistry mechanisms, where accurate prediction of the reactant mixing characteristics provides a basis for accurate examination of the

Fig. 16. Example of *W* velocity variation with time within the nozzle exit plane for the higher

In conclusion, the results provided here predict that the bulk gas flow velocity components associated with vortex rotation will broaden the lineshape for quantum transitions of species within the gas flow. The studies performed here show that if the bulk gas velocity influence is ignored when fitting the lineshape and the vortices are present, the effect will be a divergence between the reduced Doppler and Lorentzian values of the temperature. The results presented here indicate that this difference will be as high as 25 K, with some

models are now frequently employed in

and *k*,

spatial scales that characterizes turbulence.

turbulence descriptions such as the *k*,

predictions from the chemistry model.

density grid over the period of the simulation.

**4. Summary and conclusions** 

The author would like to thank Drs. Kevin Hewett and John McCord for fruitful discussions regarding the methods used in experimental lineshape measurement and temperature extraction for COIL devices. Mr. Jeremy Stanford generated the computational grids used in the simulations discussed here. Computational resources were provided by the DoD High Performance Computing Modernization Challenge Allocation at the Aeronautical Systems Center, Army Research Laboratory, Engineering Research and Development Center, and Naval Oceanographic Major Shared Resource Centers.

#### **6. References**


**9** 

*Romania* 

**Fluid-Structure Interaction** 

Stoia-Djeska Marius and Safta Carmen-Anca *POLITEHNICA University of Bucharest, Bucharest* 

The fluid-structure interaction (FSI) phenomenon is the result of the interactions of multiple continuum fields. The fluid (gas or/and liquid) forces act on a neighboring elastic solid, which is deformed and thus influences the flow of the adhering fluid. Due to the deformation of the solid, both the fluid velocity and the fluid domain change. Usually, this interaction takes place in the presence of other mechanical fields and as important examples we mention here the external body forces and acoustic fields. The mathematical modeling of the FSI phenomenon requires mathematical models for each field implied in and also the coupling mechanisms between the fields. For instance, one can imagine a mathematical model of a FSI phenomenon where the fluid flow field is modeled with the unsteady Navier-Stokes equations, the solid elastic field with the Navier equations while for the acoustic field the inhomogeneous wave equation based on Lighthill's analogy can be used and the external body forces are due to gravity. The boundary conditions for the fluid, solid and acoustic models represent here the coupling mechanism between the fields. However, even everything looks simple and clear the accurate and efficient solution of FSI problems is still a highly complicated task and an open area of research. Furthermore, using complicated models as those enumerated above the only way to obtain a solution to a specific FSI problem is the numerical one. The Computational Fluid Dynamics (CFD), Computational Structural Dynamics (CSD), Finite Element Method (FEM) provide us with specific **mathematical models** and **numerical techniques** that can be coupled to build a numerical solver for a FSI problem. However, even for simple and without engineering relevance problems the computational effort required to solve a specific problem becomes huge if the **mathematical modeling of the physical phenomenon** is not carefully done. As in many other multiple field interaction problems the **level of accuracy** of the mathematical model is crucial for both the efficiency and quality of the results. Further, we notice also that the

development of in-house codes for such purposes is a difficult task.

Today we know from the experience earned from a century of practice that different levels of approximation can be used successfully to build up a simplified mathematical model for a FSI problem. This can be done by simplifying the flow equations and by choosing engineering models for the elastic structure. With appropriate coupling mechanisms one can thus obtain mathematical models which are well-matched for practical purposes, starting from the design and up to advanced (flow and/or structural dynamics) control purposes. The use of different simplifying hypothesis (about the dimensionality of the flow and structural models, for instance) requires from the engineer a deep understanding of the multifield interaction phenomenon so that by neglecting different components of the

**1. Introduction** 

