Line | Count | Source (jump to first uncovered line) |
1 | /* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */ | |
2 | ||
3 | /* | |
4 | Copyright (C) 2014 Fabien Le Floc'h | |
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 analytichestonengine.hpp | |
21 | \brief analytic Heston expansion engine | |
22 | */ | |
23 | ||
24 | #include <ql/pricingengines/blackformula.hpp> | |
25 | #include <ql/instruments/payoffs.hpp> | |
26 | #include <ql/pricingengines/vanilla/hestonexpansionengine.hpp> | |
27 | ||
28 | #if defined(QL_PATCH_MSVC) | |
29 | #pragma warning(disable: 4180) | |
30 | #endif | |
31 | ||
32 | using std::exp; | |
33 | using std::pow; | |
34 | using std::log; | |
35 | using std::sqrt; | |
36 | ||
37 | namespace QuantLib { | |
38 | ||
39 | HestonExpansionEngine::HestonExpansionEngine( | |
40 | const ext::shared_ptr<HestonModel>& model, | |
41 | HestonExpansionFormula formula) | |
42 | 0 | : GenericModelEngine<HestonModel, |
43 | 0 | VanillaOption::arguments, |
44 | 0 | VanillaOption::results>(model), |
45 | 0 | formula_(formula) { |
46 | 0 | } |
47 | ||
48 | void HestonExpansionEngine::calculate() const | |
49 | 0 | { |
50 | // this is a european option pricer | |
51 | 0 | QL_REQUIRE(arguments_.exercise->type() == Exercise::European, |
52 | 0 | "not an European option"); |
53 | ||
54 | // plain vanilla | |
55 | 0 | ext::shared_ptr<PlainVanillaPayoff> payoff = |
56 | 0 | ext::dynamic_pointer_cast<PlainVanillaPayoff>(arguments_.payoff); |
57 | 0 | QL_REQUIRE(payoff, "non plain vanilla payoff given"); |
58 | ||
59 | 0 | const ext::shared_ptr<HestonProcess>& process = model_->process(); |
60 | ||
61 | 0 | const Real riskFreeDiscount = process->riskFreeRate()->discount( |
62 | 0 | arguments_.exercise->lastDate()); |
63 | 0 | const Real dividendDiscount = process->dividendYield()->discount( |
64 | 0 | arguments_.exercise->lastDate()); |
65 | ||
66 | 0 | const Real spotPrice = process->s0()->value(); |
67 | 0 | QL_REQUIRE(spotPrice > 0.0, "negative or null underlying given"); |
68 | ||
69 | 0 | const Real strikePrice = payoff->strike(); |
70 | 0 | const Real term = process->time(arguments_.exercise->lastDate()); |
71 | ||
72 | //possible optimization: | |
73 | // if term=lastTerm & model=lastModel & formula=lastApprox, reuse approx. | |
74 | 0 | const Real forward = spotPrice*dividendDiscount/riskFreeDiscount; |
75 | 0 | Real vol; |
76 | 0 | switch(formula_) { |
77 | 0 | case LPP2: { |
78 | 0 | LPP2HestonExpansion expansion(model_->kappa(), model_->theta(), |
79 | 0 | model_->sigma(), model_->v0(), |
80 | 0 | model_->rho(), term); |
81 | 0 | vol = expansion.impliedVolatility(strikePrice, forward); |
82 | 0 | break; |
83 | 0 | } |
84 | 0 | case LPP3: { |
85 | 0 | LPP3HestonExpansion expansion(model_->kappa(), model_->theta(), |
86 | 0 | model_->sigma(), model_->v0(), |
87 | 0 | model_->rho(), term); |
88 | 0 | vol = expansion.impliedVolatility(strikePrice, forward); |
89 | 0 | break; |
90 | 0 | } |
91 | 0 | case Forde: { |
92 | 0 | FordeHestonExpansion expansion(model_->kappa(), model_->theta(), |
93 | 0 | model_->sigma(), model_->v0(), |
94 | 0 | model_->rho(), term); |
95 | 0 | vol = expansion.impliedVolatility(strikePrice, forward); |
96 | 0 | break; |
97 | 0 | } |
98 | 0 | default: |
99 | 0 | QL_FAIL("unknown expansion formula"); |
100 | 0 | } |
101 | 0 | const Real price = blackFormula(payoff, forward, vol*sqrt(term), |
102 | 0 | riskFreeDiscount, 0); |
103 | 0 | results_.value = price; |
104 | 0 | } |
105 | ||
106 | LPP2HestonExpansion::LPP2HestonExpansion(Real kappa, Real theta, Real sigma, | |
107 | 0 | Real v0, Real rho, Real term) { |
108 | 0 | ekt = exp(kappa*term); |
109 | 0 | e2kt = ekt*ekt; |
110 | 0 | e3kt = e2kt*ekt; |
111 | 0 | e4kt = e2kt*e2kt; |
112 | 0 | coeffs[0] = z0(term, kappa, theta, sigma, v0, rho); |
113 | 0 | coeffs[1] = z1(term, kappa, theta, sigma, v0, rho); |
114 | 0 | coeffs[2] = z2(term, kappa, theta, sigma, v0, rho); |
115 | 0 | } |
116 | ||
117 | 0 | static Real fastpow(Real x, int y) { |
118 | 0 | return pow(x,y); |
119 | 0 | } |
120 | ||
121 | Real LPP2HestonExpansion::impliedVolatility(const Real strike, | |
122 | 0 | const Real forward) const { |
123 | 0 | Real x = log(strike/forward); |
124 | 0 | Real vol = coeffs[0]+x*(coeffs[1]+(x*coeffs[2])); |
125 | 0 | return std::max(1e-8,vol); |
126 | 0 | } |
127 | ||
128 | Real LPP2HestonExpansion::z0(Real t, Real kappa, Real theta, | |
129 | 0 | Real delta, Real y, Real rho) const { |
130 | 0 | return (4*pow(delta,2)*kappa*(-theta - 4*ekt*(theta + kappa*t*(theta - y)) + |
131 | 0 | e2kt*((5 - 2*kappa*t)*theta - 2*y) + 2*y)* |
132 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y) + |
133 | 0 | 128*ekt*pow(kappa,3)* |
134 | 0 | pow((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y,2) + |
135 | 0 | 32*delta*ekt*pow(kappa,2)*rho* |
136 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y)* |
137 | 0 | ((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
138 | 0 | (-1 + ekt - kappa*t)*y) + |
139 | 0 | pow(delta,2)*ekt*pow(rho,2)* |
140 | 0 | (-theta + kappa*t*theta + (theta - y)/ekt + y)* |
141 | 0 | pow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
142 | 0 | (-1 + ekt - kappa*t)*y,2) + |
143 | 0 | (48*pow(delta,2)*e2kt*pow(kappa,2)*pow(rho,2)* |
144 | 0 | pow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
145 | 0 | (-1 + ekt - kappa*t)*y,2))/ |
146 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y) - |
147 | 0 | pow(delta,2)*pow(rho,2)*((1 + ekt*(-1 + kappa*t))*theta + |
148 | 0 | (-1 + ekt)*y)*pow((2 + kappa*t + ekt*(-2 + kappa*t))* |
149 | 0 | theta + (-1 + ekt - kappa*t)*y,2) + |
150 | 0 | 2*pow(delta,2)*kappa*((1 + ekt*(-1 + kappa*t))*theta + |
151 | 0 | (-1 + ekt)*y)*(theta - 2*y + |
152 | 0 | e2kt*(-5*theta + 2*kappa*t*theta + 2*y + |
153 | 0 | 8*pow(rho,2)*((-3 + kappa*t)*theta + y)) + |
154 | 0 | 4*ekt*(theta + kappa*t*theta - kappa*t*y + |
155 | 0 | pow(rho,2)*((6 + kappa*t*(4 + kappa*t))*theta - (2 + kappa*t*(2 + kappa*t))*y))) - |
156 | 0 | (8*pow(delta,2)*pow(kappa,2)*((1 + ekt*(-1 + kappa*t))*theta + |
157 | 0 | (-1 + ekt)*y)*(theta - 2*y + |
158 | 0 | e2kt*(-5*theta + 2*kappa*t*theta + 2*y + |
159 | 0 | 8*pow(rho,2)*((-3 + kappa*t)*theta + y)) + |
160 | 0 | 4*ekt*(theta + kappa*t*theta - kappa*t*y + |
161 | 0 | pow(rho,2)*((6 + kappa*t*(4 + kappa*t))*theta - (2 + kappa*t*(2 + kappa*t))*y)))) |
162 | 0 | /(-theta + kappa*t*theta + (theta - y)/ekt + y))/ |
163 | 0 | (128.*e3kt*pow(kappa,5)*pow(t,2)* |
164 | 0 | pow((-theta + kappa*t*theta + (theta - y)/ekt + y)/(kappa*t),1.5)); |
165 | 0 | } |
166 | ||
167 | Real LPP2HestonExpansion::z1(Real t, Real kappa, Real theta, | |
168 | 0 | Real delta, Real y, Real rho) const { |
169 | 0 | return (delta*rho*(-(delta*fastpow(-1 + ekt,2)*rho*(4*theta - y)*y) + |
170 | 0 | 2*ekt*fastpow(kappa,3)*fastpow(t,2)*theta* |
171 | 0 | ((2 + 2*ekt + delta*rho*t)*theta - (2 + delta*rho*t)*y) - |
172 | 0 | 2*(-1 + ekt)*kappa*(2*theta - y)* |
173 | 0 | ((-1 + ekt)*(-2 + delta*rho*t)*theta + |
174 | 0 | (-2 + 2*ekt + delta*rho*t)*y) + |
175 | 0 | fastpow(kappa,2)*t*((-1 + ekt)* |
176 | 0 | (-4 + delta*rho*t + ekt*(-12 + delta*rho*t))*fastpow(theta,2) + |
177 | 0 | 2*(-4 + 4*e2kt + delta*rho*t + 3*delta*ekt*rho*t)*theta* |
178 | 0 | y - (-4 + delta*rho*t + 2*ekt*(2 + delta*rho*t))*fastpow(y,2))))/ |
179 | 0 | (8.*fastpow(kappa,2)*t*sqrt((-theta + kappa*t*theta + (theta - y)/ekt + y)/ |
180 | 0 | (kappa*t))*fastpow((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y, |
181 | 0 | 2)); |
182 | 0 | } |
183 | ||
184 | Real LPP2HestonExpansion::z2(Real t, Real kappa, Real theta, | |
185 | 0 | Real delta, Real y, Real rho) const { |
186 | 0 | return (fastpow(delta,2)*sqrt((-theta + kappa*t*theta + (theta - y)/ekt + y)/(kappa*t))* |
187 | 0 | (-12*fastpow(rho,2)*fastpow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
188 | 0 | (-1 + ekt - kappa*t)*y,2) + |
189 | 0 | (-theta + kappa*t*theta + (theta - y)/ekt + y)* |
190 | 0 | (theta - 2*y + e2kt* |
191 | 0 | (-5*theta + 2*kappa*t*theta + 2*y + 8*fastpow(rho,2)*((-3 + kappa*t)*theta + y)) + |
192 | 0 | 4*ekt*(theta + kappa*t*theta - kappa*t*y + |
193 | 0 | fastpow(rho,2)*((6 + kappa*t*(4 + kappa*t))*theta - (2 + kappa*t*(2 + kappa*t))*y)))) |
194 | 0 | )/(16.*e2kt*fastpow(-theta + kappa*t*theta + (theta - y)/ekt + y, |
195 | 0 | 4)); |
196 | 0 | } |
197 | ||
198 | FordeHestonExpansion::FordeHestonExpansion(Real kappa, Real theta, | |
199 | Real sigma, Real v0, | |
200 | 0 | Real rho, Real term) { |
201 | 0 | Real v0Sqrt = sqrt(v0); |
202 | 0 | Real rhoBarSquare = 1 - rho * rho; |
203 | 0 | Real sigma00 = v0Sqrt; |
204 | 0 | Real sigma01 = v0Sqrt * (rho * sigma / (4 * v0)); //term in x |
205 | 0 | Real sigma02 = v0Sqrt * ((1 - 5 * rho * rho / 2) / 24 * sigma * sigma/ (v0 * v0)); //term in x*x |
206 | 0 | Real a00 = -sigma * sigma / 12 * (1 - rho * rho / 4) + v0 * rho * sigma / 4 + kappa / 2 * (theta - v0); |
207 | 0 | Real a01 = rho * sigma / (24 * v0) * (sigma * sigma * rhoBarSquare - 2 * kappa * (theta + v0) + v0 * rho * sigma); //term in x |
208 | 0 | Real a02 = (176 * sigma * sigma - 480 * kappa * theta - 712 * rho * rho * sigma * sigma + 521 * rho * rho * rho * rho * sigma * sigma + 40 * sigma * rho * rho * rho * v0 + 1040 * kappa * theta * rho * rho - 80 * v0 * kappa * rho * rho) * sigma * sigma / (v0 * v0 * 7680); |
209 | 0 | coeffs[0] = sigma00*sigma00+a00*term; |
210 | 0 | coeffs[1] = sigma00*sigma01*2+a01*term; |
211 | 0 | coeffs[2] = sigma00*sigma02*2+sigma01*sigma01+a02*term; |
212 | 0 | coeffs[3] = sigma01*sigma02*2; |
213 | 0 | coeffs[4] = sigma02*sigma02; |
214 | 0 | } |
215 | ||
216 | Real FordeHestonExpansion::impliedVolatility(const Real strike, | |
217 | 0 | const Real forward) const { |
218 | 0 | Real x = log(strike/forward); |
219 | 0 | Real var = coeffs[0]+x*(coeffs[1]+(x*(coeffs[2]+x*(coeffs[3]+(x*coeffs[4]))))); |
220 | 0 | var = std::max(1e-8,var); |
221 | 0 | return sqrt(var); |
222 | 0 | } |
223 | ||
224 | Real LPP3HestonExpansion::z0(Real t, Real kappa, Real theta, | |
225 | 0 | Real delta, Real y, Real rho) const { |
226 | 0 | return (96*pow(delta,2)*ekt*pow(kappa,3)* |
227 | 0 | (-theta - 4*ekt*(theta + kappa*t*(theta - y)) + |
228 | 0 | e2kt*((5 - 2*kappa*t)*theta - 2*y) + 2*y)* |
229 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y) + |
230 | 0 | 3072*e2kt*pow(kappa,5)* |
231 | 0 | pow((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y,2) + |
232 | 0 | 96*pow(delta,3)*ekt*pow(kappa,2)*rho* |
233 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y)* |
234 | 0 | (-2*theta - kappa*t*theta - 2*ekt*(2 + kappa*t)* |
235 | 0 | (2*theta + kappa*t*(theta - y)) + e2kt*((10 - 3*kappa*t)*theta - 3*y) + |
236 | 0 | 3*y + 2*kappa*t*y) + 768*delta*e2kt*pow(kappa,4)*rho* |
237 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y)* |
238 | 0 | ((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
239 | 0 | (-1 + ekt - kappa*t)*y) + |
240 | 0 | 6*pow(delta,3)*kappa*rho*(-theta - 4*ekt*(theta + kappa*t*(theta - y)) + |
241 | 0 | e2kt*((5 - 2*kappa*t)*theta - 2*y) + 2*y)* |
242 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y)* |
243 | 0 | ((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
244 | 0 | (-1 + ekt - kappa*t)*y) + |
245 | 0 | 24*pow(delta,2)*e2kt*pow(kappa,2)*pow(rho,2)* |
246 | 0 | (-theta + kappa*t*theta + (theta - y)/ekt + y)* |
247 | 0 | pow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
248 | 0 | (-1 + ekt - kappa*t)*y,2) + |
249 | 0 | (1152*pow(delta,2)*e3kt*pow(kappa,4)*pow(rho,2)* |
250 | 0 | pow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
251 | 0 | (-1 + ekt - kappa*t)*y,2))/ |
252 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y) - |
253 | 0 | 24*pow(delta,2)*ekt*pow(kappa,2)*pow(rho,2)* |
254 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y)* |
255 | 0 | pow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
256 | 0 | (-1 + ekt - kappa*t)*y,2) + |
257 | 0 | 80*pow(delta,3)*ekt*kappa*pow(rho,3)* |
258 | 0 | pow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
259 | 0 | (-1 + ekt - kappa*t)*y,3) + |
260 | 0 | pow(delta,3)*ekt*pow(rho,3)* |
261 | 0 | (-theta + kappa*t*theta + (theta - y)/ekt + y)* |
262 | 0 | pow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
263 | 0 | (-1 + ekt - kappa*t)*y,3) - |
264 | 0 | (1440*pow(delta,3)*e3kt*pow(kappa,3)*pow(rho,3)* |
265 | 0 | pow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
266 | 0 | (-1 + ekt - kappa*t)*y,3))/ |
267 | 0 | pow((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y,2) - |
268 | 0 | (528*pow(delta,3)*e2kt*pow(kappa,2)*pow(rho,3)* |
269 | 0 | pow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
270 | 0 | (-1 + ekt - kappa*t)*y,3))/ |
271 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y) - |
272 | 0 | 3*pow(delta,3)*pow(rho,3)*((1 + ekt*(-1 + kappa*t))*theta + |
273 | 0 | (-1 + ekt)*y)*pow((2 + kappa*t + ekt*(-2 + kappa*t))* |
274 | 0 | theta + (-1 + ekt - kappa*t)*y,3) + |
275 | 0 | 384*pow(delta,3)*e2kt*pow(kappa,3)*rho* |
276 | 0 | ((2 + kappa*t + 2*ekt*pow(2 + kappa*t,2) + |
277 | 0 | e2kt*(-10 + 3*kappa*t))*theta + |
278 | 0 | (-3 + 3*e2kt - 2*kappa*t - 2*ekt*kappa*t*(2 + kappa*t))*y) - |
279 | 0 | (576*pow(delta,3)*e2kt*pow(kappa,3)*rho* |
280 | 0 | ((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
281 | 0 | (-1 + ekt - kappa*t)*y)* |
282 | 0 | ((1 + e2kt*(-5 + 2*kappa*t + 4*pow(rho,2)*(-3 + kappa*t)) + |
283 | 0 | 2*ekt*(2 + 2*kappa*t + |
284 | 0 | pow(rho,2)*(6 + 4*kappa*t + pow(kappa,2)*pow(t,2))))*theta + |
285 | 0 | 2*(-1 + e2kt*(1 + 2*pow(rho,2)) - |
286 | 0 | ekt*(2*kappa*t + |
287 | 0 | pow(rho,2)*(2 + 2*kappa*t + pow(kappa,2)*pow(t,2))))*y))/ |
288 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y) + |
289 | 0 | pow(delta,3)*rho*((1 + ekt*(-1 + kappa*t))*theta + |
290 | 0 | (-1 + ekt)*y)*((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
291 | 0 | (-1 + ekt - kappa*t)*y)* |
292 | 0 | (theta*(12*ekt*pow(kappa,3)*pow(rho,2)*pow(t,2) + |
293 | 0 | 8*pow(-1 + ekt,2)*pow(rho,2)*theta - |
294 | 0 | (-1 + ekt)*kappa* |
295 | 0 | (3 + 8*pow(rho,2)*t*theta + ekt*(15 + 8*pow(rho,2)*(9 + t*theta))) |
296 | 0 | + 2*pow(kappa,2)*t*(pow(rho,2)*t*theta + |
297 | 0 | 2*ekt*(3 + pow(rho,2)*(12 + t*theta)) + |
298 | 0 | e2kt*(3 + pow(rho,2)*(12 + t*theta)))) - |
299 | 0 | 2*(6*ekt*pow(kappa,3)*pow(rho,2)*pow(t,2) + |
300 | 0 | 4*pow(-1 + ekt,2)*pow(rho,2)*theta + |
301 | 0 | 2*pow(kappa,2)*t*(pow(rho,2)*t*theta + |
302 | 0 | ekt*(3 + pow(rho,2)*(6 + t*theta))) - |
303 | 0 | (-1 + ekt)*kappa* |
304 | 0 | (3 + 6*pow(rho,2)*t*theta + ekt*(3 + 2*pow(rho,2)*(6 + t*theta))))* |
305 | 0 | y + 2*pow(rho,2)*pow(1 - ekt + kappa*t,2)*pow(y,2)) - |
306 | 0 | (40*pow(delta,3)*kappa*rho*((1 + ekt*(-1 + kappa*t))*theta + |
307 | 0 | (-1 + ekt)*y)*((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
308 | 0 | (-1 + ekt - kappa*t)*y)* |
309 | 0 | (theta*(12*ekt*pow(kappa,3)*pow(rho,2)*pow(t,2) + |
310 | 0 | 8*pow(-1 + ekt,2)*pow(rho,2)*theta - |
311 | 0 | (-1 + ekt)*kappa* |
312 | 0 | (3 + 8*pow(rho,2)*t*theta + |
313 | 0 | ekt*(15 + 8*pow(rho,2)*(9 + t*theta))) + |
314 | 0 | 2*pow(kappa,2)*t*(pow(rho,2)*t*theta + |
315 | 0 | 2*ekt*(3 + pow(rho,2)*(12 + t*theta)) + |
316 | 0 | e2kt*(3 + pow(rho,2)*(12 + t*theta)))) - |
317 | 0 | 2*(6*ekt*pow(kappa,3)*pow(rho,2)*pow(t,2) + |
318 | 0 | 4*pow(-1 + ekt,2)*pow(rho,2)*theta + |
319 | 0 | 2*pow(kappa,2)*t*(pow(rho,2)*t*theta + |
320 | 0 | ekt*(3 + pow(rho,2)*(6 + t*theta))) - |
321 | 0 | (-1 + ekt)*kappa* |
322 | 0 | (3 + 6*pow(rho,2)*t*theta + ekt*(3 + 2*pow(rho,2)*(6 + t*theta))) |
323 | 0 | )*y + 2*pow(rho,2)*pow(1 - ekt + kappa*t,2)*pow(y,2)))/ |
324 | 0 | (-theta + kappa*t*theta + (theta - y)/ekt + y) - |
325 | 0 | 12*pow(delta,3)*kappa*rho*((1 + ekt*(-1 + kappa*t))*theta + |
326 | 0 | (-1 + ekt)*y)*(2*theta + kappa*t*theta - y - kappa*t*y + |
327 | 0 | ekt*((-2 + kappa*t)*theta + y))* |
328 | 0 | (theta - 2*y + e2kt* |
329 | 0 | (-5*theta + 2*kappa*t*theta + 2*y + 4*pow(rho,2)*((-3 + kappa*t)*theta + y)) + |
330 | 0 | 2*ekt*(2*(theta + kappa*t*(theta - y)) + |
331 | 0 | pow(rho,2)*((6 + kappa*t*(4 + kappa*t))*theta - (2 + kappa*t*(2 + kappa*t))*y))) + |
332 | 0 | (288*pow(delta,3)*pow(kappa,2)*rho* |
333 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y)* |
334 | 0 | (2*theta + kappa*t*theta - y - kappa*t*y + ekt*((-2 + kappa*t)*theta + y))* |
335 | 0 | (theta - 2*y + e2kt* |
336 | 0 | (-5*theta + 2*kappa*t*theta + 2*y + 4*pow(rho,2)*((-3 + kappa*t)*theta + y)) + |
337 | 0 | 2*ekt*(2*(theta + kappa*t*(theta - y)) + |
338 | 0 | pow(rho,2)*((6 + kappa*t*(4 + kappa*t))*theta - (2 + kappa*t*(2 + kappa*t))*y)))) |
339 | 0 | /(-theta + kappa*t*theta + (theta - y)/ekt + y) + |
340 | 0 | 48*pow(delta,2)*ekt*pow(kappa,3)* |
341 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y)* |
342 | 0 | (theta - 2*y + e2kt* |
343 | 0 | (-5*theta + 2*kappa*t*theta + 2*y + 8*pow(rho,2)*((-3 + kappa*t)*theta + y)) + |
344 | 0 | 4*ekt*(theta + kappa*t*theta - kappa*t*y + |
345 | 0 | pow(rho,2)*((6 + kappa*t*(4 + kappa*t))*theta - (2 + kappa*t*(2 + kappa*t))*y))) - |
346 | 0 | (192*pow(delta,2)*ekt*pow(kappa,4)* |
347 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y)* |
348 | 0 | (theta - 2*y + e2kt* |
349 | 0 | (-5*theta + 2*kappa*t*theta + 2*y + 8*pow(rho,2)*((-3 + kappa*t)*theta + y)) + |
350 | 0 | 4*ekt*(theta + kappa*t*theta - kappa*t*y + |
351 | 0 | pow(rho,2)*((6 + kappa*t*(4 + kappa*t))*theta - (2 + kappa*t*(2 + kappa*t))*y)))) |
352 | 0 | /(-theta + kappa*t*theta + (theta - y)/ekt + y) + |
353 | 0 | 3*pow(delta,3)*kappa*rho*((1 + ekt*(-1 + kappa*t))*theta + |
354 | 0 | (-1 + ekt)*y)*((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
355 | 0 | (-1 + ekt - kappa*t)*y)* |
356 | 0 | (theta - 2*y + e2kt* |
357 | 0 | (-5*theta + 2*kappa*t*theta + 2*y + 8*pow(rho,2)*((-3 + kappa*t)*theta + y)) + |
358 | 0 | 4*ekt*(theta + kappa*t*theta - kappa*t*y + |
359 | 0 | pow(rho,2)*((6 + kappa*t*(4 + kappa*t))*theta - (2 + kappa*t*(2 + kappa*t))*y))) - |
360 | 0 | (12*pow(delta,3)*pow(kappa,2)*rho* |
361 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y)* |
362 | 0 | ((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
363 | 0 | (-1 + ekt - kappa*t)*y)* |
364 | 0 | (theta - 2*y + e2kt* |
365 | 0 | (-5*theta + 2*kappa*t*theta + 2*y + 8*pow(rho,2)*((-3 + kappa*t)*theta + y)) + |
366 | 0 | 4*ekt*(theta + kappa*t*theta - kappa*t*y + |
367 | 0 | pow(rho,2)*((6 + kappa*t*(4 + kappa*t))*theta - (2 + kappa*t*(2 + kappa*t))*y)))) |
368 | 0 | /(-theta + kappa*t*theta + (theta - y)/ekt + y) + |
369 | 0 | 4*pow(delta,3)*kappa*rho*((1 + ekt*(-1 + kappa*t))*theta + |
370 | 0 | (-1 + ekt)*y)*(3*(theta - 2*y)*((2 + kappa*t)*theta - (1 + kappa*t)*y) + |
371 | 0 | 3*ekt*(6*pow(theta,2) + theta*y - 2*pow(y,2) + |
372 | 0 | kappa*(13*t*pow(theta,2) + theta*(8 - 18*t*y) + 4*y*(-3 + t*y)) + |
373 | 0 | 4*pow(kappa,2)*t*(theta + t*pow(theta,2) - 2*t*theta*y + y*(-2 + t*y))) + |
374 | 0 | 3*e3kt*(10*pow(theta,2) + |
375 | 0 | 2*pow(kappa,2)*t*theta*(6 + 8*pow(rho,2) + t*theta) - 9*theta*y + 2*pow(y,2) + |
376 | 0 | kappa*(-9*t*pow(theta,2) + 4*(3 + 4*pow(rho,2))*y + |
377 | 0 | theta*(-40 - 64*pow(rho,2) + 4*t*y))) + |
378 | 0 | e2kt*(-54*pow(theta,2) + |
379 | 0 | 8*pow(kappa,4)*pow(rho,2)*pow(t,3)*(theta - y) + 39*theta*y - 6*pow(y,2) + |
380 | 0 | 24*pow(kappa,3)*pow(t,2)*(theta + 2*pow(rho,2)*theta - (1 + pow(rho,2))*y) + |
381 | 0 | 6*pow(kappa,2)*t*(3*t*pow(theta,2) - 8*(1 + pow(rho,2))*y + |
382 | 0 | theta*(16 + 24*pow(rho,2) - 3*t*y)) - |
383 | 0 | 3*kappa*(5*t*pow(theta,2) + 2*y*(8*pow(rho,2) + 3*t*y) - |
384 | 0 | theta*(32 + 64*pow(rho,2) + 17*t*y)))) - |
385 | 0 | (48*pow(delta,3)*pow(kappa,2)*rho* |
386 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y)* |
387 | 0 | (3*(theta - 2*y)*((2 + kappa*t)*theta - (1 + kappa*t)*y) + |
388 | 0 | 3*ekt*(6*pow(theta,2) + theta*y - 2*pow(y,2) + |
389 | 0 | kappa*(13*t*pow(theta,2) + theta*(8 - 18*t*y) + 4*y*(-3 + t*y)) + |
390 | 0 | 4*pow(kappa,2)*t*(theta + t*pow(theta,2) - 2*t*theta*y + y*(-2 + t*y))) + |
391 | 0 | 3*e3kt*(10*pow(theta,2) + |
392 | 0 | 2*pow(kappa,2)*t*theta*(6 + 8*pow(rho,2) + t*theta) - 9*theta*y + |
393 | 0 | 2*pow(y,2) + kappa*(-9*t*pow(theta,2) + 4*(3 + 4*pow(rho,2))*y + |
394 | 0 | theta*(-40 - 64*pow(rho,2) + 4*t*y))) + |
395 | 0 | e2kt*(-54*pow(theta,2) + |
396 | 0 | 8*pow(kappa,4)*pow(rho,2)*pow(t,3)*(theta - y) + 39*theta*y - 6*pow(y,2) + |
397 | 0 | 24*pow(kappa,3)*pow(t,2)* |
398 | 0 | (theta + 2*pow(rho,2)*theta - (1 + pow(rho,2))*y) + |
399 | 0 | 6*pow(kappa,2)*t*(3*t*pow(theta,2) - 8*(1 + pow(rho,2))*y + |
400 | 0 | theta*(16 + 24*pow(rho,2) - 3*t*y)) - |
401 | 0 | 3*kappa*(5*t*pow(theta,2) + 2*y*(8*pow(rho,2) + 3*t*y) - |
402 | 0 | theta*(32 + 64*pow(rho,2) + 17*t*y)))))/ |
403 | 0 | (-theta + kappa*t*theta + (theta - y)/ekt + y) + |
404 | 0 | (240*pow(delta,3)*e2kt*pow(kappa,2)*rho* |
405 | 0 | ((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
406 | 0 | (-1 + ekt - kappa*t)*y)* |
407 | 0 | (12*ekt*pow(kappa,3)*pow(rho,2)*pow(t,2)*(theta - y) + |
408 | 0 | 2*pow(-1 + ekt,2)*pow(rho,2)*pow(-2*theta + y,2) - |
409 | 0 | (-1 + ekt)*kappa* |
410 | 0 | (8*(1 + ekt)*pow(rho,2)*t*pow(theta,2) + |
411 | 0 | 2*y*(-3 - 3*ekt*(1 + 4*pow(rho,2)) + 2*pow(rho,2)*t*y) + |
412 | 0 | theta*(3 - 12*pow(rho,2)*t*y + ekt*(15 + pow(rho,2)*(72 - 4*t*y))) |
413 | 0 | ) + 2*pow(kappa,2)*t*(e2kt*theta* |
414 | 0 | (3 + pow(rho,2)*(12 + t*theta)) + pow(rho,2)*t*pow(theta - y,2) + |
415 | 0 | 2*ekt*(pow(rho,2)*t*pow(theta,2) - 3*(y + 2*pow(rho,2)*y) + |
416 | 0 | theta*(3 + pow(rho,2)*(12 - t*y))))))/ |
417 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y))/ |
418 | 0 | (3072.*e4kt*pow(kappa,7)*pow(t,2)* |
419 | 0 | pow((-theta + kappa*t*theta + (theta - y)/ekt + y)/(kappa*t),1.5)); |
420 | 0 | } |
421 | ||
422 | Real LPP3HestonExpansion::z1(Real t, Real kappa, Real theta, | |
423 | 0 | Real delta, Real y, Real rho) const { |
424 | 0 | return (delta*(768*e2kt*pow(kappa,4)*rho* |
425 | 0 | ((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
426 | 0 | (-1 + ekt - kappa*t)*y) - |
427 | 0 | (576*delta*e2kt*pow(kappa,3)*pow(rho,2)* |
428 | 0 | pow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
429 | 0 | (-1 + ekt - kappa*t)*y,2))/ |
430 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y) - |
431 | 0 | 10*pow(delta,2)*pow(rho,3)*pow((2 + kappa*t + ekt*(-2 + kappa*t))* |
432 | 0 | theta + (-1 + ekt - kappa*t)*y,3) + |
433 | 0 | (6*pow(delta,2)*kappa*pow(rho,3)* |
434 | 0 | pow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
435 | 0 | (-1 + ekt - kappa*t)*y,3))/ |
436 | 0 | (-theta + kappa*t*theta + (theta - y)/ekt + y) - |
437 | 0 | (3360*pow(delta,2)*e3kt*pow(kappa,3)*pow(rho,3)* |
438 | 0 | pow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
439 | 0 | (-1 + ekt - kappa*t)*y,3))/ |
440 | 0 | pow((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y,3) - |
441 | 0 | (288*pow(delta,2)*e2kt*pow(kappa,2)*pow(rho,3)* |
442 | 0 | pow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
443 | 0 | (-1 + ekt - kappa*t)*y,3))/ |
444 | 0 | pow((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y,2) + |
445 | 0 | (234*pow(delta,2)*ekt*kappa*pow(rho,3)* |
446 | 0 | pow((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
447 | 0 | (-1 + ekt - kappa*t)*y,3))/ |
448 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y) - |
449 | 0 | 96*delta*ekt*pow(kappa,3)* |
450 | 0 | ((1 + 4*ekt*(1 + kappa*t) + e2kt*(-5 + 2*kappa*t))*theta + |
451 | 0 | 2*(-1 + e2kt - 2*ekt*kappa*t)*y) - |
452 | 0 | 12*pow(delta,2)*kappa*rho*((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
453 | 0 | (-1 + ekt - kappa*t)*y)* |
454 | 0 | ((1 + 4*ekt*(1 + kappa*t) + e2kt*(-5 + 2*kappa*t))*theta + |
455 | 0 | 2*(-1 + e2kt - 2*ekt*kappa*t)*y) - |
456 | 0 | 192*pow(delta,2)*ekt*pow(kappa,2)*rho* |
457 | 0 | ((2 + kappa*t + 2*ekt*pow(2 + kappa*t,2) + |
458 | 0 | e2kt*(-10 + 3*kappa*t))*theta + |
459 | 0 | (-3 + 3*e2kt - 2*kappa*t - 2*ekt*kappa*t*(2 + kappa*t))*y) |
460 | 0 | - (12*pow(delta,2)*ekt*pow(kappa,2)*rho* |
461 | 0 | ((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
462 | 0 | (-1 + ekt - kappa*t)*y)* |
463 | 0 | ((1 + e2kt*(-5 + 2*kappa*t + 8*pow(rho,2)*(-3 + kappa*t)) + |
464 | 0 | 4*ekt*(1 + kappa*t + |
465 | 0 | pow(rho,2)*(6 + 4*kappa*t + pow(kappa,2)*pow(t,2))))*theta + |
466 | 0 | 2*(-1 + e2kt*(1 + 4*pow(rho,2)) - |
467 | 0 | 2*ekt*(kappa*t + |
468 | 0 | pow(rho,2)*(2 + 2*kappa*t + pow(kappa,2)*pow(t,2))))*y))/ |
469 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y) + |
470 | 0 | (576*pow(delta,2)*ekt*pow(kappa,2)*rho* |
471 | 0 | ((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
472 | 0 | (-1 + ekt - kappa*t)*y)* |
473 | 0 | ((1 + e2kt*(-5 + 2*kappa*t + 4*pow(rho,2)*(-3 + kappa*t)) + |
474 | 0 | 2*ekt*(2 + 2*kappa*t + |
475 | 0 | pow(rho,2)*(6 + 4*kappa*t + pow(kappa,2)*pow(t,2))))*theta + |
476 | 0 | 2*(-1 + e2kt*(1 + 2*pow(rho,2)) - |
477 | 0 | ekt*(2*kappa*t + |
478 | 0 | pow(rho,2)*(2 + 2*kappa*t + pow(kappa,2)*pow(t,2))))*y))/ |
479 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y) + |
480 | 0 | (5*pow(delta,2)*rho*((1 + ekt*(-1 + kappa*t))*theta + |
481 | 0 | (-1 + ekt)*y)* |
482 | 0 | ((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
483 | 0 | (-1 + ekt - kappa*t)*y)* |
484 | 0 | (theta*(12*ekt*pow(kappa,3)*pow(rho,2)*pow(t,2) + |
485 | 0 | 8*pow(-1 + ekt,2)*pow(rho,2)*theta - |
486 | 0 | (-1 + ekt)*kappa* |
487 | 0 | (3 + 8*pow(rho,2)*t*theta + |
488 | 0 | ekt*(15 + 8*pow(rho,2)*(9 + t*theta))) + |
489 | 0 | 2*pow(kappa,2)*t*(pow(rho,2)*t*theta + |
490 | 0 | 2*ekt*(3 + pow(rho,2)*(12 + t*theta)) + |
491 | 0 | e2kt*(3 + pow(rho,2)*(12 + t*theta)))) - |
492 | 0 | 2*(6*ekt*pow(kappa,3)*pow(rho,2)*pow(t,2) + |
493 | 0 | 4*pow(-1 + ekt,2)*pow(rho,2)*theta + |
494 | 0 | 2*pow(kappa,2)*t*(pow(rho,2)*t*theta + |
495 | 0 | ekt*(3 + pow(rho,2)*(6 + t*theta))) - |
496 | 0 | (-1 + ekt)*kappa* |
497 | 0 | (3 + 6*pow(rho,2)*t*theta + |
498 | 0 | ekt*(3 + 2*pow(rho,2)*(6 + t*theta))))*y + |
499 | 0 | 2*pow(rho,2)*pow(1 - ekt + kappa*t,2)*pow(y,2)))/ |
500 | 0 | (ekt*(-theta + kappa*t*theta + (theta - y)/ekt + y)) - |
501 | 0 | (48*pow(delta,2)*kappa*rho*((1 + ekt*(-1 + kappa*t))*theta + |
502 | 0 | (-1 + ekt)*y)* |
503 | 0 | (2*theta + kappa*t*theta - y - kappa*t*y + |
504 | 0 | ekt*((-2 + kappa*t)*theta + y))* |
505 | 0 | (theta - 2*y + e2kt* |
506 | 0 | (-5*theta + 2*kappa*t*theta + 2*y + 4*pow(rho,2)*((-3 + kappa*t)*theta + y)) + |
507 | 0 | 2*ekt*(2*(theta + kappa*t*(theta - y)) + |
508 | 0 | pow(rho,2)*((6 + kappa*t*(4 + kappa*t))*theta - (2 + kappa*t*(2 + kappa*t))*y)) |
509 | 0 | ))/(ekt*(-theta + kappa*t*theta + (theta - y)/ekt + y)) + |
510 | 0 | (96*delta*pow(kappa,3)*((1 + ekt*(-1 + kappa*t))*theta + |
511 | 0 | (-1 + ekt)*y)* |
512 | 0 | (theta - 2*y + e2kt* |
513 | 0 | (-5*theta + 2*kappa*t*theta + 2*y + 8*pow(rho,2)*((-3 + kappa*t)*theta + y)) + |
514 | 0 | 4*ekt*(theta + kappa*t*theta - kappa*t*y + |
515 | 0 | pow(rho,2)*((6 + kappa*t*(4 + kappa*t))*theta - (2 + kappa*t*(2 + kappa*t))*y)) |
516 | 0 | ))/(-theta + kappa*t*theta + (theta - y)/ekt + y) + |
517 | 0 | (9*pow(delta,2)*kappa*rho*((1 + ekt*(-1 + kappa*t))*theta + |
518 | 0 | (-1 + ekt)*y)* |
519 | 0 | ((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
520 | 0 | (-1 + ekt - kappa*t)*y)* |
521 | 0 | (theta - 2*y + e2kt* |
522 | 0 | (-5*theta + 2*kappa*t*theta + 2*y + 8*pow(rho,2)*((-3 + kappa*t)*theta + y)) + |
523 | 0 | 4*ekt*(theta + kappa*t*theta - kappa*t*y + |
524 | 0 | pow(rho,2)*((6 + kappa*t*(4 + kappa*t))*theta - (2 + kappa*t*(2 + kappa*t))*y)) |
525 | 0 | ))/(ekt*(-theta + kappa*t*theta + (theta - y)/ekt + y)) - |
526 | 0 | (48*pow(delta,2)*ekt*pow(kappa,2)*rho* |
527 | 0 | (3*(theta - 2*y)*((2 + kappa*t)*theta - (1 + kappa*t)*y) + |
528 | 0 | 3*ekt*(6*pow(theta,2) + theta*y - 2*pow(y,2) + |
529 | 0 | kappa*(13*t*pow(theta,2) + theta*(8 - 18*t*y) + 4*y*(-3 + t*y)) + |
530 | 0 | 4*pow(kappa,2)*t*(theta + t*pow(theta,2) - 2*t*theta*y + y*(-2 + t*y))) + |
531 | 0 | 3*e3kt*(10*pow(theta,2) + |
532 | 0 | 2*pow(kappa,2)*t*theta*(6 + 8*pow(rho,2) + t*theta) - 9*theta*y + |
533 | 0 | 2*pow(y,2) + kappa*(-9*t*pow(theta,2) + 4*(3 + 4*pow(rho,2))*y + |
534 | 0 | theta*(-40 - 64*pow(rho,2) + 4*t*y))) + |
535 | 0 | e2kt*(-54*pow(theta,2) + |
536 | 0 | 8*pow(kappa,4)*pow(rho,2)*pow(t,3)*(theta - y) + 39*theta*y - |
537 | 0 | 6*pow(y,2) + 24*pow(kappa,3)*pow(t,2)* |
538 | 0 | (theta + 2*pow(rho,2)*theta - (1 + pow(rho,2))*y) + |
539 | 0 | 6*pow(kappa,2)*t*(3*t*pow(theta,2) - 8*(1 + pow(rho,2))*y + |
540 | 0 | theta*(16 + 24*pow(rho,2) - 3*t*y)) - |
541 | 0 | 3*kappa*(5*t*pow(theta,2) + 2*y*(8*pow(rho,2) + 3*t*y) - |
542 | 0 | theta*(32 + 64*pow(rho,2) + 17*t*y)))))/ |
543 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y) + |
544 | 0 | (12*pow(delta,2)*kappa*rho*((1 + ekt*(-1 + kappa*t))*theta + |
545 | 0 | (-1 + ekt)*y)* |
546 | 0 | (3*(theta - 2*y)*((2 + kappa*t)*theta - (1 + kappa*t)*y) + |
547 | 0 | 3*ekt*(6*pow(theta,2) + theta*y - 2*pow(y,2) + |
548 | 0 | kappa*(13*t*pow(theta,2) + theta*(8 - 18*t*y) + 4*y*(-3 + t*y)) + |
549 | 0 | 4*pow(kappa,2)*t*(theta + t*pow(theta,2) - 2*t*theta*y + y*(-2 + t*y))) + |
550 | 0 | 3*e3kt*(10*pow(theta,2) + |
551 | 0 | 2*pow(kappa,2)*t*theta*(6 + 8*pow(rho,2) + t*theta) - 9*theta*y + |
552 | 0 | 2*pow(y,2) + kappa*(-9*t*pow(theta,2) + 4*(3 + 4*pow(rho,2))*y + |
553 | 0 | theta*(-40 - 64*pow(rho,2) + 4*t*y))) + |
554 | 0 | e2kt*(-54*pow(theta,2) + |
555 | 0 | 8*pow(kappa,4)*pow(rho,2)*pow(t,3)*(theta - y) + 39*theta*y - |
556 | 0 | 6*pow(y,2) + 24*pow(kappa,3)*pow(t,2)* |
557 | 0 | (theta + 2*pow(rho,2)*theta - (1 + pow(rho,2))*y) + |
558 | 0 | 6*pow(kappa,2)*t*(3*t*pow(theta,2) - 8*(1 + pow(rho,2))*y + |
559 | 0 | theta*(16 + 24*pow(rho,2) - 3*t*y)) - |
560 | 0 | 3*kappa*(5*t*pow(theta,2) + 2*y*(8*pow(rho,2) + 3*t*y) - |
561 | 0 | theta*(32 + 64*pow(rho,2) + 17*t*y)))))/ |
562 | 0 | (ekt*(-theta + kappa*t*theta + (theta - y)/ekt + y)) + |
563 | 0 | (240*pow(delta,2)*e2kt*pow(kappa,2)*rho* |
564 | 0 | ((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
565 | 0 | (-1 + ekt - kappa*t)*y)* |
566 | 0 | (12*ekt*pow(kappa,3)*pow(rho,2)*pow(t,2)*(theta - y) + |
567 | 0 | 2*pow(-1 + ekt,2)*pow(rho,2)*pow(-2*theta + y,2) - |
568 | 0 | (-1 + ekt)*kappa* |
569 | 0 | (8*(1 + ekt)*pow(rho,2)*t*pow(theta,2) + |
570 | 0 | 2*y*(-3 - 3*ekt*(1 + 4*pow(rho,2)) + 2*pow(rho,2)*t*y) + |
571 | 0 | theta*(3 - 12*pow(rho,2)*t*y + |
572 | 0 | ekt*(15 + pow(rho,2)*(72 - 4*t*y)))) + |
573 | 0 | 2*pow(kappa,2)*t*(e2kt*theta*(3 + pow(rho,2)*(12 + t*theta)) + |
574 | 0 | pow(rho,2)*t*pow(theta - y,2) + |
575 | 0 | 2*ekt*(pow(rho,2)*t*pow(theta,2) - 3*(y + 2*pow(rho,2)*y) + |
576 | 0 | theta*(3 + pow(rho,2)*(12 - t*y))))))/ |
577 | 0 | pow((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y,2) - |
578 | 0 | (120*pow(delta,2)*ekt*kappa*rho* |
579 | 0 | ((2 + kappa*t + ekt*(-2 + kappa*t))*theta + |
580 | 0 | (-1 + ekt - kappa*t)*y)* |
581 | 0 | (12*ekt*pow(kappa,3)*pow(rho,2)*pow(t,2)*(theta - y) + |
582 | 0 | 2*pow(-1 + ekt,2)*pow(rho,2)*pow(-2*theta + y,2) - |
583 | 0 | (-1 + ekt)*kappa* |
584 | 0 | (8*(1 + ekt)*pow(rho,2)*t*pow(theta,2) + |
585 | 0 | 2*y*(-3 - 3*ekt*(1 + 4*pow(rho,2)) + 2*pow(rho,2)*t*y) + |
586 | 0 | theta*(3 - 12*pow(rho,2)*t*y + |
587 | 0 | ekt*(15 + pow(rho,2)*(72 - 4*t*y)))) + |
588 | 0 | 2*pow(kappa,2)*t*(e2kt*theta*(3 + pow(rho,2)*(12 + t*theta)) + |
589 | 0 | pow(rho,2)*t*pow(theta - y,2) + |
590 | 0 | 2*ekt*(pow(rho,2)*t*pow(theta,2) - 3*(y + 2*pow(rho,2)*y) + |
591 | 0 | theta*(3 + pow(rho,2)*(12 - t*y))))))/ |
592 | 0 | ((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y)))/ |
593 | 0 | (1536.*e3kt*pow(kappa,6)*pow(t,2)* |
594 | 0 | pow((-theta + kappa*t*theta + (theta - y)/ekt + y)/(kappa*t),1.5)); |
595 | 0 | } |
596 | ||
597 | Real LPP3HestonExpansion::z2(Real t, Real kappa, Real theta, | |
598 | 0 | Real delta, Real y, Real rho) const{ |
599 | 0 | return (pow(delta,2)*(8*e3kt*pow(kappa,5)*pow(rho,2)*pow(t,4)*(2 + delta*rho*t)* |
600 | 0 | pow(theta,2)*(theta - y) - delta*pow(-1 + ekt,3)*rho* |
601 | 0 | (2*(-1 + ekt*(-5 + 24*pow(rho,2)))*pow(theta,3) + |
602 | 0 | (7 + ekt*(3 + 56*pow(rho,2)))*pow(theta,2)*y - |
603 | 0 | 3*(1 + ekt*(-3 + 8*pow(rho,2)))*theta*pow(y,2) + |
604 | 0 | 2*(-1 + ekt*(-1 + 2*pow(rho,2)))*pow(y,3)) - |
605 | 0 | pow(-1 + ekt,2)*kappa* |
606 | 0 | ((-4 + delta*rho*t - 8*ekt* |
607 | 0 | (2 - 12*pow(rho,2) - 4*delta*rho*t + 25*delta*pow(rho,3)*t) + |
608 | 0 | e2kt*(20 - 96*pow(rho,2) + 3*delta*rho*t + 56*delta*pow(rho,3)*t) |
609 | 0 | )*pow(theta,3) - 2*(-8 + 2*delta*rho*t + |
610 | 0 | e2kt*(24 - 80*pow(rho,2) - 9*delta*rho*t + |
611 | 0 | 24*delta*pow(rho,3)*t) - |
612 | 0 | 4*ekt*(4 - 20*pow(rho,2) - 10*delta*rho*t + 39*delta*pow(rho,3)*t) |
613 | 0 | )*pow(theta,2)*y + (5*(-4 + delta*rho*t) + |
614 | 0 | ekt*(-16 + 80*pow(rho,2) + 57*delta*rho*t - |
615 | 0 | 140*delta*pow(rho,3)*t) + |
616 | 0 | 2*e2kt*(18 - 40*pow(rho,2) - 3*delta*rho*t + |
617 | 0 | 6*delta*pow(rho,3)*t))*theta*pow(y,2) + |
618 | 0 | 2*(4 + e2kt*(-4 + 8*pow(rho,2)) - delta*rho*t + |
619 | 0 | ekt*rho*(-8*rho - 7*delta*t + 14*delta*pow(rho,2)*t))*pow(y,3)) + |
620 | 0 | ekt*(-1 + ekt)*pow(kappa,2)*t* |
621 | 0 | ((-24 + 128*pow(rho,2) + 9*delta*rho*t - 144*delta*pow(rho,3)*t - |
622 | 0 | 4*ekt*(6 - 8*pow(rho,2) - 9*delta*rho*t + 6*delta*pow(rho,3)*t) + |
623 | 0 | e2kt*(48 - 160*pow(rho,2) - 9*delta*rho*t + |
624 | 0 | 24*delta*pow(rho,3)*t))*pow(theta,3) - |
625 | 0 | (-72 + 320*pow(rho,2) + 27*delta*rho*t - 360*delta*pow(rho,3)*t - |
626 | 0 | ekt*rho*(160*rho - 81*delta*t + 348*delta*pow(rho,2)*t) + |
627 | 0 | 2*e2kt*(36 - 80*pow(rho,2) - 3*delta*rho*t + |
628 | 0 | 6*delta*pow(rho,3)*t))*pow(theta,2)*y - |
629 | 0 | 2*(32 - 128*pow(rho,2) + 12*e2kt*(-1 + 2*pow(rho,2)) - |
630 | 0 | 15*delta*rho*t + 144*delta*pow(rho,3)*t + |
631 | 0 | 2*ekt*(-10 + 52*pow(rho,2) - 13*delta*rho*t + |
632 | 0 | 58*delta*pow(rho,3)*t))*theta*pow(y,2) + |
633 | 0 | 4*(4 - 16*pow(rho,2) - 3*delta*rho*t + 18*delta*pow(rho,3)*t + |
634 | 0 | ekt*(-4 + 16*pow(rho,2) - 2*delta*rho*t + 11*delta*pow(rho,3)*t))* |
635 | 0 | pow(y,3)) - 4*e2kt*pow(kappa,4)*pow(t,3)*theta* |
636 | 0 | (2*e2kt*(-1 + 2*pow(rho,2))*pow(theta,2) + |
637 | 0 | pow(rho,2)*(4 + 13*delta*rho*t)*pow(theta - y,2) + |
638 | 0 | ekt*((-4 + 16*pow(rho,2) - 2*delta*rho*t + 9*delta*pow(rho,3)*t)* |
639 | 0 | pow(theta,2) + (4 - 32*pow(rho,2) + 2*delta*rho*t - 19*delta*pow(rho,3)*t)* |
640 | 0 | theta*y + 4*pow(rho,2)*(2 + delta*rho*t)*pow(y,2))) - |
641 | 0 | 2*ekt*pow(kappa,3)*pow(t,2)* |
642 | 0 | (-4*pow(rho,2)*(-4 + 3*delta*rho*t)*pow(theta - y,3) + |
643 | 0 | e3kt*pow(theta,2)* |
644 | 0 | ((18 - 40*pow(rho,2) - delta*rho*t + 2*delta*pow(rho,3)*t)*theta + |
645 | 0 | 12*(-1 + 2*pow(rho,2))*y) + |
646 | 0 | 2*ekt*((-9 + 36*pow(rho,2) + 19*delta*pow(rho,3)*t)*pow(theta,3) + |
647 | 0 | 2*(9 - 30*pow(rho,2) + 7*delta*pow(rho,3)*t)*pow(theta,2)*y + |
648 | 0 | (-8 + 20*pow(rho,2) + delta*rho*t - 46*delta*pow(rho,3)*t)*theta*pow(y,2) + |
649 | 0 | pow(rho,2)*(4 + 13*delta*rho*t)*pow(y,3)) + |
650 | 0 | e2kt*(8*theta*y*(-3*theta + 2*y) + |
651 | 0 | delta*rho*t*theta*(7*pow(theta,2) - 23*theta*y + 8*pow(y,2)) - |
652 | 0 | 8*pow(rho,2)*(6*pow(theta,3) - 18*pow(theta,2)*y + 11*theta*pow(y,2) - |
653 | 0 | pow(y,3)) + 4*delta*pow(rho,3)*t* |
654 | 0 | (-13*pow(theta,3) + 31*pow(theta,2)*y - 14*theta*pow(y,2) + pow(y,3))))))/ |
655 | 0 | (64.*pow(kappa,2)*t*sqrt((-theta + kappa*t*theta + (theta - y)/ekt + y)/ |
656 | 0 | (kappa*t))*pow((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y, |
657 | 0 | 4)); |
658 | 0 | } |
659 | ||
660 | Real LPP3HestonExpansion::z3(Real t, Real kappa, Real theta, | |
661 | 0 | Real delta, Real y, Real rho) const { |
662 | 0 | return (pow(delta,3)*ekt*rho*((-15*(2 + kappa*t) + |
663 | 0 | 3*e4kt*(50 - 79*kappa*t + 35*pow(kappa,2)*pow(t,2) - |
664 | 0 | 6*pow(kappa,3)*pow(t,3) + |
665 | 0 | 8*pow(rho,2)*(-18 + 15*kappa*t - 6*pow(kappa,2)*pow(t,2) + |
666 | 0 | pow(kappa,3)*pow(t,3))) + |
667 | 0 | ekt*(-3*(20 + 86*kappa*t + 29*pow(kappa,2)*pow(t,2)) + |
668 | 0 | pow(rho,2)*(432 + 936*kappa*t + 552*pow(kappa,2)*pow(t,2) + |
669 | 0 | 92*pow(kappa,3)*pow(t,3))) + |
670 | 0 | e2kt*(360 + 324*kappa*t - 261*pow(kappa,2)*pow(t,2) - |
671 | 0 | 48*pow(kappa,3)*pow(t,3) - |
672 | 0 | 4*pow(rho,2)*(324 + 378*kappa*t - 12*pow(kappa,2)*pow(t,2) - |
673 | 0 | 2*pow(kappa,3)*pow(t,3) + 23*pow(kappa,4)*pow(t,4))) + |
674 | 0 | e3kt*(3*(-140 + 62*kappa*t + 81*pow(kappa,2)*pow(t,2) - |
675 | 0 | 38*pow(kappa,3)*pow(t,3) + 8*pow(kappa,4)*pow(t,4)) + |
676 | 0 | 4*pow(rho,2)*(324 + 54*kappa*t - 114*pow(kappa,2)*pow(t,2) + |
677 | 0 | 77*pow(kappa,3)*pow(t,3) - 19*pow(kappa,4)*pow(t,4) + |
678 | 0 | 2*pow(kappa,5)*pow(t,5))))*pow(theta,3) + |
679 | 0 | (15*(7 + 4*kappa*t) + 3*e4kt* |
680 | 0 | (-79 + 70*kappa*t - 18*pow(kappa,2)*pow(t,2) + |
681 | 0 | 24*pow(rho,2)*(5 - 4*kappa*t + pow(kappa,2)*pow(t,2))) - |
682 | 0 | 3*ekt*(26 - 200*kappa*t - 87*pow(kappa,2)*pow(t,2) + |
683 | 0 | 4*pow(rho,2)*(30 + 142*kappa*t + 115*pow(kappa,2)*pow(t,2) + |
684 | 0 | 23*pow(kappa,3)*pow(t,3))) + |
685 | 0 | 2*e2kt*(3*(-66 - 195*kappa*t + 63*pow(kappa,2)*pow(t,2) + |
686 | 0 | 16*pow(kappa,3)*pow(t,3)) + |
687 | 0 | 4*pow(rho,2)*(135 + 390*kappa*t - 9*pow(kappa,2)*pow(t,2) - |
688 | 0 | 48*pow(kappa,3)*pow(t,3) + 23*pow(kappa,4)*pow(t,4))) + |
689 | 0 | e3kt*(606 + 300*kappa*t - 585*pow(kappa,2)*pow(t,2) + |
690 | 0 | 210*pow(kappa,3)*pow(t,3) - 24*pow(kappa,4)*pow(t,4) - |
691 | 0 | 4*pow(rho,2)*(270 + 282*kappa*t - 345*pow(kappa,2)*pow(t,2) + |
692 | 0 | 153*pow(kappa,3)*pow(t,3) - 29*pow(kappa,4)*pow(t,4) + |
693 | 0 | 2*pow(kappa,5)*pow(t,5))))*pow(theta,2)*y + |
694 | 0 | (-93 - 75*kappa*t + 3*e4kt* |
695 | 0 | (35 - 18*kappa*t + 24*pow(rho,2)*(-2 + kappa*t)) + |
696 | 0 | 3*ekt*(58 - 123*kappa*t - 86*pow(kappa,2)*pow(t,2) + |
697 | 0 | 4*pow(rho,2)*(12 + 80*kappa*t + 92*pow(kappa,2)*pow(t,2) + |
698 | 0 | 23*pow(kappa,3)*pow(t,3))) + |
699 | 0 | e3kt*(-3*(74 + 137*kappa*t - 100*pow(kappa,2)*pow(t,2) + |
700 | 0 | 16*pow(kappa,3)*pow(t,3)) - |
701 | 0 | 16*pow(rho,2)*(-27 - 51*kappa*t + 45*pow(kappa,2)*pow(t,2) - |
702 | 0 | 12*pow(kappa,3)*pow(t,3) + pow(kappa,4)*pow(t,4))) + |
703 | 0 | e2kt*(36 + 909*kappa*t - 42*pow(kappa,2)*pow(t,2) - |
704 | 0 | 60*pow(kappa,3)*pow(t,3) - |
705 | 0 | 4*pow(rho,2)*(108 + 462*kappa*t + 96*pow(kappa,2)*pow(t,2) - |
706 | 0 | 117*pow(kappa,3)*pow(t,3) + 23*pow(kappa,4)*pow(t,4))))*theta*pow(y,2) |
707 | 0 | + 2*(9 + 3*e4kt*(-3 + 4*pow(rho,2)) + 15*kappa*t + |
708 | 0 | e2kt*(-3*kappa*t*(33 + 10*kappa*t) + |
709 | 0 | pow(rho,2)*(36 + 192*kappa*t + 96*pow(kappa,2)*pow(t,2) - |
710 | 0 | 46*pow(kappa,3)*pow(t,3))) + |
711 | 0 | e3kt*(18 + 57*kappa*t - 12*pow(kappa,2)*pow(t,2) - |
712 | 0 | 2*pow(rho,2)*(18 + 48*kappa*t - 21*pow(kappa,2)*pow(t,2) + |
713 | 0 | 2*pow(kappa,3)*pow(t,3))) + |
714 | 0 | ekt*(3*(-6 + 9*kappa*t + 14*pow(kappa,2)*pow(t,2)) - |
715 | 0 | 2*pow(rho,2)*(6 + 48*kappa*t + 69*pow(kappa,2)*pow(t,2) + |
716 | 0 | 23*pow(kappa,3)*pow(t,3))))*pow(y,3)))/ |
717 | 0 | (96.*kappa*t*sqrt((-theta + kappa*t*theta + (theta - y)/ekt + y)/(kappa*t))* |
718 | 0 | pow((1 + ekt*(-1 + kappa*t))*theta + (-1 + ekt)*y,5)); |
719 | 0 | } |
720 | ||
721 | LPP3HestonExpansion::LPP3HestonExpansion(Real kappa, Real theta, Real sigma, | |
722 | 0 | Real v0, Real rho, Real term) { |
723 | 0 | ekt = exp(kappa*term); |
724 | 0 | e2kt = ekt*ekt; |
725 | 0 | e3kt = e2kt*ekt; |
726 | 0 | e4kt = e2kt*e2kt; |
727 | 0 | coeffs[0] = z0(term, kappa, theta, sigma, v0, rho); |
728 | 0 | coeffs[1] = z1(term, kappa, theta, sigma, v0, rho); |
729 | 0 | coeffs[2] = z2(term, kappa, theta, sigma, v0, rho); |
730 | 0 | coeffs[3] = z3(term, kappa, theta, sigma, v0, rho); |
731 | 0 | } |
732 | ||
733 | Real LPP3HestonExpansion::impliedVolatility(const Real strike, | |
734 | 0 | const Real forward) const { |
735 | 0 | Real x = log(strike/forward); |
736 | 0 | Real vol = coeffs[0]+x*(coeffs[1]+x*(coeffs[2]+x*(coeffs[3]))); |
737 | 0 | return std::max(1e-8,vol); |
738 | 0 | } |
739 | } |