Numerical instability related to hankel order

[arturgower]

This is accentuated when particles are small and there is at least two particles far apart. This primarily affects the internal waves, which makes verifying the boundary conditions difficult.

This does not happen in the master version! The way that internal waves are calculated are identical in master and version-0.2, so I think it is related to solving the scattering_matrix.

import StaticArrays: SVector
using MultipleScattering

ω = 0.1

medium = Acoustic(1.,1.,2)
sound = Acoustic(medium.ρ, 4. + 0.0im,2)
p1 = Particle(sound,Circle([-5.0,0.0], .2))
p2 = Particle(sound,Circle([5.0,0.0], .2))

source = TwoDimAcousticPlanarSource(medium, SVector(0.0,0.0), SVector(1.0,0.0), 1.)

particles = [p1, p2]
sim = FrequencySimulation(medium, particles, source)

widths = -5.4:0.01:-4.4
x_vec = [ SVector(x,.0) for x in widths]
result = run(sim, ω, x_vec; basis_order = 7)

using Plots; pyplot()
plot(widths, abs.(field(result)[:]), ylim = (0.9,2.))

But by lowering the Hankel order by one

result = run(sim, ω, x_vec; basis_order = 6)
plot(widths, abs.(field(result)[:]))

Note that as the density of the particle is the same as the background, the field should be smooth everywhere. This instability goes away if we increase the frequency or bring the particles closer together. However, it does not go away by just adding more particle between them. For example:

particles = [Particle(sound,Circle([x,0.0], .2)) for x = -5.:0.6:5 ]
sim = FrequencySimulation(medium, particles, source)

widths = -5.4:0.01:0.4
x_vec = [ SVector(x,.0) for x in widths]
result = run(sim, ω, x_vec; basis_order = 7)
plot(widths, abs.(field(result)[:]))

[arturgower]

This shows the formula we should use now

[jondea]

Great find, do you want me to do anything?

[arturgower]

Thanks man! But it's solved now and all boundary conditions seem to be working fine now. My fault for not simplifying the equation for the scattering coefficients...

[jondea]

It's interesting to note that now the T-matrices are only required in the calculation of the scattering matrix.