diff --git a/src/ellip.jl b/src/ellip.jl index 4bb43548..65c9cad9 100644 --- a/src/ellip.jl +++ b/src/ellip.jl @@ -36,6 +36,9 @@ As suggested in this paper, the domain is restricted to ``(-\infty,1]``. ellipk(m::Real) = _ellipk(float(m)) function _ellipk(m::Float64) + isnan(m) && return NaN + (m === -Inf) && return 0.0 + flag_is_m_neg = false if m < 0.0 x = m / (m-1) #dealing with negative args @@ -49,7 +52,7 @@ function _ellipk(m::Float64) return Float64(halfπ) elseif x == 1.0 - return Inf + return flag_is_m_neg ? 0.0 : Inf elseif x > 1.0 throw(DomainError(m, "`m` must lie between -Inf and 1 ---- Domain: (-Inf,1.0]")) @@ -207,6 +210,9 @@ As suggested in this paper, the domain is restricted to ``(-\infty,1]``. ellipe(m::Real) = _ellipe(float(m)) function _ellipe(m::Float64) + isnan(m) && return NaN + (m === -Inf) && return Inf + flag_is_m_neg = false if m < 0.0 x = m / (m-1) #dealing with negative args @@ -219,7 +225,7 @@ function _ellipe(m::Float64) if x == 0.0 return Float64(halfπ) elseif x == 1.0 - return 1.0 + return flag_is_m_neg ? Inf : 1.0 elseif x > 1.0 throw(DomainError(m,"`m` must lie between -Inf and 1 ---- Domain : (-inf,1.0]")) diff --git a/test/ellip.jl b/test/ellip.jl index 89033d2a..b22c4fd7 100644 --- a/test/ellip.jl +++ b/test/ellip.jl @@ -17,6 +17,9 @@ @test ellipk(1.0) == Inf @test ellipk(Float16(0.92)) ≈ 2.683551406315229344 rtol=2*eps(Float16) @test ellipk(Float32(0.92)) ≈ 2.683551406315229344 rtol=2*eps(Float32) + @test isnan(ellipk(NaN)) + @test ellipk(-Inf) == 0.0 + @test ellipk(-1e30) ≈ 0.0 atol=1e-13 @test_throws MethodError ellipk(BigFloat(0.5)) @test_throws DomainError ellipk(1.1) end @@ -38,6 +41,9 @@ @test ellipe(1.0) ≈ 1.00 rtol=2*eps() @test ellipe(Float16(0.865)) ≈ 1.1322436887003925 rtol=2*eps(Float16) @test ellipe(Float32(0.865)) ≈ 1.1322436887003925 rtol=2*eps(Float32) + @test isnan(ellipe(NaN)) + @test ellipe(-Inf) == Inf + @test ellipe(-1e16) > ellipe(-1e15) @test_throws MethodError ellipe(BigFloat(-1)) @test_throws DomainError ellipe(1.2) end