24#include <ql/math/functional.hpp>
25#include <ql/methods/finitedifferences/meshers/fdmmesher.hpp>
26#include <ql/methods/finitedifferences/operators/fdmlinearopiterator.hpp>
27#include <ql/methods/finitedifferences/operators/fdmlinearoplayout.hpp>
28#include <ql/methods/finitedifferences/utilities/fdmhestongreensfct.hpp>
29#include <ql/methods/finitedifferences/utilities/squarerootprocessrndcalculator.hpp>
30#include <ql/processes/hestonprocess.hpp>
36 ext::shared_ptr<HestonProcess> process,
39 : l0_(l0), mesher_(
std::move(mesher)), process_(
std::move(process)), trafoType_(trafoType_) {}
48 const Real x0 = std::log(s0) + (r-q-0.5*v0*
l0_*
l0_)*t;
56 for (
const auto& iter : *
mesher_->layout()) {
60 : std::exp(
mesher_->location(iter, 1));
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);
76 retVal =
process_->pdf(x, v, t, 1e-4);
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)
86 - 2*rho*(x-x0)*(v-z0)/(sd_x*sd_v))
91 QL_FAIL(
"unknown algorithm");
101 retVal*=std::pow(v, 1.0 - 2*kappa*theta/(sigma*sigma));
104 QL_FAIL(
"unknown transformation type");
107 p[iter.index()] = retVal;
1-D array used in linear algebra.
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
Real pdf(Real v, Time t) const override
Real Time
continuous quantity with 1-year units