{'fraction_non_alphanumeric': 0.052799450414642526, 'fraction_numerical': 0.022915746601894108, 'mean_word_length': 3.8212286097310937, 'pattern_counts': {'":': 0, '<': 14, '<?xml version=': 0, '>': 15, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 8, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Silence is bliss in a Platonic relationship

We describe an effective active cloaking strategy for the scalar Helmholtz equation in three dimensions where multipole active sources are located at the vertices of the Platonic solids. A"silent zone"is created interior to the imaginary Platonic solid and only the incident field remains in a defined region exterior to the silent zone and active source configuration. This distribution of the sources ensures that the implementation of the cloaking strategy is extremely efficient. In particular, once the multipole source amplitudes required at a single source location are determined, the other source amplitudes can be calculated by simple post-processing involving multiplication of the multipole source vector by a rotation matrix. The general nature of the problem means that the technique is relevant to any scalar wave field, including both acoustics and electromagnetism.

Over the last two decades significant excitement has been generated around the idea of cloaking objects, i.e. rendering them invisible to incident wave fields [1][2][3][4][5]. Cloaking strategies are either passive or active, the former requiring the fabrication of metamaterials that are able to manipulate the wave, steering it around the object regardless of the form of the incident field [6][7][8][9]. In the context of active cloaking, sources are employed to suppress a field within a certain domain [2,10,11]. Indeed, using active sources to modify fields has been of interest in science, technology and engineering for many decades: Paul Lueg's patent in 1936 on the topic of antisound illustrates the use of active sources to suppress noise in specific regions of space [12]. Anti-sound and anti-vibration are huge areas of research and a plethora of techniques have been developed in order to address numerous problems in acoustic engineering [13][14][15]. An illustration of the relation between various methods in acoustics was provided by Cheer [16]. Specific choices of sources can ensure quiet zones and illusions [11,17] and recent work has optimised control sources to reduce the scattered field in passive scenarios and cases with flow for the cloaking of specific objects [18][19][20][21].
Over the last decade or so, interest has centred on active cloaking methods that are independent of the object to be cloaked. Miller's method of cloaking measured particle motion close to the surface of the cloaking region and simultaneously excited necessary surface sources where each source amplitude depends on the measurements at all sensing points [2]. This approach is rather limited as an active cloaking method because it cannot provide a relationship between the incident field and source amplitudes. Guevara-Vasquez et al. addressed this problem [10], describing an active source method that could use multipoles to create silent zones and simultaneously ensure that the active field is zero outside the silent zone. The integral equation formulation of the problem was employed, and subsequently converted to a linear system of equations for the source amplitudes. These amplitudes are linear functions of the incident field and it was shown, by construction, that active cloaking could be realised in two dimensions. An explicit form of the relations between the sources and incident field was provided in [22]. In particular multipole sources were used, Miller's cloak was reproduced and numerical results were compared with SVD solutions of the linearised system [10,23].
Further progress was made on the form of the source coefficients in [24], where it was demonstrated that the integral representation of Vasquez et al. [10] could be reduced to closed-form explicit formulas, which bypassed the requirement for a numerical solution. Analytical expressions were provided for general incidence and specifically for plane-wave incidence. This approach was then extended to two-dimensional (vector) elastodynamics [25] and to the thin-plate equation [26]. The attractive nature of this approach to active cloaking is that the source coefficients are independent of the object to be cloaked inside the silent zone. This is a benefit over approaches that are object-dependent where alternative approaches have to be taken in resonant regimes [27]. The latter approach however, does not suffer from large amplitudes in the vicinity of the active source regions.
Although there has been extensive work regarding the active manipulation of sound in three-dimensions in general settings, see e.g. [28][29][30][31], the active cloaking approach that is independent of the object to be cloaked has thus far predominantly been conducted in the twodimensional setting, with the sole exception of the work of Guevara-Vasquez et al [32] who described the extension of the methods in [10,23] for the scalar Helmholtz equation in three dimensions. Expressions were derived for the active field in terms of an integral of the incident field and subsequently illustrated with four active sources.
Here we deduce integral expressions for active cloaking in three dimensions. We provide expressions for the active source coefficients and introduce a new methodology associated with distributions of the active sources on the vertices of the Platonic solids. This approach is fast and efficient to implement due to the sources residing on the circumsphere of the Platonic solid. Regardless Figure 1. Schematic of the three-dimensional active exterior cloaking problem. The incident plane wave ui propagating in the direction k impinges upon an arbitrary object (the orange sphere) and is scattered in the form of the field us. The amplitudes of the L active sources (the smaller spheres) are determined such that they produce an active field u d which cancels out the incident field ui in a volume containing the object and thus nullifies the scattered field us.
of the form of the incident field, once the form of multipole source coefficients has been deduced for just one of the active sources, the remaining source amplitudes can be calculated merely by solving a matrix equation in terms of the determined multipole source vector, thus providing a clear mechanism and complete exposition of three-dimensional active cloaking for the scalar Helmholtz equation.
We consider active exterior cloaking for time-harmonic waves (with dependence e −iωt where ω is the angular frequency and t is time) governed by the three-dimensional homogeneous Helmholtz equation (∇ 2 + k 2 )u(x) = 0 where k = ω/c is the wavenumber with c the speed of wave. The scalar wave field of interest is u(x), which in the context of acoustics is the velocity potential at the position x = (x, y, z).
As illustrated in Fig. 1 for the case of an incident field u i (x) that is planar, we introduce the active field u d (x) generated by L multipole sources. The field subsequently scattered from an object interior to the active field is denoted by u s (x). These fields are conveniently described in terms of spherical waves: where U m n (x) = j n (k|x|)Y m n ( x) and V m n (x) = h (1) n (k|x|)Y m n ( x) are incoming and outgoing spherical waves respectively [33], noting that j n (k|x|) is the spherical Bessel function of the first kind and h (1) n (k|x|) is the spherical Hankel function of the first kind. Further, Y m n ( x) = A m n P m n (cos θ)e imϕ is the normalized spherical harmonic function and where hats denote unit vectors and P m n (cos θ) is the associated Legendre function with order m and degree n in terms of the polar angle θ and the azimuthal angle ϕ.
Referring to Fig. 1, we write k = k k where k = (sin θ i cos ϕ i , sin θ i sin ϕ i , cos θ i ) is the unit propagating vector. For plane wave incidence, the coefficients Q nm = 4πi n Y m n ( k), thus depending only on the incident wave angles θ i , ϕ i with the overline denoting complex conjugate. In order to achieve active cloaking, source coefficients q ,nm are sought such that for some closed domain C surrounded by the active sources at x / ∈ C where = 1, 2, · · · , L, we have u i (x) + u d (x) = 0 for x ∈ C and u d (x) → 0 as |x| → ∞. We require by the first condition that the active field u d interferes destructively with the incoming wave u i such that the total field vanishes in the region C. While it leads to the nullification of wave scattering from any object in C, we also stipulate by the second condition that the radiation of u d itself to the far field is minimized. This property will leave minimal evidence of the cloak to be detected.
Consider now that the active sources are located at the vertices of an imaginary Platonic solid. This arrangement will limit the total number of sources L to five values [34] as illustrated in Fig. 2(a) -(e). The geometry ensures that the sources reside on the circumsphere of the Platonic solids such that |x | = x 0 for all in each case, with x 0 an arbitrary constant. We further set x 1 = (0, 0, −x 0 ) in every case such that the active source with = 1 is always at the lowest position of the distribution in terms of the z coordinate. To locate the remaining sources, we note that a Platonic solid consisting of q p-sided regular polygonal faces around each vertex can be characterized by a set of two indices (p, q). (For example, a regular tetrahedron with three equilateral triangles around each vertex has the indices (3,3).) It is useful to define the length of each side of a Platonic solid in terms of x 0 as sx 0 . We will show in Part 2 of Supplementary Material that the position of each active source can be fully determined using knowledge of (p, q) and s, whose values for all the five cases are listed in Table I. In particular, since the source distribution is q-fold rotationally symmetric about the z axis, if we assign the indices = 2, · · · , q + 1 to the q vertices located immediately above the source = 1 in the counterclockwise direction, then their position vectors x are given by where z 2 = (s 2 − 2)/2, r 2 = 1 − z 2 2 and ϕ = 2( − 2)π/q. Further details about the indexing of sources with 1 < ≤ L can be found also in Part 2 of Supplementary Material. Table I. The geometric properties of each Platonic distribution of active sources and the volumes of the respective cloaked regions V = |C| (corrected to four decimal places) when the sphere radius a takes the lower limit in (10).
To determine the source amplitudes q ,nm and the precise geometric shape of the cloaked region C, we employ the procedure in [32] as our starting point and then use the approach in [24] to obtain the convenient explicit expressions (see Part 1 of Supplementary Material) in the form where we recall that Q ts are the coefficients associated with the expansion of the incident field (1) and In (7), S sµ tν (x ) (defined in (27) - (30) of Supplementary Material) is a coefficient depending on the position vector x and the quantities D νn , I mµ nν (x , ∂C ) take the respective forms where prime denotes the derivative with respect to the argument and a is any value within the range which is derived in Part 2 of Supplementary Material. In (9), the surface integral I mµ nν is performed over the face ∂C parameterized by the vector y − x for y ∈ ∂C . (To the authors' knowledge, an analytic form of I mµ nν is not available at this point. However, the integral can be simplified such that its numerical evaluation becomes significantly less expensive, which will be discussed in Part 4 of Supplementary Material.) The formulations (6) -(10) hold under two conditions: 1) the domain C is completely bounded by the union of faces ∂C where = 1, 2, · · · , L; 2) |y − x | = ax 0 such that ∂C is a surface belonging to a sphere centred at x with radius ax 0 . The cloaked region C therefore consists of the domain interior to the union of surfaces formed by identical imaginary spheres of radius ax 0 located at the vertices. Although active cloaking for a plane wave u i (x) = e ik·x is considered in particular, the source amplitude in the form (6) -(10) applies to a general incident wave since the expansion coefficient Q nm suffices to fully describe the nature of u i by (1). Nonetheless, the case of plane wave incidence where Q nm = 4πi n Y m n ( k) admits a more compact form of q ,nm as (6) is reducible to Note that for computational purpose, the infinite series in (11) needs to be truncated to a finite order. We discuss further in Part 3 of Supplementary Material how to choose this truncation parameter such that the source amplitude converges with a prescribed level of accuracy. We show in Fig. 2(a) -(e) the respective cloaked region C inside the Platonic solid for each source distribution, with a taken as the lower bound in (10). (In this case, the volume V = |C| of the domain becomes the maximum.) In Fig. 2(f), a cross section through the plane z = 0 for Fig. 2(b) is also illustrated. Note that in every case each surface of integration ∂C is delimited by a total of q identical circular arcs, where q depends on L as indicated in Table I. For non-Platonic distributions of sources these arcs would not be identical.
We will first evaluate the surface integral I mµ nν (x 1 , ∂C 1 ) for the active source = 1. A parametric form for the corresponding spherical face ∂C 1 can be conveniently derived using a local spherical coordinate system centred at the position x 1 (see Part 2 of Supplementary Material). Since ∂C 1 possesses a q-fold rotational symmetry about the z axis, it can be subdivided into q congruent segments, each being the region bounded by where and ϕ t , z 2 are defined as in (5). The entire surface ∂C 1 is given by the union of segments parameterized by (12) -(15) over t = 2, · · · , q + 1. Since all surfaces of integration ∂C have an identical geometric shape within each Platonic distribution of sources, this parametric form describes not only the face ∂C 1 but also the rest as long as a suitable coordinate transformation is applied to account for their different orientations. We may also compute the volumes of the cloaked region C using the divergence theorem on (12) - (15). In Table I we list the volume V for all the five cases illustrated in Fig. 2. It is observed that amongst the Platonic source distributions, the largest cloaked region is attained when L = 12.
Given the regularity and rotational symmetry of the Platonic source distribution, we can work out the amplitudes of all sources using only the knowledge of the explicit formula (11) and the surface integral I mµ nν (x 1 , ∂C 1 ) instead of performing the integration for all the L active sources. To elucidate the approach, it is useful to write the coefficient q ,nm as q ,nm ( k, x , ∂C ) to indic-ate its dependence on the propagating vector k and the two geometric parameters. For notational brevity, we also denote the surface integral as I by suppressing the indices and arguments in (9). Define a new coordinate system x = (x , y , z ) = R( v, Θ)x where R( v, Θ) is a rotation matrix with v the unit vector representing the rotation axis and Θ the rotation angle. As illustrated in Fig. 3, the idea is to rotate the original frame x such that in the rotated frame x , the th active source takes up the bottommost location replacing the source = 1 and the corresponding bounding surface ∂C has the same orientation as ∂C 1 does in the original frame x. Under such transformation, the source amplitude in (11) becomes where k , x and ∂C are the new forms of the propagating vector, position vector and integration surface in the space x respectively and Q νµ = 4πi ν Y µ ν ( k ). Note that by construction, x = x 1 and ∂C = ∂C 1 . As the forms of D νn and I suggest in (8) and (9), they are independent of k. The rotational symmetry of the Platonic solids and our choice of R thus ensure that they remain rotationally invariant across all . In particular, we have D νn = D νn and I = I 1 . The source amplitude q ,nm ( k, x , ∂C ), is thus transformed to q 1,nm ( k , x 1 , ∂C 1 ) in the rotated frame x . While the former has to be calculated by integrating over the faces ∂C for all , the latter can be evaluated by simply integrating over ∂C 1 parameterized by (12) - (15) and replacing k by k . If we express the source coefficients as a vector q ,n = {q ,nm } n m=−n , then it can be shown that the system of linear equations D n (γ, β, α)q ,n ( k, x , ∂C ) = q 1,n ( k , x 1 , ∂C 1 ) (17) holds where D n (γ, β, α) (defined in (122) -(124) of Supplementary Material) is the Wigner D-matrix [35] with dimensions (2n + 1) × (2n + 1) and γ, β, α are the Euler angles [36] of the matrix R( v, Θ) such that with e y , e z the unit vectors along the y, z axes. By solving (17) for n = 0, 1, · · · , N with N a positive integer , we can retrieve the full set of amplitudes q ,nm ( k, x , ∂C ) for the th active source up to the N th order multipole. The configuration of the cloaking device can be completed by repeating this procedure with the corresponding k and Euler angles for different sources. The benefit of employing active sources distributed on the vertices of the Platonic solids is therefore that once one set of source coefficients is determined for a given location, those for others follow by simple post-processing (a) (b) The reduction of the scattered field from an obstacle by the active exterior cloaking method discussed is illustrated in Fig. 4 with L = 20 active sources distributed at the vertices of an imaginary regular dodecahedron with source distance kx 0 = 5π. Here we consider the problem in the context of acoustics but the principles are similar for other scalar waves. The incident plane wave is propagating in the positive x direction with θ i = π/2, ϕ i = 0. The scattering object inside the cloaked region C is a sound-soft sphere with radius kA = 3(1 − a)kx 0 ≈ 0.8579π where a = (2 √ 6/3)/[2 sin(π/3)]. Note that here a is taken as the minimum radius permissible of the imaginary spheres bounding C for L = 4, which means that the scattering sphere has a radius three times that of the inscribed sphere of C when L = 4. In both subplots the real part of the total wave field u on the cross section z = 0 is shown. In Fig. 4(a), the cloaking devices are inactive and a prominent scattering pattern including distorted wavefronts and a shadow region behind the sphere is observed. In Fig. 4(b), the cloaking devices are activated and the multipole order of each active source is taken as N = 10. The series expansions for the source coefficients in (11) are truncated such that the active field produced by each source is within 1% relative error. The sphere now resides in a region with literally zero wave amplitude. The straighter wavefronts and the absence of the shadow region in the figure indicate that the incident wave is scattered only slightly by the sphere inside the quiet zone and the sources radiate little into the exterior of the silent region, demonstrating the effectiveness of the cloak. Note that the wave field diverges within the small neighbourhoods centred at the active point sources. In practice these large fields are confined within the finite-sized sources. For visualization purpose, we have cropped the excessively large parts to a sufficiently small value.
While it is seen that our cloaking method can effectively suppress scattering on a local cross section, the cloaking effect can be assessed globally by examining the total power W radiated by the entire system into the far field. We may look at the quantity σ defined by which, in the context of acoustic wave, is the ratio between the sound power detected in the far field after (when u d = 0) and before (when u d = 0) the cloak is activated. Detailed expressions for W (u d = 0) and W (u d = 0) are derived in Part 6 of Supplementary Material. In Fig. 5, we plot the sound power level (SWL) defined by SWL = 10 log σ against the nondimensionalized radius kA of the sound-soft sphere inside the silent zone for the five source distributions depicted in Fig. 2. The modal order N , the sphere radius A and the direction of incident wave remain the same as those in Fig. 4. The SWL is simulated for wavenumbers in the range π/2 ≤ kx 0 ≤ 6π at an interval of π/2 (or 0.0858π ≤ kA ≤ 1.0294π at an interval of around 0.0858π). Note that at each observation point, the truncation parameter of the corresponding source amplitude is chosen to achieve a relative error less than 1% in the active field generated by each source. For our cloaking approach to work we require σ < 1 and thus 10 log σ < 0. It is observed that in all cases of L the SWL remains negative for the range of kA studied here with a more significant reduction attained for lower frequencies. As  the number of sources L increases the cloaking effect improves in general with a maximum reduction of about 70 decibels achieved when L = 20 at kA ≈ 0.1716π. Further study about how changes in different parameters affect the performance of the cloaking system can be found also in Part 6 of Supplementary Material. In summary, we have formulated an efficient threedimensional active exterior cloaking strategy for the scalar Helmholtz equation, which employs multipolar sources distributed at the vertices of the Platonic solids to create a domain with zero total wave amplitude. This approach minimizes the incident wave impinging upon the object inside the region to suppress the scattered wave whilst simultaneously ensuring that the active source radiation is minimized. The Platonic distribution of the active sources means that we only need to determine the source amplitudes as an integral expression at one source location in terms of the incident field and those at other locations will follow from simple post-processing, by exploiting the symmetry and regularity of the Platonic solids.
This work was supported by a University of For a field u satisfying the homogeneous Helmholtz equation (∇ 2 + k 2 )u = 0 in a closed arbitrary domain C bounded by the surface ∂ C without any radiating source inside, the Kirchhoff-Helmholtz integral equation [37] states that where G(x, y) is the Green's function in threedimensional free space [38] defined by and ∇ y f (y) · n is the normal derivative of f at y ∈ ∂ C with the outward-pointing unit normal n. As illustrated in Fig. 6, we designate a region C ⊆ C such that for Applying this condition and substituting (21) into the first case of (20), we can write the active field u d in terms of the incident field u i as The expression (22) implies from a physical perspective that u d is equivalent to the superposition of a continuous source distribution on the bounding face ∂ C with amplitudes depending only on information related to u i on ∂ C. The next step is to replace this distribution with discrete point sources located at x / ∈ C where = 1, 2, · · · , L. We first partition ∂ C into the segments ∂ C such that Given two fixed spatial points r 0 , r 0 ∈ R 3 , a general position vector r − r 0 can be decomposed as r − r 0 = r − r 0 + r 0 − r 0 . Take r > and r < such that The addition theorem for spherical wavefunctions [33] states that where z n (.) is the spherical Bessel or Hankel function of any kind and S mµ nν (r < ), S mµ nν (r > ) are constants dependent on the position vectors r < , r > such that In (27), q 0 and G(n, m, ν, −µ, q) are defined as n (k|r − r 0 |), r = x, r 0 = y and r 0 = x as indicated in Fig. 6. Since |x − x | > |y − x | for y ∈ ∂ C , with the first case of (26) we can show that , Now (31) has the same basis function V m n (x − x ) as the ansatz in (2). Similarly, with the second case of (26) it can be proved that if we take z n (k|r − r 0 |) = j n (k|r − r 0 |), r = y, r 0 = 0 and r 0 = x as |y − x | < |x | . Note that the expres-sions (31) and (32) hold for all x ∈ C where |x − x | > |y − x |. The arbitrary nature of the domain C means that we can deform ∂ C onto the bounding face ∂C of the cloaked region C and take |y − x | = a x 0 for y ∈ ∂C with a x 0 a positive constant (possibly different for each ). As a result, ∂C is a continuous set of points located on a sphere with centre x and radius a x 0 . If we repeat this procedure for = 1, 2, · · · , L, then the cloaked region C will be the closed domain interior to a region D consisting of a union of L imaginary spheres S such that Substituting (31) and (32) into (23) and comparing with the ansatz for the active field u d (2), we can show that the source coefficient q ,nm can be written in terms of the expansion coefficient of the incident field Q nm in the form where q ,nm,ts Using the property that the outward unit normal n points towards x given the geometry of ∂C , we can replace the parameterizing vector y by y − x and evaluate the derivative terms as The expressions (6) -(9) can be obtained by substituting (36) into (34) - (35).
Note that the form of the source amplitude (6) -(9) can be simplified further when the incident wave u i is a plane wave with u i (y) = e ik·y , which is evident if we write with The simplified expression (11) will then follow upon the substitution of (38) into (6) - (9). As the total wave amplitude is driven to zero within C with q ,nm taking the forms (34) - (35), it remains to be shown that the active sources produce negligibly small radiation to the far field. We first notice by the ansatz (2) that the integral form of u d (x) in (22) should be defined for all x ∈ R 3 (although it is equal to −u i (x) only for x ∈ C). The derivation of the source coefficient in the form (34) involves expanding the Green's function in the integrand of (22) using the first case of the addition formula (26), which holds only for x ∈ R 3 \ D with D defined in (33). From (20), the Kirchhoff-Helmholtz integral equation is identically zero for x ∈ R 3 \ C. For multipolar sources of sufficiently high order, one can thus deduce that u d (x), in the form (2) and (35), vanishes identically for x ∈ R 3 \(C ∪ D), which is the region exterior to the imaginary spheres centred at the active sources.  (20) holds and the cloaked region C ⊆ C (the grey area). The continuous source distribution on the bounding surface ∂ C is replaced by the point active sources (the small dots) at x where = 1, 2, · · · , L. The face ∂ C is divided into different sections and particularly for y ∈ ∂ C (the green curve), we decompose the position vector x−y (the black arrow) into x − x (the red dashed arrow) and y − x (the blue dashed arrow) and apply the addition formula (26) to expand the wavefunctions as one centred at the source at x (the green dot).

Geometry of the source distribution and the cloaked region C
Given the characteristic indices (p, q) of a Platonic solid, the length of its side sx 0 and the position vector of the bottommost source x 1 = (0, 0, −x 0 ), we aim to determine the locations of the remaining sources. Since the source distribution has a q-fold rotational symmetry about the z axis by the definition of q, the sources must be distributed on a number of horizontal planes sliced through the polyhedron and the number of sources residing on each plane must be an integral multiple of q (except when the source is located on the z axis like the = 1 one). We can assign the indices 1 ≤ ≤ L to the active sources based on two rules: 1) the plane it is located on starting from the bottom; 2) in either the clockwise or counterclockwise direction for sources on the same plane. Any sources can thus be located by specifying a set of two integers (i, j). The first index i = ( − 1)/q + 1 represents the order of the plane and ranges from i = 1 to i = N L = (L − 1)/q + 1, where . is the ceiling function. The second index j = mod( − 2, q) denotes the order of the source on a particular plane and has the range 0 ≤ j ≤ q − 1, where mod(., .) is the modulo function. By simple geometric arguments we can show for all cases of L that where Note that j is not needed when = 1 or L provided that L = 4 as these two sources are positioned on the z axis. For convenience, we set j = 0 With x in the form (39) - (44), the numbering of sources on a plane goes counterclockwise starting from Quadrant I for i ≤ max ((N L + 1)/2, 2) and clockwise starting from Quadrant II for i > max ((N L + 1)/2, 2). Given the two indices i, j, we can retrieve as In Fig. 7 we illustrate the labelling of L = 20 active sources located at the vertices of an imaginary regular dodecahedron with p = 5, q = 3 based on the principles discussed above. Fig. 7(a) shows that when the position vector of the source = 1 is taken as As we have demonstrated in Part 1, the cloaked region C is the domain enclosed by the imaginary spheres S , each centred at x with radius a x 0 . We will keep a constant at a for all so that all spherical faces ∂C have an identical geometric shape. For our formulation of the source amplitudes to hold, a must be within the range Note that we have a < 1 as otherwise all the space interior to the source distribution would be occupied by the spheres and C would not exist. As for the lower bound of (46), we need to ensure that adjacent spheres intersect with each other so that C remains completely bounded by the spheres. The possible intersecting point furthest away from the origin is the centroid of each flat polygonal face of the Platonic solid, which suggests that a must be at least as large as s/[2 sin(π/p)] since this is the distance between the centroid and each vertex of a p-sided regular polygon with side s.
When we evaluate the coefficient q 1,nm of the source = 1, we need to first perform the surface integral I mµ nν over the spherical face ∂C 1 . To obtain a parametric form of ∂C 1 with respect to the point x 1 , we see from Fig.  8 that ∂C 1 is bounded by a set of circular arcs, each being part of the intersecting circle (the blue solid and red dashed curve) between S 1 and S for = 2, · · · , q + 1, where the spheres S centre at the sources on the plane i = 2. The equations for S 1 and S are thus |(x, y, z) − x 0 (r 2 cos ϕ , r 2 sin ϕ , z 2 )| = ax 0 , respectively, with z 2 , r 2 given by (41), (42) and ϕ = ϕ 2,j = 2mod( − 2, q)π/q = 2( − 2)π/q by (43), (44) and the fact that − 2 ≤ q − 1 < q. Eliminating the squared terms from (47) and (48), we have (r 2 cos ϕ )x + (r 2 sin ϕ )y + (1 + z 2 )z = 0, which gives the plane that the intersecting circle between S 1 and S lies on. If we introduce the parameterization (x, y, z) = x 0 (a sin θ cos ϕ, a sin θ sin ϕ, −1 + a cos θ) and apply the identity sin 2 θ + cos 2 θ = 1, then (49) becomes where h = (1 − z 2 )/(1 + z 2 ). Solving for cos θ in terms of ϕ yields Note that (51) represents the entire circumference of the intersecting circle between S 1 and S . The plus and minus form equal each other when the square root term vanishes, which gives When a takes its minimum value as stated in (46), we can show that indeed ϕ = ϕ ± π/q in (52) by substituting the values of h and lower bound of a corresponding to the five different source distributions. Considering the regularity of the Platonic solids, we can deduce that the values of ϕ and θ defined by (52) and (53) give the centroids of the two adjacent polygonal faces that meet at both source 1 and . The intersecting circle between S 1 and S thus passes through these two centroids and they are where the two cases of (51) coincide with one another. As Fig. 8 shows for the minimum case of a, only the circular arc interior to the Platonic solid (the blue solid curve) forms part of the boundary of ∂C 1 and it is defined by smaller values of θ (as viewed from x 1 ) than the exterior one (the red dashed curve) is. Since cos θ is monotonically decreasing for θ ∈ [0, π], only the plus form of (51) represents the boundary of ∂C 1 . For values of a greater than the lower bound of (46), the intersecting circle no longer passes through the centroids. It can be proved that the two points where the plus and minus part of the circumference meet each other occur at values of ϕ such that |ϕ − ϕ | > π/q since cos(ϕ − ϕ ) can be shown to be monotonically decreasing with a. Nevertheless, taking into account the q-fold rotational symmetry of ∂C 1 about the z axis, we can show that each of the q circular arcs bounding the face ∂C 1 spans only the range |ϕ − ϕ | ≤ π/q and the plus form of (51) suffices to fully describe the boundary of ∂C 1 .
The volume V = |C| of the cloaked region C can be found using the divergence theorem [39], which states that for any vector field f , with n the outward unit normal emanated from the area element dS(x). Taking f (x) = x = (x, y, z), we have as each ∂C is identical to ∂C 1 in terms of its geometric shape. If we shift the reference point from the origin to x 1 by replacing x − x 1 with x, then (55) becomes where a 0 = ax 0 and g (ϕ) is the plus form of (51). Note that (56) is applicable to the whole range of a in (46). Figure 8. Illustration of the geometry of the face ∂C1, which is part of the surface of the sphere S1 centred at the source = 1 (the green dot), when a takes the minimum value in (46). Its boundary is formed by the intersection of S1 with the spheres S , each centred at the source = 2, · · · , q + 1 (the red dots). The blue solid part of the intersecting circle is represented by the plus form of (51) while the red dashed part is by the minus one. The two circular arcs meet with each other at the centroids of two adjacent polygonal faces.

Truncation error of the source coefficients q ,nm
In (11), we have expressed the source amplitude q ,nm as an infinite series in terms of the incident wave coefficient Q νµ . To numerically implement the cloaking approach, only a finite number of terms from the series can be retained which leads to truncation error. We will show that this error can be bounded by choosing an appropriate value of truncation order.
Recall from (23) in Part 1 the integral form of the active field u d with If we substitute the expansion formula of V 0 0 (x − y) in (31) into (57) and compare the expression with the ansatz of u d in (2), then we can write the source amplitude q ,nm as For a plane wave incidence, we have by (37) that where Q νµ = 4πi ν Y µ ν ( k). From (58) and (59) we observe that the truncation error of q ,nm originates from the truncation of the expansion of u i . If we discard the terms with ν > V in (59) and denote the truncated source amplitude and incident plane wave as q ,nm and u i respectively, then the truncation error is It has already been shown in [41] and [42] that and when V = eka 0 /2 +∆ (as |y−x | = a 0 for y ∈ ∂C ) with ∆ a nonnegative integer, (61) is bounded as The result (62) states that the truncation error in u i decays exponentially for truncation order greater than the critical value eka 0 /2 . As for the truncation error of ∇u i , we can rewrite it in a form similar to (61): where j ν (ka 0 ) is expressed as a linear combination of j ν−1 (ka 0 ) and j ν+1 (ka 0 ) in the third line [43]. Taking V = eka 0 /2 + ∆ + 1 and referring to (62), we have Note that at this value of V , (62) becomes By [44], the term U m n (y − x ) can be bounded as with Due to the rotational symmetry of the source configuration, the surface area of each piece of spherical face dC is same as that of dC 1 . The equation (53) defines the polar angle of the two end points of the circular arc cos θ = g (ϕ) with respect to x 1 for = 2, · · · , q + 1. For −π/q ≤ ϕ ≤ π/q, we thus have g (ϕ) ≥ a. Using these arguments, we can show that While it has been shown that there exists a critical value for the truncation order V such that the truncation error of the source amplitude q ,nm approaches zero for V sufficiently larger than it, the exact minimum value of V required to achieve a certain level of accuracy in computation remains to be determined as the upper bound (69) depends on the arbitrary parameter ∆. We are concerned about how an error in q ,nm will lead to discrepancy in the active field u d | generated by the th source and thus degrade the performance of the cloaking system. Define A measure to quantify this effect is the relative squared truncation error averaged over the surface of a unit sphere Ω centred at x , denoted by ε, where To establish an upper bound on ε, we can use the maximum value of |q ,nm | as a proxy for the benchmark in (71). Note that with |∇ y u i (y) · n| = |ik k · ( y − x )e ik·(y−x ) e ik·x | ≤ k (73) and (66), we have Substituting (69) and (74) into (71) and noting that N is the multipole order for each active source, we can obtain an upper bound for in the form with H n (x) = (2n + 1)(n + 1)|h (1) Given the tolerance level E, the minimum truncation order V required for the source amplitude q ,nm such that < E is thus V = eka 0 /2 + ∆ + 1, where Applying this expression of V onto the settings of Fig. 4 with ka 0 = 5π · s 2 sin(π/5) , N = 10 and E = 1% where s is the length of the side of a regular dodecahedron with x 0 = 1, we truncate the series expansion of the source coefficient in (11) at V = 16.

Computation of the surface integral I mµ nν (x1, ∂C1)
To implement the cloaking strategy discussed, we need to evaluate the surface integral which will be applied to determine the strengths of all active sources given that they are equidistant from the origin and all spherical faces ∂C have an identical geometric shape in a Platonic source distribution. By the definition of the normalized spherical harmonic function and the parametric form of ∂C 1 derived in Part 2 (and stated in (12) -(15)), we can write the integral as Now replace the azimuthal angle ϕ by ϕ + ϕ in the integration with respect to ϕ. This substitution eliminates the dependence of g (ϕ) on the index as can be seen from (14) and thus we may denote g(ϕ) = g (ϕ + ϕ ). Note that g(ϕ) is even about ϕ = 0. The surface integral then becomes Up to this point it appears that simplifying the surface integral into its analytic form is not trivial. Nevertheless, we may still look for ways to speed up the computation of the expressions (81) -(83). Note that the most computationally-expensive part is the double integral J mµ nν and given the ranges of the indices n, m, ν, µ, a direct evaluation of (81) -(83) over these ranges would entail performing the double integral for a total of N n=0 V ν=0 (2n + 1)(2ν + 1) = (N + 1) 2 (V + 1) 2 (84) times, where N is the multiple order of each active source and V is the truncation order of the source amplitude discussed in Part 3. The computation will be costly especially when higher order sources are used in the cloaking system or when the frequency of the incident wave is high which requires a larger truncation order. We are therefore seeking a simpler form of J mµ nν and thus I mµ nν (x 1 , ∂C 1 ) that can facilitate the computation.
By [40], the associated Legendre functions in (83) can be written in an alternative form where In (87), (p) τ is called the Pochhammer symbol [40] and vanishes for τ ≥ 1 − p if p is a nonpositive integer. Referring back to (86), we observe that (1 − m) τ > 0 for all τ since 1 − m must be positive by our choice of m and c m n,τ is hence well-defined for all τ . Note also that (−n) τ and thus c m n,τ vanish identically for τ ≥ n + 1, which means that the series expansion for P m n (w) in (85) is finite and terminates at τ = n. Substituting (85) into (89) and replacing w by 1 − 2w inside the integral, we have where the dependence of g on ϕ is understood and C mµ nν,τ σ = 2c −|m| n,τ c −|µ| ν,σ . In (91), B x (a, b) is the incomplete beta function [40] defined by which is convergent when a is positive [45]. This condition is again satisfied since in (91), τ +σ + |m| + |µ| 2 +1 ≥

Now substituting (91) into (82) gives
With (93), we have expressed J mµ nν (ϕ) as a linear combination of integrals in terms of the incomplete beta function, which can be readily evaluated using the corresponding built-in function of some common numerical software and thus reduces a fair portion of the computational cost.
While the sum of integrals J mµ nν now takes a simpler form which is more convenient to evaluate, the integral still needs to be performed for many times when the multipole order N or the truncation parameter V becomes large. In the following we outline how we can make use of the properties of J mµ nν to significantly reduce the number of integrations required. Inspection of its form in (93) and the definition of C mµ nν,τ σ , a mµ nν in (86), (90) respectively shows that J mµ nν is symmetric with respect to the two groups of indices (n, m) and (ν, µ) such that J mµ nν = J µm νn . Without loss of generality, we may consider only the case n ≤ ν where 0 ≤ n ≤ min(N, V ) and n ≤ ν ≤ max(N, V ). The symmetry of J mµ nν means that A further inspection of the integrand in (93) reveals that there exists different combinations of the indices τ, σ, m, µ that indeed give the same value of the integral. Take i = τ + σ, j = |m| + |µ|, k = |µ − m| and rewrite (93) as Note that i, j, k remain invariant under the swap of the indices in (94). The strategy is to group the summands in (93) by identifying the integrals that have the same value of i = τ + σ and collecting the coefficients C mµ nν,τ σ multiplied to each of them. With this approach, we can avoid repeatedly evaluating integrals that give the same result. To determine the permutations of (τ, σ), we note that n ≤ ν is already assumed and there are only three possible cases for the value of i: h Case 1: i ≤ n Note that i, j ≤ N + V since As for k, we have |µ − m| ≤ |µ| + |m| = j. The lower bound of k depends on whether j is odd or even because if the sum of two integers is odd (even), then their difference must also be odd (even). The possible values taken by k are thus k = 1, 3, 5, · · · , j − 2, j, j is odd 0, 2, 4, · · · , j − 2, j, j is even (102) or more compactly, k = mod(j, 2), mod(j, 2) + 2, · · · , j − 2, j. Therefore under the change of indices (n, ν, m, µ) → (i, j, k), the integral B ijk is now evaluated for a total of . Our method has led to around 85% decrease in the number of integrations done (along with an integrand less costly to evaluate). The reduction is more significant for larger N and V as the leading order of (84) is quartic but only cubic in (103).

Transformation approach with rotation of the coordinate system
Consider a rotation in three-dimensional space R( v, Θ) which maps the coordinate system x = (x, y, z) to x = (x , y , z ) = Rx with v = (v x , v y , v z ) the unit rotation axis vector and Θ the rotation angle (in a counterclockwise sense). We require that in the rotated frame x , 1) the th active source occupies the bottommost location with position vector x = (0, 0, −x 0 ) = x 1 ; 2) the rotated spherical face ∂C is aligned in the same orientation as ∂C 1 is in the original frame x such that ∂C = ∂C 1 . We are seeking the form of v , Θ and thus the rotation matrix R for = 1, 2, · · · , L.
We first consider the active sources located on the layers with 2 ≤ i ≤ max ((N L + 1)/2, 2). The rotation matrix R can be found using the Rodrigues' rotation formula [46]: where I is the 3×3 identity matrix; v T denotes the transpose of v and V is defined by To find the rotation axis v, we recall from Part 2 that with a taking the minimum value in (46) the vertices of the cloaked region C coincide with the centroid of each p-sided polygonal face of the Platonic solid. By the definition of p, both the polyhedron and C have a p-fold rotational symmetry about an axis through the origin and any vertex of C. Since the source = 1 is the common vertex of q congruent polygonal faces, these faces give a total of q centroids which the rotation axes may pass through. We can take to be the position vector of the centroid associated with the rotation R that satisfies the two conditions stipulated above. By the regularity and symmetry of the source distribution, v is perpendicular to the polygonal face on which the centroid is located. It is hence straightforward to show that in a Cartesian coordinate system centred at the origin, v = x 0 1 − a 2 a cos ψ 2,j , a sin ψ 2,j , − 1 − a 2 , where a is the minimum radius in (46) and ψ 2,j = ϕ 2,j + π/q with ϕ 2,j and j defined as in Part 2. The corresponding unit vector v is thus v = a cos ψ 2,j , a sin ψ 2,j , − 1 − a 2 .
When a takes values greater than a, the vertices of C no longer coincide with the centroid of each polygonal face of the Platonic solid. With (51) the polar angle of the q vertices immediately above the source = 1 is given by × cos 2 (π/q)] 1/2 / a 1 + h cos 2 (π/q) (108) in a spherical coordinate system centred at x 1 . In Part 2 it has already been shown in (52) that Substituting (109) into (108), we have Note that at a = a, we recover the relation cos θ = a in (53) which defines the centroid of each polygonal face that meets at the source = 1. If we again let the rotation axis v be the position vector of the vertex of C that is associated with the rotation of the source for = 2, · · · , q + 1, then with some algebra and the expression (110) we can prove that in a Cartesian coordinate system centred at the origin, Comparing (111) with (106), we can deduce that the vertex of C, the centroid of the polygonal face and the origin are collinear with one another, which implies that the unit vector v in the form (107) indeed applies to the whole range of a.
Since each polygonal face consists of p sides and both the first and th source are located at the vertices of this polygon, the rotation angle Θ must be an integral multiple of 2π/p. Given the form of v in (107), a clockwise rotation of the original frame x by an angle of 2(i−1)π/p about the axis v will 'apparently' move the source located on the ith layer to the bottommost location. (A concrete example for the case L = 20 will be illustrated later in Fig. 9.) Provided the convention that positive Θ denotes a counterclockwise rotation of the frame x, we take Θ = −2(i − 1)π/p. Substituting (107) and this value of Θ into (104), we have the rotation matrix R = {R ,st } 3 s,t=1 given by where stu is the Levi-Civita symbol [39] defined by stu = 1, if (s, t, u) = (1, 2, 3), (2, 3, 1) or (3, 1, 2) −1, if (s, t, u) = (2, 1, 3), (3, 2, 1) or (1,3,2) .
As for the active sources located on the layers i > max ((N L + 1)/2, 2), we note that with the orientation of each Platonic solid specified by (39) - (44) in Part 2, they possess a two-fold rotational symmetry about the y axis, except for the case L = 4. In mathematical terms, this rotation can be expressed by the matrix R( e y , t) where with e y the unit vector pointing in the positive y direction and t the angle of rotation. By setting t = π we can show that for i > max ((N L + 1)/2, 2), The relation (115) means that we can always match a source located on the plane i > max ((N L + 1)/2, 2) with one on a lower layer i ≤ max ((N L + 1)/2, 2) and their position vectors are related to each other by (115). We can make use of the results derived in the previous case and express the rotation matrix R for i > max ((N L + 1)/2, 2) as and R ( v , Θ ) is in the same form as (112) with v defined by (107) and −2(i − 1)π/p replaced by −2(N L − i)π/p.
To illustrate the method outlined above, we look at the case L = 20 in Fig. 9, which shows the different rotation axes associated with the rotation that transforms the position vector x to x 1 = (0, 0, −x 0 ) (the yellow dot at the bottom in both subfigures). Fig. 9(a) shows the case for the layers 2 ≤ i ≤ max ((N L + 1)/2, 2). For L = 20, it involves those sources with i = 2, 3, 4 (or equivalently 2 ≤ ≤ 10). Since (p, q) = (5, 3) by Table I, they form a total of three regular pentagons that meet at the first source. These sources can be classified into three groups with components of a group having the same value of j and represented by dots of the same colour. The three groups are j = 0 coloured in red, j = 1 in blue and j = 2 in green. By the numbering of the sources on the layers i = 2, 3, 4 in Fig. 7(b) -(c) in Part 2, the three groups consist of sources with = {2, 5, 8}, {3, 6, 9} and {4, 7, 10} respectively. Sources in the same group share the common rotation axis v , where is the index of any source in the group since the form of v in (107) depends only on j and is the same for all components of a group. (Note that in Fig. 9(a), the minimum value of in the group is taken.) Each of these rotation axes points in the direction (107) to-wards the centroid of the pentagon that the group is located on. The angle of rotation for each active source in the group is given by Θ = −2(i − 1)π/5 with i the order of the layer the source resides on. For example, for the group = {2, 5, 8}, we have i = {2, 3, 4}. A clockwise rotation by 2π/5, 4π/5, 6π/5 of the frame x about the axis v 2 will 'bring' these sources onto the lowest position respectively. In Fig. 9(b), the case with i > max ((N L + 1)/2, 2) is illustrated. The two-fold rotational symmetry of the source distribution about the unit vector in the y direction e y (the brown arrow) means that we can match any source with i = 5, 6, 7, 8 to another source with i = 4, 3, 2, 1 which has the same value of j by the relation (115). In the figure we focus on the case with j = 0 and depict the sources which form a pair according to (115) in the same colour. Again by the labelling of sources in Fig. 7  The rotation matrix R is related to how the form of the active field u d transforms between the original space x and the rotated one x . (The subscript of R is suppressed for convenience as its dependence on is understood.) With the incident wave propagating in the direction k, the active field u d | radiated by the th source can be posed in the original frame x as where q ,nm ( k, x , ∂C ) is the source coefficient that we are seeking. Alternatively, if we switch to the rotated frame x , then u d | is in the form where q ,nm is the source amplitude in the rotated space with k = R k. The equality between q ,nm ( k , x , ∂C ) and q 1,nm ( k , x 1 , ∂C 1 ) holds since = q 1,nm ( k , x 1 , ∂C 1 ).
Note that q 1,nm ( k , x 1 , ∂C 1 ) can be evaluated using the explicit formula (6) -(10) with the parametric form of ∂C 1 in (12) - (15) and k in place of k. (In (120), the form of the source coefficient with a general incident wave is used as the case of plane wave incidence should follow immediately.) The forms of representation (118) and (119) are equivalent with each other as they describe the same quantity u d | but only in different basis functions. The correspondence between the two expressions can be established by considering the transformation of the spherical harmonic function Y m n ( x − x ) between the two frames [47]: where D n m m (γ, β, α) is the (m , m)th entry of the Wigner D-matrix [35] with degree n defined by In (123) -(124), the summation over j is such that the factorial terms are nonnegative. The arguments (γ, β, α) are the Euler angles [36] of the rotation matrix R( v, Θ) when it is decomposed into three separate rotations about the z, y, z axes of the original space x such that R( v, Θ) = R( e z , γ)R( e y , β)R( e z , α), where R( e y , β) is already defined in (114) and Comparing elementwise between the left and right side of (125), we have with the entries of R described by either (112) or (116), depending on the value of . If we substitute (121) into (118) and compare with (119), then for m = −n, · · · , n, we have n m=−n D n m m (γ, β, α)q ,nm ( k, x , ∂C ) = q 1,nm ( k , x 1 , ∂C 1 ).
Note that (130) can be written more compactly in the matrix form (17) if we set D n = {D n m m } n m ,m=−n and q ,n = {q ,nm } n m=−n .

Wave scattering and parameter study
Consider a smooth object subject to an incident scalar wave u i and denote the ensuing scattered wave as u s . The two wave fields can be expanded in terms of the basis functions as with Q nm and a nm the incident and scattering coefficient respectively. For simplicity, we only look at the case of acoustic wave but our arguments can be easily generalized to other cases like electromagnetic or elastic wave. It has been derived by [48] in the acoustic context that where T nm,ji is a transition matrix depending on the geometry of the scattering object and the boundary condition imposed on its surface. For a perfect sphere of radius A, T nm,ji takes the form for sound-hard case (134) j n (kA) h n (kA) for sound-soft case (135) with δ hk the Kronecker delta function.
By adopting procedures similar to what is done in [49], it can be shown that W , the acoustic power scattered to the far field by the object in three dimensions averaged over a wave period, is given by with ω the angular frequency and ρ the density of the medium. We can apply (136) as a measure to assess the effectiveness of our cloaking approach discussed before by defining the quantity σ in the form where u d and a nm , defined in (2) and (3), are the active field produced by the sources and the scattering coefficient of the sphere inside the cloaked region. The term f nm is the far-field amplitude of u d such that In (137), for the case without any control, the coefficient a nm depends solely on the scattered field u s and is thus equivalent to a nm (u d = 0). However, when the cloaking devices are activated (u d = 0), u s is influenced also by u d and a nm becomes a nm (u d = 0). Note that the contribution from f nm should also be included in W (u d = 0) since we require not only a reduction in the scattered field but also minimal radiation from the active sources to the far field and (138) has the same basis function as (132) does. The quantity σ is therefore the ratio of the total power radiated by the whole system before and after control.
To ensure that the cloaking method is effective, we need σ < 1.
In the following we will find the expressions for the three coefficients in (137). Substituting Q nm = Q nm and (135) into (133) gives a nm (u d = 0) in the form when the scattering object is a sound-soft sphere. As for a nm (u d = 0), we recall from Part 1 that we apply the addition theorem for spherical wavefunctions in (26) - (30) to determine the explicit form of the source amplitudes. Similarly, to obtain the form of a nm (u d = 0), we need to first invoke the addition formula to write the active field u d in terms of the incoming wavefunction U m n (x). Since |x| < |x | on the surface of the scattering sphere, by putting z n (k|r − r 0 |) = h (1) n (k|r − r 0 |), r = x, r 0 = x and r 0 = 0 in the second case of (26) for = 1, 2, · · · , L, we can rewrite the active field from every source in a form similar to (131). However, given the complexity of the expressions in (26) - (30) and the fact that the evaluation has to be repeated for L times, each with a different position vector x , the computation of the coefficients related to the active field will be expensive. A way to circumvent this problem is to again make use of the regularity and rotational symmetry of the Platonic source distribution like what we have done in Part 5 when we determine the amplitudes of all active sources using only knowledge of that of the bottommost source and the geometry of Platonic solids. To dis-tinguish between the contributions from different active sources, we now denote the rotated space for each source as x ( ) = (x ( ) , y ( ) , z ( ) ) = R x (which is different from the notation in Part 5). In the rotated frame x ( ) , we have the position vector x ( ) = x 1 = (0, 0, −x 0 ), the spherical face ∂C ( ) = ∂C 1 and the propagating vector of the incident wave k ( ) = R k. Note that x (1) = x and k (1) = k. In the new notation, the equation of transformation between the original and rotated spaces for the spherical harmonic function (121) becomes for = 1, 2, · · · , L, where γ , β , α are the Euler angles of the rotation matrix R . The system of linear equations for the source amplitudes in (130) is now Substituting (140) and (130) into the expression of u d | in (118) for = 1, 2, · · · , L and summing the active fields over all sources, we have the total active field u d in the form where the dependence of q 1,nm on the geometric parameters x 1 and ∂C 1 is already understood. An explicit formula for q 1,nm is already available in (6) -(10) along with the parametric form for the surface of integration ∂C 1 in (12) - (15). Note that the expression of S mµ nν (r 0 − r 0 ) in (27) can be simplified for which is the pointing in the positive z direction. The spherical harmonic function in (27) now becomes where q = q 0 + 2q and the exponent µ − m must vanish since the azimuthal angle ϕ is undefined for θ = 0. By [40], the associated Legendre polynomial P p q (x) can be written as With p = 0 and x = 1, we apply the general Leibniz rule [50] on the derivative such that By (144) and (145) where q 0 ≡ |n − ν| since in (28), only the first case is possible with µ = m.
where |x ( ) | < |x 1 | = x 0 on the surface of the scattering sphere. Note that the term S mm nν (−x 1 ) is now independent of and needs to be evaluated only for −x 1 using the simplified form in (146). The Kronecker delta function in (146) also filters out the azimuthal modes with µ = m. The next step is to revert (147) to the originial space x so that it has the same basis function as the incident field u i does in (131). Referring back to how the spherical harmonic function is transformed from x to x ( ) in (140) and the order of the Euler angles in (125) in Part 5, we can simply reverse the rotation and write The method illustrated above where we make use of the property that the spherical harmonic function vanishes for nonzero azimuthal order on the z axis can also be applied to simplify the expression of the source coefficient q ,nm in (6) -(10) under a general wave incidence. As we evaluate the surface integral I mµ nν (x 1 , ∂C 1 ) and adopt the rotational approach in Part 5 to determine all source amplitudes, the term S sµ tν (x ) needs to be computed only for x 1 = −x 0 e z . Similar to what is done in (143) -(146), we have Y µ−s q0+2q (− e z ) = (−1) q0+2q [2(q 0 + 2q) + 1]/(4π)δ µs and thus S sµ tν (x 1 ) = 2 √ πδ sµ (−1) s (t+ν−q0)/2 q=0 i ν−t−(q0+2q) 2(q 0 + 2q) + 1 × j q0+2q (kx 0 )G(t, s, ν, −µ, q 0 + 2q), where q 0 ≡ |t − ν|. The term q 1,nm,ts from (7) now becomes q 1,nm,ts = ∞ ν=0 S ss tν (x 1 )D νn I ms nν (x 1 , ∂C 1 ).
In Fig. 10, we plot the sound power level SWL = 10 log σ as a function of kA, which is the dimensionless radius of the scattering sphere inside the silent region. The quantity σ is in the form (137) with expressions for a nm (u d = 0), a nm (u d = 0) and f nm given by (139), (150) and (153) respectively. The active sources consist of monopoles and dipoles only with multipole order N = 1. The settings of the plot are otherwise identical to those of Fig. 5. Comparing Fig. 10 with Fig. 5, we observe that the range of kA where there is a substantial power reduction is mostly confined to kA < π/4 in all cases of L. The reduction attained is also less than that in Fig.  5. For kA > π/4 the power is even amplified slightly when L = 12, 20. Nonetheless, the plot shows that in more practical scenarios where only monopole and dipole sources are realizable, our cloaking approach is still capable of considerably reducing the power radiated to far field for lower frequencies. Another point to note from both Fig. 5 and 10 is that while the cloaking effect gets better in general for more larger value of L, increasing the number of sources does not necessarily increase the reduction at some frequencies. This phenomenon is especially obvious when we increase from L = 6 to 8 or from 12 to 20. A possible reason is that the volume of the cloaked region indeed becomes smaller in these two cases as can be seen in Table I. In Fig. 11, we show the multipole order N required for each active source in order to achieve a minimum reduction of 10 decibels for each observation point in the wavenumber ranged π/2 ≤ kx 0 ≤ 6π (or 0.0858π ≤ kA ≤ 1.0294π) for the five different source distributions. It ap-pears that the more sources employed, the less multipoles required to reach the same level of reduction in general. The result also seems to agree with what is previously suggested by Fig. 5 and 10 as the order of multipoles needed to attain a certain degree of reduction becomes higher when the frequency increases. Figure 10. The sound power level SWL = 10 log σ versus the nondimensionalized radius kA of the sound-soft sphere inside the cloaked region for five different number of sources L and π/2 ≤ kx0 ≤ 6π at an interval of π/2 (or 0.0858π ≤ kA ≤ 1.0294π at an interval of around 0.0858π). Here a multipole order of N = 1 (i.e., monopole and dipole only) is used. The other parameters are the same as those in Fig. 5. Figure 11. The order of multipole required for the five source configurations to attain a minimum reduction of 10 decibels over the range π/2 ≤ kx0 ≤ 6π at an interval of π/2 (or 0.0858π ≤ kA ≤ 1.0294π at an interval of around 0.0858π). Note that when L = 4, the SWL fails to reach −10 decibels at kA ≈ 0.8579π, 0.9437π, 1.0294π in the range N ≤ 30. Here the values of N which can achieve a reduction closest to 10 decibels are chosen for these three cases instead.