Coverage Report

Created: 2026-08-14 07:10

next uncovered line (L), next uncovered region (R), next uncovered branch (B)
/src/quantlib/ql/pricingengines/swaption/gaussian1dfloatfloatswaptionengine.cpp
Line
Count
Source
1
/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2
3
/*
4
 Copyright (C) 2013, 2015 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
 <https://www.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/gaussian1dfloatfloatswaptionengine.hpp>
21
#include <ql/experimental/coupons/swapspreadindex.hpp> // internal
22
#include <ql/math/interpolations/cubicinterpolation.hpp>
23
#include <ql/payoff.hpp>
24
25
namespace QuantLib {
26
27
0
    void Gaussian1dFloatFloatSwaptionEngine::calculate() const {
28
29
0
        QL_REQUIRE(arguments_.settlementMethod != Settlement::ParYieldCurve,
30
0
                   "cash settled (ParYieldCurve) swaptions not priced with "
31
0
                   "Gaussian1dFloatFloatSwaptionEngine");
32
33
0
        Date settlement = model_->termStructure()->referenceDate();
34
35
0
        if (arguments_.exercise->dates().back() <=
36
0
            settlement) { // swaption is expired,
37
                          // possibly generated swap
38
                          // is not
39
                          // valued
40
0
            results_.value = 0.0;
41
0
            return;
42
0
        }
43
44
0
        rebatedExercise_ =
45
0
            ext::dynamic_pointer_cast<RebatedExercise>(arguments_.exercise);
46
47
0
        std::pair<Real, Real> result =
48
0
            npvs(settlement, 0.0, includeTodaysExercise_, true);
49
50
0
        results_.value = result.first;
51
0
        results_.additionalResults["underlyingValue"] = result.second;
52
0
    }
53
54
    Real
55
    Gaussian1dFloatFloatSwaptionEngine::underlyingNpv(const Date &expiry,
56
0
                                                      const Real y) const {
57
0
        return npvs(expiry, y, true).second;
58
0
    }
59
60
0
    Swap::Type Gaussian1dFloatFloatSwaptionEngine::underlyingType() const {
61
0
        return arguments_.swap->type();
62
0
    }
63
64
    // NOLINTNEXTLINE(readability-const-return-type)
65
0
    const Date Gaussian1dFloatFloatSwaptionEngine::underlyingLastDate() const {
66
0
        Date l1 = arguments_.leg1PayDates.back();
67
0
        Date l2 = arguments_.leg2PayDates.back();
68
0
        return l2 >= l1 ? l2 : l1;
69
0
    }
70
71
    // NOLINTNEXTLINE(readability-const-return-type)
72
0
    const Array Gaussian1dFloatFloatSwaptionEngine::initialGuess(const Date &expiry) const {
73
74
0
        Size idx1 =
75
0
            std::upper_bound(arguments_.leg1ResetDates.begin(),
76
0
                             arguments_.leg1ResetDates.end(), expiry - 1) -
77
0
            arguments_.leg1ResetDates.begin();
78
79
        // very simple initial guess
80
        // check guess for nominal and weighted maturity !
81
82
0
        Real nominalSum1 = 0.0;
83
0
        for (Size i = idx1; i < arguments_.leg1ResetDates.size(); i++) {
84
0
            nominalSum1 += arguments_.nominal1[i];
85
0
        }
86
0
        Real nominalAvg1 = nominalSum1 /
87
0
            (arguments_.leg1ResetDates.size() - idx1);
88
0
        Real weightedMaturity1 = 0.0;
89
0
        for (Size i = idx1; i < arguments_.leg1ResetDates.size(); i++) {
90
0
            weightedMaturity1 +=
91
0
                arguments_.leg1AccrualTimes[i] * arguments_.nominal1[i];
92
0
        }
93
0
        weightedMaturity1 /= nominalAvg1;
94
95
0
        return {
96
0
            nominalAvg1,
97
0
            weightedMaturity1,
98
0
            0.03 // ???
99
0
        };
100
0
    }
101
102
    // calculate npv and underlying npv as of expiry date
103
    std::pair<Real, Real>
104
    Gaussian1dFloatFloatSwaptionEngine::npvs(const Date& expiry,
105
                                             const Real y,
106
                                             const bool includeExerciseOnExpiry,
107
0
                                             const bool considerProbabilities) const {
108
109
        // pricing
110
111
        // event dates are coupon fixing dates and exercise dates
112
        // we explicitly estimate cms and also libor coupons (although
113
        // the latter could be calculated analytically) to make the code
114
        // simpler
115
116
0
        std::vector<Date> events;
117
0
        events.insert(events.end(), arguments_.exercise->dates().begin(),
118
0
                      arguments_.exercise->dates().end());
119
0
        events.insert(events.end(), arguments_.leg1FixingDates.begin(),
120
0
                      arguments_.leg1FixingDates.end());
121
0
        events.insert(events.end(), arguments_.leg2FixingDates.begin(),
122
0
                      arguments_.leg2FixingDates.end());
123
0
        std::sort(events.begin(), events.end());
124
125
0
        auto it = std::unique(events.begin(), events.end());
126
0
        events.resize(std::distance(events.begin(), it));
127
128
        // only events on or after expiry are of interest by definition of the
129
        // deal part that is exericsed into.
130
131
0
        auto filit = std::upper_bound(events.begin(), events.end(),
132
0
                                      expiry - (includeExerciseOnExpiry ? 1 : 0));
133
0
        events.erase(events.begin(), filit);
134
135
0
        int idx = events.size() - 1;
136
137
0
        FloatFloatSwap swap = *arguments_.swap;
138
0
        Option::Type type =
139
0
            arguments_.type == Swap::Payer ? Option::Call : Option::Put;
140
141
0
        Array npv0(2 * integrationPoints_ + 1, 0.0),
142
0
            npv1(2 * integrationPoints_ + 1, 0.0); // arrays for npvs of the
143
                                                   // option
144
0
        Array npv0a(2 * integrationPoints_ + 1, 0.0),
145
0
            npv1a(2 * integrationPoints_ + 1, 0.0); // arrays for npvs of the
146
                                                    // underlying
147
0
        Array z = model_->yGrid(stddevs_, integrationPoints_);
148
0
        Array p(z.size(), 0.0), pa(z.size(), 0.0);
149
150
        // for probability computation
151
0
        std::vector<Array> npvp0, npvp1;
152
        // how many active exercise dates are there ?
153
0
        Size noEx =  arguments_.exercise->dates().size() -
154
0
            (std::upper_bound(arguments_.exercise->dates().begin(),
155
0
                         arguments_.exercise->dates().end(),
156
0
                         expiry - (includeExerciseOnExpiry ? 1 : 0)) -
157
0
             arguments_.exercise->dates().begin());
158
0
        Size exIdx = noEx; // current exercise index
159
0
        if (considerProbabilities && probabilities_ != None) {
160
0
            for (Size i = 0; i < noEx+1 ; ++i) {
161
0
                Array npvTmp0(2 * integrationPoints_ + 1, 0.0);
162
0
                Array npvTmp1(2 * integrationPoints_ + 1, 0.0);
163
0
                npvp0.push_back(npvTmp0);
164
0
                npvp1.push_back(npvTmp1);
165
0
            }
166
0
        }
167
        // end probability computation
168
169
0
        Date event1 = Date(), event0;
170
0
        Time event1Time = Null<Real>(), event0Time;
171
172
0
        ext::shared_ptr<IborIndex> ibor1 =
173
0
            ext::dynamic_pointer_cast<IborIndex>(arguments_.index1);
174
0
        ext::shared_ptr<SwapIndex> cms1 =
175
0
            ext::dynamic_pointer_cast<SwapIndex>(arguments_.index1);
176
0
        ext::shared_ptr<SwapSpreadIndex> cmsspread1 =
177
0
            ext::dynamic_pointer_cast<SwapSpreadIndex>(arguments_.index1);
178
0
        ext::shared_ptr<IborIndex> ibor2 =
179
0
            ext::dynamic_pointer_cast<IborIndex>(arguments_.index2);
180
0
        ext::shared_ptr<SwapIndex> cms2 =
181
0
            ext::dynamic_pointer_cast<SwapIndex>(arguments_.index2);
182
0
        ext::shared_ptr<SwapSpreadIndex> cmsspread2 =
183
0
            ext::dynamic_pointer_cast<SwapSpreadIndex>(arguments_.index2);
184
185
0
        QL_REQUIRE(ibor1 != nullptr || cms1 != nullptr || cmsspread1 != nullptr,
186
0
                   "index1 must be ibor or swap or swap spread index");
187
0
        QL_REQUIRE(ibor2 != nullptr || cms2 != nullptr || cmsspread2 != nullptr,
188
0
                   "index2 must be ibor or swap or swap spread index");
189
190
0
        do {
191
192
            // we are at event0 date, which can be a structured coupon fixing
193
            // date or an exercise date or both.
194
195
0
            bool isEventDate = true;
196
0
            if (idx == -1) {
197
0
                event0 = expiry;
198
0
                isEventDate = false;
199
0
            } else {
200
0
                event0 = events[idx];
201
0
                if (event0 == expiry)
202
0
                    idx = -1; // avoid double roll back if expiry equal to
203
                              // earliest event date
204
0
            }
205
206
0
            bool isExercise =
207
0
                std::find(arguments_.exercise->dates().begin(), arguments_.exercise->dates().end(),
208
0
                          event0) != arguments_.exercise->dates().end();
209
210
0
            bool isLeg1Fixing =
211
0
                std::find(arguments_.leg1FixingDates.begin(), arguments_.leg1FixingDates.end(),
212
0
                          event0) != arguments_.leg1FixingDates.end();
213
214
0
            bool isLeg2Fixing =
215
0
                std::find(arguments_.leg2FixingDates.begin(), arguments_.leg2FixingDates.end(),
216
0
                          event0) != arguments_.leg2FixingDates.end();
217
218
0
            event0Time = std::max(
219
0
                model_->termStructure()->timeFromReference(event0), 0.0);
220
221
            // todo add openmp support later on (as in gaussian1dswaptionengine)
222
223
0
            for (Size k = 0; k < (event0 > expiry ? npv0.size() : 1); k++) {
224
225
                // roll back
226
227
0
                Real price = 0.0, pricea = 0.0;
228
0
                if (event1Time != Null<Real>()) {
229
0
                    Real zSpreadDf = oas_.empty()
230
0
                                         ? Real(1.0)
231
0
                                         : std::exp(-oas_->value() *
232
0
                                                    (event1Time - event0Time));
233
0
                    Array yg =
234
0
                        model_->yGrid(stddevs_, integrationPoints_, event1Time,
235
0
                                      event0Time, event0 > expiry ? z[k] : y);
236
0
                    CubicInterpolation payoff0(
237
0
                        z.begin(), z.end(), npv1.begin(),
238
0
                        CubicInterpolation::Spline, true,
239
0
                        CubicInterpolation::Lagrange, 0.0,
240
0
                        CubicInterpolation::Lagrange, 0.0);
241
0
                    CubicInterpolation payoff0a(
242
0
                        z.begin(), z.end(), npv1a.begin(),
243
0
                        CubicInterpolation::Spline, true,
244
0
                        CubicInterpolation::Lagrange, 0.0,
245
0
                        CubicInterpolation::Lagrange, 0.0);
246
0
                    for (Size i = 0; i < yg.size(); i++) {
247
0
                        p[i] = payoff0(yg[i], true);
248
0
                        pa[i] = payoff0a(yg[i], true);
249
0
                    }
250
0
                    CubicInterpolation payoff1(
251
0
                        z.begin(), z.end(), p.begin(),
252
0
                        CubicInterpolation::Spline, true,
253
0
                        CubicInterpolation::Lagrange, 0.0,
254
0
                        CubicInterpolation::Lagrange, 0.0);
255
0
                    CubicInterpolation payoff1a(
256
0
                        z.begin(), z.end(), pa.begin(),
257
0
                        CubicInterpolation::Spline, true,
258
0
                        CubicInterpolation::Lagrange, 0.0,
259
0
                        CubicInterpolation::Lagrange, 0.0);
260
0
                    for (Size i = 0; i < z.size() - 1; i++) {
261
0
                        price += Gaussian1dModel::gaussianShiftedPolynomialIntegral(
262
0
                                     0.0, payoff1.cCoefficients()[i],
263
0
                                     payoff1.bCoefficients()[i],
264
0
                                     payoff1.aCoefficients()[i], p[i], z[i],
265
0
                                     z[i], z[i + 1]) *
266
0
                                 zSpreadDf;
267
0
                        pricea += Gaussian1dModel::gaussianShiftedPolynomialIntegral(
268
0
                                      0.0, payoff1a.cCoefficients()[i],
269
0
                                      payoff1a.bCoefficients()[i],
270
0
                                      payoff1a.aCoefficients()[i], pa[i], z[i],
271
0
                                      z[i], z[i + 1]) *
272
0
                                  zSpreadDf;
273
0
                    }
274
0
                    if (extrapolatePayoff_) {
275
0
                        if (flatPayoffExtrapolation_) {
276
0
                            price +=
277
0
                                Gaussian1dModel::gaussianShiftedPolynomialIntegral(
278
0
                                    0.0, 0.0, 0.0, 0.0, p[z.size() - 2],
279
0
                                    z[z.size() - 2], z[z.size() - 1], 100.0) *
280
0
                                zSpreadDf;
281
0
                            price += Gaussian1dModel::gaussianShiftedPolynomialIntegral(
282
0
                                         0.0, 0.0, 0.0, 0.0, p[0], z[0], -100.0,
283
0
                                         z[0]) *
284
0
                                     zSpreadDf;
285
0
                            pricea +=
286
0
                                Gaussian1dModel::gaussianShiftedPolynomialIntegral(
287
0
                                    0.0, 0.0, 0.0, 0.0, pa[z.size() - 2],
288
0
                                    z[z.size() - 2], z[z.size() - 1], 100.0) *
289
0
                                zSpreadDf;
290
0
                            pricea += Gaussian1dModel::gaussianShiftedPolynomialIntegral(
291
0
                                          0.0, 0.0, 0.0, 0.0, pa[0], z[0],
292
0
                                          -100.0, z[0]) *
293
0
                                      zSpreadDf;
294
0
                        } else {
295
0
                            if (type == Option::Call)
296
0
                                price +=
297
0
                                    Gaussian1dModel::gaussianShiftedPolynomialIntegral(
298
0
                                        0.0,
299
0
                                        payoff1.cCoefficients()[z.size() - 2],
300
0
                                        payoff1.bCoefficients()[z.size() - 2],
301
0
                                        payoff1.aCoefficients()[z.size() - 2],
302
0
                                        p[z.size() - 2], z[z.size() - 2],
303
0
                                        z[z.size() - 1], 100.0) *
304
0
                                    zSpreadDf;
305
0
                            if (type == Option::Put)
306
0
                                price +=
307
0
                                    Gaussian1dModel::gaussianShiftedPolynomialIntegral(
308
0
                                        0.0, payoff1.cCoefficients()[0],
309
0
                                        payoff1.bCoefficients()[0],
310
0
                                        payoff1.aCoefficients()[0], p[0], z[0],
311
0
                                        -100.0, z[0]) *
312
0
                                    zSpreadDf;
313
0
                            if (type == Option::Call)
314
0
                                pricea +=
315
0
                                    Gaussian1dModel::gaussianShiftedPolynomialIntegral(
316
0
                                        0.0,
317
0
                                        payoff1a.cCoefficients()[z.size() - 2],
318
0
                                        payoff1a.bCoefficients()[z.size() - 2],
319
0
                                        payoff1a.aCoefficients()[z.size() - 2],
320
0
                                        pa[z.size() - 2], z[z.size() - 2],
321
0
                                        z[z.size() - 1], 100.0) *
322
0
                                    zSpreadDf;
323
0
                            if (type == Option::Put)
324
0
                                pricea +=
325
0
                                    Gaussian1dModel::gaussianShiftedPolynomialIntegral(
326
0
                                        0.0, payoff1a.cCoefficients()[0],
327
0
                                        payoff1a.bCoefficients()[0],
328
0
                                        payoff1a.aCoefficients()[0], pa[0],
329
0
                                        z[0], -100.0, z[0]) *
330
0
                                    zSpreadDf;
331
0
                        }
332
0
                    }
333
0
                }
334
335
0
                npv0[k] = price;
336
0
                npv0a[k] = pricea;
337
338
                // for probability computation
339
0
                if (considerProbabilities && probabilities_ != None) {
340
0
                    for (Size m = 0; m < npvp0.size(); m++) {
341
0
                        Real price = 0.0;
342
0
                        if (event1Time != Null<Real>()) {
343
0
                            Real zSpreadDf =
344
0
                                oas_.empty()
345
0
                                    ? Real(1.0)
346
0
                                    : std::exp(-oas_->value() *
347
0
                                               (event1Time - event0Time));
348
0
                            Array yg = model_->yGrid(
349
0
                                stddevs_, integrationPoints_, event1Time,
350
0
                                event0Time, event0 > expiry ? z[k] : 0.0);
351
0
                            CubicInterpolation payoff0(
352
0
                                z.begin(), z.end(), npvp1[m].begin(),
353
0
                                CubicInterpolation::Spline, true,
354
0
                                CubicInterpolation::Lagrange, 0.0,
355
0
                                CubicInterpolation::Lagrange, 0.0);
356
0
                            for (Size i = 0; i < yg.size(); i++) {
357
0
                                p[i] = payoff0(yg[i], true);
358
0
                            }
359
0
                            CubicInterpolation payoff1(
360
0
                                z.begin(), z.end(), p.begin(),
361
0
                                CubicInterpolation::Spline, true,
362
0
                                CubicInterpolation::Lagrange, 0.0,
363
0
                                CubicInterpolation::Lagrange, 0.0);
364
0
                            for (Size i = 0; i < z.size() - 1; i++) {
365
0
                                price +=
366
0
                                    Gaussian1dModel::gaussianShiftedPolynomialIntegral(
367
0
                                        0.0, payoff1.cCoefficients()[i],
368
0
                                        payoff1.bCoefficients()[i],
369
0
                                        payoff1.aCoefficients()[i], p[i], z[i],
370
0
                                        z[i], z[i + 1]) *
371
0
                                    zSpreadDf;
372
0
                            }
373
0
                            if (extrapolatePayoff_) {
374
0
                                if (flatPayoffExtrapolation_) {
375
0
                                    price +=
376
0
                                        Gaussian1dModel::gaussianShiftedPolynomialIntegral(
377
0
                                                  0.0, 0.0, 0.0, 0.0,
378
0
                                                  p[z.size() - 2],
379
0
                                                  z[z.size() - 2],
380
0
                                                  z[z.size() - 1], 100.0) *
381
0
                                        zSpreadDf;
382
0
                                    price +=
383
0
                                        Gaussian1dModel::gaussianShiftedPolynomialIntegral(
384
0
                                                  0.0, 0.0, 0.0, 0.0, p[0],
385
0
                                                  z[0], -100.0, z[0]) *
386
0
                                        zSpreadDf;
387
0
                                } else {
388
0
                                    if (type == Option::Call)
389
0
                                        price +=
390
0
                                            Gaussian1dModel::gaussianShiftedPolynomialIntegral(
391
0
                                                      0.0,
392
0
                                                      payoff1.cCoefficients()
393
0
                                                          [z.size() - 2],
394
0
                                                      payoff1.bCoefficients()
395
0
                                                          [z.size() - 2],
396
0
                                                      payoff1.aCoefficients()
397
0
                                                          [z.size() - 2],
398
0
                                                      p[z.size() - 2],
399
0
                                                      z[z.size() - 2],
400
0
                                                      z[z.size() - 1], 100.0) *
401
0
                                            zSpreadDf;
402
0
                                    if (type == Option::Put)
403
0
                                        price +=
404
0
                                            Gaussian1dModel::gaussianShiftedPolynomialIntegral(
405
0
                                                      0.0,
406
0
                                                      payoff1
407
0
                                                          .cCoefficients()[0],
408
0
                                                      payoff1
409
0
                                                          .bCoefficients()[0],
410
0
                                                      payoff1
411
0
                                                          .aCoefficients()[0],
412
0
                                                      p[0], z[0], -100.0,
413
0
                                                      z[0]) *
414
0
                                            zSpreadDf;
415
0
                                }
416
0
                            }
417
0
                        }
418
419
0
                        npvp0[m][k] = price;
420
0
                    }
421
0
                }
422
                // end probability computation
423
424
                // event date calculations
425
426
0
                if (isEventDate) {
427
428
0
                    Real zk = event0 > expiry ? z[k] : y;
429
430
0
                    if (isLeg1Fixing) { // if event is a fixing date and
431
                                        // exercise date,
432
                        // the coupon is part of the exercise into right (by
433
                        // definition)
434
0
                        Size j = std::find(arguments_.leg1FixingDates.begin(),
435
0
                                           arguments_.leg1FixingDates.end(),
436
0
                                           event0) -
437
0
                                 arguments_.leg1FixingDates.begin();
438
0
                        Real zSpreadDf =
439
0
                            oas_.empty()
440
0
                                ? Real(1.0)
441
0
                                : std::exp(
442
0
                                      -oas_->value() *
443
0
                                      (model_->termStructure()
444
0
                                           ->dayCounter()
445
0
                                           .yearFraction(
446
0
                                                event0,
447
0
                                                arguments_.leg1PayDates[j])));
448
0
                        bool done = false;
449
0
                        do {
450
0
                            Real amount;
451
0
                            if (arguments_.leg1IsRedemptionFlow[j]) {
452
0
                                amount = arguments_.leg1Coupons[j];
453
0
                            } else {
454
0
                                Real estFixing = 0.0;
455
0
                                if (ibor1 != nullptr) {
456
0
                                    estFixing = model_->forwardRate(
457
0
                                        arguments_.leg1FixingDates[j], event0,
458
0
                                        zk, ibor1);
459
0
                                }
460
0
                                if (cms1 != nullptr) {
461
0
                                    estFixing = model_->swapRate(
462
0
                                        arguments_.leg1FixingDates[j],
463
0
                                        cms1->tenor(), event0, zk, cms1);
464
0
                                }
465
0
                                if (cmsspread1 != nullptr)
466
0
                                    estFixing =
467
0
                                        cmsspread1->gearing1() *
468
0
                                            model_->swapRate(
469
0
                                                arguments_.leg1FixingDates[j],
470
0
                                                cmsspread1->swapIndex1()
471
0
                                                    ->tenor(),
472
0
                                                event0, zk,
473
0
                                                cmsspread1->swapIndex1()) +
474
0
                                        cmsspread1->gearing2() *
475
0
                                            model_->swapRate(
476
0
                                                arguments_.leg1FixingDates[j],
477
0
                                                cmsspread1->swapIndex2()
478
0
                                                    ->tenor(),
479
0
                                                event0, zk,
480
0
                                                cmsspread1->swapIndex2());
481
0
                                Real rate =
482
0
                                    arguments_.leg1Spreads[j] +
483
0
                                    arguments_.leg1Gearings[j] * estFixing;
484
0
                                if (arguments_.leg1CappedRates[j] !=
485
0
                                    Null<Real>())
486
0
                                    rate = std::min(
487
0
                                        arguments_.leg1CappedRates[j], rate);
488
0
                                if (arguments_.leg1FlooredRates[j] !=
489
0
                                    Null<Real>())
490
0
                                    rate = std::max(
491
0
                                        arguments_.leg1FlooredRates[j], rate);
492
0
                                amount = rate * arguments_.nominal1[j] *
493
0
                                         arguments_.leg1AccrualTimes[j];
494
0
                            }
495
496
0
                            npv0a[k] -=
497
0
                                amount *
498
0
                                model_->zerobond(arguments_.leg1PayDates[j],
499
0
                                                 event0, zk, discountCurve_) /
500
0
                                model_->numeraire(event0Time, zk,
501
0
                                                  discountCurve_) *
502
0
                                zSpreadDf;
503
504
0
                            if (j < arguments_.leg1FixingDates.size() - 1) {
505
0
                                j++;
506
0
                                done =
507
0
                                    (event0 != arguments_.leg1FixingDates[j]);
508
0
                            } else
509
0
                                done = true;
510
511
0
                        } while (!done);
512
0
                    }
513
514
0
                    if (isLeg2Fixing) { // if event is a fixing date and
515
                                        // exercise date,
516
                        // the coupon is part of the exercise into right (by
517
                        // definition)
518
0
                        Size j = std::find(arguments_.leg2FixingDates.begin(),
519
0
                                           arguments_.leg2FixingDates.end(),
520
0
                                           event0) -
521
0
                                 arguments_.leg2FixingDates.begin();
522
0
                        Real zSpreadDf =
523
0
                            oas_.empty()
524
0
                                ? Real(1.0)
525
0
                                : std::exp(
526
0
                                      -oas_->value() *
527
0
                                      (model_->termStructure()
528
0
                                           ->dayCounter()
529
0
                                           .yearFraction(
530
0
                                                event0,
531
0
                                                arguments_.leg2PayDates[j])));
532
0
                        bool done;
533
0
                        do {
534
0
                            Real amount;
535
0
                            if (arguments_.leg2IsRedemptionFlow[j]) {
536
0
                                amount = arguments_.leg2Coupons[j];
537
0
                            } else {
538
0
                                Real estFixing = 0.0;
539
0
                                if (ibor2 != nullptr)
540
0
                                    estFixing = model_->forwardRate(arguments_.leg2FixingDates[j],event0,zk,ibor2);
541
0
                                if (cms2 != nullptr)
542
0
                                    estFixing = model_->swapRate(arguments_.leg2FixingDates[j],cms2->tenor(),event0,zk,cms2);
543
0
                                if (cmsspread2 != nullptr)
544
0
                                    estFixing =
545
0
                                        cmsspread2->gearing1() *
546
0
                                            model_->swapRate(
547
0
                                                arguments_.leg2FixingDates[j],
548
0
                                                cmsspread2->swapIndex1()
549
0
                                                    ->tenor(),
550
0
                                                event0, zk,
551
0
                                                cmsspread2->swapIndex1()) +
552
0
                                        cmsspread2->gearing2() *
553
0
                                            model_->swapRate(
554
0
                                                arguments_.leg2FixingDates[j],
555
0
                                                cmsspread2->swapIndex2()
556
0
                                                    ->tenor(),
557
0
                                                event0, zk,
558
0
                                                cmsspread2->swapIndex2());
559
0
                                Real rate =
560
0
                                    arguments_.leg2Spreads[j] +
561
0
                                    arguments_.leg2Gearings[j] * estFixing;
562
0
                                if (arguments_.leg2CappedRates[j] !=
563
0
                                    Null<Real>())
564
0
                                    rate = std::min(
565
0
                                        arguments_.leg2CappedRates[j], rate);
566
0
                                if (arguments_.leg2FlooredRates[j] !=
567
0
                                    Null<Real>())
568
0
                                    rate = std::max(
569
0
                                        arguments_.leg2FlooredRates[j], rate);
570
0
                                amount = rate * arguments_.nominal2[j] *
571
0
                                         arguments_.leg2AccrualTimes[j];
572
0
                            }
573
574
0
                            npv0a[k] +=
575
0
                                amount *
576
0
                                model_->zerobond(arguments_.leg2PayDates[j],
577
0
                                                 event0, zk, discountCurve_) /
578
0
                                model_->numeraire(event0Time, zk,
579
0
                                                  discountCurve_) *
580
0
                                zSpreadDf;
581
0
                            if (j < arguments_.leg2FixingDates.size() - 1) {
582
0
                                j++;
583
0
                                done =
584
0
                                    (event0 != arguments_.leg2FixingDates[j]);
585
0
                            } else
586
0
                                done = true;
587
588
0
                        } while (!done);
589
0
                    }
590
591
0
                    if (isExercise) {
592
0
                        Size j = std::find(arguments_.exercise->dates().begin(),
593
0
                                           arguments_.exercise->dates().end(),
594
0
                                           event0) -
595
0
                                 arguments_.exercise->dates().begin();
596
0
                        Real rebate = 0.0;
597
0
                        Real zSpreadDf = 1.0;
598
0
                        Date rebateDate = event0;
599
0
                        if (rebatedExercise_ != nullptr) {
600
0
                            rebate = rebatedExercise_->rebate(j);
601
0
                            rebateDate = rebatedExercise_->rebatePaymentDate(j);
602
0
                            zSpreadDf =
603
0
                                oas_.empty()
604
0
                                    ? Real(1.0)
605
0
                                    : std::exp(-oas_->value() *
606
0
                                               (model_->termStructure()
607
0
                                                    ->dayCounter()
608
0
                                                    .yearFraction(event0,
609
0
                                                                  rebateDate)));
610
0
                        }
611
0
                        Real exerciseValue =
612
0
                            (type == Option::Call ? 1.0 : -1.0) * npv0a[k] +
613
0
                            rebate * model_->zerobond(rebateDate, event0) *
614
0
                                zSpreadDf / model_->numeraire(event0Time, zk,
615
0
                                                              discountCurve_);
616
617
0
                        if (considerProbabilities && probabilities_ != None) {
618
0
                            if (exIdx == noEx) {
619
                                // if true we are at the latest date,
620
                                // so we init
621
                                // the no call probability
622
0
                                npvp0.back()[k] =
623
0
                                    probabilities_ == Naive
624
0
                                        ? Real(1.0)
625
0
                                        : 1.0 / (model_->zerobond(
626
0
                                                     event0Time, 0.0, 0.0,
627
0
                                                     discountCurve_) *
628
0
                                                 model_->numeraire(
629
0
                                                     event0, z[k],
630
0
                                                     discountCurve_));
631
0
                            }
632
0
                            if (exerciseValue >= npv0[k]) {
633
0
                                npvp0[exIdx-1][k] =
634
0
                                    probabilities_ == Naive
635
0
                                        ? Real(1.0)
636
0
                                        : 1.0 / (model_->zerobond(
637
0
                                                     event0Time, 0.0, 0.0,
638
0
                                                     discountCurve_) *
639
0
                                                 model_->numeraire(
640
0
                                                     event0Time, z[k],
641
0
                                                     discountCurve_));
642
0
                                for (Size ii = exIdx; ii < noEx+1; ++ii)
643
0
                                    npvp0[ii][k] = 0.0;
644
0
                            }
645
0
                        }
646
                        // end probability computation
647
648
0
                        npv0[k] = std::max(npv0[k], exerciseValue);
649
0
                    }
650
0
                }
651
0
            }
652
653
0
            if(isExercise)
654
0
                --exIdx;
655
656
0
            npv1.swap(npv0);
657
0
            npv1a.swap(npv0a);
658
659
            // for probability computation
660
0
            if(considerProbabilities && probabilities_ != None) {
661
0
                for(Size i=0;i<npvp0.size();++i) {
662
0
                    npvp1[i].swap(npvp0[i]);
663
0
                }
664
0
            }
665
            // end probability computation
666
667
0
            event1 = event0;
668
0
            event1Time = event0Time;
669
670
0
        } while (--idx >= -1);
671
672
0
        std::pair<Real, Real> res(
673
0
            npv1[0] * model_->numeraire(event1Time, y, discountCurve_),
674
0
            npv1a[0] * model_->numeraire(event1Time, y, discountCurve_) *
675
0
                (type == Option::Call ? 1.0 : -1.0));
676
677
        // for probability computation
678
0
        if (considerProbabilities && probabilities_ != None) {
679
0
            std::vector<Real> prob(noEx+1);
680
0
            for (Size i = 0; i < noEx+1; i++) {
681
0
                prob[i] = npvp1[i][0] *
682
0
                          (probabilities_ == Naive
683
0
                               ? 1.0
684
0
                               : model_->numeraire(0.0, 0.0, discountCurve_));
685
0
            }
686
0
            results_.additionalResults["probabilities"] = prob;
687
0
        }
688
        // end probability computation
689
690
0
        return res;
691
0
    }
692
}