RungeKutta4.hpp

Go to the documentation of this file.

#pragma ident "$Id: RungeKutta4.hpp 1957 2009-06-19 14:39:39Z afarris $"



//============================================================================
//
//  This file is part of GPSTk, the GPS Toolkit.
//
//  The GPSTk is free software; you can redistribute it and/or modify
//  it under the terms of the GNU Lesser General Public License as published
//  by the Free Software Foundation; either version 2.1 of the License, or
//  any later version.
//
//  The GPSTk 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 Lesser General Public License for more details.
//
//  You should have received a copy of the GNU Lesser General Public
//  License along with GPSTk; if not, write to the Free Software Foundation,
//  Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
//  
//  Copyright 2004, The University of Texas at Austin
//
//============================================================================


/*
 * @file RungeKutta4.hpp
 * RungeKutta4 class. Provides several versions of the Runge Kutta integration
 * The algorithms are based on Ch.16 of "Numerical Recipes in C" but the
 * implementation is entirely independent, i.e., this isn't the source 
 * distributed with the text.
 */

#ifndef GPSTK_RUNGEKUTTA4_H
#define GPSTK_RUNGEKUTTA4_H

#include "Matrix.hpp"

namespace gpstk
{


   class RungeKutta4 
   {
   public:
      RungeKutta4(const Matrix<double>& initialState,
                  double initialTime=0,
                  double timeEpsilon=1e-18)
            : currentState(initialState), currentTime(initialTime), 
              teps(timeEpsilon), M(initialState.rows()), N(initialState.cols()),
              k1(M,N), k2(M,N), k3(M,N), k4(M,N), yn(M,N), tempy(M,N)
         { }     

      void integrateTo (double nextTime, 
                        double stepSize = 0);

      void integrateTo (double nextTime,
                        Matrix<double>& error,
                        double stepSize = 0);

      virtual gpstk::Matrix<double>& 
      derivative(long double time,
                 const gpstk::Matrix<double>& inState,
                 gpstk::Matrix<double>& inStateDot) = 0;

      double getTime(void) 
      { return currentTime; }

      const Matrix<double>& getState(void) 
      { return currentState; }

   protected:
      
      double currentTime;

      gpstk::Matrix<double> currentState;

      double teps;   //< Precision for time calculations and comparisons
      int M;         //< Number of rows in the state
      int N;         //< Number of columns in the state

   private:

         RungeKutta4(const RungeKutta4& cloneDonor);

         RungeKutta4& operator= (const RungeKutta4& right);

         Matrix<double> k1, k2, k3, k4, yn, tempy;

   }; // end class RungeKutta4


} // end namespace gpstk

#endif

---