Formulations (Under Development)
The governing equations of the hydro-acoustic wave include the continuity
equation
(1)\[ \begin{align}\begin{aligned}\newcommand{\dpd}[3][]{\mathinner{
\dfrac{\partial{^{#1}}#2}{\partial{#3^{#1}}}
}}\\\dpd{\rho}{t} + \sum_{i=1}^3 \dpd{\rho v_i}{x_i} = 0\end{aligned}\end{align} \]
and the momentum equations
(2)\[\dpd{\rho v_j}{t}
+ \sum_{i=1}^3 \dpd{(\rho v_iv_j + \delta_{ij}p)}{x_i}
= \dpd{}{x_j}\left(\lambda\sum_{k=1}^3\dpd{v_k}{x_k}\right)
+ \sum_{i=1}^3 \dpd{}{x_i}
\left[\mu \left( \dpd{v_i}{x_j}+\dpd{v_j}{x_i} \right)\right],
\quad j = 1, 2, 3\]
where \(\rho\) is the density, \(v_1, v_2,\) and \(v_3\) the
Cartesian component of the velocity, \(p\) the pressure,
\(\delta_{ij}, i, j = 1, 2, 3\) the Kronecker delta, \(\lambda\) the
second viscosity coefficien, \(\mu\) the dynamic viscosity coefficient,
\(t\) the time, and \(x_1, x_2\), and \(x_3\) the Cartesian
coordinate axes. Newtonian fluid is assumed.
The above four equations in (1) and (2) have five
independent variables \(\rho, p, v_1, v_2\), and \(v_3\), and hence are
not closed without a constitutive relation. In the bulk
package,
the constitutive relation (or the equation of state) of choice is
\[K = \rho\dpd{p}{\rho}\]
where \(K\) is a constant and the bulk modulus. We chose to use the
density \(\rho\) as the independent variable, and integrate the equation of
state to be
(3)\[p = p_0 + K \ln\frac{\rho}{\rho_0}\]
where \(p_0\) and \(\rho_0\) are constants. Substituting (3)
into (2) gives
(4)\[\dpd{\rho v_j}{t} + \sum_{i=1}^3\dpd{}{x_i}
\left[\rho v_iv_j
+ \delta_{ij}\left(p_0 + K\ln\frac{\rho}{\rho_0}\right) \right]
= \sum_{i=1}^3 \dpd{}{x_i}
\left[\delta_{ij} \lambda \sum_{k=1}^3 \dpd{v_k}{x_k}
+ \mu \left( \dpd{v_i}{x_j}+\dpd{v_j}{x_i} \right)\right],
\quad j = 1, 2, 3\]
Jacobian Matrices
We proceed to analyze the advective part of the governing equations
(1) and (4). Define the conservation variables
(5)\[ \begin{align}\begin{aligned}\newcommand{\bvec}[1]{\mathbf{#1}}
\newcommand{\defeq}{\buildrel{\text{def}}\over{=}}\\\begin{split}\bvec{u} \defeq \left(\begin{array}{c}
\rho \\ \rho v_1 \\ \rho v_2 \\ \rho v_3
\end{array}\right)\end{split}\end{aligned}\end{align} \]
and flux functions
(6)\[\begin{split}\bvec{f}^{(1)} \defeq \left(\begin{array}{c}
\rho v_1 \\
\rho v_1^2 + K\ln\frac{\rho}{\rho_0} + p_0 \\
\rho v_1v_2 \\ \rho v_1v_3
\end{array}\right), \quad
\bvec{f}^{(2)} \defeq \left(\begin{array}{c}
\rho v_2 \\ \rho v_1v_2 \\
\rho v_2^2 + K\ln\frac{\rho}{\rho_0} + p_0 \\
\rho v_2v_3
\end{array}\right), \quad
\bvec{f}^{(3)} \defeq \left(\begin{array}{c}
\rho v_3 \\ \rho v_1v_3 \\ \rho v_2v_3 \\
\rho v_3^2 + K\ln\frac{\rho}{\rho_0} + p_0
\end{array}\right)\end{split}\]
Aided by the definition of conservation variables in (5), the flux
functions defined in (6) can be rewritten with \(u_1, u_2, u_3\),
and \(u_4\)
(7)\[\begin{split}\bvec{f}^{(1)} = \left(\begin{array}{c}
u_2 \\
\frac{u_2^2}{u_1} + K\ln\frac{u_1}{\rho_0} + p_0 \\
\frac{u_2u_3}{u_1} \\ \frac{u_2u_4}{u_1}
\end{array}\right), \quad
\bvec{f}^{(2)} = \left(\begin{array}{c}
u_3 \\ \frac{u_2u_3}{u_1} \\
\frac{u_3^2}{u_1} + K\ln\frac{u_1}{\rho_0} + p_0 \\
\frac{u_3u_4}{u_1}
\end{array}\right), \quad
\bvec{f}^{(3)} = \left(\begin{array}{c}
u_4 \\ \frac{u_2u_4}{u_1} \\ \frac{u_3u_4}{u_1} \\
\frac{u_4^2}{u_1} + K\ln\frac{u_1}{\rho_0} + p_0
\end{array}\right)\end{split}\]
By using (5) and (7), the left-hand-side of the governing
equations can be cast into the conservative form
(8)\[\dpd{\bvec{u}}{t} + \sum_{i=1}^3\dpd{\bvec{f}^{(i)}}{x_i} = 0\]
Aided by using the chain rule, (8) can be rewritten in the
quasi-linear form
(9)\[\dpd{\bvec{u}}{t} + \sum_{i=1}^3\mathrm{A}^{(i)}\dpd{\bvec{u}}{x_i} = 0\]
where the Jacobian matrices \(\mathrm{A}^{(1)}, \mathrm{A}^{(2)}\), and
\(\mathrm{A}^{(3)}\) are defined as
(10)\[ \begin{align}\begin{aligned}\newcommand{\pd}[3][]{
\tfrac{\partial{^{#1}}#2}{\partial{#3^{#1}}}
}\\\begin{split}\mathrm{A}^{(i)} \defeq \left(\begin{array}{cccc}
\pd{f_1^{(i)}}{u_1} & \pd{f_1^{(i)}}{u_2} &
\pd{f_1^{(i)}}{u_3} & \pd{f_1^{(i)}}{u_4} \\
\pd{f_2^{(i)}}{u_1} & \pd{f_2^{(i)}}{u_2} &
\pd{f_2^{(i)}}{u_3} & \pd{f_2^{(i)}}{u_4} \\
\pd{f_3^{(i)}}{u_1} & \pd{f_3^{(i)}}{u_2} &
\pd{f_3^{(i)}}{u_3} & \pd{f_3^{(i)}}{u_4} \\
\pd{f_4^{(i)}}{u_1} & \pd{f_4^{(i)}}{u_2} &
\pd{f_4^{(i)}}{u_3} & \pd{f_4^{(i)}}{u_4}
\end{array}\right), \quad i = 1, 2, 3\end{split}\end{aligned}\end{align} \]
Aided by using (7), the Jacobian matrices defined in (10) can
be written out as
(11)\[\begin{split}\mathrm{A}^{(1)} &= \left(\begin{array}{cccc}
0 & 1 & 0 & 0 \\
-\frac{u_2^2}{u_1^2} + \frac{K}{u_1} & 2\frac{u_2}{u_1} & 0 & 0 \\
-\frac{u_2u_3}{u_1^2} & \frac{u_3}{u_1} & \frac{u_2}{u_1} & 0 \\
-\frac{u_2u_4}{u_1^2} & \frac{u_4}{u_1} & 0 & \frac{u_2}{u_1}
\end{array}\right), \quad
\mathrm{A}^{(2)} = \left(\begin{array}{cccc}
0 & 0 & 1 & 0 \\
-\frac{u_2u_3}{u_1^2} & \frac{u_3}{u_1} & \frac{u_2}{u_1} & 0 \\
-\frac{u_3^2}{u_1^2} + \frac{K}{u_1} & 0 & 2\frac{u_3}{u_1} & 0 \\
-\frac{u_3u_4}{u_1^2} & 0 & \frac{u_4}{u_1} & \frac{u_3}{u_1}
\end{array}\right), \\
\mathrm{A}^{(3)} &= \left(\begin{array}{cccc}
0 & 0 & 0 & 1 \\
-\frac{u_2u_4}{u_1^2} & \frac{u_4}{u_1} & 0 & \frac{u_2}{u_1} \\
-\frac{u_3u_4}{u_1^2} & 0 & \frac{u_4}{u_1} & \frac{u_3}{u_1} \\
-\frac{u_4^2}{u_1^2} & 0 & 0 & 2\frac{u_4}{u_1}
\end{array}\right)\end{split}\]
Hyperbolicity
Hyperbolicity is a prerequisite for us to apply the space-time CESE method to a
system of first-order PDEs. For the governing equations, (1) and
(4), to be hyperbolic, the linear combination of the three
Jacobian matrices of their quasi-linear for must be diagonalizable. The
spectrum of the linear combination must be all real, too [Warming75],
[Chen12].
To facilitate the analysis, we chose to use the non-conservative version of the
governing equations. Define the non-conservative variables
(12)\[\begin{split}\tilde{\bvec{u}} \defeq \left(\begin{array}{c}
\rho \\ v_1 \\ v_2 \\ v_3
\end{array}\right) =
\left(\begin{array}{c}
u_1 \\ \frac{u_2}{u_1} \\ \frac{u_3}{u_1} \\ \frac{u_4}{u_1}
\end{array}\right)\end{split}\]
Aided by using (12) and (5), the transformation between the
conservative variables and the non-conservative variables can be done with the
transformation Jacobian defined as
(13)\[\begin{split}\mathrm{P} \defeq \dpd{\tilde{\bvec{u}}}{\bvec{u}} =
\left(\begin{array}{cccc}
\pd{\tilde{u}_1}{u_1} & \pd{\tilde{u}_1}{u_2} &
\pd{\tilde{u}_1}{u_3} & \pd{\tilde{u}_1}{u_4} \\
\pd{\tilde{u}_2}{u_1} & \pd{\tilde{u}_2}{u_2} &
\pd{\tilde{u}_2}{u_3} & \pd{\tilde{u}_2}{u_4} \\
\pd{\tilde{u}_3}{u_1} & \pd{\tilde{u}_3}{u_2} &
\pd{\tilde{u}_3}{u_3} & \pd{\tilde{u}_3}{u_4} \\
\pd{\tilde{u}_4}{u_1} & \pd{\tilde{u}_4}{u_2} &
\pd{\tilde{u}_4}{u_3} & \pd{\tilde{u}_4}{u_4}
\end{array}\right) = \left(\begin{array}{cccc}
1 & 0 & 0 & 0 \\
-\frac{u_2}{u_1^2} & \frac{1}{u_1} & 0 & 0 \\
-\frac{u_3}{u_1^2} & 0 & \frac{1}{u_1} & 0 \\
-\frac{u_4}{u_1^2} & 0 & 0 & \frac{1}{u_1}
\end{array}\right) = \left(\begin{array}{cccc}
1 & 0 & 0 & 0 \\
-\frac{v_1}{\rho} & \frac{1}{\rho} & 0 & 0 \\
-\frac{v_2}{\rho} & 0 & \frac{1}{\rho} & 0 \\
-\frac{v_3}{\rho} & 0 & 0 & \frac{1}{\rho}
\end{array}\right)\end{split}\]
Aided by writing (5) as
\[\begin{split}\bvec{u} = \left(\begin{array}{c}
\tilde{u}_1 \\
\tilde{u}_1\tilde{u}_2 \\ \tilde{u}_1\tilde{u}_3 \\ \tilde{u}_1\tilde{u}_4
\end{array}\right)\end{split}\]
the inverse matrix of \(\mathrm{P}\) can be obtained
(14)\[\begin{split}\mathrm{P}^{-1} = \dpd{\bvec{u}}{\tilde{\bvec{u}}} =
\left(\begin{array}{cccc}
1 & 0 & 0 & 0 \\
\tilde{u}_2 & \tilde{u}_1 & 0 & 0 \\
\tilde{u}_3 & 0 & \tilde{u}_1 & 0 \\
\tilde{u}_4 & 0 & 0 & \tilde{u}_1
\end{array}\right) = \left(\begin{array}{cccc}
1 & 0 & 0 & 0 \\
v_1 & \rho & 0 & 0 \\
v_2 & 0 & \rho & 0 \\
v_3 & 0 & 0 & \rho
\end{array}\right)\end{split}\]
and \(\mathrm{P}^{-1}\mathrm{P} = \mathrm{P}\mathrm{P}^{-1} =
\mathrm{I}_{4\times4}\).
The transformation matrix \(\mathrm{P}\) can be used to cast the
conservative equations, (9), into non-conservative ones.
Pre-multiplying \(\pd{\tilde{\bvec{u}}}{\bvec{u}}\) to (9) gives
(15)\[\dpd{\tilde{\bvec{u}}}{t} +
\sum_{i=1}^3\tilde{\mathrm{A}}^{(i)}\dpd{\tilde{\bvec{u}}}{x_i} = 0\]
where
(16)\[\tilde{\mathrm{A}}^{(1)} \defeq
\mathrm{P}\mathrm{A}^{(1)}\mathrm{P}^{-1}, \quad
\tilde{\mathrm{A}}^{(2)} \defeq
\mathrm{P}\mathrm{A}^{(2)}\mathrm{P}^{-1}, \quad
\tilde{\mathrm{A}}^{(3)} \defeq
\mathrm{P}\mathrm{A}^{(3)}\mathrm{P}^{-1}\]
To help obtaining the expression of \(\tilde{\mathrm{A}}^{(1)},
\tilde{\mathrm{A}}^{(2)}\), and \(\tilde{\mathrm{A}}^{(3)}\), we substitute
(5) into (11) and get
(17)\[\begin{split}\mathrm{A}^{(1)} &= \left(\begin{array}{cccc}
0 & 1 & 0 & 0 \\
-v_1^2 + \frac{K}{\rho} & 2v_1 & 0 & 0 \\
-v_1v_2 & v_2 & v_1 & 0 \\
-v_1v_3 & v_3 & 0 & v_1
\end{array}\right), \quad
\mathrm{A}^{(2)} = \left(\begin{array}{cccc}
0 & 0 & 1 & 0 \\
-v_1v_2 & v_2 & v_1 & 0 \\
-v_2^2 + \frac{K}{\rho} & 0 & 2v_2 & 0 \\
-v_2v_3 & 0 & v_3 & v_2
\end{array}\right), \\
\mathrm{A}^{(3)} &= \left(\begin{array}{cccc}
0 & 0 & 0 & 1 \\
-v_1v_3 & v_3 & 0 & v_1 \\
-v_2v_3 & 0 & v_3 & v_2 \\
-v_3^2 + \frac{K}{\rho} & 0 & 0 & 2v_3
\end{array}\right)\end{split}\]
The Jacobian matrices in (16) can be spelled out with the
expressions in (13), (14), and (17)
(18)\[\begin{split}\tilde{\mathrm{A}}^{(1)} = \left(\begin{array}{cccc}
v_1 & \rho & 0 & 0 \\
\frac{K}{\rho^2} & v_1 & 0 & 0 \\
0 & 0 & v_1 & 0 \\
0 & 0 & 0 & v_1
\end{array}\right), \quad
\tilde{\mathrm{A}}^{(2)} = \left(\begin{array}{cccc}
v_2 & 0 & \rho & 0 \\
0 & v_2 & 0 & 0 \\
\frac{K}{\rho^2} & 0 & v_2 & 0 \\
0 & 0 & 0 & v_2
\end{array}\right), \quad
\tilde{\mathrm{A}}^{(3)} = \left(\begin{array}{cccc}
v_3 & 0 & 0 & \rho \\
0 & v_3 & 0 & 0 \\
0 & 0 & v_3 & 0 \\
\frac{K}{\rho^2} & 0 & 0 & v_3
\end{array}\right)\end{split}\]
(15) is hyperbolic where the linear combination of its Jacobian
matrices \(\tilde{\mathrm{A}}^{(1)}\), \(\tilde{\mathrm{A}}^{(2)}\),
and \(\tilde{\mathrm{A}}^{(3)}\)
(19)\[\begin{split}\tilde{\mathrm{R}} \defeq \sum_{i=1}^3 k_i\tilde{\mathrm{A}}^{(i)}
= \left(\begin{array}{cccc}
\sum_{i=1}^3 k_iv_i & k_1\rho & k_2\rho & k_3\rho \\
\frac{k_1K}{\rho^2} & \sum_{i=1}^3 k_iv_i & 0 & 0 \\
\frac{k_2K}{\rho^2} & 0 & \sum_{i=1}^3 k_iv_i & 0 \\
\frac{k_3K}{\rho^2} & 0 & 0 & \sum_{i=1}^3 k_iv_i
\end{array}\right)\end{split}\]
where \(k_1, k_2\), and \(k_3\) are real and bounded.
The linearly combined Jacobian matrix can be used to rewrite the
three-dimensional governing equation, (15), into one-dimensional
space
(20)\[\dpd{\tilde{\bvec{u}}}{t} + \tilde{\mathrm{R}}\dpd{\tilde{\bvec{u}}}{y} = 0\]
where
(21)\[y \defeq \sum_{i=1}^3 k_ix_i\]
and aided by (19) and the chain rule
\[\sum_{i=1}^3 \tilde{\mathrm{A}}^{(i)}
\dpd{\tilde{\bvec{u}}}{x_i} =
\sum_{i=1}^3 \tilde{\mathrm{A}}^{(i)}
\dpd{\tilde{\bvec{u}}}{y} \dpd{y}{x_i} =
\sum_{i=1}^3 k_{i} \tilde{\mathrm{A}}^{(i)}
\dpd{\tilde{\bvec{u}}}{y} =
\tilde{\mathrm{R}}\dpd{\tilde{\bvec{u}}}{y}\]
The eigenvalues of \(\tilde{\mathrm{R}}\) can be found by solving the
polynomial equation \(\det(\tilde{\mathrm{R}} -
\lambda\mathrm{I}_{4\times4}) = 0\) for \(\lambda\), and are
(22)\[\lambda_{1, 2, 3, 4} = \phi, \phi,
\phi+\sqrt{\frac{K\psi}{\rho}},
\phi-\sqrt{\frac{K\psi}{\rho}}\]
where \(\phi \defeq \sum_{i=1}^3 k_iv_i\), and \(\psi \defeq
\sum_{i=1}^3 k_i^2\). The corresponding eigenvector matrix is
(23)\[\begin{split}\mathrm{T} = \left(\begin{array}{cccc}
0 & 0 &
\sqrt{\frac{\rho^3\psi}{K}} & \sqrt{\frac{\rho^3\psi}{K}} \\
k_3 & 0 & k_1 & -k_1 \\
0 & k_3 & k_2 & -k_2 \\
-k_1 & -k_2 & k_3 & -k_3
\end{array}\right)\end{split}\]
The left eigenvector matrix is
(24)\[\begin{split}\mathrm{T}^{-1} = \left(\begin{array}{cccc}
0 & -\frac{k_1^2-\psi}{k_3\psi} & -\frac{k_1k_2}{k_3\psi} & -\frac{k_1}{\psi} \\
0 & -\frac{k_1k_2}{k_3\psi} & -\frac{k_2^2-\psi}{k_3\psi} & -\frac{k_2}{\psi} \\
\frac{1}{2\sqrt{\frac{\rho^3\psi}{K}}} &
\frac{k_1}{2\psi} & \frac{k_2}{2\psi} & \frac{k_3}{2\psi} \\
\frac{1}{2\sqrt{\frac{\rho^3\psi}{K}}} &
-\frac{k_1}{2\psi} & -\frac{k_2}{2\psi} & -\frac{k_3}{2\psi}
\end{array}\right)\end{split}\]
Riemann Invariants
Aided by (23) and (24), \(\tilde{\mathrm{R}}\)
can be diagonalized as
(25)\[\begin{split}\mathrm{\Lambda} \defeq \left(\begin{array}{cccc}
\lambda_1 & 0 & 0 & 0 \\
0 & \lambda_2 & 0 & 0 \\
0 & 0 & \lambda_3 & 0 \\
0 & 0 & 0 & \lambda_4
\end{array}\right) = \left(\begin{array}{cccc}
\phi & 0 & 0 & 0 \\
0 & \phi & 0 & 0 \\
0 & 0 & \phi + \sqrt{\frac{K\psi}{\rho}} & 0 \\
0 & 0 & 0 & \phi - \sqrt{\frac{K\psi}{\rho}}
\end{array}\right) = \mathrm{T}^{-1}\tilde{\mathrm{R}}\mathrm{T}\end{split}\]
(25) defines the eigenvalue matrix of \(\tilde{\mathrm{R}}\).
Aach element in the diagonal of the eigenvalue matrix is the eigenvalue listed
in (22). Pre-multiplying (20) with
\(\mathrm{T}^{-1}\) gives
\[\mathrm{T}^{-1}\dpd{\tilde{\bvec{u}}}{t}
+ \mathrm{\Lambda}\mathrm{T}^{-1}\dpd{\tilde{\bvec{u}}}{y} = 0\]
Define the characteristic variables
(26)\[\begin{split}\hat{\bvec{u}} \defeq \left(\begin{array}{c}
-\frac{k_1^2-\psi}{k_3\psi}v_1 - \frac{k_1k_2}{k_3\psi}v_2 - \frac{k_1}{\psi}v_3 \\
-\frac{k_1k_2}{k_3\psi}v_1 - \frac{k_2^2-\psi}{k_3\psi}v_2 - \frac{k_2}{\psi}v_3 \\
-\sqrt{\frac{K}{\rho\psi}} + \frac{k_1}{2\psi}v_1 + \frac{k_2}{2\psi}v_2 + \frac{k_3}{2\psi}v_3 \\
-\sqrt{\frac{K}{\rho\psi}} - \frac{k_1}{2\psi}v_1 - \frac{k_2}{2\psi}v_2 - \frac{k_3}{2\psi}v_3
\end{array}\right)\end{split}\]
such that
\[\mathrm{T}^{-1} = \dpd{\hat{\bvec{u}}}{\tilde{\bvec{u}}}\]
Then aided by the chain rule, we can write
(27)\[\dpd{\hat{\bvec{u}}}{t} + \mathrm{\Lambda}\dpd{\hat{\bvec{u}}}{y} = 0\]
The components of \(\hat{\bvec{u}}\) defined in (26) are the
Riemann invariants.
Diffusion Term Treatment
The momentum equation (4) contains the diffusion term
\[\sum_{i=1}^3 \dpd{}{x_i}
\left[\delta_{ij} \lambda \sum_{k=1}^3 \dpd{v_k}{x_k}
+ \mu \left( \dpd{v_i}{x_j}+\dpd{v_j}{x_i} \right)\right],
\quad j = 1, 2, 3\]
Define
(28)\[\begin{split}\bvec{g}^{(1)} \defeq \left(\begin{array}{c}
0 \\
\lambda\sum_{k=1}^3\dpd{v_k}{x_k} + 2\mu\dpd{v_1}{x_1} \\
\mu(\dpd{v_1}{x_2} + \dpd{v_2}{x_1}) \\
\mu(\dpd{v_1}{x_3} + \dpd{v_3}{x_1})
\end{array}\right), \quad
\bvec{g}^{(2)} \defeq \left(\begin{array}{c}
0 \\
\mu(\dpd{v_2}{x_1} + \dpd{v_1}{x_2}) \\
\lambda\sum_{k=1}^3\dpd{v_k}{x_k} + 2\mu\dpd{v_2}{x_2} \\
\mu(\dpd{v_2}{x_3} + \dpd{v_3}{x_2})
\end{array}\right), \quad
\bvec{g}^{(3)} \defeq \left(\begin{array}{c}
0 \\
\mu(\dpd{v_3}{x_1} + \dpd{v_1}{x_3}) \\
\mu(\dpd{v_3}{x_2} + \dpd{v_2}{x_3}) \\
\lambda\sum_{k=1}^3\dpd{v_k}{x_k} + 2\mu\dpd{v_3}{x_3}
\end{array}\right)\end{split}\]