QuantLib: a free/open-source library for quantitative finance
fully annotated source code - version 1.34
Loading...
Searching...
No Matches
cevrndcalculator.cpp
Go to the documentation of this file.
1/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2
3/*
4 Copyright (C) 2018 Klaus Spanderen
5
6 This file is part of QuantLib, a free-software/open-source library
7 for financial quantitative analysts and developers - http://quantlib.org/
8
9 QuantLib is free software: you can redistribute it and/or modify it
10 under the terms of the QuantLib license. You should have received a
11 copy of the license along with this program; if not, please email
12 <quantlib-dev@lists.sf.net>. The license is also available online at
13 <http://quantlib.org/license.shtml>.
14
15 This program is distributed in the hope that it will be useful, but WITHOUT
16 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17 FOR A PARTICULAR PURPOSE. See the license for more details.
18*/
19
20/*! \file cevrndcalculator.cpp */
21
22#include <ql/errors.hpp>
27#include <boost/math/special_functions/gamma.hpp>
28#include <boost/math/distributions/non_central_chi_squared.hpp>
29
30namespace QuantLib {
31
33 : f0_(f0),
34 alpha_(alpha),
35 beta_(beta),
36 delta_((1.0-2.0*beta)/(1.0-beta)),
37 x0_(X(f0)) {
38 QL_REQUIRE(beta != 1.0, "beta can not be one");
39 }
40
42 if (delta_ < 2.0)
43 return 1.0-boost::math::gamma_p(-0.5*delta_+1.0,x0_/(2.0*t));
44 else
45 return 0.0;
46 }
47
49 return std::pow(f, 2.0*(1.0-beta_))/squared(alpha_*(1.0-beta_));
50 }
51
53 return std::pow(x*squared(alpha_*(1.0-beta_)),
54 1.0/(2.0*(1.0-beta_)));
55 }
56
58 const Real y = X(f);
59
60 if (delta_ < 2.0) {
61 return boost::math::pdf(
62 boost::math::non_central_chi_squared_distribution<Real>(
63 4.0-delta_, y/t), x0_/t)/t * 2.0*(1.0-beta_)*y/f;
64 }
65 else {
66 return boost::math::pdf(
67 boost::math::non_central_chi_squared_distribution<Real>(
68 delta_, x0_/t), y/t)/t * 2.0*(beta_-1.0)*y/f;
69 }
70 }
71
73 const Real y = X(f);
74
75 if (delta_ < 2.0)
76 return 1.0 - boost::math::cdf(
77 boost::math::non_central_chi_squared_distribution<Real>(
78 2.0-delta_, y/t), x0_/t);
79 else
80 return 1.0 - boost::math::cdf(
81 boost::math::non_central_chi_squared_distribution<Real>(
82 delta_, x0_/t), y/t);
83 }
84
86 const Real a = x0_/t;
87 const Real b = 2.0 - delta_;
88
89 c = std::max(c, -0.45*b);
90
91 const Real h = 1 - 2*(b+c)*(b+3*c)/(3*squared(b+2*c));
92 const Real p = (b+2*c)/squared(b+c);
93 const Real m = (h-1)*(1-3*h);
94
95 const Real u = (std::pow(a/(b+c), h) - (1 + h*p*(h-1-0.5*(2-h)*m*p)))/
96 (h*std::sqrt(2*p)*(1+0.5*m*p));
97
98 return u - x;
99 }
100
102 if (delta_ < 2.0) {
103 if (f0_ < QL_EPSILON || q < massAtZero(t))
104 return 0.0;
105
106 const Real x = InverseCumulativeNormal()(1-q);
107
108 auto cdfApprox = [&](Real _c){ return sankaranApprox(_c, t, x); };
109
110 const Real y0 = X(f0_)/t;
111
112 try {
113 Brent brent;
114 brent.setMaxEvaluations(20);
115 const Real guess =
116 invX(brent.solve(cdfApprox, 1e-8, y0, 0.02*y0) * t);
117
118 return InvCDFHelper(this, guess, 1e-8, 100).inverseCDF(q, t);
119 }
120 catch (...) {
121 return InvCDFHelper(this, f0_, 1e-8, 100).inverseCDF(q, t);
122 }
123 }
124 else {
125 const Real x = t * boost::math::quantile(
126 boost::math::non_central_chi_squared_distribution<Real>(
127 delta_, x0_/t), 1-q);
128 return invX(x);
129 }
130 }
131}
Brent 1-D solver.
risk neutral density calculator for the constant elasticity of variance (CEV) model
Brent 1-D solver
Definition: brent.hpp:37
CEVRNDCalculator(Real f0, Real alpha, Real beta)
Real pdf(Real f, Time t) const override
Real massAtZero(Time t) const
Real invcdf(Real q, Time t) const override
Real cdf(Real f, Time t) const override
Real sankaranApprox(Real f, Time t, Real x) const
Inverse cumulative normal distribution function.
void setMaxEvaluations(Size evaluations)
Definition: solver1d.hpp:238
Real solve(const F &f, Real accuracy, Real guess, Real step) const
Definition: solver1d.hpp:84
const DefaultType & t
Classes and functions for error handling.
#define QL_REQUIRE(condition, message)
throw an error if the given pre-condition is not verified
Definition: errors.hpp:117
ext::function< Real(Real)> b
#define QL_EPSILON
Definition: qldefines.hpp:178
Real Time
continuous quantity with 1-year units
Definition: types.hpp:62
QL_REAL Real
real number
Definition: types.hpp:50
functionals and combinators not included in the STL
Definition: any.hpp:35
T squared(T x)
Definition: functional.hpp:37
normal, cumulative and inverse cumulative distributions
ext::shared_ptr< YieldTermStructure > q
Real beta
Definition: sabr.cpp:200
Real alpha
Definition: sabr.cpp:200