Reputation: 21
I want to compute the volume integral of a function f(x,y,z) over a cylinder whose base has a stadium-like shape (it is a rectangle with semi-circles at the ends).
Although I can evaluate f(x,y,z) for every point (x,y,z) in space and that the shape itself is relatively simple, I couldn't find a way to write down analitically the limits of integration for x, y and z (or for any other common coordinates such as spherical or cylindrical) given that the origin is outside the volume (and I can't just move the origin for a more suitable position because of other reasons).
Therefore, I really think the only way is to do it numerically. It's the first time I face this isssue, but I would like to do it in Python. I searched for Python modules and the standard one seems to be scipy.integrate (I guess Numpy also has similar functions), but the functions require the limits of integration, which I can't specify. Is there a way to overcome this problem using scipy.integrate?
Upvotes: 1
Views: 1511
Reputation: 58721
A robust way for integrating function over arbitrary shapes is to first create a mesh of the that domain, for example made of tetrahedra, then integrate the function over all tetrahedra and sum the results.
There are various Python packages that can help you with that (some of them by me).
For meshing, see https://stackoverflow.com/a/37797164/353337.
For integration over tetrehedra, check quadpy.
Here an example:
import pygmsh
import numpy as np
import quadpy
with pygmsh.occ.Geometry() as geom:
union = geom.boolean_union([
geom.add_rectangle((0.0, 0.0, 0.0), 1.0, 1.0),
geom.add_disk((0.0, 0.5), 0.5),
geom.add_disk((1.0, 0.5), 0.5),
])
geom.extrude(union, [0, 0, 1])
mesh = geom.generate_mesh()
tetra = mesh.get_cells_type("tetra")
t = np.moveaxis(mesh.points[tetra], 0, 1)
scheme = quadpy.t3.get_good_scheme(5)
val = scheme.integrate(
lambda x: np.exp(x[0]),
t
)
print(np.sum(val))
3.3393111615779754
Upvotes: 2