Write a single, efficient, general program, in MATLAB, Fortran, C that compares the performance of the three numerical methods:
Apply these methods to the inhomogeneous, linear Elliptic PDE (LEPDE):
0 = F(x,y,u,u_{x},u_{y},u_{xx},u_{yy}) = A_{1}*u_{xx} + A_{2}*u_{yy} + B_{1}(x,y)*u_{x} + B_{2}(x,y)*u_{y} +C(x,y)*u + D(x,y);
Caveat Usor: You are responsible for checking the correctness of the formulas, given here and for any modifications needed to transform them to useful code.
Detailed Method Specifications:
The SOR and CGM methods are all derived from the Euler's explicit method which in the linear algebra formulation (with arrays marked with an M for matrix to avoid confusion with the model) is
AM*XM = BM; XM = {U_{i,j} | interior, row order: for j = 1 : m_{2}, i = 1 : m_{1}};
XM_{k} = XM^{(k)}; XM^{(k)}_{i+(j-1)*m1} = U_{i,j}^{(k)}; AM = LM+UM-DM; BM = {RHS vector, including boundary conditions};
h_{q} = (b_{q} - a_{q})/(m_{q}+1); for q = 1 : 2,
1. Successive Over Relaxation (SOR) Method:
Let XM_{0} = -w*DM^{-1}*BM to start.
For k = 0 : k_{max} While (||RM_{k}||_{2} > tol*||BM||_{2}) Do:
XM_{k+1} = XM_{k} - w*DM^{-1}*RM_{k}; %%SOR
RM_{k} = BM-LM*XM_{k+1}-(UM-DM)*XM_{k} ;
(DM-w*LM)*XM_{k+1} = ((1-w)*DM + w*UM)*XM_{k}-w*BM ;
U_{i,j}^{(k+1)} = c_{0,0,i,j}*U_{i,j}^{(k)} + c_{-1,0,i,j}*U_{i-1,j}^{(k+1)} + c_{0,-1,i,j}*U_{i,j-1}^{(k+1)} + c_{1,0,i,j}*U_{i+1,j}^{(k)} + c_{0,1,i,j}*U_{i,j+1}^{(k)} + c_{0,i,j} ,
2. Alternating Directions Implicit (ADI) Method:
The ADI method is based on operator splitting by dimension for a higher order approximation to replace the Crank-Nicolson Implicit Method CNIM in higher dimension to preserve the tridiagonal advantage in using the Thomas algorithm (See your class notes or Morton and Mayer, the second opinion text). The method specified here is a modification of the PRADI variant of ADI for the variable coefficients and extra drift-source terms here. Use the same initial condition X_{0} as for the SOR method. Let the ADI acceleration parameter be symbolically given by
p = padi = 0.5*Delta_t*Adh2max; where Adh2max = max_i[A_{i}/h_{i}^{2}] ;
R_{i} = (A_{i}/h_{i}^{2})/Adh2max;
ka_{max} = 0.5*k_{max};
for k=0 : ka_{max}
While (||U^{(k+1)}-U^{(k)}||_{2} >
tol*||U^{(k)}||_{2}) Do:
%%% X-ADI:
(1-p*R_{1}*(delta_{1}^{2}
+ 0.5*(B_{1,i,j}/A_{1,i,j})*h_{1}*Delta_{1} + 0.5*(C_{i,j}/A_{1,i,j})*h_{1}^{2}*I))
[U_{i,j}^{(k+0.5)}]
=(1+p*R_{2}*(delta_{2}^{2}
+ 0.5*(B_{2,i,j}/A_{2,i,j})*h_{2}*Delta_{2} + 0.5*(C_{i,j}/A_{2,i,j})*h_{2}^{2}*I))
[U_{i,j}^{(k)}] + 2*p*D_{i,j}/Adh2max ;
%%% Y-ADI:
(1-p*R_{2}*(delta_{2}^{2}
+ 0.5*(B_{2,i,j}/A_{2,i,j})*h_{2}*Delta_{2} + 0.5*(C_{i,j}/A_{2,i,j})*h_{2}^{2}*I))
[U_{i,j}^{(k+1)}]
=(1+p*R_{1}*(delta_{1}^{2}
+ 0.5*(B_{1,i,j}/A_{1,i,j})*h_{1}*Delta_{1} + 0.5*(C_{i,j}/A_{1,i,j})*h_{1}^{2}*I))
[U_{i,j}^{(k)}] + 2*p*D_{i,j}/Adh2max ;
3. Conjugate Gradient Method (CGM):
The CGM method is best modified to account for non-symmetric matrices A and for program efficiency.
For this method, use one of the following professional codes given in MATLAB from the NeTLib software archive (the class pseudo-code only gives the major ideas):
function [x, error, iter, flag] = cgs(A, x, b, M, max_it, tol)
function [x, error, iter, flag] = bicg(A, x, b, M, max_it, tol)
Test Case Specifications:
Your program should be fairly general and not made just for the particular values or functions here, but you can use the following values as a test case as long as your program is not specifically tailored to it. where a_{1} < x < b_{1}, a_{2} < y < b_{2}. On
A_{1} = exp(1.5); A_{2} = exp(1.25);
B_{1}(x,y) = exp(-1.3)*(x-a_{1})/sqrt(1+x^{2}+y^{2});
B_{2}(x,y)=exp(-1.1)*(y-a_{2})/sqrt(1+x^{2}+y^{2});
C(x,y) = -exp(1.5)/sqrt(1+x^{2}+y^{2});
D(x,y) = -exp(-2.5)/sqrt(1+x^{2}+y^{2});
u(x,a_{2})
= 6*exp(0.5)*(x-a_{1})/(b_{1}-a_{1}) ;
u(a_{1},y)
= 3*exp(0.5)*(y-a_{2})/(b_{2}-a_{2}) ;
u(x,b_{2})
= exp(0.5)*(3+(x-a_{1})/(b_{1}-a_{1})) ;
u(b_{1},y)
= 2*exp(0.5)*(3-(y-a_{2})/(b_{2}-a_{2})) ;
Output Specifications:
Efficiency Hints:
MCS 571 HomePage: http://www.math.uic.edu/~hanson/mcs571/
Email Comments or Questions to hanson A T uic edu