Example 2: A diffusion model for Call Options

In Beckers (1980) there is a model for a Stochastic Differential Equation of option pricing. The idea is a "constant elasticity of variance diffusion (or CEV) class"

 dS = \mu S dt + \sigma S^{\alpha/2}dW,\; 0 \le \alpha \lt 2.
The Black-Scholes model is the limiting case \alpha \to 2. A numerical solution of this diffusion model yields the price of a call option. Various values of the strike price K , time values, \sigma and power coefficient \alpha are used to evaluate the option price at values of the underlying price. The sets of parameters in the computation are:

  1. power \alpha = {2.0, 1.0, 0.0}
  2. strike price  K = {15.0, 20.0, 25.0}
  3. volatility \sigma = {0.2, 0.3, 0.4}
  4. times until expiration = {1/12, 4 /12, 7 /12}
  5. underlying prices = {19.0, 20.0, 21.0}
  6. interest rate r = 0.05
  7. x_{\min} = 0, x_{\max} = 60
  8. nx = 121, \, n = 3 \times nx = 363
With this model the Feynman-Kac differential equation is defined by identifying: The payoff function is the "vanilla option", p(x) = \max(x-K, 0).

import com.imsl.math.*;
import java.util.*;

public class FeynmanKacEx2 {

    public static void main(String args[]) throws Exception {
        // Compute Constant Elasticity of Variance Model for Vanilla Call
        // The set of strike prices
        double[] strikePrices = {15.0, 20.0, 25.0};
        // The set of sigma values
        double[] sigma = {0.2, 0.3, 0.4};
        // The set of model diffusion powers
        double[] alpha = {2.0, 1.0, 0.0};
        // Time values for the options
        int nt = 3;
        double[] tGrid = {1.0 / 12.0, 4.0 / 12.0, 7.0 / 12.0};
        // Values of the underlying where evaluations are made
        double[] evaluationPoints = {19.0, 20.0, 21.0};
        // Value of the interest rate and continuous dividend
        double r = 0.05, dividend = 0.0;
        // Values of the min and max underlying values modeled
        double xMin = 0.0, xMax = 60.0;
        // Define parameters for the integration step. */
        int nxGrid = 121;
        int ntGrid = 3;
        int nv = 3;
        int nint = nxGrid - 1, n = 3 * nxGrid;
        double[] xGrid = new double[nxGrid];
        double dx;
        // Number of left/right boundary conditions
        int nlbcd = 3, nrbcd = 3;
        // Time dependency
        boolean[] timeDependence = new boolean[7];
        double[] userData = new double[6];
        int j;
        // Define equally-spaced grid of points for the underlying price
        dx = (xMax - xMin) / ((double) nint);
        xGrid[0] = xMin;
        xGrid[nxGrid - 1] = xMax;
        for (int i = 2; i <= nxGrid - 1; i++) {
            xGrid[i - 1] = xGrid[i - 2] + dx;
        }

        System.out.printf("%2sConstant Elasticity of Variance Model"
                + " for Vanilla Call%n", " ");
        System.out.printf(Locale.ENGLISH,
                "%7sInterest Rate:%7.3f   Continuous Dividend:%7.3f%n",
                " ", r, dividend);
        System.out.printf(Locale.ENGLISH,
                "%7sMinimum and Maximum Prices of Underlying:%7.2f%7.2f%n",
                " ", xMin, xMax);
        System.out.printf("%7sNumber of equally spaced spline knots:%4d%n",
                " ", nxGrid - 1);
        System.out.printf("%7sNumber of unknowns:%4d%n%n", " ", n);
        System.out.printf(Locale.ENGLISH,
                "%7sTime in Years Prior to Expiration: %7.4f%7.4f%7.4f%n",
                " ", tGrid[0], tGrid[1], tGrid[2]);
        System.out.printf(Locale.ENGLISH,
                "%7sOption valued at Underlying Prices:%7.2f%7.2f%7.2f%n%n",
                " ", evaluationPoints[0], evaluationPoints[1],
                evaluationPoints[2]);

        for (int i1 = 1; i1 <= 3; i1++) /* Loop over power */ {
            for (int i2 = 1; i2 <= 3; i2++) /* Loop over volatility */ {
                for (int i3 = 1; i3 <= 3; i3++) /* Loop over strike price */ {
                    // Pass data through into evaluation methods.
                    userData[0] = strikePrices[i3 - 1];
                    userData[1] = xMax;
                    userData[2] = sigma[i2 - 1];
                    userData[3] = alpha[i1 - 1];
                    userData[4] = r;
                    userData[5] = dividend;

                    FeynmanKac callOption = new FeynmanKac(
                            new MyCoefficients(userData));

                    // Right boundary condition is time-dependent
                    timeDependence[5] = true;
                    callOption.setTimeDependence(timeDependence);
                    callOption.computeCoefficients(nlbcd, nrbcd,
                            new MyBoundaries(userData), xGrid, tGrid);
                    double[][] optionCoefficients
                            = callOption.getSplineCoefficients();
                    double[][] splineValuesOption = new double[ntGrid][];

                    // Evaluate solution at vector evaluationPoints, at each time
                    // value prior to expiration.
                    for (int i = 0; i < ntGrid; i++) {
                        splineValuesOption[i] = callOption.getSplineValue(
                                evaluationPoints, optionCoefficients[i + 1], 0);
                    }

                    System.out.printf(Locale.ENGLISH, "%2sStrike=%5.2f, Sigma="
                            + "%5.2f, Alpha=%5.2f%n", " ", strikePrices[i3 - 1],
                            sigma[i2 - 1], alpha[i1 - 1]);
                    for (int i = 0; i < nv; i++) {
                        System.out.printf("%23sCall Option Values%2s",
                                " ", " ");
                        for (j = 0; j < nt; j++) {
                            System.out.printf(Locale.ENGLISH, "%7.4f ",
                                    splineValuesOption[j][i]);
                        }
                        System.out.printf("%n");
                    }
                    System.out.printf("%n");
                }
            }
        }
    }

    static class MyCoefficients implements FeynmanKac.PdeCoefficients {

        final double zero = 0.0, half = 0.5;
        private double sigma, strikePrice, interestRate;
        private double alpha, dividend;

        public MyCoefficients(double[] myData) {
            this.strikePrice = myData[0];
            this.sigma = myData[2];
            this.alpha = myData[3];
            this.interestRate = myData[4];
            this.dividend = myData[5];
        }

        // The coefficient sigma(x)
        public double sigma(double x, double t) {
            return (sigma * Math.pow(x, alpha * half));
        }

        // The coefficient derivative d(sigma) / dx
        public double sigmaPrime(double x, double t) {
            return (half * alpha * sigma * Math.pow(x, alpha * half - 1.0));
        }

        // The coefficient mu(x)
        public double mu(double x, double t) {
            return ((interestRate - dividend) * x);
        }

        // The coefficient kappa(x)
        public double kappa(double x, double t) {
            return interestRate;
        }
    }

    static class MyBoundaries implements FeynmanKac.Boundaries {

        private double xMax, df, interestRate, strikePrice;

        public MyBoundaries(double[] myData) {
            this.strikePrice = myData[0];
            this.xMax = myData[1];
            this.interestRate = myData[4];
        }

        public void leftBoundaries(double time, double[][] bndCoeffs) {
            bndCoeffs[0][0] = 1.0;
            bndCoeffs[0][1] = 0.0;
            bndCoeffs[0][2] = 0.0;
            bndCoeffs[0][3] = 0.0;
            bndCoeffs[1][0] = 0.0;
            bndCoeffs[1][1] = 1.0;
            bndCoeffs[1][2] = 0.0;
            bndCoeffs[1][3] = 0.0;
            bndCoeffs[2][0] = 0.0;
            bndCoeffs[2][1] = 0.0;
            bndCoeffs[2][2] = 1.0;
            bndCoeffs[2][3] = 0.0;
        }

        public void rightBoundaries(double time, double[][] bndCoeffs) {
            df = Math.exp(interestRate * time);
            bndCoeffs[0][0] = 1.0;
            bndCoeffs[0][1] = 0.0;
            bndCoeffs[0][2] = 0.0;
            bndCoeffs[0][3] = xMax - df * strikePrice;
            bndCoeffs[1][0] = 0.0;
            bndCoeffs[1][1] = 1.0;
            bndCoeffs[1][2] = 0.0;
            bndCoeffs[1][3] = 1.0;
            bndCoeffs[2][0] = 0.0;
            bndCoeffs[2][1] = 0.0;
            bndCoeffs[2][2] = 1.0;
            bndCoeffs[2][3] = 0.0;
        }

        public double terminal(double x) {
            double zero = 0.0;
            // The payoff function
            double value = Math.max(x - strikePrice, zero);
            return value;
        }
    }
}

Output

  Constant Elasticity of Variance Model for Vanilla Call
       Interest Rate:  0.050   Continuous Dividend:  0.000
       Minimum and Maximum Prices of Underlying:   0.00  60.00
       Number of equally spaced spline knots: 120
       Number of unknowns: 363

       Time in Years Prior to Expiration:  0.0833 0.3333 0.5833
       Option valued at Underlying Prices:  19.00  20.00  21.00

  Strike=15.00, Sigma= 0.20, Alpha= 2.00
                       Call Option Values   4.0624  4.2576  4.4734 
                       Call Option Values   5.0624  5.2505  5.4492 
                       Call Option Values   6.0624  6.2486  6.4386 

  Strike=20.00, Sigma= 0.20, Alpha= 2.00
                       Call Option Values   0.1312  0.5951  0.9693 
                       Call Option Values   0.5024  1.0880  1.5093 
                       Call Option Values   1.1980  1.7478  2.1745 

  Strike=25.00, Sigma= 0.20, Alpha= 2.00
                       Call Option Values   0.0000  0.0111  0.0751 
                       Call Option Values   0.0000  0.0376  0.1630 
                       Call Option Values   0.0006  0.1036  0.3150 

  Strike=15.00, Sigma= 0.30, Alpha= 2.00
                       Call Option Values   4.0636  4.3405  4.6627 
                       Call Option Values   5.0625  5.2951  5.5794 
                       Call Option Values   6.0624  6.2712  6.5248 

  Strike=20.00, Sigma= 0.30, Alpha= 2.00
                       Call Option Values   0.3107  1.0261  1.5479 
                       Call Option Values   0.7317  1.5404  2.0999 
                       Call Option Values   1.3762  2.1674  2.7362 

  Strike=25.00, Sigma= 0.30, Alpha= 2.00
                       Call Option Values   0.0005  0.1124  0.3564 
                       Call Option Values   0.0035  0.2184  0.5565 
                       Call Option Values   0.0184  0.3869  0.8230 

  Strike=15.00, Sigma= 0.40, Alpha= 2.00
                       Call Option Values   4.0755  4.5143  4.9673 
                       Call Option Values   5.0660  5.4210  5.8328 
                       Call Option Values   6.0633  6.3588  6.7306 

  Strike=20.00, Sigma= 0.40, Alpha= 2.00
                       Call Option Values   0.5109  1.4625  2.1260 
                       Call Option Values   0.9611  1.9934  2.6915 
                       Call Option Values   1.5807  2.6088  3.3202 

  Strike=25.00, Sigma= 0.40, Alpha= 2.00
                       Call Option Values   0.0081  0.3302  0.7795 
                       Call Option Values   0.0287  0.5178  1.0656 
                       Call Option Values   0.0820  0.7690  1.4097 

  Strike=15.00, Sigma= 0.20, Alpha= 1.00
                       Call Option Values   4.0624  4.2479  4.4312 
                       Call Option Values   5.0624  5.2479  5.4312 
                       Call Option Values   6.0624  6.2479  6.4312 

  Strike=20.00, Sigma= 0.20, Alpha= 1.00
                       Call Option Values   0.0000  0.0219  0.1051 
                       Call Option Values   0.1497  0.4107  0.6484 
                       Call Option Values   1.0832  1.3314  1.5773 

  Strike=25.00, Sigma= 0.20, Alpha= 1.00
                       Call Option Values  -0.0000 -0.0000  0.0000 
                       Call Option Values  -0.0000 -0.0000  0.0000 
                       Call Option Values  -0.0000 -0.0000  0.0000 

  Strike=15.00, Sigma= 0.30, Alpha= 1.00
                       Call Option Values   4.0624  4.2479  4.4312 
                       Call Option Values   5.0624  5.2479  5.4312 
                       Call Option Values   6.0624  6.2479  6.4312 

  Strike=20.00, Sigma= 0.30, Alpha= 1.00
                       Call Option Values   0.0010  0.0786  0.2208 
                       Call Option Values   0.1993  0.4997  0.7539 
                       Call Option Values   1.0835  1.3444  1.6022 

  Strike=25.00, Sigma= 0.30, Alpha= 1.00
                       Call Option Values  -0.0000  0.0000  0.0000 
                       Call Option Values  -0.0000  0.0000  0.0000 
                       Call Option Values  -0.0000  0.0000  0.0004 

  Strike=15.00, Sigma= 0.40, Alpha= 1.00
                       Call Option Values   4.0624  4.2479  4.4312 
                       Call Option Values   5.0624  5.2479  5.4312 
                       Call Option Values   6.0624  6.2479  6.4312 

  Strike=20.00, Sigma= 0.40, Alpha= 1.00
                       Call Option Values   0.0072  0.1540  0.3446 
                       Call Option Values   0.2498  0.5950  0.8728 
                       Call Option Values   1.0868  1.3795  1.6586 

  Strike=25.00, Sigma= 0.40, Alpha= 1.00
                       Call Option Values   0.0000  0.0000  0.0000 
                       Call Option Values   0.0000  0.0000  0.0005 
                       Call Option Values   0.0000  0.0002  0.0057 

  Strike=15.00, Sigma= 0.20, Alpha= 0.00
                       Call Option Values   4.0624  4.2479  4.4311 
                       Call Option Values   5.0624  5.2479  5.4312 
                       Call Option Values   6.0624  6.2479  6.4312 

  Strike=20.00, Sigma= 0.20, Alpha= 0.00
                       Call Option Values   0.0001  0.0001  0.0002 
                       Call Option Values   0.0816  0.3316  0.5748 
                       Call Option Values   1.0817  1.3308  1.5748 

  Strike=25.00, Sigma= 0.20, Alpha= 0.00
                       Call Option Values   0.0000 -0.0000 -0.0000 
                       Call Option Values   0.0000 -0.0000 -0.0000 
                       Call Option Values  -0.0000  0.0000 -0.0000 

  Strike=15.00, Sigma= 0.30, Alpha= 0.00
                       Call Option Values   4.0624  4.2479  4.4312 
                       Call Option Values   5.0624  5.2479  5.4312 
                       Call Option Values   6.0624  6.2479  6.4312 

  Strike=20.00, Sigma= 0.30, Alpha= 0.00
                       Call Option Values   0.0000 -0.0000  0.0026 
                       Call Option Values   0.0894  0.3326  0.5753 
                       Call Option Values   1.0826  1.3306  1.5749 

  Strike=25.00, Sigma= 0.30, Alpha= 0.00
                       Call Option Values   0.0000 -0.0000  0.0000 
                       Call Option Values   0.0000 -0.0000  0.0000 
                       Call Option Values   0.0000 -0.0000 -0.0000 

  Strike=15.00, Sigma= 0.40, Alpha= 0.00
                       Call Option Values   4.0624  4.2479  4.4312 
                       Call Option Values   5.0624  5.2479  5.4312 
                       Call Option Values   6.0624  6.2479  6.4312 

  Strike=20.00, Sigma= 0.40, Alpha= 0.00
                       Call Option Values  -0.0000  0.0001  0.0108 
                       Call Option Values   0.0985  0.3383  0.5781 
                       Call Option Values   1.0830  1.3306  1.5749 

  Strike=25.00, Sigma= 0.40, Alpha= 0.00
                       Call Option Values   0.0000  0.0000  0.0000 
                       Call Option Values   0.0000 -0.0000  0.0000 
                       Call Option Values   0.0000 -0.0000  0.0000 

Link to Java source.