function nma_P2DDIRSOR_S() %script solves 2D poisson on unit square, zero boundary conditions, SOR method % file name: nma_P2DDIRSOR_S.m % % This script solves the 2D poisson on unit square for zero % boundary conditions using SOR method % by Nasser M. Abbasi % %---------- CONFIGURATION PARAMETERS --------------------------- h = 1/40; % mesh spacing tolerance = 10^-7; % tolerance for convergence using relative residual % define the RHS, or f below %myforce = @(X,Y) -exp( -(X-0.25).^2 - (Y-0.6).^2 ); % RHS of PDE myforce = @(X,Y) exp(X) .* sin(Y); %RHS DOPLOTS = false; %set to false if do not want to see plots w = 2*(1-pi*h); %optimal omega %---------- END CONFIGURATION PARAMETERS ------------------------ [X,Y] = meshgrid(0:h:1, 0:h:1); f = myforce(X,Y); % evaluate f(x,y) on grid n = size(f,1); % total number of grid points on one side normf = norm(f,2); % save ||f|| for use later % initialize space (grid) for residual calculation and for solution resid = zeros(n,n); u = resid; unew = u; done = false; % flag set to true in loop below when it converges k = 1; % initialize iteration counter while not(done) for i = 2 : n-1 for j = 2 : n-1 unew(i,j) = (w/4) * ( unew(i-1,j) + u(i+1,j) + unew(i,j-1)+... u(i,j+1) - h^2*f(i,j)) + (1-w)*u(i,j); end end u = unew; for i = 2 : n-1 for j = 2 : n-1 resid(i,j)= f(i,j) - 1/h^2 * ( u(i-1,j) + u(i+1,j) + ... u(i,j-1) + u(i,j+1) - 4*u(i,j) ); end end if ( norm(resid,2) / normf)