
r     = 0.05;
sigma = 0.2;
T     = 2.5;
S     = 4;
K     = 3;


eps = 1e-20;

delta = european_call(r,sigma,T,S,K,'delta');
V     = european_call(r,sigma,T,S+i*eps,K,'value');

disp(sprintf('delta (analytic) = %g',delta))
disp(sprintf('delta (CVT)      = %g',imag(V)/eps))

gamma = european_call(r,sigma,T,S,K,'gamma');
delta = european_call(r,sigma,T,S+i*eps,K,'delta');

disp(sprintf('gamma (analytic) = %g',gamma))
disp(sprintf('gamma (CVT)      = %g',imag(delta)/eps))

vega  = european_call(r,sigma,T,S,K,'vega');
V     = european_call(r,sigma+i*eps,T,S,K,'value');

disp(sprintf('vega (analytic)  = %g',vega))
disp(sprintf('vega (CVT)       = %g',imag(V)/eps))


%
% plot price versus S(0)
%

S = linspace(0,5);
K = 3;

V = european_call(r,sigma,T,S,K,'value');

figure(1)
plot(S,V)
xlabel('S'); ylabel('V')

%
% plot price versus strike
%

S = 2.5;
K = linspace(1,5);

V = european_call(r,sigma,T,S,K,'value');

figure(2)
plot(K,V)
xlabel('K'); ylabel('V')

