1 /* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2 
3 /*
4  Copyright (C) 2015 Johannes Göttker-Schnetmann
5  Copyright (C) 2015 Klaus Spanderen
6 
7  This file is part of QuantLib, a free-software/open-source library
8  for financial quantitative analysts and developers - http://quantlib.org/
9 
10  QuantLib is free software: you can redistribute it and/or modify it
11  under the terms of the QuantLib license.  You should have received a
12  copy of the license along with this program; if not, please email
13  <quantlib-dev@lists.sf.net>. The license is also available online at
14  <http://quantlib.org/license.shtml>.
15 
16  This program is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
18  FOR A PARTICULAR PURPOSE.  See the license for more details.
19 */
20 
21 /*! \file hestonblackvolsurface.hpp
22     \brief Black volatility surface back by Heston model
23 */
24 
25 #include <ql/math/functional.hpp>
26 #include <ql/math/solvers1d/brent.hpp>
27 #include <ql/time/calendars/nullcalendar.hpp>
28 #include <ql/pricingengines/blackformula.hpp>
29 #include <ql/termstructures/volatility/equityfx/hestonblackvolsurface.hpp>
30 #include <ql/functional.hpp>
31 #include <boost/scoped_ptr.hpp>
32 #include <limits>
33 
34 namespace QuantLib {
35 
36     namespace {
blackValue(Option::Type optionType,Real strike,Real forward,Real maturity,Volatility vol,Real discount,Real npv)37         Real blackValue(Option::Type optionType, Real strike,
38                         Real forward, Real maturity,
39                         Volatility vol, Real discount, Real npv) {
40 
41             return blackFormula(optionType, strike, forward,
42                                 std::max(0.0, vol)*std::sqrt(maturity),
43                                 discount)-npv;
44         }
45     }
46 
HestonBlackVolSurface(const Handle<HestonModel> & hestonModel,const AnalyticHestonEngine::ComplexLogFormula cpxLogFormula,const AnalyticHestonEngine::Integration & integration)47     HestonBlackVolSurface::HestonBlackVolSurface(
48         const Handle<HestonModel>& hestonModel,
49         const AnalyticHestonEngine::ComplexLogFormula cpxLogFormula,
50         const AnalyticHestonEngine::Integration& integration)
51     : BlackVolTermStructure(
52           hestonModel->process()->riskFreeRate()->referenceDate(),
53           NullCalendar(),
54           Following,
55           hestonModel->process()->riskFreeRate()->dayCounter()),
56       hestonModel_(hestonModel),
57       cpxLogFormula_(cpxLogFormula),
58       integration_(integration) {
59         registerWith(hestonModel_);
60     }
61 
dayCounter() const62     DayCounter HestonBlackVolSurface::dayCounter() const {
63         return hestonModel_->process()->riskFreeRate()->dayCounter();
64     }
maxDate() const65     Date HestonBlackVolSurface::maxDate() const {
66         return Date::maxDate();
67     }
minStrike() const68     Real HestonBlackVolSurface::minStrike() const {
69         return 0.0;
70     }
maxStrike() const71     Real HestonBlackVolSurface::maxStrike() const {
72         return std::numeric_limits<Real>::max();
73     }
74 
blackVarianceImpl(Time t,Real strike) const75     Real HestonBlackVolSurface::blackVarianceImpl(Time t, Real strike) const {
76         return square<Real>()(blackVolImpl(t, strike))*t;
77     }
78 
blackVolImpl(Time t,Real strike) const79     Volatility HestonBlackVolSurface::blackVolImpl(Time t, Real strike) const {
80         using namespace ext::placeholders;
81         const ext::shared_ptr<HestonProcess> process = hestonModel_->process();
82 
83         const DiscountFactor df = process->riskFreeRate()->discount(t, true);
84         const DiscountFactor div = process->dividendYield()->discount(t, true);
85         const Real spotPrice = process->s0()->value();
86 
87         const Real fwd = spotPrice
88             * process->dividendYield()->discount(t, true)
89             / process->riskFreeRate()->discount(t, true);
90 
91 
92         const PlainVanillaPayoff payoff(
93             fwd > strike ? Option::Put : Option::Call, strike);
94 
95         const Real kappa = hestonModel_->kappa();
96         const Real theta = hestonModel_->theta();
97         const Real rho   = hestonModel_->rho();
98         const Real sigma = hestonModel_->sigma();
99         const Real v0    = hestonModel_->v0();
100 
101         const boost::scoped_ptr<AnalyticHestonEngine> hestonEngine(
102             new AnalyticHestonEngine(
103                 hestonModel_.currentLink(), cpxLogFormula_, integration_));
104 
105         Real npv;
106         Size evaluations;
107 
108         AnalyticHestonEngine::doCalculation(
109             df, div, spotPrice, strike, t,
110             kappa, theta, sigma, v0, rho,
111             payoff, integration_, cpxLogFormula_,
112             hestonEngine.get(), npv, evaluations);
113 
114         if (npv <= 0.0) return std::sqrt(theta);
115 
116         Brent solver;
117         solver.setMaxEvaluations(10000);
118         const Volatility guess = std::sqrt(theta);
119         const Real accuracy = std::numeric_limits<Real>::epsilon();
120 
121         const ext::function<Real(Real)> f = ext::bind(
122             &blackValue, payoff.optionType(), strike, fwd, t, _1, df, npv);
123 
124         return solver.solve(f, accuracy, guess, 0.01);
125     }
126 }
127