Reciprocity and representation theorems for flux- and field-normalised decomposed wave fields

TL;DR

This paper develops reciprocity and representation theorems for flux- and field-normalised decomposed wave fields, providing new integral formulas useful for reflection imaging and multiple scattering analysis in wave propagation.

Contribution

It systematically analyzes decomposition operators and derives reciprocity and representation theorems for decomposed wave fields, including Kirchhoff-Helmholtz integrals.

Findings

Derived reciprocity theorems for decomposed wave fields.

Established representation theorems including Kirchhoff-Helmholtz integrals.

Applied results to reflection imaging and multiple scattering scenarios.

Abstract

We consider wave propagation problems in which there is a preferred direction of propagation. To account for propagation in preferred directions, the wave equation is decomposed into a set of coupled equations for waves that propagate in opposite directions along the preferred axis. This decomposition is not unique. We discuss flux-normalised and field-normalised decomposition in a systematic way, analyse the symmetry properties of the decomposition operators and use these symmetry properties to derive reciprocity theorems for the decomposed wave fields, for both types of normalisation. Based on the field-normalised reciprocity theorems, we derive representation theorems for decomposed wave fields. In particular we derive double- and single-sided Kirchhoff-Helmholtz integrals for forward and backward propagation of decomposed wave fields. The single-sided Kirchhoff-Helmholtz integrals for backward propagation of field-normalised decomposed wave fields find applications in reflection imaging, accounting for multiple scattering.

Equations

$$-i\omega \alpha P+{\partial }_j{Q}_j=B,$$

$$-i\omega \beta Q_j+{\partial }_{j}P=C_j.$$

$$\displaystyle\beta\partial_{j}\Bigl{(}\frac{1}{\beta}\partial_{j}P\Bigr{)}+k^{2}P=\beta\partial_{j}\Bigl{(}\frac{1}{\beta}C_{j}\Bigr{)}+i\omega\beta B,$$

$$k^2=\alpha \beta {\omega }^2.$$

$$\displaystyle\partial_{3}{\bf q}={{\mbox{\boldmath${\cal A}$}}}{\bf q}+{\bf d},$$

$$\mathbf{q}=\left(\begin{array}{c}P\\ Q_3\end{array}\right),\mathbf{d}=\left(\begin{array}{c}{C}_3\\ B_0\end{array}\right),$$

$$B_0=B+\frac{1}{i\omega }{\partial }_{\nu }\frac{1}{\beta }{C}_{\nu }$$

$$\displaystyle\quad{{\mbox{\boldmath${\cal A}$}}}=\begin{pmatrix}0&{\cal A}_{12}\\ {\cal A}_{21}&0\end{pmatrix},$$

$${\mathcal{A}}_{12}=i\omega \beta ,$$

$${\mathcal{A}}_{21}=i\omega \alpha -\frac{1}{i\omega }{\partial }_{\nu }\frac{1}{\beta }{\partial }_{\nu }.$$

$$j=\frac{1}{4}\{P^{\ast }{Q}_3+Q_3^{\ast }P\}.$$

$${\int }_{\mathbb{S}}(\mathcal{U}f{)}^{t}g\mathrm{d}{\mathbf{x}}_{\mathrm{L}}={\int }_{\mathbb{S}}{f}^t({\mathcal{U}}^{t}g)\mathrm{d}{\mathbf{x}}_{\mathrm{L}}.$$

$${\int }_{\mathbb{S}}(\mathcal{U}f{)}^{†}g\mathrm{d}{\mathbf{x}}_{\mathrm{L}}={\int }_{\mathbb{S}}{f}^{†}({\mathcal{U}}^{†}g)\mathrm{d}{\mathbf{x}}_{\mathrm{L}}.$$

$$\displaystyle{{\mbox{\boldmath${\cal A}$}}}^{t}{\bf N}$$

$$\displaystyle{{\mbox{\boldmath${\cal A}$}}}^{\dagger}{\bf K}$$

$$\mathbf{N}=\left(\begin{array}{cc}0& 1\\ -1& 0\end{array}\right),\mathbf{K}=\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right).$$

$$j=\frac{1}{4}{\mathbf{q}}^{†}\mathbf{K}\mathbf{q}.$$

$$\displaystyle\int_{\mathbb{D}}\Bigl{(}\bigl{(}({{\mbox{\boldmath${\cal A}$}}}{\bf q}_{A})^{t}+{\bf d}_{A}^{t}\bigr{)}{\bf N}{\bf q}_{B}+{\bf q}_{A}^{t}{\bf N}\bigl{(}{{\mbox{\boldmath${\cal A}$}}}{\bf q}_{B}+{\bf d}_{B}\bigr{)}\Bigr{)}{\rm d}{\bf x}=\int_{\partial\mathbb{D}}{\bf q}_{A}^{t}{\bf N}{\bf q}_{B}n_{3}{\rm d}{\bf x}.$$

$$\displaystyle\int_{\mathbb{D}}\bigl{(}{\bf d}_{A}^{t}{\bf N}{\bf q}_{B}+{\bf q}_{A}^{t}{\bf N}{\bf d}_{B}\bigr{)}{\rm d}{\bf x}=\int_{{\partial\mathbb{D}}_{0,1}}{\bf q}_{A}^{t}{\bf N}{\bf q}_{B}n_{3}{\rm d}{\bf x}_{\rm L}.$$

$$\displaystyle\int_{\mathbb{D}}(-B_{A}P_{B}+P_{A}B_{B}){\rm d}{\bf x}=\int_{{\partial\mathbb{D}}_{0,1}}\frac{1}{i\omega\beta}\bigl{(}P_{A}\partial_{3}P_{B}-(\partial_{3}P_{A})P_{B}\bigr{)}n_{3}{\rm d}{\bf x}_{\rm L}.$$

$$\displaystyle\int_{\mathbb{D}}\bigl{(}{\bf d}_{A}^{\dagger}{\bf K}{\bf q}_{B}+{\bf q}_{A}^{\dagger}{\bf K}{\bf d}_{B}\bigr{)}{\rm d}{\bf x}=\int_{{\partial\mathbb{D}}_{0,1}}{\bf q}_{A}^{\dagger}{\bf K}{\bf q}_{B}n_{3}{\rm d}{\bf x}_{\rm L}.$$

$$\displaystyle\int_{\mathbb{D}}(B_{A}^{*}P_{B}+P_{A}^{*}B_{B}){\rm d}{\bf x}=\int_{{\partial\mathbb{D}}_{0,1}}\frac{1}{i\omega\beta}\bigl{(}P_{A}^{*}\partial_{3}P_{B}-(\partial_{3}P_{A})^{*}P_{B}\bigr{)}n_{3}{\rm d}{\bf x}_{\rm L}.$$

$$\displaystyle\frac{1}{4}\int_{\mathbb{D}}\bigl{(}{\bf d}^{\dagger}{\bf K}{\bf q}+{\bf q}^{\dagger}{\bf K}{\bf d}\bigr{)}{\rm d}{\bf x}=\frac{1}{4}\int_{{\partial\mathbb{D}}_{0,1}}{\bf q}^{\dagger}{\bf K}{\bf q}n_{3}{\rm d}{\bf x}_{\rm L}$$

$$\displaystyle\frac{1}{4}\int_{\mathbb{D}}(B^{*}P+P^{*}B){\rm d}{\bf x}=\frac{1}{4}\int_{{\partial\mathbb{D}}_{0,1}}\frac{1}{i\omega\beta}\bigl{(}P^{*}\partial_{3}P-(\partial_{3}P)^{*}P\bigr{)}n_{3}{\rm d}{\bf x}_{\rm L},$$

$${\mathcal{H}}_2$$

$${\mathcal{H}}_2=k_s^2+{\partial }_{\nu }{\partial }_{\nu },$$

$$k_s^2=k^2-\frac{3({\partial }_{\nu }\beta )({\partial }_{\nu }\beta )}{4{\beta }^2}+\frac{({\partial }_{\nu }{\partial }_{\nu }\beta )}{2\beta }.$$

$$\displaystyle{{\mbox{\boldmath${\cal A}$}}}=\begin{pmatrix}0&i\omega\beta\\ -\frac{1}{i\omega\sqrt{\beta}}{\cal H}_{2}\frac{1}{\sqrt{\beta}}&0\end{pmatrix}.$$

$$\displaystyle{{\mbox{\boldmath${\cal A}$}}}={{\mbox{\boldmath${\cal L}$}}}{{\mbox{\boldmath${\cal H}$}}}{{\mbox{\boldmath${\cal L}$}}}^{-1},$$

$$\mathcal{H}$$

$$\mathcal{L}$$

$$\displaystyle{{\mbox{\boldmath${\cal L}$}}}^{-1}$$

$$\omega \beta$$

Full text

I Introduction

In many wave propagation problems it is possible to define a preferred direction of propagation. For example, in ocean acoustics, waves propagate primarily in the horizontal direction in an acoustic wave guide, bounded by the water surface and the ocean bottom. Similarly, in communication engineering, microwaves or optical waves propagate as beams through electromagnetic or optical wave guides. Wave propagation in preferred directions is not restricted to wave guides. For example, in geophysical reflection imaging applications, seismic or electromagnetic waves propagate mainly in the vertical direction (downward and upward) through a laterally unbounded medium.

To account for propagation in preferred directions, the wave equation for the full wave field can be decomposed into a set of coupled equations for waves that propagate in opposite directions along the preferred axis (for example leftward and rightward in ocean acoustics, or downward and upward in reflection imaging). In the literature on electromagnetic wave propagation these oppositely propagating waves are often called “bidirectional beams” (Hoekstra97OQE, ; Stralen97OQE, ) whereas in the acoustic literature they are usually called “one-way wave fields” (Claerbout71GEO, ; Berkhout82Book, ; McCoy87WM, ; Holberg88GP, ; Fishman93RS, ). In this paper we use the latter terminology.

There is a vast amount of literature on the analytical and numerical aspects of one-way wave propagation (Fishman84JMP, ; Fishman84JMP2, ; Halpern88JASA, ; Weston89JMP, ; Fishman91WM, ; Grimbergen98GEO2, ). A discussion of this is beyond the scope of this paper. Instead, we concentrate on the choice of the decomposition operator and the consequences for reciprocity and representation theorems.

Decomposition of a wave field into one-way wave fields is not unique. In particular, the amplitudes of the one-way wave fields can be scaled in different ways. In this paper we distinguish between so-called “flux-normalised” and “field-normalised” one-way wave fields. The square of the amplitude of a flux-normalised one-way wave field is by definition the power-flux density (or, for quantum-mechanical waves, the probability-flux density) in the direction of preference. Field-normalised one-way wave fields, on the other hand, are scaled such that the sum of the two oppositely propagating components equals the full wave field. These two forms of normalisation have been briefly analysed by De Hoop (Hoop92PHD, ; Hoop96JMP, ). From this analysis it appeared that the operators for flux-normalised decomposition exhibit more symmetry than the operators for field-normalised decomposition. Exploiting the symmetry of the flux-normalised decomposition operators, the author derived reciprocity and representation theorems for flux-normalised one-way wave fields (Wapenaar96GJI1, ; Wapenaar96GJI2, ).

The first aim of this paper is to discuss flux-normalised versus field-normalised decomposition in a systematic way. In particular, it will be shown that reciprocity theorems for field-normalised one-way wave fields can be derived in a similar way as those for flux-normalised one-way wave fields, even though the operators for field-normalised decomposition exhibit less symmetry.

The second aim is to discuss representation theorems for field-normalised one-way wave fields in a systematic way. This discussion includes links to “classical” Kirchhoff-Helmholtz integrals for one-way wave fields as well as to recent single-sided representations for backward propagation, used for example in Marchenko imaging (Wapenaar2014GEO, ). Despite the links to earlier results, the discussed representations are more general. An advantage of the representations for field-normalised one-way wave fields is that a straightforward summation of the one-way wave fields gives the full wave field.

We restrict the discussion to scalar wave fields. In section II we formulate a unified scalar wave equation for acoustic waves, horizontally polarised shear waves, transverse electric and transverse magnetic EM waves and, finally, quantum-mechanical waves. Next, we reformulate the unified wave equation into a matrix-vector form, discuss symmetry properties of the operator matrix and use this to derive reciprocity theorems in matrix-vector form. In section III we decompose the matrix-vector wave equation into a coupled system of equations for oppositely propagating one-way wave fields. We separately consider flux-normalisation and field-normalisation and derive reciprocity theorems for one-way wave fields, using both normalisations. In section IV we extensively discuss representation theorems for field-normalised one-way wave fields and indicate applications. We end with conclusions in section V.

II Unified wave equation and its symmetry properties

II.1 Unified scalar wave equation

Using a unified notation, wave propagation in a lossless medium (or, for quantum-mechanical waves, in a lossless potential) is governed by the following two equations in the space-frequency domain

[TABLE]

Here $i$ is the imaginary unit and $\omega$ the angular frequency (in this paper we consider positive frequencies only). Operator $\partial_{j}$ stands for the spatial differential operator $\partial/\partial x_{j}$ and Einstein’s summation convention applies to repeated subscripts. $P({\bf x},\omega)$ and $Q_{j}({\bf x},\omega)$ are space- and frequency-dependent wave field quantities, $\alpha({\bf x})$ and $\beta({\bf x})$ are real-valued space-dependent parameters, and $B({\bf x},\omega)$ and $C_{j}({\bf x},\omega)$ are space- and frequency-dependent source distributions. Parameters $\alpha$ and $\beta$ are both assumed to be positive, hence, metamaterials are not considered in this paper. All quantities are specified in Table 1 for different wave phenomena and are discussed in more detail below. As indicated in the first column of Table 1, we consider 3D and 2D wave problems. For the 3D situation, ${\bf x}=(x_{1},x_{2},x_{3})$ is the 3D Cartesian coordinate vector and lower-case Latin subscripts take on the values 1, 2 and 3. For the 2D situation, ${\bf x}=(x_{1},x_{3})$ is the 2D Cartesian coordinate vector and lower-case Latin subscripts take on the values 1 and 3 only.

Table 1: Quantities in unified equations (1) and (2).

$P$ $Q_{j}$ $\alpha$ $\beta$ $B$ $C_{j}$

1. Acoustic waves (3D) $p$ $v_{j}$ $\kappa$ $\rho$ $q$ $f_{j}$

2. SH waves (2D) $v_{2}$ $-\tau_{2j}$ $\rho$ $\frac{1}{\mu}$ $f_{2}$ $2h_{2j}$

3. TE waves (2D) $E_{2}$ $-\epsilon_{2jk}H_{k}$ $\varepsilon$ $\mu$ $-J_{2}^{\rm e}$ $\epsilon_{2jk}J_{k}^{\rm m}$

4. TM waves (2D) $H_{2}$ $\epsilon_{2jk}E_{k}$ $\mu$ $\varepsilon$ $-J_{2}^{\rm m}$ $-\epsilon_{2jk}J_{k}^{\rm e}$

5. Quantum waves (3D) $\Psi$ $\frac{2\hbar}{mi}\partial_{j}\Psi$ $4-\frac{4V}{\hbar\omega}$ $\frac{m}{2\hbar\omega}$

The unified boundary conditions at an interface between two media with different parameters state that $P$ and $n_{j}Q_{j}$ are continuous over the interface. Here $n_{j}$ represents the components of the normal vector ${\bf n}=(n_{1},n_{2},n_{3})$ at the interface for the 3D situation, or ${\bf n}=(n_{1},n_{3})$ for the 2D situation.

We discuss the quantities in Table 1 in more detail. The quantities in row 1, associated to 3D acoustic wave propagation in a lossless fluid medium, are acoustic pressure $p({\bf x},\omega)$, particle velocity $v_{j}({\bf x},\omega)$, compressibility $\kappa({\bf x})$, mass density $\rho({\bf x})$, volume-injection rate density $q({\bf x},\omega)$ and external force density $f_{j}({\bf x},\omega)$. For 2D horizontally polarised shear waves in a lossless solid medium, we have in row 2 horizontal particle velocity $v_{2}({\bf x},\omega)$, shear stress $\tau_{2j}({\bf x},\omega)$, mass density $\rho({\bf x})$, shear modulus $\mu({\bf x})$, external horizontal force density $f_{2}({\bf x},\omega)$ and external shear deformation rate density $h_{2j}({\bf x},\omega)$. Rows 3 and 4 contain the quantities for 2D electromagnetic wave propagation, with TE standing for transverse electric and TM for transverse magnetic. The quantities are electric field strength $E_{k}({\bf x},\omega)$, magnetic field strength $H_{k}({\bf x},\omega)$, permittivity $\varepsilon({\bf x})$, permeability $\mu({\bf x})$, external electric current density $J_{k}^{\rm e}({\bf x},\omega)$ and external magnetic current density $J_{k}^{\rm m}({\bf x},\omega)$. Furthermore, $\epsilon_{ijk}$ is the alternating tensor (or Levi-Civita tensor), with $\epsilon_{123}=\epsilon_{312}=\epsilon_{231}=1$, $\epsilon_{213}=\epsilon_{321}=\epsilon_{132}=-1$, and all other components being zero. In row 5, the quantities related to 3D quantum-mechanical wave propagation are wave function $\Psi({\bf x},\omega)$, potential $V({\bf x})$, particle mass $m$ and $\hbar=h/2\pi$, with $h$ Planck’s constant.

By eliminating $Q_{j}$ from equations (1) and (2) we obtain the unified scalar wave equation

[TABLE]

with wave number $k$ defined via

[TABLE]

II.2 Unified wave equation in matrix-vector form

We define a configuration with a preferred direction and reorganise equations (1) and (2) accordingly.

Consider a 3D spatial domain $\mathbb{D}$, enclosed by surface ${\partial\mathbb{D}}$. This surface consists of two planar surfaces ${\partial\mathbb{D}}_{0}$ and ${\partial\mathbb{D}}_{1}$ perpendicular to the $x_{3}$-axis and a cylindrical surface ${\partial\mathbb{D}}_{\rm cyl}$ with its axis parallel to the $x_{3}$-axis, see Figure 1. The surfaces ${\partial\mathbb{D}}_{0}$ and ${\partial\mathbb{D}}_{1}$ are situated at $x_{3}=x_{3,0}$ and $x_{3}=x_{3,1}$, respectively, with $x_{3,1}>x_{3,0}$. In general these surfaces do not coincide with physical boundaries. Surface $\mathbb{S}$ in Figure 1 is a cross-section of $\mathbb{D}$ at arbitrary $x_{3}$. The parameters $\alpha({\bf x})$ and $\beta({\bf x})$ are piecewise continuous smoothly varying functions in $\mathbb{D}$, with discontinuous jumps only at interfaces that are perpendicular to the $x_{3}$-axis (hence, $P$ and $Q_{3}$ are continuous over the interfaces). In the lateral direction the domain $\mathbb{D}$ can be bounded or unbounded. When $\mathbb{D}$ is laterally bounded, the configuration in Figure 1 represents a wave guide. For this situation we assume that homogeneous Dirichlet or Neumann boundary conditions apply, i.e., $P=Q_{3}=0$ or $n_{\nu}\partial_{\nu}P=n_{\nu}\partial_{\nu}Q_{3}=0$ at ${\partial\mathbb{D}}_{\rm cyl}$, where lower-case Greek subscripts take on the values 1 and 2. When $\mathbb{D}$ is laterally unbounded (for example for reflection imaging applications), the cylindrical surface ${\partial\mathbb{D}}_{\rm cyl}$ has an infinite radius and we assume that $P$ and $Q_{3}$ have “sufficient decay” at infinity. For the 2D situation, the configuration is a cross-section of the 3D situation for $x_{2}=0$ and lower-case Greek subscripts take on the value 1 only.

We reorganise equations (1) and (2) into a matrix-vector wave equation which acknowledges the $x_{3}$-direction as the direction of preference. By eliminating the lateral components $Q_{1}$ and $Q_{2}$ (or, for 2D wave problems, the lateral component $Q_{1}$), we obtain (Corones75JMAA, ; Ursin83GEO, ; Fishman84JMP, ; Wapenaar89Book, ; Hoop96JMP, )

[TABLE]

where wave vector ${\bf q}$ and source vector ${\bf d}$ are defined as

[TABLE]

with

[TABLE]

and operator matrix ${\cal A}$ defined as

[TABLE]

with

[TABLE]

The notation in the right-hand side of equations (7) and (10) should be understood in the sense that differential operators act on all factors to the right of it. Hence, operator $\partial_{\nu}\frac{1}{\beta}\partial_{\nu}$, applied via equation (5) to $P$, stands for $\partial_{\nu}(\frac{1}{\beta}\partial_{\nu}P)$.

Note that the quantities contained in the wave vector ${\bf q}$ are continuous over interfaces perpendicular to the $x_{3}$-axis. Moreover, these quantities constitute the power-flux density (or, for quantum-mechanical waves, the probability-flux density) in the $x_{3}$-direction via

[TABLE]

where the asterisk denotes complex conjugation.

II.3 Symmetry properties of the operator matrix

We discuss the symmetry properties of the operator matrix ${\cal A}$. First, consider a general operator ${\cal U}$ (which can be a scalar or a matrix), containing space-dependent parameters and differential operators $\partial_{\nu}$. We introduce the transpose operator ${\cal U}^{t}$ via the following integral property

[TABLE]

Here ${\bf x}_{\rm L}$ is the lateral coordinate vector, with ${\bf x}_{\rm L}=(x_{1},x_{2})$ for 3D and ${\bf x}_{\rm L}=x_{1}$ for 2D wave problems. $\mathbb{S}$ denotes an integration surface perpendicular to the $x_{3}$-axis at arbitrary $x_{3}$, with edge $\partial\mathbb{S}$, see Figure 1. The quantities $f({\bf x}_{\rm L})$ and $g({\bf x}_{\rm L})$ are space-dependent test functions (scalars or vectors). When these functions are vectors, $f^{t}$ is the transpose of $f$; when they are scalars, $f^{t}$ is equal to $f$. When $\mathbb{S}$ is bounded, homogeneous Dirichlet or Neumann conditions are imposed at $\partial\mathbb{S}$. When $\mathbb{S}$ is unbounded, $\partial\mathbb{S}$ has an infinite radius and $f({\bf x}_{\rm L})$ and $g({\bf x}_{\rm L})$ are assumed to have sufficient decay along $\mathbb{S}$ towards infinity. Operator ${\cal U}$ is said to be symmetric when ${\cal U}^{t}={\cal U}$ and skew-symmetric when ${\cal U}^{t}=-{\cal U}$. For the special case that ${\cal U}=\partial_{\nu}$, equation (12) implies $\partial_{\nu}^{t}=-\partial_{\nu}$ (via integration by parts along $\mathbb{S}$). Hence, operator $\partial_{\nu}$ is skew-symmetric.

We introduce the adjoint operator ${\cal U}^{\dagger}$ (i.e., the complex conjugate transpose of ${\cal U}$) via the integral property

[TABLE]

When the test functions are vectors, $f^{\dagger}$ is the complex conjugate transpose of $f$; when they are scalars, $f^{\dagger}$ is the complex conjugate of $f$. Operator ${\cal U}$ is said to be Hermitian (or self-adjoint) when ${\cal U}^{\dagger}={\cal U}$ and skew-Hermitian when ${\cal U}^{\dagger}=-{\cal U}$. For the operators ${\cal A}_{12}$ and ${\cal A}_{21}$, defined in equations (9) and (10), we find ${\cal A}_{12}^{t}={\cal A}_{12}$, ${\cal A}_{21}^{t}={\cal A}_{21}$, ${\cal A}_{12}^{\dagger}=-{\cal A}_{12}$ and ${\cal A}_{21}^{\dagger}=-{\cal A}_{21}$. Hence, operators ${\cal A}_{12}$ and ${\cal A}_{21}$ are symmetric and skew-Hermitian. With these relations, we find for the operator matrix ${\cal A}$

[TABLE]

with

[TABLE]

Note that, using the expressions for ${\bf q}$ and ${\bf K}$ in equations (6) and (16), we can rewrite equation (11) for the power-flux density (or, for quantum-mechanical waves, the probability-flux density) as

[TABLE]

II.4 Reciprocity theorems

We derive reciprocity theorems between two independent solutions of wave equation (5) for the configuration of Figure 1. We consider two states $A$ and $B$, characterised by wave vectors ${\bf q}_{A}({\bf x},\omega)$ and ${\bf q}_{B}({\bf x},\omega)$, obeying wave equation (5), with source vectors ${\bf d}_{A}({\bf x},\omega)$ and ${\bf d}_{B}({\bf x},\omega)$. In domain $\mathbb{D}$, the parameters $\alpha$ and $\beta$, and hence the matrix operator ${\cal A}$, are chosen the same in the two states (outside ${\partial\mathbb{D}}$ they may be different in the two states). Consider the quantity \partial_{3}\bigl{(}{\bf q}_{A}^{t}{\bf N}{\bf q}_{B}\bigr{)} in domain $\mathbb{D}$. Applying the product rule for differentiation, using equation (5) for both states, integrating the result over $\mathbb{D}$ and applying the theorem of Gauss yields

[TABLE]

Here $n_{3}$ is the component parallel to the $x_{3}$-axis of the outward pointing normal vector on ${\partial\mathbb{D}}$, with $n_{3}=-1$ at ${\partial\mathbb{D}}_{0}$, $n_{3}=+1$ at ${\partial\mathbb{D}}_{1}$ and $n_{3}=0$ at ${\partial\mathbb{D}}_{\rm cyl}$, see Figure 1. In the following the integral on the right-hand side is restricted to the horizontal surfaces ${\partial\mathbb{D}}_{0}$ and ${\partial\mathbb{D}}_{1}$, which together are denoted by ${{\partial\mathbb{D}}_{0,1}}$. The integral on the left-hand side can be written as $\int_{\mathbb{D}}(\cdots){\rm d}{\bf x}=\int_{x_{3,0}}^{x_{3,1}}{\rm d}x_{3}\int_{\mathbb{S}}(\cdots){\rm d}{\bf x}_{\rm L}$. Using equation (12) for the integral along $\mathbb{S}$ and symmetry property (14), it follows that the two terms in equation (18) containing operator ${\cal A}$ cancel each other. Hence, we are left with

[TABLE]

This is a convolution-type reciprocity theorem (Fokkema93Book, ; Hoop95Book, ; Achenbach2003Book, ), because products like ${\bf q}_{A}^{t}({\bf x},\omega){\bf N}{\bf q}_{B}({\bf x},\omega)$ in the frequency domain correspond to convolutions in the time domain. A more familiar form is obtained by substituting the expressions for ${\bf q}$, ${\bf d}$ and ${\bf N}$ (equations 6 and 16), choosing $C_{j}=0$ and using equation (2) to eliminate $Q_{3}$, which gives

[TABLE]

Next, consider the quantity \partial_{3}\bigl{(}{\bf q}_{A}^{\dagger}{\bf K}{\bf q}_{B}\bigr{)} in domain $\mathbb{D}$. Following the same steps as above, using equations (13) and (15) instead of (12) and (14), we obtain

[TABLE]

This is a correlation-type reciprocity theorem (Bojarski83JASA, ), because products like ${\bf q}_{A}^{\dagger}({\bf x},\omega){\bf K}{\bf q}_{B}({\bf x},\omega)$ in the frequency domain correspond to correlations in the time domain. Substituting the expressions for ${\bf q}$, ${\bf d}$ and ${\bf K}$ and choosing $C_{j}=0$ yields the more familiar form

[TABLE]

We obtain a special case by choosing states $A$ and $B$ identical. Dropping the subscripts $A$ and $B$ in equations (21) and (22) and multiplying both sides of these equations by $\frac{1}{4}$, gives

[TABLE]

and

[TABLE]

respectively. These equations quantify conservation of power (or, for quantum-mechanical waves, probability).

III Decomposed wave equation and its symmetry properties

III.1 General decomposition of the matrix-vector wave equation

To facilitate the decomposition of the matrix-vector wave equation (equation 5), we recast the operator matrix ${\cal A}$ into a somewhat different form. To this end we introduce an operator ${\cal H}_{2}$, according to

[TABLE]

with operator ${\cal A}_{21}$ defined in equation (10) and wavenumber $k$ in equation (4). Operator ${\cal H}_{2}$ can be rewritten as a Helmholtz operator (Wapenaar89Book, ; Hoop92PHD, )

[TABLE]

with the scaled wavenumber $k_{s}$ defined as (Brekhovskikh60Book, )

[TABLE]

Note that ${\cal H}_{2}^{t}={\cal H}_{2}$ and ${\cal H}_{2}^{\dagger}={\cal H}_{2}$, hence operator ${\cal H}_{2}$ is symmetric and self-adjoint and its spectrum is real-valued (with positive and negative eigenvalues). Using equation (25), we rewrite operator matrix ${\cal A}$, defined in equation (8), as

[TABLE]

Next, we decompose this operator matrix as follows

[TABLE]

with

[TABLE]

Operators ${\cal H}_{1}$, ${\cal L}_{1}$ and ${\cal L}_{2}$ are pseudo-differential operators (Fishman84JMP, ; Fishman93RS, ; Wapenaar89Book, ; Hoop92PHD, ; Hoop94WM, ; Wapenaar96GJI1, ; Fishman2000JMP, ; Thomson2015GJI1, ; Thomson2015GJI2, ). The decomposition expressed by equation (29) is not unique, hence, different choices for operators ${\cal H}_{1}$, ${\cal L}_{1}$ and ${\cal L}_{2}$ are possible. We discuss two of these choices in detail in the next two sections. Here we derive some general relations that are independent of these choices.

By substituting equations (28) and (30) $-$ (32) into equation (29) we obtain the following relations

[TABLE]

We introduce a decomposed field vector ${\bf p}$ and a decomposed source vector ${\bf s}$ via

[TABLE]

where

[TABLE]

Substitution of equations (29), (35) and (36) into the matrix-vector wave equation (5) yields

[TABLE]

Substituting equations (30) $-$ (32) and (37) into equation (38) gives

[TABLE]

This is a system of coupled one-way wave equations. From the first term on the right-hand side it follows that the one-way wave fields $P^{+}$ and $P^{-}$ propagate in the positive and negative $x_{3}$-direction, respectively. The second term on the right-hand side accounts for coupling between $P^{+}$ and $P^{-}$. The last term on the right-hand side contains sources $S^{+}$ and $S^{-}$ which emit waves in the positive and negative $x_{3}$-direction, respectively.

We conclude this section by substituting equations (35) and (36) into equations (19), (21) and (23). Using equations (12) and (13) for the integration along the lateral coordinates this yields

[TABLE]

These equations form the basis for reciprocity theorems for the decomposed field and source vectors ${\bf p}$ and ${\bf s}$ in the next two sections.

III.2 Flux-normalised decomposition and reciprocity theorems

The first choice of operators ${\cal H}_{1}$, ${\cal L}_{1}$ and ${\cal L}_{2}$ obeying equations (33) and (34) is (Hoop92PHD, ; Hoop96JMP, ; Wapenaar96GJI1, )

[TABLE]

Operator ${\cal H}_{1}$, which is the square root of the Helmholtz operator ${\cal H}_{2}$, is commonly known as the square-root operator (Claerbout71GEO, ; Berkhout82Book, ; Fishman84JMP, ). Like the Helmholtz operator ${\cal H}_{2}$, the square-root operator ${\cal H}_{1}$ is a symmetric operator (Wapenaar96GJI1, ), hence ${\cal H}_{1}^{t}={\cal H}_{1}$. For the adjoint square-root operator we have ${\cal H}_{1}^{\dagger}=({\cal H}_{1}^{t})^{*}={\cal H}_{1}^{*}$. The spectrum of ${\cal H}_{1}$ is real-valued for propagating waves and imaginary-valued for evanescent waves. Hence, unlike the Helmholtz operator, the square-root operator is not self-adjoint. If we neglect evanescent waves, we may approximate the adjoint square-root operator as ${\cal H}_{1}^{\dagger}\approx{\cal H}_{1}$. Similar relations hold for the square root of the square-root operator and its inverse, hence \bigl{(}{\cal H}_{1}^{\pm\frac{1}{2}}\bigr{)}^{t}={\cal H}_{1}^{\pm\frac{1}{2}} and, neglecting evanescent waves, \bigl{(}{\cal H}_{1}^{\pm\frac{1}{2}}\bigr{)}^{\dagger}\approx{\cal H}_{1}^{\pm\frac{1}{2}}. From here onward we replace $\approx$ by $=$ when the only approximation is the negligence of evanescent waves. Using these symmetry relations for ${\cal H}_{1}$ and equations (16), (31), (44) and (45), we obtain

[TABLE]

and, neglecting evanescent waves,

[TABLE]

with

[TABLE]

Hence, equations (40) $-$ (42) simplify to

[TABLE]

By substituting the expressions for ${\bf p}$, ${\bf s}$, ${\bf N}$ and ${\bf J}$ (equations 37, 16 and 48), we obtain

[TABLE]

Note that, since the right-hand side of equation (54) is equal to the right-hand side of equation (24), it quantifies the power flux (or the probability-flux for quantum-mechanical waves) through the surface ${{\partial\mathbb{D}}_{0,1}}$. Therefore we call $P^{+}$ and $P^{-}$ flux-normalised one-way wave fields. Consequently, equations (52) and (53) are reciprocity theorems of the convolution type and correlation type, respectively, for flux-normalised one-way wave fields. These theorems have been derived previously (Wapenaar96GJI1, ) and have found applications in advanced wave field imaging methods for active and passive data (Wapenaar2004GJI, ; Kumar2006JGR, ; Fan2006GEO, ; Weglein2006GEO, ; Slob2009IEEE, ; Ruigrok2012GRL, ; Wapenaar2014JASA, ; Slob2014GEO, ; Ravasi2016GJI, ; Ravasi2017GEO, ; Elison2018GJI, ; Staring2018GEO, ).

III.3 Field-normalised decomposition and reciprocity theorems

The second choice of operators ${\cal H}_{1}$, ${\cal L}_{1}$ and ${\cal L}_{2}$ obeying equations (33) and (34) is (Wapenaar89Book, )

[TABLE]

Only the Helmholtz operator ${\cal H}_{2}$ is the same as in the previous section (it is defined in equation 26). The operators ${\cal H}_{1}$, ${\cal L}_{1}$ and ${\cal L}_{2}$ are different from those in the previous section, but for convenience we use the same symbols. Using {\bf q}={{\mbox{\boldmath{\cal L}}}}{\bf p} (equation 35) and equations (6), (31), (37) and (56), we find

[TABLE]

hence, $P^{+}$ and $P^{-}$ have the same physical dimension as the full field variable $P$ (which is defined in Table 1 for different wave phenomena). Therefore we call $P^{+}$ and $P^{-}$ field-normalised one-way wave fields (for convenience we use the same symbols as in the previous section).

The square-root operator ${\cal H}_{2}^{\frac{1}{2}}$ is symmetric, but ${\cal H}_{1}$ defined in equation (55) is not. From this equation it easily follows that ${\cal H}_{1}$, premultiplied by $\beta^{-1}$ is symmetric, hence

[TABLE]

and, neglecting evanescent waves,

[TABLE]

Using these symmetry relations for $\frac{1}{\beta}{\cal H}_{1}$ and equations (16), (31), (56) and (57), we obtain

[TABLE]

and, neglecting evanescent waves,

[TABLE]

Using this in equations (40) and (41) yields

[TABLE]

By substituting the expressions for ${\bf p}$, ${\bf s}$, ${\bf N}$ and ${\bf J}$ (equations 37, 16 and 48), using equations (12) and (13), we obtain

[TABLE]

We aim to remove the operator ${\cal H}_{1}$ from these equations. From equation (39) and (56) we obtain

[TABLE]

with ${\cal L}_{2}$ defined in equation (57). Assuming that in state $A$ the derivatives in the $x_{3}$-direction of the parameters $\alpha$ and $\beta$ at ${{{\partial\mathbb{D}}_{0,1}}}$ vanish and there are no sources at ${{\partial\mathbb{D}}_{0,1}}$, we find from equations (67) and (68)

[TABLE]

Below we use this to remove ${\cal H}_{1}$ from the right-hand sides of equations (65) and (66). Next, we aim to remove ${\cal H}_{1}$ from the left-hand sides of these equations. From {\bf s}={{\mbox{\boldmath{\cal L}}}}^{-1}{\bf d} (equation 36) and equations (6), (32), (37), (56) and (57), we find

[TABLE]

or

[TABLE]

We define new decomposed sources $B_{0}^{+}$ and $B_{0}^{-}$, according to

[TABLE]

Using equations (69) and (72) in the right- and left-hand sides of equations (65) and (66), we obtain

[TABLE]

Equations (LABEL:eq94) and (LABEL:eq95) are reciprocity theorems of the convolution type and correlation type, respectively, for field-normalised one-way wave fields. These theorems are modifications of previously obtained results (Berkhout89GEO, ; Wapenaar89GEO, ). The main modification is that we applied decomposition at both sides of the equations instead of at the right-hand sides only. Moreover, in the present derivation the condition for the validity of equation (69) is only imposed for state $A$. In the next section we use equations (LABEL:eq94) and (LABEL:eq95) to derive representation theorems for field-normalised one-way wave fields and we indicate applications.

IV Field-normalised representation theorems

IV.1 Green’s functions

Representation theorems are obtained by substituting Green’s functions in reciprocity theorems. Our aim is to introduce one-way Green’s functions, to be used in the reciprocity theorems for field-normalised one-way wave fields (equations LABEL:eq94 and LABEL:eq95). First, we introduce the full Green’s function $G({\bf x},{\bf x}_{A},\omega)$ as a solution of the unified wave equation (3) for a unit monopole point source at ${\bf x}_{A}$, with $B({\bf x},\omega)=\delta({\bf x}-{\bf x}_{A})$ and $C_{j}({\bf x},\omega)=0$. Hence,

[TABLE]

As boundary condition we impose the radiation condition (i.e., outward propagating waves at infinity). Next, we introduce one-way Green’s function as solutions of the coupled one-way equations (67) and (68) for a unit monopole point source at ${\bf x}_{A}$. Hence, we choose again $B({\bf x},\omega)=\delta({\bf x}-{\bf x}_{A})$ and $C_{j}({\bf x},\omega)=0$. Using equations (72) and (7), we define decomposed sources as $B_{0}^{\pm}=B^{\pm}=B=\pm 2{\cal L}_{2}S^{\pm}$, with ${\cal L}_{2}$ defined in equation (57), or

[TABLE]

We consider two sets of one-way Green’s functions. For the first set we choose a point source $S^{+}({\bf x},\omega)=\frac{1}{2}{\cal L}_{2}^{-1}B^{+}({\bf x},\omega)$, with $B^{+}({\bf x},\omega)=\delta({\bf x}-{\bf x}_{A})$, which emits waves from ${\bf x}_{A}$ in the positive $x_{3}$-direction, and we set $S^{-}({\bf x},\omega)$ equal to zero. Hence, for this first set, one-way equations (67) and (68) become

[TABLE]

Here $G^{\pm,+}$ stands for $G^{\pm,+}({\bf x},{\bf x}_{A},\omega)$. The second superscript ($+$) indicates that the source at ${\bf x}_{A}$ emits waves in the positive $x_{3}$-direction. The first superscript ($\pm$) denotes the propagation direction at ${\bf x}$. For the second set of one-way Green’s functions we choose a point source $S^{-}({\bf x},\omega)=-\frac{1}{2}{\cal L}_{2}^{-1}B^{-}({\bf x},\omega)$, with $B^{-}({\bf x},\omega)=\delta({\bf x}-{\bf x}_{A})$, which emits waves from ${\bf x}_{A}$ in the negative $x_{3}$-direction, and we set $S^{+}({\bf x},\omega)$ equal to zero. Hence, for this second set, one-way equations (67) and (68) become

[TABLE]

Here $G^{\pm,-}$ stands for $G^{\pm,-}({\bf x},{\bf x}_{A},\omega)$, with the second superscript ($-$) indicating that the source at ${\bf x}_{A}$ emits waves in the negative $x_{3}$-direction. Like for the full Green’s function $G({\bf x},{\bf x}_{A},\omega)$, we impose radiation conditions for both sets of one-way Green’s functions.

To find a relation between the full Green’s function and the one-way Green’s functions, we evaluate $\beta\partial_{3}\frac{1}{\beta}\partial_{3}(G^{+,+}+G^{-,+}+G^{+,-}+G^{-,-})$ using equations (77) $-$ (80), (25), (55) and (57). This gives equation (75), with $G$ replaced by $G^{+,+}+G^{-,+}+G^{+,-}+G^{-,-}$. Since the full Green’s function and the one-way Green’s functions obey the same radiation conditions, we thus find

[TABLE]

This very simple relation is a consequence of the field-normalised decomposition, introduced in section III.3.

IV.2 Source-receiver reciprocity

We derive source-receiver reciprocity relations for the field-normalised one-way Green’s functions introduced in the previous section. To this end we make use of the reciprocity theorem of the convolution type for field-normalised one-way wave fields, equation (LABEL:eq94). This theorem was derived for the configuration of Figure 1, assuming that in domain $\mathbb{D}$, the parameters $\alpha$ and $\beta$ are the same in the two states (see section II.4). Outside $\mathbb{D}$ these parameters may be different in the two states. For the Green’s state we choose the parameters for $x_{3}\leq x_{3,0}$ and for $x_{3}\geq x_{3,1}$ independent of the $x_{3}$-coordinate, according to $\alpha({\bf x}_{\rm L})$ and $\beta({\bf x}_{\rm L})$. Hence, if we let the Green’s state (with a point source at ${\bf x}_{A}$ in $\mathbb{D}$) take the role of state $A$, then the condition for the validity of equation (69) is fulfilled. Moreover, the Green’s functions are purely outward propagating at ${{\partial\mathbb{D}}_{0,1}}$ (because outside $\mathbb{D}$ no scattering occurs along the $x_{3}$-coordinate). Hence, $G^{+,\pm}({\bf x},{\bf x}_{A},\omega)=0$ at ${\partial\mathbb{D}}_{0}$ and $G^{-,\pm}({\bf x},{\bf x}_{A},\omega)=0$ at ${\partial\mathbb{D}}_{1}$. We let a second Green’s state (with a point source at ${\bf x}_{B}$ in $\mathbb{D}$ and the same parameters $\alpha$ and $\beta$ as in state $A$, inside as well as outside $\mathbb{D}$) take the role of state $B$. Hence, $G^{+,\pm}({\bf x},{\bf x}_{B},\omega)=0$ at ${\partial\mathbb{D}}_{0}$ and $G^{-,\pm}({\bf x},{\bf x}_{B},\omega)=0$ at ${\partial\mathbb{D}}_{1}$. With only outward propagating waves at ${{\partial\mathbb{D}}_{0,1}}$, the surface integral on the right-hand side of equation (LABEL:eq94) vanishes. Hence, taking into account that $B_{0}^{\pm}=B^{\pm}$ (since $C_{j}=0$), equation (LABEL:eq94) simplifies to

[TABLE]

First, we consider sources emitting waves in the positive $x_{3}$-direction in both Green’s states, hence $B_{A}^{+}=\delta({\bf x}-{\bf x}_{A})$, $B_{A}^{-}=0$, $P_{A}^{\pm}=G^{\pm,+}({\bf x},{\bf x}_{A},\omega)$, $B_{B}^{+}=\delta({\bf x}-{\bf x}_{B})$, $B_{B}^{-}=0$ and $P_{B}^{\pm}=G^{\pm,+}({\bf x},{\bf x}_{B},\omega)$. Substituting this into equation (82) yields

[TABLE]

see Figure 2(a). Next, we replace the source in state $B$ by one emitting waves in the negative $x_{3}$-direction, hence $B_{B}^{+}=0$, $B_{B}^{-}=\delta({\bf x}-{\bf x}_{B})$ and $P_{B}^{\pm}=G^{\pm,-}({\bf x},{\bf x}_{B},\omega)$. This gives

[TABLE]

see Figure 2(b). By replacing also the source in state $A$ by one emitting waves in the negative $x_{3}$-direction, according to $B_{A}^{+}=0$, $B_{A}^{-}=\delta({\bf x}-{\bf x}_{A})$ and $P_{A}^{\pm}=G^{\pm,-}({\bf x},{\bf x}_{A},\omega)$, we obtain

[TABLE]

see Figure 2(c). Finally, changing the source in state $B$ back to the one emitting waves in the positive $x_{3}$-direction yields

[TABLE]

see Figure 2(d).

Source-receiver reciprocity relations similar to equations (83) $-$ (86) were previously derived for flux-normalised one-way Green’s functions (Wapenaar96GJI2, ), except that two of those relations involve a change of sign when interchanging the source and the receiver. The absence of sign changes in equations (83) $-$ (86) is due to the definition of $B_{0}^{\pm}$ in equation (72). Moreover, unlike the flux-normalised reciprocity relations, the field-normalised source-receiver reciprocity relations of equations (83) $-$ (86) have a very straightforward relation with the well-known source-receiver reciprocity relation for the full Green’s function. By separately summing the left- and right-hand sides of equations (83) $-$ (86) and using equation (81), we simply obtain

[TABLE]

IV.3 Kirchhoff-Helmholtz integrals for forward propagation

We derive Kirchhoff-Helmholtz integrals of the convolution type for field-normalised one-way wave fields. For state $B$ we consider the decomposed actual field, with sources only outside $\mathbb{D}$, hence, $B_{0,B}^{\pm}=0$ in $\mathbb{D}$ and $P_{B}^{\pm}=P^{\pm}({\bf x},\omega)$. The parameters $\alpha$ and $\beta$ are the actual parameters inside as well as outside $\mathbb{D}$. For state $A$ we choose the Green’s state with a unit point source at ${\bf x}_{A}$ in $\mathbb{D}$. The parameters $\alpha$ and $\beta$ in $\mathbb{D}$ are the same as those in state $B$, but for $x_{3}\leq x_{3,0}$ and for $x_{3}\geq x_{3,1}$ they are chosen independent of the $x_{3}$-coordinate. Hence, the condition for the validity of equation (69) is fulfilled. First, we consider a source in state $A$ which emits waves in the positive $x_{3}$-direction, hence $B_{A}^{+}=\delta({\bf x}-{\bf x}_{A})$, $B_{A}^{-}=0$ and $P_{A}^{\pm}=G^{\pm,+}({\bf x},{\bf x}_{A},\omega)$. Substituting all this into equation (LABEL:eq94) (with $B_{0,A}^{\pm}=B_{A}^{\pm}$) gives

[TABLE]

Next, we replace the source in state $A$ by one which emits waves in the negative $x_{3}$-direction, hence $B_{A}^{+}=0$, $B_{A}^{-}=\delta({\bf x}-{\bf x}_{A})$ and $P_{A}^{\pm}=G^{\pm,-}({\bf x},{\bf x}_{A},\omega)$. Equation (LABEL:eq94) thus gives

[TABLE]

Recall that ${{\partial\mathbb{D}}_{0,1}}$ consists of ${\partial\mathbb{D}}_{0}$ (with $n_{3}=-1$) and ${\partial\mathbb{D}}_{1}$ (with $n_{3}=+1$), see Figure 1. Since $G^{+,\pm}({\bf x},{\bf x}_{A},\omega)=0$ at ${\partial\mathbb{D}}_{0}$ and $G^{-,\pm}({\bf x},{\bf x}_{A},\omega)=0$ at ${\partial\mathbb{D}}_{1}$ (because outside $\mathbb{D}$ no scattering occurs along the $x_{3}$-coordinate in state $A$), the first term under the integral in equations (88) and (89) gives a contribution only at ${\partial\mathbb{D}}_{1}$ and the second term only at ${\partial\mathbb{D}}_{0}$. Hence,

[TABLE]

Note that there is no contribution from $P^{-}({\bf x},\omega)$ at ${\partial\mathbb{D}}_{0}$ nor from $P^{+}({\bf x},\omega)$ at ${\partial\mathbb{D}}_{1}$, see Figure 3.

We conclude this section by considering a special case. Suppose the source of the actual field (state $B$) is located at ${\bf x}_{B}$ in the half-space $x_{3}<x_{3,0}$. Then, by taking $x_{3,1}\to\infty$, the field $P^{-}$ at ${\partial\mathbb{D}}_{1}$ vanishes. This leaves the single-sided representation

[TABLE]

Note that we included the source coordinate vector ${\bf x}_{B}$ in the argument list of $P^{\pm}({\bf x}_{A},{\bf x}_{B},\omega)$. This representation is an extension of a previously derived result (Berkhout89GEO, ), in which the fields were decomposed at ${\partial\mathbb{D}}_{0}$ but not at ${\bf x}_{A}$. It describes forward propagation of the one-way field $P^{+}({\bf x},{\bf x}_{B},\omega)$ from the surface ${\partial\mathbb{D}}_{0}$ to ${\bf x}_{A}$ (with ${\bf x}_{A}$ and ${\bf x}_{B}$ defined at opposite sides of ${\partial\mathbb{D}}_{0}$). In the following two sections we discuss representations for backward propagation of one-way wave fields.

IV.4 Kirchhoff-Helmholtz integrals for backward propagation (double-sided)

We derive Kirchhoff-Helmholtz integrals of the correlation type for field-normalised one-way wave fields. For state $B$ we consider the decomposed actual field, with a point source at ${\bf x}_{B}$ and source spectrum $s(\omega)$. The parameters $\alpha$ and $\beta$ are the actual parameters inside as well as outside $\mathbb{D}$. For state $A$ we choose the Green’s state with a unit point source at ${\bf x}_{A}$ in $\mathbb{D}$. The parameters $\alpha$ and $\beta$ in $\mathbb{D}$ are the same as those in state $B$, but for $x_{3}\leq x_{3,0}$ and for $x_{3}\geq x_{3,1}$ they are chosen independent of the $x_{3}$-coordinate. Hence, the condition for the validity of equation (69) is fulfilled. First, we consider sources emitting waves in the positive $x_{3}$-direction in both states, hence $B_{A}^{+}=\delta({\bf x}-{\bf x}_{A})$, $B_{A}^{-}=0$, $P_{A}^{\pm}=G^{\pm,+}({\bf x},{\bf x}_{A},\omega)$, $B_{B}^{+}=\delta({\bf x}-{\bf x}_{B})s(\omega)$, $B_{B}^{-}=0$ and $P_{B}^{\pm}=P^{\pm,+}({\bf x},{\bf x}_{B},\omega)$. Substituting this into equation (LABEL:eq95) (with $B_{0,A}^{\pm}=B_{A}^{\pm}$ and $B_{0,B}^{\pm}=B_{B}^{\pm}$) gives

[TABLE]

where $\chi$ is the characteristic function of the domain $\mathbb{D}$. It is defined as

[TABLE]

Since $G^{+,+}({\bf x},{\bf x}_{A},\omega)=0$ at ${\partial\mathbb{D}}_{0}$ and $G^{-,+}({\bf x},{\bf x}_{A},\omega)=0$ at ${\partial\mathbb{D}}_{1}$ (because outside $\mathbb{D}$ no scattering occurs along the $x_{3}$-coordinate in state $A$), the first term under the integral in equation (92) gives a contribution only at ${\partial\mathbb{D}}_{1}$ and the second term only at ${\partial\mathbb{D}}_{0}$. Hence,

[TABLE]

Next, we replace the source in state $B$ by one emitting waves in the negative $x_{3}$-direction, hence $B_{B}^{+}=0$, $B_{B}^{-}=\delta({\bf x}-{\bf x}_{B})s(\omega)$ and $P_{B}^{\pm}=P^{\pm,-}({\bf x},{\bf x}_{B},\omega)$. This gives

[TABLE]

By replacing also the source in state $A$ by one emitting waves in the negative $x_{3}$-direction, according to $B_{A}^{+}=0$, $B_{A}^{-}=\delta({\bf x}-{\bf x}_{A})$, $P_{A}^{\pm}=G^{\pm,-}({\bf x},{\bf x}_{A},\omega)$, we obtain

[TABLE]

Finally, changing the source in state $B$ back to the one emitting waves in the positive $x_{3}$-direction yields

[TABLE]

Equation (96) is an extension of a previously derived result (Wapenaar89GEO, ), in which the fields were decomposed at ${{\partial\mathbb{D}}_{0,1}}$ but not at ${\bf x}_{A}$ and ${\bf x}_{B}$. Equations (94), (95) and (97) are further variations. Equation (97) is visualised in Figure 4. Together, these equations describe backward propagation of the one-way wave fields $P^{-,\pm}({\bf x},{\bf x}_{B},\omega)$ from ${\partial\mathbb{D}}_{0}$ and $P^{+,\pm}({\bf x},{\bf x}_{B},\omega)$ from ${\partial\mathbb{D}}_{1}$ to ${\bf x}_{A}$. Except for some special cases, the integrals along ${\partial\mathbb{D}}_{1}$ do not vanish by taking $x_{3,1}\to\infty$. Hence, unlike the forward propagation representation (90), the double-sided backward propagation representations (94) $-$ (97) in general do not simplify to single-sided representations. In the next section we discuss an alternative method to derive single-sided representations for backward propagation.

We conclude this section by considering a special case. Suppose that in state $B$ the parameters $\alpha$ and $\beta$ are the same as in state $A$ not only in $\mathbb{D}$ but also outside $\mathbb{D}$. Then $P^{\pm,\pm}({\bf x},{\bf x}_{B},\omega)=G^{\pm,\pm}({\bf x},{\bf x}_{B},\omega)s(\omega)$ for all ${\bf x}$. Substituting this into representations (94) $-$ (97), summing the left- and right-hand sides of these representations separately and dividing both sides by $s(\omega)$, using equations (81) and (87) and assuming that ${\bf x}_{B}$ is located in $\mathbb{D}$, we obtain

[TABLE]

where the so-called homogeneous Green’s function $G_{\rm h}({\bf x}_{A},{\bf x}_{B},\omega)$ is defined as

[TABLE]

(with $\Re$ denoting the real part), and where $G^{\pm}({\bf x},{\bf x}_{A},\omega)=G^{\pm,+}({\bf x},{\bf x}_{A},\omega)+G^{\pm,-}({\bf x},{\bf x}_{A},\omega)$ (and a similar expression for $G^{\pm}({\bf x},{\bf x}_{B},\omega)$). Equation (IV.4) is akin to the well-known representation for the homogeneous Green’s function (Porter70JOSA, ; Oristaglio89IP, ), but with decomposed Green’s functions under the integrals. The simple relation between representations (94) $-$ (97) on the one hand and the homogeneous Green’s function representation (IV.4) on the other hand is a consequence of the field-normalised decomposition, introduced in section III.3.

IV.5 Kirchhoff-Helmholtz integrals for backward propagation (single-sided)

The complex-conjugated Green’s functions $\{\partial_{3}G^{\pm,\pm}({\bf x},{\bf x}_{A},\omega)\}^{*}$ under the integrals in equations (94) $-$ (97) can be seen as focusing functions, which focus the wave fields $P^{\pm,\pm}({\bf x},{\bf x}_{B},\omega)$ onto a focal point ${\bf x}_{A}$. However, this focusing process requires that these wave fields are available at two boundaries ${\partial\mathbb{D}}_{0}$ and ${\partial\mathbb{D}}_{1}$, enclosing the focal point ${\bf x}_{A}$. Here we discuss single-sided field-normalised focusing functions $f_{1}^{\pm}({\bf x},{\bf x}_{A},\omega)$ and we use these in modifications of reciprocity theorems (LABEL:eq94) and (LABEL:eq95) to derive single-sided Kirchhoff-Helmholtz integrals for backward propagation.

We start by defining a new domain $\mathbb{D}_{A}$, enclosed by two surfaces ${\partial\mathbb{D}}_{0}$ and ${\partial\mathbb{D}}_{A}$ perpendicular to the $x_{3}$-axis at $x_{3}=x_{3,0}$ and $x_{3}=x_{3,A}$, respectively, with $x_{3,A}>x_{3,0}$, see Figure 5. Hence, ${\partial\mathbb{D}}_{A}$ is chosen such that it contains the focal point ${\bf x}_{A}$. The two surfaces ${\partial\mathbb{D}}_{0}$ and ${\partial\mathbb{D}}_{A}$ are together denoted by ${{\partial\mathbb{D}}_{0,A}}$. The focusing functions $f_{1}^{\pm}({\bf x},{\bf x}_{A},\omega)$, which will play the role of state $A$ in the reciprocity theorems, obey the one-way wave equations (67) and (68) (but without the source terms $S^{\pm}$), with parameters $\alpha$ and $\beta$ in $\mathbb{D}_{A}$ equal to those in the actual state $B$, and independent of the $x_{3}$-coordinate for $x_{3}\leq x_{3,0}$ and for $x_{3}\geq x_{3,A}$. Hence, the condition for the validity of equation (69) is fulfilled. Analogous to equation (58), the field-normalised focusing functions $f_{1}^{\pm}({\bf x},{\bf x}_{A},\omega)$ are related to the full focusing function $f_{1}({\bf x},{\bf x}_{A},\omega)$, according to

[TABLE]

The focusing function $f_{1}^{+}({\bf x},{\bf x}_{A},\omega)$ is incident to the domain $\mathbb{D}_{A}$ from the half-space $x_{3}<x_{3,0}$ (see Figure 5). It propagates and scatters in the inhomogeneous domain $\mathbb{D}_{A}$, focuses at ${\bf x}_{A}$ on surface ${\partial\mathbb{D}}_{A}$ and continues as $f_{1}^{+}({\bf x},{\bf x}_{A},\omega)$ in the half-space $x_{3}>x_{3,A}$. The back-scattered field leaves $\mathbb{D}_{A}$ via surface ${\partial\mathbb{D}}_{0}$ and continues as $f_{1}^{-}({\bf x},{\bf x}_{A},\omega)$ in the half-space $x_{3}<x_{3,0}$. The focusing conditions at the focal plane ${\partial\mathbb{D}}_{A}$ are (Wapenaar2014GEO, )

[TABLE]

Here ${\bf x}_{{\rm L},A}$ denotes the lateral coordinates of ${\bf x}_{A}$. The operators $\partial_{3}$ and the factor $\frac{1}{2}i\omega\beta({\bf x}_{A})$ are not necessary to define the focusing conditions but are chosen for later convenience. To avoid instability, evanescent waves are excluded from the focusing functions. This implies that the delta function in equation (101) should be interpreted as a spatially band-limited delta function. Note that the sifting property of the delta function, $h({\bf x}_{{\rm L},A})=\int_{\mathbb{S}}\delta({\bf x}_{\rm L}-{\bf x}_{{\rm L},A})h({\bf x}_{\rm L}){\rm d}{\bf x}_{\rm L}$, remains valid for a spatially band-limited delta function, assuming $h({\bf x}_{\rm L})$ is also spatially band-limited.

We now derive single-sided Kirchhoff-Helmholtz integrals for backward propagation. We consider the reciprocity theorems for field-normalised one-way wave fields (equations LABEL:eq94 and LABEL:eq95), with $\mathbb{D}$ and ${{\partial\mathbb{D}}_{0,1}}$ replaced by $\mathbb{D}_{A}$ and ${{\partial\mathbb{D}}_{0,A}}$, respectively. For state $A$ we consider the focusing functions discussed above, hence, $B_{A}^{+}({\bf x},\omega)=B_{A}^{-}({\bf x},\omega)=0$ and $P_{A}^{\pm}({\bf x},\omega)=f_{1}^{\pm}({\bf x},{\bf x}_{A},\omega)$. For state $B$ we consider the decomposed actual field, with a point source at ${\bf x}_{B}$ in the half-space $x_{3}>x_{3,0}$ and source spectrum $s(\omega)$. The parameters $\alpha$ and $\beta$ in state $B$ are the actual parameters inside as well as outside ${{\partial\mathbb{D}}_{0,A}}$. First, we consider a source in state $B$ which emits waves in the positive $x_{3}$-direction, hence $B_{B}^{+}({\bf x},\omega)=\delta({\bf x}-{\bf x}_{B})s(\omega)$, $B_{B}^{-}({\bf x},\omega)=0$ and $P_{B}^{\pm}({\bf x},\omega)=P^{\pm,+}({\bf x},{\bf x}_{B},\omega)$. Substituting all this into equations (LABEL:eq94) and (LABEL:eq95) (with $B_{0}^{\pm}=B^{\pm}$), using equations (101) and (102) in the integrals along ${\partial\mathbb{D}}_{A}$, gives

[TABLE]

and

[TABLE]

where $\chi_{A}$ is the characteristic function of the domain $\mathbb{D}_{A}$. It is defined by equation (93), with $\mathbb{D}$ and ${{\partial\mathbb{D}}_{0,1}}$ replaced by $\mathbb{D}_{A}$ and ${{\partial\mathbb{D}}_{0,A}}$, respectively. Next, we replace the source in state $B$ by one which emits waves in the negative $x_{3}$-direction, hence $B_{B}^{+}({\bf x},\omega)=0$, $B_{B}^{-}({\bf x},\omega)=\delta({\bf x}-{\bf x}_{B})s(\omega)$ and $P_{B}^{\pm}({\bf x},\omega)=P^{\pm,-}({\bf x},{\bf x}_{B},\omega)$. This gives

[TABLE]

and

[TABLE]

Equations (103) $-$ (106) are single-sided representations for backward propagation of the one-way wave fields $P^{\pm,\pm}({\bf x},{\bf x}_{B},\omega)$ from ${\partial\mathbb{D}}_{0}$ to ${\bf x}_{A}$. Similar results have been previously obtained (Wapenaar2017GP2, ; Neut2017JASA, ), but without decomposition at ${\bf x}_{B}$. An advantage of these equations over equations (94) $-$ (97) is that the backward propagated fields $P^{\pm,\pm}({\bf x}_{A},{\bf x}_{B},\omega)$ are expressed entirely in terms of integrals along the surface ${\partial\mathbb{D}}_{0}$.

Single-sided representations containing the field-normalised focusing functions $f_{1}^{\pm}({\bf x},{\bf x}_{A},\omega)$ find applications for example in reflection imaging methods which account for multiple scattering. In these methods, the focusing functions are retrieved from the reflection response at the surface ${\partial\mathbb{D}}_{0}$, using the Marchenko method (Wapenaar2014GEO, ; Singh2017GEO2, ; Costa2018GJI, ; Wapenaar2019SE, ).

We conclude this section by considering a special case. Suppose that in state $B$ the parameters $\alpha$ and $\beta$ are the same as in state $A$ throughout space. Then $P^{\pm,\pm}({\bf x},{\bf x}_{B},\omega)=G^{\pm,\pm}({\bf x},{\bf x}_{B},\omega)s(\omega)$ for all ${\bf x}$. Moreover, $P^{+,\pm}({\bf x},{\bf x}_{B},\omega)=0$ for ${\bf x}$ at ${\partial\mathbb{D}}_{0}$. Substituting this into representations (103) $-$ (106), summing the left- and right-hand sides of these representations separately, dividing both sides by $s(\omega)$ and using equation (100), we obtain

[TABLE]

(with $\Im$ denoting the imaginary part), where $G^{-}({\bf x},{\bf x}_{B},\omega)=G^{-,+}({\bf x},{\bf x}_{B},\omega)+G^{-,-}({\bf x},{\bf x}_{B},\omega)$. Taking the real part of both sides gives

[TABLE]

where $G_{\rm h}({\bf x}_{A},{\bf x}_{B},\omega)$ is the homogeneous Green’s function, defined in equation (99). Unlike in equation (IV.4), here the homogeneous Green’s function is represented by a single integral along the surface ${\partial\mathbb{D}}_{0}$, containing field normalised one-way focusing and Green’s functions.

V Conclusions

We have considered flux-normalised and field-normalised decomposition of scalar wave fields into coupled one-way wave fields. The operators for field-normalised decomposition exhibit less symmetry than those for flux-normalised decomposition. Nevertheless, we have shown that reciprocity theorems can be derived for field-normalised one-way wave fields in a similar way as those for flux-normalised one-way wave fields. An additional condition for the reciprocity theorems for field-normalised one-way wave fields is that in one of the states the derivatives in the $x_{3}$-direction of the parameters $\alpha$ and $\beta$ vanish at the boundary of the considered domain. This condition is easily fulfilled when one of the states is a Green’s function or a focusing function, for which the parameters $\alpha$ and $\beta$ can be freely chosen at and outside the boundary of the domain.

We have used the reciprocity theorems for field-normalised one-way wave fields as a starting point for deriving representation theorems for field-normalised one-way wave fields in a systematic way. We obtained representations for forward and for backward propagation of one-way wave fields. These representations account for multiple scattering. Whereas the Kirchhoff-Helmholtz integrals for forward propagation can be easily transformed into single-sided representations, this transformation is less straightforward for the Kirchhoff-Helmholtz integrals for backward propagation. By replacing the Green’s functions by focusing functions we obtained single-sided representations for backward propagation of field-normalised one-way wave fields. These representations are particularly useful to retrieve wave fields in the interior of a domain in situations where measurements can be carried out only at a single surface. An important application is reflection imaging, accounting for multiple scattering.

Data Availability

No datasets have been used for this study.

Conflicts of Interest

The author declares that he has no competing interests.

Acknowledgements

I thank an anonymous reviewer for the constructive review, which helped to improve the paper. This work has received funding from the European Union’s Horizon 2020 research and innovation programme: European Research Council (grant agreement 742703).