  The implementation of discontinuous functions occurs in many of today's
state-of-the-art partial differential equation solvers. However, in finite
element methods, this poses an inherent difficulty: efficient quadrature rules
available when integrating functions whose discontinuity falls in the element's
interior are for low order degree polynomials, not easily extended to higher
order degree polynomials, and cover a restricted set of geometries. Many
approaches to this issue have been developed in recent years. Among them one of
the most elegant and versatile is the equivalent polynomial technique. This
method replaces the discontinuous function with a polynomial, allowing
integration to occur over the entire domain rather than integrating over
complex subdomains. Although eliminating the issues involved with discontinuous
function integration, the equivalent polynomial tactic introduces its problems.
The exact subdomain integration requires a machinery that quickly grows in
complexity when increasing the polynomial degree and the geometry dimension,
restricting its applicability to lower order degree finite element families.
The current work eliminates this issue. We provide algebraic expressions to
exactly evaluate the subdomain integral of any degree polynomial on parent
finite element shapes cut by a planar interface. These formulas also apply to
the exact evaluation of the embedded interface integral. We provide recursive
algorithms that avoid overflow in computer arithmetic for standard finite
element geometries: triangle, square, cube, tetrahedron, and prism, along with
a hypercube of arbitrary dimensions.
