QuantLib: a free/open-source library for quantitative finance
Fully annotated sources - version 1.32
Loading...
Searching...
No Matches
gaussian1dnonstandardswaptionengine.cpp
1/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2
3/*
4 Copyright (C) 2013 Peter Caspers
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#include <ql/pricingengines/swaption/gaussian1dnonstandardswaptionengine.hpp>
21#include <ql/rebatedexercise.hpp>
22#include <ql/time/daycounters/actualactual.hpp>
23#include <ql/quotes/simplequote.hpp>
24#include <ql/math/interpolations/cubicinterpolation.hpp>
25#include <ql/payoff.hpp>
26
27using std::exp;
28
29namespace QuantLib {
30
31 Real
32 Gaussian1dNonstandardSwaptionEngine::underlyingNpv(const Date &expiry,
33 const Real y) const {
34
35 // determine the indices on both legs representing the cashflows that
36 // are part of the exercise into right
37
38 Size fixedIdx =
39 std::upper_bound(arguments_.fixedResetDates.begin(),
40 arguments_.fixedResetDates.end(), expiry - 1) -
41 arguments_.fixedResetDates.begin();
42 Size floatingIdx =
43 std::upper_bound(arguments_.floatingResetDates.begin(),
44 arguments_.floatingResetDates.end(), expiry - 1) -
45 arguments_.floatingResetDates.begin();
46
47 // calculate the npv of these cashflows conditional on y at expiry
48
49 Real type = (Real)arguments_.type;
50
51 Real npv = 0.0;
52 for (Size i = fixedIdx; i < arguments_.fixedResetDates.size(); i++) {
53 npv -=
54 arguments_.fixedCoupons[i] *
55 model_->zerobond(arguments_.fixedPayDates[i], expiry, y,
56 discountCurve_) *
57 (oas_.empty()
58 ? Real(1.0)
59 : exp(-oas_->value() *
60 model_->termStructure()->dayCounter().yearFraction(
61 expiry, arguments_.fixedPayDates[i])));
62 }
63
64 for (Size i = floatingIdx; i < arguments_.floatingResetDates.size();
65 i++) {
66 Real amount;
67 if (!arguments_.floatingIsRedemptionFlow[i])
68 amount = (arguments_.floatingGearings[i] *
69 model_->forwardRate(
70 arguments_.floatingFixingDates[i], expiry, y,
71 arguments_.swap->iborIndex()) +
72 arguments_.floatingSpreads[i]) *
73 arguments_.floatingAccrualTimes[i] *
74 arguments_.floatingNominal[i];
75 else
76 amount = arguments_.floatingCoupons[i];
77 npv +=
78 amount * model_->zerobond(arguments_.floatingPayDates[i],
79 expiry, y, discountCurve_) *
80 (oas_.empty()
81 ? Real(1.0)
82 : exp(-oas_->value() *
83 model_->termStructure()->dayCounter().yearFraction(
84 expiry, arguments_.floatingPayDates[i])));
85 }
86
87 return type * npv;
88 }
89
90 Swap::Type Gaussian1dNonstandardSwaptionEngine::underlyingType() const {
91 return arguments_.swap->type();
92 }
93
94 // NOLINTNEXTLINE(readability-const-return-type)
95 const Date Gaussian1dNonstandardSwaptionEngine::underlyingLastDate() const {
96 return arguments_.fixedPayDates.back();
97 }
98
99 // NOLINTNEXTLINE(readability-const-return-type)
100 const Array Gaussian1dNonstandardSwaptionEngine::initialGuess(const Date &expiry) const {
101
102 Size fixedIdx =
103 std::upper_bound(arguments_.fixedResetDates.begin(),
104 arguments_.fixedResetDates.end(), expiry - 1) -
105 arguments_.fixedResetDates.begin();
106
107 Array initial(3);
108 Real nominalSum = 0.0, weightedRate = 0.0, ind = 0.0;
109 for (Size i = fixedIdx; i < arguments_.fixedResetDates.size(); i++) {
110 nominalSum += arguments_.fixedNominal[i];
111 Real rate = arguments_.fixedRate[i];
112 if (close(rate, 0.0))
113 rate = 0.03; // this value is at least better than zero
114 weightedRate += arguments_.fixedNominal[i] * rate;
115 if (arguments_.fixedNominal[i] > 1E-8) // exclude zero nominal periods
116 ind += 1.0;
117 }
118 Real nominalAvg = nominalSum / ind;
119
120 QL_REQUIRE(nominalSum > 0.0,
121 "sum of nominals on fixed leg must be positive ("
122 << nominalSum << ")");
123
124 weightedRate /= nominalSum;
125 initial[0] = nominalAvg;
126 initial[1] =
127 model_->termStructure()->timeFromReference(underlyingLastDate()) -
128 model_->termStructure()->timeFromReference(expiry);
129 initial[2] = weightedRate;
130
131 return initial;
132 }
133
134 void Gaussian1dNonstandardSwaptionEngine::calculate() const {
135
136 QL_REQUIRE(arguments_.settlementMethod != Settlement::ParYieldCurve,
137 "cash settled (ParYieldCurve) swaptions not priced with "
138 "Gaussian1dNonstandardSwaptionEngine");
139
140 Date settlement = model_->termStructure()->referenceDate();
141
142 if (arguments_.exercise->dates().back() <=
143 settlement) { // swaption is expired, possibly generated swap is not
144 // valued
145 results_.value = 0.0;
146 return;
147 }
148
149 ext::shared_ptr<RebatedExercise> rebatedExercise =
150 ext::dynamic_pointer_cast<RebatedExercise>(arguments_.exercise);
151
152 int idx = arguments_.exercise->dates().size() - 1;
153 int minIdxAlive = static_cast<int>(
154 std::upper_bound(arguments_.exercise->dates().begin(),
155 arguments_.exercise->dates().end(), settlement) -
156 arguments_.exercise->dates().begin());
157
158 NonstandardSwap swap = *arguments_.swap;
159 Option::Type type =
160 arguments_.type == Swap::Payer ? Option::Call : Option::Put;
161
162 Array npv0(2 * integrationPoints_ + 1, 0.0),
163 npv1(2 * integrationPoints_ + 1, 0.0);
164 Array z = model_->yGrid(stddevs_, integrationPoints_);
165 Array p(z.size(), 0.0);
166
167 // for probability computation
168 std::vector<Array> npvp0, npvp1;
169 if (probabilities_ != None) {
170 for (int i = 0; i < idx - minIdxAlive + 2; ++i) {
171 Array npvTmp0(2 * integrationPoints_ + 1, 0.0);
172 Array npvTmp1(2 * integrationPoints_ + 1, 0.0);
173 npvp0.push_back(npvTmp0);
174 npvp1.push_back(npvTmp1);
175 }
176 }
177 // end probabkility computation
178
179 Date expiry1 = Null<Date>(), expiry0;
180 Time expiry1Time = Null<Real>(), expiry0Time;
181
182 do {
183
184 if (idx == minIdxAlive - 1)
185 expiry0 = settlement;
186 else
187 expiry0 = arguments_.exercise->dates()[idx];
188
189 expiry0Time = std::max(
190 model_->termStructure()->timeFromReference(expiry0), 0.0);
191
192 Size j1 =
193 std::upper_bound(arguments_.fixedResetDates.begin(),
194 arguments_.fixedResetDates.end(), expiry0 - 1) -
195 arguments_.fixedResetDates.begin();
196 Size k1 =
197 std::upper_bound(arguments_.floatingResetDates.begin(),
198 arguments_.floatingResetDates.end(), expiry0 - 1) -
199 arguments_.floatingResetDates.begin();
200
201 // todo add openmp support later on (as in gaussian1dswaptionengine)
202
203 for (Size k = 0; k < (expiry0 > settlement ? npv0.size() : 1);
204 k++) {
205
206 Real price = 0.0;
207 if (expiry1Time != Null<Real>()) {
208 Real zSpreadDf =
209 oas_.empty() ? Real(1.0)
210 : std::exp(-oas_->value() *
211 (expiry1Time - expiry0Time));
212 Array yg = model_->yGrid(stddevs_, integrationPoints_,
213 expiry1Time, expiry0Time,
214 expiry0 > settlement ? z[k] : 0.0);
215 CubicInterpolation payoff0(
216 z.begin(), z.end(), npv1.begin(),
217 CubicInterpolation::Spline, true,
218 CubicInterpolation::Lagrange, 0.0,
219 CubicInterpolation::Lagrange, 0.0);
220 for (Size i = 0; i < yg.size(); i++) {
221 p[i] = payoff0(yg[i], true);
222 }
223 CubicInterpolation payoff1(
224 z.begin(), z.end(), p.begin(),
225 CubicInterpolation::Spline, true,
226 CubicInterpolation::Lagrange, 0.0,
227 CubicInterpolation::Lagrange, 0.0);
228 for (Size i = 0; i < z.size() - 1; i++) {
229 price += Gaussian1dModel::gaussianShiftedPolynomialIntegral(
230 0.0, payoff1.cCoefficients()[i],
231 payoff1.bCoefficients()[i],
232 payoff1.aCoefficients()[i], p[i], z[i],
233 z[i], z[i + 1]) *
234 zSpreadDf;
235 }
236 if (extrapolatePayoff_) {
237 if (flatPayoffExtrapolation_) {
238 price +=
239 Gaussian1dModel::gaussianShiftedPolynomialIntegral(
240 0.0, 0.0, 0.0, 0.0, p[z.size() - 2],
241 z[z.size() - 2], z[z.size() - 1], 100.0) *
242 zSpreadDf;
243 price += Gaussian1dModel::gaussianShiftedPolynomialIntegral(
244 0.0, 0.0, 0.0, 0.0, p[0], z[0], -100.0,
245 z[0]) *
246 zSpreadDf;
247 } else {
248 if (type == Option::Call)
249 price +=
250 Gaussian1dModel::gaussianShiftedPolynomialIntegral(
251 0.0,
252 payoff1.cCoefficients()[z.size() - 2],
253 payoff1.bCoefficients()[z.size() - 2],
254 payoff1.aCoefficients()[z.size() - 2],
255 p[z.size() - 2], z[z.size() - 2],
256 z[z.size() - 1], 100.0) *
257 zSpreadDf;
258 if (type == Option::Put)
259 price +=
260 Gaussian1dModel::gaussianShiftedPolynomialIntegral(
261 0.0, payoff1.cCoefficients()[0],
262 payoff1.bCoefficients()[0],
263 payoff1.aCoefficients()[0], p[0], z[0],
264 -100.0, z[0]) *
265 zSpreadDf;
266 }
267 }
268 }
269
270 npv0[k] = price;
271
272 // for probability computation
273 if (probabilities_ != None) {
274 for (Size m = 0; m < npvp0.size(); m++) {
275 Real price = 0.0;
276 if (expiry1Time != Null<Real>()) {
277 Real zSpreadDf =
278 oas_.empty()
279 ? Real(1.0)
280 : std::exp(-oas_->value() *
281 (expiry1Time - expiry0Time));
282 Array yg = model_->yGrid(
283 stddevs_, integrationPoints_, expiry1Time,
284 expiry0Time, expiry0 > settlement ? z[k] : 0.0);
285 CubicInterpolation payoff0(
286 z.begin(), z.end(), npvp1[m].begin(),
287 CubicInterpolation::Spline, true,
288 CubicInterpolation::Lagrange, 0.0,
289 CubicInterpolation::Lagrange, 0.0);
290 for (Size i = 0; i < yg.size(); i++) {
291 p[i] = payoff0(yg[i], true);
292 }
293 CubicInterpolation payoff1(
294 z.begin(), z.end(), p.begin(),
295 CubicInterpolation::Spline, true,
296 CubicInterpolation::Lagrange, 0.0,
297 CubicInterpolation::Lagrange, 0.0);
298 for (Size i = 0; i < z.size() - 1; i++) {
299 price +=
300 Gaussian1dModel::gaussianShiftedPolynomialIntegral(
301 0.0, payoff1.cCoefficients()[i],
302 payoff1.bCoefficients()[i],
303 payoff1.aCoefficients()[i], p[i], z[i],
304 z[i], z[i + 1]) *
305 zSpreadDf;
306 }
307 if (extrapolatePayoff_) {
308 if (flatPayoffExtrapolation_) {
309 price +=
310 Gaussian1dModel::gaussianShiftedPolynomialIntegral(
311 0.0, 0.0, 0.0, 0.0,
312 p[z.size() - 2],
313 z[z.size() - 2],
314 z[z.size() - 1], 100.0) *
315 zSpreadDf;
316 price +=
317 Gaussian1dModel::gaussianShiftedPolynomialIntegral(
318 0.0, 0.0, 0.0, 0.0, p[0],
319 z[0], -100.0, z[0]) *
320 zSpreadDf;
321 } else {
322 if (type == Option::Call)
323 price +=
324 Gaussian1dModel::gaussianShiftedPolynomialIntegral(
325 0.0,
326 payoff1.cCoefficients()
327 [z.size() - 2],
328 payoff1.bCoefficients()
329 [z.size() - 2],
330 payoff1.aCoefficients()
331 [z.size() - 2],
332 p[z.size() - 2],
333 z[z.size() - 2],
334 z[z.size() - 1], 100.0) *
335 zSpreadDf;
336 if (type == Option::Put)
337 price +=
338 Gaussian1dModel::gaussianShiftedPolynomialIntegral(
339 0.0,
340 payoff1
341 .cCoefficients()[0],
342 payoff1
343 .bCoefficients()[0],
344 payoff1
345 .aCoefficients()[0],
346 p[0], z[0], -100.0,
347 z[0]) *
348 zSpreadDf;
349 }
350 }
351 }
352
353 npvp0[m][k] = price;
354 }
355 }
356 // end probability computation
357
358 if (expiry0 > settlement) {
359 Real floatingLegNpv = 0.0;
360 for (Size l = k1; l < arguments_.floatingCoupons.size();
361 l++) {
362 Real zSpreadDf =
363 oas_.empty()
364 ? Real(1.0)
365 : std::exp(
366 -oas_->value() *
367 (model_->termStructure()
368 ->dayCounter()
369 .yearFraction(
370 expiry0,
371 arguments_
372 .floatingPayDates[l])));
373 Real amount;
374 if (arguments_.floatingIsRedemptionFlow[l])
375 amount = arguments_.floatingCoupons[l];
376 else
377 amount = arguments_.floatingNominal[l] *
378 arguments_.floatingAccrualTimes[l] *
379 (arguments_.floatingGearings[l] *
380 model_->forwardRate(
381 arguments_.floatingFixingDates[l],
382 expiry0, z[k],
383 arguments_.swap->iborIndex()) +
384 arguments_.floatingSpreads[l]);
385 floatingLegNpv +=
386 amount *
387 model_->zerobond(arguments_.floatingPayDates[l],
388 expiry0, z[k], discountCurve_) *
389 zSpreadDf;
390 }
391 Real fixedLegNpv = 0.0;
392 for (Size l = j1; l < arguments_.fixedCoupons.size(); l++) {
393 Real zSpreadDf =
394 oas_.empty()
395 ? Real(1.0)
396 : std::exp(
397 -oas_->value() *
398 (model_->termStructure()
399 ->dayCounter()
400 .yearFraction(
401 expiry0,
402 arguments_.fixedPayDates[l])));
403 fixedLegNpv +=
404 arguments_.fixedCoupons[l] *
405 model_->zerobond(arguments_.fixedPayDates[l],
406 expiry0, z[k], discountCurve_) *
407 zSpreadDf;
408 }
409 Real rebate = 0.0;
410 Real zSpreadDf = 1.0;
411 Date rebateDate = expiry0;
412 if (rebatedExercise != nullptr) {
413 rebate = rebatedExercise->rebate(idx);
414 rebateDate = rebatedExercise->rebatePaymentDate(idx);
415 zSpreadDf =
416 oas_.empty()
417 ? Real(1.0)
418 : std::exp(
419 -oas_->value() *
420 (model_->termStructure()
421 ->dayCounter()
422 .yearFraction(expiry0, rebateDate)));
423 }
424 Real exerciseValue =
425 ((type == Option::Call ? 1.0 : -1.0) *
426 (floatingLegNpv - fixedLegNpv) +
427 rebate * model_->zerobond(rebateDate, expiry0, z[k],
428 discountCurve_) *
429 zSpreadDf) /
430 model_->numeraire(expiry0Time, z[k], discountCurve_);
431
432 // for probability computation
433 if (probabilities_ != None) {
434 if (idx == static_cast<int>(
435 arguments_.exercise->dates().size()) -
436 1) // if true we are at the latest date,
437 // so we init
438 // the no call probability
439 npvp0.back()[k] =
440 probabilities_ == Naive
441 ? Real(1.0)
442 : 1.0 / (model_->zerobond(expiry0Time, 0.0,
443 0.0,
444 discountCurve_) *
445 model_->numeraire(expiry0, z[k],
446 discountCurve_));
447 if (exerciseValue >= npv0[k]) {
448 npvp0[idx - minIdxAlive][k] =
449 probabilities_ == Naive
450 ? Real(1.0)
451 : 1.0 /
452 (model_->zerobond(expiry0Time, 0.0,
453 0.0,
454 discountCurve_) *
455 model_->numeraire(expiry0Time, z[k],
456 discountCurve_));
457 for (Size ii = idx - minIdxAlive + 1;
458 ii < npvp0.size(); ii++)
459 npvp0[ii][k] = 0.0;
460 }
461 }
462 // end probability computation
463
464 npv0[k] = std::max(npv0[k], exerciseValue);
465 }
466 }
467
468 npv1.swap(npv0);
469
470 // for probability computation
471 if (probabilities_ != None) {
472 for (Size i = 0; i < npvp0.size(); i++)
473 npvp1[i].swap(npvp0[i]);
474 }
475 // end probability computation
476
477 expiry1 = expiry0;
478 expiry1Time = expiry0Time;
479
480 } while (--idx >= minIdxAlive - 1);
481
482 results_.value = npv1[0] * model_->numeraire(0.0, 0.0, discountCurve_);
483
484 // for probability computation
485 if (probabilities_ != None) {
486 std::vector<Real> prob(npvp0.size());
487 for (Size i = 0; i < npvp0.size(); i++) {
488 prob[i] = npvp1[i][0] *
489 (probabilities_ == Naive
490 ? 1.0
491 : model_->numeraire(0.0, 0.0, discountCurve_));
492 }
493 results_.additionalResults["probabilities"] = prob;
494 }
495 // end probability computation
496 }
497}
QuantLib::Array
1-D array used in linear algebra.
Definition: array.hpp:52
QuantLib::Array::swap
void swap(Array &) noexcept
Definition: array.hpp:546
QuantLib::Array::end
const_iterator end() const
Definition: array.hpp:511
QuantLib::Array::size
Size size() const
dimension of the array
Definition: array.hpp:495
QuantLib::Array::begin
const_iterator begin() const
Definition: array.hpp:503
QuantLib::CubicInterpolation
Cubic interpolation between discrete points.
Definition: cubicinterpolation.hpp:105
QuantLib::CubicInterpolation::cCoefficients
const std::vector< Real > & cCoefficients() const
Definition: cubicinterpolation.hpp:182
QuantLib::CubicInterpolation::bCoefficients
const std::vector< Real > & bCoefficients() const
Definition: cubicinterpolation.hpp:181
QuantLib::CubicInterpolation::aCoefficients
const std::vector< Real > & aCoefficients() const
Definition: cubicinterpolation.hpp:180
QuantLib::Date
Concrete date class.
Definition: date.hpp:125
QuantLib::Gaussian1dModel::gaussianShiftedPolynomialIntegral
static Real gaussianShiftedPolynomialIntegral(Real a, Real b, Real c, Real d, Real e, Real h, Real x0, Real x1)
Definition: gaussian1dmodel.cpp:240
QuantLib::Gaussian1dNonstandardSwaptionEngine::initialGuess
const Array initialGuess(const Date &expiry) const override
Definition: gaussian1dnonstandardswaptionengine.cpp:100
QuantLib::Gaussian1dNonstandardSwaptionEngine::underlyingNpv
Real underlyingNpv(const Date &expiry, Real y) const override
Definition: gaussian1dnonstandardswaptionengine.cpp:32
QuantLib::Null
template class providing a null value for a given type.
Definition: null.hpp:76
QuantLib::Time
Real Time
continuous quantity with 1-year units
Definition: types.hpp:62
QuantLib::Real
QL_REAL Real
real number
Definition: types.hpp:50
QuantLib::Size
std::size_t Size
size of a container
Definition: types.hpp:58
QuantLib
Definition: any.hpp:35
QuantLib::swap
void swap(Array &v, Array &w) noexcept
Definition: array.hpp:903
QuantLib::close
bool close(const Quantity &m1, const Quantity &m2, Size n)
Definition: quantity.cpp:163
Generated by Doxygen 1.9.5