1 /* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2 
3 /*
4  Copyright (C) 2014 Peter Caspers
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 zabrsmilesection.hpp
21     \brief zabr smile section
22 */
23 
24 #ifndef quantlib_zabr_smile_section_hpp
25 #define quantlib_zabr_smile_section_hpp
26 
27 #include <ql/pricingengines/blackformula.hpp>
28 #include <ql/termstructures/volatility/smilesection.hpp>
29 #include <ql/time/daycounters/actual365fixed.hpp>
30 #include <ql/experimental/volatility/zabr.hpp>
31 #include <ql/termstructures/volatility/smilesectionutils.hpp>
32 #include <vector>
33 
34 using std::exp;
35 
36 namespace QuantLib {
37 
38 // Evaluation Tags
39 
40 struct ZabrShortMaturityLognormal {};
41 struct ZabrShortMaturityNormal {};
42 struct ZabrLocalVolatility {};
43 struct ZabrFullFd {};
44 
45 template <typename Evaluation> class ZabrSmileSection : public SmileSection {
46   public:
47     ZabrSmileSection(Time timeToExpiry,
48                      Rate forward,
49                      const std::vector<Real>& zabrParameters,
50                      const std::vector<Real>& moneyness = std::vector<Real>(),
51                      Size fdRefinement = 5);
52     ZabrSmileSection(const Date& d,
53                      Rate forward,
54                      const std::vector<Real>& zabrParameters,
55                      const DayCounter& dc = Actual365Fixed(),
56                      const std::vector<Real>& moneyness = std::vector<Real>(),
57                      Size fdRefinement = 5);
58 
minStrike() const59     Real minStrike() const { return 0.0; }
maxStrike() const60     Real maxStrike() const { return QL_MAX_REAL; }
atmLevel() const61     Real atmLevel() const { return model_->forward(); }
optionPrice(Rate strike,Option::Type type=Option::Call,Real discount=1.0) const62     Real optionPrice(Rate strike, Option::Type type = Option::Call,
63                      Real discount = 1.0) const {
64         return optionPrice(strike, type, discount, Evaluation());
65     }
66 
model()67     ext::shared_ptr<ZabrModel> model() { return model_; }
68 
69   protected:
volatilityImpl(Rate strike) const70     Volatility volatilityImpl(Rate strike) const {
71         return volatilityImpl(strike, Evaluation());
72     }
73 
74   private:
init(const std::vector<Real> & moneyness)75     void init(const std::vector<Real> &moneyness) {
76         init(moneyness, Evaluation());
77         init2(Evaluation());
78         init3(Evaluation());
79     }
80     void init(const std::vector<Real> &moneyness, ZabrShortMaturityLognormal);
81     void init(const std::vector<Real> &moneyness, ZabrShortMaturityNormal);
82     void init(const std::vector<Real> &moneyness, ZabrLocalVolatility);
83     void init(const std::vector<Real> &moneyness, ZabrFullFd);
84     void init2(ZabrShortMaturityLognormal);
85     void init2(ZabrShortMaturityNormal);
86     void init2(ZabrLocalVolatility);
87     void init2(ZabrFullFd);
88     void init3(ZabrShortMaturityLognormal);
89     void init3(ZabrShortMaturityNormal);
90     void init3(ZabrLocalVolatility);
91     void init3(ZabrFullFd);
92     Real optionPrice(Rate strike, Option::Type type, Real discount,
93                      ZabrShortMaturityLognormal) const;
94     Real optionPrice(Rate strike, Option::Type type, Real discount,
95                      ZabrShortMaturityNormal) const;
96     Real optionPrice(Rate strike, Option::Type type, Real discount,
97                      ZabrLocalVolatility) const;
98     Real optionPrice(Rate strike, Option::Type type, Real discount,
99                      ZabrFullFd) const;
100     Volatility volatilityImpl(Rate strike, ZabrShortMaturityLognormal) const;
101     Volatility volatilityImpl(Rate strike, ZabrShortMaturityNormal) const;
102     Volatility volatilityImpl(Rate strike, ZabrLocalVolatility) const;
103     Volatility volatilityImpl(Rate strike, ZabrFullFd) const;
104     ext::shared_ptr<ZabrModel> model_;
105     Evaluation evaluation_;
106     Rate forward_;
107     std::vector<Real> params_;
108     const Size fdRefinement_;
109     std::vector<Real> strikes_, callPrices_;
110     ext::shared_ptr<Interpolation> callPriceFct_;
111     Real a_, b_;
112 };
113 
114 template <typename Evaluation>
ZabrSmileSection(Time timeToExpiry,Rate forward,const std::vector<Real> & zabrParams,const std::vector<Real> & moneyness,const Size fdRefinement)115 ZabrSmileSection<Evaluation>::ZabrSmileSection(
116     Time timeToExpiry, Rate forward, const std::vector<Real> &zabrParams,
117     const std::vector<Real> &moneyness, const Size fdRefinement)
118     : SmileSection(timeToExpiry, DayCounter()), forward_(forward),
119       params_(zabrParams), fdRefinement_(fdRefinement) {
120     init(moneyness);
121 }
122 
123 template <typename Evaluation>
ZabrSmileSection(const Date & d,Rate forward,const std::vector<Real> & zabrParams,const DayCounter & dc,const std::vector<Real> & moneyness,const Size fdRefinement)124 ZabrSmileSection<Evaluation>::ZabrSmileSection(
125     const Date &d, Rate forward, const std::vector<Real> &zabrParams,
126     const DayCounter &dc, const std::vector<Real> &moneyness,
127     const Size fdRefinement)
128     : SmileSection(d, dc, Date()), forward_(forward), params_(zabrParams),
129       fdRefinement_(fdRefinement) {
130     init(moneyness);
131 }
132 
133 template <typename Evaluation>
init(const std::vector<Real> &,ZabrShortMaturityLognormal)134 void ZabrSmileSection<Evaluation>::init(const std::vector<Real> &,
135                                         ZabrShortMaturityLognormal) {
136 
137     model_ = ext::make_shared<ZabrModel>(
138         exerciseTime(), forward_, params_[0], params_[1],
139                       params_[2], params_[3], params_[4]);
140 }
141 
142 template <typename Evaluation>
init(const std::vector<Real> & a,ZabrShortMaturityNormal)143 void ZabrSmileSection<Evaluation>::init(const std::vector<Real> &a,
144                                         ZabrShortMaturityNormal) {
145     init(a, ZabrShortMaturityLognormal());
146 }
147 
148 template <typename Evaluation>
init(const std::vector<Real> & moneyness,ZabrLocalVolatility)149 void ZabrSmileSection<Evaluation>::init(const std::vector<Real> &moneyness,
150                                         ZabrLocalVolatility) {
151 
152     QL_REQUIRE(params_.size() >= 5,
153                "zabr expects 5 parameters (alpha,beta,nu,rho,gamma) but ("
154                    << params_.size() << ") given");
155 
156     model_ = ext::make_shared<ZabrModel>(
157         exerciseTime(), forward_, params_[0], params_[1],
158                       params_[2], params_[3], params_[4]);
159 
160     // set up strike grid for local vol or full fd flavour of this section
161     // this is shared with SmileSectionUtils - unify later ?
162     static const Real defaultMoney[] = {
163         0.0, 0.01, 0.05, 0.10, 0.25, 0.40, 0.50, 0.60, 0.70, 0.80, 0.90,
164         1.0, 1.25, 1.5,  1.75, 2.0,  5.0,  7.5,  10.0, 15.0, 20.0};
165     std::vector<Real> tmp;
166     if (moneyness.empty())
167         tmp = std::vector<Real>(defaultMoney, defaultMoney + 21);
168     else
169         tmp = std::vector<Real>(moneyness);
170 
171     strikes_.clear(); // should not be necessary, anyway
172     Real lastF = 0.0;
173     bool firstStrike = true;
174     for (Size i = 0; i < tmp.size(); i++) {
175         Real f = tmp[i] * forward_;
176         if (f > 0.0) {
177             if (!firstStrike) {
178                 for (Size j = 1; j <= fdRefinement_; ++j) {
179                     strikes_.push_back(lastF +
180                                        ((double)j) * (f - lastF) /
181                                            (fdRefinement_ + 1));
182                 }
183             }
184             firstStrike = false;
185             lastF = f;
186             strikes_.push_back(f);
187         }
188     }
189 }
190 
191 template <typename Evaluation>
init(const std::vector<Real> & moneyness,ZabrFullFd)192 void ZabrSmileSection<Evaluation>::init(const std::vector<Real> &moneyness,
193                                         ZabrFullFd) {
194     init(moneyness, ZabrLocalVolatility());
195 }
196 
197 template <typename Evaluation>
init2(ZabrShortMaturityLognormal)198 void ZabrSmileSection<Evaluation>::init2(ZabrShortMaturityLognormal) {}
199 
200 template <typename Evaluation>
init2(ZabrShortMaturityNormal)201 void ZabrSmileSection<Evaluation>::init2(ZabrShortMaturityNormal) {}
202 
203 template <typename Evaluation>
init2(ZabrLocalVolatility)204 void ZabrSmileSection<Evaluation>::init2(ZabrLocalVolatility) {
205     callPrices_ = model_->fdPrice(strikes_);
206 }
207 
208 template <typename Evaluation>
init2(ZabrFullFd)209 void ZabrSmileSection<Evaluation>::init2(ZabrFullFd) {
210     callPrices_.resize(strikes_.size());
211 #pragma omp parallel for
212     for (long i = 0; i < (long)strikes_.size(); i++) {
213         callPrices_[i] = model_->fullFdPrice(strikes_[i]);
214     }
215 }
216 
217 template <typename Evaluation>
init3(ZabrShortMaturityLognormal)218 void ZabrSmileSection<Evaluation>::init3(ZabrShortMaturityLognormal) {}
219 
220 template <typename Evaluation>
init3(ZabrShortMaturityNormal)221 void ZabrSmileSection<Evaluation>::init3(ZabrShortMaturityNormal) {}
222 
223 template <typename Evaluation>
init3(ZabrLocalVolatility)224 void ZabrSmileSection<Evaluation>::init3(ZabrLocalVolatility) {
225     strikes_.insert(strikes_.begin(), 0.0);
226     callPrices_.insert(callPrices_.begin(), forward_);
227 
228     callPriceFct_ = ext::shared_ptr<Interpolation>(new CubicInterpolation(
229         strikes_.begin(), strikes_.end(), callPrices_.begin(),
230         CubicInterpolation::Spline, true, CubicInterpolation::SecondDerivative,
231         0.0, CubicInterpolation::SecondDerivative, 0.0));
232     // callPriceFct_ =
233     //     ext::shared_ptr<Interpolation>(new LinearInterpolation(
234     //         strikes_.begin(), strikes_.end(), callPrices_.begin()));
235 
236     callPriceFct_->enableExtrapolation();
237 
238     // on the right side we extrapolate exponetially (because spline
239     // does not make sense)
240     // we precompute the necessary parameters here
241     static const Real eps = 1E-5; // gap for first derivative computation
242 
243     Real c0 = (*callPriceFct_)(strikes_.back());
244     Real c0p = ((*callPriceFct_)(strikes_.back() - eps) - c0) / eps;
245 
246     a_ = c0p / c0;
247     b_ = std::log(c0) + a_ * strikes_.back();
248 }
249 
250 template <typename Evaluation>
init3(ZabrFullFd)251 void ZabrSmileSection<Evaluation>::init3(ZabrFullFd) {
252     init3(ZabrLocalVolatility());
253 }
254 
255 template <typename Evaluation>
256 Real
optionPrice(Real strike,Option::Type type,Real discount,ZabrShortMaturityLognormal) const257 ZabrSmileSection<Evaluation>::optionPrice(Real strike, Option::Type type,
258                                           Real discount,
259                                           ZabrShortMaturityLognormal) const {
260     return SmileSection::optionPrice(strike, type, discount);
261 }
262 
263 template <typename Evaluation>
optionPrice(Real strike,Option::Type type,Real discount,ZabrShortMaturityNormal) const264 Real ZabrSmileSection<Evaluation>::optionPrice(Real strike, Option::Type type,
265                                                Real discount,
266                                                ZabrShortMaturityNormal) const {
267     return bachelierBlackFormula(
268         type, strike, forward_,
269         model_->normalVolatility(strike) * std::sqrt(exerciseTime()), discount);
270 }
271 
272 template <typename Evaluation>
optionPrice(Rate strike,Option::Type type,Real discount,ZabrLocalVolatility) const273 Real ZabrSmileSection<Evaluation>::optionPrice(Rate strike, Option::Type type,
274                                                Real discount,
275                                                ZabrLocalVolatility) const {
276     Real call = strike <= strikes_.back() ? (*callPriceFct_)(strike)
277                                           : exp(-a_ * strike + b_);
278     if (type == Option::Call)
279         return call * discount;
280     else
281         return (call - (forward_ - strike)) * discount;
282 }
283 
284 template <typename Evaluation>
optionPrice(Rate strike,Option::Type type,Real discount,ZabrFullFd) const285 Real ZabrSmileSection<Evaluation>::optionPrice(Rate strike, Option::Type type,
286                                                Real discount,
287                                                ZabrFullFd) const {
288     return optionPrice(strike, type, discount, ZabrLocalVolatility());
289 }
290 
291 template <typename Evaluation>
292 Real
volatilityImpl(Rate strike,ZabrShortMaturityLognormal) const293 ZabrSmileSection<Evaluation>::volatilityImpl(Rate strike,
294                                              ZabrShortMaturityLognormal) const {
295     strike = std::max(1E-6, strike);
296     return model_->lognormalVolatility(strike);
297 }
298 
299 template <typename Evaluation>
300 Real
volatilityImpl(Rate strike,ZabrShortMaturityNormal) const301 ZabrSmileSection<Evaluation>::volatilityImpl(Rate strike,
302                                              ZabrShortMaturityNormal) const {
303     Real impliedVol = 0.0;
304     try {
305         Option::Type type;
306         if (strike >= model_->forward())
307             type = Option::Call;
308         else
309             type = Option::Put;
310         impliedVol =
311             blackFormulaImpliedStdDev(type, strike, model_->forward(),
312                                       optionPrice(strike, type, 1.0), 1.0) /
313             std::sqrt(exerciseTime());
314     } catch (...) {
315     }
316     return impliedVol;
317 }
318 
319 template <typename Evaluation>
volatilityImpl(Rate strike,ZabrLocalVolatility) const320 Real ZabrSmileSection<Evaluation>::volatilityImpl(Rate strike,
321                                                   ZabrLocalVolatility) const {
322     return volatilityImpl(strike, ZabrShortMaturityNormal());
323 }
324 
325 template <typename Evaluation>
volatilityImpl(Rate strike,ZabrFullFd) const326 Real ZabrSmileSection<Evaluation>::volatilityImpl(Rate strike,
327                                                   ZabrFullFd) const {
328     return volatilityImpl(strike, ZabrShortMaturityNormal());
329 }
330 }
331 
332 #endif
333