CommodityAveragePriceOptionMonteCarloEngine Class Reference

#include <qle/pricingengines/commodityapoengine.hpp>

Inheritance diagram for CommodityAveragePriceOptionMonteCarloEngine:

Collaboration diagram for CommodityAveragePriceOptionMonteCarloEngine:

Public Member Functions
	CommodityAveragePriceOptionMonteCarloEngine (const QuantLib::Handle< QuantLib::YieldTermStructure > &discountCurve, const QuantLib::Handle< QuantExt::BlackScholesModelWrapper > &model, QuantLib::Size samples, QuantLib::Real beta=0.0, const QuantLib::Size seed=42)
	CommodityAveragePriceOptionMonteCarloEngine (const QuantLib::Handle< QuantLib::YieldTermStructure > &discountCurve, const QuantLib::Handle< QuantLib::BlackVolTermStructure > &vol, QuantLib::Size samples, QuantLib::Real beta=0.0, const QuantLib::Size seed=42)
void	calculate () const override
Public Member Functions inherited from CommodityAveragePriceOptionBaseEngine
	CommodityAveragePriceOptionBaseEngine (const QuantLib::Handle< QuantLib::YieldTermStructure > &discountCurve, const QuantLib::Handle< QuantExt::BlackScholesModelWrapper > &model, QuantLib::Real beta=0.0)
	CommodityAveragePriceOptionBaseEngine (const QuantLib::Handle< QuantLib::YieldTermStructure > &discountCurve, const QuantLib::Handle< QuantLib::BlackVolTermStructure > &vol, QuantLib::Real beta=0.0)

Private Member Functions
void	calculateSpot () const
	Calculations when underlying swap references a commodity spot price. More...
void	calculateFuture () const
	Calculations when underlying swap references a commodity spot price. More...
void	setupFuture (std::vector< QuantLib::Real > &outVolatilities, QuantLib::Matrix &outSqrtCorr, std::vector< QuantLib::Real > &outPrices, std::vector< QuantLib::Size > &futureIndex, QuantLib::Real strike) const
std::vector< QuantLib::Real >	timegrid (std::vector< QuantLib::Date > &outDates) const

Private Attributes
QuantLib::Size	samples_
QuantLib::Size	seed_

Additional Inherited Members
Protected Member Functions inherited from CommodityAveragePriceOptionBaseEngine
QuantLib::Real	rho (const QuantLib::Date &ed_1, const QuantLib::Date &ed_2) const
	Return the correlation between two future expiry dates ed_1 and ed_2. More...
bool	isModelDependent () const
bool	barrierTriggered (const Real price, const bool logPrice) const
bool	alive (const bool barrierTriggered) const
Protected Attributes inherited from CommodityAveragePriceOptionBaseEngine
QuantLib::Handle< QuantLib::YieldTermStructure >	discountCurve_
QuantLib::Handle< QuantLib::BlackVolTermStructure >	volStructure_
QuantLib::Real	beta_
QuantLib::Real	logBarrier_

Detailed Description

Commodity APO Monte Carlo Engine Monte Carlo implementation of the APO payoff Reference: Iain Clark, Commodity Option Pricing, Wiley, section 2.7.4, equations (2.118) and (2.126)

Definition at line 124 of file commodityapoengine.hpp.

Constructor & Destructor Documentation

CommodityAveragePriceOptionMonteCarloEngine() [1/2]

CommodityAveragePriceOptionMonteCarloEngine	(	const QuantLib::Handle< QuantLib::YieldTermStructure > &	discountCurve,
		const QuantLib::Handle< QuantExt::BlackScholesModelWrapper > &	model,
		QuantLib::Size	samples,
		QuantLib::Real	beta = 0.0,
		const QuantLib::Size	seed = 42
	)

Definition at line 126 of file commodityapoengine.hpp.

        : CommodityAveragePriceOptionBaseEngine(discountCurve, model, beta), samples_(samples), seed_(seed) {}

CommodityAveragePriceOptionMonteCarloEngine() [2/2]

CommodityAveragePriceOptionMonteCarloEngine	(	const QuantLib::Handle< QuantLib::YieldTermStructure > &	discountCurve,
		const QuantLib::Handle< QuantLib::BlackVolTermStructure > &	vol,
		QuantLib::Size	samples,
		QuantLib::Real	beta = 0.0,
		const QuantLib::Size	seed = 42
	)

Definition at line 133 of file commodityapoengine.hpp.

        : CommodityAveragePriceOptionBaseEngine(discountCurve, vol, beta), samples_(samples), seed_(seed) {}

Member Function Documentation

calculate()

void calculate ( ) const

override

Definition at line 318 of file commodityapoengine.cpp.

                                                                  {
    if (isModelDependent()) {
        // Switch implementation depending on whether underlying swap references a spot or future price
        if (arguments_.flow->useFuturePrice()) {
            calculateFuture();
        } else {
            calculateSpot();
        }
    }
}

Here is the call graph for this function:

calculateSpot()

void calculateSpot ( ) const

private

Calculations when underlying swap references a commodity spot price.

Definition at line 329 of file commodityapoengine.cpp.

                                                                      {
 
    // Discount factor to the APO payment date
    Real discount = discountCurve_->discount(arguments_.flow->date());
 
    // Put call indicator
    Real omega = arguments_.type == Option::Call ? 1.0 : -1.0;
 
    // Variable to hold the payoff
    Real payoff = 0.0;
 
    // Vector of timesteps from today = t_0 out to last pricing date t_n
    // i.e. {t_1 - t_0, t_2 - t_1,..., t_n - t_{n-1}}
    vector<Date> dates;
    vector<Real> dt = timegrid(dates);
 
    // On each Monte Carlo sample, we must generate the spot price process path over n time steps.
    LowDiscrepancy::rsg_type rsg = LowDiscrepancy::make_sequence_generator(dt.size(), seed_);
 
    // We will read the volatility off the surface at the effective strike
    // We will only call this method when the effectiveStrike > 0 but will check anyway
    Real effectiveStrike = arguments_.effectiveStrike - arguments_.accrued;
    QL_REQUIRE(effectiveStrike > 0.0, "calculateSpot: expected effectiveStrike to be positive");
 
    // Precalculate:
    // exp(-1/2 forward variance over dt): \exp \left(-\frac{1}{2} \int_{t_{i-1}}^{t_i} \sigma^2(u) du \right)
    // forward std dev over dt: \exp \left(\sqrt{\int_{t_{i-1}}^{t_i} \sigma^2(u) du} \right)
    // forward ratio: \exp \left(\int_{t_{i-1}}^{t_i} r(u) - r_f(u) \, du \right) = \frac{F(0, t_i)}{F(0, t_{i-1})}
    //                first period is multiplied by S(0)
    // factors array is exp(-1/2 forward variance over dt) x forward ratio
    Array expHalfFwdVar(dt.size(), 0.0);
    Array fwdStdDev(dt.size(), 0.0);
    Array fwdRatio(dt.size(), 0.0);
    Time t = 0.0;
    for (Size i = 0; i < dt.size(); i++) {
        t += dt[i];
        expHalfFwdVar[i] = volStructure_->blackForwardVariance(t - dt[i], t, effectiveStrike);
        fwdStdDev[i] = sqrt(expHalfFwdVar[i]);
        expHalfFwdVar[i] = exp(-expHalfFwdVar[i] / 2.0);
        Real fxRate{1.};
        if(arguments_.flow->fxIndex())
            fxRate = arguments_.flow->fxIndex()->fixing(dates[i+1]);
        fwdRatio[i] = fxRate * arguments_.flow->index()->fixing(dates[i + 1]);
        if (i > 0) {
            if(arguments_.flow->fxIndex())
                fxRate = arguments_.flow->fxIndex()->fixing(dates[i]);
            fwdRatio[i] /= (fxRate * arguments_.flow->index()->fixing(dates[i]));
        }
    }
    Array factors = expHalfFwdVar * fwdRatio;
 
    // Loop over each sample
    auto m = arguments_.flow->indices().size();
    for (Size k = 0; k < samples_; k++) {
 
        // Sequence is n independent standard normal random variables
        vector<Real> sequence = rsg.nextSequence().value;
 
        // Calculate payoff on this sample - price holds the evolved price
        Real price = 0.0;
        Real samplePayoff = 0.0;
        Array path(sequence.begin(), sequence.end());
        path = Exp(path * fwdStdDev) * factors;
        bool barrierTriggered = false;
        for (Size i = 0; i < dt.size(); i++) {
 
            // Update price
            price = i == 0 ? path[i] : price * path[i];
 
            // Update sum of the spot prices on the pricing dates after today
            samplePayoff += price;
 
            // check barrier
            if (arguments_.barrierStyle == Exercise::American)
                barrierTriggered = barrierTriggered || this->barrierTriggered(price, false);
        }
        // Average price on this sample
        samplePayoff /= m;
 
        // Finally, the payoff on this sample
        samplePayoff = max(omega * (samplePayoff - effectiveStrike), 0.0);
 
        // account for barrier
        if (arguments_.barrierStyle == Exercise::European)
            barrierTriggered = this->barrierTriggered(price, false);
 
        if (!alive(barrierTriggered))
            samplePayoff = 0.0;
 
        // Update the final average
        // Note: payoff * k / (k + 1) - left to right precedence => double division.
        payoff = k == 0 ? samplePayoff : payoff * k / (k + 1) + samplePayoff / (k + 1);
    }
 
    // Populate the result value
    results_.value = arguments_.quantity * arguments_.flow->gearing() * payoff * discount;
}

Here is the call graph for this function:

Here is the caller graph for this function:

calculateFuture()

void calculateFuture ( ) const

private

Calculations when underlying swap references a commodity spot price.

Definition at line 427 of file commodityapoengine.cpp.

                                                                        {
 
    // this method uses checkBarrier() on log prices, therefore we have to init the log barrier level
    if (arguments_.barrierLevel != Null<Real>())
        logBarrier_ = std::log(arguments_.barrierLevel);
 
    // Discount factor to the APO payment date
    Real discount = discountCurve_->discount(arguments_.flow->date());
 
    // Put call indicator
    Real omega = arguments_.type == Option::Call ? 1.0 : -1.0;
 
    // Variable to hold the payoff
    Real payoff = 0.0;
 
    // We will read the volatility off the surface the effective strike
    // We will only call this method when the effectiveStrike > 0 but will check anyway
    Real effectiveStrike = arguments_.effectiveStrike - arguments_.accrued;
    QL_REQUIRE(effectiveStrike > 0.0, "calculateFuture: expected effectiveStrike to be positive");
 
    // Unique future expiry dates i.e. contracts, their volatilities and the correlation matrix between them
    Matrix sqrtCorr;
    vector<Real> vols;
    vector<Real> prices;
    vector<Size> futureIndex;
    setupFuture(vols, sqrtCorr, prices, futureIndex, effectiveStrike);
 
    // Vector of timesteps from today = t_0 out to last pricing date t_n
    // i.e. {t_1 - t_0, t_2 - t_1,..., t_n - t_{n-1}}. Don't need the dates here.
    vector<Date> dates;
    vector<Real> dt = timegrid(dates);
 
    // On each Monte Carlo sample, we must generate the paths for N (size of vols) future contracts where
    // each path has n time steps. We will represent the paths with an N x n matrix. First step is to fill the
    // matrix with N x n _independent_ standard normal variables. Then correlate the N variables in each column
    // using the sqrtCorr matrix and then fill each entry F_{i, j} in the matrix with the value of the i-th
    // future price process at timestep j. Note, we will possibly simulate contracts past their expiries but not
    // use the price in the APO rate averaging.
    LowDiscrepancy::rsg_type rsg = LowDiscrepancy::make_sequence_generator(vols.size() * dt.size(), seed_);
 
    // Precalculate exp(-0.5 \sigma_i^2 \delta t_j) and std dev = sqrt(\delta t_j) \sigma_i
    Matrix drifts(vols.size(), dt.size(), 0.0);
    Matrix stdDev(vols.size(), dt.size(), 0.0);
    Array logPrices(vols.size());
    for (Size i = 0; i < drifts.rows(); i++) {
        logPrices[i] = std::log(prices[i]);
        for (Size j = 0; j < drifts.columns(); j++) {
            drifts[i][j] = -vols[i] * vols[i] * dt[j] / 2.0;
            stdDev[i][j] = vols[i] * sqrt(dt[j]);
        }
    }
 
    // Loop over each sample
    auto m = arguments_.flow->indices().size();
    Matrix paths(vols.size(), dt.size());
    for (Size k = 0; k < samples_; ++k) {
 
        // Sequence is N x n independent standard normal random variables
        const vector<Real>& sequence = rsg.nextSequence().value;
 
        // paths initially holds independent standard normal random variables
        std::copy(sequence.begin(), sequence.end(), paths.begin());
 
        // Correlate the random variables in each column
        paths = sqrtCorr * paths;
 
        // Fill the paths
        for (Size i = 0; i < paths.rows(); ++i) {
            for (Size j = 0; j < paths.columns(); ++j) {
                paths[i][j] = (j == 0 ? logPrices[i] : paths[i][j - 1]) + drifts[i][j] + stdDev[i][j] * paths[i][j];
            }
        }
 
        // Calculate the sum of the commodity future prices on the pricing dates after today
        Real samplePayoff = 0.0;
        bool barrierTriggered = false;
        Real price = 0.0;
        for (Size j = 0; j < dt.size(); ++j) {
            price = paths[futureIndex[j]][j];
            if (arguments_.barrierStyle == Exercise::American)
                barrierTriggered = barrierTriggered || this->barrierTriggered(price, true);
            samplePayoff += std::exp(price);
        }
 
        // Average price on this sample
        samplePayoff /= m;
 
        // Finally, the payoff on this sample
        samplePayoff = max(omega * (samplePayoff - effectiveStrike), 0.0);
 
        // account for barrier
        if (arguments_.barrierStyle == Exercise::European)
            barrierTriggered = this->barrierTriggered(price, true);
 
        if (!alive(barrierTriggered))
            samplePayoff = 0.0;
 
        // Update the final average
        // Note: payoff * k / (k + 1) - left to right precedence => double division.
        payoff = k == 0 ? samplePayoff : payoff * k / (k + 1) + samplePayoff / (k + 1);
    }
 
    // Populate the result value
    results_.value = arguments_.quantity * arguments_.flow->gearing() * payoff * discount;
}

Here is the call graph for this function:

Here is the caller graph for this function:

setupFuture()

void setupFuture ( std::vector< QuantLib::Real > &  outVolatilities, QuantLib::Matrix &  outSqrtCorr, std::vector< QuantLib::Real > &  outPrices, std::vector< QuantLib::Size > &  futureIndex, QuantLib::Real  strike  ) const

private

Prepare data for APO calculation. The outVolatilities parameter will be populated with separate future contract volatilities taking into account the strike level. The number of elements of outVolatilities gives the number, N, of future contracts involved in the non-accrued portion of the APO. The matrix outSqrtCorr is populated with the square root of the correlation matrix between the future contracts. The outPrices vector will be populated with the current future price values. The futureIndex is populated with the index of the future to be used on each timestep in the simulation.

Definition at line 533 of file commodityapoengine.cpp.

                                                                                 {
 
    // Ensure that vectors are empty
    outVolatilities.clear();
    prices.clear();
    futureIndex.clear();
 
    // Will hold the result
    map<Date, Real> result;
 
    // Note that here we make the simplifying assumption that the volatility can be read from the volatility
    // term structure at the future contract's expiry date. In most cases, if the volatility term structure is built
    // from options on futures, the option contract expiry will be a number of days before the future contract expiry
    // and we should really read off the term structure at that date. Also populate a temp vector containing the key
    // dates for use in the loop below where we populate the sqrt correlation matrix.
 
    // Initialise result with expiry date keys that are still live in the APO
    Date today = Settings::instance().evaluationDate();
    set<Date> expiryDates;
    for (const auto& p : arguments_.flow->indices()) {
        if (p.first > today) {
            Date expiry = p.second->expiryDate();
            // If expiry has not been encountered yet
            if (expiryDates.insert(expiry).second) {
                outVolatilities.push_back(volStructure_->blackVol(expiry, strike));
                Real fxRate{1.};
                if(arguments_.flow->fxIndex())
                    fxRate = arguments_.flow->fxIndex()->fixing(expiry);
                prices.push_back(fxRate*p.second->fixing(today));//check if today should not be p.first
            }
            futureIndex.push_back(expiryDates.size() - 1);
        }
    }
 
    // Populate the outSqrtCorr matrix
    vector<Date> vExpiryDates(expiryDates.begin(), expiryDates.end());
    outSqrtCorr = Matrix(vExpiryDates.size(), vExpiryDates.size(), 1.0);
    for (Size i = 0; i < vExpiryDates.size(); i++) {
        for (Size j = 0; j < i; j++) {
            outSqrtCorr[i][j] = outSqrtCorr[j][i] = rho(vExpiryDates[i], vExpiryDates[j]);
        }
    }
    outSqrtCorr = pseudoSqrt(outSqrtCorr, SalvagingAlgorithm::None);
}

Here is the call graph for this function:

Here is the caller graph for this function:

timegrid()

vector< Real > timegrid ( std::vector< QuantLib::Date > &  outDates ) const

private

Return the \(n\) timesteps from today, \(t_0\), up to \(t_n\) where \(n > 0\). Note that each \(t_i\) corresponds to a pricing date \(d_i\) that is after today. The method returns the vector of time deltas \(t_i - t_{i-1}\) for \(i=1,\ldots,n\) and populates the vector outDates with the dates \(d_0, d_1,...,d_n\). Note that the size of outDates is one larger than the size of the return vector.

Definition at line 580 of file commodityapoengine.cpp.

                                                                                               {
 
    // Initialise result with expiry date keys that are still live in the APO
    // i.e. {t_1 - t_0, t_2 - t_0,..., t_n - t_0} where t_0 corresponds to today
    vector<Real> times;
    outDates.clear();
    Date today = Settings::instance().evaluationDate();
    outDates.push_back(today);
    for (const auto& p : arguments_.flow->indices()) {
        if (p.first > today) {
            outDates.push_back(p.first);
            times.push_back(volStructure_->timeFromReference(p.first));
        }
    }
 
    // Populate time deltas i.e. {t_1 - t_0, t_2 - t_1,..., t_n - t_{n-1}}
    vector<Real> dt(times.size());
    adjacent_difference(times.begin(), times.end(), dt.begin());
 
    return dt;
}

Here is the caller graph for this function:

Member Data Documentation

samples_

QuantLib::Size samples_

private

Definition at line 166 of file commodityapoengine.hpp.

seed_

QuantLib::Size seed_

private

Definition at line 167 of file commodityapoengine.hpp.