Black and Scholes formula

Here is an implementation of the Black and Scholes formula for simple and digital European options in F#:

//pricing.fs file

#light

open System
open Options

// Abramowitz and Stegun formula for the
// standard normal cumulative distribution
let normalCDF x =
 let b1 =  0.319381530
 let b2 = -0.356563782
 let b3 =  1.781477937
 let b4 = -1.821255978
 let b5 =  1.330274429
 let p  =  0.2316419
 let c  =  0.39894228
 let RationalApproximation t =
  c * exp( -x * x / 2.0 ) *
  t *( t *( t * ( t * ( t * b5 + b4 ) + b3 ) + b2 ) + b1 )
 if (x >= 0.0)
 then 1.0 - RationalApproximation(1.0 / (1.0 + p * x))
 else RationalApproximation(1.0 / (1.0 - p * x))

//Black and Scholes model parameters
type BlackScholesParameters =
{ S0: float;     //current underlying price
  sigma: float;  //volatility
  r: float;}     //risk free rate


let blackSholesPricing myOption parameters =
  match myOption with
 | SimpleOption (payoffType, K, T) ->
     let d1 = (Math.Log(parameters.S0/K)
       + (parameters.r + 0.5*parameters.sigma**2.0)*T)
       / (parameters.sigma * Math.Sqrt(T))
     let d2 = d1 - parameters.sigma * Math.Sqrt(T)
     match payoffType with
    | Call -> parameters.S0 * normalCDF(d1)
           - K * Math.Exp(-parameters.r * T) * normalCDF(d2)
    | Put  -> K * Math.Exp(-parameters.r * T) * normalCDF(-d2)
           - parameters.S0 * normalCDF(-d1)

 | DigitalOption (payoffType, K, T) ->
     let d2 = (Math.Log(parameters.S0/K)
          + (parameters.r + 0.5*parameters.sigma**2.0)*T)
          / (parameters.sigma * Math.Sqrt(T))
          - parameters.sigma * Math.Sqrt(T)
     match payoffType with
    | Call -> Math.Exp(-parameters.r * T) * normalCDF(d2)
    | Put  -> Math.Exp(-parameters.r * T) * normalCDF(-d2)

We can now price an option using Black and Scholes formula, and test the convergence of our Monte Carlo pricer:

//exe.fs file

#light

open System
open Options
open Pricing
open MonteCarlo

let pricingParameters = {S0=100.0; sigma=0.19; r=0.06}

let K = 102.5 // Strike
let T = 0.3   // Maturity in years
let option = SimpleOption(Call, K, T)

//BS price
let priceBS = blackSholesPricing option pricingParameters

let monteCarlo = simulPayoffs option pricingParameters.S0
             pricingParameters.r T pricingParameters.sigma

//MC price (10,000 simulations)
let priceMC1 = mean (monteCarlo 10000)
//MC price (100,000 simulations)
let priceMC2 = mean (monteCarlo 100000)
//MC price (1,000,000 simulations)
let priceMC3 = mean (monteCarlo 1000000)

//MC relative errors
let error1 = (priceMC1 - priceBS) / priceBS
let error2 = (priceMC2 - priceBS) / priceBS
let error3 = (priceMC3 - priceBS) / priceBS

printf "Black and Scholes price: %f\n" priceBS
printf "Monte Carlo (10,000 simulations) price: %f,
error = %f\n" priceMC1 error1
printf "Monte Carlo (100,000 simulations) price: %f,
error = %f\n" priceMC2 error2
printf "Monte Carlo (1,000,000 simulations) price: %f,
error = %f\n" priceMC3 error3

Output:

Black and Scholes price: 3.836593
MC (10,000 simulations) price: 3.784829, error = -0.013492
MC (100,000 simulations) price: 3.883840, error = 0.012315
MC (1,000,000 simulations) price: 3.833737, error = -0.000744