Skip to content

Instantly share code, notes, and snippets.

@Fs02
Created February 2, 2017 19:02
Show Gist options
  • Select an option

  • Save Fs02/44ed2948bcab49ff0e5e40e23ed13355 to your computer and use it in GitHub Desktop.

Select an option

Save Fs02/44ed2948bcab49ff0e5e40e23ed13355 to your computer and use it in GitHub Desktop.
Test Levy Random
// Example program
#include <iostream>
#include <string>
#include <random>
#include <chrono>
// https://www.quantstart.com/articles/Statistical-Distributions-in-C
double inv_cdf(const double& quantile) {
// This is the Beasley-Springer-Moro algorithm which can
// be found in Glasserman [2004]. We won't go into the
// details here, so have a look at the reference for more info
static double a[4] = { 2.50662823884,
-18.61500062529,
41.39119773534,
-25.44106049637};
static double b[4] = { -8.47351093090,
23.08336743743,
-21.06224101826,
3.13082909833};
static double c[9] = {0.3374754822726147,
0.9761690190917186,
0.1607979714918209,
0.0276438810333863,
0.0038405729373609,
0.0003951896511919,
0.0000321767881768,
0.0000002888167364,
0.0000003960315187};
if (quantile >= 0.5 && quantile <= 0.92)
{
double num = 0.0;
double denom = 1.0;
for (int i=0; i<4; i++)
{
num += a[i] * pow((quantile - 0.5), 2*i + 1);
denom += b[i] * pow((quantile - 0.5), 2*i);
}
return num/denom;
}
else if (quantile > 0.92 && quantile < 1)
{
double num = 0.0;
for (int i=0; i<9; i++)
{
num += c[i] * pow((log(-log(1-quantile))), i);
}
return num;
}
else
{
return -1.0*inv_cdf(1-quantile);
}
}
int main()
{
int seed = std::chrono::system_clock::now().time_since_epoch().count();
std::default_random_engine generator(seed);
std::uniform_real_distribution<double> uniform_real(0.0, 1.0);
const double a = 0.0;
const double b = 2.0;
for (std::size_t i = 0; i < 10; ++i)
std::cout << a + b/std::pow(inv_cdf(1-uniform_real(generator)/2.0),2) << std::endl;
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment