1. Introduction
There is need to understand the laminar mixing and combustion that commonly occur within turbulent eddies. These laminar flamelet sub-domains experience significant strain. Some important work has been done here but typically for counterflows or simple vortex structures in two dimensions or with axisymmetry and often with a constant-density approximation. See Liñán (Reference Liñán1974), Marble (Reference Marble and Casci1985), Karagozian & Marble (Reference Karagozian and Marble1986), Cetegen & Sirignano (Reference Cetegen and Sirignano1988, Reference Cetegen and Sirignano1990), Peters (Reference Peters2000) and Pierce & Moin (Reference Pierce and Moin2004). Linan and Peters focused on the counterflow configuration. Liñán & Williams (Reference Liñán and Williams1993) examined ignition for non-premixed, unsteady, two-dimensional counterflow with variable density. Karagozian and Marble examined a three-dimensional flow with radial inward velocity, axial jetting and a vortex centred on the axis. The flame sheet wrapped around the axis due to the vorticity. Pierce and Moin modified the counterflow configuration by fixing domain size and forcing flux to zero at the boundaries.
These models are built around the postulate that the flamelets are always non-premixed (i.e. diffusion) flames and subject to flow strain. Nguyen, Popov & Sirignano (Reference Nguyen, Popov and Sirignano2018) and Nguyen & Sirignano (Reference Nguyen and Sirignano2018) employed the Pierce–Moin flamelet approach in the simulation of a single-injector rocket engine. They showed the importance of flamelets subject to high strain rates. However, contradictions occurred in that both premixed flames and non-premixed flames appeared in the predictions. In fact, they report multi-branched flames; in particular, the combination is often seen of a fuel-lean premixed-flame branch with a branch consisting of a merged diffusion flame and fuel-rich premixed flame.
Experiments and asymptotic analyses (Hamins, Thridandam & Seshadri Reference Hamins, Thridandam and Seshadri1985) showed that a partially premixed fuel-lean flame and a diffusion flame may co-exist in a counterflow with opposing streams of heptane vapour and methane–oxygen–nitrogen mixture. Thus, a need exists for flamelet theory to address both premixed and non-premixed flames. Recently, Rajamanickam et al. (Reference Rajamanickam, Coenen, Sanchez and Williams2019) has provided an interesting three-dimensional triple-flame analysis. Sirignano (Reference Sirignano2019b) has provided a counterflow analysis with three-dimensional strain and shown the possibility for a variety of flame configurations to exist depending on the compositions of the inflowing streams: (i) three flames including fuel-lean partially premixed, non-premixed (i.e. diffusion controlled) and fuel-rich partially premixed; (ii) non-premixed and fuel-rich partially premixed; (iii) fuel-lean partially premixed and non-premixed; (iv) non-premixed; and (v) premixed. López-Cámara, Jordà Juanós & Sirignano (Reference López-Cámara, Jordà Juanós and Sirignano2019) has extended the counterflow analysis to consider detailed kinetics for methane–oxygen detailed chemical kinetics.
There is a strong need to study mixing and combustion in three-dimensional flows with both imposed normal strain and shear strain (and therein imposed vorticity with global circulation). It is well known that the vorticity vector will tend to align with the direction of tensile (i.e. extensional) normal strain in a flow, which leads us to choose a certain three-dimensional configuration that combines the mixing layer and the counterflow. We also expect a material interface to align to be normal to the direction of the compressive normal strain. See Nomura & Elghobashi (Reference Nomura and Elghobashi1992), Nomura & Elghobashi (Reference Nomura and Elghobashi1993), Boratav, Elghobashi & Zhong (Reference Boratav, Elghobashi and Zhong1996) and Boratav, Elghobashi & Zhong (Reference Boratav, Elghobashi and Zhong1998). In this work, we extend flamelet theory in two significant aspects: the inclusion of both premixed and non-premixed-flame structures and the extension to three-dimensional fields with both shear and normal strains. We address a steady three-dimensional mixing-layer flow with primary flow component 
$u$ in the 
$x$-direction and an imposed counterflow with 
$v$ and 
$w$ velocity components in the 
$y$ and 
$z$ directions, respectively. The particular flow configuration considered here is sketched in figure 1. The monotonic profile of 
$u(\kern-0.4pt y)$ at fixed 
$x$ and 
$z$ is shown with the features of a traditional mixing layer. It also shows the imposed compressive normal strain is in the 
$y$-direction and the commensurate expansive normal strain is in the 
$z$-direction. So, convergence of streamline projections occurs in the 
$x$–
$y$ plane with divergence of streamline projections in the other two planes.
Figure 1. Sketch of mixing-layer flow with imposed counterflow.
Combustion, variable density and variable properties are examined. The classical counterflow treatment by Peters (Reference Peters2000) has two opposing streams, one of fuel or fuel plus a chemically inert gas and the other of oxidizer or oxidizer plus an inert gas. Our computations address that situation where a single diffusion flame exists. We also provide here some background analytical considerations for situations where the inflowing streams from 
$y(\infty )$ and 
$y(-\infty )$ (each of which also now has a parallel component of forced velocity in the 
$x$-direction) may consist of either only one reactant or a combustible mixture of fuel and oxidizer, thereby allowing another flame besides the simple diffusion flame to co-exist. Propane and oxygen are specifically considered with one-step, Westbrook & Dryer (Reference Westbrook and Dryer1984) kinetics; however, the qualitative conclusions are expected to be more general. Those one-step kinetic relations were obtained by fitting to experiments for premixed flames and are expected to be less accurate for diffusion flames. Nevertheless, we may accept some error here because diffusion rather than kinetics is rate controlling for diffusion flames. The approach here expands on recent work (Sirignano Reference Sirignano2019a) which used infinite kinetics for three-dimensional counterflow diffusion flames and another work (Sirignano Reference Sirignano2019b) which used one-step kinetics.
In § 2, the analysis is presented. In sequence, the three-dimensional problem is reduced to a two-dimensional form and then, for the downstream mixing-layer flow, to a one-dimensional similar form. The system of ordinary differential equations (ODEs) is presented for the thermo-chemical variables and the velocity components. An analysis is given in § 2.5 to relate a conserved scalar to the velocity solution. The validity of the similar solution form for mixing layers with certain thin reaction zones is discussed using concepts from singular perturbation theory. The chemical kinetic model is described in a form to be used as a source term for an ODE. In § 2.8, the relevance of these new findings for flamelet theory is discussed. Then, in § 3, the findings from calculations for the mixing layer with imposed counterflow are presented for both the non-reacting flow and the flow with the diffusion flame. Conclusions are presented in § 4.
2. Analysis
The planar mixing layer has been widely used in combustion studies for development of flamelet models. In particular, a similar solution can be produced, offering the convenience of reduction to a system of ODEs to describe a multi-dimensional configuration. Here, we develop a similar solution for a specific three-dimensional configuration. Consider a mixing layer with primary flow in the 
$x$-direction and diffusion primarily in the 
$y$-direction. An oxidizer rich gas enters the mixing-layer domain at negative 
$y$ values while a fuel-rich gas enters at positive 
$y$ values. Shear strain and associated vorticity with a vector in the 
$z$-direction are obviously introduced here. Superimposed on this otherwise-planar flow is a counterflow with an incoming stream with compressive (i.e. negative) normal strain in the 
$y$-direction and an outgoing stream with expansive (i.e. positive) normal strain in the 
$z$-direction. The velocity 
$\boldsymbol {u}$ has the components 
$u, v$ and 
$w$ in the 
$x,y$ and 
$z$ directions, respectively. The interface of the two incoming, opposing transverse flows is specified at 
$y = 0$. The chosen configuration here has no 
$x$-component to the pressure gradient. The classical boundary-layer approximation is made, rendering the 
$y$-component of the pressure gradient to be negligible compared to other terms. If the approaching streams have the same pressure at a distance from the interface and its viscous layer, we expect that, in a frame of reference attached to the interface, momentum balance for steady flow yields 
$\rho _{-\infty } v_{-\infty }^2 = \rho _{\infty } v_{\infty }^2$ where 
$\rho$ is fluid density. The 
$y$-directed inflowing streams in all cases bring together fluids of differing temperature 
$T$ and/or mass fractions 
$Y_m$ where each value of the index 
$m$ indicates a specific species that forms the composition of the fluid. The two streams will generally have different upstream values for 
$v, \rho , T, Y_m$ and enthalpy 
$h$. Pressure 
$p$ will be given the same values for the two inflowing streams.
The study here can be viewed as an extension of two-dimensional mixing and combustion in a shear layer to include counterflow in the transverse plane. Also, it may be viewed as an extension of the counterflow study of Liñán & Williams (Reference Liñán and Williams1993) to include imposed shear flow and vorticity.
A Newtonian fluid with the Stokes hypothesis is examined. Radiation and gravity are neglected. Fickian mass diffusion and Fourier heat conduction are considered so that all fluid properties are continuous across the interface. In the reactive case, kinetic energy and viscous dissipation will be neglected in the energy consideration because of a focus on low-Mach-number flow. However, those terms will be kept in the reactive case. The properties 
$\mu , \lambda , D, c_p$ and 
$h_{f,m}$ are the dynamic viscosity, thermal conductivity, mass diffusivity, specific heat at constant pressure and heat of formation of species 
$m$, respectively; 
$\omega _m$ is the chemical reaction rate of species 
$m$ (in reciprocal time units) for the oxidation process. The Prandtl number (
$Pr \equiv \mu c_p/\lambda$) and the Schmidt number (
$Sc \equiv \mu / (\rho D)$) will be assumed to have the same constant value. Thus, the Lewis number (
$Le = Sc/Pr$) has unitary value.
2.1. Three-dimensional formulation
The governing equations for steady three-dimensional flow with the boundary-layer approximation are given as
\begin{gather} \frac{\partial (\rho u)}{\partial x} + \frac{\partial (\rho v)}{\partial y} + \frac{\partial (\rho w)}{\partial z} = 0, \end{gather}
\begin{gather} \rho u\frac{\partial u}{\partial x} + \rho v\frac{\partial u}{\partial y} + \rho w\frac{\partial u}{\partial z} = \frac{\partial}{\partial y} \left(\mu\frac{\partial u}{\partial y}\right), \end{gather}
\begin{gather} \rho u\frac{\partial w}{\partial x} + \rho v\frac{\partial w}{\partial y} + \rho w\frac{\partial w}{\partial z}+\frac{\partial p}{\partial z} = \frac{\partial}{\partial y} \left(\mu\frac{\partial w}{\partial y}\right), \end{gather}
\begin{gather} \rho u\frac{\partial h}{\partial x} + \rho v\frac{\partial h}{\partial y} + \rho w\frac{\partial h}{\partial z} = \frac{1}{Pr}\frac{\partial}{\partial y} \left( \mu \frac{\partial h}{\partial y} \right) -\rho \varSigma_{m=1}^N h_{f,m} \omega_m, \end{gather}
\begin{gather} \rho u\frac{\partial Y_m}{\partial x} + \rho v\frac{\partial Y_m}{\partial y} + \rho w\frac{\partial Y_m}{\partial z} = \frac{1}{Pr}\frac{\partial}{\partial y} \left( \mu \frac{\partial Y_m}{\partial y} \right) + \rho \omega_m; \quad m=1, 2,\ldots, N. \end{gather}
The equation for the 
$y$-component of momentum has been replaced by the classical boundary-layer approximation that no variation of pressure in the 
$y$-direction through the mixing layer exists. Note that pressure will have a gradient in the 
$y$-direction outside of the boundary layer to create the counterflow. However, the gradient must change direction in the mixing layer passing through the zero value. Thus, the classical zero-gradient assumption using the boundary-layer approximation for a thin mixing layer is valid.
The boundary conditions have the velocity component 
$u$, the enthalpy 
$h$ and the mass fractions 
$Y_m$ going to prescribed constants as 
$y \rightarrow \infty$ and as 
$y \rightarrow - \infty$. In those same limits, 
$v$ and 
$w$ go to a counterflow configuration described as a potential flow in the 
$y$–
$z$ plane where compressive normal strain occurs in the 
$y$-direction with tensile normal strain in the 
$z$-direction. Thus, as 
$y \rightarrow \infty$
\begin{equation} \left. \begin{gathered} u \rightarrow u_{\infty}; \quad v \rightarrow - \kappa_{\infty}y;\quad w \rightarrow \kappa_{\infty}z, \\ h \rightarrow h_{\infty}; \quad Y_m \rightarrow Y_{m, \infty};\quad m=1, 2,\ldots, N , \end{gathered} \right\} \end{equation}
and as 
$y \rightarrow -\infty$
\begin{equation} \left. \begin{gathered} u \rightarrow u_{-\infty};\quad v \rightarrow - \kappa_{-\infty}y;\quad w \rightarrow \kappa_{-\infty}z, \\ h \rightarrow h_{-\infty}; \quad Y_m \rightarrow Y_{m, -\infty};\quad m=1, 2,\ldots, N . \end{gathered} \right\} \end{equation}2.2. Reduction to two dimensions and conserved scalars
Following Rajamanickam et al. (Reference Rajamanickam, Coenen, Sanchez and Williams2019) , we consider 
$w = \kappa z$, differing in that we allow 
$\kappa$ to be a function of 
$x$ and 
$y$ rather than constant. Accordingly, 
$p$ will also vary with 
$z$. In particular, since 
$p$ will not vary with 
$x$ or 
$y$, we have 
$p=p_o - (1/2)\rho _{\infty }w_{\infty }^2 = p_o - 1/2)\rho _{\infty }\kappa _{\infty }^2z^2$ where 
$p_o$ is constant. All other variables (
$u, v, \rho , h, Y_m$) will depend only on 
$x$ and 
$y$. The resulting two-dimensional system of equations follows:
\begin{gather} \frac{\partial (\rho u)}{\partial x} + \frac{\partial (\rho v)}{\partial y} + \rho \kappa = 0, \end{gather}
\begin{gather} \rho u\frac{\partial u}{\partial x} + \rho v\frac{\partial u}{\partial y} = \frac{\partial}{\partial y} \left(\mu\frac{\partial u}{\partial y}\right), \end{gather}
\begin{gather} \rho u\frac{\partial \kappa}{\partial x} + \rho v\frac{\partial \kappa}{\partial y} + \rho \kappa^2 - \rho_{\infty} \kappa_{\infty}^2 = \frac{\partial}{\partial y} \left(\mu\frac{\partial \kappa}{\partial y}\right). \end{gather}
Since the pressure does not vary with 
$y$, we must have 
$\rho _{\infty } \kappa _{\infty }^2 = \rho _{-\infty } \kappa _{-\infty }^2$. Thus, the derivates of 
$\kappa$ will go to zero at plus and minus infinity for the 
$y$ value:
\begin{gather} \rho u\frac{\partial h}{\partial x} + \rho v\frac{\partial h}{\partial y} = \frac{1}{Pr}\frac{\partial}{\partial y} \left( \mu \frac{\partial h}{\partial y} \right) -\rho \varSigma^N_{m=1}h_{f,m} \omega_m, \end{gather}
\begin{gather} \rho u\frac{\partial Y_m}{\partial x} + \rho v\frac{\partial Y_m}{\partial y} = \frac{1}{Pr}\frac{\partial}{\partial y} \left( \mu \frac{\partial Y_m}{\partial y} \right) + \rho \omega_m;\quad m=1, 2,\ldots, N. \end{gather}
Here, the sensible enthalpy 
$h =c_pT$ is based on the assumption of a calorically perfect gas. When normalized by ambient conditions, the non-dimensional values of enthalpy and temperature are identical. For simplification, we neglect the effect of species composition on specific heats and the specific gas constant. For the one-step kinetics considered here, the last term in (2.11) may be replaced by 
$- \rho Q \omega _F$ where 
$Q$ is the heating value (energy/mass) of the fuel and 
$\omega _F < 0$ is the chemical oxidation rate of the fuel. Then, 
$Y_F$ and 
$Y_O$ are the mass fractions of propane and oxygen, respectively. 
$\nu =0.275$ is the stoichiometric ratio of propane mass to oxygen mass.
The boundary conditions on (2.9)–(2.12) involve specifications of the dependent variables at 
$y=\infty$ and 
$y= - \infty$ as well as their values at an upstream value of the coordinate 
$x$. From these primitive equations for 
$h$ and 
$Y_m$, we may form equations for conserved scalars following concepts introduced by Shvab (Reference Shvab1948) and Zel'dovich (Reference Zel'dovich1949). We define two conserved scalars as
\begin{equation} \left. \begin{gathered} \alpha \equiv Y_F - \nu Y_O, \\ \beta \equiv h + \nu Y_O Q, \end{gathered} \right\} \end{equation}to obtain
\begin{equation} \left. \begin{gathered} \rho u\frac{\partial \alpha}{\partial x} + \rho v\frac{\partial \alpha}{\partial y} = \frac{1}{Pr}\frac{\partial}{\partial y} \left( \mu \frac{\partial \alpha}{\partial y} \right),\\ \rho u\frac{\partial \beta}{\partial x} + \rho v\frac{\partial \beta}{\partial y} = \frac{1}{Pr}\frac{\partial}{\partial y} \left( \mu \frac{\partial \beta}{\partial y} \right). \end{gathered} \right\} \end{equation}
As noted by Williams (Reference Williams1985), these conserved scalars essentially are types of Crocco relations which will be discussed below. Note that, for the non-reacting case, 
$h$ and 
$Y_m$ are conserved scalars satisfying this same partial differential equation.
For the higher-Mach-number case, (2.11) should include the viscous dissipation term. Then, for 
$Pr =1$, it can be combined with (2.9) to yield the equation
\begin{equation} \rho u\frac{\partial (h+u^2/2)}{\partial x} + \rho v\frac{\partial (h+u^2/2)}{\partial y} = \frac{\partial}{\partial y} \left( \mu \frac{\partial (h+u^2/2)}{\partial y} \right) -\rho \varSigma^N_{m=1}h_{f,m} \omega_m. \end{equation}
In this higher-Mach-number scenario, it is useful to re-define 
$\beta \equiv h + u^2/2 + \nu Y_O Q$ in order to satisfy (2.14).
The boundary conditions for 
$u, h, Y_m, \alpha$ and 
$\beta$ at 
$y = \infty$ and at 
$y= - \infty$ remain constant as both 
$x$ and 
$z$ vary. The boundary values 
$\kappa (\infty ) = \kappa _{\infty }$ and 
$\kappa (- \infty ) = \kappa _{-\infty }$ will vary with 
$x$. As downstream distance 
$x$ increases, the mixing-layer solution becomes less dependent on upstream inflow profiles; in fact, it is known to become independent asymptotically, depending only on the boundary conditions.
When 
$Pr = 1$, the solutions for 
$\alpha , \beta$ and, in the non-reacting case, 
$h$ and 
$Y_m$, become linear functions of 
$u$ as known through the classical Crocco integral (Crocco Reference Crocco1932). In the non-reacting case where terms of the order of the square of the Mach number are retained and 
$Pr =1$, we repeat Crocco's finding that 
$h + u^2/2$ is linear in 
$u$.
2.3. Similarity in one-dimensional form
Our development of the similarity follows a pattern originally developed for compressible boundary layer flow. However, the appearance of the 
$z$-momentum equation causes a need for new features.
Howarth (Reference Howarth1948) assumed a perfect gas with dynamic viscosity 
$\mu$ directly proportional to temperature 
$T$ and used a transformation of variable 
$\bar {y} \thicksim \int (\nu )^{-1/2}\,{\textrm {d}y}$ where 
$\nu$ is the kinematic viscosity. This leads to 
$\bar {y} \thicksim p^{1/2} \int (1/T)\, {\textrm {d}y}$. Stewartson (Reference Stewartson1949) simply stated 
$\bar {y} \thicksim \int \rho \,{\textrm {d}y}$ which does not require an assumption about the relation between temperature and viscosity. Dorodnitsyn (Reference Dorodnitsyn1942) did parallel work earlier. Lees (Reference Lees1956) generalized the transformation for situation where free stream velocity and pressure varied in the streamwise direction. The density-weighted transformation is used here to replace 
$y$ with 
$\bar {y} \equiv \int ^y_0 \rho (y') \,{\textrm {d}y}'$. The product 
$\rho \mu$ is assumed to remain constant and equal to 
$\rho _{\infty } \mu _{\infty }$ throughout the mixing layer. (This is convenient but not necessary to obtain a similar solution.) The similarity variable 
$\eta$ and other variables are defined as
\begin{equation} \left. \begin{gathered} \eta \equiv \frac{\bar{y}}{g(x)},\\ g(x) \equiv \sqrt{\frac{2 \rho_{\infty} \mu_{\infty} x}{u_{\infty}} }, \\ G(\eta) \equiv \frac{\kappa g^2}{\rho_{\infty} \mu_{\infty}} = \frac{2 \kappa x} {u_{\infty}} =\frac{2 w x}{z u_{\infty}}, \\ E \equiv \int_0^{\eta} G(\eta')\,\textrm{d} \eta'. \end{gathered} \right\} \end{equation}
Note the implication that 
$\kappa$ will vary as the reciprocal of 
$x$; 
$\eta , G$ and 
$E$ are non-dimensional here; 
$G$ is both the normalized 
$z$-component of velocity and the indicator of the imposed normal strain.
In standard fashion, a pseudo-streamfunction 
$\varPsi$ and a pseudo-
$y$-component of velocity 
$\tilde {v}$ are created mimicking an incompressible flow:
\begin{equation} \left. \begin{gathered} \varPsi = f(\eta)g(x) u_{\infty}, \\ u = \frac{\partial \varPsi}{\partial \bar{y}} = u_{\infty} \frac{\textrm{d}f}{\textrm{d}\eta}, \\ \rho_{\infty}\mu_{\infty}\tilde{v} = - \frac{\partial \varPsi}{\partial x} = \rho v + u\int_0^y\frac{\partial \rho}{\partial x}\,{\textrm{d}y}' + \sqrt{\frac{\rho_{\infty}\mu_{\infty}u_{\infty} }{2x}}E. \end{gathered} \right\} \end{equation}It follows that
\begin{equation} \tilde{v} g = \eta \frac{\textrm{d}f}{\textrm{d}\eta} - f. \end{equation}
Thereby, 
$\tilde {v}$ will vary as a function of the similarity variable multiplied by the reciprocal of the square root of 
$x$. So, 
$\tilde {v}g$ becomes a function only of the similarity variable 
$\eta$.
Generally, for two-dimensional mixing layers and boundary layers with variable density, the interest in the precise determination of 
$v$ beyond (2.17) has not been high because 
$v^2 \ll u^2$. Here, because of the imposed counterflow, we have interest in the determination of 
$vg$ as a function of 
$\eta$. The second term on the right-hand side of the second line of (2.17) requires attention in order to determine 
$v(\eta , x)$. In that term, a derivative with respect to 
$x$ is taken of an integral over 
$y$ space. We may write
\begin{equation} y = \int^{\bar{y}}_0 \frac{1}{\rho}\,\textrm{d}\bar{y}' =g(x) \tilde{I}(\eta), \end{equation}
where 
$\tilde {I}(\eta ) \equiv \int ^{\eta }_0 (1/\rho (\eta '))\,\textrm {d}\eta '$. Furthermore, taking the derivative at constant 
$y$,
\begin{equation} \left. \begin{gathered} \frac{\textrm{d}y}{\textrm{d}x} = 0 = \tilde{I}\frac{\textrm{d}g}{\textrm{d}x} +g\frac{\textrm{d}\tilde{I}}{\textrm{d}x} = \tilde{I}\frac{\textrm{d}g}{\textrm{d}x} +\left. \frac{g}{\rho}\frac{\textrm{d}\eta}{\textrm{d}x}\right|_{y = constant}, \\ \frac{\textrm{d}\bar{y}}{\textrm{d}x} = \eta \frac{\textrm{d}g}{\textrm{d}x} + g\frac{\textrm{d}\eta}{\textrm{d}x} = (\eta - \tilde{I} \rho)\frac{\textrm{d}g}{\textrm{d}x}. \end{gathered} \right\} \end{equation}
\begin{align} \rho v g \!=\! \rho_{\infty}\mu_{\infty}\left(\eta \frac{\textrm{d}f}{\textrm{d}\eta} - f\right) - u(\eta - \tilde{I} \rho)g \frac{\textrm{d}g}{\textrm{d}x}- g\sqrt{\frac{\rho_{\infty}\mu_{\infty}u_{\infty} }{2x}}E \!=\! \rho_{\infty}\mu_{\infty} \left[ \tilde{I}\rho\frac{\textrm{d}f}{\textrm{d}\eta}- f -E \right]. \end{align}
Libby & Liu (Reference Libby and Liu1968) produced the parallel analytical result for the two-dimensional boundary-layer flow but never computed 
$v$. Pruett (Reference Pruett1993) made 
$v$ calculations for a wall-layer similar solution. Kennedy & Gatski (Reference Kennedy and Gatski1994) made 
$v$ calculations for a mixing-layer similar solution. Those studies all considered a two-dimensional shear layer without imposed normal strain, i.e. 
$E =0$.
Using the perfect-gas relation and neglecting terms of order Mach number squared, we find
\begin{equation} \left(\frac{\rho}{\rho_{\infty}} \right)^{-1}= \frac{T}{T_{\infty}} = \frac{h}{h_{\infty}} = \frac{\mu}{\mu_{\infty}}. \end{equation}
For the perfect gas, 
$\tilde {I}$ may be written with 
$\tilde {h}$ as the integrand and (2.21) for 
$v$ can be modified:
\begin{equation} \rho v g = \rho_{\infty}\mu_{\infty} \left[ \frac{\int_0^{\eta}\tilde{h}(\eta')\,\textrm{d} \eta'}{\tilde{h}}\frac{\textrm{d}f}{\textrm{d}\eta}- f -E \right]. \end{equation}2.4. Formulation of ODEs
We define 
$( )' \equiv \textrm {d}( )/\textrm {d}\eta , \tilde {h} \equiv h/ h_{\infty }= \rho _{\infty }/\rho$, and 
$\tilde {\beta } \equiv \beta /h_{\infty }$. Now, the following ODEs and boundary conditions do follow:
\begin{gather} f''' + (\,f + E) f'' = 0;\quad f'(-\infty)=\frac{u_{-\infty}}{u_{\infty}}; \quad f(0) = 0; \quad f'(\infty)= 1, \end{gather}
\begin{gather} G'' + (\,f + E) G' - G^2 + G_{\infty}^2 \tilde{h} =0; \quad G(-\infty)= G_{-\infty}; \quad G(\infty)= G_{\infty}. \end{gather}
Equation (2.24a–d) differs from the classical Blasius equation because of the presence of 
$E$ which couples it to (2.25a–c). That equation may actually be interpreted as a third-order ODE since 
$G = \textrm {d}E/\textrm {d}\eta .$ In the coefficient of the first derivatives in the above and the following equations, 
$f$ describes the contribution to transverse transport due to the shear strain while 
$E$ represents the contribution to transverse transport due to the imposed normal strain. The boundary conditions on 
$G$ are chosen so that pressure will be the same with variation for 
$z$ within the two free streams. Thus, 
$\rho w^2$ will be same in the two streams for a given 
$x$ and 
$z$. If the two free streams have the same temperature (and therefore the same density), 
$G_{-\infty }= G_{\infty }$; otherwise, they must differ:
\begin{gather} \tilde{h}'' + Pr(\,f + E) \tilde{h}' = 2Pr\frac{\omega_F x}{u_{\infty}}\frac{Q}{h_{\infty}}; \quad \tilde{h}(-\infty)=\frac{h_{-\infty}}{h_{\infty}}; \quad \tilde{h}(\infty)= 1, \end{gather}
\begin{gather} Y_F'' + Pr(\,f + E) Y_F' =- 2Pr\frac{\omega_F x}{u_{\infty}}; \quad Y_F(-\infty)=Y_{F, -\infty}; \quad Y_F(\infty)= Y_{F,\infty}. \end{gather}
As written, 
$\omega _F$ is to be taken as negative when fuel is being consumed:
\begin{gather} \alpha'' + Pr(\,f + E) \alpha' = 0; \quad \alpha(-\infty)=\alpha_{-\infty}; \quad \alpha(\infty)= \alpha_{\infty}, \end{gather}
\begin{gather} \tilde{\beta}'' + Pr(\,f + E) \tilde{\beta}' = 0; \quad \tilde{\beta}(-\infty)=\tilde{h}_{-\infty}+ \frac{Q\nu Y_{O,-\infty}}{h_{\infty}}; \quad \tilde{\beta}(\infty)= 1 + \frac{Q\nu Y_{O,\infty}}{h_{\infty}}. \end{gather} Several important parameters can be identified in the equations. A Damköhler number 
$Da$ will be embedded in the chemical-rate function 
$\omega _F$ to be discussed later. Of course, the composition and temperature at 
$\eta = \infty$ and 
$\eta = - \infty$ will be influential; 
$Pr$ affects the mass and energy diffusion rate as compared to the diffusion of momentum due to viscosity and 
$G_{\infty }$ is the normalized magnitude of the imposed normal-strain rate 
$\kappa$ using 
$u_{\infty }/x$ (the reciprocal of a residence time) as the normalizing factor. The shear strain rate is estimated by 
$(u_{\infty } - u_{-\infty })/ \delta (x)$ where 
$\delta (x)$ is the mixing-layer thickness. Using the same normalization factor, the normalized shear strain is estimated as 
$(x/\delta )(1 - u_{-\infty }/u_{\infty })$. Clearly, the velocity ratio 
$u_{-\infty }/u_{\infty }$ is important. 
$\delta$ will be affected by both the shear strain and the compressive strain. For the pure counterflow without shear, 
$\delta _C = O((\nu /\kappa )^{1/2})$ where 
$\nu$ is the kinematic viscosity. For the pure shear layer without compressive normal strain, 
$\delta _S = O( x/Re_x^{1/2} ) = O((x\nu /u_{\infty })^{1/2})$. Thus, 
$\delta _S/\delta _C = O((\kappa x/u_{\infty })^{1/2} )$ and 
$G = O(( \delta _S/\delta _C )^2)$. So, high (low) values of 
$G$ will imply that the compressive (shear) strain plays a dominant role in determining layer thickness.
The physical problem has three characteristic times and, through their reciprocals, three characteristic rates: chemical reaction rate, shear-strain rate and normal-strain rate. This results in two important non-dimensional ratios of the times (or rates). Through the boundary conditions on 
$u$ and 
$G$ at plus infinity and minus infinity, the ratio of normal-strain rate to shear-strain rate is imposed. In § 2.7, a Damköhler number will be introduced; together with the ratio of velocities at plus and minus infinity, it reflects the ratio of chemical rate to shear-strain rate.
Equations (2.24a–d)–(2.29a–c) could be made more general by considering variation of 
$\rho \mu$ through the mixing layer. Terms with the first derivative in 
$\eta$ space of the normalized product 
$\rho \mu / (\rho _{\infty }\mu _{\infty })$ would appear in the equations. This secondary effect will be neglected here. In the classical two-dimensional mixing layer without imposed normal strain (i.e. 
$\kappa = 0, G =0$ and 
$E=0$), (2.25a–c) disappears while (2.24a–d) and (2.26a–c)–(2.29a–c) simplify to the well-known forms. For the case where no counterflow is applied (i.e. 
$E = 0$), (2.24a–d) has been extended to cases where 
$\textrm {d}(\rho \mu )/\textrm {d}\eta \neq 0$. Poblador-Ibanez, Davis & Sirignano (Reference Poblador-Ibanez, Davis and Sirignano2020) have examined the solution of (2.21) where 
$\textrm {d}(\rho \mu )/\textrm {d}\eta \neq 0$ and 
$E=0$. Also, they show favourable comparisons between similar solutions and two-dimensional computational results for 
$v$.
2.5. Generalized Crocco integral
Crocco (Reference Crocco1932) developed his integral solution for a compressible wall boundary layer at significant values of free-stream Mach number. Assuming 
$Pr = 1$, he showed that, knowing the solution for velocity 
$u$ and the boundary conditions for the conserved scalar, a similar solution existed with the conserved scalar 
$h_o \equiv h + u^2/2$ linear in 
$u$. Here, we show that, for 
$Pr \neq 1$, a conserved scalar solution may still be found as an integral function of the velocity derivative. In a non-reacting case at Mach number 
$M \ll 1$, 
$h_o \approx h$; thus, the static enthalpy 
$h$ may be approximated as a conserved scalar because viscous dissipation is negligible. We also show that the enthalpy, with account for viscous dissipation at any subsonic 
$M$ value and any 
$Pr$ value, can be described as an integral function of the velocity derivative.
Visualize (2.24a–d) as a first-order ODE governing 
$\,f'' = \textrm {d}u/\textrm {d}\eta$ and (2.26a–c)–(2.29a–c) as first-order ODEs governing the first derivative of the dependent variable. In the case without imposed normal strain (i.e. 
$G = 0$), the solutions for 
$f(\eta )$ and 
$u(\eta )$ will be independent of 
$Pr$ and 
$M$. (Realize that a coupling actually applies through 
$\eta (x,y)$.) Therefore 
$G$ and thus 
$E$ are coupled with enthalpy and thereby have dependencies on both 
$Pr$ and 
$M$. In that case, 
$f$ and 
$u$ will also experience the dependencies. Consider the conserved scalars 
$h$ and 
$Y_m$ in the non-reacting case and 
$\alpha$ and 
$\beta$ in any case. Clearly, a solution exists where the first derivative of the scalar is linear in 
$(\textrm {d}u/\textrm {d}\eta )^{Pr}$. In fact, since all derivatives go to zero at plus and minus infinity, we have a direct proportionality. Subsequently, we obtain the generalized Crocco integral for any conserved scalar 
$CS$:
\begin{equation} \frac{CS(\eta) - CS_{-\infty}}{CS_{\infty} - CS_{-\infty}} = \dfrac{\displaystyle\int_{-\infty}^{\eta} \left(\dfrac{\textrm{d} u}{\textrm{d}\eta}\right)^{Pr}\,\textrm{d}\zeta}{\displaystyle\int_{-\infty}^{\infty} \left(\dfrac{\textrm{d}u}{\textrm{d}\eta}\right)^{Pr}\,\textrm{d}\eta} = \dfrac{\displaystyle\int_{-\infty}^{\eta} \left(\,f''(\zeta)\right)^{Pr}\,\textrm{d}\zeta}{\displaystyle\int_{-\infty}^{\infty} \left(\,f''(\eta)\right)^{Pr}\,\textrm{d}\eta} = \frac{J(\eta)}{J(\infty)}, \end{equation}
where solution of (2.24a–d) as a first-order ODE for 
$f''$ yields
\begin{equation} J(\eta) \equiv \int^{\eta}_{-\infty} \left(\,f''(\zeta)\right)^{Pr}\,\textrm{d}\zeta =\int_{-\infty}^{\eta} e^{-Pr I(\zeta)}\,\textrm{d}\zeta; \quad I(\eta) \equiv \int_{0}^{\eta} [f + E ]\,\textrm{d}\zeta. \end{equation}
Equations (2.24a–d) and (2.25a–c) must still be integrated to determine the functions 
$f$ and 
$E$ which appear in the integrand of (2.31a,b). The integration will couple with other equations since 
$\tilde {h}$ appears in (2.25a–c). For the non-reacting case, the functions 
$\tilde {h}, Y_F, \alpha$ and 
$\beta$ may be determined using the generalized Crocco integral in (2.30). For the reacting case, the integral may be used to evaluate 
$\alpha$ and 
$\beta$ but either 
$\tilde {h}$ or 
$Y_F$ must be solved from (2.26a–c) or (2.27a–c).
Obviously, in cases where 
$Pr$ and Schmidt number 
$Sc$ differ, 
$Sc$ should appear in the exponent for mass fraction (for the non-reacting case only) and 
$\alpha$ solutions. In such a case, (2.26a–c) must be re-formulated since it is now based on 
$Pr = Sc$.
It is possible to generalize the Crocco integral also in the case where Mach-number squared terms are not neglected. Here, under the boundary-layer approximation the dissipation term 
$\mu (\partial u/\partial y)^2$ is added to the right-hand sides of (2.4) and (2.11). Consequently, the term 
$-Pr (u_{\infty }^2/h_{\infty })(\,f'')^2$ must be added to the right-hand side of (2.26a–c). The non-reacting case will be considered so that 
$\omega _F =0$ in that equation. Again, the second-order ODEs will be treated as first order for determining the first derivative. Then, a simple integration follows for determining 
$\tilde {h}$ from knowledge of its first derivative. The result is
\begin{equation} \left. \begin{aligned} \tilde{h} & = \frac{h_{-\infty}}{h_{\infty}} + C\frac{\displaystyle\int_{-\infty}^{\eta}(\,f'')^{Pr}\,\textrm{d}\zeta}{ \displaystyle\int_{-\infty}^{\infty}(\,f'')^{Pr}\,\textrm{d}\zeta } - Pr\frac{u_{\infty}^2}{h_{\infty}}\int_{-\infty}^{\eta}(\,f'')^{Pr}\left[\int_{-\infty}^{\zeta} (\,f'')^{2 - Pr}\,\textrm{d} \theta\right]\,\textrm{d} \zeta,\\ & = \frac{h_{-\infty}}{h_{\infty}} + C \frac{J(\eta)}{J(\infty)} - Pr \varOmega^2\frac{u_{\infty}^2}{h_{\infty}}\int_{-\infty}^{\eta}\int_{-\infty}^{\zeta} e^{-PrI(\zeta) -(2-Pr)I(\theta)}\,\textrm{d}\theta\,\textrm{d}\zeta,\\ C & \equiv 1 - h_{-\infty}/h_{\infty} +(Pr u_{\infty}^2/h_{\infty})\int_{-\infty}^{\infty}(\,f'')^{Pr} \left[\int_{-\infty}^{\eta}(\,f'')^{2 - Pr}\,\textrm{d}\zeta \right]\,\textrm{d}\eta,\\ & = 1 - h_{-\infty}/h_{\infty} + (Pr \varOmega^2u_{\infty}^2/h_{\infty})\int_{-\infty}^{\infty} \int_{-\infty}^{\eta} e^{-Pr I(\eta)-(2 - Pr)I(\zeta)}\,\textrm{d}\zeta\,\textrm{d}\eta,\\ \varOmega & \equiv f''_{max} = \frac{f'(\infty) - f'(-\infty)}{\displaystyle\int_{-\infty}^{\infty} e^{-I(\eta)}\,\textrm{d}\eta }. \end{aligned} \right\} \end{equation}
Here, 
$\varOmega$ is a non-dimensional indicator of the magnitude of the maximum shear-strain rate. Upon neglect of terms of 
$O(u_{\infty }^2/h_{\infty })$, (2.32) reduces to (2.30). With imposed normal strain, 
$f$ couples to 
$G$, which in turn couples to 
$\tilde {h}$; thus, (2.32) is not closed in that case. When 
$Pr =1$, the integrals can be evaluated analytically yielding the original Crocco integral relation where 
$h + u^2/2$ becomes linear in 
$u = u_{\infty }\,f'$. For the perfect gas, the parameter 
$u_{\infty }^2/h_{\infty } = (\gamma -1)M^2$. Here, 
$M$ is the Mach number of the free stream at 
$y =\infty$ which is generally the faster stream in the calculations presented. Illingworth (Reference Illingworth1949) in a flat-plate boundary-layer analysis produced an integral similar to the form of the second line of (2.32) to evaluate enthalpy and presented some calculations for 
$Pr =0.725$. However, the connection with the velocity derivative and the nature of the result as a generalized Crocco integral was not mentioned.
The analytical results in this subsection apply with or without imposed normal strain. They can be applied to a wall boundary layer by simply replacing boundary conditions at minus infinity by the wall boundary conditions. Therefore, we have a general analytical relation for the scalar variable in terms of the Blasius solution.
2.6. Non-similar behaviour of the reaction zone
There is a challenge in defending the claim that (2.26a–c) and (2.27a–c) are in similar form because of the appearance of the factor 
$\omega _F x/ u_{\infty }$ in the source and sink terms; 
$\omega _F$ will depend on pressure, temperature and mass fractions which are expected to be similar (or at least near similar) while 
$x$ is clearly non-similar. We cannot expect that 
$\omega _F \sim 1/x$. A suitable explanation can be made using the concept of inner and outer solutions from singular perturbation theory (Van Dyke Reference Van Dyke1964; Kevorkian & Cole Reference Kevorkian and Cole1996).
The reaction rates are negligible outside of thin reaction zones within the mixing layer. Within those zones, second derivatives become notably larger than first derivatives and (2.26a–c) and (2.27a–c) become
\begin{equation} \left. \begin{gathered} \tilde{h}'' \approx 2Pr\frac{\omega_F x}{u_{\infty}}\frac{Q}{h_{\infty}}; \\ Y_F'' \approx - 2Pr\frac{\omega_F x}{u_{\infty}}. \end{gathered} \right\} \end{equation}These relations determine our inner solutions which must be matched to the outer solutions which apply to the much larger portion of the mixing layer where reaction rate is negligible. Integration of the inner solutions over the reaction zone yields
\begin{equation} \left. \begin{gathered} \tilde{h}'|_+ - \tilde{h}'|_- \approx 2 Pr \frac{Q}{h_{\infty}} \int_-^+\frac{\omega_F x}{u_{\infty}}\,\textrm{d}\eta; \\ Y_F'|_+ - Y_F'|_- \approx - 2 Pr \int_-^+ \frac{\omega_F x}{u_{\infty}}\,\textrm{d}\eta. \end{gathered} \right\} \end{equation}The ‘plus’ and ‘minus’ subscripts and integral limits denote the far edges of the inner zone: namely, the infinity limits on the small inner scale in singular perturbation theory. These results apply for both diffusion flames and premixed flames.
First, let us discuss diffusion flames. Physically, we have a diffusion-controlled situation where peak temperature in the reaction zone varies weakly. Thus, the peak value of the reaction rate per unit volume 
$\omega _F$ varies weakly with downstream position and the reaction-zone thickness adjusts to have production and consumption rates match the diffusion rate. The gradients of the scalars in the 
$y$ space will decrease as 
$1/\sqrt {x}$ which implies no change with 
$x$ for the derivatives in 
$\eta$ space. The important point is that the width of the reaction zone measured in 
$\eta$ space narrows with increasing downstream distance 
$x$. In particular, the measure 
$\Delta y$ of the zone width behaves as 
$1/\sqrt {x}$ which implies that the ratio of reaction-zone thickness to mixing-layer thickness behaves as 
$1/x$. Consequently, 
$\Delta \eta$ for the reaction zone varies as 
$1/x$. Thus, the integrals with their integrands proportional to 
$x$ will actually not vary with 
$x$ so that the jump in first-derivative values with respect to 
$\eta$ are not dependent on 
$x$. The jumps give the important quantities of heat production and fuel-mass consumption in the reaction zone and its impact on the remainder of the mixing layer. In other words, the rate of heat and mass diffusion in or out of the reaction zone depends solely on the variable 
$\eta$, giving similarity for the bulk of the mixing layer. An important point here is that, under the boundary-layer approximation, the communication between the reaction zone at any 
$x$ position and the rest of the mixing layer occurs only in the 
$y$ direction through transport and diffusion with the latter dominant; there is no direct transfer of information from the reaction zone at one 
$x$ value to the reaction zone at another 
$x$ value. Note that the infinite-kinetics model yields consistent behaviour with our results for the derivatives outside the reaction zone.
A premixed flame which is not excessively fuel lean or fuel rich will propagate in a wave-like manner relative to the combustible mixture. The rate of fuel consumption and heat production (and the flame speed) depends on the product of diffusion rate and chemical reaction rate. Flame speed (relative to the upstream incoming fluid velocity), flame thickness and scalar gradients in 
$y$ space will not vary with 
$x$. The implication therefore is that flame speed, flame thickness, and scalar gradients in 
$\eta$ space vary as 
$1/\sqrt {x}, 1/\sqrt {x},$ and 
$\sqrt {x}$, respectively. The imposed counterflow velocity decreases as 
$1/\sqrt {x}$, and therefore the premixed flame transverse motion across the mixing layer could not be arrested at the same 
$\eta$ value for all 
$x$ values. Premixed-flame position would change with 
$x$ position, moving away from 
$\eta =0$, in both 
$y$ space and 
$\eta$ space while diffusion flame position would not change in 
$\eta$ space. Thus, the similar solution that predicts premixed or partially premixed flames of the classical form with wave-like propagation cannot be accurate and should only be trusted to support the plausibility of the premixed-flame occurrence. The premixed configuration deserves further examination with a two-dimensional analysis addressing (2.8)–(2.14). In this work, premixed flames with mixture ratios far from stoichiometric will be examined; they are more likely to be diffusion controlled with weak influence of chemical kinetics in determining a flame speed.
For reacting counterflow calculations with one-step propane–oxygen kinetics, Sirignano (Reference Sirignano2019b) found that, in a configuration with a diffusion flame and a fuel-rich flame, the latter flame did not show a tendency to propagate at speed proportional to the square root of the integrated reaction rate through the reaction zone. The change in flame speed with increasing Damköhler number (
$Da$) was notably weaker than the expected square root dependence for a premixed laminar flame. The situation was closer to one where the width of the reaction zone increased with the reciprocal of the reaction rate; thus, there was little variation in the integrated reaction rate or total consumption rate for fuel over the zone. Accordingly, the mean temperature and jumps in enthalpy gradient and mass fraction gradient across the reaction zone did not differ much depending on 
$Da$. It appears that the flame was a zone whose volume adjusted (at constant reaction rate) to accommodate the mass flux rate passing through the fuel-rich flame and moving towards the diffusion flame.
2.7. Chemical kinetic model
In a one-step chemical reaction, each species is consumed or produced at a rate in direct proportion to the rate of some other species that is produced or consumed. We will focus on propane–oxygen flows with one-step kinetics. However, results are expected to be qualitatively more general, applying to situations with more detailed kinetics and to other hydrocarbon/oxygen-or-air combination. Westbrook & Dryer (Reference Westbrook and Dryer1984) kinetics are used; they were developed for premixed flames but any error for non-premixed flames is viewed as tolerable here because diffusion would be rate controlling. The reaction rate (rate of change of fuel mass fraction) in units of reciprocal seconds is given as
\begin{equation} \omega_F = - A {\rho}^{0.75} Y_F^{0.1} Y_O^{1.65} e^{-50.237/\tilde{h}}, \end{equation}
where the ambient reference temperature is set at 300 K and density 
$\rho$ is to be given in units of kilograms per cubic metre. Here, 
$A =4.788 \times 10^8\ {(\textrm {kg}\ \textrm {m}^{-3})}^{-0.75}/s$. The dimensional reciprocal of residence time 
$u_{\infty }/x$ is used to normalize time and reaction rate. In non-dimensional terms,
\begin{equation} \left. \begin{gathered} \frac{2\omega_F x}{u_{\infty}}= - \frac{2 A {\rho_{\infty}}^{0.75}}{u_{\infty}/x}\tilde{h}^{-0.75} Y_F^{0.1} Y_O^{1.65} e^{-50.237/\tilde{h}}, \\ \frac{2\omega_F x}{u_{\infty}}= - \frac{Da}{\tilde{h}^{0.75}} Y_F^{0.1} Y_O^{1.65} e^{-50.237/\tilde{h}}. \end{gathered} \right\} \end{equation}
The above equation defines the Damköhler number 
$Da$. Furthermore, we set 
$Da \equiv K Da_{ref}$ where
\begin{equation} Da_{ref} \equiv \frac{{2A(10\ \textrm{kg}\ \textrm{m}^{-3})}^{0.75}}{(20/s)} =2.693 \times 10^6; \quad K \equiv \left[\frac{\rho_{\infty}}{10\ \textrm{kg}\ \textrm{m}^{-3}}\right]^{0.75}\frac{20/s}{u_{\infty}/x}, \end{equation}
where 
$\rho _{\infty } =10\ \textrm {kg}\ \textrm {m}^{-3}$ and 
$u_{\infty }/x =20/s$ are arbitrarily chosen as reference values for density and reaction rate, respectively. The reference value for density implies an elevated pressure. The 
$20/s$ reference value is in the middle of an interesting range for this chemical reaction. Clearly, there is no need to set pressure (or its proxy, density) and the strain rate separately for a one-step reaction. For propane and oxygen, the mass stoichiometric ratio 
$\nu =0.275$.
The non-dimensional parameter 
$K$ will increase (decrease) as 
$u_{\infty }/x$ decreases (increases) and/ or the pressure increases (decreases). 
$K =1$ is our base case and the range covered will include 
$10^{-2} \leq K\leq 3$.
2.8. The 
$\varSigma$ space
Flamelet theory (Peters Reference Peters2000; Pierce & Moin Reference Pierce and Moin2004) has evolved with the use of a conserved scalar as the independent variable replacing the 
$y$ or 
$\eta$ coordinate. Bilger (Reference Bilger1976) has emphasized the use of element-based mass fractions which become conserved scalars because chemistry does not destroy atoms but only changes molecules. Bilger refers to it as a ‘mixture fraction’. It only is useful as a replacement for 
$y$ if it remains monotonic in 
$y$. That will be always true for the steady-state, one-dimensional case; however, in the multi-dimensional case or in the one-dimensional unsteady case, monotonic behaviour of the upstream boundary conditions or initial conditions becomes a requirement. Thereby, the use of the mixture fraction is not always optimal. Sirignano (Reference Sirignano2019a) argued that any conserved scalar that was monotonic varying in the direction normal to the flame surface would suffice. In fact, it need not have physical meaning. A simple option is to use the solution 
$\varSigma$ of the following equation and boundary conditions:
\begin{equation} \left. \begin{gathered} \varSigma'' +Pr (\,f + E) \varSigma'= 0, \\ \varSigma(-\infty) = 0; \quad \varSigma(\infty) = 1. \end{gathered} \right\} \end{equation}Using the generalized Crocco integral, the result is
\begin{equation} \varSigma(\eta) = \frac{J(\eta)}{J(\infty)}; \quad J(\eta) \equiv \int_{-\infty}^{\eta} e^{-I(\eta')}\,\textrm{d}\eta'; \quad I(\eta) \equiv \int_{-\infty}^{\eta} Pr[f + E ]\,\textrm{d}\zeta. \end{equation}
Here, 
$\varSigma$ equals a normalized steady-state conserved scalar. For two examples at any 
$Pr$ value,
\begin{align} \varSigma = \frac{\alpha(\eta) - \alpha(-\infty)}{\alpha(\infty) - \alpha(-\infty)} = \frac{\beta(\eta) - \beta(-\infty)}{\beta(\infty) - \beta(-\infty)}. \end{align}
For 
$Pr =1$, the most natural choice for a shear layer is
\begin{equation} \varSigma = \frac{u(\eta) - u_{\infty}}{u_{\infty}- u_{-\infty}}. \end{equation}
For a general value of 
$Pr$, 
$\varSigma$ as a functional of 
$u$ can still be given. These results apply for the single diffusion flame as well as for a multiple-flame configuration. For the diffusion flame case with oxygen and fuel only coming in opposing flows, the mixture fraction 
$Z$ commonly used is simply the normalized 
$\alpha$ conserved scalar. Thereby, 
$Z = \varSigma$ in that case; however, in a broader set of problems, they are not always the same. 
$Z$ will not vary from zero to unity if the bounding compositions are not pure fuel and pure oxidizer.
In similar fashion to Peters (Reference Peters2000), the independent variable 
$y$ can be replaced by 
$\varSigma (y)$ in (2.27a–c). The result is
\begin{equation} \left. \begin{gathered} \chi \frac{\textrm{d}^2 Y_F}{\textrm{d} \varSigma^2}+ Pr\frac{\omega_F x}{u_{\infty} }= 0,\\ \chi \frac{\textrm{d}^2 \tilde{h}}{\textrm{d} \varSigma^2}- Pr\frac{Q}{h_{\infty}}\frac{\omega_F x}{u_{\infty} }= 0, \\ \chi \equiv \frac{1}{2} \left(\frac{\textrm{d} \varSigma}{\textrm{d} \eta}\right)^2 = \frac{1}{2 \rho^2} \left(\frac{\textrm{d} \varSigma}{\textrm{d}y}\right)^2 = \frac{1}{2} \frac{e^{-2I(\eta)}}{ J^2(\infty) }, \end{gathered} \right\} \end{equation}
where 
$\chi$ is commonly named the scalar dissipation rate. However, for laminar mixing-layer and boundary-layer flows at 
$Pr =1$, it is better described as one half of the square of the strain rate. In (2.42), 
$\eta (\varSigma ) = J^{Inv} (\varSigma J(\infty ))$ must be substituted where 
$J^{Inv}$ is the inverse function of 
$J$ which must be determined numerically or approximated.
Equation (2.24a–d) also presents 
$f'$ as a ‘conserved scalar’; thus, it will also be linear in 
$\varSigma$. In fact, the velocity field can be used to determine 
$\varSigma$:
\begin{equation} \left. \begin{gathered} \varSigma \equiv \frac{J(\eta)}{J(\infty)} = \frac{\displaystyle\int^{\eta}_{-\infty} \left(\,f''(\zeta)\right)^{Pr}\,\textrm{d}\zeta}{ \displaystyle\int^{\infty}_{-\infty} \left(\,f''(\eta)\right)^{Pr}\,\textrm{d}\eta },\\ \chi = \frac{1}{2} \frac{\left(\textrm{d}u/\textrm{d}\eta\right)^{2 Pr}}{ \left(\displaystyle\int^{\infty}_{-\infty} \left(\textrm{d}u/\textrm{d}\eta\right)^{Pr}\,\textrm{d}\eta \right)^2}. \end{gathered} \right\} \end{equation}
The ODEs for 
$Y_F$ and 
$\tilde {h}$ given in (2.42) still apply. The only parameter affecting the solutions for 
$Y_F(\varSigma )$ and 
$\tilde {h}(\varSigma )$ is the product 
$Pr\,Da= 2.693\times 10^6 PrK$. For most of the 
$\varSigma$ space, the reaction rate is negligible and a linear relation without explicit dependence on any parameter appears. However, matching with the solution in the reaction zone affects the slope of the linear relation. An explicit dependence on 
$PrK$ can only occur in the narrow reaction zone between two larger linear domains which must be matched with the linear regions.
The values of 
$\tilde {h}$ and 
$Y_F$ will have significant variation within a narrow region around each of the 
$\varSigma$ values where a reaction zone exists. On both sides of that narrow region, 
$\tilde {h}$ and 
$Y_m$ will be linear in 
$\varSigma$. Equations for the conserved scalars, 
$\alpha$ and 
$\beta$, can be created, producing homogeneous equations with linear solutions in 
$\varSigma$. There is little reason though to solve the ODEs in 
$\varSigma$ space since they couple back to 
$f(\eta )$ and 
$G(\eta )$. It is more sensible to solve (2.24a–d)–(2.29a–c) or (2.24a–d)–(2.27a–c) and (2.30) with use of (2.40) than to integrate the ODEs in 
$\varSigma$ space. Equation (2.42) is interesting but really is only useful for producing solutions with constant-density counterflow where 
$f$ becomes linear in 
$y$ and 
$E$ is not present.
2.9. Numerical method
The system of ODEs is solved numerically using a relaxation method and central differences. Typically, solution over the range 
$-3 \leq \eta \leq 3$ provides adequate fittings to the asymptotic behaviours. Most calculations have 
$K = 1, Pr = 1, u_{-\infty }/ u_{\infty } =0.25$, and 
$G_{\infty }=1$ with emphasis on the effect of variation on composition of the free streams. However, the effects of 
$K, Pr, u_{-\infty }/ u_{\infty }$, 
$G_{\infty }$ and 
$G_{-\infty }$ variations are shown as well. The ambient temperatures in the two incoming streams are generally taken to be identical here at a value of 300 K with calorically perfect-gas relations yielding density and enthalpy. Typically, except where noted, the temperatures, densities, and enthalpies of the two ambient incoming flows are identical to each other. Sirignano (Reference Sirignano2019b) has shown that the difference in ambient temperatures can have some importance on the flow. Several qualitatively different cases with regard to the compositions of the two free streams are examined.
3. Results
The results for a non-reacting case will be discussed first in the next subsection. In the following subsections, the results for a configuration with a mixing layer containing a diffusion flame will be considered followed by an examination of results for mixing layers with two and three flames. Consequences of the imposed normal strain, differences in the two free-stream velocities, the difference between kinematic viscosity and thermal diffusivity will be considered through the variations of 
$G, u_{-\infty }/u_{\infty }$ and 
$Pr$, respectively. The importance of compressibility and viscous dissipation will be considered in the non-reacting case through the variation of 
$M$ and the consideration of kinetic energy and dissipated energy. In the configurations with flames, the impact of chemical reaction rates will be examined through variation of 
$Da$.
3.1. Non-reacting mixing layer
For the non-reacting case, the effects of Mach number 
$M$ are considered. Thereby, kinetic energy and viscous dissipation are not neglected although they are for reacting flows. The impact of imposed normal strain via counterflow and the thermal diffusivity are addressed through variation of 
$G$ and 
$Pr$. Both the velocity and scalar fields are examined. The consequences of the generalized Crocco relation is presented.
Figure 2 applies to a case where the two free streams have the same thermodynamic conditions. As shown by (2.22), the normalized enthalpy in our figures also represents the normalized absolute temperature. The small increase of temperature and enthalpy within the mixing layer in figure 2(a) here is due only to viscous dissipation. The amount of energy dissipated is not large compared to the thermal energy in the free stream, even for the significant value of 
$M = 0.5$. In figure 2(b), it is shown that, here with 
$Pr =1$, the classical Crocco relation holds. Figures 3–7 have results for cases where the two ambient temperatures are distinct and the impact of viscous dissipation is not so readily seen in the plots where enthalpy is monotonic across the layer. The value of 
$M$ is determined by using the higher free-stream velocity at 
$y =\infty$ and the lower free-stream speed of sound at 
$y= \infty$; Therefore, it is the highest Mach number in the flow field. Figures 3 and 4 show variations of 
$h, h_o, f, E, u/u_{\infty }, vg/\mu _{\infty }$ and 
$G$ across the mixing layer for 
$Pr = 0.7, 1.0$ and 
$1.3$. The effect of 
$Pr$ on 
$h$ and therefore on 
$h_o$ is significant; its effect on the transverse velocities 
$v$ and 
$w$ is slight, as shown via 
$vg/\mu _{\infty }$ and 
$G$; and the effect on 
$f, E$ and 
$u/u_{\infty }$ is negligible.The classical linear relation between 
$h_o$ and 
$u$ is not followed in figure 3(c) unless 
$Pr =1$. The differences in values for 
$h_o$ appear to be 
$O(Pr -1)$ as expected from (2.32). The viscous dissipation levels remain low throughout this 
$Pr$ range.
Figure 2. Non-reacting mixing layer with imposed normal strain, effect of viscous dissipation and Crocco integral solution: (a) 
$h(\eta )/h_{\infty }$ and (b) 
$h_o/h_{\infty }$ versus 
$u/u_{\infty }$. Here 
$Pr =1.0; G_{\infty } = 1.0; u_{-\infty } / u_{\infty }= 0.25; M =0.5$.
Figure 3. Non-reacting mixing layer with imposed normal strain, effect of Prandtl number: (a) 
$h(\eta )/h_{\infty },\ (b)\ h_o(\eta )/h_{\infty }$ and (c) 
$h_o/h_{\infty }$ versus 
$u/u_{\infty }.\ \text{Here } Pr =0.7, 1.0, 1.3; G_{\infty } = 1.0; u_{-\infty } / u_{\infty }= 0.25; M =0.5$.
Figure 4. Non-reacting mixing layer with imposed normal strain, effect of Prandtl number: (a) 
$f(\eta ),\, (b)\, E(\eta ),\, (c)\, u(\eta )/u_{\infty },\, (d)\, v(\eta )g/\mu _{\infty }$ and (e) 
$G(\eta )$. Here 
$Pr =0.7, 1.0, 1.3; G_{\infty } = 1.0; u_{-\infty } / u_{\infty }= 0.25; M =0.5$.
Now, the effect of 
$M$ is shown in figure 5 for 
$Pr =1$. There is negligible effect on static enthalpy but the influence on kinetic energy and therefore on 
$h_o$ is clear. The slope of 
$h_o(u)$ increases with 
$M$.
Figure 5. Non-reacting mixing layer with imposed normal strain, effect of Mach number: (a) 
$h(\eta )/h_{\infty },\ (b)\ h_o(\eta )/h_{\infty }$ and (c) 
$h_o/h_{\infty }$ versus 
$u/u_{\infty }$. Here 
$Pr =1.0; G_{\infty } = 1.0; u_{-\infty } / u_{\infty }= 0.25; M =0.25, 0.50, 0.75$.
The rate of normal strain is varied in figures 6 and 7. As the non-dimensional rate of strain 
$G$ increases with increase in 
$G_{\infty }$, the mixing layer gets thinner and the slopes of both scalars and velocity increase. The linear relation between 
$h_o$ and 
$u$ does not change as 
$G_{\infty }$ increases for 
$Pr = 1$. For 
$Pr\neq 1$, 
$h_o$ is no longer a conserved scalar; that is, viscous dissipation would appear in an equation describing it. Thereby, the function 
$h_o(u)$ becomes influenced by 
$Pr$ and linearity can be lost. The gradients of 
$h$ and 
$h_o$ increase as 
$Pr$ increases. As 
$G_{\infty }$ increases, the values of 
$G, E$ and 
$v$ tend to increase substantially affecting layer thickness, gradients and transport rates. In the case where 
$G=0$ and 
$E = 0$, the value of 
$v$ is very small; however, it is not identically zero everywhere. Without imposed normal strain, a normal-strain rate and a non-zero 
$v$ occur because the accelerating flow for 
$\eta <0$ has converging streamlines while the decelerating flow for 
$\eta >0$ has diverging streamlines, each producing a positive 
$v$. Realize that we have arbitrarily taken 
$v=0$ at 
$y=0$ (i.e. 
$\eta =0$) to represent the dividing streamline between the two streams. The exact (small) 
$y(x)$ value for the dividing streamline can easily be found by a well-known correction that sets the 
$y$-momentum flux of the two streams in balance. The correction does not modify strain rates or mixing rates; thus, we bypass it here.
Figure 6. Effect of normal-strain rate on non-reacting mixing layer: (a) 
$h(\eta )/h_{\infty }, (b)\ h_o(\eta )/h_{\infty }$ and (c) 
$h_o/h_{\infty }$ versus 
$u/u_{\infty }$. Here 
$Pr =1.0; G_{\infty } = 0, 1.0,4.0; u_{-\infty } / u_{\infty }= 0.25; M =0.5$.
Figure 7. Effect of normal-strain rate on non-reacting mixing layer: (a) 
$f(\eta ),\ (b)\ E(\eta ), (c)\ u(\eta )/u_{\infty }, (d)\ v(\eta )g/\mu _{\infty}$ and (e) 
$G(\eta )$. Here 
$Pr = 1.0; G_{\infty } = 0, 1.0, 4.0; u_{-\infty } / u_{\infty }= 0.25; M =0.5$.
3.2. Diffusion flame
Results for the basic parameter case with a single diffusion flame are shown in figures 8–10. The higher velocity stream at 
$y = \infty$ is composed purely of propane (i.e. 
$Y_{F,\infty } =1$) while the slower stream at 
$y = -\infty$ contains pure oxygen (i.e. 
$Y_{O,-\infty } =1$). In the basic case, the parameters are 
$K =1.0 , Pr =1.0 , G_{\infty } = 1.0$ and 
$u_{-\infty } / u_{\infty }= 0.25$. The reference temperature is the 300 K value taken as boundary conditions at 
$\eta = \infty$ and 
$\eta = - \infty$. This temperature determines the ambient values for enthalpy, density and dynamic viscosity.
Figure 8. Single diffusion flame solutions: (a) enthalpy 
$h/h_{\infty }$; (b) reaction rate 
$\omega _F x/u_{\infty }$; (c,d) mass fractions 
$Y_F, \ \nu Y_O$; and (e,f) scalars 
$\alpha , \beta .$ Here 
$K =1.0; Pr =1.0; G_{\infty } = 1.0; u_{-\infty } / u_{\infty }= 0.25$.
The first figure shows from the thermochemical variables that the reaction zone stationed around 
$\eta = -0.5$ is quite narrow and orders of magnitude thinner than the mixing layer. So, clearly, diffusion is rate controlling here, supporting the use of a similar solution.
Figure 9 shows that the velocity 
$u$ is monotonic through the shear layer, as expected in the absence of an imposed pressure gradient in flow direction. The Blasius function 
$f$ is also monotonic. However, the transverse velocities do not vary monotonically with 
$\eta$ as shown in the sub-figures for 
$vg/\mu _{\infty }$ and 
$G$; the variation of density due to heat release cause overshoots of this velocity component similar to findings for the counterflow (Sirignano Reference Sirignano2019a,Reference Sirignanob). In particular, the strain rates near the reaction zone are substantially augmented in magnitude compared to those imposed by the same free-stream flows in a non-reacting case. The integral 
$E$ is monotonic by nature. Its magnitude compares with 
$f$ in its effect on the solution. In fact, in this case, the imposed normal strain represented through 
$E$ has somewhat more influence than the shear strain represented through 
$f$.
Figure 9. Single diffusion flame solutions for dynamic field variables: (a) 
$f,\ (b)\ E, (c)\ u/u_{\infty },\ (d)\ vg/\mu _{\infty } \text{ and } (e)\ G$. Here 
$K =1.0; Pr =1.0; G_{\infty } = 1.0; u_{-\infty } / u_{\infty }= 0.25$.
The behaviour of enthalpy and mass fractions of the reactants are shown in figure 10 to be linear in 
$\varSigma$ space except within the narrow reaction zones. For our perfect gas, temperature will have the qualitatively identical behaviour. Thus, in analogy to the counterflow case, the 
$\varSigma$ space provides interesting information but it can only be calculated with coupling to the system of ODEs.Here, 
$\varSigma$ could be a conserved scalar or could be formed from the velocity field via 
$J(\eta )$.
Figure 10. Single diffusion flame solutions in 
$\varSigma$ space: (a) enthalpy 
$h/h_{\infty }$; and (b,c) mass fractions 
$Y_F,\ \nu Y_O$. Here 
$K =1.0; Pr =1.0; G_{\infty } = 1.0; u_{-\infty } / u_{\infty }= 0.25$.
In order to assure that the similarity approximation is reasonable here, we reduce the 
$Da$ severalfold up to an order of magnitude. In figure 11, the results for the 
$h/h_{\infty }, Y_f$ and 
$\nu Y_O$ are displayed in a fashion that zooms with an expansion on values around the reaction zone. The reaction zone increases somewhat in size and the peak value of enthalpy decreases with decreasing 
$Da$ but the zone still remains very narrow compared to mixing-layer dimensions. Furthermore, the asymptotic behaviours at the edge of the reaction zone remain independent of 
$Da$, indicating that the integral of the reaction rate is unchanged with changing 
$Da$. Thereby, we know that the similar solution is quite satisfactory for this diffusion-flame configuration. With 
$Da$ varying by an order of magnitude, the fields for 
$u/u_{\infty }, G, f$ and 
$E$ show no differences. The plot of scalars in 
$\varSigma$-space shows no significant difference either.
Figure 11. Effect of Damköhler number on single diffusion flame solutions: (a) enthalpy 
$h/h_{\infty }$; and (b,c) mass fractions 
$Y_F,\ \nu Y_O$. Here 
$K =0.1, 0.3, 1.0; Pr =1.0; G_{\infty } = 1.0; u_{-\infty } / u_{\infty }= 0.25$.
The Prandtl number 
$Pr$ is a parameter of importance. In practical situations, its magnitude will not vary that widely; however, it has consequences and its understanding is valuable. Increases in 
$Pr$ make scalar gradients increase as shown in figure 12. If 
$Pr > 1(<1)$, the normalized scalar gradients become greater than (smaller than) the normalized velocity gradients. Accordingly, as 
$Pr$ increases from 0.7 to 1.3, the peak value of the source in the energy equation (2.26a–c) becomes larger as shown in figure 12(b) since it contains 
$Pr$ multiplied by the reaction rate. The steeper gradients need not imply increasing transport rates since the diffusivities for mass and energy decrease as 
$Pr$ increases. There is a shift of the reaction zone towards the fuel side with increased 
$Pr$. The 
$u$-component of velocity has insignificant change due to variation of 
$Pr$ but more change is seen in the transverse components indicated by the behaviour of 
$G$ indicated in figure 13.
Figure 12. Effect of Prandtl number on single diffusion flame solutions: (a) enthalpy 
$h/h_{\infty }$; (b) reaction rate 
$\omega _F x/u_{\infty }$; (c,d) mass fractions 
$Y_F,\ \nu Y_O$. Here 
$K =1.0; Pr= 0.7, 1.0, 1.3; G_{\infty } = 1.0; u_{-\infty } / u_{\infty }= 0.25$.
Figure 13. Effect of Prandtl number on single diffusion flame solutions for dynamic field variables: (a) 
$E \text{ and } (b)\ G$. Here 
$K =1.0; Pr= 0.7, 1.0, 1.3; G_{\infty } = 1.0; u_{-\infty } / u_{\infty }= 0.25$.
The transverse normal-strain rate imposed through the counterflow has major effects as shown in figures 14 and 15. Scalar variables and their transport rates are affected in major ways. The effect of adding normal strain to a shear layer can be seen by starting with 
$G_{\infty } = 0$ and increasing it in those figures. The effect of adding shear to a pure counterflow can be inferred by starting with a large value of 
$G_{\infty }$ and decreasing it. With increasing normal strain, i.e. increasing value of 
$G_{\infty }$, the mixing layer gets thinner, increasing transport rates. Heat production rate and fuel consumption rate in the diffusion flame accommodate to the increased transport rate. Note that flame location moves towards the transverse velocity stagnation plane at 
$\eta = 0$ but the velocity at the peak reaction rate actually increases with increasing strain rate. Without imposed normal strain, 
$G=0$ and 
$E=0$ throughout the flow. Thus, 
$w=0$ everywhere in that case. However, the transverse velocity component 
$v$ has non-zero values due to converging and diverging streamlines and gas expansion caused by heat release. As before, we arbitrarily set 
$v(x,0) =0$ and bypass the opportunity to correct it based on transverse momentum balance.
Figure 14. Effect of normal-strain rate on single diffusion flame solutions: (a,b) enthalpy 
$h/h_{\infty }$ in both 
$\eta$ and 
$\varSigma$ spaces; (c,d) mass fractions 
$Y_F, \ \nu Y_O$. Here 
$K =1.0; Pr= 1.0; G_{\infty } = 0, 1.0, 4.0; u_{-\infty } / u_{\infty }= 0.25$.
Figure 15. Effect of normal-strain rate on single diffusion flame solutions for dynamic field variables: (a) 
$f,\, (b)\, E,\, (c)\, u/u_{\infty },\, (d)\, vg/\mu _{\infty }$ and (e) G. Here 
$K =1.0; Pr= 1.0; G_{\infty } = 0, 1.0, 4.0; u_{-\infty } / u_{\infty }= 0.25$.
The velocity ratio 
$u_{-\infty }/u_{\infty }$ does have some impact on the velocity field as shown in figure 16. As velocity ratio increases, the shear-strain rate decreases. While the imposed normal-strain rate remains fixed in value, it gains some importance relative to the shear; in particular, some narrowing of the mixing layer occurs due to the reduction of the flow displacement effect.
Figure 16. Effect of velocity ratio on single diffusion flame solutions for enthalpy and dynamic field variables: (a) 
$h/h_{\infty},\, (b)\, f,\, (c)\, u/u_{\infty },\, (d)\, G\, \text{ and } (e)\, E.$ Here 
$K =1.0; Pr =1.0; G_{\infty } = 1.0; u_{-\infty } / u_{\infty }= 0, 0.25, 0.50$.
There are small differences displayed in figure 17 concerning the orientation of shear strain relative to the scalar gradient when fuel is moved from the higher-speed stream to the lower-speed stream with the opposite movement for the oxygen. The reaction zone always locates towards the oxygen side but slightly more so when oxygen flows in the slower stream. The gradients are always greater on the oxygen side but a little more so when the oxygen stream is the slower but accelerating stream with more convergence of the streamlines. Negligible effect of the reversal is seen on the 
$u$ velocity or 
$f$ or on the conserved scalars. Some effect is shown on the 
$G$ and 
$E$ variables.
Figure 17. Effect of reactant-stream reversal on single diffusion flame solutions for enthalpy and dynamic field variables: (a) 
$h/h_{\infty }$, (b) 
$Y_F$, (c) 
$\nu Y_O$, (d) 
$G$ and (e) 
$E$. Here 
$K =1.0; Pr =1.0; G_{\infty } = 1.0; u_{-\infty } / u_{\infty }= 0.25$. Fuel in higher-speed stream and oxygen in lower-speed stream: blue line. Oxygen in higher-speed stream and fuel in lower-speed stream: red line.
3.3. Multiple flames
If there is a combustible mixture in at least one of the two free streams and at least one reactant in the other stream, it is possible to have multiple flames: a diffusion flame with one or two partially premixed flames. The well-known triple flame exhibits that character. This is well known for triple flames where several flame branches can appear (Jordà Juanós & Sirignano Reference Jordà Juanós and Sirignano2014; Rajamanickam et al. Reference Rajamanickam, Coenen, Sanchez and Williams2019). It also is found in counterflows (López-Cámara et al. Reference López-Cámara, Jordà Juanós and Sirignano2019; Sirignano Reference Sirignano2019a; López-Cámara, Jordà Juanós & Sirignano Reference López-Cámara, Jordà Juanós and Sirignano2020). This understanding leads to the consideration here of situations where a fuel-rich combustible mixture exists in the faster stream at 
$\eta = \infty$ and a fuel-lean combustible mixture exists in the slower stream at 
$\eta = -\infty$. Figure 18 shows results with 
$Da$ varying over several orders of magnitude; specifically 
$K = 0.01, 0.1$ and 
$1.0$.
Figure 18. Multiple-flame solutions for enthalpy, mass fractions and transverse velocity: (a) 
$h/h_{\infty }$, (b) 
$Y_F$, (c) 
$\nu Y_O$ and (d) 
$v g / \mu _{\infty }$. Here 
$K =0.01, 0.1, 1.0; Pr =1.0 ; G_{\infty } = 1.0; u_{-\infty } / u_{\infty }= 0.25$.
Results show three flames: a fuel-lean premixed flame to the left, a diffusion flame in the centre, and a fuel-rich premixed flame to the right. It is noteworthy that the diffusion flame is to the left of the 
$v = 0$ surface because of asymmetry in the reaction-rate law. Excess oxygen from the left diffuses and advects through the fuel-lean premixed flame to reach the diffusion flame. Excess fuel from the right diffuses (but against advection) to reach the diffusion flame. Heat diffuses from the diffusion flame towards both premixed flames.
If the premixed flames behave as classical isolated flames, the speed of propagation should increase with 
$\sqrt {K}$ which does not occur here. Also, propagation out of the mixing layer would occur. In fact, the propagation speed here is rather low and there is evidence of some diffusion control, especially for the fuel-rich flame and somewhat for the fuel-lean flame at low 
$Da$. The dependence on 
$Da$ and the re-location of the fuel-lean flame disallows rigour in the use of the similar solution for this case. Thus, the results in figure 18 are self-contradictory and cannot be regarded as quantitatively accurate. However, they provide qualitative guidance about the behaviour and the ability of three flames to co-exist. To obtain trustworthy quantitative information, the two-dimensional system of (2.8)–(2.12) should be solved.
The above results with two combustible free streams lead to the examination of a configuration with pure oxygen in the slower stream at 
$\eta = -\infty$ and a mixture in the faster stream at 
$\eta = \infty$ which is still more fuel rich and therefore less reactive than the configuration in figure 18. These results are shown in figure 19. The reaction-zone locations for both the diffusion flame and the weaker fuel-rich premixed flame do not change as 
$Da$ is varied over one order of magnitude. The transverse velocity 
$v$ also remains unchanged at the flame locations. The fuel-rich premixed flame here exhibits no wave-like character and is clearly diffusion controlled.
Figure 19. Two-flame solutions for enthalpy and dynamic field variables: (a) 
$h/h_{\infty }$, (b) 
$\omega_{F} x/u_{\infty}$, (c) 
$Y_F$, (d) 
$\nu Y_O$ and (e) 
$v g/\mu_{\infty}$. Here 
$K =1.0; Pr =1.0; G_{\infty } = 1.0; u_{-\infty }/u_{\infty } = 0.25$.
Note that multiple-flame results here agree qualitatively with the counterflow-flame results of Sirignano (Reference Sirignano2019a) and López-Cámara et al. (Reference López-Cámara, Jordà Juanós and Sirignano2019, Reference López-Cámara, Jordà Juanós and Sirignano2020). Sirignano used one-step kinetics for propane and oxygen while Lopez et al. used detailed kinetics with methane–air chemistry.
4. Concluding remarks
A three-dimensional configuration for a mixing layer between two streams with imposed normal strain in the transverse plane has been analysed through one-dimensional similar solutions for both non-reacting and reacting configurations. That is, the two ambient streams have velocity components reflecting both parallel co-flow and counterflow. A similar system of the Navier–Stokes equations coupled with equations for scalar transport is developed and solved. Variable density, temperature and composition are considered. One free stream is pure oxygen or sometimes a fuel-lean combustible mixture while the other stream is pure propane or sometimes a fuel-rich combustible mixture.
Here, a novel feature is that three characteristic rates (and their reciprocals as times) are imposed on the flamelet model: chemical reaction rate, shear-strain rate and normal-strain rate. This can be viewed as an evolution of mixing and combustion studies for shear layers to include imposed normal strain. Or it can be described as an advancement to include imposed shear and vorticity on the common counterflow flamelet model. The results are significant; for example, the shear strain affects the likelihood of multiple-flame branches. The enhancement of the mixing and combustion rates by imposed normal strain on a shear layer can be very substantial. The results show the need for flamelet models with both shear strain and normal strain. Existing models with only normal strain cannot adequately reflect the situation in practical turbulent combustors where shear has obvious importance.
The Damköhler number 
$Da$ is an important parameter in the reacting case which has been varied over a range of one order of magnitude here. In the definition of 
$Da$ here, the reaction rate is normalized using a residence time 
$x/u_{\infty }$. The diffusion flame becomes thinner with increasing 
$Da$ but its position does not change in either the single-flame or multi-flame configurations. The premixed flames show some diffusion control, especially at lower values of 
$Da$ and less flammable inflowing mixture ratios. In some cases, we can expect at least the fuel-rich premixed flame to fit the similar solution format with no change in location as 
$Da$ varies.
The impact of the magnitude of the dimensional shear-strain rate is made solely through the velocity ratio 
$u_{-\infty }/u_{\infty }$. The imposed normal-strain rate is described through the parameter 
$G$. The increase in velocity ratio causes a decrease in the shear-strain rate while an increase in 
$G$ results in an increased normal-strain rate. The two strains have opposing effects on mixing-layer thickness and thereby on scalar gradients, mixing rates, and burning rates. Increased strain rate increases slightly both the displacement of streamlines and the width of the mixing layer. Increases in the normal-strain rate can profoundly decrease the mixing-layer width and increase transport rates.
The velocity profiles and the scalar profiles are shown to depend on the Prandtl number as well as the strain rates; higher 
$Pr$ caused steeper gradients and more narrow thermal and compositional layers. A generalization has been developed for the classical Crocco integral to address non-unitary Prandtl number.
A closed form for the transverse velocity is found to match the similar solution formulation.
For the non-reacting subsonic case, the effect of Mach number 
$M$ is examined. Viscous dissipation is shown to have a modest role in the thermal behaviour.
Flamelet theory as a closure model for turbulent combustion has been based on the tracking of two variables: a normalized conserved scalar and the strain rate; the latter may be given either directly or through a progress variable. Mixture fraction has traditionally been used for the conserved scalar. A new normalized scalar is presented; the variable 
$\varSigma$ can be built on the Shvab–Zel’dovich conserved scalars or it can be a functional of the velocity field.
In the future, studies with more refined and detailed descriptions of chemical kinetics, transport and thermophysical properties should be considered. The information, insights and guidance appearing here should provide some guidance. It must be understood, as explained in § 2.6, that, with the chemical reaction confined to a thin zone, the reduction to a similar solution in one contrived dimension (i.e. 
$\eta$) is reasonable. However, as the situation moves towards an interest in extinction or ignition with a broader flame zone or towards a more thorough examination of the multiple-flame configurations, our three-dimensional problem can, at best, be reduced only to a two-dimensional problem.
Acknowledgements
This research was supported by the Air Force Office of Scientific Research under Grant FA9550-18-1-0392 with Dr M. Birkan as the scientific officer.
Declaration of interest
The author reports no conflict of interest.
