A quadratic approximation to American prices due to Barone-Adesi and Whaley.

We now discuss an approximation to the option price of an American option on a
commodity having a continuous payout, described in Barone-Adesi and Whaley (1987)
(BAW). Their solution technique finds an approximate solution to the
differential equation

To do their approximation, BAW decomposes the American price into the European
price and the early exercise premium

Here is the early exercise premium. This will also have to satisfy the partial differential equation. To come up with an approximation BAW transformed equation (7.1) into one where the terms involving are neglible, removed these, and ended up with a standard linear homeogenous second order equation, which has a known solution.

The functional form of the approximation is

where

and solves

In implementing this, the only problem is finding the critical value .

This is the classical problem of finding a root of the equation

This is solved using Newton's algorithm for finding the root.

We start by finding a first ``seed'' value .
The next estimate of is found

At each step we need to evaluate and its derivative .

Call option

where is the Black Scholes value for commodities.

Similarly, for a put option the approximation is

and we solve for again by Newton's procedure, where now

// file approx_am_call.cc // author: Bernt Arne Oedegaard // implements the quadratic approximation of Barone-Adesi and Whaley // described in Journal of Finance, june 87. #include <cmath> #include <algorithm> #include "normdist.h" // normal distribution #include "fin_algoritms.h" // define other option pricing formulas const double ACCURACY=1.0e-6; double option_price_american_call_approximated_baw( double S, double X, double r, double b, double sigma, double time) { double sigma_sqr = sigma*sigma; double time_sqrt = sqrt(time); double nn = 2.0*b/sigma_sqr; double m = 2.0*r/sigma_sqr; double K = 1.0-exp(-r*time); double q2 = (-(nn-1)+sqrt(pow((nn-1),2.0)+(4*m/K)))*0.5; // seed value from paper double q2_inf = 0.5 * ( (-nn-1.0) + sqrt(pow((nn-1),2.0)+4.0*m)); double S_star_inf = X / (1.0 - 1.0/q2_inf); double h2 = -(b*time+2.0*sigma*time_sqrt)*(X/(S_star_inf-X)); double S_seed = X + (S_star_inf-X)*(1.0-exp(h2)); int no_iterations=0; // iterate on S to find S_star, using Newton steps double Si=S_seed; double g=1; double gprime=1.0; while ((fabs(g) > ACCURACY) && (fabs(gprime)>ACCURACY) // to avoid exploding Newton's && ( no_iterations++<500) && (Si>0.0)) { double c = option_price_european_call_payout(Si,X,r,b,sigma,time); double d1 = (log(Si/X)+(b+0.5*sigma_sqr)*time)/(sigma*time_sqrt); g=(1.0-1.0/q2)*Si-X-c+(1.0/q2)*Si*exp((b-r)*time)*N(d1); gprime=( 1.0-1.0/q2)*(1.0-exp((b-r)*time)*N(d1)) +(1.0/q2)*exp((b-r)*time)*n(d1)*(1.0/(sigma*time_sqrt)); Si=Si-(g/gprime); }; double S_star = 0; if (fabs(g)>ACCURACY) { S_star = S_seed; } // did not converge else { S_star = Si; }; double C=0; double c = option_price_european_call_payout(S,X,r,b,sigma,time); if (S>=S_star) { C=S-X; } else { double d1 = (log(S_star/X)+(b+0.5*sigma_sqr)*time)/(sigma*time_sqrt); double A2 = (1.0-exp((b-r)*time)*N(d1))* (S_star/q2); C=c+A2*pow((S/S_star),q2); }; return max(C,c); // know value will never be less than BS value }

For further details of the Barone Adesi quadratic approximation: Main source is the original paper: Barone-Adesi and Whaley (1987). The approximation is also discussed in section 14A of Hull (1993).