QuantLib: a free/open-source library for quantitative finance
fully annotated source code - version 1.34
Loading...
Searching...
No Matches
fdmhestongreensfct.cpp
Go to the documentation of this file.
1/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2
3/*
4 Copyright (C) 2014 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 fdmhestongreenfct.cpp
21 \brief Heston Fokker-Planck Green's function
22*/
23
31#include <utility>
32
33namespace QuantLib {
34
35 FdmHestonGreensFct::FdmHestonGreensFct(ext::shared_ptr<FdmMesher> mesher,
36 ext::shared_ptr<HestonProcess> process,
38 const Real l0)
39 : l0_(l0), mesher_(std::move(mesher)), process_(std::move(process)), trafoType_(trafoType_) {}
40
42
43 const Rate r = process_->riskFreeRate()->forwardRate(0, t, Continuous);
44 const Rate q = process_->dividendYield()->forwardRate(0,t, Continuous);
45
46 const Real s0 = process_->s0()->value();
47 const Real v0 = process_->v0();
48 const Real x0 = std::log(s0) + (r-q-0.5*v0*l0_*l0_)*t;
49
50 const Real rho = process_->rho();
51 const Real theta = process_->theta();
52 const Real kappa = process_->kappa();
53 const Real sigma = process_->sigma();
54
55 Array p(mesher_->layout()->size());
56 for (const auto& iter : *mesher_->layout()) {
57 const Real x = mesher_->location(iter, 0);
59 ? mesher_->location(iter, 1)
60 : std::exp(mesher_->location(iter, 1));
61
62 Real retVal;
63 switch (algorithm) {
64 case ZeroCorrelation:
65 {
66 const Real sd_x = l0_*std::sqrt(v0*t);
67 const Real p_x = M_1_SQRTPI*M_SQRT1_2/sd_x
68 * std::exp(-0.5*squared((x - x0)/sd_x));
70 v0, kappa, theta, sigma).pdf(v, t);
71
72 retVal = p_v*p_x;
73 }
74 break;
75 case SemiAnalytical:
76 retVal = process_->pdf(x, v, t, 1e-4);
77 break;
78 case Gaussian:
79 {
80 const Real sd_x = l0_*std::sqrt(v0*t);
81 const Real sd_v = sigma*std::sqrt(v0*t);
82 const Real z0 = v0 + kappa*(theta - v0)*t;
83 retVal = 1.0/(M_TWOPI*sd_x*sd_v*std::sqrt(1-rho*rho))
84 *std::exp(-( squared((x-x0)/sd_x)
85 + squared((v-z0)/sd_v)
86 - 2*rho*(x-x0)*(v-z0)/(sd_x*sd_v))
87 /(2*(1-rho*rho)) );
88 }
89 break;
90 default:
91 QL_FAIL("unknown algorithm");
92 }
93
94 switch (trafoType_) {
96 retVal*=v;
97 break;
99 break;
101 retVal*=std::pow(v, 1.0 - 2*kappa*theta/(sigma*sigma));
102 break;
103 default:
104 QL_FAIL("unknown transformation type");
105 }
106
107 p[iter.index()] = retVal;
108 }
109
110 return p;
111 }
112}
1-D array used in linear algebra.
Definition: array.hpp:52
const ext::shared_ptr< FdmMesher > mesher_
const FdmSquareRootFwdOp::TransformationType trafoType_
FdmHestonGreensFct(ext::shared_ptr< FdmMesher > mesher, ext::shared_ptr< HestonProcess > process, FdmSquareRootFwdOp::TransformationType trafoType_, Real l0=1.0)
const ext::shared_ptr< HestonProcess > process_
Array get(Time t, Algorithm algorithm) const
const DefaultType & t
#define QL_FAIL(message)
throw an error (possibly with file and line information)
Definition: errors.hpp:92
Heston Fokker-Planck Green's function.
iterator for a linear fdm operator
memory layout of a fdm linear operator
mesher for a fdm grid
const ext::shared_ptr< FdmMesher > mesher_
Real Time
continuous quantity with 1-year units
Definition: types.hpp:62
QL_REAL Real
real number
Definition: types.hpp:50
Real Rate
interest rates
Definition: types.hpp:70
Heston stochastic process.
Real kappa
Real theta
Real v0
Real rho
Real sigma
functionals and combinators not included in the STL
#define M_TWOPI
#define M_1_SQRTPI
#define M_SQRT1_2
Definition: any.hpp:35
T squared(T x)
Definition: functional.hpp:37
STL namespace.
ext::shared_ptr< YieldTermStructure > q
ext::shared_ptr< YieldTermStructure > r
ext::shared_ptr< BlackVolTermStructure > v
risk neutral terminal density calculator for the square root process