%CH19        Program for Chapter 19
%
% Up and out call option
% Evaluates Black-Scholes formula and also uses Monte Carlo

randn('state',100)

%%%%%%%%% Problem and method parameters %%%%%%%%%%%
S = 5; E = 6; sigma = 0.25; r = 0.05; T = 1; B = 9; 
Dt = 1e-3; N = T/Dt; M = 1e4;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%%%%%%%%%%%%%%%%%%%%%% Black-Scholes value %%%%%%%%%%%%%%%%%%%%%
tau = T;
power1 = -1 + (2*r)/(sigma^2);
power2 = 1 + (2*r)/(sigma^2);
d1 = (log(S/E) + (r + 0.5*sigma^2)*(tau))/(sigma*sqrt(tau));
d2 = d1 - sigma*sqrt(tau);
e1 = (log(S/B) + (r + 0.5*sigma^2)*(tau))/(sigma*sqrt(tau));
e2 = (log(S/B) + (r - 0.5*sigma^2)*(tau))/(sigma*sqrt(tau));
f1 = (log(S/B) - (r - 0.5*sigma^2)*(tau))/(sigma*sqrt(tau));
f2 = (log(S/B) - (r + 0.5*sigma^2)*(tau))/(sigma*sqrt(tau));
g1 = (log(S*E/(B^2)) - (r - 0.5*sigma^2)*(tau))/(sigma*sqrt(tau));
g2 = (log(S*E/(B^2)) - (r + 0.5*sigma^2)*(tau))/(sigma*sqrt(tau));
Nd1 = 0.5*(1+erf(d1/sqrt(2))); Nd2 = 0.5*(1+erf(d2/sqrt(2)));
Ne1 = 0.5*(1+erf(e1/sqrt(2))); Ne2 = 0.5*(1+erf(e2/sqrt(2)));
Nf1 = 0.5*(1+erf(f1/sqrt(2))); Nf2 = 0.5*(1+erf(f2/sqrt(2)));
Ng1 = 0.5*(1+erf(g1/sqrt(2))); Ng2 = 0.5*(1+erf(g2/sqrt(2)));
a = (B/S)^power1; b = (B/S)^power2;
bsval = S*(Nd1-Ne1-b*(Nf2-Ng2)) - E*exp(-r*tau)*(Nd2-Ne2-a*(Nf1-Ng1))
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

V = zeros(M,1);
for i = 1:M
    Svals = S*cumprod(exp((r-0.5*sigma^2)*Dt+sigma*sqrt(Dt)*randn(N,1)));
    Smax = max(Svals);
    if Smax < B 
       V(i) = exp(-r*T)*max(Svals(end)-E,0);
    end
end
aM = mean(V); bM = std(V);
conf = [aM - 1.96*bM/sqrt(M), aM + 1.96*bM/sqrt(M)]
