DirichletAnnulusApproxFun.jl

DirichletAnnulusApproxFun.jl

using ApproxFun
using LinearAlgebra
using Interact
using Plots

a = 1.; b=5.
Ω = a..b
# Ω = Chebyshev(a..b)
r = Fun(identity, Ω)
Δ_rad = 𝒟^2 + 1.0/r * 𝒟
L = -Δ_rad


rhos = Dict(
    "ρ(r) = 1" =>Fun(r -> 1., Ω),
    "ρ(r) = r" => r,
    "ρ(r) = 1/r" => 1/r,
    "ρ(r) = b-r" => b -r,
)
bdry = Dirichlet()    # zero potential outside annulus

ks = collect(keys(rhos))
@manipulate for key in ks
    ρ = rhos[key]
    Φ = [bdry; L] \ [[0.,0.], ρ]
    E = -𝒟*Φ
    r0 = first(roots(E))
    C = cumsum(r*ρ)
    plt = plot(xlabel="radius", title="r0=$r0")
    plot!(Φ, label="Potential")
    # plot!(E, label="Electric field")
    # plot!(C, label="Charge cdf")
    plot!(ρ, label="Charge density")
    vline!([r0], label="left current = $(C(r0)/C(b))")
    plt
end