JuliaGeometry / MeshIntegrals.jl

Numerical integration over Meshes.jl geometry domains
https://juliageometry.github.io/MeshIntegrals.jl/
MIT License
28 stars 4 forks source link
differential-geometry julia numerical-integration quadrature-methods scientific-computing

MeshIntegrals.jl

Docs-stable Docs-dev License: MIT ColPrac

Build Status codecov Coveralls Aqua QA

MeshIntegrals.jl uses differential forms to enable fast and easy numerical integration of arbitrary integrand functions over domains defined via Meshes.jl geometries. This is achieved using:

These solvers have support for integrand functions that produce scalars, vectors, and Unitful.jl Quantity types. While HCubature.jl does not natively support Quantity type integrands, this package provides a compatibility layer to enable this feature.

Usage

Basic

integral(f, geometry)

Performs a numerical integration of some integrand function f over the domain specified by geometry. The integrand function can be anything callable with a method f(::Meshes.Point). A default integration method will be automatically selected according to the geometry: GaussKronrod() for 1D, and HAdaptiveCubature() for all others.

integral(f, geometry, rule)

Performs a numerical integration of some integrand function f over the domain specified by geometry using the specified integration rule, e.g. GaussKronrod(). The integrand function can be anything callable with a method f(::Meshes.Point).

Additionally, several optional keyword arguments are defined in the API to provide additional control over the integration mechanics.

Aliases

lineintegral(f, geometry)
surfaceintegral(f, geometry)
volumeintegral(f, geometry)

Alias functions are provided for convenience. These are simply wrappers for integral that also validate that the provided geometry has the expected number of parametric dimensions. Like with integral, a rule can also optionally be specified as a third argument.

Demo

using Meshes
using MeshIntegrals
using Unitful

# Define a path that approximates a sine-wave on the xy-plane
mypath = BezierCurve(
    [Point(t*u"m", sin(t)*u"m", 0.0u"m") for t in range(-pi, pi, length=361)]
)

# Map f(::Point) -> f(x, y, z) in unitless coordinates
f(p::Meshes.Point) = f(ustrip(to(p))...)

# Integrand function in units of Ohms/meter
f(x, y, z) = (1 / sqrt(1 + cos(x)^2)) * u"Ω/m"

integral(f, mypath)
# -> Approximately 2*Pi Ω