Logkernel evaluated at ±1 (#52)

* Support logkernel evaluate at ±1

* v0.3.6

Author: dlfivefifty

Project.toml

@@ -1,7 +1,7 @@
 name = "SingularIntegrals"
 uuid = "d7440221-8b5e-42fc-909c-0567823f424a"
 authors = ["Sheehan Olver <solver@mac.com>"]
-version = "0.3.5"
+version = "0.3.6"
 
 [deps]
 ArrayLayouts = "4c555306-a7a7-4459-81d9-ec55ddd5c99a"

src/logkernel.jl

@@ -42,7 +42,7 @@ end
 """
     logkernel(P)
 
-applies the log kernel log(abs(x-t)) to the columns of a quasi matrix, i.e., `(log.(abs.(x - x')) * P)/π`
+applies the log kernel log(abs(x-t)) to the columns of a quasi matrix, i.e., `(log.(abs.(x - x')) * P)`
 """
 logkernel(P, z...) = logkernel_layout(MemoryLayout(P), P, z...)
 
@@ -145,10 +145,12 @@ function complexlogkernel_recurrence(wP::Weighted{<:Any,<:ChebyshevU})
     (2*(n .+ 2)) ./ (n .+ 3), Zeros{R}(∞), (n .+ 2) .* (n .- 1) ./ (n .* (n .+ 3))
 end
 
+
+zlog(z) = ifelse(iszero(z), z, z*log(z))
 function complexlogkernel_ics(P::Weighted{<:Any,<:Legendre}, z::Number)
     T = promote_type(eltype(P), typeof(z))
-    r0 = (1 + z)log(1 + z) - (z-1)log(z-1) - 2one(T)
-    r1 = (z+1)*r0/2 + 1 - (z+1)log(z+1)
+    r0 = zlog(1 + z) - zlog(z-1) - 2one(T)
+    r1 = (z+1)*r0/2 + 1 - zlog(z+1)
     r2 = z*r1 + 2*one(T)/3
     r0,r1,r2
 end

test/test_logkernel.jl

@@ -26,6 +26,13 @@ using ClassicalOrthogonalPolynomials: affine
         @test complexlogkernel(exp.(x), 2 + im) ≈ sum(log.((2+im) .- x) .* exp.(x))
         @test_throws ArgumentError complexlogkernel(Jacobi(0.1,0.2), 2+im)
     end
+
+    @testset "endpoints" begin
+        @test complexlogkernel(exp.(x), 1) ≈ logkernel(exp.(x), 1) ≈ -1.9569548200977875
+        @test real(complexlogkernel(exp.(x), -1+0im)) ≈ logkernel(exp.(x), -1) ≈ 0.27395419528475684
+        @test complexlogkernel(exp.(x), -1+0.0im) ≈ 0.27395419528475684+7.38400687288264im
+        @test complexlogkernel(exp.(x), -1-0.0im) ≈ 0.27395419528475684-7.38400687288264im
+    end
 end
 
 @testset "LogKernelPoint" begin
@@ -37,10 +44,12 @@ end
         @test L isa LogKernelPoint{Float64,ComplexF64,ComplexF64,Float64,ChebyshevInterval{Float64}}
         @test (L * wU)[1:5] ≈ [  -1.2919202947616695, -0.20965486677056738, 0.6799687631764493, 0.13811497572177128, -0.2289481463304956]
 
-        @test L * (wU / wU \ @.(exp(x) * sqrt(1 - x^2))) ≈ -1.4812979070884382
+        fU = wU / wU \ @.(exp(x) * sqrt(1 - x^2))
+        @test L * fU ≈ logkernel(fU, z) ≈ -1.4812979070884382
 
         wT = Weighted(ChebyshevT())
-        @test L * (wT / wT \ @.(exp(x) / sqrt(1 - x^2))) ≈ -1.9619040529776954 #mathematica
+        fT = (wT / wT \ @.(exp(x) / sqrt(1 - x^2)))
+        @test L * fT ≈ logkernel(fT, z) ≈ -1.9619040529776954 #mathematica
     end
 
     @testset "Real point" begin