1 /* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2
3 /*
4 Copyright (C) 2013 Gary Kennedy
5 Copyright (C) 2015 Peter Caspers
6 Copyright (C) 2017 Klaus Spanderen
7 Copyright (C) 2020 Marcin Rybacki
8
9 This file is part of QuantLib, a free-software/open-source library
10 for financial quantitative analysts and developers - http://quantlib.org/
11
12 QuantLib is free software: you can redistribute it and/or modify it
13 under the terms of the QuantLib license. You should have received a
14 copy of the license along with this program; if not, please email
15 <quantlib-dev@lists.sf.net>. The license is also available online at
16 <http://quantlib.org/license.shtml>.
17
18 This program is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
20 FOR A PARTICULAR PURPOSE. See the license for more details.
21 */
22
23
24 #include "blackformula.hpp"
25 #include "utilities.hpp"
26 #include <ql/pricingengines/blackformula.hpp>
27
28 #include <boost/math/special_functions/fpclassify.hpp>
29
30 using namespace QuantLib;
31 using namespace boost::unit_test_framework;
32
33
testBachelierImpliedVol()34 void BlackFormulaTest::testBachelierImpliedVol(){
35
36
37 BOOST_TEST_MESSAGE("Testing Bachelier implied vol...");
38
39 Real forward = 1.0;
40 Real bpvol = 0.01;
41 Real tte = 10.0;
42 Real stdDev = bpvol*std::sqrt(tte);
43 Option::Type optionType = Option::Call;
44 Real discount = 0.95;
45
46 Real d[] = {-3.0, -2.0, -1.0, -0.5, 0.0, 0.5, 1.0, 2.0, 3.0};
47 for(Size i=0;i<LENGTH(d);++i){
48
49
50 Real strike = forward - d[i] * bpvol * std::sqrt(tte);
51
52 Real callPrem = bachelierBlackFormula(optionType, strike, forward, stdDev, discount);
53
54 Real impliedBpVol = bachelierBlackFormulaImpliedVol(optionType,strike, forward, tte, callPrem, discount);
55
56 if (std::fabs(bpvol-impliedBpVol)>1.0e-12){
57 BOOST_ERROR("Failed, expected " << bpvol << " realised " << impliedBpVol );
58 }
59 }
60 }
61
testChambersImpliedVol()62 void BlackFormulaTest::testChambersImpliedVol() {
63
64 BOOST_TEST_MESSAGE("Testing Chambers-Nawalkha implied vol approximation...");
65
66 Option::Type types[] = {Option::Call, Option::Put};
67 Real displacements[] = {0.0000, 0.0010, 0.0050, 0.0100, 0.0200};
68 Real forwards[] = {-0.0010, 0.0000, 0.0050, 0.0100, 0.0200, 0.0500};
69 Real strikes[] = {-0.0100, -0.0050, -0.0010, 0.0000, 0.0010, 0.0050,
70 0.0100, 0.0200, 0.0500, 0.1000};
71 Real stdDevs[] = {0.10, 0.15, 0.20, 0.30, 0.50, 0.60, 0.70,
72 0.80, 1.00, 1.50, 2.00};
73 Real discounts[] = {1.00, 0.95, 0.80, 1.10};
74
75 Real tol = 5.0E-4;
76
77 for (Size i1 = 0; i1 < LENGTH(types); ++i1) {
78 for (Size i2 = 0; i2 < LENGTH(displacements); ++i2) {
79 for (Size i3 = 0; i3 < LENGTH(forwards); ++i3) {
80 for (Size i4 = 0; i4 < LENGTH(strikes); ++i4) {
81 for (Size i5 = 0; i5 < LENGTH(stdDevs); ++i5) {
82 for (Size i6 = 0; i6 < LENGTH(discounts); ++i6) {
83 if (forwards[i3] + displacements[i2] > 0.0 &&
84 strikes[i4] + displacements[i2] > 0.0) {
85 Real premium = blackFormula(
86 types[i1], strikes[i4], forwards[i3],
87 stdDevs[i5], discounts[i6],
88 displacements[i2]);
89 Real atmPremium = blackFormula(
90 types[i1], forwards[i3], forwards[i3],
91 stdDevs[i5], discounts[i6],
92 displacements[i2]);
93 Real iStdDev =
94 blackFormulaImpliedStdDevChambers(
95 types[i1], strikes[i4], forwards[i3],
96 premium, atmPremium, discounts[i6],
97 displacements[i2]);
98 Real moneyness = (strikes[i4] + displacements[i2]) /
99 (forwards[i3] + displacements[i2]);
100 if(moneyness > 1.0) moneyness = 1.0 / moneyness;
101 Real error = (iStdDev - stdDevs[i5]) / stdDevs[i5] * moneyness;
102 if(error > tol)
103 BOOST_ERROR(
104 "Failed to verify Chambers-Nawalkha "
105 "approximation for "
106 << types[i1]
107 << " displacement=" << displacements[i2]
108 << " forward=" << forwards[i3]
109 << " strike=" << strikes[i4]
110 << " discount=" << discounts[i6]
111 << " stddev=" << stdDevs[i5]
112 << " result=" << iStdDev
113 << " exceeds maximum error tolerance");
114 }
115 }
116 }
117 }
118 }
119 }
120 }
121 }
122
testRadoicicStefanicaImpliedVol()123 void BlackFormulaTest::testRadoicicStefanicaImpliedVol() {
124
125 BOOST_TEST_MESSAGE(
126 "Testing Radoicic-Stefanica implied vol approximation...");
127
128 const Time T = 1.7;
129 const Rate r = 0.1;
130 const DiscountFactor df = std::exp(-r*T);
131
132 const Real forward = 100;
133
134 const Volatility vol = 0.3;
135 const Real stdDev = vol * std::sqrt(T);
136
137 const Option::Type types[] = { Option::Call, Option::Put };
138 const Real strikes[] = {
139 50, 60, 70, 80, 90, 100, 110, 125, 150, 200, 300 };
140
141 const Real tol = 0.02;
142
143 for (Size i=0; i < LENGTH(strikes); ++i) {
144 const Real strike = strikes[i];
145 for (Size j=0; j < LENGTH(types); ++j) {
146 const Option::Type type = types[j];
147
148 const ext::shared_ptr<PlainVanillaPayoff> payoff(
149 ext::make_shared<PlainVanillaPayoff>(type, strike));
150
151 const Real marketValue = blackFormula(payoff, forward, stdDev, df);
152
153 const Real estVol = blackFormulaImpliedStdDevApproximationRS(
154 payoff, forward, marketValue, df) / std::sqrt(T);
155
156 const Real error = std::fabs(estVol - vol);
157 if (error > tol) {
158 BOOST_ERROR("Failed to verify Radoicic-Stefanica"
159 "approximation for "
160 << type
161 << "\n forward :" << forward
162 << "\n strike :" << strike
163 << "\n discount :" << df
164 << "\n implied vol :" << vol
165 << "\n result :" << estVol
166 << "\n error :" << error
167 << "\n tolerance :" << tol);
168 }
169 }
170 }
171 }
172
testRadoicicStefanicaLowerBound()173 void BlackFormulaTest::testRadoicicStefanicaLowerBound() {
174
175 BOOST_TEST_MESSAGE("Testing Radoicic-Stefanica lower bound...");
176
177 // testing lower bound plot figure 3.1 from
178 // "Tighter Bounds for Implied Volatility",
179 // J. Gatheral, I. Matic, R. Radoicic, D. Stefanica
180 // https://papers.ssrn.com/sol3/papers.cfm?abstract_id=2922742
181
182 const Real forward = 1.0;
183 const Real k = 1.2;
184
185 for (Real s=0.17; s < 2.9; s+=0.01) {
186 const Real strike = std::exp(k)*forward;
187 const Real c = blackFormula(Option::Call, strike, forward, s);
188 const Real estimate = blackFormulaImpliedStdDevApproximationRS(
189 Option::Call, strike, forward, c);
190
191 const Real error = s - estimate;
192 if (boost::math::isnan(estimate) || std::fabs(error) > 0.05) {
193 BOOST_ERROR("Failed to lower bound Radoicic-Stefanica"
194 "approximation for "
195 << "\n forward :" << forward
196 << "\n strike :" << k
197 << "\n stdDev :" << s
198 << "\n result :" << estimate
199 << "\n error :" << error);
200
201 }
202
203 if (c > 1e-6 && error < 0.0) {
204 BOOST_ERROR("Failed to verify Radoicic-Stefanica is lower bound"
205 << "\n forward :" << forward
206 << "\n strike :" << k
207 << "\n stdDev :" << s
208 << "\n result :" << estimate
209 << "\n error :" << error);
210 }
211 }
212 }
213
testImpliedVolAdaptiveSuccessiveOverRelaxation()214 void BlackFormulaTest::testImpliedVolAdaptiveSuccessiveOverRelaxation() {
215 BOOST_TEST_MESSAGE("Testing implied volatility calculation via "
216 "adaptive successive over-relaxation...");
217
218 SavedSettings backup;
219
220 const DayCounter dc = Actual365Fixed();
221 const Date today = Date(12, July, 2017);
222 Settings::instance().evaluationDate() = today;
223
224 const Date exerciseDate = today + Period(15, Months);
225 const Time exerciseTime = dc.yearFraction(today, exerciseDate);
226
227 const ext::shared_ptr<YieldTermStructure> rTS = flatRate(0.10, dc);
228 const ext::shared_ptr<YieldTermStructure> qTS = flatRate(0.06, dc);
229
230 const DiscountFactor df = rTS->discount(exerciseDate);
231
232 const Volatility vol = 0.20;
233 const Real stdDev = vol * std::sqrt(exerciseTime);
234
235 const Real s0 = 100;
236 const Real forward= s0 * qTS->discount(exerciseDate)/df;
237
238 const Option::Type types[] = { Option::Call, Option::Put };
239 const Real strikes[] = { 50, 60, 70, 80, 90, 100, 110, 125, 150, 200 };
240 const Real displacements[] = { 0, 25, 50, 100};
241
242 const Real tol = 1e-8;
243
244 for (Size i=0; i < LENGTH(strikes); ++i) {
245 const Real strike = strikes[i];
246
247 for (Size j=0; j < LENGTH(types); ++j) {
248 const Option::Type type = types[j];
249
250 const ext::shared_ptr<PlainVanillaPayoff> payoff(
251 ext::make_shared<PlainVanillaPayoff>(type, strike));
252
253 for (Size k=0; k < LENGTH(displacements); ++k) {
254
255 const Real displacement = displacements[k];
256 const Real marketValue = blackFormula(
257 payoff, forward, stdDev, df, displacement);
258
259 const Real impliedStdDev = blackFormulaImpliedStdDevLiRS(
260 payoff, forward, marketValue, df, displacement,
261 Null<Real>(), 1.0, tol, 100);
262
263 const Real error = std::fabs(impliedStdDev - stdDev);
264 if (error > 10*tol) {
265 BOOST_ERROR("Failed to calculated implied volatility"
266 " with adaptive successive over-relaxation"
267 << "\n forward :" << forward
268 << "\n strike :" << strike
269 << "\n stdDev :" << stdDev
270 << "\n displacement:" << displacement
271 << "\n result :" << impliedStdDev
272 << "\n error :" << error
273 << "\n tolerance :" << tol);
274 }
275 }
276 }
277 }
278 }
279
assertBlackFormulaForwardDerivative(Option::Type optionType,const std::vector<Real> & strikes,Real bpvol)280 void assertBlackFormulaForwardDerivative(
281 Option::Type optionType,
282 const std::vector<Real> &strikes,
283 Real bpvol)
284 {
285 Real forward = 1.0;
286 Real tte = 10.0;
287 Real stdDev = bpvol * std::sqrt(tte);
288 Real discount = 0.95;
289 Real displacement = 0.01;
290 Real bump = 0.0001;
291 Real epsilon = 1.e-10;
292 std::string type = optionType == Option::Call ? "Call" : "Put";
293
294 for (std::vector<Real>::const_iterator it = strikes.begin(); it != strikes.end(); ++it)
295 {
296 Real strike = *it;
297 Real delta = blackFormulaForwardDerivative(
298 optionType, strike, forward, stdDev, discount, displacement);
299 Real bumpedDelta = blackFormulaForwardDerivative(
300 optionType, strike, forward + bump, stdDev, discount, displacement);
301
302 Real basePremium = blackFormula(
303 optionType, strike, forward, stdDev, discount, displacement);
304 Real bumpedPremium = blackFormula(
305 optionType, strike, forward + bump, stdDev, discount, displacement);
306 Real deltaApprox = (bumpedPremium - basePremium) / bump;
307
308 /*! Based on the Mean Value Theorem, the below inequality
309 should hold for any function that is monotonic in the
310 area of the bump.
311 */
312 bool success =
313 (std::max(delta, bumpedDelta) + epsilon > deltaApprox)
314 && (deltaApprox > std::min(delta, bumpedDelta) - epsilon);
315
316 if (!success)
317 {
318 BOOST_ERROR("Failed to calculate the derivative of the"
319 " Black formula w.r.t. forward"
320 << "\n option type :" << type
321 << "\n forward :" << forward
322 << "\n strike :" << strike
323 << "\n stdDev :" << stdDev
324 << "\n displacement :" << displacement
325 << "\n analytical delta :" << delta
326 << "\n approximated delta:" << deltaApprox);
327 }
328 }
329 }
330
testBlackFormulaForwardDerivative()331 void BlackFormulaTest::testBlackFormulaForwardDerivative() {
332
333 BOOST_TEST_MESSAGE("Testing forward derivative of the Black formula...");
334
335 std::vector<Real> strikes;
336 strikes.push_back(0.1);
337 strikes.push_back(0.5);
338 strikes.push_back(1.0);
339 strikes.push_back(2.0);
340 strikes.push_back(3.0);
341 const Real vol = 0.1;
342 assertBlackFormulaForwardDerivative(Option::Call, strikes, vol);
343 assertBlackFormulaForwardDerivative(Option::Put, strikes, vol);
344 }
345
testBlackFormulaForwardDerivativeWithZeroStrike()346 void BlackFormulaTest::testBlackFormulaForwardDerivativeWithZeroStrike() {
347
348 BOOST_TEST_MESSAGE("Testing forward derivative of the Black formula "
349 "with zero strike...");
350
351 std::vector<Real> strikes;
352 strikes.push_back(0.0);
353 const Real vol = 0.1;
354 assertBlackFormulaForwardDerivative(Option::Call, strikes, vol);
355 assertBlackFormulaForwardDerivative(Option::Put, strikes, vol);
356 }
357
testBlackFormulaForwardDerivativeWithZeroVolatility()358 void BlackFormulaTest::testBlackFormulaForwardDerivativeWithZeroVolatility() {
359
360 BOOST_TEST_MESSAGE("Testing forward derivative of the Black formula "
361 "with zero volatility...");
362
363 std::vector<Real> strikes;
364 strikes.push_back(0.1);
365 strikes.push_back(0.5);
366 strikes.push_back(1.0);
367 strikes.push_back(2.0);
368 strikes.push_back(3.0);
369 const Real vol = 0.0;
370 assertBlackFormulaForwardDerivative(Option::Call, strikes, vol);
371 assertBlackFormulaForwardDerivative(Option::Put, strikes, vol);
372 }
373
assertBachelierBlackFormulaForwardDerivative(Option::Type optionType,const std::vector<Real> & strikes,Real bpvol)374 void assertBachelierBlackFormulaForwardDerivative(
375 Option::Type optionType,
376 const std::vector<Real> &strikes,
377 Real bpvol)
378 {
379 Real forward = 1.0;
380 Real tte = 10.0;
381 Real stdDev = bpvol * std::sqrt(tte);
382 Real discount = 0.95;
383 Real bump = 0.0001;
384 Real epsilon = 1.e-10;
385 std::string type = optionType == Option::Call ? "Call" : "Put";
386
387 for (std::vector<Real>::const_iterator it = strikes.begin(); it != strikes.end(); ++it)
388 {
389 Real strike = *it;
390 Real delta = bachelierBlackFormulaForwardDerivative(
391 optionType, strike, forward, stdDev, discount);
392 Real bumpedDelta = bachelierBlackFormulaForwardDerivative(
393 optionType, strike, forward + bump, stdDev, discount);
394
395 Real basePremium = bachelierBlackFormula(
396 optionType, strike, forward, stdDev, discount);
397 Real bumpedPremium = bachelierBlackFormula(
398 optionType, strike, forward + bump, stdDev, discount);
399 Real deltaApprox = (bumpedPremium - basePremium) / bump;
400
401 /*! Based on the Mean Value Theorem, the below inequality
402 should hold for any function that is monotonic in the
403 area of the bump.
404 */
405 bool success =
406 (std::max(delta, bumpedDelta) + epsilon > deltaApprox)
407 && (deltaApprox > std::min(delta, bumpedDelta) - epsilon);
408
409 if (!success)
410 {
411 BOOST_ERROR("Failed to calculate the derivative of the"
412 " Bachelier Black formula w.r.t. forward"
413 << "\n option type :" << type
414 << "\n forward :" << forward
415 << "\n strike :" << strike
416 << "\n stdDev :" << stdDev
417 << "\n analytical delta :" << delta
418 << "\n approximated delta:" << deltaApprox);
419 }
420 }
421 }
422
testBachelierBlackFormulaForwardDerivative()423 void BlackFormulaTest::testBachelierBlackFormulaForwardDerivative() {
424
425 BOOST_TEST_MESSAGE("Testing forward derivative of the "
426 "Bachelier Black formula...");
427
428 std::vector<Real> strikes;
429 strikes.push_back(-3.0);
430 strikes.push_back(-2.0);
431 strikes.push_back(-1.0);
432 strikes.push_back(-0.5);
433 strikes.push_back(0.0);
434 strikes.push_back(0.5);
435 strikes.push_back(1.0);
436 strikes.push_back(2.0);
437 strikes.push_back(3.0);
438 const Real vol = 0.001;
439 assertBachelierBlackFormulaForwardDerivative(Option::Call, strikes, vol);
440 assertBachelierBlackFormulaForwardDerivative(Option::Put, strikes, vol);
441 }
442
testBachelierBlackFormulaForwardDerivativeWithZeroVolatility()443 void BlackFormulaTest::testBachelierBlackFormulaForwardDerivativeWithZeroVolatility() {
444
445 BOOST_TEST_MESSAGE("Testing forward derivative of the Bachelier Black formula "
446 "with zero volatility...");
447
448 std::vector<Real> strikes;
449 strikes.push_back(-3.0);
450 strikes.push_back(-2.0);
451 strikes.push_back(-1.0);
452 strikes.push_back(-0.5);
453 strikes.push_back(0.0);
454 strikes.push_back(0.5);
455 strikes.push_back(1.0);
456 strikes.push_back(2.0);
457 strikes.push_back(3.0);
458 const Real vol = 0.0;
459 assertBachelierBlackFormulaForwardDerivative(Option::Call, strikes, vol);
460 assertBachelierBlackFormulaForwardDerivative(Option::Put, strikes, vol);
461 }
462
suite()463 test_suite* BlackFormulaTest::suite() {
464 test_suite* suite = BOOST_TEST_SUITE("Black formula tests");
465
466 suite->add(QUANTLIB_TEST_CASE(
467 &BlackFormulaTest::testBachelierImpliedVol));
468 suite->add(QUANTLIB_TEST_CASE(
469 &BlackFormulaTest::testChambersImpliedVol));
470 suite->add(QUANTLIB_TEST_CASE(
471 &BlackFormulaTest::testRadoicicStefanicaImpliedVol));
472 suite->add(QUANTLIB_TEST_CASE(
473 &BlackFormulaTest::testRadoicicStefanicaLowerBound));
474 suite->add(QUANTLIB_TEST_CASE(
475 &BlackFormulaTest::testImpliedVolAdaptiveSuccessiveOverRelaxation));
476 suite->add(QUANTLIB_TEST_CASE(
477 &BlackFormulaTest::testBlackFormulaForwardDerivative));
478 suite->add(QUANTLIB_TEST_CASE(
479 &BlackFormulaTest::testBlackFormulaForwardDerivativeWithZeroStrike));
480 suite->add(QUANTLIB_TEST_CASE(
481 &BlackFormulaTest::testBlackFormulaForwardDerivativeWithZeroVolatility));
482 suite->add(QUANTLIB_TEST_CASE(
483 &BlackFormulaTest::testBachelierBlackFormulaForwardDerivative));
484 suite->add(QUANTLIB_TEST_CASE(
485 &BlackFormulaTest::testBachelierBlackFormulaForwardDerivativeWithZeroVolatility));
486
487 return suite;
488 }
489