A simple Julia package to perform trapezoidal integration over common Julia arrays.
the package is now registered on Julia Registry, so it can be added as follows
import Pkg
Pkg.pkg"add Trapz"
using Trapz
vx=range(0,1,length=100)
vy=range(0,2,length=200)
vz=range(0,3,length=300)
M=[x^2+y^2+z^2 for x=vx,y=vy,z=vz]
I=trapz((vx,vy,vz),M)
print("result: ",I)
result: 28.00030
using Trapz
using Printf
Base.show(io::IO, f::Float64) = @printf(io, "%1.5f", f)
function test(λ)
R=@trapz 0:0.0001:π x (sin(λ*x)/2, cos(λ*x)/2, cos(λ*x)^2/π)
print("λ = ",λ," result of integrals: ",R)
end
test(0.5)
λ = 0.50000 result of integrals: (0.99995, 1.00000, 0.50000)
test(1.0)
λ = 1.00000 result of integrals: (1.00000, 0.00005, 0.49997)
test(2.0)
λ = 2.00000 result of integrals: (0.00000, -0.00005, 0.49997)
using BenchmarkTools
@btime trapz($(vx,vy,vz),$M);
3.131 ms (4 allocations: 157.30 KiB)
@btime trapz($(:,vy, vz),$M);
3.084 ms (3 allocations: 157.20 KiB)
@btime trapz($(:,vy,:),$M);
4.090 ms (2 allocations: 234.45 KiB)
in this example we are calculating 3 multidimensional integrals simultaneously, in other words we are calculating a multidimensional (3D) integral of a vector function
using BenchmarkTools
vx=range(0,1,length=100)
vy=range(0,2,length=200)
vz=range(0,3,length=300)
function integr(vx,vy,vz)
@trapz vx x @trapz vy y @trapz vz z (x*x+y*y+z*z, x*y*z, cos(x*y)+cos(x*z)+cos(y*z))
end
@btime integr($vx,$vy,$vz)
129.633 ms (0 allocations: 0 bytes)
(28.00030, 4.50000, 9.93814)
using PyCall
np=pyimport("numpy")
timenumpy = @belapsed np.trapz(np.trapz(np.trapz($M,$vz),$vy),$vx)
timejulia = @belapsed trapz($(vx,vy,vz),$M)
how_faster=timenumpy/timejulia
print("Trapz.jl is ~ ",how_faster," times faster than numpy's trapz")
Trapz.jl is ~ 7.34493 times faster than numpy's trapz