EllipticFunction.hpp

Go to the documentation of this file.
00001 /**
00002  * \file EllipticFunction.hpp
00003  * \brief Header for GeographicLib::EllipticFunction class
00004  *
00005  * Copyright (c) Charles Karney (2008, 2009) <charles@karney.com>
00006  * and licensed under the LGPL.  For more information, see
00007  * http://geographiclib.sourceforge.net/
00008  **********************************************************************/
00009 
00010 #if !defined(GEOGRAPHICLIB_ELLIPTICFUNCTION_HPP)
00011 #define GEOGRAPHICLIB_ELLIPTICFUNCTION_HPP "$Id: EllipticFunction.hpp 6866 2010-09-11 02:15:29Z karney $"
00012 
00013 #include "GeographicLib/Constants.hpp"
00014 
00015 namespace GeographicLib {
00016 
00017   /**
00018    * \brief Elliptic functions needed for TransverseMercatorExact
00019    *
00020    * This provides the subset of elliptic functions needed for
00021    * TransverseMercatorExact.  For a given ellipsoid, only parameters \e
00022    * e<sup>2</sup> and 1 - \e e<sup>2</sup> are needed.  This class taken the
00023    * parameter as a constructor parameters and caches the values of the
00024    * required complete integrals.  A method is provided for Jacobi elliptic
00025    * functions and for the incomplete elliptic integral of the second kind in
00026    * terms of the amplitude.
00027    *
00028    * The computation of the elliptic integrals uses the algorithms given in
00029    * - B. C. Carlson,
00030    *   <a href="http://dx.doi.org/10.1007/BF02198293"> Computation of elliptic
00031    *   integrals</a>, Numerical Algorithms 10, 13&ndash;26 (1995).
00032    * .
00033    * The computation of the Jacobi elliptic functions uses the algorithm given
00034    * in
00035    * - R. Bulirsch,
00036    *   <a href="http://dx.doi.org/10.1007/BF01397975"> Numerical Calculation of
00037    *   Elliptic Integrals and Elliptic Functions</a>, Numericshe Mathematik 7,
00038    *   78&ndash;90 (1965).
00039    * .
00040    * The notation follows Abramowitz and Stegun, Chapters 16 and 17.
00041    **********************************************************************/
00042   class EllipticFunction {
00043   private:
00044     typedef Math::real real;
00045     static const real tol, tolRF, tolRD, tolRG0, tolJAC, tolJAC1;
00046     enum { num = 10 }; // Max depth required for sncndn.  Probably 5 is enough.
00047     static real RF(real x, real y, real z) throw();
00048     static real RD(real x, real y, real z) throw();
00049     static real RG0(real x, real y) throw();
00050     const real _m, _m1;
00051     mutable bool _init;
00052     mutable real _kc, _ec, _kec;
00053     bool Init() const throw();
00054   public:
00055 
00056     /**
00057      * Constructor.
00058      *
00059      * @param[in] m the parameter which must lie in [0, 1].  (No checking
00060      *   is done.)
00061      **********************************************************************/
00062     explicit EllipticFunction(real m) throw();
00063 
00064     /**
00065      * @return the parameter \e m.
00066      **********************************************************************/
00067     Math::real m() const throw() { return _m; }
00068 
00069     /**
00070      * @return the complementary parameter \e m' = (1 - \e m).
00071      **********************************************************************/
00072     Math::real m1() const throw() { return _m1; }
00073 
00074     /**
00075      * @return the complete integral of first kind, \e K(\e m).
00076      **********************************************************************/
00077     Math::real K() const throw() { _init || Init(); return _kc; }
00078 
00079     /**
00080      * @return the complete integral of second kind, \e E(\e m).
00081      **********************************************************************/
00082     Math::real E() const throw() { _init || Init(); return _ec; }
00083 
00084     /**
00085      * @return the difference \e K(\e m) - \e E(\e m) (which can be computed
00086      *   directly).
00087      **********************************************************************/
00088     Math::real KE() const throw() { _init || Init(); return _kec; }
00089 
00090     /**
00091      * The Jacobi elliptic functions.
00092      *
00093      * @param[in] x the argument.
00094      * @param[out] sn sn(<i>x</i>|<i>m</i>).
00095      * @param[out] cn cn(<i>x</i>|<i>m</i>).
00096      * @param[out] dn dn(<i>x</i>|<i>m</i>).
00097      **********************************************************************/
00098     void sncndn(real x, real& sn, real& cn, real& dn) const throw();
00099 
00100     /**
00101      * The incomplete integral of the second kind.
00102      *
00103      * @param[in] phi
00104      * @return int sqrt(1 -  \e m sin<sup>2</sup>(\e phi)) \e dphi.
00105      **********************************************************************/
00106     Math::real E(real phi) const throw();
00107 
00108     /**
00109      * The incomplete integral of the second kind in terms of Jacobi elliptic
00110      * functions
00111      *
00112      * @param[in] sn
00113      * @param[in] cn
00114      * @param[in] dn
00115      * @return int dn(\e w)<sup>2</sup> \e dw (A+S 17.2.10).
00116      *
00117      * Instead of specifying the ampltiude \e phi, we provide \e sn = sin(\e
00118      * phi), \e cn = cos(\e phi), \e dn = sqrt(1 - \e m sin<sup>2</sup>(\e
00119      * phi)).
00120      **********************************************************************/
00121     Math::real E(real sn, real cn, real dn) const throw();
00122   };
00123 
00124 
00125 } // namespace GeographicLib
00126 
00127 #endif