## 5 Finite elements with embedded gradient discontinuity  
An alternative approach to model localised discontinuities within the FEM consists in embedding the discontinuity directly into the elements by adding degrees of freedom at elemental level instead of nodal level (as in the XFEM). The resulting methods appear in the literature with different names and they can be considered somehow as precursors of the XFEM (e.g. [59, 60]). Further references, mainly focusing on fracture mechanics, can be retrieved in the comparative study by [PERSON] [61]. The basic idea is to use standard FE interpolation as from Eq. (20) in all the elements but those intersected by \(\Gamma\). In these latter, the standard interpolation is enriched by suitable functions \(\hat{\psi}_{i}\) ( \(i=1,\ldots,n_{\
u}\) ), which vanish at the nodes and reproduce the desired discontinuity. The resulting interpolation can be written as  
\[T\left(\mathbf{x},t\right)=\sum_{i=1}^{n}N_{i}\left(\mathbf{x}\right)T_{i} \left(t\right)+\sum_{i=1}^{n_{\
u}}\hat{\psi}_{i}\left(\mathbf{x},t\right) \hat{a}_{i}\left(t\right), \tag{45}\]where \(\hat{a}_{i}\) are element-level parameters. The enrichment functions \(\hat{\psi_{i}}\) depend on the time through the interface-position, which, in turn, is described by the time-evolving level set function. Eq. (45) is the embedded-discontinuity version of XFEM's Eq. (31) and can be rewritten in matrix form as  
\[T\left(\mathbf{x},t\right)=\hat{\mathbf{N}}\left(\mathbf{x},t\right)\hat{ \mathbf{T}}\left(t\right), \tag{46}\]  
where the matrix \(\hat{\mathbf{N}}\left(\mathbf{x},t\right)=\left[N_{1}\left(\mathbf{x}\right),\ldots,N_{n}\left(\mathbf{x}\right)\!,\hat{\psi_{1}}\left(\mathbf{x},t\right) \!\ldots\!,\hat{\psi_{n}}\left(\mathbf{x},t\right)\right]\) collects the standard and discontinuity-embedding shape functions, and the vector \(\hat{\mathbf{T}}\left(t\right)=\left[T_{1}\left(t\right),\ldots,T_{n}\left(t \right),\hat{a}_{1}\left(t\right),\ldots,\hat{a}_{n_{\
u}}\left(t\right) \right]^{T}\) the standard nodal degrees of freedom and the additional elemental ones. The numerical model is developed by following the same steps as in Sections 4.2 and 4.3, so that the same formal expressions obtained there for the XFEM hold for the embedded-discontinuity FEM, once Eq. (31) is replaced by Eq. (45). In order to explain the procedure, a two-node element with an embedded gradient discontinuity for the one-dimensional Stefan problem (see also Section 6) is developed here. In this simple case the standard FE shape functions are linear (Figure 4a) and the enrichment of the element containing the phase interface can be obtained by a piecewise-linear hat-shaped function \(\hat{\psi}\left(x,t\right)\) (Figure 4c), whose contribution vanishes at the nodes, assumes a unit value at the discontinuity, and it is controlled by the elemental degree of freedom \(\hat{a}\). It can be observed that the definition of \(\hat{\psi}\) given above does not cover the case in which the phase interface coincides with a node. When the interface reaches the local \(i\)-th node (\(i=1,2\)), the function \(\hat{\psi}\) coincides with \(N_{i}\), in the element that last contained the interface. However, \(\hat{\psi}\) is still linearly independent from \(N_{i}\) at the global level, as visualised in Figure 5, and no particular indeterminacy arises, unless the interface reaches the domain boundary. In the last, quite particular case, the loss of linear independency may also compromise the convergence of the gradient discontinuity at the phase boundary, as observed by [PERSON] and [PERSON] [15], and, in any case, the conditioning of the system matrixes deteriorates as the phase interface gets closer to a boundary node. Figure 4: Shape functions within one element in local coordinate \(\xi\): **(a)** standard, **(b)** enriched in XFEM, and **(c)** enriched in embedded-discontinuity FEM. Figure 5: Standard FE shape functions and embedded-discontinuity function.