Gas dynamics and vorticity generation in laser-induced breakdown of air

Research has shown that the ignition characteristics of laser induced plasmas in fuel-air mixtures are influenced by the gas dynamics effects induced during the gas breakdown stage. Here, we present numerical modeling of the fluid mechanics induced by breakdown (plasma formation) from a nanosecond near-infrared (NIR) laser pulse in air. The simulations focus on the post-discharge kernel dynamics with the goal of developing a better understanding on how vorticity is generated during the kernel cooling phase. Initial conditions (ICs) of kernel shape, temperature and pressure (corresponding to the end of the laser pulse) are found from experimental Rayleigh scattering data. It is shown that this method for determining ICs is preferred versus use of the Taylor-Sedov blast wave theory as it provides a more accurate description of the starting field. Past experimental observations have revealed that the gas dynamics of nanosecond laser sparks typically lead to the formation of an asymmetric torus with a frontal lobe propagating towards the laser source. We show that the development of the asymmetric torus is governed by strong vorticity generated through baroclinic torque arising from the blast wave that forms at the kernel boundary. Initially, the blast takes the shape of the teardrop kernel but then evolves into a spherical front during the first ~10 µs because the blast wave strength varies along its circumference. This spatial variation leads to a misalignment between the pressure and density gradients and generation of vorticity by baroclinic torque. Ultimately, the observed flow-field is dictated by how the energy was initially deposited around the beam waist during breakdown. As such, one can tailor the aerodynamics induced during the cooling and recombination phase by controlling the energy deposition profile.


I. Introduction
Laser induced breakdown has been studied for a wide range of combustion applications including as an ignition source for reciprocating engines [1][2][3], gas turbines [4][5][6] and rocket engines [7,8], and as a means to provide flame stabilization in high speed diffusion jet flames [9][10][11][12]. For widely used nanosecond pulses, a typical breakdown event is governed by two processes. An initial pool of free electrons is formed by multi-photon ionization (MPI) whereby a gas molecule becomes ionized through the absorption of multiple-photons within the molecule's excitation lifetime. This is followed by the plasma growth stage via electron avalanche ionization (EAI). During EAI, the free electrons generated by MPI are accelerated in the presence of the laser's electric field (i.e., inverse bremsstrahlung absorption) and generate further gas ionization through collisions with the nearby gas molecules. During this period, the plasma rapidly reaches high temperatures (~30,000 K) and high pressures (several tens of bar) [13,14]. The initial plasma kernel has a characteristic teardrop shape with the axial length governed by the focusing conditions and beam energy. The plasma kernel grows with high velocity (~0.1-1 km/s) during the first microsecond after the pulse [15]. An interesting feature of the kernel's subsequent development is the presence of vorticity within the plasma even when the breakdown was generated in an initially quiescent gas. Vorticity in the fluid leads to the collapse of the teardrop kernel along its optical axis and, ultimately, to the formation of an asymmetric torus with a lobe feature pointing in the direction of the laser, at several tens of microseconds later. When experimentally viewed by (planar) 2-D imaging methods, the flow structure appears as a "three-lobe" pattern due to the persisting torus and frontal lobe (also referred to as "third lobe") [16][17][18]. Shown in Fig. 1 is an example of a kernel investigated as part of this work which has collapsed into an asymmetric torus at 100 µs after the laser pulse. The gas dynamics involved in the breakdown process, in particular the vorticity generation mechanism leading to the creation of the three-lobe structure, have not been sufficiently explained in the literature. The formation mechanism of the third lobe was first discussed by Bradley et al. [15] who posited that the third lobe forms due to the interaction between the rarefaction waves (which constitute the tail end of the blast wave) and the hot gas kernel present at the beam focus. The inward flow induced by the passing rarefaction during first ~5-10 µs is key to the mechanism proposed by Bradley. Similar to a Taylor instability, the gas of lower density (hotter) is accelerated more by this outward pressure gradient than the neighboring layer of gas that has higher density (colder). The resulting torque generates a pair of counter-rotating vortices-one at the upstream (laser incident) side and one at the downstream side of the kernel. This (nearly adjacent) vortex pair causes the hot gas kernel left after breakdown to collapse into a torus. The imbalance of vortex strength between the pair also generates the third (frontal) lobe. However, no actual measurements or computational fluid dynamics (CFD) simulations were performed to confirm this mechanism. Moreover, detailed analysis of the kernel topology has suggested that vorticity might be present at even earlier times, before the blast wave has separated from the kernel. In another study, Dors et al. [19] conducted a CFD analysis of a breakdown kernel from a 10 ns Nd:YAG laser (1064 nm) pulse with Epulse=25 mJ. Their study assumed that the plasma is in equilibrium and the initial kernel profile was determined from experimental observations. Interestingly, while their results predict the formation of an asymmetric toroid, the direction of propagation of the third lobe was towards the cooler side of the post-breakdown kernel which is contrary to what has been experimentally reported [20,21]. This suggests that obtaining the correct direction of fluid rotation induced by vorticity is very sensitive to the initial conditions.
It is important to note the significant role of vorticity for practical applications including the viability of flame ignition. Endo et al. [22] indicated that the plasma-driven fluid dynamics dictates the flame kernel topology. Furthermore, they suggest that flame vorticity entrains the surrounding combustible mixture which leads to an increase of the effective kernel energy in the early stages of flame development. A comparative experimental study between single-and dual-pulse laser breakdown ignition by the present authors showed that the plasma-induced vorticity is responsible for the flame kernel extinction at the lean limit [23]. Moreover, it was shown that by controlling how the energy is deposited in the gas along the optical axis, one can control the direction of fluid rotation including changing the direction of the third lobe (i.e. making it appear towards or away from the incoming beam) or even completely suppressing it. Therefore, understanding the mechanisms through which vorticity is generated in a laser spark (and tailoring it) can have a wide range of benefits such as reducing flame stretching in highly turbulent combustion devices (e.g. aero-turbines), flame stabilization in ramjet or scramjet propulsion systems, increase in combustion efficiency for stationary natural gas engines, and optimization of plasma aerodynamic flow control devices.
In this paper we use computational fluid dynamics to investigate the gas dynamic effects induced by a nanosecond laser spark in atmospheric air. Initial conditions (after the end of the incident pulse) are found from experimental studies. The temporal evolution of the flow field and the corresponding kernel dynamics are modeled under the 2-D axis-symmetric assumption. A vorticity equation analysis shows that vorticity is generated mainly by baroclinic torque due to misalignment between the pressure and density gradients within the kernel. Pressure contours show that the strength of the blast wave varies along its circumference during the first ~10 µs after the laser pulse, which results in a misalignment of the pressure gradient with respect to the density gradient (which follows the temperature gradient). We find that, ultimately, the initial shape of the breakdown kernel dictates the subsequent gas dynamics. The teardrop kernel gives rise to two zones of significant recirculation at the front (facing the laser source) and the back (away from the laser source) of the kernel that have different strengths. These two vortices govern the collapse of the laser kernel into the three-lobe structure observed experimentally. These findings provide new insight into the vorticity generation, and its implications to the flow field, relative to the current literature.

II. Mathematical Model
The fluid mechanics generated by the laser breakdown process were numerically modeled using a Riemann solver developed in-house based on the Navier-Stokes equations. At this stage, the model neglects the effects of plasma kinetics and focuses exclusively on the fluid dynamics induced by the hot gas kernel in a non-reactive air gas mixture (assumed to be 79% N2 and 21% O2). As such, the simulations presented here begin after the laser pulse has ended with initial conditions determined from experiments (see below). Future work will include kinetics, including attempting to model the ionization (laser on) phase, but it is important to note that several past studies have found that the observed gas dynamics effects are dominated by the interaction between the blast wave and the hot gas kernel after the pulse (i.e. not depending on kinetics) [13,16,24,25]. Moreover, post-discharge heating by vibrationaltranslational (V-T) relaxation reactions or by quenching of metastable states take place on a longer timescale (hundred of microseconds) compared to the timescale of blast-wave detachment and kernel collapse (~0.5-10 µs). The blast itself is created by the ultra-fast heating mechanism of electron avalanche ionization which takes place during the laser pulse. This mechanism is taken into account in the current simulations through initial conditions based on the plasma properties measured using Rayleigh/Thomson scattering at the end of the pulse. Hence, neglecting the plasma chemistry should still allow modeling and insight into how the flow field develops provided that correct profile of pressure and temperature at the end of the pulse are used in computations.
The governing equations written in cylindrical coordinates with axial symmetry (zero azimuthal velocity) are presented below in conservative form: (1) where: (2) In Eq. (2), U is the vector of conservative variables, F(U) a and G(U) a describe the advection flux (in the r and z direction respectively), F(U) d , G(U) d represent the diffusion flux, and S(U) is the source term that is the result of the coordinate system transformation from Cartesian to cylindrical. is the gas density, p is pressure, u and v are the radial (r-) and axial (z-) components of velocity, E is the total energy (the sum of kinetic and internal energy), , are the shear stress component, and qi represents the i th component of the conductive heat flux. In this coordinate system the incident laser is oriented along the z-axis. + ( ) + ( ) = ( ) + ( ) + ( ) 4 The system of conservation laws is complemented with two additional relations that describe the thermodynamic state of the system: = where the first relation in Eq. (3) is the ideal gas law (Rg is the specific gas constant), and the second is the caloric equation of state (h is the specific enthalpy). The caloric perfect gas assumption does not fully apply in the case of the high-temperature laser plasma; however, the numerical model assumes a thermally perfect gas which involves solving an implicit nonlinear equation in temperature in each computational cell and at each time step. The following relation for temperature is derived from Eq. (3): In Eq. (4) the enthalpy, h(T), is computed using the Shomate interpolation polynomials for molecular nitrogen and oxygen according to the NIST database [26]. Viscosity is modeled using Sutherland's law, assuming a Newtonian fluid with temperature-dependent viscosity, and the conductive fluxes are calculated with temperature-dependent thermal conductivities [27].
The conservation equations given in Eq. (1) are numerically integrated using a dimensional splitting approach. The advection terms are solved using Roe's flux differencing scheme with slope limiters [28] and the diffusion terms are solved using a second order centered difference scheme [29]. Finally, the resulting ODEs are integrated with an explicit 4 th order Runge-Kutta IV method.
The governing equations are solved on a structured Cartesian grid using a control volume discretization approach. A piece-wise distribution of initial data over the computational cells is employed where the value at the center of a cell represents the average of each variable over the cell. An optimum cell resolution is required in order to resolve the phenomenon of interest while minimizing the computational requirements. Experimental results indicate that kernel evolution takes place over two timescales suggesting that two different cell sizes can be used in the simulation. In the first 20 the kernel grows significantly; the temperature drops by more than half of its initial value and vorticity is generated via baroclinic torque (see Section III-D for a detailed description). Over this time interval, an optimum cell size of 15 microns was determined through a grid sensitivity analysis. At times greater than 20 , the cell size is doubled to 30 microns in order to save computational time since the gradients of pressure and temperature are much weaker and the toroidal structure has already formed. Additionally, the axial-symmetric nature of the kernel under investigation allows for modeling only half of the kernel which further reduce the computational cost. The code uses non-reflective boundary conditions at the edge of the domain. This allows the blast wave to leave the computational domain with minimal disturbances [30,31].

A. Determination of Initial Conditions
Since the current simulations do not model the optical breakdown process, the choice of the initial conditions (corresponding to the end of the laser pulse) becomes critical. These initial conditions are found from experimental studies based on breakdown from a focused Nd:YAG beam (1064 nm) with representative parameters of 10 ns pulse duration, 70 mJ energy and waist of 170 m. The initial kernel shape is based on broadband optical emission images taken immediately after the pulse has ended [32]. The pressure and temperature of the kernel at the end of the discharge are determined from 1-D Rayleigh/Thomson scattering measurements [33]. A typical scattering image is show in the left of Figure 1 [33][34][35]. The right of Fig. 2 shows the resulting 1-D gas density profiles (normalized by ambient density, 0) along the Rayleigh probe beam for different times after the end of the breakdown pulse. At all times, far from the plasma (where the gas is undisturbed) the density remains unchanged, i.e. /0=1. Moving inwards (from either side) toward the center of the domain there is a sharp increase in density due to the shock wave at the front of the blast that accompanies the breakdown process. At regions inside the blast, the density steadily drops due to progressively higher temperatures. (At early times, an increase in signal is observed right at the center of the kernel due to the presence of Thomson scattering from the free electrons; this signal contribution does not correspond to the true gas density and is excluded from determination of initial conditions.) These density profiles (right of Fig. 2) contain information on the shock wave location, blast wave strength and kernel size, at various times after the breakdown, and are thus well suited for determination of initial conditions for the simulation. The method used to determine the initial conditions (corresponding to the end of the breakdown pulsed) is as follows. We assume a teardrop shape (major axis along the laser propagation direction) for the starting kernel from our past imaging measurements [32]. We then iteratively run the CFD model for various values of peak temperature and kernel size and compare the resulting density profiles to their experimental counterparts (right of Fig.1), seeking the best agreement. The pressure profile is linked to the temperature profile using the ideal gas equation of state under the assumption that the gas density should remain unchanged during the laser energy deposition phase. This hypothesis is valid as long as the laser pulse duration is much shorter than the acoustic-time scale (where the latter is given by the radial length of kernel, divided by the speed of sound). The simulated density profiles are shown with solid lines in the right of Fig. 2. In terms of selecting the needed initial conditions (for end of laser pulse) we set a criterion that both the location and the amplitude of the peak of the blast wave should agree, between simulated and Rayleigh profiles, to better than 5%. If one considers how much the (initial condition) pressure and temperature fields can vary, while still satisfying the 5% agreement criterion, one finds that the peak temperature varies be less than ~3% (i.e. ~500 K maximum change in peak temperature of ~15,000 K). Such a small temperature change will not change the hydrodynamic effects presented here. The initial temperature and pressure profiles are shown in Figure 3.  An important contribution of this study is the method, outlined above, of using Rayleigh scattering profiles to infer the pressure and temperature of the breakdown plasma at the end of the laser pulse. Previous work in the field has typically used Taylor-Sedov blast wave theory which relies on matching the experimentally observed blast wave radii (usually from schlieren) with those predicted from the theory. However, we find that that this method grossly over-estimates the kernel's thermodynamic state at the end of the pulse. Taylor's blast wave analysis assumes a perfectly spherical shock wave front. However, this is not true for a typical laser spark that more closely resembles a cylinder than a sphere. In this case, the blast wave speed varies strongly with location along the circumference thus confounding analysis based on a simpler spherical wave. Additionally, Taylor-Sedov analysis assumes an instantaneous explosion from an infinitesimally small volume. This means that, while the Taylor-Sedov theory can be used to determine the energy released based on measurements of the shock wave radii once it becomes spherical, one cannot simply extrapolate temperature and pressure immediately after the end of the nanosecond pulse. Our analysis shows that the size of the initial kernel plays a key role in determining the correct temperature and pressure immediately after the pulse. To illustrate this point, Fig. 4 shows simulated cases for two different starting kernel configurations that result in similar blast wave radius and strength variations over time. We see that when considering both size and temperature, one can get nearly degenerate blast wave evolution for different parameter combinations, indicating that the size also plays a significant role in determining the correct thermodynamic state at the end of the laser pulse. If this size effect is ignored the kernel temperature is grossly over-estimated, with previous literature using Taylor-Sedov analysis reporting temperatures in excess of 100,000 K and pressures between 5,000-10,000 bar [15,36,37].

B. Kernel Dynamics
A comparison between experimental and synthetic schlieren images showing the kernel development after the breakdown pulse is presented in Fig. 5. The synthetic images were obtained by plotting the contours of the vertical density gradient (which would correspond to a horizontal knife edge orientation in experiments). The images show that the kernel maintains its initial teardrop shape for the first ~5 µs after the pulse. However, its volume has increased by approximately a factor of 4 during this initial stage. This is followed by compression of the kernel along the optical propagation axis between ~5-30 µs. At 30 µs after the pulse, the kernel presents two strong indentations at the leading and trailing edge. As will be shown in detail in the following section, this is caused by two significant recirculation regions that entrains the cold surrounding air from either side and push this air towards the center of the kernel. At 50 µs, the three-lobe structure (i.e., the two side lobes and the frontal lobe) is already formed. In 3-D, this structure resembles an asymmetric torus in which the plane of stagnation at the center is broken because the re-circulation at the back of the kernel is stronger. In the final stages of the kernel development, the frontal lobe continues to grow towards the laser source as more surrounding fluid is pushed through the middle of the torus. The kernel maintains this shape for several hundred of microseconds until it is completely dissipated by Rayleigh-Taylor instabilities formed inside the flow. The actual lifetime of the three-lobe structure depends on the energy absorbed inside the gas and the focusing conditions, but it is typically between 0.1-1 ms. A comparison between the evolution of maximum kernel temperature obtained experimentally by Rayleigh scattering and CFD is shown in in Fig. 6. Interestingly, the temperature decay presents three successive regions: (a) a region of sharp temperature decrease associated with fast kernel expansion at times <~1 s, (b) a region where the rate of temperature decrease slows between ~1-15 µs, and (c) a region governed by toroidal dynamics in which gradual kernel cooling is observed between 15-200 µs after the pulse. The relatively uniform temperature in region (b) is a new finding and is due to heat diffusion from the kernel being balanced by heating of the kernel center caused by the axial compression that was previously described. Temperature measurements at earlier times are not shown as Thomson scattering and pressure changes precludes temperature determination by the current Rayleigh scattering method. The comparison shows that the CFD simulation captures the temperature evolution during the toroidal development stage.

C. Blast Wave Propagation
Previous work on laser breakdown has outlined the importance of blast wave dynamics in the generation of vorticity. Blast waves are observed in many types of ultrafast heat releases from point sources. A typical blast wave front consists of a leading shock wave that compresses and heats the gas ahead of it and a trailing rarefaction wave which acts to reduce the gas density behind the wave returning it to its initial (unperturbed) state. During the interaction between the heated gas and the rarefaction, a zone of pressure gradient reversal is observed when the pressure behind the blast becomes smaller than the ambient pressure. In bomb detonations, this is typically referred as a "blast wind" which tends to suck destroyed items back towards the epicenter of the explosion. Figures 7 and 8 show the computed pressure fields for times before and after blast wave detachment respectively. Since the laser energy is deposited inside the gas on a much shorter timescale (nanoseconds) compared to the hydrodynamic scale of the laser plasma (order microseconds based on speed of sound and characteristic spark dimension), a pressure buildup inside of the plasma kernel accompanies the rise in temperature during the energy addition process. In the early stages of the kernel development, acoustic waves coalesce at the edge of the kernel thereby forming a strong blast wave. This process is accompanied by a fast decrease in pressure at the center of the kernel during the first few hundred nanoseconds. Interestingly, an increase in pressure at the center of the kernel is observed at ~300 ns after the pulse, corresponding to the formation of a secondary pressure front (see below). Despite being weaker in amplitude, this front propagates faster than that of the main blast since the density of the fluid around it is much lower than that of the far-field; this secondary pressure front then catches up (and merges) with the main blast wave at around ~1.5 µs after the pulse. The resulting blast initially takes the shape of the post-laser breakdown teardrop kernel (Fig. 3). However, as the time evolves the blast wave front becomes progressively more spherical. This phenomenon is described in more detail below and has a significant impact on vorticity generation. Once the blast wave detaches from the kernel (at time ~1 s), the pressure at the epicenter gradually decreases. At 1 µs after the pulse, the rarefaction wave that closely follows behind the shock front begins to decrease the pressure inside the gas kernel until it becomes smaller than the ambient pressure. The outward pointing pressure gradient observed in the 2.5 µs and 5 µs frames at the center of the kernel gives rise to a fluid flow directed towards the center of the kernel. This inward flow initiates the observed kernel collapse into the asymmetric torus shown in Fig. 6. Finally, the pressure at the center of the kernel recovers to atmospheric in ~10-15 µs after the pulse and the blast wave becomes fully detached and perfectly spherical. A more detailed analysis of the secondary pressure wave observed between 0.3-1.5 µs is warranted as this feature has not been previously observed to our knowledge. The wave velocity is subsonic (e.g. ~600 m/s corresponding to M=0.3 at 0.3 µs after the pulse) meaning the wave is not a secondary blast. Being subsonic, this secondary front does not appreciably compress the fluid which is why it is not visible in either schlieren or Rayleigh scattering measurements. The source of this secondary front is revealed in Fig. 9. Inspecting the velocity profiles at 0.4 µs after the pulse reveals that significant recirculation takes place inside of the kernel. Since the pressure at the epicenter is still above atmospheric there is a net acceleration towards the outside. However, the strong source of vorticity near the kernel boundary tends to recirculate the fluid back to the kernel center leading to a temporary increase in gas pressure. This is a significant finding because previous research attribute the presence of vorticity within the kernel to (only) the interaction between the rarefaction and the remaining hot gas kernel behind the blast [13], yet the present results suggest that vorticity is generated within the kernel even during the (earlier) blast wave formation process. This result has been checked for grid resolutions as small as 1 µm and it does not appear to be a numerical artefact. A comprehensive discussion of vorticity generation is given in the next section. Right: Velocity field imposed on velocity contours showing the presence of vorticity at early times during the kernel development (t=0.4 µs is shown here). The region used for the right panel is shown with the red rectangle in the left panel, i.e. the right plot is made for half of the kernel.

D. Mechanism of Vorticity Generation
Fig. 10 provides vorticity contours for various time delays during the kernel development. Positive vorticity means that the fluid is rotating counter-clockwise in the plane of the figure (color coded red in the images below) while negative vorticity indicates clockwise rotation (shown in blue). As previously noted, vorticity is generated during the first few hundred of nanoseconds after the laser pulse, well before the leading blast wave has detached itself from the main kernel. The two counter rotating vortices push the fluid radially inwards, thus causing the formation of the secondary pressure front observed between ~0.3-1.5 µs. This source of vorticity persists until the outward pressure gradient is established inside the kernel (due to formation of the rarefaction wave at the kernel boundary). As the kernel collapses, the two vortices change their direction of rotation and entrain the surrounding cold fluid and push it along the optical axis towards the center of the kernel (between ~2 and 100 s). This causes the kernel to compress faster over the optical axis than radially (compare the images for 4 and 7 s). At 7 µs, one notices indentations at the top and bottom of the kernel which grow as time progresses. The vortices gradually move towards the center of the kernel during the collapse, until they collide with each other a few tens of microseconds after the pulse (note that the rotations of the vortices are opposite in the images for 7 s versus 100 s in Fig. 10), and ultimately pass one another. When the vortices collide, the stronger vortex (originally at the bottom) is able to displace the weaker vortex (originally at the top) radially outwards and continue its upward motion. Finally, the stronger vortex, positioned at the top, continues to entrain surrounding fluid from the bottom which is expelled at the top giving rise to the frontal lobe which propagates towards the laser source. The kernel boundary (shown with dashed black lines in Fig. 10 and Fig. 11) was determined using an edge finding algorithm that tracks the location of maximum temperature gradient along the perimeter of the kernel. To better understand the mechanism(s) through which vorticity is generated, each of the terms in the vorticity equation are evaluated individually in CFD post-processing: The first term on the right-hand side describes enhancement of vorticity by stretching and tilting due to the presence of velocity gradients. This is the mechanism through which turbulent eddies transfer energy to smaller scales in turbulent flows. The second term describes the effects of flow compressibility on the vorticity generation. This term enhances vorticity through the flow-field as the flow undergoes expansion or contraction. The third term describes the effects of viscous diffusion on vorticity generation. Finally, the last term describes vorticity generation by baroclinic torque, i.e. misalignment between density and pressure gradients. Evaluating each term shows that baroclinic torque is the main mechanism through which vorticity is generated with some contribution also from the compressibility (contraction) term. Fig. 11 shows 2-D contours of the rate of vorticity from the dominant terms as well as the total rate. The images show that vorticity is generated predominantly by baroclinic torque and that the expansion term only acts as a source of enhancement/dissipation of vorticity that has already been generated in the fluid. Most of vorticity is generated during the first ~5 µs after the laser pulse with the baroclinic torque accounting for ~90% of the total vorticity generated at this stage. At ~10 µs after the pulse, the blast wave becomes fully detached and perfectly spherical thus removing the main source of baroclinicity. As time progresses, the contribution of the expansion term to the total vorticity generation increases. The kernel compression (which was initially induced by baroclinic torque) now gradually becomes the main contributor to vorticity enhancement. However, the rate of vorticity generation has dropped by two orders of magnitude: from ~10 12 1/s 2 at 4 µs down to ~10 10 1/s 2 at 10 µs. A comparison between the pressure contours presented in the previous section and Fig. 11 shows that the generation of baroclinic torque is related to the blast wave dynamics. Baroclinic torque is generated in two stages. Before the shock wave detaches from the kernel surface, baroclinic torque induces vorticity which leads to the generation of the secondary pressure front discussed earlier. Then, after the formation of the rarefaction wave at the tail of the blast, the direction of baroclinic torque is reversed. As show in Fig. 11, baroclinic torque is the main contributor to the axial collapse of the kernel during this stage. After ~10 µs the effect of baroclinic torque becomes negligible (note that the rate of vorticity generation drops by two orders of magnitude between 1-10 µs). This time delay also corresponds to the complete blast wave separation from the kernel. Figure 11: Rate of vorticity generation by baroclinic torque (top) and kernel contraction (middle) compared to the total rate of vorticity generation (bottom). The shape of the actual spark kernel is shown with dotted line contour. The images show how vorticity is generated during the kernel collapse phase (t=0.1-10 µs).
A physical explanation of the mechanism that generates baroclinic torque is provided with the aid of Fig. 12. As the shock wave is formed at the boundary of the kernel during the early stages of kernel development, it has a nonuniform pressure around its circumference. The pressure contours in the top of Fig. 12 show that the blast is stronger on the sides of the kernel where the front is more planar compared to the top and bottom where the front is more spherical. This is to be expected because the energy is dissipated more quickly in a spherical shock versus a planar one due to stretching of the wavefront during the associated increase in perimeter. A consequence of the non-uniform pressure field is that the pressure gradient at the "corner" of the kernel is slightly tilted towards the higher pressure part (side region) of the shock front. In contrast, the gradient of density simply follows the temperature gradient (which points directly towards the epicenter of the blast) but it is inverted since higher temperature leads to lower density in a region of constant pressure. The broken line in the central panel of Fig. 12 indicates the orientation of the density gradient (as determined from the left panel). The two gradients and the ensuing vorticity direction are shown together schematically in the right panel of Fig. 12. This misalignment between the density and pressure gradients at the corners of the kernel leads to the generation of vorticity by baroclinic torque. At t~1.5-2 µs after breakdown, the detachment of the blast wave from kernel establishes an inward pointing pressure gradient which is caused by the inward flow induced by the rarefaction wave. This leads to a decrease, and ultimately to a direction reversal, of the vorticity generated by baroclinic torque. The bottom part of the blast wave has an even tighter radius of curvature than the top leading to similar behavior but with larger gradient misalignment, thus making the vortex at the bottom stronger than the one at the top. As discussed previously, this will ultimately dictate the direction of propagation of the third lobe, i.e., towards the direction of the plasma forming laser for the present case where the energy deposition leads to a stronger vortex on the downstream (bottom) side. Many of these ideas can be conceptually understood by considering the blast wave dynamics of a kernel initially shaped as a long thin cylinder: In this case the blast comprises two spherical fronts from the top and bottom of the cylinder and a planar front from the sides; the planar side front will propagates faster eventually catching up to the spherical fronts and ultimately become spherical (with baroclinic torque again generated at the "corners" until the wave is fully spherical). In this scenario, a perfectly symmetric torus will be obtained (with no third lobe observed) since the vortex pair will have equal strength. Such kernel dynamics are sometimes encountered in nanosecond pin-to-pin discharges when the energy is distributed homogenously inside the electrode gap [24,38,39]. Moreover, while the present focus is on NIR laser sparks, the mechanism described here can also be used to explain the hydrodynamics observed in visible and ultra-violet laser induced plasma.

IV. Conclusions and Future Work
A numerical study of the gas dynamics induced by a laser spark was presented in this manuscript. It was found that the blast wave dynamics play a key role in the development of the experimentally observed asymmetric torus (three-lobe structure). Initially, the kernel takes the shape of a teardrop with the hotter region facing towards the laser propagation direction. This will ultimately influence how vorticity is being generated inside the kernel during its subsequent evolution. The strength of the ensuing blast wave, which initially takes the shape of the kernel, varies along its circumference during the early stages of kernel development. This leads to misalignment between the density and pressure gradient and generates vorticity by baroclinic torque. An analysis of the vorticity equation shows that baroclinic torque is the main mechanism that governs the collapse of the initial teardrop kernel into an asymmetric torus. Baroclinic torque is formed in two stages: first, before the shock wave detaches from the boundary of the kernel, baroclinic torque generates vorticity that brings fluid radially inwards and leads to the formation of a secondary pressure front within the kernel; second, several microseconds after the pulse, the rarefaction wave generates a pressure gradient reversal inside the kernel. This changes the direction of rotation induced by the baroclinic torque and leads to the axial collapse of the kernel. Past explanations of vorticity generation in the laser induced plasmas, in particular from Bradley [13], considered vorticity formation due to the (later) interaction of the rarefaction wave and kernel but not the possibility of (earlier) vorticity formation before detachment of the shock wave.
Ultimately, what dictates the direction in which the third lobe propagates is which of the two counter-rotating vortices is stronger. This depends primarily on the initial kernel topology. If the kernel is hotter at the front (which is the typical configuration for the optical breakdown process in which EAI is the dominant ionization mechanism), then the vortex at the bottom will be stronger (due to larger misalignment between the pressure and the density gradients from the tighter radius) causing entrained gas to expel a lobe towards the front side (laser incident side). From a practical point of view, the following finding is important: if one desires to control the kernel dynamics in applications such as ignition, flame stabilization or aerodynamic flow control, one needs to control how the energy is deposited around the beam waist. This is consistent with our past work using dual-pulse laser plasma showing that the lobe can present at the front or back sides of the kernel, or not at all, depending on energy deposition [23].