Calculating Greeks in C

In quantitative finance, the Greeks are the quantities representing the sensitivity of the price of derivatives such as options to a change in underlying parameters on which the value of an instrument or portfolio of financial instruments is dependent. The name is used because the most common of these sensitivities are often denoted by Greek letters. There are plenty of online resources for understanding the underlying mathematical formulas but here's a quick script

//
//  main.c
//  stat598a1
//
//  Created by Fizi Yadav on 1/20/14.
//  Copyright (c) 2014 Fizi Yadav. All rights reserved.
//

#include <stdio.h>
#include <math.h>


double getCDF(double);
double stdnorm(double);
double getd(double S, double K, double r, double v, double T, double t);
double getCallPrice(double S, double K, double r, double v, double T,double t);
double getDelta(double S, double K, double r, double v, double T, double t);
double getGamma(double S, double K, double r, double v, double T, double t);
double getVega(double S, double K, double r, double v, double T, double t);
double getTheta(double S, double K, double r, double v, double T, double t);
double getRho(double S, double K, double r, double v, double T, double t);
void printGreeks(double S, double K, double r, double v, double T, double t);


int main(int argc, const char * argv[])
{
    
    double S = 40.0;   // Option price
    double K = 45.0;   // Strike price
    double r = 0.08;    // Risk-free rate. Percentages must be provided as decimal eg: (5%) as (0.05)
    double v = 0.05;    // Volatility of the underlying. Percentages must be provided as decimal
    double T = 3.0;     // Time of maturity. Must be in years
    double t = 0.0;     // Start time
    
    printf("\nParameters:\n");
    printf("Underlying Asset Price: %lf\n", S);
    printf("Strike Price %lf\n", K);
    printf("Risk-Free Rate: %lf\n", r);
    printf("Volatility: %lf\n", v);
    printf("Time to maturity: %lf\n", T-t);
    
    printf("\nCall Price: %lf\n", getCallPrice(S, K, r, v, T, t));
    
    printGreeks(S, K, r, v, T, t);

    return 0;
}


//calculate normal CDF given x

double getCDF(double x){
    
    double b0 = 0.2316419;
    double b1 = 0.319381530;
    double b2 = -0.356563782;
    double b3 = 1.781477937;
    double b4 = -1.821255978;
    double b5 = 1.330274429;
    double t = 1/(1+b0*x);
    //double stdnorm = 0.398942*pow(2.71828, -0.5*pow(x, 2));
    
    double prod = b1*t+b2*pow(t,2)+b3*pow(t,3)+b4*pow(t,4)+b5*pow(t,5);
    
    if (x >= 0.0) {
        return (1-stdnorm(x)*prod);
    } else {
        return 1.0 - getCDF(-x);
    }
    
}

//calculate call price of a european option

double getCallPrice(double S, double K, double r, double v, double T, double t){
    
    double d1=getd(S, K, r, v, T, t);
    double d2 = d1 - v*sqrt(T-t);
    
    printf("d1: %lf\n", d1);
    printf("d2: %lf\n", d2);
    
    return (S*getCDF(d1))-(K*exp(-r*(T-t))*getCDF(d2));
    
}

//print the greeks of a european call option

void printGreeks(double S, double K, double r, double v, double T, double t){
    
    /* some sensitivities are quoted in scaled-down terms, to match the scale of likely changes in the parameters. Rho is reported divided by 100 , vega by 100 (1 vol point change), and theta by 365 or 252 (1 day decay based on either calendar days or trading days per year). */
    
    printf("\nGreeks:\n");
    
    printf("Delta: %lf\n", getDelta(S, K, r, v, T,t));
    printf("Gamma: %lf\n", getGamma(S, K, r, v, T,t));
    printf("Theta: %lf\n", getTheta(S, K, r, v, T,t)/365.0);
    printf("Vega: %lf\n", getVega(S, K, r, v, T,t)/100.0);
    printf("Rho: %lf\n", getRho(S, K, r, v, T,t)/100.0);
    
}


// Standard normal probability density function
double stdnorm(double x) {
    
    return (1.0/(pow(2*M_PI,0.5)))*exp(-0.5*x*x);
}

double getd(double S, double K, double r, double v, double T, double t){
    
    return (log(S/K)+(r+pow(v,2)*0.5)*(T-t))/(v*sqrt(T-t));
    
}



// Calculate the European call Delta

double getDelta(double S, double K, double r, double v, double T, double t) {
    return getCDF(getd(S, K, r, v, T, t));
}

// Calculate the European call Gamma

double getGamma(double S, double K, double r, double v, double T, double t) {
    return stdnorm(getd(S, K, r, v, T, t))/(S*v*sqrt(T-t));
}

// Calculate the European call Vega

double getVega(double S, double K, double r, double v, double T, double t) {
    return S*stdnorm(getd(S, K, r, v, T, t))*sqrt(T-t);
}

// Calculate the European call Theta

double getTheta(double S, double K, double r, double v, double T, double t) {
    
    double d1=getd(S, K, r, v, T, t);
    double d2 = d1 - v*sqrt(T-t);
    
    return (-(S*stdnorm(d1)*v)/(2*sqrt(T-t)))
    - (r*K*exp(-r*(T-t))*getCDF(d2));
}

// Calculate the European call Rho

double getRho(double S, double K, double r, double v, double T, double t) {
    
    return K*(T-t)*exp(-r*(T-t))*getCDF(getd(S, K, r, v, T, t)-(v*sqrt(T-t)));
}