# Efficient Discretization of the Laplacian: Application to Moving Boundary Problems
## 2 Materials and Methods
### Application to the Numerical Solution of the Stefan Problem
#### 2.5.3 STEP (ii) (a)  
We denote by \(\Gamma^{n}\) the approximation of \(\Gamma_{n\delta}\). As we have a quadrilateral mesh, \(\Gamma^{n}\) becomes a broken line. It is determined by a number of \(Nx\) points such that their horizontal coordinates are fixed and uniformly distributed. According to (13), we can approximate the evolution of \(\Gamma_{t}\) at each point \((x(t),y(t))\in\Gamma_{t}\) using the first-order approximation of (12) as  
\[ \begin{split}&(x((n+1)\delta),y((n+1)\delta))\approx(x(n\delta), y(n\delta))+\delta\cdot\partial_{t}\Gamma_{n\delta}(x(n\delta),y(n\delta))\\ &=(x(n\delta),y(n\delta))+\delta\cdot[[\partial_{v}u(x(n\delta),y(n\delta))]]\,,\end{split}  \tag{21}\]  
which also defines a time stepping. Observe that using only (21) for the given grid points results in general points on \(\Gamma^{n+1}\), which do not necessarily have \(x\)-coordinates fitting to the mesh, Figure 4.  
Before dealing with this problem, we should first compute \([[\partial_{\mathbf{v}}\mathbf{u}^{n}]]\), for which we need the gradients \(\
abla\mathbf{u}_{i}^{n}\). With these, for the jump term in (21), we have  
\[[[\partial_{\mathbf{v}}u(x(n\delta),y(n\delta))]]=(\mathbf{v}\cdot(\
abla\mathbf{u}_{0}^{n} -\
abla\mathbf{u}_{1}^{n}))\mathbf{v}.\]  
Another chief problem is that function values are given in the grid points, where the outward normal of \(\Gamma^{n}\) is not defined. Instead, to move each broken line segment of \(\Gamma^{n}\), we compute the normal vector \(\mathbf{v}\) in the midpoint of these sections. These are given as \(\frac{x_{i}+x_{i+1}}{2},i=1,2,\ldots,Nx\), see Figure 4. We have to approximate the gradients here, which is performed in the following way. We first develop a bilinear interpolation of the numerical solution on the quadrilateral with the vertices \(\mathbf{x}_{i},\mathbf{x}_{i+1},\mathbf{x}_{i+Nx+3},\mathbf{x}_{(i+1)+Nx+3}\) lying on the top of the zero level set as follows. \[\mathbf{u}_{0}^{n}(x,y)\approx a_{0}+a_{1}\mathbf{x}+a_{2}y+a_{3}xy, \tag{22}\]  
where the coefficients \(a_{j}\) are obtained by solving a linear system. Shifting the origin to \(\mathbf{x}_{i}=(x_{j},y_{k})\) and evaluating (22) here, we first observe that \(a_{0}=0\). Similarly, relating point values and the right-hand side of (22) at the points \(\mathbf{x}_{i},\mathbf{x}_{i+1},\mathbf{x}_{i+Nx+3},\mathbf{x}_{(i+1)+Nx+3}\), we obtain the equalities  
\[ \begin{cases}0=a_{1}h_{x}+a_{2}(y_{k+1}-y_{k})+a_{3}h_{x}(y_{k+1}-y_{k})\\ \mathbf{u}_{0}^{n+1}(\mathbf{x}_{j+1+Nx+3})=a_{1}h_{x}+a_{2}(y_{k+1+Nx+3}-y_{k})+a_{3}h _{x}(y_{k+1+Nx+3}-y_{k})\\ \mathbf{u}_{0}^{n}(\mathbf{x}_{j+Nx+3})=a_{2}(y_{k+Nx+3}-y_{k}),\end{cases} \]  
where the number in the bracket of the right-hand side yields the vector component. Similarly, we interpolate the numerical solution on the quadrilaterals lying on the bottom of \(\Gamma^{n}\) to obtain  
\[\mathbf{u}_{1}^{n}(x,y)\approx b_{0}+b_{1}\mathbf{x}+b_{2}y+b_{3}xy. \tag{23}\]  
Corresponding to the previous case, we obtain  
\[ \begin{cases}\mathbf{u}_{1}^{n}(\mathbf{x}_{i-(Nx+3)})=b_{2}(y_{k-(Nx+3)}-y_{k})\\ \mathbf{u}_{1}^{n}(\mathbf{x}_{i+1-(Nx+3)})=b_{1}h_{x}+b_{2}(y_{k+1-(Nx+3)}-y_{k})+b_{ 3}h_{x}(y_{k+1-(Nx+3)}-y_{k})\\ 0=b_{1}h_{x}+b_{2}(y_{k+1}-y_{k}).\end{cases} \]  
Then, by (22) and (23), we obtain  
\[\
abla\mathbf{u}_{0}^{n}(x_{j+1/2},y_{k+1/2})=(a_{1}+a_{3}y_{k+1/2},a_{2}+a_{3}x_ {j+1/2}), \tag{24}\]  
\[\
abla\mathbf{u}_{1}^{n}(x_{j+1/2},y_{k+1/2})=(b_{1}+b_{3}y_{k+1/2},b_{2}+b_{3}x_ {j+1/2}). \tag{25}\]  
Next, we move the segment according to (21). That is, we compute  
\[\mathbf{d}_{i}=\delta\cdot(\mathbf{v}_{i}\cdot(\
abla\mathbf{u}_{0}^{n}-\
abla\mathbf{u}_{1}^ {n}))\mathbf{v}_{i} \tag{26}\]  
using (24), (25) and shift the segment between \(\mathbf{x}_{i}\) and \(\mathbf{x}_{i+1}\) with \(\mathbf{d}_{i}\) in (26). In this way, applying it on two neighboring segments between \(\mathbf{x}_{i}\) and \(\mathbf{x}_{i+1}\) and between \(\mathbf{x}_{i}\) and \(\mathbf{x}_{i+1}\), we have two new vertices of \(\Gamma^{n+1}\) on the fixed vertical lines.