On the synthesis of acoustic sources with controllable near fields

TL;DR

This paper introduces a method for designing acoustic sources that can generate specific near-field patterns in free space and ocean environments, supported by theoretical analysis and numerical simulations.

Contribution

It presents a novel strategy for synthesizing acoustic sources with controllable near fields, including multiple pattern approximation and low far-field signature techniques.

Findings

Successful approximation of prescribed near-field patterns

Effective control of multiple disjoint near-field regions

Low far-field signature while maintaining near-field accuracy

Abstract

In this paper we present a strategy for the the synthesis of acoustic sources with controllable near fields in free space and finite depth homogeneous ocean environments. We first present the theoretical results at the basis of our discussion and then, to illustrate our findings we focus on the following three particular examples: 1. acoustic source approximating a prescribed field pattern in a given bounded sub- region of its near field. 2. acoustic source approximating different prescribed field patterns in given disjoint bounded near field sub-regions. 3. acoustic source approximating a prescribed back-propagating field in a given bounded near field sub-region while maintaining a very low far field signature. For each of these three examples, we discuss the optimization scheme used to approx- imate their solutions and support our claims through relevant numerical simulations.

Equations

$$i)$$

$$\vspace{0.15cm}\left\{\vspace{0.15cm}\begin{array}[]{llll}\Delta u+k^{2}u=0\mbox{ in }{\mathbb{R}}^{3}\!\setminus\!\overline{D}_{a}\vspace{0.15cm},\\ \nabla u\cdot{\bf n}=v_{n},(\mbox{ or }u=p_{b})\mbox{ on }\partial D_{a}\vspace{0.15cm},\\ \displaystyle\left(\frac{{\bf x}}{|{\bf x}|},\nabla u({\bf x})\right)\!-\!iku({\bf x})\!=\!o\left(\frac{1}{|{\bf x}|}\!\right)\!,\mbox{ as }|{\bf x}|\rightarrow\infty\mbox{ uniformly for all $\displaystyle\frac{{\bf x}}{|{\bf x}|},$ }\end{array}\right.$$

$$u\approx u_1,\text{mbox}in{D}_1,\text{mbox}andu\approx 0,\text{mbox}in{D}_2,$$

$$v_n(\mathbf{x})=$$

$$p_b(\mathbf{x})=$$

$$u(\mathbf{x})={\int }_{\partial D_{a^{\prime }}}w(y)\frac{\partial \Phi (\mathbf{x},\mathbf{y})}{\partial {\mathbf{n}}_{\mathbf{y}}}d{s}_{\mathbf{y}},\text{mbox}for\mathbf{x}\in {\mathbb{R}}^3\setminus {\overline{D}}_a,$$

$$\vspace{0.15cm}\left\{\vspace{0.15cm}\begin{array}[]{llll}\Delta u+k^{2}u=0\mbox{ in }{\mathbb{R}}_{h}^{3}\!\setminus\!\overline{D}_{a}\vspace{0.15cm},\\ \nabla u\cdot{\bf n}=v_{n},(\mbox{ or }u=p_{b}),\mbox{ on }\partial D_{a}\vspace{0.15cm},\\ u\!=\!0\mbox{ on }z\!=\!0,\;\;\displaystyle\frac{\partial u}{\partial z}\!=\!0\mbox{ on }z\!=\!h,\\ \mbox{Radiation condition at infinity uniformly when }|\tilde{{\bf x}}|\rightarrow\infty,\end{array}\right.$$

$$u\approx u_1,\text{mbox}in{D}_1,\text{mbox}andu\approx 0,\text{mbox}in{D}_2,$$

$$u(\tilde{\mathbf{x}},z)=\sum _{n=0}^{\infty }{\varphi }_n(z){\psi }_n(\tilde{\mathbf{x}}),\text{mbox}forr>R,$$

$${\varphi }_n=\sin [k(1-a_n^2{)}^{\frac{1}{2}}z],\text{mbox}with{a}_n={\left[1-\frac{(2n+1{)}^2{\pi }^2}{4k^2{h}^2}\right]}^{\frac{1}{2}},$$

$$\underset{r\to \infty }{\lim }{r}^{\frac{1}{2}}\left(\frac{\partial {\psi }_n}{\partial r}-ik{a}_n{\psi }_n\right)=0,$$

$$v_n(\mathbf{x})=$$

$$p_b(\mathbf{x})=$$

$$G(z,\zeta ,|\tilde{\mathbf{x}}-\tilde{\mathbf{y}}|)=\Phi (\mathbf{x},\mathbf{y})+{\Phi }_1(z,\zeta ,|\tilde{\mathbf{x}}-\tilde{\mathbf{y}}|),$$

$$u(\mathbf{x})={\int }_{\partial D_{a^{\prime }}}w(\mathbf{y})\frac{\partial G(x,y)}{\partial {\mathbf{n}}_y}d{s}_{\mathbf{y}},\text{mbox}for\mathbf{x}\in {\mathbb{R}}^3\setminus {\overline{D}}_a,$$

$$\mathcal{D}{w}_{\alpha }(\mathbf{x})={\eta }_1{\int }_{\partial D_{a^{\prime }}}{w}_{\alpha }(\mathbf{y})\frac{\partial \Phi (\mathbf{x},\mathbf{y})}{\partial {\nu }_{\mathbf{y}}}d{S}_{\mathbf{y}}+i{\eta }_2{\int }_{\partial D_{a^{\prime }}}{w}_{\alpha }(\mathbf{y})\Phi (\mathbf{x},\mathbf{y})d{S}_{\mathbf{y}},$$

$$F(w)=\frac{1}{\Vert f{\Vert }_{L^2(\partial D_c)}^2}\Vert \mathcal{D}w-u_1{\Vert }_{L^2(\partial W_1)}^2+\mu \Vert \mathcal{D}w{\Vert }_{L^2(\partial W_2)}^2+\alpha \Vert w{\Vert }_{L^2(\partial D_{a^{\prime }})}^2,$$

$$\mu=\left\{\begin{array}[]{lll}0,&\mbox{ if }D_{2}=\emptyset,\\ 1,&\mbox{ if }D_{2}\mbox{ is bounded },\\ \frac{1}{4\pi R^{2}},&\mbox{ if }D_{2}={\mathbb{R}}^{3}\setminus B_{R}({\bf 0}),\end{array}\right.$$

$$D_1=\{(r,\theta ,\varphi ),r\in [0.011,0.015],\theta \in [-\frac{\pi }{4},\frac{\pi }{4}],\varphi \in [\frac{3\pi }{4},\frac{5\pi }{4}]\}.$$

Full text

acoustic source approximating a prescribed field pattern in a given bounded sub-region of its near field. 2. 2.

acoustic source approximating different prescribed field patterns in given disjoint bounded near field sub-regions. 3. 3.

acoustic source approximating a prescribed back-propagating field in a given bounded near field sub-region while maintaining a very low far field signature.

For each of these three examples, we discuss the optimization scheme used to approximate their solutions and support our claims through relevant numerical simulations.

1 Introduction and main results

The problem of active control of acoustic fields is well studied in the literature with a multitude of ideas and techniques presented (see monographs [1, 2]). The main strategies for active sound control are based on the use of boundary controls or secondary sources.

Applications of sound field control ideas are very important and they include: active noise control [5] (see also the pioneer works [3, 4]), acoustic field reproduction [7, 8, 9, 10, 11] and active control of scattered sound fields [12, 13, 14, 15, 16, 17, 18, 19, 20]. A rigorous comparative analysis of the theoretical similarities and respective challenges for these three areas of applications is done in [21].

In a recent development in [22] (see also [23] for the low frequency approximation), a general analytical approach based on the theory of boundary layer potentials was proposed for the active acoustic control problem in homogeneous environments. Then, in [24], building up on [22], the authors presented a thorough two dimensional sensitivity analysis for the synthesis of time-harmonic weak radiators with controllable patterns in some exterior region and proved that such acoustic sources will be feasible only if the region of control is in the reactive near-field of the source.

The work presented in this paper uses ideas from, and is relevant to, a wide array of important research areas: acoustic wave field synthesis, inverse source problems, optimization, personal audio techniques, acoustic near field control. We are making use of the theoretical results developed in [22] and, through a Tikhonov regularization procedure (with Morozov discrepancy principle for the choice of the regularization parameter), we synthesize acoustic sources in one of the following scenarios:

$1.$ Sources approximating a given pattern in a prescribed exterior near field sub-region.

$2.$ Sources approximating a given pattern in a prescribed sub-region of their near field while having a null in a different given sub-region of their near field.

$3.$ Sources which have a very weak field in a given (sufficiently far) exterior annuli while approximating a given pattern in a prescribed sub-region of their near field.

The first type of sources are relevant for the problem of acoustic rendering [6, 7]. The second type of sources above present an interest for the problem of personal audio studied in [9, 10, 11, 12] where we assume that by superposition our strategy will imply the possibility to approximate, with one source, different given sound patterns in disjoint regions of space. For the third type of sources above, although our strategy works well for the general question of synthesis of weak acoustic radiators approximating any given pattern in the near field control region, we focused on the problem of characterizing the necessary inputs (normal velocity or pressure) on the boundary of the source so that it approximates a backward propagating plane wave in the region of control while maintaining a very weak field in the given exterior annuli. This problem is relevant for the question of acoustic shielding or cloaking since by using a similar strategy we believe we can synthesize a planar array with similar properties: having a very weak field in an exterior annuli (where enemy detection measurements are taken) while approximating a given backward propagating plane wave in a near field region in front of it. Thus, by superposition, such an array could, when paired with a time control loop for the detection of interrogating signals, annihilate through destructive interference any incoming signal in its near field region without a large signature in its far field (i.e., shielding an object located behind the array). Then, a compact volume surrounded by a similar conformal array would lead to an active cloaking device for any object located inside.

The results presented in the literature regarding pattern synthesis use arrays of secondary sources (usually approximated by point sources) to control the field in interior regions (i.e., located in the interior of the geometric convex hull of the point sources), or focus on planar rendering (i.e. control in a horizontal plane) or assume that the field to be approximated propagate away from the source to be synthesized.

In the present paper we propose a theoretical optimization strategy for the synthesis of acoustic sources which approximate different prescribed field patterns in given disjoint * exterior * regions in free space and finite depth homogeneous ocean environments. To simplify the exposition, for the numerical support we consider only the case of sources in free space and focus on the three particular cases listed above.

The paper is organized as follows: In Section 2 we present the theoretical results in two parts: first, in Subsection 2.1 we first briefly recall the theoretical results of [22] for acoustic control in free space and then, in Subsection 2.2 we discuss their extension to the problem of underwater acoustic control in the context of a constant depth homogeneous ocean environments. In Section 3, we build up on our previous results in [24] and discuss the $L^{2}$ - Tikhonov regularization with Morozov discrepancy numerical approximation for the acoustic control problem in 3D and (assuming the superposition principle) without loosing the generality present numerical simulations in the three important situations listed above: first, in Subsection 3.1 we present the synthesis of an acoustic source approximating a prescribed plane wave in a give near field sub-region; then in Subsection 3.2 we present the synthesis of an acoustic source with a null in a given sub-region of its near field and approximating an outgoing plane wave in a disjoint near field sub-region; and finally, in Subsection 3.3 we synthesize a very weak acoustic radiator (almost non-radiating source (ANR)) approximating, in a sub-region of its near field, a given backward propagating (propagating towards the source) plane wave.

2 Theoretical results

In this section we will present the theoretical results behind our optimization scheme presented below. In Section 2.1 we will recall the results of [22] developed for the free space environments (i.e., homogeneous media with no boundaries and radiating condition at infinity) and then in Section 2.2 we will present their extension to the case of finite depth homogeneous ocean environments as introduced in [25, 26] (i.e., infinite rectangular waveguide with constant depth along z direction, $z\in[h,0]$ for some $h<0$, and pressure release boundary at the water-air interface $z=0$, total reflecting boundary at the ocean bottom interface $z=h$ together with radiation condition at infinity).

We consider the source support represented by $D_{a}$, a compact region of space with smooth enough boundary (Lipschitz boundary will be enough) and as in [22], assume that $D_{1}\Subset{\mathbb{R}}^{3}$ and $D_{2}\Subset{\mathbb{R}}^{3}$, with $D_{1}\cap D_{2}=\emptyset$ and $\{D_{1}\cup D_{2}\}\cap D_{a}=\emptyset$. Without losing the generality we assume $u_{1}$ is a solution of the Helmholtz equation in neighbourhoods of $D_{1}$ and $u_{2}=0$. With the these general hypotheses, in what follows the following three geometrical situations will be considered:

[TABLE]

At this point we mention that the theoretical results of [22] hold true for any finite number of mutually disjoint regions $D_{i}$, $i\in\{1,...,n\}$ satisfying $\{\cup_{i}D_{i}\}\cap D_{a}=\emptyset$ and respectively $n$ scalar acoustic fields $u_{i}$ in $D_{i}$.

In the case of free space environments the main question is to characterize boundary inputs (normal velocity or pressure) on the boundary of the source such that the acoustic field radiated by it has the property that it approximates $u_{1}$ in $D_{1}$ and $u_{2}=0$ in region $D_{2}$ respectively (the condition on $D_{2}$ is not needed in case (2) $i)$ above).

2.1 Acoustic control in free space

In this section we will recall the result obtained in [22] in the geometrical setting described at (2): A source approximating prescribed acoustic patterns $u_{1}$ and $u_{2}$ in two given disjoint exterior regions, $D_{1}$ and respectively $D_{2}$. Mathematically this can be written as follows:

**Problem 1. **

Find normal velocity ${v}_{n}$ (or pressure $p_{b}$) on the boundary of the antenna $\partial D_{a}$ so that,

[TABLE]

(where $\Delta$ denotes the 3D Laplace operator, $\nabla$ denotes the 3D gradient, $q=o(a)$ means $\displaystyle\lim_{a\rightarrow 0}q=0$, $\partial S$ denotes the boundary of the set $S\Subset{\mathbb{R}}^{3}$ and here ${\bf n}$ denotes the exterior normal to $\partial D_{a}$) and the following approximations hold true,

[TABLE]

where the approximation in (2.3) is in the sense of smooth norms (e.g., twice differentiable functions)

As a consequence of results in [22] we have that Problem 1 above can be answered in the affirmative in all the geometrical configurations described at (2). Indeed, if $k$ is not a resonance (i.e., for all wave numbers $k$ except a discrete set [22]) we have that there exists an infinite class of smooth functions $w$ (i.e., infinitely differentiable) so that normal velocity $v_{n}$ (or pressure $p_{b}$ ) given by,

[TABLE]

(where $\rho$ denotes the density of the surrounding medium, $D_{a}^{\prime}\Subset D_{a}$ is a smooth (i.e. with $C^{2}$ boundary) compact region, ${\bf n}_{\bf y}$ denotes the exterior normal to $\partial D_{a}^{\prime}$ computed in ${\bf y}\in\partial D_{a}^{\prime}$ and $\Phi$ is the free space fundamental solution of the Helmholtz equation), will generate the required acoustic field $u$ satisfying (2.2) and (2.3).

Remark 2.1**.**

Note that the fact that $D_{a}^{\prime}$ is smooth with $D_{a}^{\prime}\Subset D_{a}$ in (2.4) or (2.5) implies that the boundary input $v_{n}$ or $p_{b}$ is smooth on $\partial D_{a}$. Moreover, this permits us to assume minimal smoothness for the boundary of the actual physical source $\partial D_{a}$ (i.e., just enough to have the exterior problem well posed and thus Lipschitz will suffice) which may be very important for some applications.

Remark 2.2**.**

Observe that the normal velocity $v_{n}$ (or pressures $p_{b}$) defined at (2.4), (or (2.5)) generate a solution $u$ of (2.2), (2.3) represented as a double layer potential defined by

[TABLE]

but, it is elementary to see how the results of [22] can be extended to obtain solutions of (2.2), (2.3) represented as single layer potentials or a linear combination between double layer and single layer potentials.

2.2 Acoustic control in homogeneous oceans of constant depth

For the case of constant depth homogeneous ocean environments we model the surrounding homogeneous media as an infinite rectangular wave-guide, with constant depth $h$, i.e, $z\in[h,0]$ (where $z$ denotes the vertical coordinate in a rectangular coordinate system and $h<0$), and assume a pressure release condition at the water-air interface, i.e., zero pressure at $z=0$, and total pressure reflection at the bottom ocean interface, zero normal pressure at $z=h$ interface [25, 26].

Let ${\mathbb{R}}_{h}^{3}=\{{\bf x}=(\tilde{{\bf x}},z)\in{\mathbb{R}}^{3},\tilde{{\bf x}}\in{\mathbb{R}}^{2},h\leq z\leq 0\}$ and consider domains $D_{1}$, $D_{2}$ and functions $u_{1}$ and $u_{2}=0$ as in Section 2.1. Assuming cylindrical coordinates and using the same notations as in (2.2), (2.3) the problem can be formulated mathematically as follows:

**Problem 2. **

Find normal velocity ${v}_{n}$ (or the pressures $p_{b}$) on the boundary of the antenna $\partial D_{a}$ so that $u$, the solution of,

[TABLE]

satisfies

[TABLE]

where as above in Section 2.1 the condition on $D_{2}$ is not needed in the case when $D_{2}=\emptyset$.

We mention that, the radiation condition at infinity in problem (2.6) is understood as in [25, 26], i.e., for the solution $u$ represented in normal mode expansion

[TABLE]

(where $r=|\tilde{{\bf x}}|$) we have that

[TABLE]

and each of the $\psi_{n}$ satisfy the following radiation condition when $r\rightarrow\infty$

[TABLE]

uniformly for $\theta\in[0,2\pi]$ where $\theta=\tan^{-1}(\frac{x_{2}}{x_{1}})$. Next, we will describe how the results in [22] can be extended for this case. Indeed, as above, assuming that $k$ is not a resonant frequency, we have that there exists an infinite class of smooth functions $w$ so that normal velocity $v_{n}$ (or pressures $p_{b}$) given by,

[TABLE]

(where $\rho$ denotes the density of the surrounding media, $D_{a}^{\prime}\Subset D_{a}$ is a smooth region, ${\bf n}_{\bf y}$ denotes the exterior normal to $\partial D_{a}^{\prime}$ computed in ${\bf y}\in\partial D_{a}^{\prime}$ and $G$ is the Green’s function associated to problem (2.6)), will generate the required acoustic field $u$ satisfying (2.6) and (2.7). Indeed, it is observed in [26] that for ${\bf y}=(\zeta,\tilde{{\bf y}})$, with $\tilde{{\bf y}}=(y_{1},y_{2})$, the Green function $G$ associate to problem (2.6) is given by,

[TABLE]

where $\Phi$ is the fundamental free space solution of Helmholtz equation and $\Phi_{1}$ above is bounded and continuous at $z=\zeta$ and $\tilde{{\bf x}}=\tilde{{\bf y}}$. Based on these onsiderations it is concluded in [26] that the double layer and the single layer operators associated with the Green’s function $G$ have the same compactness properties and satisfy the same jump relations as the classical layer potentials associate to $\Phi$. Thus, by using this together with a few elementary technical adjustments it can be proved that the results presented in [22] will extend to this case, i.e., normal velocities (or pressures) given by (2.9) (or (2.10)) will generate acoustic fields described by double layer potentials associated to $G$ and satisfying (2.6) and (2.7). Moreover, we make the observation that the expressions (2.9), (2.10) can be used in computations since the Green’s function $G$ is computed explicitly in [26].

In this regard we mention that the statement of Remark 2.1 apply to the case of finite depth homogeneous oceans as well. The following remark is similar in spirit with Remark 2.2 but we presented here for the sake of completeness.

Remark 2.3**.**

The normal velocity $v_{n}$ (or the pressure $p_{b}$) given at (2.9) (or (2.10)) generates a solution $u$ of (2.6), (2.7) represented as a double layer potential defined by

[TABLE]

but, in a similar manner as above, the results of [22] could be easily extended to obtain solutions of (2.6), (2.7) represented as single layer potentials or a linear combination between double layer and single layer potentials.

3 Optimization schemes and Numerical simulations

In this section we describe the mathematical ideas behind the optimization scheme used towards the approximation of solutions to (2.2), (2.3) and respectively (2.6), (2.7).

The $L^{2}$- optimization and sensitivity analysis for the 2D formulation of the problem (2.2), (2.3) in the case (2) $iii)$ (with $\{D_{1}\cup D_{2}\}\cap D_{a}=\emptyset$ as above), was performed in [24] where it was numerically observed that good approximation of a stable solution with minimal power budget is achieved in the reactive near field of the source, i.e., when $D_{1}$ in (2.2) is located very close to the source $D_{a}$.

Similarly as in the 2D case treated in [24], the 3D $L^{2}$- optimization scheme for problem (2.2), (2.3) is based on Tikhonov regularization with Morozov discrepancy principle. In this context, as in [22], [24] regularity results and the well posedness of the interior and exterior acoustic boundary value problem (recall that $k$ was chosen non resonant) imply that in order to achieve approximate smooth controls in $D_{1}$ and $D_{2}$ it will be sufficient to have approximate $L^{2}$ controls on the boundaries of two slightly larger sets, $W_{1}\Subset{\mathbb{R}}^{3}$ and $W_{2}\Subset{\mathbb{R}}^{3}$, i.e., with $D_{1}\Subset W_{1}$ and $D_{2}\Subset W_{2}$. From (2.4) and Remark 2.2 it follows that solutions of (2.2), (2.3) can be approximated by a linear combination of double and single layer potentials, i.e.

[TABLE]

where $\eta_{1},\eta_{2}\in\mathbb{R}$ are fixed parameters and where $w_{\alpha}$ is the Tikhonov regularization solution, i.e., minimizer of the following discrepancy functional,

[TABLE]

with the regularization parameter $\alpha$ chosen according to the Morozov Discrepancy principle (see [24] for a 2D implementation in the case (2) $iii)$ and [27, 28] for the general theoretical discussion) and the weight $\mu$ above given by

[TABLE]

with $B_{R}({\bf 0})$ denoting the ball centred in the origin with radius $R$ such that $D_{a}\cup D_{1}\Subset B_{R}({\bf 0})$.

For the numerical simulations we make use of the spherical harmonic decomposition for $w_{\alpha}$ (the density of the layer potential operators used to represent the solution (3.11)) and through the method of moments and Tikhonov regularization we approximate a solution of the problem (2.2) and (2.3) in all the geometrical situations described above at (2). In this regard in all of the numerical simulations below we assumed 30 spherical harmonic orders in the spherical harmonic decomposition of $w_{\alpha}$.

In all the simulations the fictitious domain $D_{a^{\prime}}$ appearing in our strategy is the ball centered in the origin and radius $0.01$ and the actual physical source boundary $\partial D_{a}$ must be located in the annuli $0.01<r<dist(D_{a},D_{1})$ (where $dist$ here denotes the distance between the two sets) and for all the cases considered at (2) $D_{1}$ is given by

[TABLE]

In the remainder of the paper will present numerical simulations of our strategy and, to simplify the exposition, we will focus only on problem (2.2) and (2.3). Thus, the next three sections will present our numerical simulations for the Tikhonov regularization solution corresponding to problem (2.2) and (2.3) as follows: Section 3.1 discuses the case (2) $i)$, Section 3.2 discusses the case (2) $ii)$ while Section 3.3 discusses the case (2) $iii)$.

3.1 Synthesis of a prescribed pattern in a subregion of the source near-field

In this section we present the Tikhonov regularization solution for the problem (2.2), (2.3) introduced in Section 2.1 describing the applications to the synthesis of acoustic sources approximating a given field patterns in a near field bounded region $D_{1}$, see Figure 1.

Thus, we will show the performance of the Tikhonov solution described in (3.11), (3.12), (3.13) in the case (2) $i)$. Without loss of generality we consider the case when the source to be synthesized approximates in region $D_{1}$, described at (3.14), an outgoing plane wave propagating along the negative $x$-axis, $u_{1}=e^{-ixk}$ with wave number $k=10$.

First, in Figure 2 we present a cross-sectional view of the generated field along $z=0$ in a region characterized by $(x,y)\in[-5,5]^{2}$. This plot indicates the synthesised source causality (i.e., the fact that the source field is outgoing). This fact can also be observed in the time domain simulation presented in animation 1 where the propagating time-harmonic field generated by the synthesized source is shown.

In Figure 3 we present the quality of our control results in the region of interest $D_{1}$ as required in (2.3). The left and center plots in the figure describe respectively the field generated by the source, and the plane wave to be approximated $u_{1}=e^{-ixk}$ in region $D_{1}$. The good accuracy of our approximation $O(10^{-3})$ can be observed in the right picture of Figure 3 where the relative pointwise error between $u$ (the solution of (2.2)) and $u_{1}=e^{-ixk}$ (the field to be approximated) is presented.

Figure 4 shows six a cross-sectional views of the generated field along $z=0$ in a near-field region characterized by $(x,y)\in[-0.02,0.02]^{2}$. More explicitly, in order left to right from top left to bottom right plot, we present six cross-sectional ($z=0$) time-snapshots ( for $kct=\{\frac{15}{50}\pi,\frac{16}{50}\pi,\frac{17}{50}\pi,\frac{18}{50}\pi,\frac{19}{50}\pi,\frac{20}{50}\pi\}$) of the time-harmonic field generated by the synthesized source in it’s near field region, including the region of interest $D_{1}$ (where $c$ here was used as the speed of sound in air). The color scheme in the plots is (truncated to 1 light yellow and -1 dark blue) with the antenna region (colored cyan) not included in the numerical simulations and with the black stripe representing field amplitudes of $\approx 0.6$. Following the plots in order from top left to bottom right plot it can be observed how the source works to approximate an outgoing plane wave $u_{1}=e^{-ixk}e^{-ikct}$ in region $D_{1}$ (e.g., corresponding rectilinear black strip outgoing propagating through the control region). Indeed, the plots of Figure 4 show the propagation of the generated field by focusing on the portion of the field with amplitude $\approx 0.6$ marked as a dark stripe. It can be observed how this portion of the field enters region $D_{1}$ at time $kct=\frac{17}{50}\pi$ in a nearly rectilinear shape and continues keeping the same rectilinear form (indicating plane wave character of the approximated field in the control region) outgoing throughout a neighbourhood of region $D_{1}$.

The time domain animation

animation 2 presents the cross-sectional view along $z=0$ of the time-harmonic evolution of the field generated by the synthesized source and respectively the propagating plane wave $u_{1}=e^{-ixk}e^{-ikct}$ in a near field region given by $(x,y)\in[-0.1,0.1]^{2}$. The multimedia file shows two animations: the top one describing the time propagation of the generated field and the bottom one describing the time propagation of the plane wave $u_{1}=e^{-ixk}e^{-ikct}$. The color scheme in the movies is (truncated to 1 light yellow and -1 dark blue) with the antenna region removed from the simulations (colored cyan) and with the black stripe representing amplitude values around $\approx 0.6$, and the white stripe representing amplitude values around $\approx-0.6$ respectively. Observing that there will be two black stripes and respectively two white stripes per period for the approximated plane wave one can see in the animation the good accuracy of the approximation of the outgoing plane wave $u_{1}=e^{-ixk}e^{-ikct}$ in region $D_{1}$.

Figure 5 describes the density $w_{\alpha}$ (see (3.11)) on the boundary of the fictitious domain $D_{a^{\prime}}$ as an indication of the possible complexity of the required source inputs $v_{n}$ or $p_{b}$ described at (2.4) or (2.5). In the left plot of the figure we present the density values on the surface of $D_{a^{\prime}}$ viewed in a side 3D perspective and for better visualization we show two more plots in the figure: the center plot shows the density values on the part of the surface facing region $D_{1}$ while the right plot of the figure presents the density values on the oposite part of the surface.

3.2 Synthesis of different prescribed patterns in disjoint subregions of the source near-field

In this section we present the Tikhonov regularization solution for the problem (2.2), (2.3) described in Section 2.1 describing the applications to the synthesis of acoustic sources approximating two different field patterns in two prescribed disjoint near field regions. Thus, we show next the performance of the Tikhonov solution described in (3.11), (3.12), (3.13) in the case (2) $ii)$.

As in Section 3.1, we consider the case when the synthesized source approximates in region $D_{1}$, described at (3.14), an outgoing plane wave propagating along the negative $x_{1}$-axis, $u_{1}=e^{-ixk}$ with wave number $k=10$ while having a null in region $D_{2}=\{(r,\theta,\phi),r\in[0.011,0.015],\theta\in[-\frac{\pi}{4},\frac{\pi}{4}],\phi\in[-\frac{\pi}{4},\frac{\pi}{4}]\}+(0.018,0,0)$, (i.e. region $D_{2}$ is the same as region $D_{1}$ but shifted to the right along $x_{1}$ axis 0.018 units). A sketch of the geometries are presented with Figure 6.

In Figure 7 we present a cross-sectional view of the generated field along $z=0$ in a region characterized by $(x,y)\in[-5,5]^{2}$. As above, this plot supports the claim about the source causality (i.e., the fact that the source field is outgoing). This fact can be better observed in the time domain simulation presented in

animation 3 where the propagating time-harmonic field generated by the synthesized source is shown.

In Figure 8 we present the quality of our control results in the two regions of interest $D_{1},D_{2}$ as required in (2.3). The left and center plots on the top row in the figure describe respectively the field generated by the source, and the outgoing plane wave to be approximated $u_{1}=e^{-ixk}$. The good accuracy of our approximation $O(10^{-3})$ can be observed in the right picture on the top row of Figure 8 where the relative pointwise error between $u$ (the field generated solution of (2.2)) and $u_{1}=e^{-ixk}$ (the field to be approximated) is presented. The fourth picture in Figure 8 (bottom row of the figure) presents a scattered plot with the values of the generated field in region $D_{2}$ where very small values of the field (associated with the null effect as required from the optimization procedure) can be observed.

Figure 9, shows six a cross-sectional views of the generated field along $z=0$ in a near-field region characterized by $(x,y)\in[-0.1,0.1]^{2}$. More explicitly, in order left to right and from top left to bottom right plot, we present six cross-sectional ($z=0$) time-snapshots ($kct=\{\frac{83}{50}\pi,\frac{84}{50}\pi,\frac{85}{50}\pi,\frac{86}{50}\pi,\frac{87}{50}\pi,\frac{88}{50}\pi\}$) of the time-harmonic field generated by the synthesized source in it’s near field region, including the two regions of interest $D_{1},D_{2}$. The color scheme in the plots is (truncated to 1 light yellow and -1 dark blue) with the antenna region not included in the numerical simulations and thus corresponding to zero values field (cyan color) and with the black stripe representing field amplitudes of $\approx 0.6$. Following the plots in order from top left to bottom right plot it can be observed how the source works to approximate a plane wave corresponding to a straight black strip in region $D_{1}$ to its left while maintaining a null in region $D_{2}$ to its right. Indeed, the plots of Figure 9 show the time propagation of the generated field by focusing on the portion of the field with amplitude $\approx 0.6$ marked as a dark stripe. It can be observed how this portion of the field enters region $D_{1}$ at time $kct=\frac{85}{50}\pi$ in a nearly rectilinear shape and continues keeping the same form outgoing throughout a neighbourhood of region $D_{1}$ (thus indicating the good plane wave structure of the approximated field in the control region) while in all the plots the field is approximately zero in region $D_{2}$.

The time domain animation

animation 4 enhances the message of Figure 9. The multimedia file presents the cross-sectional view along $z=0$ of the time-harmonic evolution of the field generated by the synthesized source and respectively the propagating plane wave $u_{1}=e^{-ixk}e^{-ikct}$ in a near field region given by $(x,y)\in[-0.1,0.1]^{2}$ in two simultaneous animations: the top one describing the time propagation of the generated field and the bottom one describing the time propagation of the plane wave $u_{1}=e^{-ixk}e^{-ikct}$. The color scheme in the movies is (truncated to 1 light yellow and -1 dark blue) with the antenna region removed from the simulations (colored cyan) and with the black stripe representing amplitude values around $\approx 0.6$, and the white stripe representing amplitude values around $\approx-0.6$ respectively. As above, we point out that there are two black stripes and respectively two white stripes per period for the approximated plane wave. The animation clearly shows the good accuracy of the approximation in region $D_{1}$ as well as the null in region $D_{2}$.

Figure 10 describes the density $w_{\alpha}$ (see (3.11)) on the boundary of the fictitious domain $D_{a^{\prime}}$ as an indication of the possible complexity of the required source inputs $v_{n}$ or $p_{b}$ described at (2.4) or (2.5). In the left plot in the figure we present the density values on the surface of $D_{a^{\prime}}$ viewed from a 3D side perspective and for better visualization we show two more plots in the figure; the center plot shows the density values on the part of the surface facing region $D_{1}$; the right plot of the figure presents the density values on the part of the surface facing $D_{2}$.

3.3 Almost non-radiating acoustic sources with controllable near fields

In this section we present a third application of the results discussed in Section 2.1 and study the problem in the case (2) $iii)$. As an extreme example, we show how the source synthesised by our scheme approximates an incoming plane wave in region $D_{1}$ while having a very small field beyond radius $R=10$, (i.e. region $D_{2}$ in this case is the exterior of , $B_{10}({\bf 0})$, i.e., the ball centered at the origin with radius 10).

Thus, in what follows we consider the case when the source is required to approximate $u_{1}=e^{ixk}$ (plane wave propagating towards the source in the positive $x$ direction) with wave number $k=10$ in region $D_{1}$ described at (3.14) while having a very small field in region $D_{2}={\mathbb{R}}^{3}\!\setminus\!B_{10}({\bf 0})$. A sketch of the geometries are presented in Figure 11.

Figure 12 shows a cross-sectional view of the generated field along $z=0$ in a region characterized by $(x,y)\in[-5,5]^{2}$. This plot proves the source causality (i.e., the fact that the synthesized source field is outgoing). This fact can also be observed in the time domain simulation presented in

animation 5 where the propagating time-harmonic field generated by the synthesized source is shown.

In Figure 13 we show the quality of our control results in region $D_{1}$ as required by (2.3). The left and center plots in the figure describe respectively the field generated by the source $u$, and the plane wave to be approximated $u_{1}=e^{ixk}$. The good accuracy of our approximation ( ($O(10^{-3})$) can be observed in the right plot of Figure 13 where the relative pointwise error between the synthesised field $u$ and $u_{1}=e^{ixk}$ is presented.

Figure 14 shows the fast decay in region $D_{2}$ as required by (2.3). Indeed, the left plot of the figure describes the very small values of the generated field computed on the sphere of radius 10. On the other hand, the right plot in Figure 14 describes the absolute values of $\displaystyle r\cdot\sup_{B_{r}({\bf 0})}|u|$ as a function of $r\in(10,1000)$. The asymptotic limit of this function, $O(10^{-2})$, is the supremum value of far field pattern and this once more confirms the fact that the source synthesized by our scheme is a weak radiator. In fact, we also computed the actual power radiated by this source, i.e., $P=Re(\int_{S}u^{*}{\bf v}\cdot{\bf n})$ where $u$ represents the pressure field solution of (2.2), (2.3) , $A^{*}$ denotes the complex conjugate of complex quantity $A$, and ${\bf v}\cdot{\bf n}$ denotes the normal velocity on a sphere $S$ surrounding the source $D_{a}$, and we found that it is of order $O(10^{-7})$, once more indicating a very weak radiator.

Figure 15 shows nine a cross-sectional views of the generated field along $z=0$ in a near-field region characterized by $(x,y)\in[-0.1,0.1]^{2}$. The figure describes in order left to right from top left to bottom right plot, nine cross-sectional ($z=0$) time-snapshots ($kct=\{\frac{79}{50}\pi,\frac{80}{50}\pi,\frac{81}{50}\pi,\frac{81.5}{50}\pi,\frac{81.7}{50}\pi,\frac{82}{50}\pi,\frac{82.1}{50}\pi,\frac{83}{50}\pi,\frac{84}{50}\pi\}$) of the time-harmonic field generated by the synthesized source in it’s near field region. The color scheme in the plots is (truncated to 1 light yellow and -1 dark blue) with the antenna region (coloured cyan) not included in the numerical simulations and with the black stripe representing field amplitudes of $\approx 0.6$.

Following the plots in order from top left to bottom right plot it can be observed how the source works to approximate the incoming plane wave $u_{1}=e^{ixk}e^{-ikct}$ in $D_{1}$. Indeed, the $(a),(b)$ plots show how the source works on creating a plane wave in region $D_{1}$. The $(c),(d),(e),(f)$ plots are zoomed in closer to the antenna in a region $(x,y)\in[-0.02,0.02]^{2}$ so that the approximation of the incoming plane wave is better observed. It can be seen in these plots how the portion of the fields with values $\approx 0.6$ (black stripe) enters region $D_{1}$ at time $kct=\frac{81.5}{50}\pi$ in a nearly rectilinear shape (i.e. corresponding to plane wave character) and propagates towards the source while keeping the same rectilinear throughout a neighbourhood of region $D_{1}$. In the last two plots we see how the generated field looses form near the source and propagates away from it in the far field (i.e., corresponding to a causal source).

The time domain animation

animation 6 presents the cross-sectional view along $z=0$ of the time-harmonic evolution of the field generated by the synthesized source and respectively the incoming plane wave $u_{1}=e^{ixk}e^{-ikct}$ in a near field region given by $(x,y)\in[-0.1,0.1]^{2}$. The multimedia file shows two animations: the top one describing the time propagation of the generated field and the bottom one describing the time propagation of the plane wave $u_{1}=e^{ixk}e^{-ikct}$. The color scheme in the movies is (truncated to 1 light yellow and -1 dark blue) with the antenna region removed from the simulations and with the black stripe representing amplitude values around $\approx 0.6$, and the white stripe representing amplitude values around $\approx-0.6$ respectively. As above note that there will be two black stripes and respectively two white stripes per period for the approximated plane wave. The animations show the good accuracy of the approximation of the outgoing plane wave $u_{1}=e^{-ixk}e^{-ikct}$ in region $D_{1}$.

Figure 16 shows the density $w_{\alpha}$ (see (3.11)) on the boundary of the synthesized source. In the left plot in the figure we present the density values on the surface of the source viewed from a 3D side perspective and for better visualization we show two more plots in the figure; the center plot shows the density values on the part of the surface facing region $D_{1}$; the right plot of the figure presents the density values on the opposite pole of the source.

4 Conclusions and Future Work

In the time-harmonic regime we described a unified framework where the possibility to characterize boundary source inputs for the control of acoustic fields in homogeneous infinite media (Section 2.1) or finite depth homogeneous ocean environments (Section 2.2) was theoretically established. We then presented a 3D optimization scheme (Section 3) based on the Method of Moments and Tikhonov regularization with Morozov discrepancy principle for the numerical characterization of boundary inputs (normal velocity or pressures) necessary on the boundary of the source to obtain the desired control effects. We then numerically discussed the performance of our scheme for the problem in homogeneous infinite environments in the three distinct geometrical situations described at (2). We did not show any explicit numerical simulations for the case of finite depth homogeneous ocean environments but, as shown in Section 2.2, in this case the associated Green’s function is computed explicitly in [26] and is a continuous perturbation of the free space Helmholtz fundamental solution, hence the optimization scheme presented in the paper can be adapted to this case as well and this will be a part of our future work.

We believe that the results presented in this report may be relevant to: the problem of acoustic rendering or covert communications (Section 3.1); the problem of generating personal audio spots (Section 3.2); the question of acoustic protection where, paired with a time control loop for the detection of interrogating signals, a planar (or conformal) array of source elements similar to the one described Section 3.3 is synthesized as a very weak radiator (thus unobservable to a far field measurement device) with a controllable near-field such that it nulls (by destructive interference in a region to the left of the array) an interrogating field coming from the left thus protecting a region behind the array, i.e., to its right side, (or within its convex hull if a conformal array is used).

We mention also that the actual physical source boundary $\partial D_{a}$ does not need to be very smooth in general (Lipschitz suffice) and it must exists between the fictitious domain $D_{a^{\prime}}$ and the boundary of $D_{1}$. In this regard, our investigations suggest also that, assuming more harmonics in the expansion of the density $w_{\alpha}$ we could achieve the same degree of control with region $D_{1}$ located further away from the source (giving thus more freedom in the choice of an actual physical boundary $\partial D_{a}$) or when one considers more then two regions of control. In this regard, in the spirit of [24], a detailed study of the sensitivity of the optimization scheme with respect to parameters such as, relative position of the control regions $D_{1}$ and $D_{2}$ and their distance from the source $D_{a}$, power budget and oscillatory character of the source input as well as acoustic intensity of the source will be presented in forthcoming reports.

We observe that Figure 5, Figure 10 and Figure 16 indicate the fact that, for all the cases considered in the paper, the synthesized source requires a very complex input on its boundary, i.e., with sub-areas characterized by small values and fast oscillations (e.g. , the part facing $D_{1}$ ) and other sub-areas characterized by very large values and slower oscillations (e.g., the pole opposite to $D_{1}$). In the context of linear approximation where everything can be scaled down appropriately, the results presented above in Figures 8, 10, 13 and 16 corroborated with the superposition principle suggest the possibility of approximating arbitrary given patterns (with small amplitudes) in each of the regions of interest. In this context, an investigation of similar optimization methods minimizing the discrepancy functional introduced at (3.12) but where instead of the penalty term $\alpha||w||_{L^{2}}^{2}$ one penalizes for the $TV$ or $L^{1}$ norm of the normal velocity (or pressure) ((2.4) or (2.5) for free space and (2.9) or (2.10) for homogeneous ocean environments), will be considered in a forthcoming report.

Last but not least, we believe that by using the superposition principle our strategy can be extended to the time-domain where different prescribed broadband signals could in principle be synthesised in disjoint near field regions. This, and, if needed, a direct time-domain analysis will be another important part of our future work.

Acknowledgment

This work was supported by the Office of Naval Research under the award N00014-15-1-2462.