Awesome

Trapz.jl

A simple Julia package to perform trapezoidal integration over common Julia arrays. the package is registered on Julia Registry, so it can be added as follows

import Pkg
Pkg.pkg"add Trapz"

Main Usage Example:

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

Example Usage of @trapz macro:

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)

Benchmarks

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)

Benchmarks & example for @trapz macro

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)

Comparison to Numpy

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