QuantLib › quantlib-dev

Re: [QuantLib-svn] SF.net SVN: quantlib:[18270] trunk/QuantLib

Classic

List

Threaded

2 messages Options

Luigi Ballabio

Re: [QuantLib-svn] SF.net SVN: quantlib:[18270] trunk/QuantLib

Last one :)

Is the order of convergence (or a lack of one) known when an
extrapolation is created, given h and f? If so, we could pass it to
the constructor and leave the interface for calling the formula just
as operator()(t) so that we don't have to pass the order again and
again,

We might also overload it with operator()(t,s) if we didn't pass the
order of convergence. Just a thought, no big deal anyway---it's just
that formulaUnknownConvergence looks awfully long...

Luigi

On Sun, Jun 10, 2012 at 10:56 AM, <[hidden email]> wrote:

> Revision: 18270
> http://quantlib.svn.sourceforge.net/quantlib/?rev=18270&view=rev
> Author: klausspanderen
> Date: 2012-06-10 08:56:09 +0000 (Sun, 10 Jun 2012)
> Log Message:
> -----------
> added Richardson extrapolation
>
> Modified Paths:
> --------------
> trunk/QuantLib/ql/math/Makefile.am
> trunk/QuantLib/ql/math/all.hpp
> trunk/QuantLib/test-suite/interpolations.cpp
> trunk/QuantLib/test-suite/interpolations.hpp
>
> Added Paths:
> -----------
> trunk/QuantLib/ql/math/richardsonextrapolation.cpp
> trunk/QuantLib/ql/math/richardsonextrapolation.hpp
>
> Added: trunk/QuantLib/ql/math/richardsonextrapolation.cpp
> ===================================================================
> --- trunk/QuantLib/ql/math/richardsonextrapolation.cpp (rev 0)
> +++ trunk/QuantLib/ql/math/richardsonextrapolation.cpp 2012-06-10 08:56:09 UTC (rev 18270)
> @@ -0,0 +1,73 @@
> +/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
> +
> +/*
> + Copyright (C) 2012 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
> + <[hidden email]>. The license is also available online at
> + <http://quantlib.org/license.shtml>.
> +
> + 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.
> +*/
> +
> +/*! \file richardsonextrapolation.cpp
> +*/
> +
> +#include <ql/errors.hpp>
> +#include <ql/math/solvers1d/brent.hpp>
> +#include <ql/math/richardsonextrapolation.hpp>
> +
> +#include <cmath>
> +
> +namespace QuantLib {
> + namespace {
> + class RichardsonEqn {
> + public:
> + RichardsonEqn(Real fh, Real ft, Real fs, Real t, Real s)
> + : fh_(fh), ft_(ft), fs_(fs), t_(t), s_(s) { }
> +
> + Real operator()(Real k) const {
> + return ft_ + (ft_-fh_)/(std::pow(t_, k)-1.0)
> + - ( fs_ + (fs_-fh_)/(std::pow(s_, k)-1.0));
> + }
> + private:
> + const Real fh_, ft_, fs_, t_, s_;
> + };
> +
> + }
> +
> + RichardsonExtrapolation::RichardsonExtrapolation(
> + Real h, const boost::function<Real (Real)>& f)
> + : h_(h), fh_(f(h)), f_(f) { }
> +
> +
> + Real RichardsonExtrapolation::formula(Real t, Real k) const {
> +
> + QL_REQUIRE(t > 1, "scaling factor must be greater than 1");
> + const Real tk = std::pow(t,k);
> +
> + return (tk*f_(h_/t)-fh_)/(tk-1.0);
> + }
> +
> + Real RichardsonExtrapolation::formulaUnknownConvergence(Real t, Real s)
> + const {
> + QL_REQUIRE(t > 1 && s > 1, "scaling factors must be greater than 1");
> +
> + const Real ft = f_(h_/t);
> + const Real fs = f_(h_/s);
> +
> + const Real k = Brent().solve(RichardsonEqn(fh_, ft, fs, t, s),
> + 1e-8, 0.05, 10);
> +
> + const Real ts = std::pow(s,k);
> +
> + return (ts*fs-fh_)/(ts-1.0);
> + }
> +}
>
>
> Property changes on: trunk/QuantLib/ql/math/richardsonextrapolation.cpp
> ___________________________________________________________________
> Added: svn:mime-type
> + text/plain
>
> Added: trunk/QuantLib/ql/math/richardsonextrapolation.hpp
> ===================================================================
> --- trunk/QuantLib/ql/math/richardsonextrapolation.hpp (rev 0)
> +++ trunk/QuantLib/ql/math/richardsonextrapolation.hpp 2012-06-10 08:56:09 UTC (rev 18270)
> @@ -0,0 +1,62 @@
> +/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
> +
> +/*
> + Copyright (C) 2012 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
> + <[hidden email]>. The license is also available online at
> + <http://quantlib.org/license.shtml>.
> +
> + 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.
> +*/
> +
> +/*! \file richardsonextrapolation.hpp
> +*/
> +
> +#ifndef quantlib_richardson_extrapolation_hpp
> +#define quantlib_richardson_extrapolation_hpp
> +
> +#include <ql/types.hpp>
> +
> +#include <boost/function.hpp>
> +
> +namespace QuantLib {
> +
> + //! Richardson Extrapolatio
> + /*! References:
> +
> + http://en.wikipedia.org/wiki/Richardson_extrapolation
> + */
> +
> + class RichardsonExtrapolation {
> + public:
> + RichardsonExtrapolation(Real h,
> + const boost::function<Real (Real)>& f);
> +
> + /*! Extrapolation for known order of convergence
> + \param step size for extrapolation. Is often 2
> + \param k order of convergence
> + */
> + Real formula(Real t, Real k) const;
> +
> + /*! Extrapolation for unknown order of convergence
> + \param t first step size for the extrapolation.
> + \param s second step size for the extrapolation
> + */
> + Real formulaUnknownConvergence(Real t, Real s) const;
> +
> + private:
> + const Real h_;
> + const Real fh_;
> + const boost::function<Real (Real)> f_;
> + };
> +}
> +
> +#endif
>
>
> Property changes on: trunk/QuantLib/ql/math/richardsonextrapolation.hpp
> ___________________________________________________________________
> Added: svn:mime-type
> + text/plain
>
> Modified: trunk/QuantLib/test-suite/interpolations.cpp
> ===================================================================
> --- trunk/QuantLib/test-suite/interpolations.cpp 2012-05-30 22:32:14 UTC (rev 18269)
> +++ trunk/QuantLib/test-suite/interpolations.cpp 2012-06-10 08:56:09 UTC (rev 18270)
> @@ -36,6 +36,7 @@
> #include <ql/math/integrals/simpsonintegral.hpp>
> #include <ql/math/kernelfunctions.hpp>
> #include <ql/math/functional.hpp>
> +#include <ql/math/richardsonextrapolation.hpp>
> #include <ql/math/randomnumbers/sobolrsg.hpp>
> #include <ql/math/optimization/levenbergmarquardt.hpp>
> #include <boost/assign/std/vector.hpp>
> @@ -1610,6 +1611,35 @@
> }
> }
>
> +namespace {
> + Real f(Real h) {
> + return std::pow( 1.0 + h, 1/h);
> + }
> +}
> +
> +void InterpolationTest::testRichardsonExtrapolation() {
> + BOOST_MESSAGE("Testing Richardson Extrapolation...");
> +
> + /* example taken from
> + * http://www.ipvs.uni-stuttgart.de/abteilungen/bv/lehre/
> + * lehrveranstaltungen/vorlesungen/WS0910/
> + * NSG_termine/dateien/Richardson.pdf
> + */
> +
> + const RichardsonExtrapolation extrap(0.1, f);
> +
> + const Real stepSize = 2.0;
> + const Real orderOfConvergence = 1.0;
> +
> + const Real tol = 0.00001;
> + const Real expected = 2.71285;
> + const Real calculated = extrap.formula(stepSize, orderOfConvergence);
> +
> + if (std::fabs(expected-calculated) > tol) {
> + BOOST_ERROR("failed to reproduce Richardson extrapolation");
> + }
> +}
> +
> test_suite* InterpolationTest::suite() {
> test_suite* suite = BOOST_TEST_SUITE("Interpolation tests");
>
> @@ -1636,6 +1666,8 @@
> suite->add(QUANTLIB_TEST_CASE(
> &InterpolationTest::testKernelInterpolation2D));
> suite->add(QUANTLIB_TEST_CASE(&InterpolationTest::testBicubicDerivatives));
> + suite->add(QUANTLIB_TEST_CASE(
> + &InterpolationTest::testRichardsonExtrapolation));
>
> return suite;
> }
>
> Modified: trunk/QuantLib/test-suite/interpolations.hpp
> ===================================================================
> --- trunk/QuantLib/test-suite/interpolations.hpp 2012-05-30 22:32:14 UTC (rev 18269)
> +++ trunk/QuantLib/test-suite/interpolations.hpp 2012-06-10 08:56:09 UTC (rev 18270)
> @@ -46,6 +46,8 @@
> static void testKernelInterpolation();
> static void testKernelInterpolation2D();
> static void testBicubicDerivatives();
> + static void testRichardsonExtrapolation();
> +
> static boost::unit_test_framework::test_suite* suite();
> };
>
>
> This was sent by the SourceForge.net collaborative development platform, the world's largest Open Source development site.
>
>
> ------------------------------------------------------------------------------
> Live Security Virtual Conference
> Exclusive live event will cover all the ways today's security and
> threat landscape has changed and how IT managers can respond. Discussions
> will include endpoint security, mobile security and the latest in malware
> threats. http://www.accelacomm.com/jaw/sfrnl04242012/114/50122263/
> _______________________________________________
> QuantLib-cvs mailing list
> [hidden email]
> https://lists.sourceforge.net/lists/listinfo/quantlib-cvs

------------------------------------------------------------------------------
Live Security Virtual Conference
Exclusive live event will cover all the ways today's security and
threat landscape has changed and how IT managers can respond. Discussions
will include endpoint security, mobile security and the latest in malware
threats. http://www.accelacomm.com/jaw/sfrnl04242012/114/50122263/
_______________________________________________
QuantLib-dev mailing list
[hidden email]
https://lists.sourceforge.net/lists/listinfo/quantlib-dev

Klaus Spanderen-2

Re: [QuantLib-svn] SF.net SVN: quantlib:[18270] trunk/QuantLib

Hi Luigi,

good idea, thanks for the hint. I'm going to change it.

regards

Klaus

On Wednesday, June 27, 2012 05:44:16 PM Luigi Ballabio wrote:

> Last one :)

> Is the order of convergence (or a lack of one) known when an

> extrapolation is created, given h and f? If so, we could pass it to

> the constructor and leave the interface for calling the formula just

> as operator()(t) so that we don't have to pass the order again and

> again,

> We might also overload it with operator()(t,s) if we didn't pass the

> order of convergence. Just a thought, no big deal anyway---it's just

> that formulaUnknownConvergence looks awfully long...

> Luigi