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

source code of quantlib/test-suite/blackformula.cpp