function [u,xgrid] = implicit_fd_SOR( payoff_fn, u_m_inf, u_p_inf, r, sigma, xLeft, xRight, Nx, tau_Max, M )
% 
% Solves the diffusion equation using an implicit finite-difference with the SOR iterative scheem
% 
% See the Epage 169 in Wilmotts book: The Mathematics of Financial Derivatives
% 
% Written by:
% -- 
% John L. Weatherwax                2008-06-11
% 
% email: wax@alum.mit.edu
% 
% Please send comments and especially bug reports to the
% above email address.
% 
%-----

xgrid = linspace( xLeft, xRight, Nx+1 ); 
dx    = (xRight-xLeft)/Nx; 
dt    = tau_Max/M; 

a        = dt/(dx*dx); 
n_eps    = 1e-8; 
omega    = 1.0; 
domega   = 0.05; 
oldloops = 10000; 
k        = r/(0.5*sigma^2); 

% get the initial conditions: 
tau    = 0.0; 
oldu   = payoff_fn( xgrid, tau, k ); oldu = oldu(:).'; 

u      = zeros(M+1,Nx+1);    % <- space for our results ...
u(1,:) = oldu; 
newu   = zeros(1,Nx+1); 
b      = zeros(1,Nx+1); 
for m=1:M,
  tau = m*dt; 

  b( 2:(end-1) ) = oldu( 2:(end-1) ); 
  
  % update the endpoints: 
  oldu(1)   = u_m_inf( xgrid(1),   tau, k ); 
  oldu(end) = u_p_inf( xgrid(end), tau, k ); 

  [newu,loops] = SOR_solver( oldu, b, a, omega, n_eps );
  if( loops>oldloops ) domega = -domega; end 
  omega    = omega + domega; 
  oldloops = loops; 
  % prepare for the next pass: 
  oldu = newu;   

  % store 
  u(m+1,:) = newu; 
end