# A Compact FEM Implementation for Parabolic Integro-Differential Equations in 2D
## 7 Overview of an Implementation in Matlab
* **Code for the backward-Euler scheme with the left rectangular rule, corresponding to Equation (9) and employing (7):**

```
1:procedureBELRFEMPIDE(\(c4n,n4e,dir,N,dt\))
2:\(nE\leftarrow\mathrm{NROW}(n4e)\)
3:\(nC\leftarrow\mathrm{NROW}(c4n)\)
4:\(U\leftarrow\mathrm{ZEROS}(nC,N+1)\)
5:\(FN\leftarrow\mathrm{UNIQUE}([12\ \ldots nC])\backslash\mathrm{UNIQUE}(dir)\)
6:\(c1\gets c4n(n4e(:,1),:)\)
7:\(d21\gets c4n(n4e(:,2),:)-c1\)
8:\(d31\gets c4n(n4e(:,3),:)-c1\)
9:\(area4\gets 2*(d21(:,1).*d31(:,2)-d21(:,2). *d31(:,1))\)
10:\(A\leftarrow\textsc{StiffAssemb}(n4e,nC,area4,d21,d31,nE)\)
11:\(B\leftarrow\textsc{MassAssemb}(nC,n4e^{T},\frac{1}{4}*area4)\)
12:\(U(:,1)\leftarrow\textsc{U0}(c4n)\)
13:for\(j\gets 2\) to \((N+1)\), do
14:\(U(:,j)\leftarrow\textsc{solveBELRR}(c4n,n4e,dir,dt,FN,A,B,n,U,c1,d21,d31, area4)\)
15:endfor
16:return\(U\)
17:endprocedure
```

**Auxiliary function**

```
18:proceduresolveBELRR(\(c4n,n4e,dir,dt,FN,A,B,n,U,c1,d21,d31,area4\))
2:\(r\leftarrow\mathrm{NROW}(c4n)\)
3:\(fsT\leftarrow\textsc{f}(c1+\frac{d21}{3}+d31,(n-1)*dt)\)
4:\(L\leftarrow\textsc{repmat}(dt.*area4.*fsT. *\frac{1}{12},3,1)\)
5:\(b\leftarrow\textsc{accumarray}(n4e(:),L,[r\ 1])\)
6:\(Int\gets dt^{2}*A*U(:,1)*\textsc{integral}((n-1)*dt,0)\)
7:if\(n>2\), then
8:for\(k\gets 2\) to \((n-1)\), do
9:\(Int\gets Int+dt^{2}*A*U(:,k)*\textsc{integral}((n-1)*dt,(k-1)*dt)\)
10:endfor
11:endif
12:\(b\gets b+B*U(:,n-1)+Int\)
13:\(u\leftarrow\mathrm{SPARSE}(r,1)\)
14:\(u_{d}\leftarrow\mathrm{UNIQUE}(dir)\)
15:\(u(u_{d})\leftarrow\mathrm{UD}(c4n(u_{d}. :),(n-1)*dt)\)
16:\(b\gets b-(dt*A+B)*u\)
17:\(M_{1}\gets dt*A(FN,FN)+B(FN,FN)\)
18:\(u(FN)\leftarrow\mathrm{INV}(M_{1})*b(FN)\)
19:\(solBELR\gets u\)
20:return\(solBELR\)
21:endprocedure
```

**Algorithm A4** Vectorized version of backward-Euler scheme with the left rectangular rule. * **Code for the backward-Euler scheme with the right rectangular rule, corresponding to Equation (9) and employing (8):**

```
1:procedureBERRFEMPIDE(\(c4n,n4e,dir,N,dt\))
2:\(nE\leftarrow\)NROW(\(n4e\))
3:\(nC\leftarrow\)NROW(\(c4n\))
4:\(U\leftarrow\)ZEROS(\(nC,N+1\))
5:\(FN\leftarrow\)UNIQUE(\([12\ \ldots\ nC])\backslash\)UNIQUE(\(dir\))
6:\(c1\gets c4n(n4e(:,1),:)\)
7:\(d21\gets c4n(n4e(:,2),:)-c1\)
8:\(d31\gets c4n(n4e(:,3),:)-c1\)
9:\(area4\gets 2*(d21(:,1)*d31(:,2)-d21(:,2). *d31(:,1))\)
10:\(A\leftarrow\)StiffAssemb(\(n4e,nC,area4,d21,d31,nE\))
11:\(B\leftarrow\)MassAssemb(\(nC,n4e^{T},\frac{1}{4}*area4\))
12:\(U(:,1)\leftarrow\)U0(\(c4n\))
13:for\(j\gets 2\) to \((N+1)\), do
14:\(U(:,j)\leftarrow\)solveBERRR(\(c4n,n4e,dir,dt,FN,A,B,n,U,c1,d21,d31,area4\))
15:endfor
16:return\(U\)
17:endprocedure
```

**Auxiliary function**

```
18:proceduresolveBERRR(\(c4n,n4e,dir,dt,FN,A,B,n,U,c1,d21,d31,area4\))
20:\(r\leftarrow\)NROW(\(c4n\))
21:\(fS\leftarrow\)F(\(c1+\frac{d21+d31}{3},(n-1)*dt\))
22:\(L\leftarrow\)repmat(\(dt.*area4.*fS. *\frac{1}{12},3,1\))
23:\(b\leftarrow\)accumaray(\(n4e(:),L,[r\ 1]\))
24:\(Int\gets 0\)
25:if\(n>2\), then
26:\(Int\gets dt^{2}*A*U(:,1)*\textsc{intgral}((n-1)*dt,dt)\)
27:for\(k\gets 2\) to \((n-2)\), do
28:\(Int\gets Int+dt^{2}*A*U(:,k)*\textsc{intgral}((n-1)*dt,k*dt)\)
29:endfor
30:endif
31:\(b\gets b+B*U(:,n-1)+Int\)
32:\(u\leftarrow\)SPARSE(\(r,1\))
33:\(u_{d}\leftarrow\)UNIQUE(\(dir\))
34:\(u(u_{d})\leftarrow\)UD(\(c4n(u_{d},:),(n-1)*dt\))
35:\(b\gets b-(dt*A+B+dt^{2}*A*\textsc{intgral}((n-1)*dt,(n-1)*dt))*u\)
36:\(M_{1}\gets dt*A(FN,FN)+B(FN,FN)+dt^{2}*A(FN,FN)*\textsc{intgral}((n-1)*dt,(n-1)*dt)\)
37:\(u(FN)\leftarrow\)INV(\(M_{1}\))*\(b(FN\))
38:\(solBERR\gets u\)
39:return\(solBERR\)
40:endprocedure
```

**Algorithm A5** Vectorized version of backward-Euler scheme with the left rectangular rule. ```
1:procedureCNFEMPIDE(\(c4n,n4e,dir,N,dt\))
2:\(nE\leftarrow\)NROW(\(n4e\))
3:\(nC\leftarrow\)NROW(\(c4n\))
4:\(U\leftarrow\)ZEROS(\(nC,N+1\))
5:\(FN\leftarrow\)UNIQUE(\([12\;\ldots nC]\))UNIQUE(\(dir\))
6:\(c1\gets c4n(n4e(:,1),. )\)
7:\(d21\gets c4n(n4e(:,2),. )\)\(-c1\)
8:\(d31\gets c4n(n4e(:,3),. )\)\(-c1\)
9:\(area4\gets 2*(d21(:,1),.*d31(:,2)-d21(:,2),.