Usage¶

To use Noisy Duck in a project:

import noisyduck

Example: Annular Cylindrical Duct¶

Numerical eigenvalue/eigenvector decomposition:

import noisyduck as nd
eigenvalues, eigenvectors = nd.annulus.numerical.decomposition(omega,m,r,rho,u,v,w,p,gam,filter='acoustic')

Analytical eigenvalue/eigenvector decomposition:

import noisyduck as nd
eigenvalues, eigenvectors = nd.annulus.analytical.decomposition(omega,m,mach,ri,ro,n)

See cylindrical_annulus_uniform_flow.py.

Cylindrical Annuluar Duct¶

Analytical¶

This module provides a functionality for computing the eigenvalue/eigenvector decomposition of a uniform axial flow in an annular cylindrical duct. The decomposition is based on an analytical solution of the convected wave equation for pressure, yielding the acoustic part of the eigen decomposition. The eigenvectors from this decomposition correspond specifically to the acoustic pressure disturbance.

Theory:¶

The analysis here follows that of Moinier and Giles, “Eigenmode Analysis for Turbomachinery Applications”, Journal of Propulsion and Power, Vol. 21, No. 6, 2005.

For a uniform axial mean, linear pressure perturbations satisfy the convected wave equation

\[\left( \frac{\partial}{\partial t} + M \frac{\partial}{\partial z} \right)^2 p = \nabla^2 p \quad\quad\quad \lambda < r < 1\]

where \(M\) is the Mach number of the axial mean flow and \(r\) has been normalized by the outer duct radius. Boundary conditions for a hard-walled duct imply zero radial velocity, which is imposed using the condition

\[\frac{\partial p}{\partial r} = 0 \quad\quad \text{at} \quad\quad r=\lambda,1\]

Conducting a normal mode analysis of the governing equation involves assuming a solution with a form

\[p(r,\theta,z,t) = e^{j(\omega t + m \theta + k z)}P(r)\]

This leads to a Bessel equation

\[\frac{1}{r} \frac{d}{dr} \left(r \frac{dP}{dr} \right) + \left( \mu^2 - \frac{m^2}{r^2} \right)P = 0 \quad\quad \lambda < r < 1\]

where \(\mu^2 = (\omega + M k)^2 - k^2\). The solution to the Bessel equation is

\[P(r) = a J_m(\mu r) + b Y_m(\mu r)\]

where \(J_m\) and \(Y_m\) are Bessel functions. Applying boundary conditions to the general solution gives a set of two equations

\[\begin{split}\begin{bmatrix} J'_m(\mu \lambda) & Y'_m(\mu \lambda) \\ J'_m(\mu) & Y'_m(\mu) \end{bmatrix} \begin{bmatrix} a \\ b \end{bmatrix} = 0\end{split}\]

which has nontrivial solutions as long as the determinant is zero. Solving for the zeros of the determinant give \(\mu\), at which point the quadratic equation above can be solved for the axial wavenumbers \(k\).

Example:¶

eigenvalues, eigenvectors =
noisyduck.annulus.analytical.decomposition(omega,m,mach,r,n)

noisyduck.annulus.analytical.compute_eigenvalues(omega, mach, zeros)[source]¶

This procedure compute the analytical eigenvalues for the convected wave equation. A uniform set of nodes is created to search for sign changes in the value for the eigensystem. When a sign change is detected, it is known that a zero is close. The eigensystem is then passed to a bisection routine to find the location of the zero. The location corresponds to the eigenvalue for the system.

Parameters:	omega (float) – temporal wavenumber. m (int) – circumferential wavenumber. mach (float) – Mach number. ri (float) – inner radius of a circular annulus. ro (float) – outer radius of a circular annulus. n (int) – number of eigenvalues to compute.
Returns:	An array of the first ‘n’ eigenvalues for the system.

noisyduck.annulus.analytical.compute_eigenvector(r, sigma, m, zero)[source]¶

Return the eigenvector for the system.

Parameters:	r (np.array(float)) – array of radial locations. sigma (float) – ratio of inner to outer radius, ri/ro. m (int) – circumferential wavenumber. zero (float) – a zero of the determinant of the convected wave equation
Returns:	the eigenvector associated with the input ‘m’ and ‘zero’, evaluated at radial locations defined by the input array ‘r’. Length of the return array is len(r).

noisyduck.annulus.analytical.compute_zeros(m, mach, ri, ro, n)[source]¶

This procedure compute the zeros of the determinant for the convected wave equation.

A uniform set of nodes is created to search for sign changes in the value of the function. When a sign change is detected, it is known that a zero is close. The function is then passed to a bisection routine to find the location of the zero.

Parameters:	m (int) – circumferential wavenumber. mach (float) – Mach number ri (float) – inner radius of a circular annulus. ro (float) – outer radius of a circular annulus. n (int) – number of eigenvalues to compute.
Returns:	An array of the first ‘n’ eigenvalues for the system.

noisyduck.annulus.analytical.decomposition(omega, m, mach, r, n)[source]¶

This procedure computes the analytical eigen-decomposition of the convected wave equation. The eigenvectors returned correspond specifically with acoustic pressure perturbations.

Inner and outer radii are computed using min(r) and max(r), so it is important that these end-points are included in the incoming array of radial locations.

Parameters:	omega (float) – temporal wave number. m (int) – circumferential wave number. mach (float) – Mach number. r (float) – array of radius stations, including rmin and rmax. n (int) – number of eigenvalues/eigenvectors to compute.
Returns:	a tuple containing an array of eigenvalues, and an array of eigenvectors evaluated at radial locations.
Return type:	(eigenvalues, eigenvectors)

noisyduck.annulus.analytical.eigensystem(b, m, ri, ro)[source]¶

Computes the function associated with the eigensystem of the convected wave equation. The location of the zeros for this function correspond to the eigenvalues for the convected wave equation.

The solution to the Bessel equation with boundary conditions applied yields a system of two linear equations. .. math:

A x =
\begin{bmatrix}
    J'_m(\mu \lambda)   &   Y'_m(\mu \lambda) \\
    J'_m(\mu)           &   Y'_m(\mu)
\end{bmatrix}
\begin{bmatrix}
    x_1 \\ x_2
\end{bmatrix}
= 0

This procedure evaluates the function .. math:

det(A) = f(b) = J_m(b*ri)*Y_m(b*ro)  -  J_m(b*ro)*Y_m(b*ri)

So, this procedure can be passed to another routine such as numpy to find zeros.

Parameters:	b (float) – coordinate. m (int) – circumferential wave number. ri (float) – inner radius of a circular annulus. ro (float) – outer radius of a circular annulus.

Numerical¶

This module provides a functionality for computing the numerical eigenvalue/eigenvector decomposition of a uniform axial flow in an annular cylindrical duct. The decomposition is based on a normal mode analysis of the three-dimensional linearized Euler equations, which yields an eigensystem that is discretized and solved numerically.

Theory:¶

Filtering:¶

Example:¶

eigenvalues, eigenvectors_l, eigenvectors_r =
    noisyduck.annulus.numerical.decomposition(omega,
                                              m,
                                              r,
                                              rho,
                                              vr,
                                              vt,
                                              vz,
                                              p,
                                              gam,
                                              filter='acoustic')

noisyduck.annulus.numerical.construct_numerical_eigensystem_general(omega, m, r, rho, vr, vt, vz, p, gam, perturb_omega=True)[source]¶

Constructs the numerical representation of the eigenvalue problem associated with the three-dimensional linearized euler equations subjected to a normal mode analysis.

NOTE: If perturb_omega=True, a small imaginary part is added to the temporal frequency to facilitate determining the propagation direction of eigenmodes based on the sign of the imaginary part of their eigenvalue. That is: \(\omega = \omega - 10^{-5}\omega j\). See Moinier and Giles[2].

[1] Kousen, K. A., “Eigenmodes of Ducted Flows With Radially-Dependent: Axial and Swirl Velocity Components”, NASA/CR 1999-208881, March 1999.
[2] Moinier, P., and Giles, M. B., “Eigenmode Analysis for Turbomachinery: Applications”, Journal of Propulsion and Power, Vol. 21, No. 6, November-December 2005.

Parameters:

omega (float) – temporal frequency.
m (int) – circumferential wavenumber.
r (float) – array of equally-spaced radius locations for the discretization, including end points.
rho (float) – mean density.
vr (float) – mean radial velocity.
vt (float) – mean tangential velocity.
vz (float) – mean axial velocity.
p (float) – mean pressure.
gam (float) – ratio of specific heats.
perturb_omega (bool) – If true, small imaginary part is added to the temporal frequency. Can help in determining direction of propagation.

Returns:

left-hand side of generalized eigenvalue problem, right-hand: side of generalized eigenvalue problem.

Return type:

(M, N)

noisyduck.annulus.numerical.decomposition(omega, m, r, rho, vr, vt, vz, p, gam, filter='None', alpha=1e-05, equation='general', perturb_omega=True)[source]¶

Compute the numerical eigen-decomposition of the three-dimensional linearized Euler equations for a cylindrical annulus.

Parameters:	omega (float) – temporal frequency. m (int) – circumferential wavenumber. r (float) – array of equally-spaced radius locations for the discretization, including end points. rho (float) – mean density. vr (float) – mean radial velocity. vt (float) – mean tangential velocity. vz (float) – mean axial velocity. p (float) – mean pressure. gam (float) – ratio of specific heats. filter (string, optional) – Optional filter for eigenmodes. values = [‘None’, ‘acoustic’] alpha (float, optional) – Criteria governing filtering acoustic modes. equation (string, optional) – Select from of governing equation for the decomposition. values = [‘general’, ‘radial equilibrium’] perturb_omega (bool) – If true, small imaginary part is added to the temporal frequency. Can help in determining direction of propagation.
Returns:	a tuple containing an array of eigenvalues, an array of left eigenvectors evaluated at radial locations, an array of right eigenvectors evaluated at radial locations.
Return type:	(eigenvalues, left_eigenvectors, right_eigenvectors)

Note

The eigenvectors being returned include each field \([\rho,v_r,v_t,v_z,p]\). The primitive variables can be extracted into their own eigenvectors by copying out those entries from the returned eigenvectors as:

res = len(r)
rho_eigenvectors = eigenvectors[0*res:1*res,:]
vr_eigenvectors  = eigenvectors[1*res:2*res,:]
vt_eigenvectors  = eigenvectors[2*res:3*res,:]
vz_eigenvectors  = eigenvectors[3*res:4*res,:]
p_eigenvectors   = eigenvectors[4*res:5*res,:]