maphenon.h

Go to the documentation of this file.

/////////////////////////////////////////////////////////////////////////////
///
/// @file maphenon.h
///
/// The Henon map.
/// This file contains the definition of a class that represents
/// the Henon map. The class includes a method for computing
/// a rigorous outer approximation of the image of an open rectangle.
///
/// @author Pawel Pilarczyk
///
/////////////////////////////////////////////////////////////////////////////

// Copyright (C) 1997-2010 by Pawel Pilarczyk.
//
// This file is part of my research software package.  This is free software;
// you can redistribute it and/or modify it under the terms of the GNU
// General Public License as published by the Free Software Foundation;
// either version 2 of the License, or (at your option) any later version.
//
// This software 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
// GNU General Public License for more details.
//
// You should have received a copy of the GNU General Public License along
// with this software; see the file "license.txt".  If not, write to the
// Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
// MA 02111-1307, USA.

// Started on December 22, 2008. Last revision: December 22, 2008.


#ifndef _FINRESDYN_MAPHENON_H_
#define _FINRESDYN_MAPHENON_H_


// include some standard C++ header files
#include <iostream>

// include local header files
#include "rounding.h"
#include "opensqr.h"


// --------------------------------------------------
// ----------------- the Henon map ------------------
// --------------------------------------------------

/// The Henon map.
/// This class defines the Henon map: (x,y) -> (1 + y - ax^2, bx),
/// where a = 1.4 and b = 0.3.
template <class NumType>
class tMapHenon
{
public:
        /// The type of numbers in the Henon map.
        typedef NumType NumberType;

        /// The constructor of a Henon map object.
        tMapHenon (const NumType &aNumerator = 14,
                const NumType &aDenominator = 10,
                const NumType &bNumerator = 3,
                const NumType &bDenominator = 10);

        /// Returns the dimension of the phase space.
        static int dim ();

        /// Returns an initial point for constructing a cover
        /// of the attractor. The vector provided is iterated
        /// for the given time by the map, and the result is stored
        /// back in the provided vector.
        void initialPoint (NumType *x, int iterCount) const;

        /// Returns a rigorous bound of the image of the given rectangle
        /// by the Henon map.
        template <class InputType, class ResultType>
        void image (const InputType *xMin, const InputType *xMax,
                ResultType *yMin, ResultType *yMax) const;

private:
        /// A lower bound for the parameter 'a' in the Henon map.
        const NumType aMin;

        /// An upper bound for the parameter 'a' in the Henon map.
        const NumType aMax;

        /// A lower bound for the parameter 'b' in the Henon map.
        const NumType bMin;

        /// An upper bound for the parameter 'b' in the Henon map.
        const NumType bMax;

}; /* class tMapHenon */

// --------------------------------------------------

template <class NumType>
inline tMapHenon<NumType>::tMapHenon
        (const NumType &aNumerator, const NumType &aDenominator,
        const NumType &bNumerator, const NumType &bDenominator):
        aMin (tRounding<NumType>::div_down (aNumerator, aDenominator)),
        aMax (tRounding<NumType>::div_up (aNumerator, aDenominator)),
        bMin (tRounding<NumType>::div_down (bNumerator, bDenominator)),
        bMax (tRounding<NumType>::div_up (bNumerator, bDenominator))
{
        if ((aMin < 0) || (bMin < 0))
                throw "Trying to use negative parameters of the Henon map.";
        return;
} /* tMapHenon::tMapHenon */

template <class NumType>
inline int tMapHenon<NumType>::dim ()
{
        return 2;
} /* tMapHenon::dim */

template <class NumType>
inline void tMapHenon<NumType>::initialPoint (NumType *x, int iterCount)
        const
{
        // iterate the point by the Henon map: (x,y) -> (1 + y - ax^2, bx)
        for (int i = 0; i < iterCount; ++ i)
        {
                NumType tmp = 1 + x [1] - aMin * x [0] * x [0];
                x [1] = bMin * x [0];
                x [0] = tmp;
        }
        return;
} /* tMapHenon::initialPoint */

template <class NumType>
template <class InputType, class ResultType>
void tMapHenon<NumType>::image
        (const InputType *xMin, const InputType *xMax,
        ResultType *yMin, ResultType *yMax) const
{
        // extract the rounding class to shorten the notation
        typedef tRounding<ResultType> rnd;

        // compute elements of the formulas for the image
        ResultType x2Min, x2Max;
        openSqr (xMin [0], xMax [0], x2Min, x2Max);
        ResultType ax2Min = rnd::mul_down (aMin, x2Min);
        ResultType ax2Max = rnd::mul_up (aMax, x2Max);
        ResultType y1Min = rnd::add_down (1, xMin [1]);
        ResultType y1Max = rnd::add_up (1, xMax [1]);
        ResultType bxMin = rnd::mul_down (bMin, xMin [0]);
        ResultType bxMax = rnd::mul_up (bMax, xMax [0]);

        // compute the final formulas of the image of the rectangle
        yMin [0] = rnd::sub_down (y1Min, ax2Max);
        yMax [0] = rnd::sub_up (y1Max, ax2Min);
        yMin [1] = bxMin;
        yMax [1] = bxMax;

        return;
} /* tMapHenon::image */


#endif // _FINRESDYN_MAPHENON_H_