How Can I Implement Early Exercise in a PDE Method pricing an American Put Option?

· Pricing and Hedging Tutoriual

An American option is an option which the owner can exercise at any time during its lifetime. That means the option’s value cannot drop below the exercise value, i.e. the option value V(t) of an American put option satisfies

V(t)\geq \max(K-S(t),0), \forall t. (1)

We use the above condition in the PDE solver (How Can I Price an Option with a PDE Method in Matlab?) to price an American put option with S=100K=100T=2r=0.05 and \sigma=0.4.

All posts in this series:

Implementation of the PDE Solver

Introducing the constraint (1) can be done by adding a single line of code within the solver loop:

% solve for V(t_0)
for tau=1:time_steps
  V(:,tau+1) = M\V(:,tau);
  V = max(K-S',V);

Analyzing the convergence of the solver, we will increase “refine” from 1 to 9 and then take a close look at the results. The complete solver with iteration of the refine is:

% init S with a structured mesh
refine = 9;

sigma = 0.4; r = 0.05; K = 100; T = 1;

S = [0 60:20:160 200 1000];

Sn = [];
for iterRefine=1:refine
 % refine mesh
 if iterRefine>1
   Sn(1:2:length(S)*2) = S;
   Sn(2:2:length(S)*2-1) = (S(2:length(S)) -S(1:length(S)-1))/2 + S(1:length(S)-1);
   S = Sn;

 time_steps = round(256 *2^(iterRefine-1) * T);
 dtau = T/time_steps;

 V = max(K-S',0); % V(t=T) = payoff at maturity

 % Abbreviation
 Sm = S(1:end-2);
 Si = S(2:end-1);
 Sp = S(3:end);

 m_iminus = -(sigma^2*Si.^2*dtau)./((Sp-Sm).*(Si-Sm)) + Si.*r*dtau./(Sp-Sm);

 m_i = 1 + (sigma^2*Si.^2*dtau).*(1./( (Sp-Sm).*(Sp-Si) ) +1./( (Sp-Sm).*(Si-Sm) ) ) + r*dtau;
 m_i = [1+r*dtau m_i 1]; % boundary

 m_iplus =-(sigma^2*Si.^2*dtau)./((Sp-Sm).*(Sp-Si)) - Si.*r*dtau./(Sp-Sm);

 % Same but more efficient than: M = diag(m_iminus,-1) + diag(m_i) + diag(m_iplus,1);
 M = spdiags([ [m_iminus 0 0]; m_i; [0 0 m_iplus]]', -1:1, length(m_i), length(m_i));

 for tau=1:time_steps
   V = M\V;
   V = max(K-S',V);

 % save value at S=100 for convergence study
 V_res(iterRefine,1) = V(S==100);

On my laptop, this code runs about 30 sec. Then, the option prices for S=100 for the different refine values are in V_res. An increase of the refinement doubles the nodes in the mesh for S and doubles the number of time steps. So, what happens to the error?

 refine        V_res  Difference       Factor
      1      13.1160           -            -
      2      13.5188      0.4028            -
      3      13.6265      0.1077       3.7399
      4      13.6555      0.0290       3.7185
      5      13.6637      0.0083       3.5074
      6      13.6663      0.0026       3.2000
      7      13.6671      0.0008       3.0796
      8      13.6674      0.0003       2.9366
      9      13.6675      0.0001       2.6696

In the above table, we see that increasing the “refine” by one (doubling the nodes in S mesh and doubling the time steps) reduces the error to about 1/3.

Graph illustrating V_diff

Possible further improvements

Instead of the implicit Euler time-stepping, one can increase the convergence combining implicit and explicit Euler to Crank-Nickolson. Then, using Crank-Nickolson and apply the constraint (1) in an implicit way, the increasing “refine” by one reduces the error to 1/4. That means, highly accurate results are obtained in less time (see e.g. Peter Forsyth: Quadratic convergence of a penalty method for valuing American options).


It is easy to impose the early exercise constraint in a PDE solver. The presented implementation only requires a single extra line of code compared to the European Option PDE solver.

Leave a Reply

Fill in your details below or click an icon to log in: Logo

You are commenting using your account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s

This site uses Akismet to reduce spam. Learn how your comment data is processed.