upgrade for v0.7 and improve code

Author: jlapeyre

.travis.yml

@@ -3,11 +3,9 @@ os:
   - linux
   - osx
 julia:
-  - 0.6
-#  - nightly
+  - 0.7
+  - nightly
 notifications:
   email: false
 after_success:
-  - julia -e 'Pkg.add("Documenter")'
-  - julia -e 'cd(Pkg.dir("InverseLaplace")); include(joinpath("docs", "make.jl"))'
   - julia -e 'Pkg.add("Coverage"); cd(Pkg.dir("InverseLaplace")); using Coverage; Coveralls.submit(process_folder()); Codecov.submit(process_folder())'

LICENSE.md

@@ -1,6 +1,6 @@
 The InverseLaplace.jl package is licensed under the MIT "Expat" License:
 
-> Copyright (c) 2015 2016: John Lapeyre.
+> Copyright (c) 2015-2018: John Lapeyre.
 >
 > Permission is hereby granted, free of charge, to any person obtaining
 > a copy of this software and associated documentation files (the

README.md

@@ -5,13 +5,15 @@
 [![](https://img.shields.io/badge/docs-latest-blue.svg)](https://jlapeyre.github.io/InverseLaplace.jl/latest)
 Linux, OSX: [![Build Status](https://travis-ci.org/jlapeyre/InverseLaplace.jl.svg)](https://travis-ci.org/jlapeyre/InverseLaplace.jl)   Windows: [![Build Status](https://ci.appveyor.com/api/projects/status/github/jlapeyre/InverseLaplace.jl?branch=master&svg=true)](https://ci.appveyor.com/project/jlapeyre/inverselaplace-jl)       [![Coverage Status](https://coveralls.io/repos/github/jlapeyre/InverseLaplace.jl/badge.svg?branch=master)](https://coveralls.io/github/jlapeyre/InverseLaplace.jl?branch=master) [![codecov](https://codecov.io/gh/jlapeyre/InverseLaplace.jl/branch/master/graph/badge.svg)](https://codecov.io/gh/jlapeyre/InverseLaplace.jl)
 
-
 This implements some numerical methods for computing inverse Laplace transforms in Julia.
 
+InverseLaplace v0.1.0 is the last version that supports Julia v0.6
+
 See the documentation https://jlapeyre.github.io/InverseLaplace.jl/latest .
 
 Note: A small part of `InverseLaplace.jl` depends on `Optim.jl`, which is currently
-broken on Julia v0.7. So `InverseLaplace.jl` is broken on Julia v0.7.
-
+broken on Julia v0.7. So optimization is disabled for `InverseLaplace.jl` versions
+greater than v0.1.0
 
-Note: the last version of this module supporting Julia v0.4 is tagged v0.0.2
+Note: the last version of this module supporting Julia v0.4 is tagged v0.0.2.
+the last version of this module supporting Julia v0.6 is tagged v0.1.0.

REQUIRE

@@ -1,3 +1,5 @@
-julia 0.6
-Optim
-Compat 0.17.0
+julia 0.7-beta
+SpecialFunctions
+# Optim
+AbstractFFTs
+FFTW

appveyor.yml

@@ -1,9 +1,9 @@
 environment:
   matrix:
-  - JULIAVERSION: "julialang/bin/winnt/x86/0.6/julia-0.6-latest-win32.exe"
-  - JULIAVERSION: "julialang/bin/winnt/x64/0.6/julia-0.6-latest-win64.exe"
-  # - JULIAVERSION: "julianightlies/bin/winnt/x86/julia-latest-win32.exe"
-  # - JULIAVERSION: "julianightlies/bin/winnt/x64/julia-latest-win64.exe"
+  - JULIA_URL: "https://julialang-s3.julialang.org/bin/winnt/x86/0.7/julia-0.7-latest-win32.exe"
+  - JULIA_URL: "https://julialang-s3.julialang.org/bin/winnt/x64/0.7/julia-0.7-latest-win64.exe"
+  - JULIA_URL: "https://julialangnightlies-s3.julialang.org/bin/winnt/x86/julia-latest-win32.exe"
+  - JULIA_URL: "https://julialangnightlies-s3.julialang.org/bin/winnt/x64/julia-latest-win64.exe"
 
 branches:
   only:
@@ -17,9 +17,15 @@ notifications:
     on_build_status_changed: false
 
 install:
+  - ps: "[System.Net.ServicePointManager]::SecurityProtocol = [System.Net.SecurityProtocolType]::Tls12"
+# If there's a newer build queued for the same PR, cancel this one
+  - ps: if ($env:APPVEYOR_PULL_REQUEST_NUMBER -and $env:APPVEYOR_BUILD_NUMBER -ne ((Invoke-RestMethod `
+        https://ci.appveyor.com/api/projects/$env:APPVEYOR_ACCOUNT_NAME/$env:APPVEYOR_PROJECT_SLUG/history?recordsNumber=50).builds | `
+        Where-Object pullRequestId -eq $env:APPVEYOR_PULL_REQUEST_NUMBER)[0].buildNumber) { `
+        throw "There are newer queued builds for this pull request, failing early." }
 # Download most recent Julia Windows binary
   - ps: (new-object net.webclient).DownloadFile(
-        $("http://s3.amazonaws.com/"+$env:JULIAVERSION),
+        $env:JULIA_URL,
         "C:\projects\julia-binary.exe")
 # Run installer silently, output to C:\projects\julia
   - C:\projects\julia-binary.exe /S /D=C:\projects\julia
@@ -31,4 +37,4 @@ build_script:
       Pkg.clone(pwd(), \"InverseLaplace\"); Pkg.build(\"InverseLaplace\")"
 
 test_script:
-  - C:\projects\julia\bin\julia --check-bounds=yes -e "Pkg.test(\"InverseLaplace\")"
+  - C:\projects\julia\bin\julia -e "Pkg.test(\"InverseLaplace\")"

docs/src/index.md

@@ -5,6 +5,8 @@
 The source repository is [https://github.com/jlapeyre/InverseLaplace.jl](https://github.com/jlapeyre/InverseLaplace.jl).
 
 This package provides three numerical algorithms for numerically inverting Laplace transforms.
+InverseLaplace v0.1.0 is the last version that supports Julia v0.6.
+Optimization of the Weeks method is temporarily disabled for Julia v0.7.
 
 ## Contents
 

src/InverseLaplace.jl

@@ -1,31 +1,31 @@
-VERSION >= v"0.4.0-dev+6521" && __precompile__()
+__precompile__()
 
 module InverseLaplace
 
-using Compat
-
+import SpecialFunctions
 # I have to export everything in order for Documenter.jl to find
 # the strings. A few hours of work would probably be enough to solve the problem.
 export ILt, setNterms
-export Weeks, WeeksErr, optimize, opteval, setparameters
+# export optimize, opteval # Broken at the moment
+export Weeks, WeeksErr, setparameters
 export Talbot, GWR, ILT
 export TransformPair, ILtPair, abserr, iltpair_power
 export ilt, talbot, gwr
 
-@compat abstract type AbstractILt end
+abstract type AbstractILt end
 
 """
    ft::Talbot = Talbot(func::Function, Nterms::Integer=32)
 
 return `ft`, which estimates the inverse Laplace transform of `func` with
 the fixed Talbot algorithm. `ft(t)` evaluates the transform at `t`.  You may
-want to tune `Nterms` together with `setprecision(BigFloat,x)`.
+want to tune `Nterms` together with `setprecision(BigFloat, x)`.
 
 # Example
 
 Compute the inverse transform of the transform of `cos` at argument `pi/2`.
 ```julia-repl
-julia> ft = Talbot(s -> s/(s^2+1), 80);
+julia> ft = Talbot(s -> s/(s^2 + 1), 80);
 
 julia> ft(pi/2)
 -3.5366510684573195e-5
@@ -42,12 +42,13 @@ julia> Float64(ft(big(pi)/2))
     function for complex arguments. The GWR method is, in general, less accurate and
     less stable, but does not evaluate the Laplace transform function for complex arguments.
 """
-type Talbot{T<:Base.Callable} <: AbstractILt
-    LSfunc::T  # Laplace space function
+struct Talbot{T<:Base.Callable} <: AbstractILt
+    Laplace_space_func::T  # Laplace space function
     Nterms::Int
 end
 
-Talbot(LSfunc::Base.Callable) = Talbot(LSfunc,32)
+const talbot_default_num_terms = 32
+Talbot(Laplace_space_func::Base.Callable) = Talbot(Laplace_space_func, talbot_default_num_terms)
 
 """
    ft::GWR = GWR(func::Function, Nterms::Integer=16)
@@ -59,18 +60,19 @@ the GWR algorithm. `ft(t)` evaluates the transform at `t`.
 
 Compute the inverse transform of the transform of `cos` at argument `pi/2`.
 ```
-julia> ft = GWR(s -> s/(s^2+1), 16);
+julia> ft = GWR(s -> s/(s^2 + 1), 16);
 
 julia> ft(pi/2)
 -0.001
 ```
 """
-type GWR{T<:Base.Callable} <: AbstractILt
-    LSfunc::T  # Laplace space function
+struct GWR{T<:Base.Callable} <: AbstractILt
+    Laplace_space_func::T  # Laplace space function
     Nterms::Int
 end
 
-GWR(LSfunc::Base.Callable) = GWR(LSfunc,16)
+const gwr_default_num_terms = 16
+GWR(Laplace_space_func::Base.Callable) = GWR(Laplace_space_func, gwr_default_num_terms)
 
 """
     ILT(function, Nterms=32)
@@ -80,14 +82,13 @@ This is an alias for the default `Talbot()` method.
 ILT(args...) = Talbot(args...)
 
 function Base.show(io::IO, ailt::AbstractILt)
-    print(io, string(typeof(ailt)), "(Nterms=", ailt.Nterms,')')
+    print(io, string(typeof(ailt)), "(Nterms=", ailt.Nterms, ')')
 end
 
 # TODO get rid of this in favor of above.
 # Rely on broadcasting, as well.
-
-type ILt{T<:Base.Callable, V<:Base.Callable} <: AbstractILt
-    LSfunc::T
+struct ILt{T<:Base.Callable, V<:Base.Callable} <: AbstractILt
+    Laplace_space_func::T
     iltfunc::V
     Nterms::Int
 end
@@ -109,7 +110,7 @@ of `32` may give extremely inaccurate estimates.
 ## Example
 
 ```julia
-julia> itr = ILt( s -> 1/(1+s^2), talbot);
+julia> itr = ILt( s -> 1/(1 + s^2), talbot);
 
 julia> itr([ pi/4, pi/2, 3*pi/4, -pi])
 4-element Array{Float64,1}:
@@ -119,7 +120,7 @@ julia> itr([ pi/4, pi/2, 3*pi/4, -pi])
  -3.66676e-12
 ```
 """
-ILt(func,iltfunc) = ILt(func, iltfunc, 32)
+ILt(func, iltfunc) = ILt(func, iltfunc, talbot_default_num_terms)
 
 # Make this the default
 """
@@ -128,8 +129,7 @@ ILt(func,iltfunc) = ILt(func, iltfunc, 32)
 `ilt` is an alias for the default inverse Laplace transform method `talbot`.
 """
 ilt(args...) = talbot(args...)
-ILt(func) = ILt(func, talbot, 32)
-
+ILt(func) = ILt(func, talbot, talbot_default_num_terms)
 
 """
     setNterms(ailt::AbstractILt, Nterms::Integer)
@@ -140,18 +140,19 @@ calls `ailt(t)` reflect the new value of `Nterms`.
 """
 setNterms(ailt::AbstractILt, N::Integer) = (ailt.Nterms = N)
 
-(ailt::ILt)(t) = ailt.iltfunc(ailt.LSfunc, t, ailt.Nterms)
-(ailt::ILt)(t,N) = ailt.iltfunc(ailt.LSfunc, t, N)
+## Make the ILT data types callable.
+(ailt::ILt)(t) = ailt.iltfunc(ailt.Laplace_space_func, t, ailt.Nterms)
+(ailt::ILt)(t,N) = ailt.iltfunc(ailt.Laplace_space_func, t, N)
 
-(ailt::Talbot)(t) = talbot(ailt.LSfunc, t, ailt.Nterms)
-(ailt::Talbot)(t, n::Integer) = talbot(ailt.LSfunc, t, n)
+(ailt::Talbot)(t) = talbot(ailt.Laplace_space_func, t, ailt.Nterms)
+(ailt::Talbot)(t, n::Integer) = talbot(ailt.Laplace_space_func, t, n)
 
-(ailt::GWR)(t) = gwr(ailt.LSfunc, t, ailt.Nterms)
-(ailt::GWR)(t, n::Integer) = gwr(ailt.LSfunc, t, n)
+(ailt::GWR)(t) = gwr(ailt.Laplace_space_func, t, ailt.Nterms)
+(ailt::GWR)(t, n::Integer) = gwr(ailt.Laplace_space_func, t, n)
 
 include("fixed_talbot.jl")
 include("gwr.jl")
 include("weeks.jl")
-include("test.jl")
+include("pairtest.jl")
 
 end # module

src/fixed_talbot.jl

@@ -1,10 +1,11 @@
+# Fixed Talbot method
+#
 # Abate, J. and Valkó, P.P.
 # Multi-precision Laplace transform inversion
 # International Journal for Numerical Methods in Engineering, Vol. 60 (Iss. 5-7)  2004  pp 979–993
-# Fixed Talbot method
 
 """
-    talbot(func::Function, t::AbstractFloat, M::Integer=32)
+    talbot(func::Function, t::AbstractFloat, M::Integer=talbot_default_num_terms)
 
 Evaluate the inverse Laplace transform of `func` at the point `t`. Use `M` terms in the algorithm.
 For `typeof(t)` is `Float64`, the default for `M` is `32`. For `BigFloat` the default is `64`.
@@ -14,38 +15,38 @@ If `BigFloat` precision is larger than default, try increasing `M`. `talbot is v
 # Example
 
 ```jldoctest
-julia> InverseLaplace.talbot( s -> 1/s^3,  3)
+julia> InverseLaplace.talbot( s -> 1/s^3, 3)
 4.50000000000153
 ```
 
 !!! note
     This function uses the fixed Talbot method. It evaluates `func` for complex arguments.
 """
 function talbot(func, t, M)
-    bM = convert(typeof(t),M)
-    r = (2 * bM) / (5*t)
-    term = (1//2) * exp(r*t) * func(r)
+    bM = convert(typeof(t), M)
+    r = (2 * bM) / (5 * t)
+    term = (1//2) * exp(r * t) * func(r)
     for i in 1:M-1
         theta = i * (pi/bM)
-        s = r*theta*(complex(cot(theta),one(theta)))
-        sigma = theta + (theta*cot(theta)-1)*cot(theta)
-        term += real(exp(t*s) * complex(one(t),sigma) * func(s))
+        s = r * theta * (complex(cot(theta), one(theta)))
+        sigma = theta + (theta * cot(theta) - 1) * cot(theta)
+        term += real(exp(t * s) * complex(one(t), sigma) * func(s))
     end
-    return term * 2 / (5*t)
+    return term * 2 / (5 * t)
 end
 
-talbot(func,t) = talbot(func,t,32)
-talbot(func,t::BigFloat) = talbot(func,t,64)
-talbot(func,t::Integer,args...) = talbot(func,BigFloat(t),args...)
+talbot(func, t) = talbot(func, t, talbot_default_num_terms)
+const talbot_BigFloat_default_num_terms = 64
+talbot(func, t::BigFloat) = talbot(func, t, talbot_BigFloat_default_num_terms)
+talbot(func, t::Integer, args...) = talbot(func, BigFloat(t), args...)
 # Hmm, at some point, one of these routines actually gave a Rational result. Don't recall how.
 # But, it can't be talbot.
-talbot(func,t::Rational,args...) = talbot(func,BigFloat(t),args...)
+talbot(func, t::Rational, args...) = talbot(func, BigFloat(t), args...)
 
 # Operate on an array of values of t. A single function evaluation
 # f(s) is used for all t
 # This gives more inaccurate results the further values of
 # t are from tmax
-
 """
     talbotarr(func, ta::AbstractArray, M)
 
@@ -56,31 +57,28 @@ than a vectorized application of `talbot`, but is in general, less accurate.
 """
 function talbotarr(func, t::AbstractArray, M)
     tt = typeof(t[1])
-    bM = convert(tt,M)
+    bM = convert(tt, M)
     terms = similar(t)
     tmax = maximum(t)
-    r = (2 * bM) / (5*tmax)
+    r = (2 * bM) / (5 * tmax)
     fr = (1//2) * func(r)
     for j in 1:length(terms)
-        terms[j] =  exp(r*t[j]) * fr
+        terms[j] =  exp(r * t[j]) * fr
     end
     for i in 1:M-1
         theta = i * (pi/bM)
-        s = r*theta*(complex(cot(theta),one(theta)))
-        sigma = theta + (theta*cot(theta)-1)*cot(theta)
-        fs = complex(one(tt),sigma) * func(s)
+        s = r * theta * (complex(cot(theta), one(theta)))
+        sigma = theta + (theta * cot(theta) - 1) * cot(theta)
+        fs = complex(one(tt), sigma) * func(s)
         for j in 1:length(terms)
-            terms[j] += real(exp(t[j]*s) * fs)
+            terms[j] += real(exp(t[j] * s) * fs)
         end
     end
     for j in 1:length(terms)
-        terms[j] *= 2/(5*tmax)
+        terms[j] *= 2/(5 * tmax)
     end
     return terms
 end
 
-talbotarr(func,t::AbstractArray) = talbotarr(func,t,32)
-talbotarr(func,t::Vector{BigFloat}) = talbotarr(func,t,64)
-#talbotarr(func,t::T) = talbotarr(func,t,64) where T <: AbstractArray{V} where V <: BigFloat
-
-####
+talbotarr(func, t::AbstractArray) = talbotarr(func, t, talbot_default_num_terms)
+talbotarr(func, t::Vector{BigFloat}) = talbotarr(func, t, talbot_BigFloat_default_num_terms)