Logo Search packages:      
Sourcecode: quantlib version File versions  Download package


/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */

 Copyright (C) 2005, 2006 Klaus Spanderen

 This file is part of QuantLib, a free-software/open-source library
 for financial quantitative analysts and developers - http://quantlib.org/

 QuantLib is free software: you can redistribute it and/or modify it
 under the terms of the QuantLib license.  You should have received a
 copy of the license along with this program; if not, please email
 <quantlib-dev@lists.sf.net>. The license is also available online at

 This program is distributed in the hope that it will be useful, but WITHOUT
 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
 FOR A PARTICULAR PURPOSE.  See the license for more details.

#include <ql/legacy/libormarketmodels/lfmcovarproxy.hpp>
#include <ql/math/integrals/kronrodintegral.hpp>

namespace QuantLib {
        const boost::shared_ptr<LmVolatilityModel>& volaModel,
        const boost::shared_ptr<LmCorrelationModel>& corrModel)

    : LfmCovarianceParameterization(corrModel->size(), corrModel->factors()),
      corrModel_(corrModel) {

        QL_REQUIRE(volaModel_->size() == corrModel_->size(),
            "different size for the volatility (" << volaModel_->size() <<
            ") and correlation (" << corrModel_->size() <<
            ") models");

    LfmCovarianceProxy::volatilityModel() const {
        return volaModel_;

    LfmCovarianceProxy::correlationModel() const {
        return corrModel_;

    Disposable<Matrix> LfmCovarianceProxy::diffusion(Time t,
                                                     const Array& x) const {

        Matrix pca = corrModel_->pseudoSqrt(t, x);
        Array  vol = volaModel_->volatility(t, x);
        for (Size i=0; i<size_; ++i) {
            std::transform(pca.row_begin(i), pca.row_end(i),
                           std::bind2nd(std::multiplies<Real>(), vol[i]));

        return pca;

    Disposable<Matrix> LfmCovarianceProxy::covariance(Time t,
                                                      const Array& x) const {

        Array  volatility  = volaModel_->volatility(t, x);
        Matrix correlation = corrModel_->correlation(t, x);

        Matrix tmp(size_, size_);
        for (Size i=0; i<size_; ++i) {
            for (Size j=0; j<size_; ++j) {
                tmp[i][j] = volatility[i]*correlation[i][j]*volatility[j];

        return tmp;

    class LfmCovarianceProxy::Var_Helper {
        Var_Helper(const LfmCovarianceProxy* proxy, Size i, Size j);

        Real operator()(Real t) const;
        const Size i_, j_;
        const LmVolatilityModel*  const volaModel_;
        const LmCorrelationModel* const corrModel_;

    LfmCovarianceProxy::Var_Helper::Var_Helper(const LfmCovarianceProxy* proxy,
                                               Size i, Size j)
    : i_(i),
      corrModel_(proxy->corrModel_.get()) {

    Real LfmCovarianceProxy::Var_Helper::operator()(Real t) const {
        Volatility v1, v2;

        if (i_ == j_) {
            v1 = v2 = volaModel_->volatility(i_, t);
        } else {
            v1 = volaModel_->volatility(i_, t);
            v2 = volaModel_->volatility(j_, t);

        return  v1 * corrModel_->correlation(i_, j_, t) * v2;

    Real LfmCovarianceProxy::integratedCovariance(
                               Size i, Size j, Time t, const Array& x) const {

        if (corrModel_->isTimeIndependent()) {
            try {
                // if all objects support these methods
                // thats by far the fastest way to get the
                // integrated covariance
                return corrModel_->correlation(i, j, 0.0, x)
                        * volaModel_->integratedVariance(j, i, t, x);
            catch (Error&) {
                // okay proceed with the
                // slow numerical integration routine

        QL_REQUIRE(x.empty(), "can not handle given x here");

        Real tmp=0.0;
        Var_Helper helper(this, i, j);

        GaussKronrodAdaptive integrator(1e-10, 10000);
        for (Size k=0; k<64; ++k) {
            tmp+=integrator(helper, k*t/64., (k+1)*t/64.);
        return tmp;


Generated by  Doxygen 1.6.0   Back to index