% drives the various finite difference methods from Chapter 8
% 
% 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.
% 
%-----

close all; drawnow; 
clear all; 
clear functions; 

addpath( '../Chapter3' );     % <- in analytic solution to the black-scholes solution 

% initialize some parameters of our method:
% -- these are choosen to duplicate some of Figure 8.12 (the european put) from Wilmott's 
%    book epage 174
% 
sigma  = 0.4; 
r      = 0.1;
E      = 10.0; 
k      = r/(0.5*sigma^2); 
T      = 3/12;                % <- time to expiration (in years)
t      = 0;
  
% the change of variables from financial variables to dimensionless variables is given by 
% 
% x   = log(S/E)
% tau = 0.5*sigma^2 * (T-t) 
% v   = P/E
% u   = v * exp( -alpha * x - beta * tau )
% 

% initialize the "x=log(S/E)" grid ... we will convert to financial variables after solving 
%
SRight = 16.0;                % <- +\infty 
SLeft  = 1e-9;                % <- -\infty 
xLeft  = log(SLeft/E); 
xRight = log(SRight/E); 

Nx = 1000;
dx = (xRight-xLeft)/Nx; 

a  = 0.5; 
a  = 1.0; 
a  = 5.0; 

dt = a*dx^2;
% we integrate the diffusion equation from tau=0 (t=T expiration) to tau=0.5*sigma^2 (t=0 now)
tau_Max = (0.5*sigma^2)*T; 
M       = ceil(tau_Max/dt);
fprintf('running alpha = %10.6f\n',a); 
%[u,xgrid] = implicit_fd_LU(@tran_payoff_put, @u_m_inf_put, @u_p_inf_put, r, sigma, xLeft, xRight, Nx, tau_Max, M );
[u,xgrid] = implicit_fd_SOR(@tran_payoff_put, @u_m_inf_put, @u_p_inf_put, r, sigma, xLeft, xRight, Nx, tau_Max, M );
%[u,xgrid] = implicit_fd_SOR(@tran_payoff_call, @u_m_inf_call, @u_p_inf_call, r, sigma, xLeft, xRight, Nx, tau_Max, M );
%[u,xgrid] = implicit_fd_LU(@tran_payoff_call, @u_m_inf_call, @u_p_inf_call, r, sigma, xLeft, xRight, Nx, tau_Max, M );

% solve some cash-or-nothing problems ... 
% 
if( 0 ) 
  B = 0.5*E; 
  global b; b = B/E; 
  [u,xgrid] = implicit_fd_LU(@tran_payoff_CON_put, @u_m_inf_CON_put, @u_p_inf_CON_put, r, sigma, xLeft, xRight, Nx, tau_Max, M );
  %[u,xgrid] = implicit_fd_LU(@tran_payoff_CON_call, @u_m_inf_CON_call, @u_p_inf_CON_call, r, sigma, xLeft, xRight, Nx, tau_Max, M );
end

% the change of variables from dimensionless variables to financial variables is given by 
S   = E*exp( xgrid ); 
t   = 0;                       % <- the financial time when we want to evaluate our options value 
tau = 0.5*(sigma^2)*(T-t);     % <- translates into this scaled time 

Spow = (S.^(0.5*(1-k))); 
Smat = repmat( Spow(:).', [M+1, 1] ); 
V    = (E^(0.5*(1+k))) * Smat * exp( -(1/4)*((k+1)^2)*tau ).*u; 

figure; ns=plot( S, V(end,:), 'x' ); grid on; hold on; xlabel( 'S' ); ylabel('C'); 

[C,P] = blackscholes(E,T,S,r,sigma); 
if( V(end,end)<1.0 )   % <- a heuristic as to which option we are pricing
  as=plot( S, P, '-og' ); 
  title( ['The European Put Example from Figure 8.12.  \alpha=',num2str(a)] ); 
else
  as=plot( S, C, '-og' ); 
  title( ['The European Call Example.  \alpha=',num2str(a)] ); 
end
legend( [ ns,as ], {'numerical solution', 'analytic solution'}, 'location', 'north' ); 
fn=['fig_8_12_a_',num2str(a),'.eps']; 
%saveas( gcf, ['../../WriteUp/Graphics/Chapter8',fn], 'epsc' ); 

% evaluate for the S values's found in the table: 
Sgrid = [ 0.00, 2, 4, 6, 8, 10:2:16 ]; 
Vgrid = interp1( S, V(end,:), Sgrid, 'linear', 'extrap' ); 
%Vgrid = interp1( S, V(end,:), Sgrid, 'nearest', 'extrap' );
Vgrid