Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

logkernel with weighted ChebyshevT #38

Closed
kburns opened this issue Feb 20, 2024 · 3 comments
Closed

logkernel with weighted ChebyshevT #38

kburns opened this issue Feb 20, 2024 · 3 comments

Comments

@kburns
Copy link

kburns commented Feb 20, 2024

I'm trying a simple modification of one of the readme examples, applying a log kernel to a weighted ChebyshevT function. This fails but I can't quite tell what's going on under the hood. Any pointers?

W = Weighted(ChebyshevT())
x = axes(W, 1)
f = expand(W, x -> exp(x) / sqrt(1-x^2))
log.(abs.(0.3 .- x')) * f

returns

ERROR: not implemented
Stacktrace:
  [1] error(s::String)
    @ Base ./error.jl:35
  [2] logkernel_layout(::ClassicalOrthogonalPolynomials.WeightedOPLayout{…}, P::Weighted{…}, z::Float64)
    @ SingularIntegrals ~/.julia/packages/SingularIntegrals/LB9Wj/src/logkernel.jl:46
  [3] logkernel(P::Weighted{Float64, ChebyshevT{Float64}}, z::Float64)
    @ SingularIntegrals ~/.julia/packages/SingularIntegrals/LB9Wj/src/logkernel.jl:32
  [4] logkernel_layout(LAY::LazyArrays.ApplyLayout{…}, V::QuasiArrays.ApplyQuasiVector{…}, y::Float64)
    @ SingularIntegrals ~/.julia/packages/SingularIntegrals/LB9Wj/src/logkernel.jl:51
  [5] logkernel_layout(::ContinuumArrays.ExpansionLayout{…}, A::QuasiArrays.ApplyQuasiVector{…}, dims::Float64)
    @ SingularIntegrals ~/.julia/packages/SingularIntegrals/LB9Wj/src/logkernel.jl:54
  [6] logkernel(P::QuasiArrays.ApplyQuasiVector{Float64, typeof(*), Tuple{Weighted{…}, LazyArrays.ApplyArray{…}}}, z::Float64)
    @ SingularIntegrals ~/.julia/packages/SingularIntegrals/LB9Wj/src/logkernel.jl:32
  [7] macro expansion
    @ ~/.julia/packages/SingularIntegrals/LB9Wj/src/logkernel.jl:14 [inlined]
  [8] mul
    @ ~/.julia/packages/ContinuumArrays/amUIf/src/operators.jl:34 [inlined]
  [9] *(A::QuasiArrays.BroadcastQuasiMatrix{…}, B::QuasiArrays.ApplyQuasiVector{…})
    @ QuasiArrays ~/.julia/packages/QuasiArrays/ZCn0F/src/matmul.jl:23
 [10] top-level scope
    @ REPL[108]:1
Some type information was truncated. Use `show(err)` to see complete types.
@dlfivefifty
Copy link
Member

Just means it hasn't been implemented yet. I guess I only got around to weighted Chebyshev U... There is a general approach for classical OPs but need to add it.

If you only need to evaluate the log kernel on the interval you can do this:

julia> using ClassicalOrthogonalPolynomials, SingularIntegrals

julia> W = Weighted(ChebyshevT())
Weighted(ChebyshevT())

julia> x = axes(W, 1)
Inclusion(-1.0 .. 1.0 (Chebyshev))

julia> f = expand(W, x -> exp(x) / sqrt(1-x^2))
Weighted(ChebyshevT()) * [1.2660658777520084, 1.1303182079849703, 0.27149533953407656, 0.044336849848663776, 0.0054742404420936785, 0.0005429263119139354, 4.497732295427654e-5, 3.198436462529006e-6, 1.992124804817033e-7, 1.1036771880698548e-8    ]

julia> x = axes(W, 1)
Inclusion(-1.0 .. 1.0 (Chebyshev))

julia> @time Lf = log.(abs.(x .- x')) * f;
  0.000020 seconds (7 allocations: 496 bytes)

julia> @time Lf[0.3]
  0.000008 seconds (1 allocation: 16 bytes)
-3.437622681488403

@kburns
Copy link
Author

kburns commented Feb 26, 2024

Thanks. How is (log.(abs.(x .- x')) * f)[0.3] different under the hood from (log.(abs.(0.3 .- x')) * f)?

@dlfivefifty
Copy link
Member

The former deduces the Chebyshev coefficients by applying a diagonal matrix:

julia> T = ChebyshevT(); W = Weighted(T);

julia> T \ logkernel(W)
ℵ₀×ℵ₀ LinearAlgebra.Diagonal{Float64, LazyArrays.ApplyArray{Float64, 1, typeof(vcat), Tuple{Float64, LazyArrays.BroadcastVector{Float64, typeof(/), Tuple{Float64, InfiniteArrays.InfUnitRange{Int64}}}}}} with indices OneToInf()×OneToInf():
 -2.17759                             
          -3.14159                     
                   -1.5708             
                           -1.0472     
                                      
                                     
                                      
                                      
                                      
                                      
                                        

(I don't know if I documented it but log.(abs(x .- x')) * f is a synonym for logkernel(f). )

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants