1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
|
/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
/*
Copyright (C) 2006 Ferdinando Ametrano
Copyright (C) 2006 Mario Pucci
Copyright (C) 2006 StatPro Italia srl
This file is part of QuantLib, a free-software/open-source library
for financial quantitative analysts and developers - http://quantlib.org/
QuantLib is free software: you can redistribute it and/or modify it
under the terms of the QuantLib license. You should have received a
copy of the license along with this program; if not, please email
<quantlib-dev@lists.sf.net>. The license is also available online at
<http://quantlib.org/license.shtml>.
This program is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
FOR A PARTICULAR PURPOSE. See the license for more details.
*/
#include <ql/termstructures/volatility/sabr.hpp>
#include <ql/utilities/dataformatters.hpp>
#include <ql/math/comparison.hpp>
#include <ql/errors.hpp>
namespace QuantLib {
Real unsafeSabrVolatility(Rate strike,
Rate forward,
Time expiryTime,
Real alpha,
Real beta,
Real nu,
Real rho) {
const Real oneMinusBeta = 1.0-beta;
const Real A = std::pow(forward*strike, oneMinusBeta);
const Real sqrtA= std::sqrt(A);
Real logM;
if (!close(forward, strike))
logM = std::log(forward/strike);
else {
const Real epsilon = (forward-strike)/strike;
logM = epsilon - .5 * epsilon * epsilon ;
}
const Real z = (nu/alpha)*sqrtA*logM;
const Real B = 1.0-2.0*rho*z+z*z;
const Real C = oneMinusBeta*oneMinusBeta*logM*logM;
const Real tmp = (std::sqrt(B)+z-rho)/(1.0-rho);
const Real xx = std::log(tmp);
const Real D = sqrtA*(1.0+C/24.0+C*C/1920.0);
const Real d = 1.0 + expiryTime *
(oneMinusBeta*oneMinusBeta*alpha*alpha/(24.0*A)
+ 0.25*rho*beta*nu*alpha/sqrtA
+(2.0-3.0*rho*rho)*(nu*nu/24.0));
Real multiplier;
// computations become precise enough if the square of z worth
// slightly more than the precision machine (hence the m)
static const Real m = 10;
if (std::fabs(z*z)>QL_EPSILON * m)
multiplier = z/xx;
else {
multiplier = 1.0 - 0.5*rho*z - (3.0*rho*rho-2.0)*z*z/12.0;
}
return (alpha/D)*multiplier*d;
}
void validateSabrParameters(Real alpha,
Real beta,
Real nu,
Real rho) {
QL_REQUIRE(alpha>0.0, "alpha must be positive: "
<< alpha << " not allowed");
QL_REQUIRE(beta>=0.0 && beta<=1.0, "beta must be in (0.0, 1.0): "
<< beta << " not allowed");
QL_REQUIRE(nu>=0.0, "nu must be non negative: "
<< nu << " not allowed");
QL_REQUIRE(rho*rho<1.0, "rho square must be less than one: "
<< rho << " not allowed");
}
Real sabrVolatility(Rate strike,
Rate forward,
Time expiryTime,
Real alpha,
Real beta,
Real nu,
Real rho) {
QL_REQUIRE(strike>0.0, "strike must be positive: "
<< io::rate(strike) << " not allowed");
QL_REQUIRE(forward>0.0, "at the money forward rate must be "
"positive: " << io::rate(forward) << " not allowed");
QL_REQUIRE(expiryTime>=0.0, "expiry time must be non-negative: "
<< expiryTime << " not allowed");
validateSabrParameters(alpha, beta, nu, rho);
return unsafeSabrVolatility(strike, forward, expiryTime,
alpha, beta, nu, rho);
}
}
|
