1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
|
/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
/*
Copyright (C) 2011 Klaus Spanderen
This file is part of QuantLib, a free-software/open-source library
for financial quantitative analysts and developers - http://quantlib.org/
QuantLib is free software: you can redistribute it and/or modify it
under the terms of the QuantLib license. You should have received a
copy of the license along with this program; if not, please email
<quantlib-dev@lists.sf.net>. The license is also available online at
<https://www.quantlib.org/license.shtml>.
This program is distributed in the hope that it will be useful, but WITHOUT
ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
FOR A PARTICULAR PURPOSE. See the license for more details.
*/
/*! \file gemanroncoroniprocess.cpp
\brief Geman-Roncoroni process
*/
#include <ql/math/functional.hpp>
#include <ql/processes/eulerdiscretization.hpp>
#include <ql/experimental/processes/gemanroncoroniprocess.hpp>
namespace QuantLib {
GemanRoncoroniProcess::GemanRoncoroniProcess(
Real x0,
Real alpha, Real beta,
Real gamma, Real delta,
Real eps, Real zeta, Real d,
Real k, Real tau,
Real sig2, Real a, Real b,
Real theta1, Real theta2, Real theta3,
Real psi)
: StochasticProcess1D(ext::shared_ptr<discretization>(
new EulerDiscretization)),
x0_(x0),
alpha_(alpha), beta_(beta),
gamma_(gamma), delta_(delta),
eps_(eps), zeta_(zeta), d_(d),
k_(k), tau_(tau),
sig2_(sig2), a_(a), b_(b),
theta1_(theta1), theta2_(theta2), theta3_(theta3),
psi_(psi) {
}
Real GemanRoncoroniProcess::x0() const {
return x0_;
}
Real GemanRoncoroniProcess::drift(Time t, Real x) const {
const Real mu = alpha_ + beta_*t + gamma_*std::cos(eps_+2*M_PI*t)
+ delta_*std::cos(zeta_+4*M_PI*t);
const Real muPrime = beta_ - gamma_*2*M_PI*std::sin(eps_+2*M_PI*t)
- delta_*4*M_PI*std::sin(zeta_+4*M_PI*t);
return muPrime + theta1_*(mu-x);
}
Real GemanRoncoroniProcess::diffusion(Time t, Real /*x*/) const {
return std::sqrt(sig2_ + a_*squared(std::cos(M_PI*t+b_)));
}
Real GemanRoncoroniProcess::stdDeviation(Time t0, Real /*x0*/, Time dt) const {
const Volatility sig2t = sig2_+a_*squared(std::cos(M_PI*t0+b_));
return std::sqrt(sig2t/(2*theta1_)*(1.0-std::exp(-2*theta1_*dt)));
}
Real GemanRoncoroniProcess::evolve(Time t0, Real x0,
Time dt, Real dw) const {
// random number generator for the jump part
if (!urng_) {
typedef PseudoRandom::urng_type urng_type;
urng_ = ext::make_shared<urng_type>((unsigned long)(1234UL * dw + 56789UL));
}
Array du(3);
du[0] = urng_->next().value;
du[1] = urng_->next().value;
return evolve(t0, x0, dt, dw, du);
}
Real GemanRoncoroniProcess::evolve(Time t0, Real x0, Time dt,
Real dw, const Array& du) const {
Real retVal;
const Time t = t0 + 0.5*dt;
const Real mu = alpha_ + beta_*t + gamma_*std::cos(eps_ +2*M_PI*t)
+ delta_*std::cos(zeta_+4*M_PI*t);
const Real j = -1.0/theta3_
*std::log(1.0+du[1]*(std::exp(-theta3_*psi_)-1.0));
if (x0 <= mu+d_) {
retVal = StochasticProcess1D::evolve(t, x0, dt, dw);
const Real jumpIntensity
= theta2_*(2.0/(1+std::fabs(std::sin(M_PI*(t-tau_)/k_)))-1);
const Time interarrival = -1.0/jumpIntensity*std::log(du[0]);
if (interarrival < dt) {
retVal += j;
}
}
else {
retVal = x0-j;
}
return retVal;
}
}
|
