\chapter{Simulations}

\begin{figure}[ht]%
\centering
\includegraphics[height=0.4\textheight]{DLRflames2.jpg}%
\caption{Image of a simple jet flame in the center with 22.1\% $CH_4$, 33.2\% $H_2$ and 44.7\% $N_2$. The jet velocity is \unit{42.1}{\metre \per \second} and $Re=15200$. The picture to the left shows Rayleigh-scattering images of the flame. The image to the right is from 2D-PLIF measurements.}%
\label{fig:DLRflames}%
\end{figure}

\section{Scenario description}
The TNF experimental data archive \cite{TNFdata} contains data from measurements on several different flame configurations. Most configurations are for a mixture of fuels, but in this thesis the Sandia $H_2$-A flame was chosen since it has the simplest fuel composition with 100 \%  $H_2$. The measurements of this flame where made at Sandia National Laboratories and the setup is documented in \cite{SandiaH2He},\cite{Barlow1994} and \cite{Barlow1996}. The fuel jet has a nozzle diameter of \unit{3.75}{\milli \metre} and fuel jet velocity of \unit{296}{\metre \per \second} which gives $Re=10000$ . The uncertainties are listed in \cite{SandiaH2He}, in the interval $0.5 < \xi / \xi_{st} < 2$, to be $\pm 3-4 \%$ for $T$, $\xi$, $N_2$ and $O_2$ mole and mass fractions. The uncertainty of measurement for $Y_{OH}$ is $\pm 15 \%$.

\begin{figure}%
\subfloat[Sketch of measurement locations]{
\includegraphics[height=0.5\textheight]{flamemeasure.png}%
\label{fig:Flamemeasure}%
}
\subfloat[Computational grid]{
\includegraphics[height=0.5\textheight]{initial_grid.jpg}%
\label{fig:Initial grid}%
}
\caption{\subref{fig:Flamemeasure} Sketch of the locations of the point measurements. Point measurements of $\xi$, $T$ and $Y_k$ where made at 7 different distances from the fuel outlet  \subref{fig:Initial grid}Initial grid of measuring points in the flame used for computations.}%
\label{fig:Flamesetup}
\end{figure}

Point measurements from the Sandia H2-A flame where curve fitted for easy use in the simulations. The experimental data for this flame is available from \cite{TNFdata}. The data was used to find $T(\xi)$, $u(\xi)$ and $Y_k(\xi)$. Simultaneous single shot data of the mass fractions of $H_2$,$O_2$,$H_2O$, $N_2$, $OH$, $NO$ where made with Laser Raman and Rayleigh scattering. LDV\footnote{Laser Doppler Velocimetry}-measurements where made of the velocity profile at ETH-Zürich \cite{flurylaser}. \ For the sake of clarity the $x$ denotes the position in the axial or streamwise direction of the reacting flow and $r$ denotes the radial direction.

The data set contains both Favre-averaged and Reynolds-averaged measured values of the temperature, mixture fraction and mass fractions.
The rms of the mixture fraction is defined as

\begin{equation}
\xi_{rms}=\sqrt{\widetilde{\xi''^2}}
\label{eq:rms}
\end{equation}

\begin{table}[ht]%
\caption{Sandia $H_2$-flame data \cite{SandiaH2He}}
\begin{tabular}{|l |c| r |}
\hline
Property & Value	& Uncertainty of measurement \\ \hline
Fuel 	& $100 \% H_2$ & \\
Jet velocity $U_{j}$ & $\unit{296}{\metre \per \second}$ & $\pm 1.5\%$  \\
Inner diameter of fuel nozzle, $d$	& \unit{3.75}{\milli \metre} &\\
$Re_d$ 	& $10000$ &\\
Coflow velocity 	& $\unit{1}{\metre \per \second}$ & $\pm 1.3 \%$ \\
Coflow temperature	& $295 K$& $\pm 2 K$	\\
Fuel temperature	& $296 K$ & $\pm 2 K$\\
Stoichiometric mixture fraction, $\xi_{st}$ & 0.0285 & \\
Visible flame length,$L$ & $\approx 180 D (\unit{675}{\milli \metre})$ & \\

\hline
\end{tabular}
\label{tab:H2data}
\end{table}

\section{Assumptions}

\begin{enumerate}
	

\item The irreversibility contribution from viscous dissipation of the flow has not been calculated since this contribution was found to be negligible in the literature (\cite{Nishida2002},\cite{Stanciu2001}). 

\item The axial variation in mixture fraction, $\xi(x)$ is assumed to be the same for the whole radial cross section of the flame. 

\item The chosen data set is for the Favre-averaged properties of the flow. The root mean-square (rms) values of each property is also given in the data set.

\item The pressure gradient in the flame is assumed to be negligible compared to the concentration gradient shown in Eq. \ref{eq:Y_largerthan_p}.

\item The scalar dissipation rate $\chi$ is assumed independent of the other variables.

\item The mean entropy production is computed in 2 dimensions.

\item In Eq. \ref{eq:FastTmixfrac} $Q$ is the lower heating value of hydrogen

\item When using the fast chemistry relations, only the mass fractions of $H_2$,$O_2$,$N_2$ and $H_2O$ where considered. For the curvefits the measured mass fraction of $OH$ was accounted for in the entropy production term from mass transfer Eq. \ref{eq:massentropy}, though since the highest measured mass fraction of $OH$ was $0.004$, its contribution was to total was minimal. 

\item $NO$-formation in the flame has not been considered since the mass fraction of $NO$ is $\mathcal{O}(10^{-4})$ in the part of the highest temperature region of the flame. This is where mass fraction of $NO$ reaches its maximum since the $N_2+O_2->NO$ reaction is highly temperature dependent.

\item The conductivity of the mixture is given by Eq. \ref{eq:conductivity of mixture}.

\item In the process of finding the scalar dissipation rate Eq. \ref{eq:Scalar dissipation} the following assumption was made 


\begin{equation}
\left( \frac{\partial \tilde{\xi}}{\partial x_j} \right)^2= \left(\frac{\partial \tilde{\xi}}{\partial x} \right)^2 + \left(\frac{\partial \tilde{\xi}}{\partial r} \right)^2 + \underbrace{\left(\frac{1}{r} \frac{\partial \tilde{\xi}}{\partial \theta} \right)^2}_{=0}.
\label{eq:nothetavar}
\end{equation}

\item The mixture is regarded as a mixture of ideal gases. By looking at a generalized compressibility chart, like Fig A-1 in \cite{MS}, this assumption can be checked. Hydrogen has a critical pressure of $13.0$ bar and critical temperature of $33.2$ K \cite[Table A-1]{MS}, and from the chart the ideal gas assumption seems to hold since $p/p_c$ is low and $T/T_c$ is high. 

\item The diffusion coefficient in the scalar dissipation model, Eq. \ref{eq:Scalar dissipation} and the relation for entropy production from heat transfer, Eq. \ref{eq:Qentropy} is the turbulent diffusion coefficient. An approximation for this must be made from a turbulence model. A mixing length model can be used. $nu_T$ is the value that is estimated. For a turbulent jet $Sc_T$ and $Pr_T$ are set to be 0.7. The Lewis number is assumed to be unity ($Le_T=\frac{Sc_T}{Pr_T}=1$ ). 

\item Using Prandtl's mixing length model one obtains the following relation for the turbulent exchange coefficient,$\nu_T$,
\begin{equation}
\nu_T=l^2 \left| \frac{d \overline{u}}{dr} \right|.
\label{eq:turbulent visc}
\end{equation}
Here $l$ is the mixing length

A simple zero-equation turbulence model was used in the computations. In a turbulent shear flow, $l$ is a function of the thickness, $\delta$, of the shear layer. Following the methodology in \cite[p.207]{warnatz2006combustion}, in a turbulent round jet in stagnant surroundings $\delta=0.085x$. The mixing length is then found as

\begin{equation}
l=\alpha \delta  \text{ for }  0 \le r \le \delta.
\label{eq:mixing length}
\end{equation}

where $\alpha=0.075$ determined from a variety of experiments for typical conditions. This correlation is taken from Launder and Spalding (1972)\cite{launder1972lectures} where they interpret $\delta$ as the distance from the center axis to the point where the fluid velocity equals 1 \% of the maximum velocity differencde across the shear layer.
\end{enumerate}

The balanced reaction for the combustion of hydrogen is

\begin{equation}
H_2+\frac{1}{\Phi}\frac{1}{2}(O_2+3.76 N_2) \rightarrow H_2O + 0.5\cdot 3.76 N_2
\label{eq:H2combustion}
\end{equation}


Since the temperature and mass fractions are given as piecewise continuous functions their derivative must also be prescribed for the fast chemistry assumption. The method of \cite{Stanciu2001} was used to define these derivatives.

\section{Curve fitting}
\label{sec:curve fitting}

The experimental data was curve fitted as a function of $r$ and $\xi$ in Matlab R2010a using the ``Basic fitting'' plot tool. This tool allows the user to choose a  The best fit for the property was chosen from a visual evaluation of the curve fitting polynomial and the its norm of residuals. In retrospect it would be advantageous to write a script for choosing the best suited curvefit for the data, since the process of doing this in the Matlab toolbox was quite time consuming. 

The curve fits for the Sandia $H_2$-A flame are listed in Appendix \ref{sec:curvefits}.

\section{Using the prescribed pdf method}


With the prescribed pdf method the flame is regarded as a distribution of mean mixture fractions, $\overline{\xi}$, and variances, $\overline{\xi'^2}$.
The calculations where made in Maple 14.00. Maple was chosen for the computations because it is powerful in solving symbolic functions. The Maple routine for calculating the entropy production terms from heat transfer, mass transfer and chemical reaction is attached in Appendix \ref{sec:mapleroutine}.

The values of $\overline{\xi}$ and $\overline{\xi}_{rms}$ where quite low in the flame. The mixture fraction is related to the mass fraction of fuel, and since the molar mass of $H_2$ is a lot lower than the molar mass of air this gives low mass fractions. which is relatively low because of the gave high exponents $a$ and $b$ in for the beta probability density function Eqs. \ref{eq:beta-func pdf}-\ref{eq:beta parameters}. This subsequently resulted in large intricate polynomials which could not be integrated analytically, and therefore they had to be solved by numerical integration. The pdf for two of the cells could not be found since the coefficients $a<0$ and $b<0$.

\section{Description of calculation routine}
A Maple routine was written 

In the experimental data the axial variation of the mixture fraction was only measured at the center of the flame. These measurements where curve fitted using the above mentioned method, .  The term $\tilde{\xi}(x)$ and its derivative $\\partial \tilde{\xi}/\partial x$ where considered to be the same outside the center of the flame. $\frac{\partial \tilde{\xi}}{\partial r}>>\partial \xi/\partial x$, therefore this simplification does not affect the end result significantly.  


\begin{figure}%
\centering
\includegraphics[height=0.2\textheight]{xi(x).jpg}%
\caption{Axial variation of $\tilde{\xi}$ at the centerline}%
\label{fig:xi(x)}%
\end{figure}

Thermodynamic data for the heat capacity, enthalpy and entropy of each species was calculated with coefficients from \cite{kee1990chemkin}.

The density of the mixture is modeled using the ideal gas law. Relating the density to mass fractions in the flame should be a better correlation. Another possible adjustment to the same ideal gas model is 

\begin{equation}
\rho(\xi)=(\xi \rho_2 T_2+(1-\xi) \rho_1 T_1)\frac{1}{T(\xi)}.
\label{eq:rho(xi)}
\end{equation}

$\rho_1$ is the density of the fuel stream and $\rho_2$ is the density of the oxidizer stream.

The released heat from \emph{equation reference to fastchem T and dT} was set to be the value of the lower heating value of hydrogen ($LHV_{H_2}=\unit{119950}{\kilo \joule \per \kilogram}$ in table A.21 in \cite{MS})

The following boundary conditions where imposed on the chemical reaction term

\begin{eqnarray}
\text{As } X_k \rightarrow 0 & X_k ln(X_k)\rightarrow 0 \\
\text{As } X_k \rightarrow 1 & X_k ln(X_k) \rightarrow 1 \\
\label{eq:chementropy BCs}
\end{eqnarray}

From \cite{Ertesvċg2004} the field of local entropy production can be calculated using the following relation

\begin{equation}
\dot{I}=T_{0} \int_V{\overline{\rho}\tilde{\sigma}dV}=T_{0} \int{}\int{}\int{\dot{\overline{S}}'''_{gen} r dr dx d\theta}
\label{eq:Total volume entropy}
\end{equation}


where $\overline{\rho} \tilde{\sigma}=\overline{\rho \sigma}$ denotes the sum of the production terms in the entropy equation. $dV=r dr dx d\theta$ where $r \in [0,r_{max,i}], x \in [0,L]$ and $\theta \in [0, 2 \pi]$. $r_{max,i}$ is the radial distance from the centerline of the outer measuring point of the $i$'th measured sheet of the flame. The flame is assumed to be symmetrical about the centerline and that there are no gradients in the $\theta$-direction.

\begin{equation}
\dot{I}=\sum_i{\sum_j{\dot{\overline{S}}'''_{gen} r_{ij} (\Delta r)_j (\Delta x)_i }}
\label{eq:cyl shell sum}
\end{equation}
The computation script that does the volume integral of the flame for sum of entropy is attached in Appendix C. 

Figure \ref{fig:Initial grid} shows the grid The data set contains measurements from a total of 72 points in the flame distributed over 7 different lengths from the fuel outlet. 

The calculated results where plotted in Matlab.

\begin{figure}%
\centering
\includegraphics[width=\columnwidth]{flameshell.jpg}%
\caption{Outer shell of the flame when Fig. \ref{fig:Initial grid} is rotated $\unit{360}{\degree}$ around the x-axis}%
\label{fig:flameshell}%
\end{figure}

Figure \ref{fig:flameshell} shows the outer shell of the computed flame when Fig. \ref{fig:Initial grid} is rotated around the x-axis.