# Integration of Polynomials Times Double Step Function in Quadrilateral Domains for XFEM Analysis
## 1 Introduction  
The exact numerical integration of discontinuous functions by means of a standard quadrature rule is a complex topic that has been investigated by various authors over time [1; 2; 3; 4; 5; 6; 7; 8]. This argument has been particularly relevant recently in the context of fracture mechanics, where Galerkin-based methods, such as the extended finite element method (XFEM) and the generalised finite element method (GFEM), have been used to analyse the mechanical behaviour of structures and solids embedding cracks in their continuum [9; 10; 11]. The major advantages of such methods, in opposition to the standard finite element method (FEM), are that they allow a regular discretisation and eliminate the need for remeshing [7]. These methods, in fact, use discontinuous functions (i.e., the Heaviside step function) and highly localised functions to reproduce the effect of the discontinuity in the solution field [12]. However, these functions result in discontinuous terms in the finite element stiffness matrix, leading to a non-negligible computational error when integrating using a classical quadrature rule due to the non-polynomial nature of the integrand [13]. Quadrature of terms including discontinuous and singular functions is usually achieved in this context by subdividing the elements crossed by discontinuities into quadrature subdomains [7; 12; 14], although alternative quadrature methods via adaptive quadrature [15], nonconvex polygons [16] and a regularised Heaviside step function [17; 18] have been proposed. Nonetheless, splitting the integration domain into subdomains spoils the eleganceof XFEM and somehow voids its main purpose of not requiring remeshing as quadrature subcells are to be introduced. In this context, an efficient solution has been presented by [PERSON] in [13; 19], where the concept of equivalent polynomials was introduced, which has been proven its effectiveness in exact numerical integration of polynomials times step functions over various domain shapes [13; 19; 20]. The equivalent polynomials concept has been explored by the authors, and it has been extended to the case of a double discontinuity within a single finite element as a first step towards the integration of an arbitrary number of discontinuities. This situation is not uncommon and can be encountered when scenarios such as crack branching, kinking or junction occur [11; 21; 22; 23; 24; 25; 26; 27; 28]. In this context, an efficient integration algorithm to exactly evaluate polynomials times double step function over quadrangular finite element domains is presented in this paper. The proposed algorithm has been designed to tackle the integration issues that arise in problems involving more than one discontinuity (such as crack branching or multiple fractures) within a single enriched finite element of a XFEM discretisation, providing an exact numerical solution without subdividing the domain. The significance of this work is nonetheless much deeper, aiming to provide an algorithm that could be remarkably useful not only in fracture mechanic problems in the XFEM/GFEM domain but also in computational geometry and as a mathematical tool to solve integrals over non-trivial shape domains with ease. Along with the proposed algorithm, the _double discontinuity equivalent polynomials_ (_DD_EQP_) library is herein presented (Supplementary Materials). The library is a double discontinuity version of the _equivalent polynomials_ library (_EQP_) recently developed by the authors [20] to evaluate integrals of polynomials times step function over various shapes of 2D and 3D domains of integration. _DD_EQP_ is a simple and direct application of the proposed integration algorithm; it can be used right away by users or implemented in more complex computational software or frameworks effortlessly. The structure of the paper is briefly summarised. In Section 2, a literature review about numerical methods in fracture mechanics and problems involving multiple discontinuities is presented in order to give context to the proposed algorithm.