Asian Option Implementation

Jan Vecer, Department of Statistics, Columbia University, vecer@stat.columbia.edu

This method can be used as a benchmark to compare with other methods for solving Asian options as long as you properly refer to this work and to my papers. The precision of this implementation depends on the build in Mathematica function NDSolve, and also on your choice of the cutoff points of the spacial boundary for the corresponding PDE. It is your responsibility to find the optimal cutoff points and WorkingPrecision depending on your current Mathematica version which works best for your parameters. For academic purposes only, usual disclaimer applies.

This is Mathematica implementation of my method for pricing discretely or continuously sampled Asian options as described in my paper Vecer, J. (2002) "Unified Asian Pricing", Risk, Vol. 15, No. 6, 113-116. The numerical precision of the solution is limited by Mathematica's function NDSolve applied for partial differential equations. The precision of NDSolve has been improving for more recent versions. The corresponding PDE has unlimited spacial boundary, so we have to choose the cutoff points, which also determines the precision of the outcome. The cutoff points which depend on the size of the volatility and time work reasonably well.

AsianContinuous[S_, K_, r_, σ_, T_] := 
      Module[{S1 = S, K1 = K, σ1 = σ, r1 = r, T1 = T}, 
        qcontinuous[t2_, T2_, r2_] := (1 - Exp[-r2*t2])/(r2*T2); 
        solution = 
          NDSolve[{(∂_t) u[x, t] == 
                0.5*σ1^2*((1 - Exp[-r1*t])/(r1*T1) - 
                        x)^2*(∂_(x, x))u[x, t], 
              u[x, 0] == If[x > 0, x, 0], 
              u[-1.0*σ1*T1, t] == 0, 
              u[1.5*σ1*T1, t] == 1.5*σ1*T1}, 
            u, {x, -1.0*σ1*T1, 1.5*σ1*T1}, {t, 0, 
              T1}]; (S1*
                u[qcontinuous[T1, T1, r1] - Exp[-r1*T1]*K1/S1, T1] /. 
              solution)[[1]]]

AsianDiscrete[S_, K_, r_, σ_, T_, n_] := 
      Module[{S1 = S, K1 = K, σ1 = σ, r1 = r, T1 = T, n1 = n}, 
        qdiscrete[t2_, T2_, r2_, n2_] := (1/n2)*
            Sum[Exp[-r2*T2*(n2 - k)/n2], {k, 
                Floor[n2*(T2 - t2)/T2] + 1, n2}]; 
        solution = 
          NDSolve[{(∂_t) u[x, t] == 
                0.5*σ1^2*(qdiscrete[t, T1, r1, n1] - 
                        x)^2*(∂_(x, x))u[x, t], 
              u[x, 0] == If[x > 0, x, 0], 
              u[-1.0*σ1*T1, t] == 0, 
              u[1.2*σ1*T1, t] == 1.2*σ1*T1}, 
            u, {x, -1.0*σ1*T1, 1.2*σ1*T1}, {t, 0, 
              T1}]; (S1*
                u[qdiscrete[T1, T1, r1, n1] - Exp[-r1*T1]*K1/S1, T1] /. 
              solution)[[1]]]

Here is the example from the paper which is also used in many other papers as a benchmark. I add timing for comparison, the value depends on the computer. The 7 values should be computed within a fraction of a second. The precision is typically the first 3 decimal digits, it can be further improved by choosing larger cutoff points in the PDE and by increasing WorkingPrecision option in NDSolve function.

Timing[AsianTest={AsianContinuous[1.9,2,.05,.5,1],AsianContinuous[2,2,.05,.5,1],AsianContinuous[2.1,2,.05,.5,1],
AsianContinuous[2,2,.02,.1,1],AsianContinuous[2,2,.18,.3,1],AsianContinuous[2,2,.0125,.25,2],AsianContinuous[2,2,.05,.5,2]}]

Here is another numerical example from my paper for discretely sampled Asian options.

Timing[AsianDicreteTest={AsianDiscrete[95,100,.1,.4,1,10],AsianDiscrete[95,100,.1,.4,1,25],AsianDiscrete[95,100,.1,.4,1,50],AsianDiscrete[95,100,.1,.4,1,125],AsianDiscrete[95,100,.1,.4,1,250],AsianDiscrete[95,100,.1,.4,1,500],AsianDiscrete[95,100,.1,.4,1,1000],AsianContinuous[95,100,.1,.4,1]}]

Timing[AsianDicreteTest={AsianDiscrete[100,100,.1,.4,1,10],AsianDiscrete[100,100,.1,.4,1,25],AsianDiscrete[100,100,.1,.4,1,50],AsianDiscrete[100,100,.1,.4,1,125],AsianDiscrete[100,100,.1,.4,1,250],AsianDiscrete[100,100,.1,.4,1,500],AsianDiscrete[100,100,.1,.4,1,1000],AsianContinuous[100,100,.1,.4,1]}]

Timing[AsianDicreteTest={AsianDiscrete[105,100,.1,.4,1,10],AsianDiscrete[105,100,.1,.4,1,25],AsianDiscrete[105,100,.1,.4,1,50],AsianDiscrete[105,100,.1,.4,1,125],AsianDiscrete[105,100,.1,.4,1,250],AsianDiscrete[105,100,.1,.4,1,500],AsianDiscrete[105,100,.1,.4,1,1000],AsianContinuous[105,100,.1,.4,1]}]