Adjoint characteristic decomposition of one-dimensional waves

TL;DR

This paper explores the connection between continuous and discrete adjoint characteristic decompositions in one-dimensional acoustic wave propagation, demonstrating their relationship via a similarity transformation in the Laplace domain.

Contribution

It establishes a theoretical link between continuous and discrete adjoint methods through a similarity transformation, applied to acoustic wave problems in the Laplace domain.

Findings

Adjoint characteristic decompositions are connected by a similarity transformation.

The framework is demonstrated on a one-dimensional acoustic resonator.

The approach enables future adjoint-based design of acoustic networks.

Abstract

Adjoint methods enable the accurate calculation of the sensitivities of a quantity of interest. The sensitivity is obtained by solving the adjoint system, which can be derived by continuous or discrete adjoint strategies. In acoustic wave propagation, continuous and discrete adjoint methods have been developed to compute the eigenvalue sensitivity to design parameters and passive devices (Aguilar, J. G. et al, 2017, J. Computational Physics, vol. 341, 163-181). In this short communication, it is shown that the continuous and discrete adjoint characteristic decompositions, and Riemann invariants, are connected by a similarity transformation. The results are shown in the Laplace domain. The adjoint characteristic decomposition is applied to a one-dimensional acoustic resonator, which contains a monopole source of sound. The proposed framework provides the foundation to tackle larger acoustic networks with a discrete adjoint approach, opening up new possibilities for adjoint-based design of problems that can be solved by the method of characteristics.

Figures

Equations

$$\frac{\partial u^{\prime }}{\partial t}+\frac{1}{\bar{\rho }}\frac{\partial p^{\prime }}{\partial x}$$

$$\frac{\partial p^{\prime }}{\partial t}+\gamma \bar{p}\frac{\partial u^{\prime }}{\partial x}$$

$$s\hat{u}+\frac{1}{\bar{\rho }}\frac{d\hat{p}}{dx}$$

$$s\hat{p}+\gamma \bar{p}\frac{d\hat{u}}{dx}$$

$$\mathcal{N}(s)\underline{\hat{q}}=0,$$

$$[[\hat{p}{]]}_{b^{-}}^{b^{+}}=0,[[\hat{u}{]]}_{b^{-}}^{b^{+}}=\frac{\gamma -1}{\gamma \bar{p}}\hat{\mathcal{Q}},$$

$$\left[\begin{array}{c}{\hat{p}}_i\\ {\hat{u}}_i\end{array}\right]=\left[\begin{array}{cc}\exp \left(-s\frac{x-b}{{\bar{c}}_i}\right)& \exp \left(s\frac{x-b}{{\bar{c}}_i}\right)\\ \frac{1}{{\bar{\rho }}_i{\bar{c}}_i}\exp \left(-s\frac{x-b}{{\bar{c}}_i}\right)& \frac{-1}{{\bar{\rho }}_i{\bar{c}}_i}\exp \left(s\frac{x-b}{{\bar{c}}_i}\right)\end{array}\right]\left[\begin{array}{c}{\hat{f}}_i\\ {\hat{g}}_i\end{array}\right]$$

$$\underset{\hat{\mathbf{q}}(x)}{\underbrace{\left[\begin{array}{c}{\hat{p}}_1\\ {\hat{u}}_1\\ {\hat{p}}_2\\ {\hat{u}}_2\end{array}\right]}}=\underset{\mathbf{\Phi }(x)}{\underbrace{\left[\begin{array}{cccc}\exp (s\frac{x-b}{{\bar{c}}_1})& 0& \exp (-s\frac{x-b}{{\bar{c}}_1})& 0\\ -\frac{1}{{\bar{\rho }}_1{\bar{c}}_1}\exp (s\frac{x-b}{{\bar{c}}_1})& 0& \frac{1}{{\bar{\rho }}_1{\bar{c}}_1}\exp (-s\frac{x-b}{{\bar{c}}_1})& 0\\ 0& \exp (-s\frac{x-b}{{\bar{c}}_2})& 0& \exp (s\frac{x-b}{{\bar{c}}_2})\\ 0& \frac{1}{{\bar{\rho }}_2{\bar{c}}_2}\exp (-s\frac{x-b}{{\bar{c}}_2})& 0& -\frac{1}{{\bar{\rho }}_2{\bar{c}}_2}\exp (s\frac{x-b}{{\bar{c}}_2})\end{array}\right]}}\underset{\hat{\mathbf{w}}\equiv \left[\begin{array}{c}{\hat{\mathbf{w}}}_{out}\\ {\hat{\mathbf{w}}}_{in}\end{array}\right]}{\underbrace{\left[\begin{array}{c}{\hat{g}}_1\\ {\hat{f}}_2\\ {\hat{f}}_1\\ {\hat{g}}_2\end{array}\right],}}$$

$${\hat{f}}_1=R_1{\hat{g}}_1\exp (-s{\tau }_1),{\hat{g}}_2=R_2{\hat{f}}_2\exp (-s{\tau }_2),$$

$$\underset{{\mathbf{N}}_{out}(s)}{\underbrace{\left[\begin{array}{cc}1+R_1\exp (-s{\tau }_1)& -1-R_2\exp (-s{\tau }_2)\\ (1-R_1\exp (-s{\tau }_1))\left(1+\frac{\gamma -1}{\gamma \bar{p}}n\exp (-s\tau )\right)& \frac{{\bar{c}}_2}{{\bar{c}}_1}(1-R_2\exp (-s{\tau }_2))\end{array}\right]}}\underset{{\hat{\mathbf{w}}}_{out}}{\underbrace{\left[\begin{array}{c}{\hat{g}}_1\\ {\hat{f}}_2\end{array}\right]}}=\left[\begin{array}{c}0\\ 0\end{array}\right],$$

$$\underset{{\mathbf{N}}_{in}(s)}{\underbrace{\left[\begin{array}{cc}-1-R_1^{-1}\exp (s{\tau }_1)& 1+R_2^{-1}\exp (s{\tau }_2)\\ -\frac{{\bar{c}}_1}{{\bar{c}}_2}(1-R_1^{-1}\exp (s{\tau }_1))\left(1+\frac{\gamma -1}{\gamma \bar{p}}n\exp (-s\tau )\right)& -(1-R_2^{-1}\exp (s{\tau }_2))\end{array}\right]}}\underset{{\hat{\mathbf{w}}}_{in}}{\underbrace{\left[\begin{array}{c}{\hat{f}}_1\\ {\hat{g}}_2\end{array}\right]}}=\left[\begin{array}{c}0\\ 0\end{array}\right],$$

$$\frac{\partial u^{+}}{\partial t}+{\bar{c}}^2\frac{\partial p^{+}}{\partial x}$$

$$\frac{\partial p^{+}}{\partial t}+\frac{\partial u^{+}}{\partial x}$$

$$-s^{\ast }{\hat{u}}^{+}+{\bar{c}}^2\frac{d{\hat{p}}^{+}}{dx}$$

$$-s^{\ast }{\hat{p}}^{+}+\frac{d{\hat{u}}^{+}}{dx}$$

$$\left[\begin{array}{c}{\hat{p}}_i^{+}\\ {\hat{u}}_i^{+}\end{array}\right]=\left[\begin{array}{cc}\exp \left(-s^{\ast }\frac{x-b}{{\bar{c}}_i}\right)& \exp \left(s^{\ast }\frac{x-b}{{\bar{c}}_i}\right)\\ -{\bar{c}}_i\exp \left(-s^{\ast }\frac{x-b}{{\bar{c}}_i}\right)& {\bar{c}}_i\exp \left(s^{\ast }\frac{x-b}{{\bar{c}}_i}\right)\end{array}\right]\left[\begin{array}{c}{\hat{g}}_{i,c}^{+}\\ {\hat{f}}_{i,c}^{+}\end{array}\right].$$

$$\underset{{\hat{\mathbf{q}}}_c^{+}(x)}{\underbrace{\left[\begin{array}{c}{\hat{p}}_1^{+}\\ {\hat{u}}_1^{+}\\ {\hat{p}}_2^{+}\\ {\hat{u}}_2^{+}\end{array}\right]}}=\underset{{\mathbf{\Phi }}_c^{+}(x)}{\underbrace{\left[\begin{array}{cccc}\exp (s^{\ast }\frac{x-b}{{\bar{c}}_1})& 0& \exp (-s^{\ast }\frac{x-b}{{\bar{c}}_1})& 0\\ {\bar{c}}_1\exp (s^{\ast }\frac{x-b}{{\bar{c}}_1})& 0& -{\bar{c}}_1\exp (-s^{\ast }\frac{x-b}{{\bar{c}}_1})& 0\\ 0& \exp (-s^{\ast }\frac{x-b}{{\bar{c}}_2})& 0& \exp (s^{\ast }\frac{x-b}{{\bar{c}}_2})\\ 0& -{\bar{c}}_2\exp (-s^{\ast }\frac{x-b}{{\bar{c}}_2})& 0& {\bar{c}}_2\exp (s^{\ast }\frac{x-b}{{\bar{c}}_2})\end{array}\right]}}\underset{{\hat{\mathbf{w}}}_c^{+}\equiv \left[\begin{array}{c}{\hat{\mathbf{w}}}_{c,out}^{+}\\ {\hat{\mathbf{w}}}_{c,in}^{+}\end{array}\right]}{\underbrace{\left[\begin{array}{c}{\hat{f}}_{1,c}^{+}\\ {\hat{g}}_{2,c}^{+}\\ {\hat{g}}_{1,c}^{+}\\ {\hat{f}}_{2,c}^{+}\end{array}\right],}}$$

$$[[{\hat{u}}^{+}{]]}_{b^{-}}^{b^{+}}=0,\frac{\gamma \bar{p}}{\gamma -1}[[{\hat{p}}^{+}{]]}_{b^{-}}^{b^{+}}=-n{\hat{p}}^{+}(b^{+})\exp (-s^{\ast }\tau ).$$

$${\hat{g}}_{1,c}^{+}=R_1^{+-1}{\hat{f}}_{1,c}^{+}\exp (-s^{\ast }{\tau }_1),{\hat{f}}_{2,c}^{+}=R_2^{+-1}{\hat{g}}_{2,c}^{+}\exp (-s^{\ast }{\tau }_2).$$

$$R_1^{+}=R_1^{-1\ast },R_2^{+}=R_2^{-1\ast }.$$

$$\underset{{\mathbf{N}}_{c,out}^{+}(s^{\ast })}{\underbrace{\left[\begin{array}{cc}1+R_1^{+-1}\exp (-s^{\ast }{\tau }_1)& -(1+R_2^{+-1}\exp (-s^{\ast }{\tau }_2))\left(1+\frac{\gamma -1}{\gamma \bar{p}}n\exp (-s^{\ast }\tau )\right)\\ 1-R_1^{+-1}\exp (-s^{\ast }{\tau }_1)& \frac{{\bar{c}}_2}{{\bar{c}}_1}(1-R_2^{+-1}\exp (-s^{\ast }{\tau }_2))\end{array}\right]}}\underset{{\hat{\mathbf{w}}}_{c,out}^{+}}{\underbrace{\left[\begin{array}{c}{\hat{f}}_{1,c}^{+}\\ {\hat{g}}_{2,c}^{+}\end{array}\right]}}=\left[\begin{array}{c}0\\ 0\end{array}\right].$$

$$\underset{{\mathbf{N}}_{c,in}^{+}(s^{\ast })}{\underbrace{\left[\begin{array}{cc}-1-R_1^{+}\exp (s^{\ast }{\tau }_1)& \left(1+R_2^{+}\exp (s^{\ast }{\tau }_2)\right)\left(1+\frac{\gamma -1}{\gamma \bar{p}}n\exp (-s^{\ast }\tau )\right)\\ -\frac{{\bar{c}}_1}{{\bar{c}}_2}\left(1-R_1^{+}\exp (s^{\ast }{\tau }_1)\right)& -\left(1-R_2^{+}\exp (s^{\ast }{\tau }_2)\right)\end{array}\right]}}\underset{{\hat{\mathbf{w}}}_{c,in}^{+}}{\underbrace{\left[\begin{array}{c}{\hat{g}}_{1,c}^{+}\\ {\hat{f}}_{2,c}^{+}\end{array}\right]}}=\left[\begin{array}{c}0\\ 0\end{array}\right].$$

$$\underset{{\mathbf{N}}_{out}^H(s^{\ast })}{\underbrace{\left[\begin{array}{cc}1+R_1^{\ast }\exp (-s^{\ast }{\tau }_1)& (1-R_1^{\ast }\exp (-s^{\ast }{\tau }_1))\left(1+\frac{\gamma -1}{\gamma \bar{p}}n\exp (-s^{\ast }\tau )\right)\\ -(1+R_2^{\ast }\exp (-s^{\ast }{\tau }_2))& \frac{{\bar{c}}_2}{{\bar{c}}_1}(1-R_2^{\ast }\exp (-s^{\ast }{\tau }_2))\end{array}\right]}}\underset{{\hat{\mathbf{w}}}_{d,out}^{+}}{\underbrace{\left[\begin{array}{c}{\hat{g}}_{1,d}^{+}\\ {\hat{f}}_{2,d}^{+}\end{array}\right]}}=\left[\begin{array}{c}0\\ 0\end{array}\right],$$

$$\underset{{\mathbf{N}}_{in}^H(s^{\ast })}{\underbrace{\left[\begin{array}{cc}-1-R_1^{-1\ast }\exp (s^{\ast }{\tau }_1)& -\frac{{\bar{c}}_1}{{\bar{c}}_2}(1-R_1^{-1\ast }\exp (s^{\ast }{\tau }_1))\left(1+\frac{\gamma -1}{\gamma \bar{p}}n\exp (-s^{\ast }\tau )\right)\\ 1+R_2^{-1\ast }\exp (s^{\ast }{\tau }_2)& -(1-R_2^{-1\ast }\exp (s^{\ast }{\tau }_2))\end{array}\right]}}\underset{{\hat{\mathbf{w}}}_{d,in}^{+}}{\underbrace{\left[\begin{array}{c}{\hat{f}}_{1,d}^{+}\\ {\hat{g}}_{2,d}^{+}\end{array}\right]}}=\left[\begin{array}{c}0\\ 0\end{array}\right].$$

$${\mathbf{N}}^H=\mathbf{Q\Lambda }{\mathbf{Q}}^{-1},{\mathbf{N}}_c^{+}={\mathbf{Q}}_c\mathbf{\Lambda }{\mathbf{Q}}_c^{-1},$$

$${\mathbf{N}}_c^{+}{\mathbf{P}}^{-1}={\mathbf{P}}^{-1}{\mathbf{N}}^H.$$

$${\mathbf{P}}^{-1}={\mathbf{Q}}_c{\mathbf{Q}}^{-1}.$$

$${\mathbf{P}}^{-1}=\left[\begin{array}{cc}\frac{N_{21}^H}{N_{c,21}^{+}}& \frac{N_{c,11}^{+}-N_{11}^H}{N_{c,21}^{+}}\\ 0& 1\end{array}\right],$$

$${\mathbf{P}}_{out}^{-1}$$

$${\mathbf{P}}_{in}^{-1}$$

$${\mathbf{P}}_{out}^{-1}{\hat{\mathbf{w}}}_{d,out}^{+}={\hat{\mathbf{w}}}_{c,out}^{+},{\mathbf{P}}_{in}^{-1}{\hat{\mathbf{w}}}_{d,in}^{+}={\hat{\mathbf{w}}}_{c,in}^{+}.$$

$${\hat{\mathbf{q}}}_c^{+}(x)=\underset{{\mathbf{\Phi }}_d^{+}(x)}{\underbrace{{\mathbf{\Phi }}_c^{+}(x)\left[\begin{array}{cc}{\mathbf{P}}_{out}^{-1}& \left[\begin{array}{cc}0& 0\\ 0& 0\end{array}\right]\\ \left[\begin{array}{cc}0& 0\\ 0& 0\end{array}\right]& {\mathbf{P}}_{in}^{-1}\end{array}\right]}}\left[\begin{array}{c}{\hat{\mathbf{w}}}_{d,out}^{+}\\ {\hat{\mathbf{w}}}_{d,in}^{+}\end{array}\right].$$

$$\delta \sigma (x_p)=\frac{{\hat{\mathbf{q}}}^{+\ast }(x_p)\cdot \mathbf{P}\hat{\mathbf{q}}(x_p)}{{\int }_0^L{\hat{\mathbf{q}}}^{+\ast }(x)\cdot \frac{\partial \mathbf{N}}{\partial \sigma }\hat{\mathbf{q}}(x)dx},$$

Full text

††journal: Journal of Computational Physics

1 Introduction

In acoustic wave propagation, the primitive variables, e.g., the acoustic velocity and pressure, can be decomposed by the method of characteristics. The solutions of the characteristic equations are the Riemann invariants, which are the acoustic travelling waves. On the one hand, it was shown that the discrete adjoint approach is very straightforward to implement because it only requires the complex transpose of the matrix applied to the Riemann invariants [1]. The sensitivity of the eigenvalue to the design parameters, which is the quantity of interest, is computed to machine precision. However, with the discrete adjoint approach, it is still an open question how to calculate the adjoint primitive variables, which contain the spatial information, from the adjoint Riemann invariants. On the other hand, the continuous adjoint approach enables the calculation of the adjoint primitive variables from the adjoint Riemann invariants, i.e., no spatial information is lost. The downside of the continuous adjoints in wave propagation is that the sensitivities are not guaranteed to be accurate to machine precision, and more involved mathematical derivation is required even in simple systems [see Secs. 2.4 and 5.4 in 1]. As reviewed in [2], the spatial information contained in the adjoint primitive variables is crucial for the design of acoustic resonators, for example, by (i) optimal placement of passive devices, such as acoustic dampers, to suppress loud oscillations; (ii) optimal change of the shape of the acoustic resonator to prevent loud oscillations from occurring; and (iii) computation of the sensitivity to the acoustic boundary conditions. In this paper, the mathematical connection that relates the characteristic decompositions, and Riemann invariants, of the continuous and discrete adjoint systems is found. Therefore, the spatial dependence of the solution is recovered in the discrete adjoint approach.

2 Wave propagation with a monopole source of sound

An acoustic resonator with a monopole source of sound, i.e., a flame, is considered as a prototypical model of a thermoacoustic system (Fig. 1). The main assumptions of the model are [1] (i) the ratio between the radius and length of the duct is sufficiently small such that only low-frequency longitudinal acoustics propagate (the axial coordinate is $x$); (ii) the length of the flame is sufficiently smaller than the acoustic wavelength such that the heat is released at $x=b$ (compact assumption); (iii) the acoustics evolve on a slow mean flow such that the mean-flow Mach number is negligible; (iv) the flow is isentropic, except at the flame location; and (v) the gas is ideal, $p=\rho\mathcal{R}T$, where $p$ is the pressure, $\rho$ is the density, $T$ is the temperature, and $\mathcal{R}$ is the gas constant. The speed of sound is $c=\sqrt{\gamma\mathcal{R}T}$, where $\gamma$ is the heat capacity ratio. For more details the reader may refer to [1]. The problem is governed by the nonlinear Euler equations with a model for the flame. Because the focus is on the acoustics of the resonator, which are considered as linear unsteady perturbations, a generic flow variable is decomposed as $(\cdot)=\bar{(\cdot)}+(\cdot)^{\prime}$, where $\bar{(\cdot)}\sim\mathcal{O}(1)$ is the steady mean flow component, and $(\cdot)^{\prime}\sim\mathcal{O}(\epsilon)$, with $\epsilon\rightarrow 0$, is the acoustic fluctuation. By substituting this decomposition in the nonlinear Euler equation, the mean-flow and acoustic equations are derived. On grouping the terms $\sim\mathcal{O}(1)$, the mean flow is fully characterized by algebraic equations that link the quantities upstream and downstream of the heat source. Two thermodynamic variables upstream of the heat source, e.g., $\bar{T}_{1}$ and $\bar{p}_{1}$, need to be set along with the increase in the mean-flow temperature, $\Delta\bar{T}$, induced by the heat source, such that $\bar{T}_{2}=\bar{T}_{1}+\Delta\bar{T}$. The mean-flow pressure is constant in a zero-Mach number flow, i.e., $\bar{p}_{1}=\bar{p}_{2}\equiv\bar{p}$. The acoustic equations are obtained by grouping the variables $\sim\mathcal{O}(\epsilon)$. The linearized momentum and energy equations read, respectively

[TABLE]

where $t$ is the time; $u^{\prime}$ and $p^{\prime}$ are the acoustic velocity and pressure, respectively. The equations are studied in the Laplace domain by modal decomposition $(\cdot)^{\prime}(x,t)=\hat{(\cdot)}(x)\exp(st)$, where $s$ is a complex number, yielding

[TABLE]

In compact notation, Eqs. (2a)-(2b) are represented by

[TABLE]

where $\mathcal{N}$ is the linear differential operator that encapsulates the boundary conditions, and $\underline{\hat{q}}=[\hat{u},\hat{p}]^{T}$. Equations (2a)-(2b) hold at any acoustic location except at the heat source, i.e., $x\in(0,L)-\{b\}$. At the heat-source location, $x=b$, the linearized momentum and energy equations in integral form provide the jump conditions [1]

[TABLE]

where $\llbracket\cdot\rrbracket^{b^{+}}_{b^{-}}$ is the difference between the quantity $(\cdot)$ evaluated downstream of the flame, at $b^{+}$, and upstream of the flame, at $b^{-}$. The heat-release is $\hat{\mathcal{Q}}$. The flame is assumed to be perfectly premixed, such that the heat is released after a time delay, $\tau$, after that the acoustic velocity has perturbed the flame’s base, i.e., $\hat{\mathcal{Q}}=n\hat{u}(b^{-})\exp(-s\tau)$, where $n$ is the flame index. This is similar to the $n-\tau$ model proposed by Crocco [3], which captures the essential time-delayed physics of thermoacoustic instabilities. For $x\in(0,L)-\{b\}$, the primitive variables in Eqs. (2a)-(2b) are decomposed in Riemann invariants through a characteristic decomposition as

[TABLE]

where $i=1$ for a quantity upstream of the heat source, and $i=2$ for a quantity downstream of the heat source. The Riemann invariants $\hat{f}_{i}$ are the right-travelling waves, while the Riemann invariants $\hat{g}_{i}$ are the left-travelling waves (Fig. 1).

To connect the spatial dependency of the continuous and discrete adjoints (Sec. 4), it is convenient to re-arrange Eq. (5) as

[TABLE]

where $\hat{\mathbf{w}}_{out}\equiv[\hat{g}_{1},\hat{f}_{2}]^{T}$ are the outgoing waves, and $\hat{\mathbf{w}}_{in}\equiv[\hat{f}_{1},\hat{g}_{2}]^{T}$ are the ingoing waves. The problem is closed by the acoustic boundary conditions at the inlet, $x=0$, and outlet, $x=L$, which link right and left travelling waves as

[TABLE]

where $R_{1}$ and $R_{2}$ are the acoustic reflection coefficients, which may be complex, $\tau_{1}={2b}/{\bar{c}_{1}}$ and $\tau_{2}={2(L-b)}/{\bar{c}_{2}}$. By combining the jump conditions (4) and relations (7) with the characteristic decomposition (5), the stability of the acoustic system is governed by a $2\times 2$ nonlinear eigenvalue problem

[TABLE]

where the subscript $out$ denotes the outgoing travelling waves. Alternatively, by using relations (7), the nonlinear eigenproblem reads

[TABLE]

where the subscript $in$ denotes the ingoing travelling waves. The primitive variables, and their spatial variation, can be recovered from the characteristic decomposition in Eq. (5).

3 Adjoint equations

Adjoint methods enable the accurate and efficient calculation of the sensitivity of a quantity of interest to the design parameters. The adjoint system can be obtained by two strategies, which will be referred to as continuous and discrete adjoints [e.g., 4, 1, 2].

3.1 Continuous adjoint approach

A continuous adjoint approach provides the adjoint partial differential equations with respect to a sesquilinear form. Here, the following inner product $\left[a(x,t),b(x,t)\right]\equiv\int_{0}^{T}\left(\int_{0}^{b-}a^{*}b\;dx+\int_{b+}^{L}a^{*}b\;dx\right)dt$ is employed as a sesquilinear form, where ∗ is the complex conjugate. $a$ and $b$ represent generic complex functions. The adjoint partial differential equations are defined via the Lagrange-Green identity $\left[\underline{\hat{q}}^{+},\mathcal{N}\underline{\hat{q}}\right]=\left[\mathcal{N}^{+}\underline{\hat{q}}^{+},\underline{\hat{q}}\right]$. On integrating by parts, the adjoint partial differential equations read

[TABLE]

To transform the adjoint equations in the Laplace domain, the adjoint modal transformation $\underline{q}^{+}(x,t)=\underline{\hat{q}}^{+}(x)\Psi(t)$ has to be found. By eliminating the time dependency from the Lagrange-Green identity $\left[\underline{\hat{q}}^{+}\Psi(t),\underline{\hat{q}}\exp(st)\right]=\left[\underline{\hat{q}}^{+}\Psi(t)\exp(s^{*}t),\underline{\hat{q}}\right]$, it follows that $\Psi(t)=\exp(-s^{*}t)$. Therefore, the adjoint modal transformation is $(\cdot)^{+}(x,t)$ $=$ $\hat{(\cdot)}^{+}(x)\exp(-s^{*}t)$. The change in sign of the Laplace variable physically signifies that the adjoint equations evolve backward in time, i.e., they are anticausal. The adjoint equations in the Laplace domain read

[TABLE]

for $x\in(0,L)-\{b\}$. This set of equations is hyperbolic, hence, the solution is obtained by an adjoint characteristic decomposition

[TABLE]

The subscript $c$ stands for “continuous” (adjoint). Equation (12) is re-arranged as

[TABLE]

where $\hat{\mathbf{w}}^{+}_{c,out}\equiv[\hat{f}^{+}_{1,c},\hat{g}^{+}_{2,c}]^{T}$ are the outgoing adjoint waves, and $\hat{\mathbf{w}}^{+}_{c,in}\equiv[\hat{g}^{+}_{1,c},\hat{f}^{+}_{2,c}]^{T}$ are the ingoing adjoint waves. At the heat-source location, the adjoint momentum and energy equations in integral form provide the adjoint jump conditions [1]

[TABLE]

The problem is closed by the adjoint boundary conditions at the inlet, $x=0$, and outlet, $x=L$, which link right and left adjoint waves as

[TABLE]

Zeroing the adjoint boundary conditions, which stem from the integration of the Lagrange-Green identity, yields the adjoint reflection coefficients at the inlet, $x=0$, and outlet, $x=L$, respectively [1]

[TABLE]

The substitution of the adjoint characteristic decomposition (12) in the adjoint jump conditions (14) yields the continuous adjoint eigenproblem

[TABLE]

Alternatively, by using relations (15), the continuous adjoint eigenproblem reads

[TABLE]

Section 2.4 of [1] shows the full derivation of the continuous adjoint equations, jump conditions and boundary conditions.

3.2 Discrete adjoint approach

In the discrete adjoints, the linearized system (or algorithm) of the numerically discretized differential equations is transposed to spawn the adjoint system. In one-dimensional wave propagation, the discrete adjoint approach is straightforward: The adjoint eigenproblem is obtained by taking the conjugate transpose of the direct matrix (8)

[TABLE]

where H denotes the complex conjugate transpose, and the subscript $d$ stands for “discrete” (adjoint). Alternatively, from (9)

[TABLE]

In the discrete adjoint approach, the adjoint characteristic decomposition is not known a priori, therefore the adjoint primitive variables, and their spatial variation, cannot be calculated. Because the matrices from the discrete and continuous adjoint approaches are not equal to each other, i.e., $\mathbf{N}^{H}\not=\mathbf{N}_{c}^{+}$ (see Eqs. (17)-(20)), the continuous adjoint characteristic decomposition111whose spatial dependency is provided by $\boldsymbol{\Phi}_{c}^{+}(x)$ in Eq. (13) cannot be applied to the discrete adjoint Riemann invariants (19)-(20) to compute the adjoint primitive variables. The problem is schematically described in Fig. 2 and can be concisely stated as

Problem: What is the adjoint characteristic decomposition in the discrete adjoint approach?

Section 4 proposes a practical method to obtain the missing spatial information in the discrete adjoint approach.

4 The connection between continuous and discrete adjoints

An eigenvalue, $\sigma$, is the solution of the transcendental equation $\det(\mathbf{N}(\sigma))=0$. Although not necessary, defective systems [5] are not considered here for brevity, which implies that the ranks of $\mathbf{N}^{H}$ and $\mathbf{N}_{c}^{+}$ are equal to 1. By inspection of Eqs. (17)-(20), the matrices $\mathbf{N}^{H}$ and $\mathbf{N}_{c}^{+}$ have the same trascendental characteristic equation, i.e., $\det(\mathbf{N}^{H}(s^{*}))=\det(\mathbf{N}^{+}_{c}(s^{*}))$, therefore, they have the same eigenvalue, $\sigma^{*}$, i.e., $\det(\mathbf{N}^{H}(\sigma^{*}))=\det(\mathbf{N}^{+}_{c}(\sigma^{*}))=0$. Thus, matrices $\mathbf{N}^{H}(\sigma^{*})$ and $\mathbf{N}^{+}_{c}(\sigma^{*})$ are similar, which is the key property that enables the connection between continuous and discrete adjoint approaches. The adjoint matrices are being evaluated at the eigenvalue $\sigma^{*}$ from now on and the dependency on it will be dropped for brevity. To find the similarity transformation, the adjoint matrices are eigendecomposed222Any other decomposition that provides a basis would work, e.g., the singular value decomposition. as

[TABLE]

The similarity transformation (21) can be numerically calculated for problems with more degrees of freedom. However, in the $2\times 2$ matrix under investigation, an analytical expression can be derived: $\mathbf{\Lambda}=\begin{bmatrix}\lambda_{+}&0\\ 0&\lambda_{-}\end{bmatrix}$, where $\lambda_{-}=0$ because matrices (21) have rank 1. Therefore, the second eigenvalue reads $\lambda_{+}=\textrm{Tr}(\mathbf{N}^{H})=\textrm{Tr}(\mathbf{N}_{c}^{+})$. By definition of similarity, $\exists\;\mathbf{P}\in\mathbb{C}^{2\times 2}$ with its inverse such that

[TABLE]

Equation (22) shows that the continuous and discrete adjoint matrices are the same operator, which is represented in two different bases. Therefore $\mathbf{Q}_{c}^{-1}\mathbf{N}_{c}^{+}\mathbf{Q}_{c}=\mathbf{Q}^{-1}\mathbf{N}^{H}\mathbf{Q}\implies\mathbf{N}_{c}^{+}\mathbf{Q}_{c}\mathbf{Q}^{-1}=\mathbf{Q}_{c}\mathbf{Q}^{-1}\mathbf{N}^{H}$. By comparison with (22), the similarity transformation can be expressed as

[TABLE]

On working through the algebra and simplifying, it is possible to show that the similarity transformation reduces to a simple expression

[TABLE]

which is a similarity transformation for $2\times 2$ matrices that have rank 1 when ${N}^{+}_{c,21}\neq 0$. The similarity transformation (24) explicitly reads

[TABLE]

The relations between the discrete and continuous adjoint Riemann invariants read

[TABLE]

This brings us to the main result of this paper.

Result: The adjoint characteristic decomposition is given by

[TABLE]

With this connection, which is numerically verified in Fig. 3, the first-order change of the eigenvalue, $\delta\sigma$, due to a perturbation $\mathbf{P}\delta(x-x_{p})$, where $\delta(x-x_{p})$ is the Dirac delta distribution centred at $x_{p}$, can be calculated from the discrete adjoint approach as [6, 2]

[TABLE]

which is useful in sensitivity analysis.

4.1 Acoustics with a mean flow

Adding flow inhomogeneities and mean flow effects, for which the continuous adjoint formulation becomes cumbersome, requires no conceptual modification to the similarity transformation. The full set of equations can be found in [e.g., see Sec. 5.4 in ref. 1] and [7], which are not reported here for brevity. The problem with a mean flow has three primitive variables, the pressure, $\hat{p}$, the velocity, $\hat{u}$, and the density, $\hat{\rho}$. Although an analytical derivation is still possible, the $3\times 3$ similarity matrix in Eq. (23) is evaluated numerically in this case. Figure 4 shows the spatial dependence of the primitive variables and their adjoint obtained with the adjoint characteristic decomposition (28).

5 Conclusions

A systematic method is proposed to connect the continuous and discrete adjoint characteristic decompositions and adjoint Riemann invariants in wave propagation. The key mathematical observation is that the two adjoint operators are similar. To keep the mathematics simple, the method is shown in a low-order model of an acoustic resonator with a monopole source of sound. Because larger acoustic networks are composed of a connection of acoustic resonators, the adjoint characteristic decomposition can be scaled up, for example, in the design of stable gas turbine combustors [e.g., 8, 9, 10]. With the proposed adjoint characteristic decomposition, the adjoint primitive variables can be obtained from the adjoint Riemann invariants of the discrete adjoint approach. Thus, the missing spatial information of the discrete adjoint approach is found. The analytical derivation is numerically verified for a case without mean flow and a case with a mean flow. Such spatial information enables the optimal passive control of acoustic oscillations by exploiting the discrete adjoint approach, which is versatile, straightforward to implement and provides sensitivities that are accurate to machine precision. This paper provides the foundation to tackle larger acoustic networks. The connection between continuous and discrete adjoints can be used in adjoint-based design of other problems that can be solved by the method of characteristics.

Acknowledgements

Financial support from the Royal Academy of Engineering Research Fellowships Scheme is gratefully acknowledged. The author thanks J. Aguilar for fruitful discussions.