## 4 Review of XFEM for Stefan problem
### Phase-interface temperature condition  
Eq. (38) represents the discrete form of the problem of heat conduction as from Eqs (1-3). The interface temperature condition expressed by Eq. (4) must be also recasted into a discrete form and enforced. This can be achieved by following two different approaches: (a) a general one in which the enforcement methods (see below) are applied directly to Eq. (4), a weak form is derived, and, then, discretised; and (b) a simplified approach in which Eq. (4) is firstly discretised and then the constraint enforcement methods are applied to the discrete form. The simplified approach (b) is followed here. This can be successfully used for the one-dimensional cases presented in the examples of this paper. However, the enforcement of interface conditions in methods derived from the PUM, such as the XFEM, is an intensively-researched field, due to the stability issues arising in higher dimensional cases (e.g. [55, 56]). Eq. (4) is discretised by selecting some points, finite in number, on the interface, namely \(\mathbf{x}_{h}^{*(p)}\), \(h=1,\ldots,n_{n}^{(p)}\), at time step \(p\). By using Eq. (33) one obtains  
\[g_{h}\left(\overline{\mathbf{T}}^{(p)}\right)=T\left(\mathbf{x}_{h}^{*(p)},t ^{(p)}\right)-T_{\mathrm{t}_{\mathrm{t}_{\mathrm{t}_{\mathrm{t}_{\mathrm{t}_{ \mathrm{t}}}}}}}=\overline{\mathbf{N}}\left(\mathbf{x}_{h}^{*(p)},p\Delta t \right)\overline{\mathbf{T}}^{(p)}-T_{\mathrm{t}_{\mathrm{t}_{\mathrm{t}_{ \mathrm{t}}}}}=0\,. \tag{40}\]  
The points where the interface condition is imposed can be individuated as the intersections of the interface with the edges of the finite element mesh. At a given time step, Eq. (40) represents a linear constraint, which can be rewritten in matrix form as  
\[\mathbf{g}^{(p)}=\mathbf{G}^{(p)}\overline{\mathbf{T}}^{(p)}-\mathbf{c}= \mathbf{0}\,, \tag{41}\]  
where the all components of the vector \(\mathbf{c}\) are equal to \(T_{\mathrm{t}_{\mathrm{t}_{\mathrm{t}_{\mathrm{t}}}}}\) and the \(h\)-th row of the matrix \(\mathbf{G}^{(p)}\) is \(\overline{\mathbf{N}}\left(\mathbf{x}_{h}^{*(p)},p\Delta t\right)\), which can be determined from the knowledge of the level set. Two general numerical methods for enforcing the constraints are considered here: (i) the penalty method and (ii) the Lagrange multiplier method. Although these methods can also be applied to non-conservative problems [57], they are presented here in their conservative version. The idea within _penalty method_ is to add to the potential, of which Eq. (38) represent the stationary condition, an artificial term \(\frac{1}{2}b\mathbf{g}^{T}\mathbf{g}\), where \(b\) is a positive constant, in general dimensional, called penalty parameter. If \(b\)  is sufficiently large, the minimum of the modified potential cannot be achieved without \(\left|\mathbf{g}\right|\) being very small, i.e. without (approximatively) satisfying the constraints. The stationary condition of the additional potential provide the penalty term \(b\mathbf{G}^{T}\left(\mathbf{G}\mathbf{\overline{T}}-\mathbf{c}\right)\) to be added to the right-hand side of Eq. (38), so that the temperature update step becomes  
\[\left(\mathbf{\mathbf{\tilde{K}}}^{(p)}+b\mathbf{K}_{\mathrm{pen}}^{(p)} \right)\mathbf{\overline{T}}^{(p)}=\mathbf{\tilde{q}}^{(p)}+b\mathbf{q}_{ \mathrm{pen}}^{(p)}\,, \tag{42}\]  
where  
\[\mathbf{K}_{\mathrm{pen}}^{(p)} = \mathbf{G}^{(p)T}\mathbf{G}^{(p)}, \tag{43}\] \[\mathbf{q}_{\mathrm{pen}}^{(p)} = \mathbf{G}^{(p)T}\mathbf{c}.\]  
The solution generated with the penalty method satisfies the condition over the interface only approximately. The choice of the parameter \(b\) is crucial, since a value 'too small' will lead to unacceptable looseness of the constraints, whereas a value 'too large' will deteriorate the conditioning of the linear problem. Examples of this are provided in Section 6.4, where a suitable dimensionless version of \(b\) is defined in order to perform problem-independent analyses.