Changes in progress, added calculations for FinDiff greeks

src/black_scholes.c

@@ -5,20 +5,22 @@
 #include <math.h>
 #include <time.h>
 
-double bsm(double spot, double rfr, double vol, double strike,
+
-		struct tm expiry, struct tm value, int type)
+void bsm(struct Option *opt)
 {
 	double d1, d2, tte, price;
 	time_t value_date, expiry_date;
 
-	expiry_date = mktime(&expiry);
+	expiry_date = mktime(opt->expiry_date);
-	value_date = mktime(&value);
+	value_date = mktime(opt->value_date);
 	tte = difftime(expiry_date, value_date) / (60 * 60 * 24 * 365);
-	d1 = log(spot / strike) + tte * (rfr + pow(vol, 2) / 2)
-		/ (vol * pow(tte, 0.5));
-	d2 = d1 - (vol * pow(tte, 0.5));
 
-	price = spot * normalcdf(d1 * type) * type -
+	d1 = (log(opt->spot / opt->strike) + tte * (opt->rfr + pow(opt->vol, 2) / 2))
-		strike * exp(-rfr * tte) * normalcdf(d2 * type) * type;
+		/ (opt->vol * pow(tte, 0.5));
-	return price;
+	d2 = d1 - (opt->vol * pow(tte, 0.5));
+
+	price = opt->spot * normalcdf(d1 * opt->type) * opt->type - opt->strike *
+		exp(-opt->rfr * tte) * normalcdf(d2 * opt->type) * opt->type;
+	opt->fv = price;
+	opt->vega = opt->spot * exp(-opt->rfr * tte) * pow(tte, 0.5) * normalpdf(d1);
 }

src/black_scholes.h

@@ -1,14 +1,7 @@
 #ifndef BLACK_SCHOLES_H_
 #define BLACK_SCHOLES_H_
 
-struct tm;
+struct Option;
-
+void bsm(struct Option *option);
-double gbm(double spot, double rfr, double vol, double strike, struct tm expiry,
-		struct tm value, int type, int sims);
-
-double normalcdf(double value);
-
-double bsm(double spot, double rfr, double vol, double strike, struct tm expiry,
-		struct tm value, int type);
 
 #endif

src/gbm_mc.c

@@ -17,26 +17,43 @@ double gbm_simulation(double spot, double rfr, double vol, double tte, double ra
 	return spot * exp(drift + stoch);
 }
 
-double gbm(double spot, double rfr, double vol, double strike, struct tm expiry,
+void gbm(struct Option *opt)
-		struct tm value, int type, int sims)
 {
-	double tte, expiry_date, value_date, level, price, rand;
+	double tte, expiry_date, value_date, level, rand;
-	double results = 0;
+	double delta_shift = 0, vega_shift = 0, theta_shift = 0, rho_shift = 0;
+	double price, delta, vega, theta, rho;
+	double base = 0;
 	int i;
 
-	if (sims < 1) sims = 1;
+	if (opt->sims < 1) opt->sims = 1;
 
-	expiry_date = mktime(&expiry);
+	expiry_date = mktime(opt->expiry_date);
-	value_date = mktime(&value);
+	value_date = mktime(opt->value_date);
 	tte = difftime(expiry_date, value_date) / (60 * 60 * 24 * 365);
 
-	for (i=0; i<sims; i++) {
+	for (i=0; i<opt->sims; i++) {
 		rand = gaussrand();
-		level = gbm_simulation(spot, rfr, vol, tte, rand);
+		level = gbm_simulation(opt->spot, opt->rfr, opt->vol, tte, rand);
-		results += max((level - strike) * type, 0);
+		base += max((level - opt->strike) * opt->type, 0);
+		level = gbm_simulation(opt->spot + 0.01, opt->rfr, opt->vol, tte, rand);
+		delta_shift += max((level - opt->strike) * opt->type, 0);
+		level = gbm_simulation(opt->spot, opt->rfr, opt->vol + 0.01, tte, rand);
+		vega_shift += max((level - opt->strike) * opt->type, 0);
+		level = gbm_simulation(opt->spot, opt->rfr, opt->vol, tte + 1/365, rand);
+		theta_shift += max((level - opt->strike) * opt->type, 0);
+		level = gbm_simulation(opt->spot, opt->rfr + 0.0001, opt->vol, tte, rand);
+		rho_shift += max((level - opt->strike) * opt->type, 0);
 	}
 
-	price = results / sims * 1 / pow((1 + rfr), tte);
+	price = base / opt->sims * 1 / pow((1 + opt->rfr), tte);
+	delta = delta_shift / opt->sims * 1 / pow((1 + opt->rfr), tte);
+	vega = vega_shift / opt->sims * 1 / pow((1 + opt->rfr), tte);
+	theta = theta_shift / opt->sims * 1 / pow((1 + opt->rfr), tte);
+	rho = rho_shift / opt->sims * 1 / pow((1 + opt->rfr), tte);
 
-	return price;
+	opt->fv = price;
+	opt->delta = (delta - price) * 100 * opt->spot;
+	opt->vega = (vega - price) * opt->spot;
+	opt->theta = (theta - price) * opt->spot;
+	opt->rho = (rho - price) * opt->spot;
 }

src/gbm_mc.h

@@ -3,8 +3,9 @@
 
 #include <time.h>
 
+struct Option;
 double gbm_simulation(double spot, double rfr, double vol, double tte, double rand);
 
-double gbm(double spot, double rfr, double vol, double strike, struct tm expiry, struct tm value, int type, int sims);
+void gbm(struct Option *opt);
 
 #endif

src/opt-pricer.c

@@ -2,6 +2,7 @@
 
 #include "black_scholes.h"
 #include "gbm_mc.h"
+#include "utils.h"
 
 #ifdef __MINGW32__
 #include "strptime.h"
@@ -72,13 +73,13 @@ int print_help(void)
 int main(int argc, char *argv[])
 {
 	double spot = 0, strike = 0, rfr = 0, vol = 0, sims = 1000;
-	double bs_price, mc_price;
+	char expiry_date[11], value_date[11], buffer[1024];
-	char expiry_date[11], buffer[512];
 	int opt, option_index = 0, type = 0;
-	struct tm expiry, value;
+	struct tm value, expiry;
+	struct Option bs_opt, mc_opt;
 
 	extern char *optarg;
-    extern int getopt();
+	extern int getopt_long();
 
 	memset(&expiry, 0, sizeof(expiry));
 	memset(&value, 0, sizeof(value));
@@ -132,19 +133,56 @@ int main(int argc, char *argv[])
 		return 1;
 	}
 
+	bs_opt.spot = spot;
+	bs_opt.strike = strike;
+	bs_opt.expiry_date = &expiry;
+	bs_opt.value_date = &value;
+	bs_opt.rfr = rfr;
+	bs_opt.vol = vol;
+	bs_opt.type = type;
+
+	mc_opt.spot = spot;
+	mc_opt.strike = strike;
+	mc_opt.expiry_date = &expiry;
+	mc_opt.value_date = &value;
+	mc_opt.rfr = rfr;
+	mc_opt.vol = vol;
+	mc_opt.type = type;
+	mc_opt.sims = sims;
+
 	strftime(expiry_date, 11, "%Y-%m-%d", &expiry);
-	bs_price = bsm(spot, rfr, vol, strike, expiry, value, type);
+	strftime(value_date, 11, "%Y-%m-%d", &value);
-	mc_price = gbm(spot, rfr, vol, strike, expiry, value, type, sims);
+	bsm(&bs_opt);
-	sprintf(&buffer[0], "Arguments:\n"
+	gbm(&mc_opt);
-			"spot: %f\n"
+	sprintf(&buffer[0],
-			"strike: %f\n"
+		"\nValuation date: %s\n\n"
-			"rfr: %f\n"
+		"	              | BS Analytic | BS Monte Carlo |\n"
-			"vol: %f\n"
+		"	 ---------------------------------------------\n"
-			"sims: %f\n"
+		"	 |Type:       | %10.1i  | %13.1i  |\n"
-			"exp: %s\n"
+		"	 |Spot:       | %10.2f  | %13.2f  |\n"
-			"BS price: %f\n"
+		"	 |Strike:     | %10.2f  | %13.2f  |\n"
-			"MC price: %f\n",
+		"	 |Risk-free:  | %10.2f%% | %13.2f%% |\n"
-			spot, strike, rfr, vol, sims, expiry_date, bs_price, mc_price);
+		"	 |Implied Vol:| %10.2f%% | %13.2f%% |\n"
+		/*
+		"	 |sims:     | %4.2f| %4.2f|\n"
+		"	 |exp:      | %4.2s| %4.2f|\n"
+		*/
+		"	 |Fair value: | %10.2f  | %13.2f  |\n"
+		"	 |Delta:      | %10.2f  | %13.2f  |\n"
+		"	 |Vega:       | %10.2f  | %13.2f  |\n"
+		"	 |Theta:      |   %8.4f  |   %11.4f  |\n"
+		"	 |Rho :       |   %8.4f  |   %11.4f  |\n\n",
+		value_date,
+		bs_opt.type, mc_opt.type,
+		bs_opt.spot, mc_opt.spot,
+		bs_opt.strike, mc_opt.strike,
+		bs_opt.rfr*100, mc_opt.rfr*100,
+		bs_opt.vol*100, mc_opt.vol*100,
+		bs_opt.fv, mc_opt.fv,
+		bs_opt.delta, mc_opt.delta,
+		bs_opt.vega, mc_opt.vega,
+		bs_opt.theta, mc_opt.theta,
+		bs_opt.rho, mc_opt.rho);
 	printf("%s", buffer);
 	return 0;
 }

src/utils.c

@@ -1,5 +1,10 @@
 #include <math.h>
 #include <stdlib.h>
+#include <time.h>
+
+
+#define M_PI	3.14159265358979323846264338327950288
+
 
 double normalcdf(double z)
 {
@@ -23,6 +28,13 @@ double normalcdf(double z)
 	return 0.5 * (1 + sign * y);
 }
 
+
+double normalpdf(double z)
+{
+	return exp(-0.5 * pow(z, 2)) / pow(M_PI * 2, 0.5);
+}
+
+
 double gaussrand()
 {
 	/* Marsaglia and Bray, ``A Convenient Method for Generating Normal Variables'' */

src/utils.h

@@ -3,6 +3,30 @@
 
 double normalcdf(double z);
 
+double normalpdf(double z);
+
 double gaussrand();
 
+struct Option {
+	/* option details */
+	int type;
+	double strike;
+	struct tm *expiry_date;
+
+	/* market data */
+	struct tm *value_date;
+	double spot;
+	double rfr;
+	double vol;
+	int sims;
+
+	/* fv and greeks */
+	double fv;
+	double delta;
+	double vega;
+	double rho;
+	double gamma;
+	double theta;
+};
+
 #endif