22#include <ql/math/distributions/normaldistribution.hpp>
23#include <ql/math/integrals/segmentintegral.hpp>
24#include <ql/math/solvers1d/brent.hpp>
25#include <ql/models/shortrate/twofactormodels/g2.hpp>
26#include <ql/pricingengines/blackformula.hpp>
34 a_(arguments_[0]), sigma_(arguments_[1]), b_(arguments_[2]),
35 eta_(arguments_[3]), rho_(arguments_[4]) {
48 ext::shared_ptr<TwoFactorModel::ShortRateDynamics>
G2::dynamics()
const {
49 return ext::shared_ptr<ShortRateDynamics>(
new
60 Real temp = 1.0 - std::exp(-(
a()+
b())*t);
61 Real temp1 = 1.0 - std::exp(-
a()*(s-t));
62 Real temp2 = 1.0 - std::exp(-
b()*(s-t));
68 0.5*sigma2*temp1*temp1*(1.0 - std::exp(-2.0*
a()*t))/a3 +
69 0.5*eta2*temp2*temp2*(1.0 - std::exp(-2.0*
b()*t))/b3 +
72 return std::sqrt(
value);
76 return A(t,T) * std::exp(-
B(
a(),(T-t))*x-
B(
b(),(T-t))*y);
80 Time bondMaturity)
const {
90 Real expat = std::exp(-
a()*t);
91 Real expbt = std::exp(-
b()*t);
94 Real valuex = cx*cx*(t + (2.0*expat-0.5*expat*expat-1.5)/
a());
95 Real valuey = cy*cy*(t + (2.0*expbt-0.5*expbt*expbt-1.5)/
b());
98 - (expat*expbt-1.0)/(
a()+
b()));
99 return valuex + valuey +
value;
104 std::exp(0.5*(
V(T-t) -
V(T) +
V(t)));
108 return (1.0 - std::exp(-x*t))/x;
111 class G2::SwaptionPricingFunction {
113 SwaptionPricingFunction(
Real a,
120 std::vector<Time> payTimes,
123 : a_(
a), sigma_(
sigma), b_(
b), eta_(
eta), rho_(
rho), w_(w), T_(start),
124 t_(
std::move(payTimes)), rate_(fixedRate), size_(t_.size()), A_(size_), Ba_(size_),
128 sigmax_ = sigma_*std::sqrt(0.5*(1.0-std::exp(-2.0*a_*T_))/a_);
129 sigmay_ = eta_*std::sqrt(0.5*(1.0-std::exp(-2.0*b_*T_))/b_);
130 rhoxy_ = rho_*eta_*sigma_*(1.0 - std::exp(-(a_+b_)*T_))/
131 ((a_+b_)*sigmax_*sigmay_);
133 Real temp = sigma_*sigma_/(a_*a_);
134 mux_ = -((temp+rho_*sigma_*eta_/(a_*b_))*(1.0 - std::exp(-
a*T_)) -
135 0.5*temp*(1.0 - std::exp(-2.0*a_*T_)) -
136 rho_*sigma_*eta_/(b_*(a_+b_))*
137 (1.0- std::exp(-(b_+a_)*T_)));
139 temp = eta_*eta_/(b_*b_);
140 muy_ = -((temp+rho_*sigma_*eta_/(a_*b_))*(1.0 - std::exp(-
b*T_)) -
141 0.5*temp*(1.0 - std::exp(-2.0*b_*T_)) -
142 rho_*sigma_*eta_/(a_*(a_+b_))*
143 (1.0- std::exp(-(b_+a_)*T_)));
145 for (
Size i=0; i<size_; i++) {
146 A_[i] = model.
A(T_, t_[i]);
147 Ba_[i] = model.
B(a_, t_[i]-T_);
148 Bb_[i] = model.
B(b_, t_[i]-T_);
152 Real mux()
const {
return mux_; }
153 Real sigmax()
const {
return sigmax_; }
155 CumulativeNormalDistribution phi;
156 Real temp = (x - mux_)/sigmax_;
157 Real txy = std::sqrt(1.0 - rhoxy_*rhoxy_);
161 for (i=0; i<size_; i++) {
162 Real tau = (i==0 ? t_[0] - T_ : t_[i] - t_[i-1]);
163 Real c = (i==size_-1 ?
Real(1.0+rate_*tau) : rate_*tau);
164 lambda[i] = c*A_[i]*std::exp(-Ba_[i]*x);
167 SolvingFunction function(lambda, Bb_) ;
169 s1d.setMaxEvaluations(1000);
170 Real searchBound = std::max(10.0*sigmay_, 1.0);
171 Real yb = s1d.solve(function, 1e-6, 0.00, -searchBound, searchBound);
173 Real h1 = (yb - muy_)/(sigmay_*txy) -
174 rhoxy_*(x - mux_)/(sigmax_*txy);
178 for (i=0; i<size_; i++) {
180 Bb_[i]*sigmay_*std::sqrt(1.0-rhoxy_*rhoxy_);
181 Real kappa = - Bb_[i] *
182 (muy_ - 0.5*txy*txy*sigmay_*sigmay_*Bb_[i] +
183 rhoxy_*sigmay_*(x-mux_)/sigmax_);
184 value -= lambda[i] *std::exp(kappa)*phi(-w_*h2);
187 return std::exp(-0.5*temp*temp)*
value/
188 (sigmax_*std::sqrt(2.0*M_PI));
193 class SolvingFunction {
195 SolvingFunction(
const Array& lambda,
const Array& Bb)
196 : lambda_(lambda), Bb_(Bb) {}
199 for (
Size i=0; i<lambda_.size(); i++) {
200 value -= lambda_[i]*std::exp(-Bb_[i]*y);
205 const Array& lambda_;
209 Real a_, sigma_, b_, eta_, rho_, w_;
211 std::vector<Time> t_;
215 Real mux_, muy_, sigmax_, sigmay_, rhoxy_;
222 "non-constant nominals are not supported yet");
227 arguments.floatingResetDates[0]);
230 std::vector<Time> fixedPayTimes(arguments.fixedPayDates.size());
231 for (
Size i=0; i<fixedPayTimes.size(); ++i)
234 arguments.fixedPayDates[i]);
236 SwaptionPricingFunction function(
a(),
sigma(),
b(),
eta(),
rho(),
241 Real upper = function.mux() + range*function.sigmax();
242 Real lower = function.mux() - range*function.sigmax();
244 return arguments.nominal * w *
termStructure()->discount(start) *
245 integrator(function, lower, upper);
Constraint imposing all arguments to be in [low,high]
Real value(const Array ¶ms, const std::vector< ext::shared_ptr< CalibrationHelper > > &)
Standard constant parameter .
Time yearFraction(const Date &, const Date &, const Date &refPeriodStart=Date(), const Date &refPeriodEnd=Date()) const
Returns the period between two dates as a fraction of year.
Analytical term-structure fitting parameter .
Two-additive-factor gaussian model class.
Real discountBond(Time now, Time maturity, Array factors) const override
G2(const Handle< YieldTermStructure > &termStructure, Real a=0.1, Real sigma=0.01, Real b=0.1, Real eta=0.01, Real rho=-0.75)
Real swaption(const Swaption::arguments &arguments, Rate fixedRate, Real range, Size intervals) const
Real sigmaP(Time t, Time s) const
Real discountBondOption(Option::Type type, Real strike, Time maturity, Time bondMaturity) const override
void generateArguments() override
Real A(Time t, Time T) const
ext::shared_ptr< ShortRateDynamics > dynamics() const override
Returns the short-rate dynamics.
Real B(Real x, Time t) const
Shared handle to an observable.
template class providing a null value for a given type.
std::pair< iterator, bool > registerWith(const ext::shared_ptr< Observable > &)
Constraint imposing positivity to all arguments
Integral of a one-dimensional function.
Arguments for swaption calculation
Term-structure consistent model class.
const Handle< YieldTermStructure > & termStructure() const
Abstract base-class for two-factor models.
Real Time
continuous quantity with 1-year units
std::size_t Size
size of a container
Real blackFormula(Option::Type optionType, Real strike, Real forward, Real stdDev, Real discount, Real displacement)