1 /* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2 
3 /*
4  Copyright (C) 2011 Klaus Spanderen
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 /*! \file expouwithjumpsprocess.cpp
21     \brief Ornstein Uhlenbeck process plus exp jumps (Kluge Model)
22 */
23 
24 #include <ql/experimental/processes/extouwithjumpsprocess.hpp>
25 #include <ql/experimental/processes/extendedornsteinuhlenbeckprocess.hpp>
26 
27 namespace QuantLib {
28 
ExtOUWithJumpsProcess(const ext::shared_ptr<ExtendedOrnsteinUhlenbeckProcess> & process,Real Y0,Real beta,Real jumpIntensity,Real eta)29     ExtOUWithJumpsProcess::ExtOUWithJumpsProcess(
30             const ext::shared_ptr<ExtendedOrnsteinUhlenbeckProcess>& process,
31             Real Y0, Real beta, Real jumpIntensity, Real eta)
32     : Y0_(Y0),
33       beta_(beta),
34       jumpIntensity_(jumpIntensity), eta_(eta),
35       ouProcess_(process) { }
36 
size() const37     Size ExtOUWithJumpsProcess::size() const {
38         return 2;
39     }
factors() const40     Size ExtOUWithJumpsProcess::factors() const {
41         return 3;
42     }
43     ext::shared_ptr<ExtendedOrnsteinUhlenbeckProcess>
getExtendedOrnsteinUhlenbeckProcess() const44         ExtOUWithJumpsProcess::getExtendedOrnsteinUhlenbeckProcess() const {
45         return ouProcess_;
46     }
beta() const47     Real ExtOUWithJumpsProcess::beta() const {
48         return beta_;
49     }
jumpIntensity() const50     Real ExtOUWithJumpsProcess::jumpIntensity() const {
51         return jumpIntensity_;
52     }
eta() const53     Real ExtOUWithJumpsProcess::eta() const {
54         return eta_;
55     }
56 
initialValues() const57     Disposable<Array> ExtOUWithJumpsProcess::initialValues() const {
58         Array retVal(2);
59         retVal[0] = ouProcess_->x0();
60         retVal[1] = Y0_;
61 
62         return retVal;
63     }
64 
drift(Time t,const Array & x) const65     Disposable<Array> ExtOUWithJumpsProcess::drift(Time t, const Array& x) const {
66         Array retVal(2);
67         retVal[0] = ouProcess_->drift(t, x[0]);
68         retVal[1] = -beta_*x[1];
69 
70         return retVal;
71     }
72 
73     Disposable<Matrix>
diffusion(Time t,const Array & x) const74     ExtOUWithJumpsProcess::diffusion(Time t, const Array& x) const {
75         Matrix retVal(2, 2, 0.0);
76         retVal[0][0] = ouProcess_->diffusion(t, x[0]);
77 
78         return retVal;
79     }
80 
evolve(Time t0,const Array & x0,Time dt,const Array & dw) const81     Disposable<Array> ExtOUWithJumpsProcess::evolve(
82         Time t0, const Array& x0, Time dt, const Array& dw) const {
83 
84         Array retVal(2);
85         retVal[0] = ouProcess_->evolve(t0, x0[0], dt, dw[0]);
86         retVal[1] = x0[1]*std::exp(-beta_*dt);
87 
88         const Real u1 = std::max(QL_EPSILON, std::min(cumNormalDist_(dw[1]),
89                                                       1.0-QL_EPSILON));
90 
91         const Time interarrival = -1.0/jumpIntensity_*std::log(u1);
92         if (interarrival < dt) {
93             const Real u2 = std::max(QL_EPSILON, std::min(cumNormalDist_(dw[2]),
94                                                           1.0-QL_EPSILON));
95             const Real jumpSize = -1.0/eta_*std::log(u2);
96             retVal[1] += jumpSize;
97         }
98         return retVal;
99     }
100 }
101