1/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2/*
3 Copyright (C) 2014, 2015, 2018 Peter Caspers
4
5 This file is part of QuantLib, a free-software/open-source library
6 for financial quantitative analysts and developers - http://quantlib.org/
7
8 QuantLib is free software: you can redistribute it and/or modify it
9 under the terms of the QuantLib license. You should have received a
10 copy of the license along with this program; if not, please email
11 <quantlib-dev@lists.sf.net>. The license is also available online at
12 <http://quantlib.org/license.shtml>.
13
14
15 This program is distributed in the hope that it will be useful, but
16 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
17 or FITNESS FOR A PARTICULAR PURPOSE. See the license for more details.
18*/
19
20/*! \file lognormalcmsspreadpricer.cpp
21*/
22
23#include <ql/experimental/coupons/cmsspreadcoupon.hpp>
24#include <ql/experimental/coupons/lognormalcmsspreadpricer.hpp>
25#include <ql/math/integrals/kronrodintegral.hpp>
26#include <ql/pricingengines/blackformula.hpp>
27#include <ql/termstructures/volatility/swaption/swaptionvolcube.hpp>
28#include <ql/optional.hpp>
29#include <utility>
30
31using std::sqrt;
32
33namespace QuantLib {
34
35 class LognormalCmsSpreadPricer::integrand_f {
36 const LognormalCmsSpreadPricer* pricer;
37 public:
38 explicit integrand_f(const LognormalCmsSpreadPricer* pricer)
39 : pricer(pricer) {}
40 Real operator()(Real x) const {
41 return pricer->integrand(x);
42 }
43 };
44
45 LognormalCmsSpreadPricer::LognormalCmsSpreadPricer(
46 const ext::shared_ptr<CmsCouponPricer>& cmsPricer,
47 const Handle<Quote>& correlation,
48 Handle<YieldTermStructure> couponDiscountCurve,
49 const Size integrationPoints,
50 const ext::optional<VolatilityType>& volatilityType,
51 const Real shift1,
52 const Real shift2)
53 : CmsSpreadCouponPricer(correlation), cmsPricer_(cmsPricer),
54 couponDiscountCurve_(std::move(couponDiscountCurve)) {
55
56 registerWith(h: correlation);
57 if (!couponDiscountCurve_.empty())
58 registerWith(h: couponDiscountCurve_);
59 registerWith(h: cmsPricer_);
60
61 QL_REQUIRE(integrationPoints >= 4,
62 "at least 4 integration points should be used ("
63 << integrationPoints << ")");
64 integrator_ =
65 ext::make_shared<GaussHermiteIntegration>(args: integrationPoints);
66
67 cnd_ = ext::make_shared<CumulativeNormalDistribution>(args: 0.0, args: 1.0);
68
69 if (!volatilityType) {
70 QL_REQUIRE(shift1 == Null<Real>() && shift2 == Null<Real>(),
71 "if volatility type is inherited, no shifts should be "
72 "specified");
73 inheritedVolatilityType_ = true;
74 volType_ = cmsPricer->swaptionVolatility()->volatilityType();
75 } else {
76 shift1_ = shift1 == Null<Real>() ? 0.0 : shift1;
77 shift2_ = shift2 == Null<Real>() ? 0.0 : shift2;
78 inheritedVolatilityType_ = false;
79 volType_ = *volatilityType;
80 }
81 }
82
83 Real LognormalCmsSpreadPricer::integrand(const Real x) const {
84
85 // this is Brigo, 13.16.2 with x = v / sqrt(2)
86
87 Real v = M_SQRT2 * x;
88 Real h =
89 k_ - b_ * s2_ * std::exp(x: (m2_ - 0.5 * v2_ * v2_) * fixingTime_ +
90 v2_ * std::sqrt(x: fixingTime_) * v);
91 Real phi1, phi2;
92 phi1 = (*cnd_)(
93 phi_ * (std::log(x: a_ * s1_ / h) +
94 (m1_ + (0.5 - rho_ * rho_) * v1_ * v1_) * fixingTime_ +
95 rho_ * v1_ * std::sqrt(x: fixingTime_) * v) /
96 (v1_ * std::sqrt(x: fixingTime_ * (1.0 - rho_ * rho_))));
97 phi2 = (*cnd_)(
98 phi_ * (std::log(x: a_ * s1_ / h) +
99 (m1_ - 0.5 * v1_ * v1_) * fixingTime_ +
100 rho_ * v1_ * std::sqrt(x: fixingTime_) * v) /
101 (v1_ * std::sqrt(x: fixingTime_ * (1.0 - rho_ * rho_))));
102 Real f = a_ * phi_ * s1_ *
103 std::exp(x: m1_ * fixingTime_ -
104 0.5 * rho_ * rho_ * v1_ * v1_ * fixingTime_ +
105 rho_ * v1_ * std::sqrt(x: fixingTime_) * v) *
106 phi1 -
107 phi_ * h * phi2;
108 return std::exp(x: -x * x) * f;
109 }
110
111 Real LognormalCmsSpreadPricer::integrand_normal(const Real x) const {
112
113 // this is http://ssrn.com/abstract=2686998, 3.20 with x = s / sqrt(2)
114
115 Real s = M_SQRT2 * x;
116
117 Real beta =
118 phi_ *
119 (gearing1_ * adjustedRate1_ + gearing2_ * adjustedRate2_ - k_ +
120 std::sqrt(x: fixingTime_) *
121 (rho_ * gearing1_ * vol1_ + gearing2_ * vol2_) * s);
122 Real f =
123 close_enough(x: alpha_, y: 0.0)
124 ? Real(std::max(a: beta, b: 0.0))
125 : psi_ * alpha_ / (M_SQRTPI * M_SQRT2) *
126 std::exp(x: -beta * beta / (2.0 * alpha_ * alpha_)) +
127 beta * (1.0 - (*cnd_)(-psi_ * beta / alpha_));
128 return std::exp(x: -x * x) * f;
129 }
130
131 void
132 LognormalCmsSpreadPricer::initialize(const FloatingRateCoupon &coupon) {
133
134 coupon_ = dynamic_cast<const CmsSpreadCoupon *>(&coupon);
135 QL_REQUIRE(coupon_, "CMS spread coupon needed");
136 index_ = coupon_->swapSpreadIndex();
137 gearing_ = coupon_->gearing();
138 spread_ = coupon_->spread();
139
140 fixingDate_ = coupon_->fixingDate();
141 paymentDate_ = coupon_->date();
142
143 // if no coupon discount curve is given just use the discounting curve
144 // from the _first_ swap index.
145 // for rate calculation this curve cancels out in the computation, so
146 // e.g. the discounting
147 // swap engine will produce correct results, even if the
148 // couponDiscountCurve is not set here.
149 // only the price member function in this class will be dependent on the
150 // coupon discount curve.
151
152 today_ = QuantLib::Settings::instance().evaluationDate();
153
154 if (couponDiscountCurve_.empty())
155 couponDiscountCurve_ =
156 index_->swapIndex1()->exogenousDiscount()
157 ? index_->swapIndex1()->discountingTermStructure()
158 : index_->swapIndex1()->forwardingTermStructure();
159
160 discount_ = paymentDate_ > couponDiscountCurve_->referenceDate()
161 ? couponDiscountCurve_->discount(d: paymentDate_)
162 : 1.0;
163
164 spreadLegValue_ = spread_ * coupon_->accrualPeriod() * discount_;
165
166 gearing1_ = index_->gearing1();
167 gearing2_ = index_->gearing2();
168
169 QL_REQUIRE(gearing1_ > 0.0 && gearing2_ < 0.0,
170 "gearing1 (" << gearing1_
171 << ") should be positive while gearing2 ("
172 << gearing2_ << ") should be negative");
173
174 c1_ = ext::make_shared<CmsCoupon>(
175 args: coupon_->date(), args: coupon_->nominal(), args: coupon_->accrualStartDate(),
176 args: coupon_->accrualEndDate(), args: coupon_->fixingDays(),
177 args: index_->swapIndex1(), args: 1.0, args: 0.0, args: coupon_->referencePeriodStart(),
178 args: coupon_->referencePeriodEnd(), args: coupon_->dayCounter(),
179 args: coupon_->isInArrears());
180
181 c2_ = ext::make_shared<CmsCoupon>(
182 args: coupon_->date(), args: coupon_->nominal(), args: coupon_->accrualStartDate(),
183 args: coupon_->accrualEndDate(), args: coupon_->fixingDays(),
184 args: index_->swapIndex2(), args: 1.0, args: 0.0, args: coupon_->referencePeriodStart(),
185 args: coupon_->referencePeriodEnd(), args: coupon_->dayCounter(),
186 args: coupon_->isInArrears());
187
188 c1_->setPricer(cmsPricer_);
189 c2_->setPricer(cmsPricer_);
190
191 if (fixingDate_ > today_) {
192
193 fixingTime_ = cmsPricer_->swaptionVolatility()->timeFromReference(
194 d: fixingDate_);
195
196 swapRate1_ = c1_->indexFixing();
197 swapRate2_ = c2_->indexFixing();
198
199 adjustedRate1_ = c1_->adjustedFixing();
200 adjustedRate2_ = c2_->adjustedFixing();
201
202 ext::shared_ptr<SwaptionVolatilityStructure> swvol =
203 *cmsPricer_->swaptionVolatility();
204 ext::shared_ptr<SwaptionVolatilityCube> swcub =
205 ext::dynamic_pointer_cast<SwaptionVolatilityCube>(r: swvol);
206
207 if(inheritedVolatilityType_ && volType_ == ShiftedLognormal) {
208 shift1_ =
209 swvol->shift(optionDate: fixingDate_, swapTenor: index_->swapIndex1()->tenor());
210 shift2_ =
211 swvol->shift(optionDate: fixingDate_, swapTenor: index_->swapIndex2()->tenor());
212 }
213
214 if (swcub == nullptr) {
215 // not a cube, just an atm surface given, so we can
216 // not easily convert volatilities and just forbid it
217 QL_REQUIRE(inheritedVolatilityType_,
218 "if only an atm surface is given, the volatility "
219 "type must be inherited");
220 vol1_ = swvol->volatility(
221 optionDate: fixingDate_, swapTenor: index_->swapIndex1()->tenor(), strike: swapRate1_);
222 vol2_ = swvol->volatility(
223 optionDate: fixingDate_, swapTenor: index_->swapIndex2()->tenor(), strike: swapRate2_);
224 } else {
225 vol1_ = swcub->smileSection(optionDate: fixingDate_,
226 swapTenor: index_->swapIndex1()->tenor())
227 ->volatility(strike: swapRate1_, type: volType_, shift: shift1_);
228 vol2_ = swcub->smileSection(optionDate: fixingDate_,
229 swapTenor: index_->swapIndex2()->tenor())
230 ->volatility(strike: swapRate2_, type: volType_, shift: shift2_);
231 }
232
233 if(volType_ == ShiftedLognormal) {
234 mu1_ = 1.0 / fixingTime_ * std::log(x: (adjustedRate1_ + shift1_) /
235 (swapRate1_ + shift1_));
236 mu2_ = 1.0 / fixingTime_ * std::log(x: (adjustedRate2_ + shift2_) /
237 (swapRate2_ + shift2_));
238 }
239 // for the normal volatility case we do not need the drifts
240 // but rather use adjusted rates directly in the integrand
241
242 rho_ = std::max(a: std::min(a: correlation()->value(), b: 0.9999),
243 b: -0.9999); // avoid division by zero in integrand
244 } else {
245 // fixing is in the past or today
246 adjustedRate1_ = c1_->indexFixing();
247 adjustedRate2_ = c2_->indexFixing();
248 }
249 }
250
251 Real LognormalCmsSpreadPricer::optionletPrice(Option::Type optionType,
252 Real strike) const {
253 // this method is only called for future fixings
254 optionType_ = optionType;
255 phi_ = optionType == Option::Call ? 1.0 : -1.0;
256 Real res = 0.0;
257 if (volType_ == ShiftedLognormal) {
258 // (shifted) lognormal volatility
259 if (strike >= 0.0) {
260 a_ = gearing1_;
261 b_ = gearing2_;
262 s1_ = swapRate1_ + shift1_;
263 s2_ = swapRate2_ + shift2_;
264 m1_ = mu1_;
265 m2_ = mu2_;
266 v1_ = vol1_;
267 v2_ = vol2_;
268 k_ = strike + gearing1_ * shift1_ + gearing2_ * shift2_;
269 } else {
270 a_ = -gearing2_;
271 b_ = -gearing1_;
272 s1_ = swapRate2_ + shift1_;
273 s2_ = swapRate1_ + shift2_;
274 m1_ = mu2_;
275 m2_ = mu1_;
276 v1_ = vol2_;
277 v2_ = vol1_;
278 k_ = -strike - gearing1_ * shift1_ - gearing2_ * shift2_;
279 res += phi_ * (gearing1_ * adjustedRate1_ +
280 gearing2_ * adjustedRate2_ - strike);
281 }
282 res +=
283 1.0 / M_SQRTPI * (*integrator_)(integrand_f(this));
284 } else {
285 // normal volatility
286 Real forward = gearing1_ * adjustedRate1_ +
287 gearing2_ * adjustedRate2_;
288 Real stddev =
289 std::sqrt(x: fixingTime_ *
290 (gearing1_ * gearing1_ * vol1_ * vol1_ +
291 gearing2_ * gearing2_ * vol2_ * vol2_ +
292 2.0 * gearing1_ * gearing2_ * rho_ * vol1_ * vol2_));
293 res =
294 bachelierBlackFormula(optionType: optionType_, strike, forward, stdDev: stddev, discount: 1.0);
295 }
296 return res * discount_ * coupon_->accrualPeriod();
297 }
298
299 Rate LognormalCmsSpreadPricer::swapletRate() const {
300 return swapletPrice() / (coupon_->accrualPeriod() * discount_);
301 }
302
303 Real LognormalCmsSpreadPricer::capletPrice(Rate effectiveCap) const {
304 // caplet is equivalent to call option on fixing
305 if (fixingDate_ <= today_) {
306 // the fixing is determined
307 const Rate Rs = std::max(
308 a: coupon_->index()->fixing(fixingDate: fixingDate_) - effectiveCap, b: 0.);
309 Rate price = gearing_ * Rs * coupon_->accrualPeriod() * discount_;
310 return price;
311 } else {
312 Real capletPrice = optionletPrice(optionType: Option::Call, strike: effectiveCap);
313 return gearing_ * capletPrice;
314 }
315 }
316
317 Rate LognormalCmsSpreadPricer::capletRate(Rate effectiveCap) const {
318 return capletPrice(effectiveCap) /
319 (coupon_->accrualPeriod() * discount_);
320 }
321
322 Real LognormalCmsSpreadPricer::floorletPrice(Rate effectiveFloor) const {
323 // floorlet is equivalent to put option on fixing
324 if (fixingDate_ <= today_) {
325 // the fixing is determined
326 const Rate Rs = std::max(
327 a: effectiveFloor - coupon_->index()->fixing(fixingDate: fixingDate_), b: 0.);
328 Rate price = gearing_ * Rs * coupon_->accrualPeriod() * discount_;
329 return price;
330 } else {
331 Real floorletPrice = optionletPrice(optionType: Option::Put, strike: effectiveFloor);
332 return gearing_ * floorletPrice;
333 }
334 }
335
336 Rate LognormalCmsSpreadPricer::floorletRate(Rate effectiveFloor) const {
337 return floorletPrice(effectiveFloor) /
338 (coupon_->accrualPeriod() * discount_);
339 }
340
341 Real LognormalCmsSpreadPricer::swapletPrice() const {
342 return gearing_ * coupon_->accrualPeriod() * discount_ *
343 (gearing1_ * adjustedRate1_ + gearing2_ * adjustedRate2_) +
344 spreadLegValue_;
345 }
346}
347

source code of quantlib/ql/experimental/coupons/lognormalcmsspreadpricer.cpp