Compute Integrals
We start with setting the Feel++ environment and loading the Feel++ library.
Set the Feel++ environment with local repository
import feelpp
import sys
app = feelpp.Environment(["myapp"],config=feelpp.localRepository(""))
1. Integral of a function
\(\int_{\Omega} f(x) dx\) is the integral of the function \(f\) over the domain \(\Omega\). The domain \(\Omega\) is a subset of \(\mathbb{R}^n\) and \(f\) is a function defined on \(\Omega\). The domain \(\Omega\) is discretized as a mesh \(\mathcal{T}\) and the function \(f\) may be approximated by a finite element function \(\hat{f}\). The integral is approximated by the sum of the integrals over the elements of the mesh \(\mathcal{T}\). The integral over each element is computed using a quadrature formula.
We first load a mesh
Load a mesh \(\mathcal{T}\)
geo=feelpp.download( "github:{repo:feelpp,path:feelpp/quickstart/laplacian/cases/feelpp2d/feelpp2d.geo}", worldComm=app.worldCommPtr() )[0]
print("geo file: {}".format(geo))
mesh = feelpp.load(feelpp.mesh(dim=2,realdim=2), geo, 0.1)
Then we can define functions \(f\) and compute the integral of \(f\) over the mesh \(\mathcal{T}\).
Compute integrals
from feelpp.integrate import integrate
i1 = integrate(range=feelpp.elements(mesh),expr="sin(x+y):x:y")
i2 = integrate(range=feelpp.boundaryfaces(mesh),expr="x*nx+y*ny+z*nz:x:y:z:nx:ny:nz")
i3 = integrate(range=feelpp.markedelements(mesh, "marker"), expr="1")
print("i1 = {}, i2 = {}, i3 = {}".format(i1,i2,i3))
Results
i1 = [3.24282386], i2 = [149.465], i3 = [0.]