[go: up one dir, main page]
1/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2
3/*
4 Copyright (C) 2006 StatPro Italia srl
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 "brownianbridge.hpp"
21#include "utilities.hpp"
22#include <ql/methods/montecarlo/brownianbridge.hpp>
23#include <ql/methods/montecarlo/pathgenerator.hpp>
24#include <ql/math/randomnumbers/sobolrsg.hpp>
25#include <ql/math/randomnumbers/inversecumulativersg.hpp>
26#include <ql/math/statistics/sequencestatistics.hpp>
27#include <ql/processes/blackscholesprocess.hpp>
28#include <ql/termstructures/yield/flatforward.hpp>
29#include <ql/termstructures/volatility/equityfx/blackconstantvol.hpp>
30#include <ql/time/calendars/nullcalendar.hpp>
31
32using namespace QuantLib;
33using namespace boost::unit_test_framework;
34
35namespace {
36
37 template <class ForwardIterator1, class ForwardIterator2>
38 Real maxDiff(ForwardIterator1 begin1, ForwardIterator1 end1,
39 ForwardIterator2 begin2) {
40 Real diff = 0.0;
41 while (begin1 != end1) {
42 diff = std::max(diff, std::fabs(*begin1 - *begin2));
43 ++begin1; ++begin2;
44 }
45 return diff;
46 }
47
48 template <class ForwardIterator1, class ForwardIterator2>
49 Real maxRelDiff(ForwardIterator1 begin1, ForwardIterator1 end1,
50 ForwardIterator2 begin2) {
51 Real diff = 0.0;
52 while (begin1 != end1) {
53 diff = std::max(diff, std::fabs((*begin1 - *begin2)/(*begin2)));
54 ++begin1; ++begin2;
55 }
56 return diff;
57 }
58
59}
60
61
62void BrownianBridgeTest::testVariates() {
63 BOOST_TEST_MESSAGE("Testing Brownian-bridge variates...");
64
65 std::vector<Time> times = {0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 2.0, 5.0};
66
67 Size N = times.size();
68
69 Size samples = 262143;
70 unsigned long seed = 42;
71 SobolRsg sobol(N, seed);
72 InverseCumulativeRsg<SobolRsg,InverseCumulativeNormal> generator(sobol);
73
74 BrownianBridge bridge(times);
75
76 SequenceStatistics stats1(N);
77 SequenceStatistics stats2(N);
78
79 std::vector<Real> temp(N);
80
81 for (Size i=0; i<samples; ++i) {
82 const std::vector<Real>& sample = generator.nextSequence().value;
83
84 bridge.transform(begin: sample.begin(), end: sample.end(), output: temp.begin());
85 stats1.add(begin: temp.begin(), end: temp.end());
86
87 temp[0] = temp[0]*std::sqrt(x: times[0]);
88 for (Size j=1; j<N; ++j)
89 temp[j] = temp[j-1] + temp[j]*std::sqrt(x: times[j]-times[j-1]);
90 stats2.add(begin: temp.begin(), end: temp.end());
91 }
92
93 // normalized single variates
94 std::vector<Real> expectedMean(N, 0.0);
95 Matrix expectedCovariance(N, N, 0.0);
96 for (Size i=0; i<N; i++)
97 expectedCovariance[i][i] = 1.0;
98
99#ifndef __FAST_MATH__
100 Real meanTolerance = 1.0e-16;
101#else
102 Real meanTolerance = 1.0e-14;
103#endif
104 Real covTolerance = 2.5e-4;
105
106 std::vector<Real> mean = stats1.mean();
107 Matrix covariance = stats1.covariance();
108
109 Real maxMeanError = maxDiff(begin1: mean.begin(), end1: mean.end(),
110 begin2: expectedMean.begin());
111 Real maxCovError = maxDiff(begin1: covariance.begin(), end1: covariance.end(),
112 begin2: expectedCovariance.begin());
113
114 if (maxMeanError > meanTolerance) {
115 Array calculated(N), expected(N);
116 std::copy(first: mean.begin(), last: mean.end(), result: calculated.begin());
117 std::copy(first: expectedMean.begin(), last: expectedMean.end(), result: expected.begin());
118 BOOST_ERROR("failed to reproduce expected mean values"
119 << "\n calculated: " << calculated
120 << "\n expected: " << expected
121 << "\n max error: " << maxMeanError);
122 }
123
124 if (maxCovError > covTolerance) {
125 BOOST_ERROR("failed to reproduce expected covariance\n"
126 << " calculated:\n" << covariance
127 << " expected:\n" << expectedCovariance
128 << " max error: " << maxCovError);
129 }
130
131 // denormalized sums along the path
132 expectedMean = std::vector<Real>(N, 0.0);
133 expectedCovariance = Matrix(N, N);
134 for (Size i=0; i<N; ++i)
135 for (Size j=i; j<N; ++j)
136 expectedCovariance[i][j] = expectedCovariance[j][i] = times[i];
137
138 covTolerance = 6.0e-4;
139
140 mean = stats2.mean();
141 covariance = stats2.covariance();
142
143 maxMeanError = maxDiff(begin1: mean.begin(), end1: mean.end(),
144 begin2: expectedMean.begin());
145 maxCovError = maxDiff(begin1: covariance.begin(), end1: covariance.end(),
146 begin2: expectedCovariance.begin());
147
148 if (maxMeanError > meanTolerance) {
149 Array calculated(N), expected(N);
150 std::copy(first: mean.begin(), last: mean.end(), result: calculated.begin());
151 std::copy(first: expectedMean.begin(), last: expectedMean.end(), result: expected.begin());
152 BOOST_ERROR("failed to reproduce expected mean values"
153 << "\n calculated: " << calculated
154 << "\n expected: " << expected
155 << "\n max error: " << maxMeanError);
156 }
157
158 if (maxCovError > covTolerance) {
159 BOOST_ERROR("failed to reproduce expected covariance\n"
160 << " calculated:\n" << covariance
161 << " expected:\n" << expectedCovariance
162 << " max error: " << maxCovError);
163 }
164}
165
166
167void BrownianBridgeTest::testPathGeneration() {
168 BOOST_TEST_MESSAGE("Testing Brownian-bridge path generation...");
169
170 std::vector<Time> times = {0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 2.0, 5.0, 7.0, 9.0, 10.0};
171
172 TimeGrid grid(times.begin(), times.end());
173
174 Size N = times.size();
175
176 Size samples = 131071;
177 unsigned long seed = 42;
178 SobolRsg sobol(N, seed);
179 InverseCumulativeRsg<SobolRsg,InverseCumulativeNormal> gsg(sobol);
180
181 Date today = Settings::instance().evaluationDate();
182 Handle<Quote> x0(ext::shared_ptr<Quote>(new SimpleQuote(100.0)));
183 Handle<YieldTermStructure> r(ext::shared_ptr<YieldTermStructure>(
184 new FlatForward(today,0.06,Actual365Fixed())));
185 Handle<YieldTermStructure> q(ext::shared_ptr<YieldTermStructure>(
186 new FlatForward(today,0.03,Actual365Fixed())));
187 Handle<BlackVolTermStructure> sigma(
188 ext::shared_ptr<BlackVolTermStructure>(
189 new BlackConstantVol(today, NullCalendar(), 0.20,Actual365Fixed())));
190
191 ext::shared_ptr<StochasticProcess1D> process(
192 new BlackScholesMertonProcess(x0, q, r, sigma));
193
194
195 PathGenerator<InverseCumulativeRsg<SobolRsg,InverseCumulativeNormal> >
196 generator1(process, grid, gsg, false);
197 PathGenerator<InverseCumulativeRsg<SobolRsg,InverseCumulativeNormal> >
198 generator2(process, grid, gsg, true);
199
200 SequenceStatistics stats1(N);
201 SequenceStatistics stats2(N);
202
203 std::vector<Real> temp(N);
204
205 for (Size i=0; i<samples; ++i) {
206 const Path& path1 = generator1.next().value;
207 std::copy(first: path1.begin()+1, last: path1.end(), result: temp.begin());
208 stats1.add(begin: temp.begin(), end: temp.end());
209
210 const Path& path2 = generator2.next().value;
211 std::copy(first: path2.begin()+1, last: path2.end(), result: temp.begin());
212 stats2.add(begin: temp.begin(), end: temp.end());
213 }
214
215 std::vector<Real> expectedMean = stats1.mean();
216 Matrix expectedCovariance = stats1.covariance();
217
218 std::vector<Real> mean = stats2.mean();
219 Matrix covariance = stats2.covariance();
220
221 Real meanTolerance = 3.0e-5;
222 Real covTolerance = 3.0e-3;
223
224 Real maxMeanError = maxRelDiff(begin1: mean.begin(), end1: mean.end(),
225 begin2: expectedMean.begin());
226 Real maxCovError = maxRelDiff(begin1: covariance.begin(), end1: covariance.end(),
227 begin2: expectedCovariance.begin());
228
229 if (maxMeanError > meanTolerance) {
230 Array calculated(N), expected(N);
231 std::copy(first: mean.begin(), last: mean.end(), result: calculated.begin());
232 std::copy(first: expectedMean.begin(), last: expectedMean.end(), result: expected.begin());
233 BOOST_ERROR("failed to reproduce expected mean values"
234 << "\n calculated: " << calculated
235 << "\n expected: " << expected
236 << "\n max error: " << maxMeanError);
237 }
238
239 if (maxCovError > covTolerance) {
240 BOOST_ERROR("failed to reproduce expected covariance\n"
241 << " calculated:\n" << covariance
242 << " expected:\n" << expectedCovariance
243 << " max error: " << maxCovError);
244 }
245}
246
247test_suite* BrownianBridgeTest::suite() {
248 auto* suite = BOOST_TEST_SUITE("Brownian bridge tests");
249 suite->add(QUANTLIB_TEST_CASE(&BrownianBridgeTest::testVariates));
250 suite->add(QUANTLIB_TEST_CASE(&BrownianBridgeTest::testPathGeneration));
251 return suite;
252}
253
254

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