1. INTRODUCTION

All practical high-power laser systems require a final optical element for collimating, focusing, steering, or combining. Except for the case of beam combining, the final optical element must withstand the full beam fluence and intensity without damage or causing significant aberrations. For high-peak-power laser systems, either laser-induced breakdown will immediately damage conventional (solid-state) optics or Kerr self-focusing will cause beam collapse within the optical element leading to laser-induced breakdown [1]. High-average-power lasers can cause thermal expansion of optical elements which introduce power- and time-dependent optical aberrations [2,3].

One solution to these issues in high-peak- or high-average-power lasers is to use larger diameter optics. This drops fluence and intensity but at the expense of system size, cost, and maneuverability. For high-average-power lasers, thermal management in conventional optics becomes more difficult as lens size increases. Temporal spreading of energy, such as in chirped-pulse amplification, can also drop peak powers but only for pulses with sufficient bandwidth and this is not effective if a final optical element is used for focusing and steering of the compressed laser pulse.

To overcome these issues, it has been proposed to use a gas vortex to create a lens [4]. A gas-based lens has a number of advantages over conventional optics. First, gas densities are about ${10}^{-3}$ that of conventional optics. This increases the Kerr self-focusing critical power from typical values of megawatts to gigawatts reducing the risk of ionization. If ionization does occur, the lower density results in a proportionate decrease in the plasma generation. Independent of the overall density, the ionization potentials of gases are often larger than that of solids which helps suppress the onset of laser-induced breakdown. For example, helium and nitrogen have ionization potentials of 24.6 eV and 14.7 eV, while ${\mathrm{SiO}}_2$ and BK7 are 9 eV and 4.5 eV, respectively. The large gas flow necessary to form density gradients makes thermal management and the risk of thermally induced optical aberrations less of a concern because the gas is being replenished on 100 μs time scales. Finally, if laser-induced breakdown does occur, the result is not the permanent damage of an optical element but brief plasma generation and then replacement with fresh gas. The gas vortex lens allows operation of high-power lasers above the damage thresholds of conventional optics and, additionally, its self-healing design allows operation near gas breakdown thresholds without risk to the optical element. Specific applications include using the gas vortex lens in a beam expander for a high-energy laser weapon system, a beam expander to non-destructively end laser-plasma filaments [5], or, after ionization, as a tunable, focusing, plasma lens for petawatt-class lasers [6].

A refractive index gradient suitable for lensing can be generated by several mechanisms: thermal gradients, gas-composition gradients, conductive-advective cooling, and centrifugal potentials. The use of thermal [7–9] and gas-composition [10,11] gradients were some of the first mechanisms reported for generating gas-based lenses. Thermal gas lenses consist of gas flow through a heated tube. The lower density (hotter gas) near the tube walls creates a positive lens. Later work extended this concept by rotating the tube about its axis so that convective flows would not distort the lens [12–15]. Gas-composition gradients are created by flowing one gas through a porous tube while a second gas (with a different refractive index) diffuses into the tube. Gas-composition gradients can also be created with spatially dependent ionization of a gas, as seen in inductive [16,17] and dielectric-barrier discharge [18,19] plasma lenses. A colliding-jet flow [20] generates density gradients when gas jets collide along the optical axis and form a density minimum on-axis via conductive-advection cooling. Offsetting the gas jets was shown to stabilize the density profile with a rotational flow [20].

The proposed device design of this work is similar to our previous work [4]. The design is seen in Fig. 1 where two, offset, colliding jets (inlets shown in yellow) are enclosed by a cylindrical chamber with two outlets (shown in blue) along the optical axis ($z$ axis). This transversely constrains the flow of the colliding jets and the resulting rotation stabilizes it against turbulence [20]. Unlike the colliding jets, the gas density gradient is created by a centrifugal potential from the rotating flow. Recent experiments have shown the current design to be able to produce a defocusing gas lens [4] but a detailed understanding of the flow structure remained elusive.

 

Fig. 1. Interior walls of the gas lens are shown in gray. The inlets and outlets are shown in yellow and blue, respectively. The optical axis corresponds to the $z$ axis. Characteristic streamlines of the flow are shown in blue for fluid parcels that originate in the $z>0$ half of the top inlet. The flow has mirror symmetry across the $z=0$ plane and two-fold rotational symmetry around the $z$ axis.

Download Full Size | PDF

Previous works have explored the concept of using a gas vortex for an optical lens [21–25]. The earlier publications report on the dependence of the focal length on the axial pressure and a theoretical description of the focal length based on an adiabatic flow [22]. A subsequent series of experimentally based publications detail the focal length dependence on the vortex tube length and the formation of an annular laser beam [23–25]. But there has yet to be full three-dimensional fluid simulations detailing the gas vortex lens and flow structure resulting optical properties.

In this paper, we will show steady-state three-dimensional fluid dynamics simulations of the gas vortex lens and time-dependent fluid simulations of the transient mass flow rate. The Gladstone–Dale relation is used to calculate the spatially dependent refractive index from the steady-state mass density. Ray optics simulations are carried out to calculate the optical phase shift and Zernike spectrum of the gas vortex lens. The dependence of optical properties on mass flow rate are shown.

2. GAS FLOW PROPERTIES

This paper will cover the specific geometry shown in Fig. 1. All discussions of dimensions will pertain to the interior walls of the gas lens. Experimentally, the bounding structure is 3D printed [4]. The gas lens body is composed of a cylindrical region oriented along the $z$ axis. The body radius and length are 3 and 1 mm. Centered on either side of the body are cylindrical outlets with radius of 0.75 mm and length of 0.5 mm. There are two rectangular inlets with the incoming flow oriented along the $y$ axis. The inlets’ cross-sectional dimensions are 1 mm by 3 mm to match the radius and thickness of the lens body. The axis of each inlet is offset from the $y$ axis by 1.5 mm in the positive and negative x-directions. The top and bottom faces of the inlet region are each positioned above and below the xz-plane by 3.5 mm, such that the resulting flow rotates around the $z$ axis with negative angular velocity.

The origin of our coordinate system is centered in the cylindrical body of the gas lens. Several locations in the gas lens will be referred to by name. A cross-sectional diagram with labels is shown in Fig. 2. The regions are defined as follows: “inlet” is the surface defined by $0\mathrm{mm}<x<3\mathrm{mm}$, $y=3.5\mathrm{mm}$, and $|z|<0.5\mathrm{mm}$; “lens body” is the cylindrical volume defined by $r=\sqrt{x^2+y^2}<3\mathrm{mm}$ and $|z|<0.5\mathrm{mm}$; “outlet” is the cylindrical volume connecting the lens body to the outlet chamber and is defined as $r<0.75\mathrm{mm}$ and $0.5\mathrm{mm}<|z|<1\mathrm{mm}$; “outlet entrances” are the circular surfaces defined as $r<0.75\mathrm{mm}$ and $|z|=0.5\mathrm{mm}$; “outlet exit” is the circular surface defined as $r<0.75\mathrm{mm}$ and $|z|=1\mathrm{mm}$; and “outlet chamber” is the volume that the gas exits into and is defined by $1\mathrm{mm}<|z|$ for all x and y. The “top inlet” is the one located in the $x>0$ and $y>0$ domain.

 

Fig. 2. Cross section of the gas lens in the $x=+0$ plane. The regions referred to by name are labeled. Named volumes and surfaces are shown with a unique color or line style, respectively.

Download Full Size | PDF

The gas vortex is driven by a fixed mass flow rate $\dot{m}$ of nitrogen gas divided evenly across the two inlets. The outlet chambers are meant to approximate flow into open atmosphere. Simulations are done with nitrogen gas exclusively and we do not attempt to model the mixing that occurs when the pure nitrogen gas vortex flow meets open atmosphere as occurs in experiments [4]. This approximation is reasonable considering that the fractional change in refractive index is of the order of $1\times {10}^{-5}$. See Appendix A for details on the equations solved, specific boundary conditions, and additional computational details.

As the gas enters the lens body through the top inlet, it transitions from a uniform flow in the negative y-direction to a rotating flow. The remainder of this section is a summary of typical characteristics of the gas flow through the lens.

Radial line outs of the velocity and number density are shown in Fig. 3 for a mass flow rate of $0.65\mathrm{g}{\mathrm{s}}^{-1}$. The flow is shown in cylindrical coordinates $(r,\varphi ,z)$ at different z-slices starting at $z=0$ in the center of the lens body, moving to $z=0.5\mathrm{mm}$ at the outlet entrance, to $z=1\mathrm{mm}$ at the outlet exit, and then in the outlet chamber out to $z=1.5\mathrm{mm}$. The angular dependence is not shown because the flow is largely axisymmetric.

 

Fig. 3. (a) Angular velocity, (b) radial velocity, (c) axial velocity, and (d) number density of the gas at different z locations. The position $z=0$ is the mid-plane of the gas lens, $z=0.5\mathrm{mm}$ is the outlet entrance, and $z=1\mathrm{mm}$ is the outlet exit. The mass flow rate is $0.65\mathrm{g}{\mathrm{s}}^{-1}$. Line outs were taken along the cylindrical radius $r$, specifically, where $y=0$. The black, dashed line is the average number density from $z=0$ to 1.5 mm.

Download Full Size | PDF

Between $z=0\mathrm{mm}$ and 0.5 mm, the gas enters the lens body at the inlets with a uniform flow in the y-direction. This is seen in the approximately $1/r$-dependence of the angular velocity $\omega$ at $z=0$ in Fig. 3(a). There is a net inward radial flow and the angular velocity is at its maximum. The resulting number density profile, Fig. 3(d), is roughly parabolic with a minimum below the background number density of the outlet chamber. The angular velocity and density behavior is characteristic of a flow such as an adiabatic compressible Rankine vortex or $n=2$ Vatistas and Aboelkassem vortex [26–28].

At $z=0.5\mathrm{mm}$, the angular velocity flattens out for $r<0.2\mathrm{mm}$, but otherwise remains unchanged from $z=0\mathrm{mm}$. This suggests the development of a solid body rotation near the axis due to viscous damping of rotational shear. At this location, the axial velocity at $x=y=0$ reaches its maximum value and develops a parabolic radial profile with a maximum near the walls of the outlet. The radial velocity, Fig. 3(b), transitions from inward to outward flow near $z=0.5\mathrm{mm}$. The radial profile of the number density begins flattening out and the effect of the outlet wall can be seen at $r=0.75\mathrm{mm}$.

Between $z=0.5\mathrm{mm}$ and $z=1\mathrm{mm}$, the angular velocity is observed continuing to flatten out due to viscous damping of rotational shear until it is almost uniform across the outlet. The axial velocity at $x=y=0$, Fig. 3(c), decreases along the outlet while just inside the walls it increases until it reaches its peak value near $z=1\mathrm{mm}$. The increase in the axial velocity is largely due to the pressure difference between the interior and exterior of the outlet. The effect of the viscous boundary layer along the interior wall of the outlet can be seen in the line outs of Figs. 3(a) and 3(c) as the velocity quickly goes to zero near $r=0.75\mathrm{mm}$. The radial velocity does not change significantly and remains an outward flow. The gas density continues to fill in near the axis while decreasing near the walls. At $z=1\mathrm{mm}$ the gas has nearly reached the uniform density of the outlet chamber.

A. Flow Characteristics from Conservation of Mass and Momentum

While the gas flow has viscous boundary layers, turbulence, and generally complex flow structure, some insight can be gained by treating the flow with the following simplifying assumptions: (1) the inlet mass density ${\rho }_i$, pressure $P_i$, and velocity $v_i$ are uniform across the inlet; (2) the flow is isentropic; (3) the outlet exit mass density ${\rho }_o$, pressure $P_o$, axial velocity $v_{o,z}$, and angular velocity $\omega$ are uniform across the outlet; (4) the z-component of the angular momentum $L_{i,o}$, mass, and enthalpy are conserved between the inlet and the outlet; and (5) the outlet density and pressure are equal to that of the outlet chamber. Several of these are directly supported by simulation. The full set of conditions of assumptions (1) and (5) are well supported. The assumption that the angular velocity is uniform across the outlet is also well supported by simulations. Simulations show that the torque on the walls is about a tenth of the angular momentum flux at the inlet supporting the assumption that the angular momentum is conserved.

For a given mass flow rate $\dot{m}$, the inlet velocity is $v_i=\dot{m}/({\rho }_i{A}_i)$, where the inlet area $A_i=2R_i{w}_i$ is the total area of both inlets, $R_i$ is the radius of the lens body, and $w_i$ is the inlet width. The offset of the inlets means the gas enters the device with angular momentum around the $z$ axis. By using assumption (1), the total angular momentum flowing into the device through both inlets in the time $\mathrm{\Delta }t$ is $L_i=\dot{m}{R}_i{v}_i\mathrm{\Delta }t/2={\dot{m}}^2{R}_i\mathrm{\Delta }t/(2{\rho }_i{A}_i)$. All of the gas that flows into the inlets must exit the outlets. The average axial velocity at the outlet is

The inlet mass density ${\rho }_i$ can be estimated as a consequence of the Bernoulli equation. The change in enthalpy is $\mathrm{\Delta }h=v_o^2/2+P_o/{\rho }_o-v_i^2/2-P_i/{\rho }_i$, where the total outlet velocity has an axial and angular component $v_o^2\approx v_{o,z}^2+v_{o,\varphi }^2$, and $P_{i,o}$ is the inlet (outlet) static pressure. The angular component of the velocity is proportional to radius $v_{\varphi }=\omega r$ and requires that the evaluation of the change in enthalpy be done in a radially average sense. The average of the outlet cross section is defined as $⟨·⟩=(2/R_o^2){\int }_0^{R_o}(·)r\mathrm{d}r$. Assuming conditions 1, 2, and 3, the mass density and pressure are related between the inlet and outlet with $P_i/P_o={({\rho }_i/{\rho }_o)}^{\gamma }$, where $\gamma$ is the adiabatic index. The radially averaged change in enthalpy is given by

 

Fig. 4. Estimated area-averaged axial velocity [solid black line from Eq. (1) with ${\rho }_o={\rho }_{\mathrm{atm}}$] and total velocity [dashed red line from Eq. (3)] at the gas lens outlet are shown as functions of the mass flow rate through the device. The density averaged axial velocity (black cross) and total velocity (red plus) from CFD++ simulations are shown for comparison. The gray line marks the speed of sound in the outlet chamber.

Download Full Size | PDF

The uniformity of $v_{o,z}$, assumption (3), is necessary to derive Eq. (2), specifically, $⟨v_z^2⟩={⟨v_z⟩}^2$, which is a consequence of assumption (3). While this approximation is not well justified from Fig. 3(c) alone, if $v_z\propto r^2$, then the uniformity approximation underestimates the second term of Eq. (2) by about 30%. This has a negligible impact on Fig. 4.

B. Mass Flow Rate

The steady-state CFD++ simulations use a constant mass flow rate boundary condition for the inlets. We use the time-dependent SPARC simulations from a previous work [4] to guide the selection of physically significant mass flow rates and provide estimates of the transient time scales. SPARC solves the time-dependent continuity equations for density, momentum, and thermal energy [29]. Details are found in [20,29]. The following are key characteristics of SPARC that are distinct from CFD++: (1) simulations are time-dependent; (2) the turbulence model is only appropriate for laminar to transitional flows; (3) only has inviscid wall boundary conditions; (4) uses a structured mesh, which results in a staircase-representation of curved surfaces; and (5) there are no inflow boundary conditions for constant mass flow rates.

The SPARC simulation is set up with the gas at rest an initial uniform number density and temperature of $2.436\times {10}^{19}{\mathrm{m}}^{-3}$ and 300 K. Attached to each inlet region is a reservoir of gas with high density. Within 1.1 mm of the inlet the number density is $5.6\times {10}^{19}{\mathrm{m}}^{-3}$ and it increases to $2.8\times {10}^{20}{\mathrm{m}}^{-3}$ farther away.

Figure 5 shows the total mass flow rate through both inlets and outlets. The characteristic distance from inlet to outlet is $\sim 4\mathrm{mm}$, which sets the time scale of pressure imbalances across the device to propagate, roughly, $4\mathrm{mm}/340\mathrm{m}{\mathrm{s}}^{-1}=12\mathrm{\mu s}$. By about 65 μs the inlet and outlet mass flow rates have stabilized between 0.5 and $1\mathrm{g}{\mathrm{s}}^{-1}$ (shown as the gray band in Fig. 5). We use these mass flow rates to guide those of the steady-state simulations for comparison.

 

Fig. 5. Total mass flow rate of the inlets (black) and outlets (red). The gray band shows the region with mass flow rates between 0.5 and $1\mathrm{g}{\mathrm{s}}^{-1}$.

Download Full Size | PDF

3. GAS LENS OPTICAL PROPERTIES

The optical properties, specifically, the focal length, astigmatism, and spherical aberration, of the gas lens are calculated using the SeaRay simulation framework [6,30] with the steady-state gas density profiles from CFD++. We utilize the ray-tracing feature of SeaRay to calculate the phase shift caused by the gas lens on a 800 nm wavelength laser beam, such as is typical in ultrashort pulse laser systems.

The gas density profile which does the lensing has a radius of $R=0.75\mathrm{mm}$ set by the radius of the outlet. As seen in Fig. 3(d), the gas density reaches atmospheric conditions close to the outlet exit. Therefore, the thickness of the gas lens is approximately the distance between both outlets, $\sim 2\mathrm{mm}$, with some residual gas outside the lens. Simulations show that the density quickly approaches that of the outlet chamber after the outlet exit. The gas density from the fluid simulations is on an unstructured mesh and only for a quarter-domain ($x>0$ and $z>0$), specifically, within the outlet wall $r<R$ and out to where the density reaches atmospheric conditions $0\mathrm{mm}<z<3\mathrm{mm}$. The density is interpolated onto a regular Cartesian mesh and symmetry is used to determine the density on the full domain of the lens $r<R$ and $|z|<3\mathrm{mm}$. SeaRay was used to calculate a spatially dependent refractive index $n(r,z)$ from the density. The refractivity for nitrogen gas at a wavelength of 800 nm and mass density of ${\rho }^{\star }=1.249\mathrm{kg}{\mathrm{m}}^{-3}$ is $n^{\star }-1=2.9624\times {10}^{-4}$ [31]. From the Gladstone–Dale relation, the density-dependent refractive index is $n(r,z)=1+(n^{\star }-1)\rho (r,z)/{\rho }^{\star }$.

A collimated (flat phase fronts), Gaussian beam is simulated entering the gas density profile and the phase is calculated directly after the gas lens. The phase fronts directly after the gas lens allows for comparison with an ideal spherical lens (see Appendix C). From the phase, the Zernike coefficients (see Appendix B) are calculated and translated into meaningful quantities, such as focal length or spherical aberration. Figure 6 shows that the Zernike spectrum for a mass flow rate of $1\mathrm{g}{\mathrm{s}}^{-1}$. The spectrum is labeled by the single-index $j$ specified by the Optical Society of America (OSA) [32] and details are given in Appendix B. The defocus mode $j=4$ and spherical aberration $j=12$ are the most important terms. The astigmatism is represented by $j=3,5$ and is small with the difference in focal lengths being less than 8%. The piston mode $j=0$ is excluded from Fig. 6 because it is insignificant to optical propagation.

 

Fig. 6. Zernike spectrum of the optical phase shift induced by a gas vortex lens with mass flow rate of $1\mathrm{g}{\mathrm{s}}^{-1}$. The OSA single-index $j$ is used [32]. The defocus and spherical aberration modes are $j=4$ and $j=12$.

Download Full Size | PDF

Figure 7 shows how the focal length of the diverging gas lens decreases with increasing mass flow rate. This is because with a higher mass flow rate, there is a larger angular momentum entering the inlet and therefore large angular velocity. The larger centrifugal forces create larger density gradients. For a mass flow rate of $1\mathrm{g}{\mathrm{s}}^{-1}$ the focal length is 0.601 m. This is consistent with experimental observations [4] which were made at 532 nm because, while the Zernike coefficients are inversely proportional to wavelength, the focal length’s wavelength dependence is from the refractivity, which only varies by $\sim 1\%$ between 800 nm and 532 nm.

 

Fig. 7. Focal length of the diverging lens as a function of mass flow rate.

Download Full Size | PDF

The main optical aberration is spherical. For the focal lengths and lens radii of interest ($f\sim 1\mathrm{m}$ and $R\sim 1\mathrm{mm}$), an ideal spherical lens has a very small contribution from the $Z_4^0(r,\varphi )$ ($j=12$) Zernike polynomial. At 800 nm wavelengths, the relative contribution of the spherical aberration term to the defocusing is about $1\times {10}^{-7}$. Therefore, any significant $r^4$-dependence to the phase is largely a deviation from an ideal spherical lens. This deviation results in paraxial and marginal rays that focus at different points. Figure 8 shows the difference between the focal lengths of the marginal and paraxial rays to the focal length of best focus. For the mass flow rate of $0.65\mathrm{g}{\mathrm{s}}^{-1}$, the spherical aberration is $\mathrm{\Delta }{f}_s=4.8\mathrm{m}$ (see Appendix C for definition). For comparison, focal length in this case is $-1.8\mathrm{m}$. Specifically, the virtual point of best focus is located 1.8 m before the gas lens. The paraxial rays have a virtual focus between the point of best focus and the gas lens. The marginal rays have a real focus after the gas lens. The axially averaged, number density, shown as the black dashed line of Fig. 3(d), begins decreasing for $r>0.7\mathrm{mm}$. Rays outside this radius will experience a weak focusing effect which contributes to the large spherical aberration. The large aberration at $0.2\mathrm{g}{\mathrm{s}}^{-1}$ is due to a weaker positive focusing effect pushing the focus of the marginal rays farther after the gas lens. Practically, the impact of the large spherical aberration can be minimized by not utilizing the full lens diameter.

 

Fig. 8. Fractional difference in focal length between paraxial and margin rays is shown as a function of mass flow rate.

Download Full Size | PDF

4. CONCLUSIONS

In conclusion, we have shown via steady-state simulations that a rotating gas flow can be used to create a compact negative lens. The focal length can be controlled by the mass flow rate through the inlets. A larger mass flow rate results in a faster flow, larger angular momentum, and shorter focal length. Spherical aberrations are the dominant deviation from an ideal lens for the geometry and flow rates investigated. Future work is needed to analyze the temporal stability of the lens’ optical properties, determine the dominant mechanism (specifically geometric features) for controlling the optical aberrations, and to better understand the limits of obtainable focal lengths.

APPENDIX A: CFD++ SIMULATIONS

We used the software package CFD++ 17.1.1 [33], which was developed and is maintained by Metacomp Technologies.

1. Simulation Geometry

For modeling purposes, the outlet boundaries conditions must be moved away from the physical outlet of the gas lens. Therefore, outlet chambers are connected to both outlets shown in Fig. 1. The outlet chamber is connected to the lens outlet with a conic frustum with half-angle of 80 deg. The outlet chamber length was always equal to twice its radius.

The geometry has a twofold rotational symmetry around the $z$ axis and mirror symmetry across the xy-plane and is shown in Fig. 9. These symmetries are reflected by solving in a volume restricted to the positive-x, positive-z quadrant.

 

Fig. 9. Simulation domain is one quadrant of the physical domain and is the union of the gas vortex lens and an outlet chamber.

Download Full Size | PDF

2. Equation Solved

The Navier–Stokes fluid equations [34] were solved using a compressible (non-preconditioned) density-based method:

(A1)$$\frac{\partial }{\partial t}\rho +\frac{\partial }{\partial x_k}(\rho v_k)=0,$$

(A2)$$\frac{\partial }{\partial t}(\rho v_i)+\frac{\partial }{\partial x_k}(p{\delta }_{ik}+\rho v_i{v}_k-{\sigma }_{ik}^{\prime })=0,$$

(A3)$$\frac{\partial }{\partial t}E+\frac{\partial }{\partial x_k}(v_k(E+p)-v_i{\sigma }_{ik}^{\prime }-k_i\frac{\partial }{\partial x_i}T)=0,$$

where $p$ is the pressure, $\rho$ is the mass density, $v_k$ is the $k$th component of the velocity, ${\sigma }_{ik}^{\prime }$ is the viscous stress tensor, $E$ is the total energy density, $T$ is the temperature, and $k_i$ is the thermal conductivity. The viscous stress tensor is defined as

(A4)$${\sigma }_{ik}^{\prime }=\mu (\frac{\partial v_i}{\partial x_k}+\frac{\partial v_k}{\partial x_i}-\frac{2}{3}{\delta }_{ik}\frac{\partial v_l}{\partial x_l})-\zeta {\delta }_{ik}\frac{\partial v_l}{\partial x_l},$$

with $\mu$, $\zeta$ being the coefficients of viscosity. CFD++ assumes that $\zeta =0$. The equation of state was that of a compressible prefect gas with the pressure given by $p=\rho R_{s}T$, where $R_s=c_p-c_v$, $R_s=k_B/m_{\mathrm{gas}}$ is the specific gas constant, $k_B=1.3807\times {10}^{-23}\mathrm{J}{\mathrm{K}}^{-1}$ is the Boltzmann constant, and $m_{\mathrm{gas}}$ is the molecular mass. The total energy density is given by $E=\rho v^2/2+p/(\gamma -1)$, where $\rho v^2$ is the dynamic pressure, $p/(\gamma -1)$ is the internal energy density of an ideal gas, $\gamma =c_p/c_v$ is the adiabatic index, and $c_p$ ($c_v$) is the specific heat at constant pressure (volume). The first viscosity coefficient and the thermal conductivity are modeled by fits to Sutherland’s law. The viscosity is given by

(A5)$$\frac{\mu }{{\mu }_0}={\left(\frac{T}{T_{0\mu }}\right)}^{3/2}\frac{T_{0\mu }+S_{\mu }}{T+S_{\mu }},$$

where ${\mu }_0$, $T_{0\mu }$ are the reference viscosity and temperature and $S_{\mu }$ is the viscosity Sutherland constant. The thermal conductivity is given by

(A6)$$\frac{k}{k_0}={\left(\frac{T}{T_{0k}}\right)}^{3/2}\frac{T_{0k}+S_k}{T+S_k},$$

where $k_0$, $T_{0k}$ are the reference viscosity and temperature and $S_k$ is the conductivity Sutherland constant.

For nitrogen gas, the Sutherland law parameters for viscosity and conductivity are ${\mu }_0=1.656\times {10}^{-5}\mathrm{kg}{\mathrm{m}}^{-1}{\mathrm{s}}^{-1}$, $T_{0\mu }=273.11\mathrm{K}$, $S_{\mu }=104.7\mathrm{K}$, $k_0=2.407\times {10}^{-2}\mathrm{W}{\mathrm{m}}^{-1}{\mathrm{K}}^{-1}$, $T_{0k}=273.11\mathrm{K}$, and $S_k=178.1\mathrm{K}$. At 293 K, the viscosity is $1.74813\times {10}^{-5}\mathrm{kg}{\mathrm{m}}^{-1}\mathrm{s}$. The molecular mass is $m_{\mathrm{gas}}=4.653\times {10}^{-26}\mathrm{kg}$. For a diatomic gas, the specific heat at constant pressure (volume) is $c_p=7R_s/2$ ($c_v=5R_s/2$) and the adiabatic index is $\gamma =1.4$.

3. Boundary Conditions

The boundary condition across the inlet (brown in Fig. 9) is “inflow-mass flow rate and temperature,” which assumes a fixed mass flow rate $\dot{m}$, constant temperature $T_{\mathrm{in}}=293\mathrm{K}$, turbulence kinetic energy $k_{\mathrm{in}}$, and dissipation ${ϵ}_{\mathrm{in}}$. At this boundary, the following condition is met: $\dot{m}={\int }_{\mathrm{inlet}}\rho \mathbf{v}·\mathrm{d}\mathbf{a}$, where $\rho$ is calculated from the equation of state $p=\rho RT$ using $T=T_{\mathrm{in}}$. The flow is assumed to be normal to the inlet surface. The pressure far from the lens outlet must return to atmospheric conditions, specifically, pressure $P_{\mathrm{atm}}=101325\mathrm{Pa}$, temperature $T_{\mathrm{atm}}=293\mathrm{K}$, and mass density ${\rho }_{\mathrm{atm}}=1.17\mathrm{kg}{\mathrm{m}}^{-3}$. The outlet surface (orange in Fig. 9) uses the “outflow-back pressure imposition” boundary condition from CFD++ with the reference pressure being $P_{\mathrm{atm}}=101325\mathrm{Pa}$. The xy-plane defines the mirror plane for the lens (purple in Fig. 9) and the boundary condition on this surface reflects that symmetry. The yz-plane divides the simulation domain such that it has twofold rotational symmetry around the $z$ axis. As a result, a zonal, “simple flow-through” boundary condition is used to match the $y<0$, yz-half-plane (red in Fig. 9) with a $\pi$-rotation around the $z$ axis to the $y>0$, yz-half-plane (green in Fig. 9). Along the walls of the device (blue in Fig. 9), a viscous (non-slip) boundary condition is used. Heat transfer at the walls is assumed to be adiabatic (zero heat flux). Wall functions are used to model the boundary layer that results from the non-slip boundary condition. CFD++ has the option to use a standard Launder–Spalding [35] or an advanced two-layer approach. All simulations shown were carried out with the advanced two-layer model using non-equilibrium wall types. This approach improved the convergence characteristic without changing the nature of the flow.

4. Turbulence Model

The estimated Reynolds number suggests that the flow is in the transitional regime between laminar and turbulent flow. This transition occurs for Poiseuille’s law for Reynolds numbers of $1\times {10}^3$ to $7.5\times {10}^4$ [36]. The Reynolds number Re is estimated by $\mathrm{Re}\sim \rho uL/\mu =\dot{m}L/(\mu A)$, where $\rho$, $u$, $L$, $\mu$, $\dot{m}$, and $A$ are the characteristic mass density, velocity, length, viscosity, mass flow rate, and cross-sectional area. Consider an approximate mass flow rate of $0.65\mathrm{g}{\mathrm{s}}^{-1}$. The inlet area is $6\times {10}^{-6}{\mathrm{m}}^2$, the characteristic length is the diameter of the lens body $L\sim 6\mathrm{mm}$, and the resulting Reynolds number is roughly $4\times {10}^4$. Near the outlet, the area is $3.5\times {10}^{-6}{\mathrm{m}}^2$ and the length is the diameter of the outlet $L\sim 1.5\mathrm{mm}$ and the resulting Reynolds number is roughly $2\times {10}^4$. It is clear from these Reynolds numbers that a turbulent flow must be modeled. A two-equation, realizable k-epsilon model [37] is used for modeling turbulence. While turbulence is not significant everywhere, it helped overall to model convergence by damping oscillatory flows with turbulence diffusion.

5. Numerical Settings

The steady-state simulations are carried out with a point-implicit time integration scheme with a localized Courant condition. The Courant number $C$ varied between 50 and 100 but is set using the “CFD++ Numerics Wizard,” which initializes numerical parameters based on the characteristic Mach speed. The peak Mach number occurs near the lens outlet. A simple estimate can be had by assuming atmospheric density near the outlet. The Mach number is then given by $M\approx \dot{m}/(c_0{\rho }_{\mathrm{atm}}{A}_o)$, where $c_0=\sqrt{\gamma R_s{T}_{\mathrm{atm}}}=349\mathrm{m}{\mathrm{s}}^{-1}$ is the speed of sound of nitrogen, ${\rho }_{\mathrm{atm}}=1.17\mathrm{kg}{\mathrm{m}}^{-3}$ is the mass density, and $A\approx 3.5\times {10}^{-6}{\mathrm{m}}^2$ is the outlet area. For $1\mathrm{g}{\mathrm{s}}^{-1}$, the Mach number is about 0.7, which is in the transonic regime. The Courant number is ramped from 1 to 75 over the first 100 steps. For this method each cell has its own time step based off the Courant number, cell volume, and maximum modulus eigenvalue of the cell. In regions where the mesh has rapid variation in cell size it is important to use the time-step spatial smoothing. This limits the possibility of loss of conservation between cells.

6. Mesh Generation in CUBIT

An unstructured, uniform, tetrahedral mesh is generated using the TetMesh scheme in Cubit 15.2. Cubit is a software package from Sandia National Laboratory for three-dimensional geometry and mesh generation. Cubit’s TetMesh scheme employs a third-party library, MeshGem, to carry out mesh generation. The approximate sizes of mesh elements with the lens body and outlet are 0.2 mm, 0.1 mm, or 0.05 mm. The smallest mesh elements are always located on the wall of the outlet where the sharpest flow gradients are found or along the $z$ axis, which is essential for lensing properties. The mesh size on the boundary of the outlet chamber is 1 mm and there is a smooth variation to the outlet exit.

A boundary layer mesh is grown along the walls with an initial layer thickness of 0.01 mm, growth rate of 1.2, and total thickness of approximately 0.1 mm. An example of the mesh elements on the yz-plane near the center of the lens body and outlet are shown in Fig. 10.

 

Fig. 10. Example of mesh ($\sim 50\mathrm{\mu m}$) near the outlet on the yz-plane.

Download Full Size | PDF

7. Residual, Mesh, and Geometry Convergences

Numerical convergence was based off the sum of all residuals (excluding turbulence) dropping by 5 orders of magnitude. Typically, all residuals have dropped by at least 4 orders of magnitude.

The mesh size around the outlet is varied from 0.2 to 0.05 mm to ensure that the density profile does not depend on the mesh. Figure 11 shows the percent difference between a 0.1 and 0.05 mm mesh size in the axially averaged density. The density is averaged over the domain of the lens body and outlet, specifically, $|z|<1\mathrm{mm}$.

 

Fig. 11. Percent difference between radial density profiles between the 100 and 50 μm simulations.

Download Full Size | PDF

Finally, the gas density and velocity are insensitive to when the outlet chamber radius is at least 8 mm. This is necessary to accurately model gas expansion into an open domain.

APPENDIX B: ZERNIKE POLYNOMIALS

Decomposition of the optical phase into an expansion of Zernike polynomials provides a useful way to characterize the optical properties of a lens when the aberrations are limited. A Zernike expansion is a complete set on the unit circle $\rho \in [0,1]$ and $\theta \in [0,2\pi )$ and is defined by

(B1)$$Z_n^m(\rho ,\theta )=\{\begin{array}{cc}{R}_n^{|m|}(\rho )\sin (|m|\theta ),& m<0\\ R_n^{|m|}(\rho )\cos (|m|\theta ),& m\ge 0\end{array},$$

with $R_n^m(\rho )$ being the radial dependence, $n\ge 0$ is the radial index, and $-n\le m\le n$ is the azimuthal index. The radial polynomial is defined as

(B2)$$R_n^m(\rho )=\sum _{k=0}^{(n-m)/2}\frac{{(-1)}^k(n-k)!}{k![(n+m)/2-k]![(n-m)/2-k]!}{\rho }^{n-2k},$$

for when $n-m$ is even, and $R_n^m(\rho )=0$ otherwise. A table of the first 15 Zernike polynomials and their descriptions is provided in Table 1.

Table 1. Functional Form of Zernike Polynomials on the Unit Circle, $\rho \in [0,1]$ and $\theta \in [0,2\pi )$, Their OSA Single-Index $j$, and Description

View Table

The Zernike polynomials satisfy the orthogonality condition

(B3)$${\int }_0^1\mathrm{d}\rho \rho {\int }_0^{2\pi }\mathrm{d}\theta Z_n^m(\rho ,\theta )Z_{n^{\prime }}^{m^{\prime }}(\rho ,\theta )=\frac{{ϵ}_m\pi }{2n+2}{\delta }_{n,n^{\prime }}{\delta }_{m,m^{\prime }},$$

where $n-m$ and $n^{\prime }-m^{\prime }$ are both even, ${ϵ}_m=2$ for $m=0$ and ${ϵ}_m=1$ otherwise, and ${\delta }_{n,n^{\prime }}$ is the Kronecker delta function. The decomposition of the phase $\mathrm{\varphi }(\rho ,\theta )$ on the unit circle is

(B4)$$\mathrm{\varphi }(\rho ,\theta )=\sum _{n=0}^{\infty }\sum _{\genfrac{}{}{0ex}{}{m=-n,}{n-m\mathrm{even}}}^n{s}_{n,m}{Z}_n^m(\rho ,\theta ),$$

where $s_{n,m}$ is the Zernike coefficient. The Zernike coefficient can be interpreted as that polynomial peak contribution to the phase shift at the edge of the unit circle because $R_n^m(1)=1$. The Zernike coefficients can be calculated as follows:

(B5)$$s_{n,m}=\frac{2n+2}{{ϵ}_m\pi }{\int }_0^1\mathrm{d}\rho \rho {\int }_0^{2\pi }\mathrm{d}\varphi \mathrm{\varphi }(\rho ,\theta )Z_n^m(\rho ,\varphi ).$$

For the purpose of plotting the Zernike spectrum, a single-index scheme is used. The single-index $j=[n(n+2)+m]/2$set forth as a standard by the OSA uniquely specifics each Zernike polynomial instead of $n$ and $m$ [32].

APPENDIX C: SPHERICAL LENSES WITH SMALL ABERRATIONS

An ideal spherical lens is one that creates spherical fronts of constant phase. The phase for a spherical wave is given by $\mathrm{\varphi }=k\sqrt{r^2+z^2}$, where $k$ is the wavenumber, and $\sqrt{r^2+z^2}$ is the distance from the focal point. After a defocusing ideal spherical lens, the phase is such that the distance from the virtual focus would be the focal length $f$. The phase fronts would be given by $\mathrm{\varphi }(\rho )=kf\sqrt{1+{\rho }^2{R}^2/f^2}$, where $\rho =r/R$ is the radius normalized to the radius of the lens $R$.

The Zernike coefficients for an ideal lens can be calculated analytically but practically for the case of $R\ll f$ a series expansion in small $R$ is all that is needed. To lowest order in the lens radius, the first three non-zero Zernike coefficients for the ideal lens are $s_{0,0}=-kf-k{R}^2{(4f)}^{-1}(1-R^2/(6f^2)+O^6(R/f))$ (piston), $s_{2,0}=-k{R}^2{(4f)}^{-1}(1-R^2/(4f^2)+O^6(R/f))$ (defocus), and $s_{4,0}=k{R}^4/(48f^3)+O^6(R/f)$ (spherical aberration). To lowest order in $R$, $s_{2,0}$ provides an estimate of the focal length of the lens using $f=-k{R}^2/(4s_{2,0})$, where negative (positive) $f$ means it is a defocusing (focusing) lens. For $k=2\pi /800\mathrm{nm}$, $R=0.75\mathrm{mm}$, and $s_{2,0}=1\mathrm{rad}$, then $f=-1.1\mathrm{m}$. The expansion for small $R/f$ is appropriate given that $R/f\sim -7\times {10}^{-4}$.

The longitudinal spherical aberration can also be estimated from the Zernike coefficients. The longitudinal spherical aberration is the distance between the focal planes of rays close to the optical axis (paraxial rays) and rays close to the edge of the lens (marginal rays). Consider the phase $\mathrm{\varphi }=s_{2,0}{Z}_2^0(\rho ,\theta )+s_{4,0}{Z}_4^0(\rho ,\theta )$, where the defocusing coefficient $s_{2,0}=-k{R}^2/(4f_0)$ is that of the ideal spherical lens. The radial wavenumber is $k_r=\partial \mathrm{\varphi }/\partial r$, which when combined with the magnitude of the wavenumber $k$ tells us the angle $\varphi$ that the ray makes with the optical axis $\sin \varphi =k_r/k$. For a nearly spherical wave, the ratio of the radius at which the ray exits the lens $r$ to focal length $f(r)$ is also $\sin \varphi =r/|f(r)|$. The sign of $f(r)$ can be resolved by noticing that for a defocusing wave, $f(r)<0$ but $k_r>0$. Therefore, the radius-dependent focal length is $f(r)=-kr/k_r(r)=-(k{R}^2/4)/(s_{2,0}+3(2r^2/R^2-1)s_{4,0})$. If we evaluate $f(r)$ for paraxial rays, $f(0)/f_0=1+3{\alpha }_s/(1-3{\alpha }_s)$, and marginal rays, $f(R)/f_0=1-3{\alpha }_s/(1+3{\alpha }_s)$, where ${\alpha }_s=s_{4,0}/s_{2,0}$. For $1/3>{\alpha }_s>0$, the paraxial rays have a longer focal length than the marginal rays. For $-1/3<{\alpha }_s<0$, the paraxial rays have a shorter focal length than the marginal rays. The change in focal length between the paraxial and margin rays (longitudinal spherical aberration) is $\mathrm{\Delta }{f}_s/f_0=(f(R)-f(0))/f_0=-6{\alpha }_s/(1-9{\alpha }_s^2)$.

Funding

Office of Naval Research (ONR) (N0001417WX01788).

Acknowledgment

Resources of the Department of Defense High Performance Computing and Modernization Program (HPCMP) were used in this work. We would like to thank the generous support from Metacomp’s support staff in troubleshooting issues with CFD++ model setup and convergence.

REFERENCES

1. H. Bercegol, P. R. Bouchut, L. Lamaignère, B. Le Garrec, and G. Razé, “The impact of laser damage on the lifetime of optical components in fusion lasers,” Proc. SPIE 5273, 312–325 (2003). [CrossRef]  

2. J. Peñano, P. Sprangle, A. Ting, R. Fischer, B. Hafizi, and P. Serafim, “Optical quality of high-power laser beams in lenses,” J. Opt. Soc. Am. B 26, 503–510 (2009). [CrossRef]  

3. P. Sprangle, B. Hafizi, A. Ting, and R. Fischer, “High-power lasers for directed-energy applications,” Appl. Opt. 54, F201–F209 (2015). [CrossRef]  

4. D. Kaganovich, L. Johnson, D. Gordon, A. Mamonau, and B. Hafizi, “Lensing properties of rotational gas flow,” Appl. Opt. 57, 9392–9396 (2018). [CrossRef]  

5. W. Liu, J.-F. Gravel, F. Théberge, A. Becker, and S. Chin, “Background reservoir: its crucial role for long-distance propagation of femtosecond laser pulses in air,” Appl. Phys. B 80, 857–860 (2005). [CrossRef]  

6. D. Gordon, A. Stamm, B. Hafizi, L. Johnson, D. Kaganovich, R. Hubbard, A. Richardson, and D. Zhigunov, “Ideal form of optical plasma lenses,” Phys. Plasmas 25, 063101 (2018). [CrossRef]  

7. D. Marcuse and S. Miller, “Analysis of a tubular gas lens,” Bell Syst. Tech. J. 43, 1759–1782 (1964). [CrossRef]  

8. A. Beck, “Thermal gas lens measurements,” Bell Syst. Tech. J. 43, 1818–1820 (1964). [CrossRef]  

9. D. Berreman, “A lens or light guide using convectively distorted thermal gradients in gases,” Bell Syst. Tech. J. 43, 1469–1475 (1964). [CrossRef]  

10. A. Beck, “Gas mixture lens measurements,” Bell Syst. Tech. J. 43, 1821–1825 (1964). [CrossRef]  

11. D. Berreman, “A gas lens using unlike, counter-flowing gases,” Bell Syst. Tech. J. 43, 1476–1479 (1964). [CrossRef]  

12. M. Notcutt, M. Michaelis, P. Cunningham, and J. Waltham, “Spinning pipe gas lens,” Opt. Laser Technol. 20, 243–250 (1988). [CrossRef]  

13. M. Michaelis, M. Notcutt, and P. Cunningham, “Drilling by gas lens focused laser,” Opt. Commun. 59, 369–374 (1986). [CrossRef]  

14. M. Michaelis, C. Dempers, M. Kosch, A. Prause, M. Notcutt, P. Cunningham, and J. Waltham, “A gas-lens telescope,” Nature 353, 547–548 (1991). [CrossRef]  

15. C. Mafusire, A. Forbes, G. Snedden, and M. Michaelis, “The spinning pipe gas lens revisited,” South African J. Sci. 104, 260–264 (2008).

16. J. Urzay, M. Mortazavi, and A. Mani, “Optical degradation of inductively-coupled plasma lenses by conversion of electromagnetic energy into unsteady flows,” Tech. Rep. (Center for Turbulence Research, 2013).

17. M. Mortazavi, J. Urzay, and A. Mani, “Computational hydrodynamics and optical performance of inductively-coupled plasma adaptive lenses,” Phys. Plasmas 22, 062110 (2015). [CrossRef]  

18. B. W. Neiswander, E. Matlis, and T. C. Corke, “Plasma lens for optical path difference control,” AIAA J. 50, 123–130 (2012). [CrossRef]  

19. B. W. Neiswander, E. Matlis, and T. C. Corke, “Geometric optimization of a cylindrical plasma adaptive optics lens,” AIAA J. 51, 657–664 (2013). [CrossRef]  

20. D. Kaganovich, D. Gordon, M. Helle, and A. Ting, “Shaping gas jet plasma density profile by laser generated shock waves,” J. Appl. Phys. 116, 013304 (2014). [CrossRef]  

21. R. J. McConnel, “Method and apparatus for refracting a laser beam,” U.S. patent 4,402,574 (September 6, 1983).

22. P. Mikheev, V. Nikolaev, S. Shalaginov, and A. Shepelenko, “Investigation of the defocusing properties of a vortical gas flow,” J. Eng. Phys. 59, 1375–1379 (1990). [CrossRef]  

23. B. Vilenchits, A. Zhdanovskii, and D. Umreiko, “Focusing properties of a vortex gas flow,” J. Appl. Spectrosc. 62, 153–155 (1995). [CrossRef]  

24. B. Vilenchits, A. Zhdanovskii, N. Lemesh, and L. Senchuk, “The effect of the scale factor on the focusing properties of a vortex gas flow,” J. Appl. Spectrosc. 64, 399–402 (1997). [CrossRef]  

25. B. Vilenchits, A. Zhdanovskii, and V. Popov, “Influence of a vortex gas flow on an axial laser beam,” J. Appl. Spectrosc. 72, 59–63 (2005). [CrossRef]  

26. Y. Aboelkassem and G. H. Vatistas, “New model for compressible vortices,” J. Fluids Eng. 129, 1073–1079 (2007). [CrossRef]  

27. G. H. Vatistas, V. Kozel, and W. Mih, “A simpler model for concentrated vortices,” Exp. Fluids 11, 73–76 (1991). [CrossRef]  

28. H. Katanoda and M. H. Bin Yusof, “Energy separation mechanism in uni-flow vortex tube using compressible vortex flow,” in International Conference on Fluid Mechanics, Heat Transfer and Thermodynamics (2014), Vol. 2014, pp. 1252–1255.

29. D. Gordon, P. Sprangle, S. Slinker, R. Fernsler, and M. Lampe, “SPARC-a simulation model for electrical charges,” Tech. Rep. NRL/MR/6790–06-8974 (U.S. Naval Research Lab, 2006).

30. D. Gordon, “SeaRay,” 2018, https://github.com/USNavalResearchLaboratory/SeaRay. Commit: d4f58d5ed5310aadafccff958c4b927bec0da93c.

31. E. R. Peck and B. N. Khanna, “Dispersion of nitrogen,” J. Opt. Soc. Am. 56, 1059–1063 (1966). Equation (4) of this reference is used. [CrossRef]  

32. L. Thibos, R. A. Applegate, J. T. Schwiegerling, and R. Webb, “Standards for reporting the optical aberrations of eyes,” in Vision Science and Its Applications, OSA Technical Digest (Optical Society of America, 2000), paper SuC1.

33. S. Chakravarthy, “A unified-grid finite volume formulation for computational fluid dynamics,” Int. J. Numer. Methods Fluids 31, 309–323 (1999). [CrossRef]  

34. L. Landau and E. Lifshitz, Fluid Mechanics, Vol. 6 of Course of Theoretical Physics (Pergamon, 1959).

35. B. E. Launder, “On the computation of convective heat transfer in complex turbulent flows,” J. Heat Transfer 110, 1112–1128 (1988). [CrossRef]  

36. E. U. Condon and H. Odishaw, Handbook of Physics (McGraw-Hill, 1958).

37. U. Goldberg, O. Peroomian, and S. Chakravarthy, “A wall-distance-free k-ε model with enhanced near-wall treatment,” J. Fluids Eng. 120, 457–462 (1998). [CrossRef]