bisect(S0,k,r,T,sigma,D1,t1)

function S_star = bisect(S0,k,r,T,sigma,D1,t1)

%bisectional method of finding the value of S_Star to be used in
%rogewhaley.m file.
%%Inputs
%S0 = Current Stock price 
%D1 = final dividend amount
%k = strike price
%sigma = Volatility
%t1: final dividend date
%T: Time to maturity in years.
%r: risk free rate

%Author: Sivakumar Batthala
%MBA candidate
%Chicago Graduate School of Business
%University of chicago
%Date:02/23/2005
%Please email [e-mail] for any clarifications or errors.

[blscall,p] = bsprice(S0, k, r, T, sigma);
S_high = S0;
temp = blscall - S_high -D1 + k;
ACCURACY = 1e-6;
S_low = 0;

while ((temp>0) && (S_high<=1e10) ) 
      S_high = S_high * 2;
      c = bsprice(S_high, k, r, T-t1, sigma);
      temp = c-S_high-D1+k;
end
      
if (S_high>1e10) %exercise not optimal
   c = bsprice(S_high, k, r, T-t1, sigma);
end
 
S_star = 0.5 * S_high; % // now find S_star that solves c=S_star-D+k
   c = bsprice(S_star, k, r, T-t1, sigma);
   test = c-S_star-D1+k;
   
   while (abs(test) > 0 && (S_high-S_low)>ACCURACY  )
      if (test<0.0)
          S_high = S_star; 
      else
          S_low = S_star;
      end
          S_star = 0.5 * (S_high + S_low);
      c = bsprice(S_star,k,r,T-t1,sigma);
      test = c-S_star-D1+k;
    end
   S_star = S_star;