Commit 500882a9 authored by Yu Pan's avatar Yu Pan
Browse files

Explain advection and diffusion term separately

parent ff4ffb9f
......@@ -3,59 +3,101 @@
In this chapter, we walk the reader through our 2D and 3D compressible Navier-Stokes Solver (CompressibleFlowSolver).
\section{Fundamental Theories of CompressibleFlowSolver}
The governing equation systems include continuity equation, momentum equations and energy equation. Write in conservative form
The governing equation systems include continuity equation, momentum equations and energy equation. Written in conservative form
\begin{equation}\label{eq1}
\frac{\partial U}{\partial t}+
\frac{\partial F_{i}(U)}{\partial x_{i}}
=\frac{\partial G_{i}(U,\nabla U)}{\partial x_{i}}
\frac{d \textbf{U}}{d t}+
\frac{\partial \textbf{F}_{i}(U)}{\partial x_{i}}+
\frac{\partial \textbf{G}_{i}(U,\nabla U)}{\partial x_{i}}
=0
\end{equation}
where $U=[\rho,\rho u_{1} \hdots \rho u_{d},E]^{T}$
where $\textbf{U}=[\rho,\rho u_{1} \hdots \rho u_{d},E]^{T}$ is the vector of conservative variables.
Inviscid flux F in i direction is
And the inviscid flux $\textbf{F}$ in i direction is
\begin{equation}
F_{i}=
\textbf{F}_{i}=
\begin{bmatrix}
\rho u_{i}\\
\rho u_{1}u_{i}+\delta_{1,i}P\\
\rho u_{1}u_{i}+\delta_{1,i}p\\
\vdots\\
\rho u_{d}u_{i}+\delta_{d,i}P\\
u_{i}(E+P)
\rho u_{d}u_{i}+\delta_{d,i}p\\
u_{i}(E+p)
\end{bmatrix}
\end{equation}
where total energy $E=\rho (c_{v}T+u_{i}u_{i}/2)$, pressure $P=(\gamma-1)(E-\rho u_{i}u_{i}/2)$
where total energy $E=\rho (c_{v}T+u_{i}u_{i}/2)$, pressure $p=(\gamma-1)(E-\rho u_{i}u_{i}/2)$
Viscous flux G in i direction is
\begin{equation}
G_{i}=
The viscous flux $\textbf{G}$ in i direction is
\begin{equation}\label{eq7}
\textbf{G}_{i}=
\begin{bmatrix}
0\\
\tau_{i1}\\
-\tau_{i1}\\
\vdots\\
\tau_{id}\\
\sum\limits_{j=1}^{d}{u_{j}\tau_{ij}}-q_{i}
-\tau_{id}\\
q_{i}-\sum\limits_{j=1}^{d}{u_{j}\tau_{ij}}
\end{bmatrix}
\end{equation}
where viscous term $\tau_{ij}=\mu(\frac{\partial u_{i}}{\partial x_{j}}+\frac{\partial u_{j}}{\partial x_{i}}-\frac{2}{3}\nabla \cdot\textbf{u}\delta_{ij})$ and heat flux $\textbf{q}=-\kappa \nabla T$
where the viscous tensor $\tau_{ij}=\mu(\frac{\partial u_{i}}{\partial x_{j}}+\frac{\partial u_{j}}{\partial x_{i}}-\frac{2}{3}\frac{\partial u_{i}}{\partial x_{i}}\delta_{ij})$ and heat flux $q_{i}=-\kappa \frac{\partial T_{i}}{\partial x_{i}}$
In Nektar++ current compressible flow solver, DG spatial discretization is in default. Multiply a test function $\phi$, integrate inside an element $K\in \Omega_{h}$ and sum over all the elements.
\begin{equation}
\sum\limits_{K\in \Omega_{h}}{\int_{K}{\frac{\partial U}{\partial t}\phi dx}}+I_{inv}=I_{vis}
In Nektar++, the compressible flow solver spatial discretization is dicountinous Galerkin. Multiply a test function $\phi$, integrate inside an element $K\in \Omega_{h}$, where $\Omega_{h}$ is the total computational region.
\begin{equation}\label{eq4}
\int_{K}{\frac{d \textbf{U}^{e}}{d t}\phi^{e} dx}+I_{inv}+I_{vis}=0
\end{equation}
where
\begin{equation}\label{eq2}
I_{inv}= \sum\limits_{K\in \Omega_{h}}{\int_{\partial K}{\sum\limits_{i=1}^{d}{\widetilde{F}_{i}(U)n_{i}\phi}}dS}-
\sum\limits_{K\in \Omega_{h}}{\int_{K}{\sum\limits_{i=1}^{d}{F_{i}(U)\frac{\partial \phi}{\partial x_{i}}dx}}}
I_{inv}={\int_{\partial K}{\sum\limits_{i=1}^{d}{\widetilde{\textbf{F}^{e}_{i}}(U)n_{i}\phi^{e}}}dS}-
{\int_{K}{\sum\limits_{i=1}^{d}{\textbf{F}^{e}_{i}(U)\frac{\partial \phi^{e}}{\partial x_{i}}dx}}}
\end{equation}
\begin{equation}\label{eq3}
I_{inv}= \sum\limits_{K\in \Omega_{h}}{\int_{\partial K}{\sum\limits_{i=1}^{d}{\widetilde{G}_{i}(U,\nabla U)n_{i}\phi}}dS}-
\sum\limits_{K\in \Omega_{h}}{\int_{K}{\sum\limits_{i=1}^{d}{G_{i}(U,\nabla U)\frac{\partial \phi}{\partial x_{i}}dx}}}
I_{inv}= {\int_{\partial K}{\sum\limits_{i=1}^{d}{\widetilde{\textbf{G}}^{e}_{i}(U,\nabla U)n_{i}\phi^{e}}}dS}-
{\int_{K}{\sum\limits_{i=1}^{d}{\textbf{G}^{e}_{i}(U,\nabla U)\frac{\partial \phi^{e}}{\partial x_{i}}dx}}}
\end{equation}
Different flux terms such as Equation \eqref{eq2} and \eqref{eq3} and source terms are treated separated and transferred through member variables as the sketch \ref{fig5} shows.
Various advection flux $\widetilde{\textbf{F}}$ and diffusion flux $\widetilde{\textbf{G}}(U,\nabla U)$ are supported in Nektar++ CompressibleSolver. And each specific type of flux such as LDG (DiffusionLDGNS) and Interior Penalty (DiffusionIP) diffusion fluxes inherit from their general flux class (Diffison).
\paragraph{Advection flux $\textbf{F}_{i}$: }
the following section takes AdvectionWeakDG.cpp as an example to demonstrate how the advection flux is calculated in Nektar++. Some matrices index can refer to \cite{KaSh05}
Separate the advection flux from Equation \eqref{eq4}.
\begin{equation}
\int_{K}{\frac{d \textbf{U}^{e}}{d t}\phi^{e} dx}+I_{inv}=0
\end{equation}
Different flux terms such as Equation \eqref{eq2} and \eqref{eq3} and source terms are treated separated and transferred through member variables as the sketch \ref{fig1} shows.
Various advection flux $\widetilde{F}$ and diffsion flux $\widetilde{G}(U,\nabla U)$ are supported in Nektar++ CompressibleSolver. And each specific flux such as LDG diffusion flux (LDGNS) inherits from its general flux type (Diffison).
Transform Equation \eqref{eq4} into elemental matrix form.
\begin{equation}\label{eq5}
(B^{e})^{T}W^{e}B^{e}\frac{d \hat{\textbf{U}}^{e}}{d t}
-\sum\limits_{i=1}^{d}{(D^{e}_{x_{i}} B^{e})^{T}}W^{e}B^{e}\hat{\textbf{F}}^{e}_{i}+
b^{e}
=0
\end{equation}
where the integration of elemental trace flux $b^{e}=\sum\limits_{i=1}^{d}{B^{e}\widetilde{\textbf{F}}^{e}_{i}}n_{i}$ and the Riemann flux $\widetilde{\textbf{F}_{i}}=\textbf{F}^{e}_{i}(U^{+},U^{-})$.
Other matrix index explanations are as Table \ref{table2}
\begin {table}[H]
\caption {Item explanation} \label{table2}
\begin{center}
\scalebox{0.9}[1.]{
\begin{tabular}{| c|c|c|c|}
\hline
$\phi^{e}$ & $\frac{\partial \phi^{e}}{\partial x_{i}}$& Diagonal weight matrix & Coefficients\\
\hline
$B^{e}$ &$(D^{e}_{x_{i}} B^{e})^{T}$& $W^{e}$& $\hat{U}^{e}$,$\hat{F}^{e}$\\
\hline
\end{tabular}
}
\end{center}
\end{table}
Write Equation \eqref{eq5} in semi-discretization form.
\begin{equation}\label{eq8}
\frac{d \hat{\textbf{U}}^{e}}{d t}
=((B^{e})^{T}W^{e}B^{e})^{-1}(\sum\limits_{i=1}^{d}{(D^{e}_{x_{i}} B^{e})^{T}}W^{e}B^{e}\hat{\textbf{F}}^{e}_{i}-
b^{e})
\end{equation}
where mass matrix $M=(B^{e})^{T}W^{e}B^{e}$.
Usually inside one specific flux type, the codes includes the following simliar process to caculate the flux. Take AdvectionWeakDG.cpp as an example.
\begin {table}[!h]
Table \ref{table1} lists some functions that the process to calculate the advection term.
\begin {table}[h]
\caption {AdvectionWeakDG.cpp} \label{table1}
\begin{center}
\scalebox{0.9}[1.]{
......@@ -63,24 +105,110 @@ Usually inside one specific flux type, the codes includes the following simliar
\hline
Variable/Function name & Physical meaning \\
\hline
m$\_$fluxVector & Advection volume flux: $F_{i}$\\
fluxVector & Advection volume flux: $\textbf{F}_{i}$\\
\hline
numflux & Advection numerical flux at trace: $\widetilde{F}_{i}$\\
numflux & Advection (Riemann) numerical flux at trace: $\widetilde{\textbf{F}}_{i}$\\
\hline
IProductWRTDerivBase & Volume flux integration: $\sum\limits_{K\in \Omega_{h}}{\int_{K}{\sum\limits_{i=1}^{d}{F_{i}(U)\frac{\partial \phi}{\partial x_{i}}dx}}}$ \\
IProductWRTDerivBase & Volume flux integration: $\sum\limits_{i=1}^{d}{(D^{e}_{x_{i}} B^{e})^{T}}W^{e}B^{e}\hat{\textbf{F}}^{e}_{i}$ \\
\hline
AddTraceIntegral & Add Surface flux into volume flux interation: $+\sum\limits_{K\in \Omega_{h}}{\int_{\partial K}{\sum\limits_{i=1}^{d}{\widetilde{F}_{i}(U)n_{i}\phi}}dS}$\\
AddTraceIntegral & Add surface flux into volume flux interation: $+b^{e}$\\
\hline
MultiplyByElmImvMass& Multiply the inverse of mass matrix to get the flux coefficients\\
MultiplyByElmImvMass& Multiply the inverse of mass matrix to get the flux coefficients $ \hat{\textbf{U}}$\\
\hline
\end{tabular}
}
\end{center}
\end{table}
\paragraph{Diffusion flux $\textbf{G}_{i}$:}
the following section takes DiffusionLDGNS.cpp as an example to demonstrate how the diffusion flux is calculated in Nektar++. Some matrices index can refer to \cite{KaSh05}.
Separate the diffusion flux from Equation \eqref{eq4}.
\begin{equation}\label{eq6}
\int_{K}{\frac{d \textbf{U}^{e}}{d t}\phi^{e} dx}+I_{vis}=0
\end{equation}
As the resolution of diffusion flux is to solve a $2^{nd}$ PDE, it needs two equations to resolve the problem. LDG separates the resolution into two $1^{st}$ PDE.
From Equation \eqref{eq7}, diffusion derivatives $\textbf{G}_{i}(P,\nabla P)$ so that can reduce the computation using primitive variables $\textbf{P}=[u_{1}\hdots u_{d}, T]^{T}$. Therefore, the first step is to calculate the derivatives $\textbf{Q}_{i}=\frac{\partial P}{\partial x_{i}}$
\begin{equation}\label{eq9}
\int_{K}{\textbf{Q}^{e}_{i} \phi^{e} dx}=
{\int_{\partial K}{\widetilde{\textbf{P}}^{e}n_{i}\phi^{e}}dS}-
\int_{K}{\textbf{P}^{e}\frac{\partial \phi^{e}}{\partial x_{i}}dx}
\end{equation}
Transform Equation \eqref{eq9} into elemental matrix form.
\begin{equation}\label{eq10}
\hat{\textbf{Q}}_{i}^{e}
=((B^{e})^{T}W^{e}B^{e})^{-1}
((D^{e}_{x_{i}} B^{e})^{T}W^{e}B^{e}\hat{\textbf{P}}^{e}-
b_{1,i}^{e})
\end{equation}
where the integration of elemental trace flux $b_{1,i}^{e}=B^{e}\widetilde{\textbf{P}}n_{i}$ and the LDG numerical flux $\widetilde{\textbf{P}}=\textbf{P}^{e}_{i}(P^{+},P^{-})$.
When we get the diffusion derivatives $\textbf{G}_{i}$, the process to solve Equation \eqref{eq6} is simliar to the process of Equation \eqref{eq8}.
\begin{equation}\label{eq8}
\frac{d \hat{\textbf{U}}^{e}}{d t}
=((B^{e})^{T}W^{e}B^{e})^{-1}(\sum\limits_{i=1}^{d}{(D^{e}_{x_{i}} B^{e})^{T}}W^{e}B^{e}\hat{\textbf{G}}^{e}_{i}-
b_{2}^{e})
\end{equation}
where the integration of elemental trace flux $b_{2}^{e}=\sum\limits_{i=1}^{d}{B^{e}\widetilde{\textbf{G}}^{e}_{i}}n_{i}$ and the LDG numerical flux $\widetilde{\textbf{G}_{i}}=\textbf{G}^{e}_{i}(G^{+},G^{-})$.
LDG numerical fluxes in Equations \eqref{eq9} and \eqref{eq8} satisfy
\begin{equation}
\tilde{\textbf{P}}=\{\textbf{P}\}+C_{12}\cdot [\textbf{P} n]
\end{equation}
where $C_{12}$ is a coefficient vector.
\begin{equation}
\tilde{\textbf{G}}=\{\textbf{G}\}-C_{12} [\textbf{G} \cdot n]
\end{equation}
where $[\quad]$ represents the jump while $\{\quad \}$ represents the average. In Nektar++, $C_{12}=\frac{1}{2}$, $C_{11}$=0
Table \ref{table3} lists some functions that the process to calculate the diffusion term.
\begin {table}[h]
\caption {DiffusionLDGNS.cpp} \label{table3}
\begin{center}
\scalebox{0.9}[1.]{
\begin{tabular}{ | c | c|}
\hline
Variable/Function name & Physical meaning \\
\hline
v\_NumericalFluxO1 & Numerical flux $\tilde{\textbf{P}}^{e}$ \\
\hline
First IProductWRTDerivBase & Volume flux $(D^{e}_{x_{i}} B^{e})^{T}W^{e}B^{e}\hat{\textbf{P}}^{e}$\\
\hline
First AddTraceIntegral& Add numerical flux integration $+b_{1,i}^{e}$ \\
\hline
MultiplyByElmImvMass& Multiply the inverse of mass matrix to get the flux coefficients $\textbf{Q}$\\
\hline
m\_fluxVectorNS & Diffusion volume flux: $\textbf{G}_{i}$\\
\hline
v\_NumericalFluxO2 & Diffusion numerical flux : $\widetilde{\textbf{G}}_{i}$\\
\hline
Second IProductWRTDerivBase & Volume flux integration: $\sum\limits_{i=1}^{d}{(D^{e}_{x_{i}} B^{e})^{T}}W^{e}B^{e}\hat{\textbf{G}}^{e}_{i}$ \\
\hline
Second AddTraceIntegral & Add surface flux into volume flux interation: $+b_{2}^{e}$\\
\hline
MultiplyByElmImvMass& Multiply the inverse of mass matrix to get the flux coefficients $\hat{\textbf{U}}$\\
\hline
\end{tabular}
}
\end{center}
\end{table}
\clearpage
\section{Data Structure of CompressibleFlowSolver}
\begin{figure}
\caption{CompressibleFlowSolver DataStructure}
\caption{CompressibleFlowSolver DataStructure}\label{fig5}
\centering
\begin{turn}{90}
\includestandalone[width=1.2\linewidth]{DataStructure}
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment