Eigenvalues problem in ApproxFun


#1

Is ApproxFun able to compute the eigenstates of a 1dimensional linear differential operator in the way Chebfun does? I assume not yet.
I am interested in the Schroedinger equation.

michele


#2

There’s some very basic support:

using ApproxFun
S = Chebyshev()
B = dirichlet()
D = Derivative(S)
x = Fun()
λ,V=eigs([B;D^2 + x^2],200)

Here, n=200 is a discretization size, and the eigenvalues/vectors are pruned to only return those which have converged.

Here all but the last operator is interpreted as a boundary condition. This syntax will likely be changed, and I’m hoping in the future to replace eigs built out of discrete matrices with a simultaneous inverse iteration (built on top of \).


#3

Thank you so much! Here are the first three states of the harmonic oscillator

computed with the following code:

S = Chebyshev(Interval(-5,5))
B = dirichlet()
D = Derivative(S)
x = Fun(identity, Interval(-5, 5))
λ, V = eigs([B;-D^2/2 + x^2/2], 200)
p = sortperm(λ)
plot(V[p[1]])
plot!(V[p[2]])
plot!(V[p[3]])

and their energies (they should be 1/2, 3/2 and 5/2)

(0.500000000076718,1.5000000036715844,2.500000084018815)

(The eigenstates are not normalized, but you are surely aware of this.)


#4

Sorry to revive this old conversation but I’m trying to do somrthing similar to @mzaffalon but I guess something has changed with ApproxFun so the code used above doesn’t seem to work anymore. Trivially, dirichlet() is now Dirichlet, but even with that fixed, I am getting an error message I don’t quite understand when trying to recreate these results.

using ApproxFun
S = Chebyshev(Interval(-5,5))
B = Dirichlet()
D = Derivative(S)
x = Fun(identity, Interval(-5, 5))
λ, V = eigs([B;-D^2/2 + x^2/2], 200)

produces the error

Implement Conversion from ApproxFun.Chebyshev{ApproxFun.Segment{Float64},Float64} to ApproxFun.ArraySpace{ApproxFun.Space,1,ApproxFun.Point{Float64},Any}
defaultConversion(::ApproxFun.Chebyshev{ApproxFun.Segment{Float64},Float64}, ::ApproxFun.ArraySpace{ApproxFun.Space,1,ApproxFun.Point{Float64},Any}) at Conversion.jl:35
#eigs#725 at eigs.jl:11 [inlined]
eigs(::ApproxFun.InterlaceOperator{Float64,1,ApproxFun.Chebyshev{ApproxFun.Segment{Float64},Float64},ApproxFun.ArraySpace{ApproxFun.Space,1,ApproxFun.Point{Float64},Any},ApproxFun.CachedIterator{Tuple{Int64,Int64},ApproxFun.BlockInterlacer{Tuple{ApproxFun.Repeated{Bool}}},Tuple{Int64,Int64,Tuple{Void},Tuple{Int64}}},ApproxFun.CachedIterator{Tuple{Int64,Int64},ApproxFun.BlockInterlacer{Tuple{Array{Int64,1},ApproxFun.Repeated{Bool}}},Tuple{Int64,Int64,Tuple{Int64,Void},Tuple{Int64,Int64}}},Tuple{ApproxFun.Infinity{Bool},ApproxFun.Infinity{Bool}}}, ::Int64) at eigs.jl:8
include_string(::String, ::String) at loading.jl:522
include_string(::String, ::String, ::Int64) at eval.jl:30
include_string(::Module, ::String, ::String, ::Int64, ::Vararg{Int64,N} where N) at eval.jl:34
(::Atom.##100#105{String,Int64,String})() at eval.jl:75
withpath(::Atom.##100#105{String,Int64,String}, ::String) at utils.jl:30
withpath(::Function, ::String) at eval.jl:38
hideprompt(::Atom.##99#104{String,Int64,String}) at repl.jl:66
macro expansion at eval.jl:73 [inlined]
(::Atom.##98#103{Dict{String,Any}})() at task.jl:80

Does anyone have any ideas on how to do this with the current version of ApproxFun?


#5

I think I changed the override to eigs(Bcs, A, n), that is, give the bcs seperate from the operator. I’ll check it works in the morning.


#6

Yes, eigs will accept arguments of that form. But in this case, it throws an error because size(Dirichlet(),1) == ∞.


#7

A work around was posted here https://github.com/JuliaApproximation/ApproxFun.jl/issues/564#issuecomment-373555416


#8

Thanks very much!

Sorry to pester you but I was looking at some examples used in chebfun and it seems that they’re able to get eigenvalues for all sorts of funny potentials, most of which don’t seem to work for me with ApproxFun.

For instance, for the absolute value potential,

S = Chebyshev(Interval(-10,10))
B = Dirichlet(S)
D = Derivative(S)
V = Fun(x -> abs(x), Interval(-10, 10))
λ, ψ = eigs(B, -D^2/2 + V, 200)

returns no eigenvalues or eigenvectors even though they must exist.
I also tried various ways of making abs(x) better behaved like

V = Fun(x -> sqrt(x^2), Interval(-10, 10))

and

V = Fun(x -> sqrt(x^2 +0.1), Interval(-10, 10))

but to no avail. Do you have any ideas about what could cause this?


#9

Don’t use the constructor for functions with singularities, and specify a tolerance for eigs. This works pretty well form me:

x = Fun(-20 .. 20)
V = abs(x)
D = Derivative()
S = space(V)
B = [Dirichlet(S); continuity(S,0:1)]
λ, v = eigs(B, -D^2 + V, 500,tolerance=1E-10)

p = plot(V; legend=false)
    for k=1:20
        plot!(real(v[k]/norm(v[k]) + λ[k]))
    end
    p


#10

Ah, I see! Maybe I’ll spend some time with this and shoot you a pull request to add something about this to the documentation.


#11

:+1: also think about adding an example to https://github.com/JuliaApproximation/ApproxFunExamples


#12

Just an FYI, the above code doesn’t seem to work on the current tagged version of Approxfun but works fine on the Master branch.