1 /* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2 
3 /*
4  Copyright (C) 2009 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 fdmblackscholesmesher.cpp
21     \brief 1-d mesher for the Black-Scholes process (in ln(S))
22 */
23 
24 #include <ql/processes/blackscholesprocess.hpp>
25 #include <ql/termstructures/yieldtermstructure.hpp>
26 #include <ql/termstructures/yield/quantotermstructure.hpp>
27 #include <ql/termstructures/volatility/equityfx/blackconstantvol.hpp>
28 #include <ql/math/distributions/normaldistribution.hpp>
29 #include <ql/methods/finitedifferences/utilities/fdmquantohelper.hpp>
30 #include <ql/methods/finitedifferences/meshers/uniform1dmesher.hpp>
31 #include <ql/methods/finitedifferences/meshers/concentrating1dmesher.hpp>
32 #include <ql/methods/finitedifferences/meshers/fdmblackscholesmesher.hpp>
33 
34 namespace QuantLib {
35 
FdmBlackScholesMesher(Size size,const ext::shared_ptr<GeneralizedBlackScholesProcess> & process,Time maturity,Real strike,Real xMinConstraint,Real xMaxConstraint,Real eps,Real scaleFactor,const std::pair<Real,Real> & cPoint,const DividendSchedule & dividendSchedule,const ext::shared_ptr<FdmQuantoHelper> & fdmQuantoHelper,Real spotAdjustment)36     FdmBlackScholesMesher::FdmBlackScholesMesher(
37         Size size,
38         const ext::shared_ptr<GeneralizedBlackScholesProcess>& process,
39         Time maturity, Real strike,
40         Real xMinConstraint, Real xMaxConstraint,
41         Real eps, Real scaleFactor,
42         const std::pair<Real, Real>& cPoint,
43         const DividendSchedule& dividendSchedule,
44         const ext::shared_ptr<FdmQuantoHelper>& fdmQuantoHelper,
45         Real spotAdjustment)
46     : Fdm1dMesher(size) {
47 
48         const Real S = process->x0();
49         QL_REQUIRE(S > 0.0, "negative or null underlying given");
50 
51         std::vector<std::pair<Time, Real> > intermediateSteps;
52         for (Size i=0; i < dividendSchedule.size(); ++i) {
53             const Time t = process->time(dividendSchedule[i]->date());
54             if (t <= maturity && t >= 0.0)
55                 intermediateSteps.push_back(
56                     std::make_pair(
57                         process->time(dividendSchedule[i]->date()),
58                         dividendSchedule[i]->amount()
59                     ) );
60         }
61 
62         const Size intermediateTimeSteps = std::max<Size>(2, Size(24.0*maturity));
63         for (Size i=0; i < intermediateTimeSteps; ++i)
64             intermediateSteps.push_back(
65                 std::make_pair((i+1)*(maturity/intermediateTimeSteps), 0.0));
66 
67         std::sort(intermediateSteps.begin(), intermediateSteps.end());
68 
69         const Handle<YieldTermStructure> rTS = process->riskFreeRate();
70 
71         const Handle<YieldTermStructure> qTS =
72             (fdmQuantoHelper) != 0 ?
73                 Handle<YieldTermStructure>(ext::make_shared<QuantoTermStructure>(
74                     process->dividendYield(), process->riskFreeRate(),
75                     Handle<YieldTermStructure>(fdmQuantoHelper->fTS_), process->blackVolatility(),
76                     strike, Handle<BlackVolTermStructure>(fdmQuantoHelper->fxVolTS_),
77                     fdmQuantoHelper->exchRateATMlevel_, fdmQuantoHelper->equityFxCorrelation_)) :
78                 process->dividendYield();
79 
80         Time lastDivTime = 0.0;
81         Real fwd = S + spotAdjustment;
82         Real mi = fwd, ma = fwd;
83 
84         for (Size i=0; i < intermediateSteps.size(); ++i) {
85             const Time divTime = intermediateSteps[i].first;
86             const Real divAmount = intermediateSteps[i].second;
87 
88             fwd = fwd / rTS->discount(divTime) * rTS->discount(lastDivTime)
89                       * qTS->discount(divTime) / qTS->discount(lastDivTime);
90 
91             mi  = std::min(mi, fwd); ma = std::max(ma, fwd);
92 
93             fwd-= divAmount;
94 
95             mi  = std::min(mi, fwd); ma = std::max(ma, fwd);
96 
97             lastDivTime = divTime;
98         }
99 
100         // Set the grid boundaries
101         const Real normInvEps = InverseCumulativeNormal()(1-eps);
102         const Real sigmaSqrtT
103             = process->blackVolatility()->blackVol(maturity, strike)
104                                                         *std::sqrt(maturity);
105 
106         Real xMin = std::log(mi) - sigmaSqrtT*normInvEps*scaleFactor;
107         Real xMax = std::log(ma) + sigmaSqrtT*normInvEps*scaleFactor;
108 
109         if (xMinConstraint != Null<Real>()) {
110             xMin = xMinConstraint;
111         }
112         if (xMaxConstraint != Null<Real>()) {
113             xMax = xMaxConstraint;
114         }
115 
116         ext::shared_ptr<Fdm1dMesher> helper;
117         if (   cPoint.first != Null<Real>()
118             && std::log(cPoint.first) >=xMin && std::log(cPoint.first) <=xMax) {
119 
120             helper = ext::shared_ptr<Fdm1dMesher>(
121                 new Concentrating1dMesher(xMin, xMax, size,
122                     std::pair<Real,Real>(std::log(cPoint.first),
123                                          cPoint.second)));
124         }
125         else {
126             helper = ext::shared_ptr<Fdm1dMesher>(
127                                         new Uniform1dMesher(xMin, xMax, size));
128 
129         }
130 
131         locations_ = helper->locations();
132         for (Size i=0; i < locations_.size(); ++i) {
133             dplus_[i]  = helper->dplus(i);
134             dminus_[i] = helper->dminus(i);
135         }
136     }
137 
138     ext::shared_ptr<GeneralizedBlackScholesProcess>
processHelper(const Handle<Quote> & s0,const Handle<YieldTermStructure> & rTS,const Handle<YieldTermStructure> & qTS,Volatility vol)139     FdmBlackScholesMesher::processHelper(const Handle<Quote>& s0,
140                                          const Handle<YieldTermStructure>& rTS,
141                                          const Handle<YieldTermStructure>& qTS,
142                                          Volatility vol) {
143 
144         return ext::make_shared<GeneralizedBlackScholesProcess>(
145 
146                 s0, qTS, rTS,
147                 Handle<BlackVolTermStructure>(
148                     ext::shared_ptr<BlackVolTermStructure>(
149                         new BlackConstantVol(rTS->referenceDate(),
150                                              Calendar(),
151                                              vol,
152                                              rTS->dayCounter()))));
153     }
154 }
155 
156