Then it will be a good exercise for you (and I will not be surprised, if for your prof too).
Here is a sample script to calculate a simple barrier option in BS-world, try to understand it and play around with it to see how many paths you need to make the simulation results converge
(Convergence as such is yet not a proof of the solution correctness the absence of convergence is a good disprove).
#calculates the price of the Barrier Call Option
#params: annualized vola, annualized interest rate, time to maturity, strike,
# Barrier, number of bank days in year, number of simulated paths
#seed for RNG (to make the simulation reproducible)
barrierOption <- function( sigma, r, T, So, K, B, N_BISDAYS, N_SIM, seed)
{
set.seed( seed );
N_DAYS = T * N_BISDAYS #trading days
vd = sigma / sqrt( N_BISDAYS )
drift = exp( (r - 0.5*sigma*sigma) * seq(1, N_DAYS)/N_BISDAYS )
payOffs = array(0, dim=N_SIM)
for (sz in 1:N_SIM)
{
path = array(0, dim=N_DAYS)
pathBM = vd*rnorm(N_DAYS)
path = So * exp( cumsum(pathBM) ) * drift
KOs = which( path <= B)
if( length(KOs) > 0 ) #Barrier hit (at least once)
payOffs[sz] = 0
else
payOffs[sz] = max( (path[N_DAYS] - K), 0 )
}
discountedPayOffs = payOffs / exp(r*T)
optionPrice = mean(discountedPayOffs)
return( optionPrice )
}
#Seed-dependent results
barrierOption (sigma=0.3, r=0.01, T=5, So=10, K=20, B=3, N_BISDAYS=242, N_SIM=1e3, seed=12345)
barrierOption (sigma=0.3, r=0.01, T=5, So=10, K=20, B=3, N_BISDAYS=242, N_SIM=1e3, seed=67890)
#Seed-independent results
barrierOption (sigma=0.3, r=0.01, T=5, So=10, K=20, B=3, N_BISDAYS=242, N_SIM=1e7, seed=12345)
barrierOption (sigma=0.3, r=0.01, T=5, So=10, K=20, B=3, N_BISDAYS=242, N_SIM=1e7, seed=67890)
#tricks to accelerate computation in R (compilation and multi-threading)
library(compiler)
enableJIT(3) #maximum optimization of JIT compiler
boFast <- cmpfun(barrierOption) #Compiled function runs faster
install.packages("snow")
install.packages("doSNOW")
install.packages("foreach")
require(snow)
require(doSNOW)
require(foreach)
N_CORES = 10 #number of CPU cores, change accordingly
cl = makeCluster(rep("localhost",N_CORES), type="SOCK")
registerDoSNOW(cl)
res <- foreach(k=1:N_CORES) %dopar% boFast(0.3, 0.01, 5, 10, 20, 3, 242, 1e6, 12345)
print(paste("Startwert ist 12345, Optionpreis ist ", mean(as.numeric(res))))
res <- foreach(k=1:N_CORES) %dopar% boFast(0.3, 0.01, 5, 10, 20, 3, 242, 1e6, 67890)
print(paste("Startwert ist 67890, Optionpreis ist ", mean(as.numeric(res))))
stopCluster(cl)