# Homework: Change this python code to R code

from math import exp, sqrt, log

from scipy.stats import norm

import datetime

def find_vol(target_value, call_put, S, K, T, r):

MAX_ITERATIONS = 100

PRECISION = 1.0e-5

m = 1.0

sigma = 0.5

for i in range(0, MAX_ITERATIONS):

price = bs_price(call_put, S, K, T, r, sigma)

vega = bs_vega(call_put, S, K, T, r, sigma)

price = price

diff = target_value - price # our root

print (i, sigma, diff)

if (abs(diff) < PRECISION):

return sigma

sigma = sigma + m*diff/vega # f(x) / f'(x)

# value wasn't found, return best guess so far

return sigma

n = norm.pdf

N = norm.cdf

def bs_price(cp_flag,S,K,T,r,v,q=0.0):

F = S*exp((r-q)*T)

d1 = (log(F/K)+(v*v/2.)*T)/(v*sqrt(T))

d2 = d1-v*sqrt(T)

if cp_flag == 'c':

price = exp(-r*T)*(F*N(d1)-K*N(d2))

else:

price = exp(-r*T)*(K*N(-d2)-F*N(-d1))

return price

def bs_vega(cp_flag,S,K,T,r,v,q=0.0):

F = S*exp((r-q)*T)

d1 = (log(F/K)+(v*v/2.)*T)/(v*sqrt(T))

return exp(-r*T)* F * sqrt(T)*n(d1)

V_market = 17.5

K = 585

T = (datetime.date(2014,10,18) - datetime.date(2014,9,8)).days / 365.

S = 586.08

r = 0.0002

cp = 'c' # call option

implied_vol = find_vol(V_market, cp, S, K, T, r)

print ('Implied vol: %.2f%%' % (implied_vol * 100))

print ('Market price = %.2f' % V_market)

print ('Model price = %.2f' % bs_price(cp, S, K, T, r, implied_vol))