diff --git a/.travis.yml b/.travis.yml index 4423298..2a89581 100644 --- a/.travis.yml +++ b/.travis.yml @@ -3,6 +3,5 @@ julia: - nightly - 1.4 - 1.5 - - 1.6 os: - linux diff --git a/Project.toml b/Project.toml index 161009d..15a46c5 100644 --- a/Project.toml +++ b/Project.toml @@ -1,7 +1,7 @@ name = "FunctionIntegrator" uuid = "7685536e-2581-4f83-bef1-2ba363c9cb91" authors = ["Brenton Horne "] -version = "0.5.1" +version = "0.6.0" [deps] FastGaussQuadrature = "442a2c76-b920-505d-bb47-c5924d526838" diff --git a/src/Rectangle.jl b/src/Rectangle.jl index e81a60a..758ea86 100644 --- a/src/Rectangle.jl +++ b/src/Rectangle.jl @@ -13,7 +13,7 @@ function rectangle_rule_left(f::Function, N::Number, a::Number, b::Number) y = 0; x = a; for i=1:N - y += h*f(x) + y += h*f(x); x += h; end return y @@ -34,7 +34,7 @@ function rectangle_rule_midpoint(f::Function, N::Number, a::Number, b::Number) y = 0; x = a+h/2; for i=1:N - y += h*f(x) + y += h*f(x); x += h; end return y @@ -55,7 +55,7 @@ function rectangle_rule_right(f::Function, N::Number, a::Number, b::Number) y = 0; x = a+h; for i=1:N - y += h*f(x) + y += h*f(x); x += h; end return y diff --git a/src/Simpson.jl b/src/Simpson.jl index ed15d48..a18a348 100644 --- a/src/Simpson.jl +++ b/src/Simpson.jl @@ -1,11 +1,5 @@ -function stepwise_simpsons(f::Function, h::Number, x::Number, i::Integer, N::Number) - if i == 1 || i == N + 1 - return h / 3 * f(x) - elseif (i % 2) == 1 - return 2 * h / 3 * f(x) - else - return 4 * h / 3 * f(x) - end +function stepwise_simpsons(f::Function, h::Number, x::Number) + return h/6 * (f(x) + 4*f(x+h/2) + f(x+h)); end function stepwise_simpsons38(f::Function, h::Number, x::Number, i::Integer, N::Number) @@ -28,15 +22,12 @@ uses [Simpson's rule](https://en.wikipedia.org/wiki/Simpson%27s_rule) to approxi """ function simpsons_rule(f::Function, N::Number, a::Number, b::Number) N = convert(Int64, N); - iseven(N) || error("N must be even in order for Simspon's rule to work properly.") h = (b-a)/N; y = 0; x = a; - for i=1:N+1 - y = y + stepwise_simpsons(f, h, x, i, N); - if i < N+1 - x += h; - end + for i=1:N + y += stepwise_simpsons(f, h, x); + x += h; end return y end diff --git a/src/Trapezoidal.jl b/src/Trapezoidal.jl index 37639c5..a34dcca 100644 --- a/src/Trapezoidal.jl +++ b/src/Trapezoidal.jl @@ -1,9 +1,5 @@ function stepwise_trapezoidal(f, h, x, i, N) - if i == 1 || i == N + 1 - return h / 2 * f(x) - else - return h * f(x) - end + return h/2 * (f(x) + f(x+h)); end """ @@ -20,11 +16,9 @@ function trapezoidal_rule(f::Function, N::Number, a::Number, b::Number) h = (b-a)/N; y = 0; x = a; - for i=1:N+1 - y = y + stepwise_trapezoidal(f, h, x, i, N); - if i < N+1 - x += h; - end + for i=1:N + y += stepwise_trapezoidal(f, h, x, i, N); + x += h; end return y end diff --git a/test/airyai.jl b/test/airyai.jl index aadbe64..4ac3f48 100644 --- a/test/airyai.jl +++ b/test/airyai.jl @@ -30,7 +30,7 @@ printstyled("Performing the Airy Ai(x) test, where Ai(x) is integrated on the se printstyled("Running: rombergs_method\n"; color = :magenta) @time @test rombergs_method(x -> airyai(x), 9, 0, 100) ≈ 1.0/3.0 printstyled("Running: simpsons_rule\n"; color = :magenta) - @time @test simpsons_rule(x -> airyai(x), 2512, 0, 100) ≈ 1.0/3.0 + @time @test simpsons_rule(x -> airyai(x), 1256, 0, 100) ≈ 1.0/3.0 printstyled("Running: simpsons38_rule\n"; color = :magenta) @time @test simpsons38_rule(x -> airyai(x), 3075, 0, 100) ≈ 1.0/3.0 printstyled("Running: trapezoidal_rule\n"; color = :magenta) diff --git a/test/besselj.jl b/test/besselj.jl index fe81769..e79f42f 100644 --- a/test/besselj.jl +++ b/test/besselj.jl @@ -32,7 +32,7 @@ printstyled("Performing the besselj test, where BesselJ_1(2) is approximated and printstyled("Running: rombergs_method\n"; color = :magenta) @time @test rombergs_method(besselj_integrand, 5, 0, pi/2) ≈ besselj(1,2) printstyled("Running: simpsons_rule\n"; color = :magenta) - @time @test simpsons_rule(besselj_integrand, 6, 0, pi/2) ≈ besselj(1,2) + @time @test simpsons_rule(besselj_integrand, 3, 0, pi/2) ≈ besselj(1,2) printstyled("Running: simpsons38_rule\n"; color = :magenta) @time @test simpsons38_rule(besselj_integrand, 9, 0, pi/2) ≈ besselj(1,2) printstyled("Running: trapezoidal_rule\n"; color = :magenta) diff --git a/test/cos_cot_integral.jl b/test/cos_cot_integral.jl index 14a0a01..f2a00e1 100644 --- a/test/cos_cot_integral.jl +++ b/test/cos_cot_integral.jl @@ -34,7 +34,7 @@ printstyled("Integrating cos^2(x)/(1+cot(x)) from 0 to pi/2 and comparing the re printstyled("Running: rombergs_method\n"; color = :magenta) @time @test rombergs_method(cos_cot_fn, 4, 0, pi/2) ≈ 0.25 printstyled("Running: simpsons_rule\n"; color = :magenta) - @time @test simpsons_rule(cos_cot_fn, 78, 0, pi/2) ≈ 0.25 + @time @test simpsons_rule(cos_cot_fn, 39, 0, pi/2) ≈ 0.25 printstyled("Running: simpsons38_rule\n"; color = :magenta) @time @test simpsons38_rule(cos_cot_fn, 96, 0, pi/2) ≈ 0.25 printstyled("Running: trapezoidal_rule\n"; color = :magenta) diff --git a/test/cosine.jl b/test/cosine.jl index cc6d61c..467979f 100644 --- a/test/cosine.jl +++ b/test/cosine.jl @@ -28,7 +28,7 @@ printstyled("Integrating cosine from 0 to pi/2 and comparing the result to the a printstyled("Running: rombergs_method\n"; color = :magenta) @time @test rombergs_method(x -> cos(x), 3, 0, pi/2) ≈ 1 printstyled("Running: simpsons_rule\n"; color = :magenta) - @time @test simpsons_rule(x -> cos(x), 40, 0, pi/2) ≈ 1 + @time @test simpsons_rule(x -> cos(x), 20, 0, pi/2) ≈ 1 printstyled("Running: simpsons38_rule\n"; color = :magenta) @time @test simpsons38_rule(x -> cos(x), 48, 0, pi/2) ≈ 1 printstyled("Running: trapezoidal_rule\n"; color = :magenta) diff --git a/test/expnx2datan.jl b/test/expnx2datan.jl index a65e1fe..ebe8632 100644 --- a/test/expnx2datan.jl +++ b/test/expnx2datan.jl @@ -34,7 +34,7 @@ printstyled("Integrating e^(-x^2)/(1+x^2) on the infinite domain [-inf, inf], or printstyled("Running: rombergs_method\n"; color = :magenta) @time @test rombergs_method(x -> exp(-x^2)/(x^2+1), 12, -100, 100) ≈ sol_9 printstyled("Running: simpsons_rule\n"; color = :magenta) - @time @test simpsons_rule(x -> exp(-x^2)/(x^2+1), 1240, -100, 100) ≈ sol_9 + @time @test simpsons_rule(x -> exp(-x^2)/(x^2+1), 620, -100, 100) ≈ sol_9 printstyled("Running: simpsons38_rule\n"; color = :magenta) @time @test simpsons38_rule(x -> exp(-x^2)/(x^2+1), 1767, -100, 100) ≈ sol_9 printstyled("Running: trapezoidal_rule\n"; color = :magenta) diff --git a/test/gaussian.jl b/test/gaussian.jl index 445190f..08d3bd8 100644 --- a/test/gaussian.jl +++ b/test/gaussian.jl @@ -34,7 +34,7 @@ printstyled("Integrate exp(-x^2) from minus infinity to positive infinity and co printstyled("Running: rombergs_method\n"; color = :magenta) @time @test rombergs_method(x -> exp(-x^2), 12, -100, 100) ≈ sqrt(pi) printstyled("Running: simpsons_rule\n"; color = :magenta) - @time @test simpsons_rule(x -> exp(-x^2), 536, -100, 100) ≈ sqrt(pi) + @time @test simpsons_rule(x -> exp(-x^2), 268, -100, 100) ≈ sqrt(pi) printstyled("Running: simpsons38_rule\n"; color = :magenta) @time @test simpsons38_rule(x -> exp(-x^2), 780, -100, 100) ≈ sqrt(pi) printstyled("Running: trapezoidal_rule\n"; color = :magenta) diff --git a/test/logarithm.jl b/test/logarithm.jl index e7b92ea..a1fbd34 100644 --- a/test/logarithm.jl +++ b/test/logarithm.jl @@ -28,7 +28,7 @@ printstyled("Integrating 1/x from 1 to e and comparing the result to the analyti printstyled("Running: rombergs_method\n"; color = :magenta) @time @test rombergs_method(x -> 1/x, 5, 1, exp(1)) ≈ 1 printstyled("Running: simpsons_rule\n"; color = :magenta) - @time @test simpsons_rule(x -> x^(-1), 70, 1, exp(1)) ≈ 1 + @time @test simpsons_rule(x -> x^(-1), 34, 1, exp(1)) ≈ 1 printstyled("Running: simpsons38_rule\n"; color = :magenta) @time @test simpsons38_rule(x -> x^(-1), 84, 1, exp(1)) ≈ 1 printstyled("Running: trapezoidal_rule\n"; color = :magenta) diff --git a/test/logoverx.jl b/test/logoverx.jl index cd25ff2..7c0d658 100644 --- a/test/logoverx.jl +++ b/test/logoverx.jl @@ -27,7 +27,7 @@ printstyled("Integrating log(x)/x from 1 to e and comparing the result to the an printstyled("Running: rombergs_method\n"; color = :magenta) @time @test rombergs_method(x -> log(x)/x, 5, 1, exp(1)) ≈ 0.5 printstyled("Running: simpsons_rule\n"; color = :magenta) - @time @test simpsons_rule(x -> log(x)/x, 100, 1, exp(1)) ≈ 0.5 + @time @test simpsons_rule(x -> log(x)/x, 50, 1, exp(1)) ≈ 0.5 printstyled("Running: simpsons38_rule\n"; color = :magenta) @time @test simpsons38_rule(x -> log(x)/x, 114, 1, exp(1)) ≈ 0.5 printstyled("Running: trapezoidal_rule\n"; color = :magenta) diff --git a/test/modbessel0.jl b/test/modbessel0.jl index b897a38..91d387f 100644 --- a/test/modbessel0.jl +++ b/test/modbessel0.jl @@ -36,7 +36,7 @@ printstyled("Approximating the modified bessel function I_1(1) and comparing it printstyled("Running: rombergs_method\n"; color = :magenta) @time @test rombergs_method(dmodified_bessel_1, 5, 0, pi/2) ≈ besseli(x,1) printstyled("Running: simpsons_rule\n"; color = :magenta) - @time @test simpsons_rule(dmodified_bessel_1, 6, 0, pi/2) ≈ besseli(x,1) + @time @test simpsons_rule(dmodified_bessel_1, 3, 0, pi/2) ≈ besseli(x,1) printstyled("Running: simpsons38_rule\n"; color = :magenta) @time @test simpsons38_rule(dmodified_bessel_1, 9, 0, pi/2) ≈ besseli(x,1) printstyled("Running: trapezoidal_rule\n"; color = :magenta) diff --git a/test/simppen.jl b/test/simppen.jl index ecf25a3..b2e04b1 100644 --- a/test/simppen.jl +++ b/test/simppen.jl @@ -30,7 +30,7 @@ printstyled("Integrating 1/sqrt(-19.6 sin(x)) from -pi to 0 and comparing the re printstyled("Running: rombergs_method on [-pi+1e-8, -1e-8]. Only a rough approximation can be realistically achieved with this function, due to the singularities.\n"; color = :magenta) @time @test abs(rombergs_method(x -> (-19.6*sin(x))^(-0.5), 24, -pi+1e-8, -1e-8) - ellipk(1/2)/sqrt(2.45)) < 1e-4 printstyled("Running: simpsons_rule on [-pi+1e-8, -1e-8]. Only a rough approximation can be realistically achieved with this function, due to the singularities.\n"; color = :magenta) - @time @test abs(simpsons_rule(x -> (-19.6*sin(x))^(-0.5), 1e8, -pi+1e-8, -1e-8) - ellipk(1/2)/sqrt(2.45)) < 1e-4 + @time @test abs(simpsons_rule(x -> (-19.6*sin(x))^(-0.5), 5e7, -pi+1e-8, -1e-8) - ellipk(1/2)/sqrt(2.45)) < 1e-4 printstyled("Running: simpsons38_rule on [-pi+1e-8, -1e-8]. Only a rough approximation can be realistically achieved with this function, due to the singularities.\n"; color = :magenta) @time @test abs(simpsons38_rule(x -> (-19.6*sin(x))^(-0.5), 1e8+2, -pi+1e-8, -1e-8) - ellipk(1/2)/sqrt(2.45)) < 1e-4 printstyled("Running: trapezoidal_rule on [-pi+1e-8, -1e-8]. Only a rough approximation can be realistically achieved with this function, due to the singularities.\n"; color = :magenta) diff --git a/test/sinexpox.jl b/test/sinexpox.jl index 931245e..d0cc99e 100644 --- a/test/sinexpox.jl +++ b/test/sinexpox.jl @@ -40,7 +40,7 @@ printstyled("Integrating sin(x^2)e^(-x)/x from 0 to infinity, with the approxima printstyled("Running: rombergs_method\n"; color = :magenta) @time @test rombergs_method(sinexpox, 12, 0, 100) ≈ sol_11 printstyled("Running: simpsons_rule\n"; color = :magenta) - @time @test simpsons_rule(sinexpox, 4202, 0, 100) ≈ sol_11 + @time @test simpsons_rule(sinexpox, 2101, 0, 100) ≈ sol_11 printstyled("Running: simpsons38_rule\n"; color = :magenta) @time @test simpsons38_rule(sinexpox, 5148, 0, 100) ≈ sol_11 printstyled("Running: trapezoidal_rule\n"; color = :magenta) diff --git a/test/sinxx.jl b/test/sinxx.jl index 1b0bc1b..b0ed751 100644 --- a/test/sinxx.jl +++ b/test/sinxx.jl @@ -36,7 +36,7 @@ printstyled("Integrating sin(x)/x from 0 to 100 and comparing it to the exact re printstyled("Running: rombergs_method\n"; color = :magenta) @time @test rombergs_method(sinxx, 9, 0, 100) ≈ sol_8 printstyled("Running: simpsons_rule\n"; color = :magenta) - @time @test simpsons_rule(sinxx, 678, 0, 100) ≈ sol_8 + @time @test simpsons_rule(sinxx, 339, 0, 100) ≈ sol_8 printstyled("Running: simpsons38_rule\n"; color = :magenta) @time @test simpsons38_rule(sinxx, 831, 0, 100) ≈ sol_8 printstyled("Running: trapezoidal_rule\n"; color = :magenta) diff --git a/test/test_7.jl b/test/test_7.jl index e135088..611122d 100644 --- a/test/test_7.jl +++ b/test/test_7.jl @@ -35,7 +35,7 @@ printstyled("Integrating (x^3+1)/(x^4 (x+1)(x^2+1)) from 1 to e and comparing th printstyled("Running: rombergs_method\n"; color = :magenta) @time @test rombergs_method(partfrac, 6, 1, exp(1)) ≈ sol_7 printstyled("Running: simpsons_rule\n"; color = :magenta) - @time @test simpsons_rule(partfrac, 200, 1, exp(1)) ≈ sol_7 + @time @test simpsons_rule(partfrac, 100, 1, exp(1)) ≈ sol_7 printstyled("Running: simpsons38_rule\n"; color = :magenta) @time @test simpsons38_rule(partfrac, 234, 1, exp(1)) ≈ sol_7 printstyled("Running: trapezoidal_rule\n"; color = :magenta)