This example evaluates the price of a European Option using two payoff strategies: Cash-or-Nothing and Vertical Spread. In the first case the payoff function is
using System; using Imsl.Math; public class FeynmanKacEx3 { public static void Main(String[] args) { int i, j; int nxGrid = 61; int ntGrid = 2; int nv = 12; // The strike price double strikePrice = 10.0; // The spread value double spreadValue = 15.0; // The Bet for the Cash-or-Nothing Call double bet = 2.0; // The sigma value double sigma = 0.4; // Time values for the options double[] tGrid = { 0.25, 0.5 }; // Values of the underlying where evaluations are made double[] evaluationPoints = new double[nv]; // Value of the interest rate and continuous dividend double r = 0.1, dividend = 0.0; // Values of the min and max underlying values modeled double xMin = 0.0, xMax = 30.0; // Define parameters for the integration step. 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; // Structure for the evaluation routines. int iData; double[] rData = new double[7]; bool[] timeDependence = new bool[7]; // Define an equally-spaced grid of points for the // underlying price dx = (xMax - xMin) / ((double)(nint)); xGrid[0] = xMin; xGrid[nxGrid - 1] = xMax; for (i = 2; i <= nxGrid - 1; i++) { xGrid[i - 1] = xGrid[i - 2] + dx; } for (i = 1; i <= nv; i++) { evaluationPoints[i - 1] = 2.0 + (i - 1) * 2.0; } rData[0] = strikePrice; rData[1] = bet; rData[2] = spreadValue; rData[3] = xMax; rData[4] = sigma; rData[5] = r; rData[6] = dividend; // Flag the difference in payoff functions // 1 for the Bet, 2 for the Vertical Spread // In both cases, no time dependencies for // the coefficients and the left boundary // conditions timeDependence[5] = true; iData = 1; FeynmanKac betOption = new FeynmanKac(new MyCoefficients(rData)); betOption.SetTimeDependence(timeDependence); betOption.ComputeCoefficients(nlbcd, nrbcd, new MyBoundaries(rData, iData), xGrid, tGrid); double[,] betOptionCoefficients = betOption.GetSplineCoefficients(); double[][] splineValuesBetOption = new double[ntGrid][]; double[] betOptionTimeCoeffs = new double[n]; iData = 2; FeynmanKac spreadOption = new FeynmanKac(new MyCoefficients(rData)); spreadOption.SetTimeDependence(timeDependence); spreadOption.ComputeCoefficients(nlbcd, nrbcd, new MyBoundaries(rData, iData), xGrid, tGrid); double[,] spreadOptionCoefficients = spreadOption.GetSplineCoefficients(); double[][] splineValuesSpreadOption = new double[ntGrid][]; double[] spreadOptionTimeCoeffs = new double[n]; // Evaluate solutions at vector evaluationPoints, at each time value // prior to expiration. for (i = 0; i < ntGrid; i++) { for (j = 0; j < n; j++) { betOptionTimeCoeffs[j] = betOptionCoefficients[i + 1, j]; spreadOptionTimeCoeffs[j] = spreadOptionCoefficients[i + 1, j]; } splineValuesBetOption[i] = betOption.GetSplineValue( evaluationPoints, betOptionTimeCoeffs, 0); splineValuesSpreadOption[i] = spreadOption.GetSplineValue( evaluationPoints, spreadOptionTimeCoeffs, 0); } Console.Out.WriteLine(" European Option Value for A Bet"); Console.Out.WriteLine( " and a Vertical Spread, 3 and 6 Months prior to Expiry"); Console.Out.WriteLine( " Number of equally spaced spline knots:{0,4:d}", nxGrid); Console.Out.WriteLine(" Number of unknowns:{0,4:d}", n); Console.Out.WriteLine( " Strike={0,5:f2}, Sigma={1,5:f2}, Interest Rate={2,5:f2}", strikePrice, sigma, r); Console.Out.WriteLine(" Bet={0,5:f2}, Spread Value={1,5:f2}", bet, spreadValue); Console.Out.WriteLine(); Console.Out.WriteLine( " Underlying A Bet Vertical Spread"); for (i = 0; i < nv; i++) { Console.Out.WriteLine( " {0,9:f4}{1,9:f4}{2,9:f4}{3,9:f4}{4,9:f4}", evaluationPoints[i], splineValuesBetOption[0][i], splineValuesBetOption[1][i], splineValuesSpreadOption[0][i], splineValuesSpreadOption[1][i]); } } // These routines define the coefficients, payoff, boundary conditions // and forcing term for American and European Options. class MyCoefficients : FeynmanKac.IPdeCoefficients { const double zero = 0.0; double sigma, strikePrice, interestRate; double spread, bet, dividend; public MyCoefficients(double[] rData) { this.strikePrice = rData[0]; this.bet = rData[1]; this.spread = rData[2]; this.sigma = rData[4]; this.interestRate = rData[5]; this.dividend = rData[6]; } // The coefficient sigma(x) public double Sigma(double x, double t) { return (sigma * x); } // The coefficient derivative d(sigma) / dx public double SigmaPrime(double x, double t) { return sigma; } // 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; } } class MyBoundaries : FeynmanKac.IBoundaries { double strikePrice, spread, bet, interestRate, df; int dataInt; public MyBoundaries(double[] rData, int iData) { this.strikePrice = rData[0]; this.bet = rData[1]; this.spread = rData[2]; this.interestRate = rData[5]; this.dataInt = iData; } 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; return; } public void RightBoundaries(double time, double[,] bndCoeffs) { // This is the discount factor using the risk-free // interest rate df = Math.Exp(interestRate * time); // Use flag passed to decide on boundary condition switch (dataInt) { case 1: bndCoeffs[0, 0] = 1.0; bndCoeffs[0, 1] = 0.0; bndCoeffs[0, 2] = 0.0; bndCoeffs[0, 3] = bet * df; break; case 2: bndCoeffs[0, 0] = 1.0; bndCoeffs[0, 1] = 0.0; bndCoeffs[0, 2] = 0.0; bndCoeffs[0, 3] = (spread - strikePrice) * df; break; } 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; return; } public double Terminal(double x) { const double zero = 0.0; double val = 0.0; switch (dataInt) { // The payoff function - Use flag passed to decide which case 1: // After reaching the strike price the payoff jumps // from zero to the bet value. val = zero; if (x > strikePrice) { val = bet; } break; case 2: /* Function is zero up to strike price. Then linear between strike price and spread. Then has constant value Spread-Strike Price after the value Spread. */ val = Math.Max(x - strikePrice, zero) - Math.Max(x - spread, zero); break; } return val; } } }
European Option Value for A Bet and a Vertical Spread, 3 and 6 Months prior to Expiry Number of equally spaced spline knots: 61 Number of unknowns: 183 Strike=10.00, Sigma= 0.40, Interest Rate= 0.10 Bet= 2.00, Spread Value=15.00 Underlying A Bet Vertical Spread 2.0000 0.0000 0.0000 0.0000 0.0000 4.0000 0.0000 0.0013 0.0000 0.0005 6.0000 0.0112 0.0729 0.0038 0.0452 8.0000 0.2686 0.4291 0.1486 0.3833 10.0000 0.9948 0.9781 0.8898 1.1907 12.0000 1.6103 1.4301 2.1904 2.2267 14.0000 1.8650 1.6926 3.4267 3.1567 16.0000 1.9335 1.8171 4.2274 3.8282 18.0000 1.9477 1.8696 4.6261 4.2499 20.0000 1.9501 1.8902 4.7903 4.4918 22.0000 1.9505 1.8979 4.8493 4.6222 24.0000 1.9506 1.9008 4.8685 4.6901Link to C# source.