$\textbf{Alternating-Direction-Implicit Method}$
For my final project for AMATH 581, I will be implementing the alternating-direction-implicit method (ADI). The method is outlined here.
This method can be used for any elliptic PDE. For this project, I will consider the Laplace Equation
$$u_{xx} + u_{yy} = 0$$
Using 2nd order difference scheme for both $u_{xx}$ and $u_{yy}$, choosing $\Delta x = \Delta y = h$, and $N$ points in each direction, we have the following:
$$u(x - h,y) + u(x + h,y) + u(x, y - h) + u(x, y + h) - 4u(x,y) =0$$
I discretize the grid so that $x = 0, h, 2h,...$ correspond to $x_i$ and $y = 0, h, 2h,...$ correspond to $y_j$, and $u(x_i, y_j) = u_{i,j}$:
$1)$
\begin{equation}
u_{i - 1,j} + u_{i+1,j} + u_{i, j-1} + u_{i, j+1} - 4u_{i,j} = 0
\end{equation}
If all the terms with index $j$ on the left hand side, we have:
$2)$
$$u_{i - 1,j} - 4u_{i,j} + u_{i+1,j} = -u_{i, j-1} -u_{i, j+1}$$
Likewise, rearranging $1)$ so the terms with index $i$ are on the left hand side:
$3)$
$$u_{i,j-1} - 4u_{i,j} + u_{i,j+1} = -u_{i-1, j} -u_{i+1, j}$$
$\textbf{Algorithm:}$
Every iteration $m$ of the ADI method has two stages. In the first stage, we will use equation $2)$ to update the grid, and for the second stage, we will use equation $3)$ to update the grid. The first update using equation $2)$ is:
$4)$
$$u_{i - 1,j}^{m+1/2} - 4u_{i,j}^{m+1/2} + u_{i+1,j}^{m+1/2} = -u_{i, j-1}^{m} -u_{i, j+1}^{m}$$
Here, the superscipts represent the iteration of the solver. For the first stage, represented by the $m+1/2$ superscript, we set up a system of equations using $4)$ for each of the interior points along the line horizontal line $y = y_j$ and solve for $u_{1,j}^{m+1/2}, u_{2,j}^{m+1/2},..., u_{N - 1,j}^{m+1/2}$. We do this for each $y_j$.
We then compute the second stage, represented by the $m+1$ superscript, with equation $3)$:
$5)$
$$u_{i,j-1}^{m+1} - 4u_{i,j}^{m+1} + u_{i,j+1}^{m+1} = -u_{i-1, j}^{m+1/2} -u_{i+1, j}^{m+1/2}$$
Now we use the equation above for each point along the horizontal line of $x = x_i$ and solve for $u_{i,1}^{m+1},u_{i,2}^{m+1},...,u_{i,N-1}^{m+1}$. We do this for each $x_i$.
Every time we complete stage one or two, we solve a relatively small $N \space \text{by} N$ system of equations $A \vec{x} = \vec {b}$, where the matrix $A$ has the following form:
$$A =
\begin{bmatrix}
-4 & 1\\
1 & -4 & 1\\\
& 1& -4 & 1\\
& & & ...
\end{bmatrix}
$$
Although we solve this small $N$ by $N$ $2N$ system $2N$ times per iteration, this method can be advantageous much more advantageous than solving for all points simultaneously, which would involve solving a banded system that is $N^2$ by $N^2$. For example, problems which have a very high number of grid points will be unfeasible to solve directly and instead solving many smaller tridiagonal systems can be much faster. Comparing ADI to other iterative methods, like Jacobi iteration or Gauss Seidel, which iteratively solve for each point in the grid, ADI solves iteratively solves for each line of points. In this way, ADI lies somewhere between a iterative solver and a direct solver.
