epoc32/include/stdapis/boost/math/quaternion.hpp
branchSymbian2
changeset 2 2fe1408b6811
     1.1 --- /dev/null	Thu Jan 01 00:00:00 1970 +0000
     1.2 +++ b/epoc32/include/stdapis/boost/math/quaternion.hpp	Tue Mar 16 16:12:26 2010 +0000
     1.3 @@ -0,0 +1,1924 @@
     1.4 +//  boost quaternion.hpp header file
     1.5 +
     1.6 +//  (C) Copyright Hubert Holin 2001.
     1.7 +//  Distributed under the Boost Software License, Version 1.0. (See
     1.8 +//  accompanying file LICENSE_1_0.txt or copy at
     1.9 +//  http://www.boost.org/LICENSE_1_0.txt)
    1.10 +
    1.11 +// See http://www.boost.org for updates, documentation, and revision history.
    1.12 +
    1.13 +#ifndef BOOST_QUATERNION_HPP
    1.14 +#define BOOST_QUATERNION_HPP
    1.15 +
    1.16 +
    1.17 +#include <complex>
    1.18 +#include <iosfwd>                                    // for the "<<" and ">>" operators
    1.19 +#include <sstream>                                    // for the "<<" operator
    1.20 +
    1.21 +#include <boost/config.hpp> // for BOOST_NO_STD_LOCALE
    1.22 +#include <boost/detail/workaround.hpp>
    1.23 +#ifndef    BOOST_NO_STD_LOCALE
    1.24 +    #include <locale>                                    // for the "<<" operator
    1.25 +#endif /* BOOST_NO_STD_LOCALE */
    1.26 +
    1.27 +#include <valarray>
    1.28 +
    1.29 +
    1.30 +
    1.31 +#include <boost/math/special_functions/sinc.hpp>    // for the Sinus cardinal
    1.32 +#include <boost/math/special_functions/sinhc.hpp>    // for the Hyperbolic Sinus cardinal
    1.33 +
    1.34 +
    1.35 +namespace boost
    1.36 +{
    1.37 +    namespace math
    1.38 +    {
    1.39 +#if     BOOST_WORKAROUND(__GNUC__, < 3)
    1.40 +        // gcc 2.95.x uses expression templates for valarray calculations, but
    1.41 +        // the result is not conforming. We need BOOST_GET_VALARRAY to get an
    1.42 +        // actual valarray result when we need to call a member function
    1.43 +    #define    BOOST_GET_VALARRAY(T,x)    ::std::valarray<T>(x)
    1.44 +        // gcc 2.95.x has an "std::ios" class that is similar to 
    1.45 +        // "std::ios_base", so we just use a #define
    1.46 +    #define    BOOST_IOS_BASE    ::std::ios
    1.47 +        // gcc 2.x ignores function scope using declarations,
    1.48 +        // put them in the scope of the enclosing namespace instead:
    1.49 +        using    ::std::valarray;
    1.50 +        using    ::std::sqrt;
    1.51 +        using    ::std::cos;
    1.52 +        using    ::std::sin;
    1.53 +        using    ::std::exp;
    1.54 +        using    ::std::cosh;
    1.55 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
    1.56 +
    1.57 +#define    BOOST_QUATERNION_ACCESSOR_GENERATOR(type)                    \
    1.58 +            type                    real() const                        \
    1.59 +            {                                                           \
    1.60 +                return(a);                                              \
    1.61 +            }                                                           \
    1.62 +                                                                        \
    1.63 +            quaternion<type>        unreal() const                      \
    1.64 +            {                                                           \
    1.65 +                return(quaternion<type>(static_cast<type>(0),b,c,d));   \
    1.66 +            }                                                           \
    1.67 +                                                                        \
    1.68 +            type                    R_component_1() const               \
    1.69 +            {                                                           \
    1.70 +                return(a);                                              \
    1.71 +            }                                                           \
    1.72 +                                                                        \
    1.73 +            type                    R_component_2() const               \
    1.74 +            {                                                           \
    1.75 +                return(b);                                              \
    1.76 +            }                                                           \
    1.77 +                                                                        \
    1.78 +            type                    R_component_3() const               \
    1.79 +            {                                                           \
    1.80 +                return(c);                                              \
    1.81 +            }                                                           \
    1.82 +                                                                        \
    1.83 +            type                    R_component_4() const               \
    1.84 +            {                                                           \
    1.85 +                return(d);                                              \
    1.86 +            }                                                           \
    1.87 +                                                                        \
    1.88 +            ::std::complex<type>    C_component_1() const               \
    1.89 +            {                                                           \
    1.90 +                return(::std::complex<type>(a,b));                      \
    1.91 +            }                                                           \
    1.92 +                                                                        \
    1.93 +            ::std::complex<type>    C_component_2() const               \
    1.94 +            {                                                           \
    1.95 +                return(::std::complex<type>(c,d));                      \
    1.96 +            }
    1.97 +        
    1.98 +        
    1.99 +#define    BOOST_QUATERNION_MEMBER_ASSIGNMENT_GENERATOR(type)                               \
   1.100 +            template<typename X>                                                            \
   1.101 +            quaternion<type> &        operator = (quaternion<X> const  & a_affecter)        \
   1.102 +            {                                                                               \
   1.103 +                a = static_cast<type>(a_affecter.R_component_1());                          \
   1.104 +                b = static_cast<type>(a_affecter.R_component_2());                          \
   1.105 +                c = static_cast<type>(a_affecter.R_component_3());                          \
   1.106 +                d = static_cast<type>(a_affecter.R_component_4());                          \
   1.107 +                                                                                            \
   1.108 +                return(*this);                                                              \
   1.109 +            }                                                                               \
   1.110 +                                                                                            \
   1.111 +            quaternion<type> &        operator = (quaternion<type> const & a_affecter)      \
   1.112 +            {                                                                               \
   1.113 +                a = a_affecter.a;                                                           \
   1.114 +                b = a_affecter.b;                                                           \
   1.115 +                c = a_affecter.c;                                                           \
   1.116 +                d = a_affecter.d;                                                           \
   1.117 +                                                                                            \
   1.118 +                return(*this);                                                              \
   1.119 +            }                                                                               \
   1.120 +                                                                                            \
   1.121 +            quaternion<type> &        operator = (type const & a_affecter)                  \
   1.122 +            {                                                                               \
   1.123 +                a = a_affecter;                                                             \
   1.124 +                                                                                            \
   1.125 +                b = c = d = static_cast<type>(0);                                           \
   1.126 +                                                                                            \
   1.127 +                return(*this);                                                              \
   1.128 +            }                                                                               \
   1.129 +                                                                                            \
   1.130 +            quaternion<type> &        operator = (::std::complex<type> const & a_affecter)  \
   1.131 +            {                                                                               \
   1.132 +                a = a_affecter.real();                                                      \
   1.133 +                b = a_affecter.imag();                                                      \
   1.134 +                                                                                            \
   1.135 +                c = d = static_cast<type>(0);                                               \
   1.136 +                                                                                            \
   1.137 +                return(*this);                                                              \
   1.138 +            }
   1.139 +        
   1.140 +        
   1.141 +#define    BOOST_QUATERNION_MEMBER_DATA_GENERATOR(type)       \
   1.142 +            type    a;                                        \
   1.143 +            type    b;                                        \
   1.144 +            type    c;                                        \
   1.145 +            type    d;
   1.146 +        
   1.147 +        
   1.148 +        template<typename T>
   1.149 +        class quaternion
   1.150 +        {
   1.151 +        public:
   1.152 +            
   1.153 +            typedef T value_type;
   1.154 +            
   1.155 +            
   1.156 +            // constructor for H seen as R^4
   1.157 +            // (also default constructor)
   1.158 +            
   1.159 +            explicit            quaternion( T const & requested_a = T(),
   1.160 +                                            T const & requested_b = T(),
   1.161 +                                            T const & requested_c = T(),
   1.162 +                                            T const & requested_d = T())
   1.163 +            :   a(requested_a),
   1.164 +                b(requested_b),
   1.165 +                c(requested_c),
   1.166 +                d(requested_d)
   1.167 +            {
   1.168 +                // nothing to do!
   1.169 +            }
   1.170 +            
   1.171 +            
   1.172 +            // constructor for H seen as C^2
   1.173 +                
   1.174 +            explicit            quaternion( ::std::complex<T> const & z0,
   1.175 +                                            ::std::complex<T> const & z1 = ::std::complex<T>())
   1.176 +            :   a(z0.real()),
   1.177 +                b(z0.imag()),
   1.178 +                c(z1.real()),
   1.179 +                d(z1.imag())
   1.180 +            {
   1.181 +                // nothing to do!
   1.182 +            }
   1.183 +            
   1.184 +            
   1.185 +            // UNtemplated copy constructor
   1.186 +            // (this is taken care of by the compiler itself)
   1.187 +            
   1.188 +            
   1.189 +            // templated copy constructor
   1.190 +            
   1.191 +            template<typename X>
   1.192 +            explicit            quaternion(quaternion<X> const & a_recopier)
   1.193 +            :   a(static_cast<T>(a_recopier.R_component_1())),
   1.194 +                b(static_cast<T>(a_recopier.R_component_2())),
   1.195 +                c(static_cast<T>(a_recopier.R_component_3())),
   1.196 +                d(static_cast<T>(a_recopier.R_component_4()))
   1.197 +            {
   1.198 +                // nothing to do!
   1.199 +            }
   1.200 +            
   1.201 +            
   1.202 +            // destructor
   1.203 +            // (this is taken care of by the compiler itself)
   1.204 +            
   1.205 +            
   1.206 +            // accessors
   1.207 +            //
   1.208 +            // Note:    Like complex number, quaternions do have a meaningful notion of "real part",
   1.209 +            //            but unlike them there is no meaningful notion of "imaginary part".
   1.210 +            //            Instead there is an "unreal part" which itself is a quaternion, and usually
   1.211 +            //            nothing simpler (as opposed to the complex number case).
   1.212 +            //            However, for practicallity, there are accessors for the other components
   1.213 +            //            (these are necessary for the templated copy constructor, for instance).
   1.214 +            
   1.215 +            BOOST_QUATERNION_ACCESSOR_GENERATOR(T)
   1.216 +            
   1.217 +            // assignment operators
   1.218 +            
   1.219 +            BOOST_QUATERNION_MEMBER_ASSIGNMENT_GENERATOR(T)
   1.220 +            
   1.221 +            // other assignment-related operators
   1.222 +            //
   1.223 +            // NOTE:    Quaternion multiplication is *NOT* commutative;
   1.224 +            //            symbolically, "q *= rhs;" means "q = q * rhs;"
   1.225 +            //            and "q /= rhs;" means "q = q * inverse_of(rhs);"
   1.226 +            
   1.227 +            quaternion<T> &        operator += (T const & rhs)
   1.228 +            {
   1.229 +                T    at = a + rhs;    // exception guard
   1.230 +                
   1.231 +                a = at;
   1.232 +                
   1.233 +                return(*this);
   1.234 +            }
   1.235 +            
   1.236 +            
   1.237 +            quaternion<T> &        operator += (::std::complex<T> const & rhs)
   1.238 +            {
   1.239 +                T    at = a + rhs.real();    // exception guard
   1.240 +                T    bt = b + rhs.imag();    // exception guard
   1.241 +                
   1.242 +                a = at; 
   1.243 +                b = bt;
   1.244 +                
   1.245 +                return(*this);
   1.246 +            }
   1.247 +            
   1.248 +            
   1.249 +            template<typename X>
   1.250 +            quaternion<T> &        operator += (quaternion<X> const & rhs)
   1.251 +            {
   1.252 +                T    at = a + static_cast<T>(rhs.R_component_1());    // exception guard
   1.253 +                T    bt = b + static_cast<T>(rhs.R_component_2());    // exception guard
   1.254 +                T    ct = c + static_cast<T>(rhs.R_component_3());    // exception guard
   1.255 +                T    dt = d + static_cast<T>(rhs.R_component_4());    // exception guard
   1.256 +                
   1.257 +                a = at;
   1.258 +                b = bt;
   1.259 +                c = ct;
   1.260 +                d = dt;
   1.261 +                
   1.262 +                return(*this);
   1.263 +            }
   1.264 +            
   1.265 +            
   1.266 +            
   1.267 +            quaternion<T> &        operator -= (T const & rhs)
   1.268 +            {
   1.269 +                T    at = a - rhs;    // exception guard
   1.270 +                
   1.271 +                a = at;
   1.272 +                
   1.273 +                return(*this);
   1.274 +            }
   1.275 +            
   1.276 +            
   1.277 +            quaternion<T> &        operator -= (::std::complex<T> const & rhs)
   1.278 +            {
   1.279 +                T    at = a - rhs.real();    // exception guard
   1.280 +                T    bt = b - rhs.imag();    // exception guard
   1.281 +                
   1.282 +                a = at;
   1.283 +                b = bt;
   1.284 +                
   1.285 +                return(*this);
   1.286 +            }
   1.287 +            
   1.288 +            
   1.289 +            template<typename X>
   1.290 +            quaternion<T> &        operator -= (quaternion<X> const & rhs)
   1.291 +            {
   1.292 +                T    at = a - static_cast<T>(rhs.R_component_1());    // exception guard
   1.293 +                T    bt = b - static_cast<T>(rhs.R_component_2());    // exception guard
   1.294 +                T    ct = c - static_cast<T>(rhs.R_component_3());    // exception guard
   1.295 +                T    dt = d - static_cast<T>(rhs.R_component_4());    // exception guard
   1.296 +                
   1.297 +                a = at;
   1.298 +                b = bt;
   1.299 +                c = ct;
   1.300 +                d = dt;
   1.301 +                
   1.302 +                return(*this);
   1.303 +            }
   1.304 +            
   1.305 +            
   1.306 +            quaternion<T> &        operator *= (T const & rhs)
   1.307 +            {
   1.308 +                T    at = a * rhs;    // exception guard
   1.309 +                T    bt = b * rhs;    // exception guard
   1.310 +                T    ct = c * rhs;    // exception guard
   1.311 +                T    dt = d * rhs;    // exception guard
   1.312 +                
   1.313 +                a = at;
   1.314 +                b = bt;
   1.315 +                c = ct;
   1.316 +                d = dt;
   1.317 +                
   1.318 +                return(*this);
   1.319 +            }
   1.320 +            
   1.321 +            
   1.322 +            quaternion<T> &        operator *= (::std::complex<T> const & rhs)
   1.323 +            {
   1.324 +                T    ar = rhs.real();
   1.325 +                T    br = rhs.imag();
   1.326 +                
   1.327 +                T    at = +a*ar-b*br;
   1.328 +                T    bt = +a*br+b*ar;
   1.329 +                T    ct = +c*ar+d*br;
   1.330 +                T    dt = -c*br+d*ar;
   1.331 +                
   1.332 +                a = at;
   1.333 +                b = bt;
   1.334 +                c = ct;
   1.335 +                d = dt;
   1.336 +                
   1.337 +                return(*this);
   1.338 +            }
   1.339 +            
   1.340 +            
   1.341 +            template<typename X>
   1.342 +            quaternion<T> &        operator *= (quaternion<X> const & rhs)
   1.343 +            {
   1.344 +                T    ar = static_cast<T>(rhs.R_component_1());
   1.345 +                T    br = static_cast<T>(rhs.R_component_2());
   1.346 +                T    cr = static_cast<T>(rhs.R_component_3());
   1.347 +                T    dr = static_cast<T>(rhs.R_component_4());
   1.348 +                
   1.349 +                T    at = +a*ar-b*br-c*cr-d*dr;
   1.350 +                T    bt = +a*br+b*ar+c*dr-d*cr;    //(a*br+ar*b)+(c*dr-cr*d);
   1.351 +                T    ct = +a*cr-b*dr+c*ar+d*br;    //(a*cr+ar*c)+(d*br-dr*b);
   1.352 +                T    dt = +a*dr+b*cr-c*br+d*ar;    //(a*dr+ar*d)+(b*cr-br*c);
   1.353 +                
   1.354 +                a = at;
   1.355 +                b = bt;
   1.356 +                c = ct;
   1.357 +                d = dt;
   1.358 +                
   1.359 +                return(*this);
   1.360 +            }
   1.361 +            
   1.362 +            
   1.363 +            
   1.364 +            quaternion<T> &        operator /= (T const & rhs)
   1.365 +            {
   1.366 +                T    at = a / rhs;    // exception guard
   1.367 +                T    bt = b / rhs;    // exception guard
   1.368 +                T    ct = c / rhs;    // exception guard
   1.369 +                T    dt = d / rhs;    // exception guard
   1.370 +                
   1.371 +                a = at;
   1.372 +                b = bt;
   1.373 +                c = ct;
   1.374 +                d = dt;
   1.375 +                
   1.376 +                return(*this);
   1.377 +            }
   1.378 +            
   1.379 +            
   1.380 +            quaternion<T> &        operator /= (::std::complex<T> const & rhs)
   1.381 +            {
   1.382 +                T    ar = rhs.real();
   1.383 +                T    br = rhs.imag();
   1.384 +                
   1.385 +                T    denominator = ar*ar+br*br;
   1.386 +                
   1.387 +                T    at = (+a*ar+b*br)/denominator;    //(a*ar+b*br)/denominator;
   1.388 +                T    bt = (-a*br+b*ar)/denominator;    //(ar*b-a*br)/denominator;
   1.389 +                T    ct = (+c*ar-d*br)/denominator;    //(ar*c-d*br)/denominator;
   1.390 +                T    dt = (+c*br+d*ar)/denominator;    //(ar*d+br*c)/denominator;
   1.391 +                
   1.392 +                a = at;
   1.393 +                b = bt;
   1.394 +                c = ct;
   1.395 +                d = dt;
   1.396 +                
   1.397 +                return(*this);
   1.398 +            }
   1.399 +            
   1.400 +            
   1.401 +            template<typename X>
   1.402 +            quaternion<T> &        operator /= (quaternion<X> const & rhs)
   1.403 +            {
   1.404 +                T    ar = static_cast<T>(rhs.R_component_1());
   1.405 +                T    br = static_cast<T>(rhs.R_component_2());
   1.406 +                T    cr = static_cast<T>(rhs.R_component_3());
   1.407 +                T    dr = static_cast<T>(rhs.R_component_4());
   1.408 +                
   1.409 +                T    denominator = ar*ar+br*br+cr*cr+dr*dr;
   1.410 +                
   1.411 +                T    at = (+a*ar+b*br+c*cr+d*dr)/denominator;    //(a*ar+b*br+c*cr+d*dr)/denominator;
   1.412 +                T    bt = (-a*br+b*ar-c*dr+d*cr)/denominator;    //((ar*b-a*br)+(cr*d-c*dr))/denominator;
   1.413 +                T    ct = (-a*cr+b*dr+c*ar-d*br)/denominator;    //((ar*c-a*cr)+(dr*b-d*br))/denominator;
   1.414 +                T    dt = (-a*dr-b*cr+c*br+d*ar)/denominator;    //((ar*d-a*dr)+(br*c-b*cr))/denominator;
   1.415 +                
   1.416 +                a = at;
   1.417 +                b = bt;
   1.418 +                c = ct;
   1.419 +                d = dt;
   1.420 +                
   1.421 +                return(*this);
   1.422 +            }
   1.423 +            
   1.424 +            
   1.425 +        protected:
   1.426 +            
   1.427 +            BOOST_QUATERNION_MEMBER_DATA_GENERATOR(T)
   1.428 +            
   1.429 +            
   1.430 +        private:
   1.431 +            
   1.432 +        };
   1.433 +        
   1.434 +        
   1.435 +        // declaration of quaternion specialization
   1.436 +        
   1.437 +        template<>    class quaternion<float>;
   1.438 +        template<>    class quaternion<double>;
   1.439 +        template<>    class quaternion<long double>;
   1.440 +        
   1.441 +        
   1.442 +        // helper templates for converting copy constructors (declaration)
   1.443 +        
   1.444 +        namespace detail
   1.445 +        {
   1.446 +            
   1.447 +            template<   typename T,
   1.448 +                        typename U
   1.449 +                    >
   1.450 +            quaternion<T>    quaternion_type_converter(quaternion<U> const & rhs);
   1.451 +        }
   1.452 +        
   1.453 +        
   1.454 +        // implementation of quaternion specialization
   1.455 +        
   1.456 +        
   1.457 +#define    BOOST_QUATERNION_CONSTRUCTOR_GENERATOR(type)                                                 \
   1.458 +            explicit            quaternion( type const & requested_a = static_cast<type>(0),            \
   1.459 +                                            type const & requested_b = static_cast<type>(0),            \
   1.460 +                                            type const & requested_c = static_cast<type>(0),            \
   1.461 +                                            type const & requested_d = static_cast<type>(0))            \
   1.462 +            :   a(requested_a),                                                                         \
   1.463 +                b(requested_b),                                                                         \
   1.464 +                c(requested_c),                                                                         \
   1.465 +                d(requested_d)                                                                          \
   1.466 +            {                                                                                           \
   1.467 +            }                                                                                           \
   1.468 +                                                                                                        \
   1.469 +            explicit            quaternion( ::std::complex<type> const & z0,                            \
   1.470 +                                            ::std::complex<type> const & z1 = ::std::complex<type>())   \
   1.471 +            :   a(z0.real()),                                                                           \
   1.472 +                b(z0.imag()),                                                                           \
   1.473 +                c(z1.real()),                                                                           \
   1.474 +                d(z1.imag())                                                                            \
   1.475 +            {                                                                                           \
   1.476 +            }
   1.477 +        
   1.478 +        
   1.479 +#define    BOOST_QUATERNION_MEMBER_ADD_GENERATOR_1(type)             \
   1.480 +            quaternion<type> &        operator += (type const & rhs) \
   1.481 +            {                                                        \
   1.482 +                a += rhs;                                            \
   1.483 +                                                                     \
   1.484 +                return(*this);                                       \
   1.485 +            }
   1.486 +    
   1.487 +#define    BOOST_QUATERNION_MEMBER_ADD_GENERATOR_2(type)                             \
   1.488 +            quaternion<type> &        operator += (::std::complex<type> const & rhs) \
   1.489 +            {                                                                        \
   1.490 +                a += rhs.real();                                                     \
   1.491 +                b += rhs.imag();                                                     \
   1.492 +                                                                                     \
   1.493 +                return(*this);                                                       \
   1.494 +            }
   1.495 +    
   1.496 +#define    BOOST_QUATERNION_MEMBER_ADD_GENERATOR_3(type)                      \
   1.497 +            template<typename X>                                              \
   1.498 +            quaternion<type> &        operator += (quaternion<X> const & rhs) \
   1.499 +            {                                                                 \
   1.500 +                a += static_cast<type>(rhs.R_component_1());                  \
   1.501 +                b += static_cast<type>(rhs.R_component_2());                  \
   1.502 +                c += static_cast<type>(rhs.R_component_3());                  \
   1.503 +                d += static_cast<type>(rhs.R_component_4());                  \
   1.504 +                                                                              \
   1.505 +                return(*this);                                                \
   1.506 +            }
   1.507 +    
   1.508 +#define    BOOST_QUATERNION_MEMBER_SUB_GENERATOR_1(type)             \
   1.509 +            quaternion<type> &        operator -= (type const & rhs) \
   1.510 +            {                                                        \
   1.511 +                a -= rhs;                                            \
   1.512 +                                                                     \
   1.513 +                return(*this);                                       \
   1.514 +            }
   1.515 +    
   1.516 +#define    BOOST_QUATERNION_MEMBER_SUB_GENERATOR_2(type)                             \
   1.517 +            quaternion<type> &        operator -= (::std::complex<type> const & rhs) \
   1.518 +            {                                                                        \
   1.519 +                a -= rhs.real();                                                     \
   1.520 +                b -= rhs.imag();                                                     \
   1.521 +                                                                                     \
   1.522 +                return(*this);                                                       \
   1.523 +            }
   1.524 +    
   1.525 +#define    BOOST_QUATERNION_MEMBER_SUB_GENERATOR_3(type)                      \
   1.526 +            template<typename X>                                              \
   1.527 +            quaternion<type> &        operator -= (quaternion<X> const & rhs) \
   1.528 +            {                                                                 \
   1.529 +                a -= static_cast<type>(rhs.R_component_1());                  \
   1.530 +                b -= static_cast<type>(rhs.R_component_2());                  \
   1.531 +                c -= static_cast<type>(rhs.R_component_3());                  \
   1.532 +                d -= static_cast<type>(rhs.R_component_4());                  \
   1.533 +                                                                              \
   1.534 +                return(*this);                                                \
   1.535 +            }
   1.536 +    
   1.537 +#define    BOOST_QUATERNION_MEMBER_MUL_GENERATOR_1(type)             \
   1.538 +            quaternion<type> &        operator *= (type const & rhs) \
   1.539 +            {                                                        \
   1.540 +                a *= rhs;                                            \
   1.541 +                b *= rhs;                                            \
   1.542 +                c *= rhs;                                            \
   1.543 +                d *= rhs;                                            \
   1.544 +                                                                     \
   1.545 +                return(*this);                                       \
   1.546 +            }
   1.547 +    
   1.548 +#define    BOOST_QUATERNION_MEMBER_MUL_GENERATOR_2(type)                             \
   1.549 +            quaternion<type> &        operator *= (::std::complex<type> const & rhs) \
   1.550 +            {                                                                        \
   1.551 +                type    ar = rhs.real();                                             \
   1.552 +                type    br = rhs.imag();                                             \
   1.553 +                                                                                     \
   1.554 +                type    at = +a*ar-b*br;                                             \
   1.555 +                type    bt = +a*br+b*ar;                                             \
   1.556 +                type    ct = +c*ar+d*br;                                             \
   1.557 +                type    dt = -c*br+d*ar;                                             \
   1.558 +                                                                                     \
   1.559 +                a = at;                                                              \
   1.560 +                b = bt;                                                              \
   1.561 +                c = ct;                                                              \
   1.562 +                d = dt;                                                              \
   1.563 +                                                                                     \
   1.564 +                return(*this);                                                       \
   1.565 +            }
   1.566 +    
   1.567 +#define    BOOST_QUATERNION_MEMBER_MUL_GENERATOR_3(type)                      \
   1.568 +            template<typename X>                                              \
   1.569 +            quaternion<type> &        operator *= (quaternion<X> const & rhs) \
   1.570 +            {                                                                 \
   1.571 +                type    ar = static_cast<type>(rhs.R_component_1());          \
   1.572 +                type    br = static_cast<type>(rhs.R_component_2());          \
   1.573 +                type    cr = static_cast<type>(rhs.R_component_3());          \
   1.574 +                type    dr = static_cast<type>(rhs.R_component_4());          \
   1.575 +                                                                              \
   1.576 +                type    at = +a*ar-b*br-c*cr-d*dr;                            \
   1.577 +                type    bt = +a*br+b*ar+c*dr-d*cr;                            \
   1.578 +                type    ct = +a*cr-b*dr+c*ar+d*br;                            \
   1.579 +                type    dt = +a*dr+b*cr-c*br+d*ar;                            \
   1.580 +                                                                              \
   1.581 +                a = at;                                                       \
   1.582 +                b = bt;                                                       \
   1.583 +                c = ct;                                                       \
   1.584 +                d = dt;                                                       \
   1.585 +                                                                              \
   1.586 +                return(*this);                                                \
   1.587 +            }
   1.588 +    
   1.589 +// There is quite a lot of repetition in the code below. This is intentional.
   1.590 +// The last conditional block is the normal form, and the others merely
   1.591 +// consist of workarounds for various compiler deficiencies. Hopefuly, when
   1.592 +// more compilers are conformant and we can retire support for those that are
   1.593 +// not, we will be able to remove the clutter. This is makes the situation
   1.594 +// (painfully) explicit.
   1.595 +    
   1.596 +#define    BOOST_QUATERNION_MEMBER_DIV_GENERATOR_1(type)             \
   1.597 +            quaternion<type> &        operator /= (type const & rhs) \
   1.598 +            {                                                        \
   1.599 +                a /= rhs;                                            \
   1.600 +                b /= rhs;                                            \
   1.601 +                c /= rhs;                                            \
   1.602 +                d /= rhs;                                            \
   1.603 +                                                                     \
   1.604 +                return(*this);                                       \
   1.605 +            }
   1.606 +
   1.607 +#if defined(__GNUC__) && (__GNUC__ < 3)
   1.608 +    #define    BOOST_QUATERNION_MEMBER_DIV_GENERATOR_2(type)                                            \
   1.609 +            quaternion<type> &        operator /= (::std::complex<type> const & rhs)                    \
   1.610 +            {                                                                                           \
   1.611 +                using    ::std::valarray;                                                               \
   1.612 +                                                                                                        \
   1.613 +                valarray<type>    tr(2);                                                                \
   1.614 +                                                                                                        \
   1.615 +                tr[0] = rhs.real();                                                                     \
   1.616 +                tr[1] = rhs.imag();                                                                     \
   1.617 +                                                                                                        \
   1.618 +                type            mixam = (BOOST_GET_VALARRAY(type,static_cast<type>(1)/abs(tr)).max)();  \
   1.619 +                                                                                                        \
   1.620 +                tr *= mixam;                                                                            \
   1.621 +                                                                                                        \
   1.622 +                valarray<type>    tt(4);                                                                \
   1.623 +                                                                                                        \
   1.624 +                tt[0] = +a*tr[0]+b*tr[1];                                                               \
   1.625 +                tt[1] = -a*tr[1]+b*tr[0];                                                               \
   1.626 +                tt[2] = +c*tr[0]-d*tr[1];                                                               \
   1.627 +                tt[3] = +c*tr[1]+d*tr[0];                                                               \
   1.628 +                                                                                                        \
   1.629 +                tr *= tr;                                                                               \
   1.630 +                                                                                                        \
   1.631 +                tt *= (mixam/tr.sum());                                                                 \
   1.632 +                                                                                                        \
   1.633 +                a = tt[0];                                                                              \
   1.634 +                b = tt[1];                                                                              \
   1.635 +                c = tt[2];                                                                              \
   1.636 +                d = tt[3];                                                                              \
   1.637 +                                                                                                        \
   1.638 +                return(*this);                                                                          \
   1.639 +            }
   1.640 +#elif    defined(BOOST_NO_ARGUMENT_DEPENDENT_LOOKUP)
   1.641 +    #define    BOOST_QUATERNION_MEMBER_DIV_GENERATOR_2(type)                         \
   1.642 +            quaternion<type> &        operator /= (::std::complex<type> const & rhs) \
   1.643 +            {                                                                        \
   1.644 +                using    ::std::valarray;                                            \
   1.645 +                using    ::std::abs;                                                 \
   1.646 +                                                                                     \
   1.647 +                valarray<type>    tr(2);                                             \
   1.648 +                                                                                     \
   1.649 +                tr[0] = rhs.real();                                                  \
   1.650 +                tr[1] = rhs.imag();                                                  \
   1.651 +                                                                                     \
   1.652 +                type            mixam = static_cast<type>(1)/(abs(tr).max)();        \
   1.653 +                                                                                     \
   1.654 +                tr *= mixam;                                                         \
   1.655 +                                                                                     \
   1.656 +                valarray<type>    tt(4);                                             \
   1.657 +                                                                                     \
   1.658 +                tt[0] = +a*tr[0]+b*tr[1];                                            \
   1.659 +                tt[1] = -a*tr[1]+b*tr[0];                                            \
   1.660 +                tt[2] = +c*tr[0]-d*tr[1];                                            \
   1.661 +                tt[3] = +c*tr[1]+d*tr[0];                                            \
   1.662 +                                                                                     \
   1.663 +                tr *= tr;                                                            \
   1.664 +                                                                                     \
   1.665 +                tt *= (mixam/tr.sum());                                              \
   1.666 +                                                                                     \
   1.667 +                a = tt[0];                                                           \
   1.668 +                b = tt[1];                                                           \
   1.669 +                c = tt[2];                                                           \
   1.670 +                d = tt[3];                                                           \
   1.671 +                                                                                     \
   1.672 +                return(*this);                                                       \
   1.673 +            }
   1.674 +#else
   1.675 +    #define    BOOST_QUATERNION_MEMBER_DIV_GENERATOR_2(type)                         \
   1.676 +            quaternion<type> &        operator /= (::std::complex<type> const & rhs) \
   1.677 +            {                                                                        \
   1.678 +                using    ::std::valarray;                                            \
   1.679 +                                                                                     \
   1.680 +                valarray<type>    tr(2);                                             \
   1.681 +                                                                                     \
   1.682 +                tr[0] = rhs.real();                                                  \
   1.683 +                tr[1] = rhs.imag();                                                  \
   1.684 +                                                                                     \
   1.685 +                type            mixam = static_cast<type>(1)/(abs(tr).max)();        \
   1.686 +                                                                                     \
   1.687 +                tr *= mixam;                                                         \
   1.688 +                                                                                     \
   1.689 +                valarray<type>    tt(4);                                             \
   1.690 +                                                                                     \
   1.691 +                tt[0] = +a*tr[0]+b*tr[1];                                            \
   1.692 +                tt[1] = -a*tr[1]+b*tr[0];                                            \
   1.693 +                tt[2] = +c*tr[0]-d*tr[1];                                            \
   1.694 +                tt[3] = +c*tr[1]+d*tr[0];                                            \
   1.695 +                                                                                     \
   1.696 +                tr *= tr;                                                            \
   1.697 +                                                                                     \
   1.698 +                tt *= (mixam/tr.sum());                                              \
   1.699 +                                                                                     \
   1.700 +                a = tt[0];                                                           \
   1.701 +                b = tt[1];                                                           \
   1.702 +                c = tt[2];                                                           \
   1.703 +                d = tt[3];                                                           \
   1.704 +                                                                                     \
   1.705 +                return(*this);                                                       \
   1.706 +            }
   1.707 +#endif /* defined(__GNUC__) && (__GNUC__ < 3) */ /* BOOST_NO_ARGUMENT_DEPENDENT_LOOKUP */
   1.708 +    
   1.709 +#if defined(__GNUC__) && (__GNUC__ < 3)
   1.710 +    #define    BOOST_QUATERNION_MEMBER_DIV_GENERATOR_3(type)                                            \
   1.711 +            template<typename X>                                                                        \
   1.712 +            quaternion<type> &        operator /= (quaternion<X> const & rhs)                           \
   1.713 +            {                                                                                           \
   1.714 +                using    ::std::valarray;                                                               \
   1.715 +                                                                                                        \
   1.716 +                valarray<type>    tr(4);                                                                \
   1.717 +                                                                                                        \
   1.718 +                tr[0] = static_cast<type>(rhs.R_component_1());                                         \
   1.719 +                tr[1] = static_cast<type>(rhs.R_component_2());                                         \
   1.720 +                tr[2] = static_cast<type>(rhs.R_component_3());                                         \
   1.721 +                tr[3] = static_cast<type>(rhs.R_component_4());                                         \
   1.722 +                                                                                                        \
   1.723 +                type            mixam = (BOOST_GET_VALARRAY(type,static_cast<type>(1)/abs(tr)).max)();  \
   1.724 +                                                                                                        \
   1.725 +                tr *= mixam;                                                                            \
   1.726 +                                                                                                        \
   1.727 +                valarray<type>    tt(4);                                                                \
   1.728 +                                                                                                        \
   1.729 +                tt[0] = +a*tr[0]+b*tr[1]+c*tr[2]+d*tr[3];                                               \
   1.730 +                tt[1] = -a*tr[1]+b*tr[0]-c*tr[3]+d*tr[2];                                               \
   1.731 +                tt[2] = -a*tr[2]+b*tr[3]+c*tr[0]-d*tr[1];                                               \
   1.732 +                tt[3] = -a*tr[3]-b*tr[2]+c*tr[1]+d*tr[0];                                               \
   1.733 +                                                                                                        \
   1.734 +                tr *= tr;                                                                               \
   1.735 +                                                                                                        \
   1.736 +                tt *= (mixam/tr.sum());                                                                 \
   1.737 +                                                                                                        \
   1.738 +                a = tt[0];                                                                              \
   1.739 +                b = tt[1];                                                                              \
   1.740 +                c = tt[2];                                                                              \
   1.741 +                d = tt[3];                                                                              \
   1.742 +                                                                                                        \
   1.743 +                return(*this);                                                                          \
   1.744 +            }
   1.745 +#elif    defined(BOOST_NO_ARGUMENT_DEPENDENT_LOOKUP)
   1.746 +    #define    BOOST_QUATERNION_MEMBER_DIV_GENERATOR_3(type)                  \
   1.747 +            template<typename X>                                              \
   1.748 +            quaternion<type> &        operator /= (quaternion<X> const & rhs) \
   1.749 +            {                                                                 \
   1.750 +                using    ::std::valarray;                                     \
   1.751 +                using    ::std::abs;                                          \
   1.752 +                                                                              \
   1.753 +                valarray<type>    tr(4);                                      \
   1.754 +                                                                              \
   1.755 +                tr[0] = static_cast<type>(rhs.R_component_1());               \
   1.756 +                tr[1] = static_cast<type>(rhs.R_component_2());               \
   1.757 +                tr[2] = static_cast<type>(rhs.R_component_3());               \
   1.758 +                tr[3] = static_cast<type>(rhs.R_component_4());               \
   1.759 +                                                                              \
   1.760 +                type            mixam = static_cast<type>(1)/(abs(tr).max)(); \
   1.761 +                                                                              \
   1.762 +                tr *= mixam;                                                  \
   1.763 +                                                                              \
   1.764 +                valarray<type>    tt(4);                                      \
   1.765 +                                                                              \
   1.766 +                tt[0] = +a*tr[0]+b*tr[1]+c*tr[2]+d*tr[3];                     \
   1.767 +                tt[1] = -a*tr[1]+b*tr[0]-c*tr[3]+d*tr[2];                     \
   1.768 +                tt[2] = -a*tr[2]+b*tr[3]+c*tr[0]-d*tr[1];                     \
   1.769 +                tt[3] = -a*tr[3]-b*tr[2]+c*tr[1]+d*tr[0];                     \
   1.770 +                                                                              \
   1.771 +                tr *= tr;                                                     \
   1.772 +                                                                              \
   1.773 +                tt *= (mixam/tr.sum());                                       \
   1.774 +                                                                              \
   1.775 +                a = tt[0];                                                    \
   1.776 +                b = tt[1];                                                    \
   1.777 +                c = tt[2];                                                    \
   1.778 +                d = tt[3];                                                    \
   1.779 +                                                                              \
   1.780 +                return(*this);                                                \
   1.781 +            }
   1.782 +#else
   1.783 +    #define    BOOST_QUATERNION_MEMBER_DIV_GENERATOR_3(type)                  \
   1.784 +            template<typename X>                                              \
   1.785 +            quaternion<type> &        operator /= (quaternion<X> const & rhs) \
   1.786 +            {                                                                 \
   1.787 +                using    ::std::valarray;                                     \
   1.788 +                                                                              \
   1.789 +                valarray<type>    tr(4);                                      \
   1.790 +                                                                              \
   1.791 +                tr[0] = static_cast<type>(rhs.R_component_1());               \
   1.792 +                tr[1] = static_cast<type>(rhs.R_component_2());               \
   1.793 +                tr[2] = static_cast<type>(rhs.R_component_3());               \
   1.794 +                tr[3] = static_cast<type>(rhs.R_component_4());               \
   1.795 +                                                                              \
   1.796 +                type            mixam = static_cast<type>(1)/(abs(tr).max)(); \
   1.797 +                                                                              \
   1.798 +                tr *= mixam;                                                  \
   1.799 +                                                                              \
   1.800 +                valarray<type>    tt(4);                                      \
   1.801 +                                                                              \
   1.802 +                tt[0] = +a*tr[0]+b*tr[1]+c*tr[2]+d*tr[3];                     \
   1.803 +                tt[1] = -a*tr[1]+b*tr[0]-c*tr[3]+d*tr[2];                     \
   1.804 +                tt[2] = -a*tr[2]+b*tr[3]+c*tr[0]-d*tr[1];                     \
   1.805 +                tt[3] = -a*tr[3]-b*tr[2]+c*tr[1]+d*tr[0];                     \
   1.806 +                                                                              \
   1.807 +                tr *= tr;                                                     \
   1.808 +                                                                              \
   1.809 +                tt *= (mixam/tr.sum());                                       \
   1.810 +                                                                              \
   1.811 +                a = tt[0];                                                    \
   1.812 +                b = tt[1];                                                    \
   1.813 +                c = tt[2];                                                    \
   1.814 +                d = tt[3];                                                    \
   1.815 +                                                                              \
   1.816 +                return(*this);                                                \
   1.817 +            }
   1.818 +#endif /* defined(__GNUC__) && (__GNUC__ < 3) */ /* BOOST_NO_ARGUMENT_DEPENDENT_LOOKUP */
   1.819 +    
   1.820 +#define    BOOST_QUATERNION_MEMBER_ADD_GENERATOR(type)   \
   1.821 +        BOOST_QUATERNION_MEMBER_ADD_GENERATOR_1(type)    \
   1.822 +        BOOST_QUATERNION_MEMBER_ADD_GENERATOR_2(type)    \
   1.823 +        BOOST_QUATERNION_MEMBER_ADD_GENERATOR_3(type)
   1.824 +        
   1.825 +#define    BOOST_QUATERNION_MEMBER_SUB_GENERATOR(type)   \
   1.826 +        BOOST_QUATERNION_MEMBER_SUB_GENERATOR_1(type)    \
   1.827 +        BOOST_QUATERNION_MEMBER_SUB_GENERATOR_2(type)    \
   1.828 +        BOOST_QUATERNION_MEMBER_SUB_GENERATOR_3(type)
   1.829 +        
   1.830 +#define    BOOST_QUATERNION_MEMBER_MUL_GENERATOR(type)   \
   1.831 +        BOOST_QUATERNION_MEMBER_MUL_GENERATOR_1(type)    \
   1.832 +        BOOST_QUATERNION_MEMBER_MUL_GENERATOR_2(type)    \
   1.833 +        BOOST_QUATERNION_MEMBER_MUL_GENERATOR_3(type)
   1.834 +        
   1.835 +#define    BOOST_QUATERNION_MEMBER_DIV_GENERATOR(type)   \
   1.836 +        BOOST_QUATERNION_MEMBER_DIV_GENERATOR_1(type)    \
   1.837 +        BOOST_QUATERNION_MEMBER_DIV_GENERATOR_2(type)    \
   1.838 +        BOOST_QUATERNION_MEMBER_DIV_GENERATOR_3(type)
   1.839 +        
   1.840 +#define    BOOST_QUATERNION_MEMBER_ALGEBRAIC_GENERATOR(type)   \
   1.841 +        BOOST_QUATERNION_MEMBER_ADD_GENERATOR(type)            \
   1.842 +        BOOST_QUATERNION_MEMBER_SUB_GENERATOR(type)            \
   1.843 +        BOOST_QUATERNION_MEMBER_MUL_GENERATOR(type)            \
   1.844 +        BOOST_QUATERNION_MEMBER_DIV_GENERATOR(type)
   1.845 +        
   1.846 +        
   1.847 +        template<>
   1.848 +        class quaternion<float>
   1.849 +        {
   1.850 +        public:
   1.851 +            
   1.852 +            typedef float value_type;
   1.853 +            
   1.854 +            BOOST_QUATERNION_CONSTRUCTOR_GENERATOR(float)
   1.855 +            
   1.856 +            // UNtemplated copy constructor
   1.857 +            // (this is taken care of by the compiler itself)
   1.858 +            
   1.859 +            // explicit copy constructors (precision-loosing converters)
   1.860 +            
   1.861 +            explicit            quaternion(quaternion<double> const & a_recopier)
   1.862 +            {
   1.863 +                *this = detail::quaternion_type_converter<float, double>(a_recopier);
   1.864 +            }
   1.865 +            
   1.866 +            explicit            quaternion(quaternion<long double> const & a_recopier)
   1.867 +            {
   1.868 +                *this = detail::quaternion_type_converter<float, long double>(a_recopier);
   1.869 +            }
   1.870 +            
   1.871 +            // destructor
   1.872 +            // (this is taken care of by the compiler itself)
   1.873 +            
   1.874 +            // accessors
   1.875 +            //
   1.876 +            // Note:    Like complex number, quaternions do have a meaningful notion of "real part",
   1.877 +            //            but unlike them there is no meaningful notion of "imaginary part".
   1.878 +            //            Instead there is an "unreal part" which itself is a quaternion, and usually
   1.879 +            //            nothing simpler (as opposed to the complex number case).
   1.880 +            //            However, for practicallity, there are accessors for the other components
   1.881 +            //            (these are necessary for the templated copy constructor, for instance).
   1.882 +            
   1.883 +            BOOST_QUATERNION_ACCESSOR_GENERATOR(float)
   1.884 +            
   1.885 +            // assignment operators
   1.886 +            
   1.887 +            BOOST_QUATERNION_MEMBER_ASSIGNMENT_GENERATOR(float)
   1.888 +            
   1.889 +            // other assignment-related operators
   1.890 +            //
   1.891 +            // NOTE:    Quaternion multiplication is *NOT* commutative;
   1.892 +            //            symbolically, "q *= rhs;" means "q = q * rhs;"
   1.893 +            //            and "q /= rhs;" means "q = q * inverse_of(rhs);"
   1.894 +            
   1.895 +            BOOST_QUATERNION_MEMBER_ALGEBRAIC_GENERATOR(float)
   1.896 +            
   1.897 +            
   1.898 +        protected:
   1.899 +            
   1.900 +            BOOST_QUATERNION_MEMBER_DATA_GENERATOR(float)
   1.901 +            
   1.902 +            
   1.903 +        private:
   1.904 +            
   1.905 +        };
   1.906 +        
   1.907 +        
   1.908 +        template<>
   1.909 +        class quaternion<double>
   1.910 +        {
   1.911 +        public:
   1.912 +            
   1.913 +            typedef double value_type;
   1.914 +            
   1.915 +            BOOST_QUATERNION_CONSTRUCTOR_GENERATOR(double)
   1.916 +            
   1.917 +            // UNtemplated copy constructor
   1.918 +            // (this is taken care of by the compiler itself)
   1.919 +            
   1.920 +            // converting copy constructor
   1.921 +            
   1.922 +            explicit                quaternion(quaternion<float> const & a_recopier)
   1.923 +            {
   1.924 +                *this = detail::quaternion_type_converter<double, float>(a_recopier);
   1.925 +            }
   1.926 +            
   1.927 +            // explicit copy constructors (precision-loosing converters)
   1.928 +            
   1.929 +            explicit                quaternion(quaternion<long double> const & a_recopier)
   1.930 +            {
   1.931 +                *this = detail::quaternion_type_converter<double, long double>(a_recopier);
   1.932 +            }
   1.933 +            
   1.934 +            // destructor
   1.935 +            // (this is taken care of by the compiler itself)
   1.936 +            
   1.937 +            // accessors
   1.938 +            //
   1.939 +            // Note:    Like complex number, quaternions do have a meaningful notion of "real part",
   1.940 +            //            but unlike them there is no meaningful notion of "imaginary part".
   1.941 +            //            Instead there is an "unreal part" which itself is a quaternion, and usually
   1.942 +            //            nothing simpler (as opposed to the complex number case).
   1.943 +            //            However, for practicallity, there are accessors for the other components
   1.944 +            //            (these are necessary for the templated copy constructor, for instance).
   1.945 +            
   1.946 +            BOOST_QUATERNION_ACCESSOR_GENERATOR(double)
   1.947 +            
   1.948 +            // assignment operators
   1.949 +            
   1.950 +            BOOST_QUATERNION_MEMBER_ASSIGNMENT_GENERATOR(double)
   1.951 +            
   1.952 +            // other assignment-related operators
   1.953 +            //
   1.954 +            // NOTE:    Quaternion multiplication is *NOT* commutative;
   1.955 +            //            symbolically, "q *= rhs;" means "q = q * rhs;"
   1.956 +            //            and "q /= rhs;" means "q = q * inverse_of(rhs);"
   1.957 +            
   1.958 +            BOOST_QUATERNION_MEMBER_ALGEBRAIC_GENERATOR(double)
   1.959 +            
   1.960 +            
   1.961 +        protected:
   1.962 +            
   1.963 +            BOOST_QUATERNION_MEMBER_DATA_GENERATOR(double)
   1.964 +            
   1.965 +            
   1.966 +        private:
   1.967 +            
   1.968 +        };
   1.969 +        
   1.970 +        
   1.971 +        template<>
   1.972 +        class quaternion<long double>
   1.973 +        {
   1.974 +        public:
   1.975 +            
   1.976 +            typedef long double value_type;
   1.977 +            
   1.978 +            BOOST_QUATERNION_CONSTRUCTOR_GENERATOR(long double)
   1.979 +            
   1.980 +            // UNtemplated copy constructor
   1.981 +            // (this is taken care of by the compiler itself)
   1.982 +            
   1.983 +            // converting copy constructors
   1.984 +            
   1.985 +            explicit                    quaternion(quaternion<float> const & a_recopier)
   1.986 +            {
   1.987 +                *this = detail::quaternion_type_converter<long double, float>(a_recopier);
   1.988 +            }
   1.989 +            
   1.990 +            explicit                    quaternion(quaternion<double> const & a_recopier)
   1.991 +            {
   1.992 +                *this = detail::quaternion_type_converter<long double, double>(a_recopier);
   1.993 +            }
   1.994 +            
   1.995 +            // destructor
   1.996 +            // (this is taken care of by the compiler itself)
   1.997 +            
   1.998 +            // accessors
   1.999 +            //
  1.1000 +            // Note:    Like complex number, quaternions do have a meaningful notion of "real part",
  1.1001 +            //            but unlike them there is no meaningful notion of "imaginary part".
  1.1002 +            //            Instead there is an "unreal part" which itself is a quaternion, and usually
  1.1003 +            //            nothing simpler (as opposed to the complex number case).
  1.1004 +            //            However, for practicallity, there are accessors for the other components
  1.1005 +            //            (these are necessary for the templated copy constructor, for instance).
  1.1006 +            
  1.1007 +            BOOST_QUATERNION_ACCESSOR_GENERATOR(long double)
  1.1008 +            
  1.1009 +            // assignment operators
  1.1010 +            
  1.1011 +            BOOST_QUATERNION_MEMBER_ASSIGNMENT_GENERATOR(long double)
  1.1012 +            
  1.1013 +            // other assignment-related operators
  1.1014 +            //
  1.1015 +            // NOTE:    Quaternion multiplication is *NOT* commutative;
  1.1016 +            //            symbolically, "q *= rhs;" means "q = q * rhs;"
  1.1017 +            //            and "q /= rhs;" means "q = q * inverse_of(rhs);"
  1.1018 +            
  1.1019 +            BOOST_QUATERNION_MEMBER_ALGEBRAIC_GENERATOR(long double)
  1.1020 +            
  1.1021 +            
  1.1022 +        protected:
  1.1023 +            
  1.1024 +            BOOST_QUATERNION_MEMBER_DATA_GENERATOR(long double)
  1.1025 +            
  1.1026 +            
  1.1027 +        private:
  1.1028 +            
  1.1029 +        };
  1.1030 +        
  1.1031 +        
  1.1032 +#undef    BOOST_QUATERNION_MEMBER_ALGEBRAIC_GENERATOR
  1.1033 +#undef    BOOST_QUATERNION_MEMBER_ADD_GENERATOR
  1.1034 +#undef    BOOST_QUATERNION_MEMBER_SUB_GENERATOR
  1.1035 +#undef    BOOST_QUATERNION_MEMBER_MUL_GENERATOR
  1.1036 +#undef    BOOST_QUATERNION_MEMBER_DIV_GENERATOR
  1.1037 +#undef    BOOST_QUATERNION_MEMBER_ADD_GENERATOR_1
  1.1038 +#undef    BOOST_QUATERNION_MEMBER_ADD_GENERATOR_2
  1.1039 +#undef    BOOST_QUATERNION_MEMBER_ADD_GENERATOR_3
  1.1040 +#undef    BOOST_QUATERNION_MEMBER_SUB_GENERATOR_1
  1.1041 +#undef    BOOST_QUATERNION_MEMBER_SUB_GENERATOR_2
  1.1042 +#undef    BOOST_QUATERNION_MEMBER_SUB_GENERATOR_3
  1.1043 +#undef    BOOST_QUATERNION_MEMBER_MUL_GENERATOR_1
  1.1044 +#undef    BOOST_QUATERNION_MEMBER_MUL_GENERATOR_2
  1.1045 +#undef    BOOST_QUATERNION_MEMBER_MUL_GENERATOR_3
  1.1046 +#undef    BOOST_QUATERNION_MEMBER_DIV_GENERATOR_1
  1.1047 +#undef    BOOST_QUATERNION_MEMBER_DIV_GENERATOR_2
  1.1048 +#undef    BOOST_QUATERNION_MEMBER_DIV_GENERATOR_3
  1.1049 +        
  1.1050 +#undef    BOOST_QUATERNION_CONSTRUCTOR_GENERATOR
  1.1051 +        
  1.1052 +        
  1.1053 +#undef    BOOST_QUATERNION_MEMBER_ASSIGNMENT_GENERATOR
  1.1054 +        
  1.1055 +#undef    BOOST_QUATERNION_MEMBER_DATA_GENERATOR
  1.1056 +        
  1.1057 +#undef    BOOST_QUATERNION_ACCESSOR_GENERATOR
  1.1058 +        
  1.1059 +        
  1.1060 +        // operators
  1.1061 +        
  1.1062 +#define    BOOST_QUATERNION_OPERATOR_GENERATOR_BODY(op)      \
  1.1063 +        {                                                    \
  1.1064 +            quaternion<T>    res(lhs);                       \
  1.1065 +            res op##= rhs;                                   \
  1.1066 +            return(res);                                     \
  1.1067 +        }
  1.1068 +        
  1.1069 +#define    BOOST_QUATERNION_OPERATOR_GENERATOR_1_L(op)                                                  \
  1.1070 +        template<typename T>                                                                            \
  1.1071 +        inline quaternion<T>    operator op (T const & lhs, quaternion<T> const & rhs)                  \
  1.1072 +        BOOST_QUATERNION_OPERATOR_GENERATOR_BODY(op)
  1.1073 +        
  1.1074 +#define    BOOST_QUATERNION_OPERATOR_GENERATOR_1_R(op)                                                  \
  1.1075 +        template<typename T>                                                                            \
  1.1076 +        inline quaternion<T>    operator op (quaternion<T> const & lhs, T const & rhs)                  \
  1.1077 +        BOOST_QUATERNION_OPERATOR_GENERATOR_BODY(op)
  1.1078 +        
  1.1079 +#define    BOOST_QUATERNION_OPERATOR_GENERATOR_2_L(op)                                                  \
  1.1080 +        template<typename T>                                                                            \
  1.1081 +        inline quaternion<T>    operator op (::std::complex<T> const & lhs, quaternion<T> const & rhs)  \
  1.1082 +        BOOST_QUATERNION_OPERATOR_GENERATOR_BODY(op)
  1.1083 +        
  1.1084 +#define    BOOST_QUATERNION_OPERATOR_GENERATOR_2_R(op)                                                  \
  1.1085 +        template<typename T>                                                                            \
  1.1086 +        inline quaternion<T>    operator op (quaternion<T> const & lhs, ::std::complex<T> const & rhs)  \
  1.1087 +        BOOST_QUATERNION_OPERATOR_GENERATOR_BODY(op)
  1.1088 +        
  1.1089 +#define    BOOST_QUATERNION_OPERATOR_GENERATOR_3(op)                                                    \
  1.1090 +        template<typename T>                                                                            \
  1.1091 +        inline quaternion<T>    operator op (quaternion<T> const & lhs, quaternion<T> const & rhs)      \
  1.1092 +        BOOST_QUATERNION_OPERATOR_GENERATOR_BODY(op)
  1.1093 +        
  1.1094 +#define    BOOST_QUATERNION_OPERATOR_GENERATOR(op)     \
  1.1095 +        BOOST_QUATERNION_OPERATOR_GENERATOR_1_L(op)    \
  1.1096 +        BOOST_QUATERNION_OPERATOR_GENERATOR_1_R(op)    \
  1.1097 +        BOOST_QUATERNION_OPERATOR_GENERATOR_2_L(op)    \
  1.1098 +        BOOST_QUATERNION_OPERATOR_GENERATOR_2_R(op)    \
  1.1099 +        BOOST_QUATERNION_OPERATOR_GENERATOR_3(op)
  1.1100 +        
  1.1101 +        
  1.1102 +        BOOST_QUATERNION_OPERATOR_GENERATOR(+)
  1.1103 +        BOOST_QUATERNION_OPERATOR_GENERATOR(-)
  1.1104 +        BOOST_QUATERNION_OPERATOR_GENERATOR(*)
  1.1105 +        BOOST_QUATERNION_OPERATOR_GENERATOR(/)
  1.1106 +
  1.1107 +
  1.1108 +#undef    BOOST_QUATERNION_OPERATOR_GENERATOR
  1.1109 +        
  1.1110 +#undef    BOOST_QUATERNION_OPERATOR_GENERATOR_1_L
  1.1111 +#undef    BOOST_QUATERNION_OPERATOR_GENERATOR_1_R
  1.1112 +#undef    BOOST_QUATERNION_OPERATOR_GENERATOR_2_L
  1.1113 +#undef    BOOST_QUATERNION_OPERATOR_GENERATOR_2_R
  1.1114 +#undef    BOOST_QUATERNION_OPERATOR_GENERATOR_3
  1.1115 +
  1.1116 +#undef    BOOST_QUATERNION_OPERATOR_GENERATOR_BODY
  1.1117 +        
  1.1118 +        
  1.1119 +        template<typename T>
  1.1120 +        inline quaternion<T>                    operator + (quaternion<T> const & q)
  1.1121 +        {
  1.1122 +            return(q);
  1.1123 +        }
  1.1124 +        
  1.1125 +        
  1.1126 +        template<typename T>
  1.1127 +        inline quaternion<T>                    operator - (quaternion<T> const & q)
  1.1128 +        {
  1.1129 +            return(quaternion<T>(-q.R_component_1(),-q.R_component_2(),-q.R_component_3(),-q.R_component_4()));
  1.1130 +        }
  1.1131 +        
  1.1132 +        
  1.1133 +        template<typename T>
  1.1134 +        inline bool                                operator == (T const & lhs, quaternion<T> const & rhs)
  1.1135 +        {
  1.1136 +            return    (
  1.1137 +                        (rhs.R_component_1() == lhs)&&
  1.1138 +                        (rhs.R_component_2() == static_cast<T>(0))&&
  1.1139 +                        (rhs.R_component_3() == static_cast<T>(0))&&
  1.1140 +                        (rhs.R_component_4() == static_cast<T>(0))
  1.1141 +                    );
  1.1142 +        }
  1.1143 +        
  1.1144 +        
  1.1145 +        template<typename T>
  1.1146 +        inline bool                                operator == (quaternion<T> const & lhs, T const & rhs)
  1.1147 +        {
  1.1148 +            return    (
  1.1149 +                        (lhs.R_component_1() == rhs)&&
  1.1150 +                        (lhs.R_component_2() == static_cast<T>(0))&&
  1.1151 +                        (lhs.R_component_3() == static_cast<T>(0))&&
  1.1152 +                        (lhs.R_component_4() == static_cast<T>(0))
  1.1153 +                    );
  1.1154 +        }
  1.1155 +        
  1.1156 +        
  1.1157 +        template<typename T>
  1.1158 +        inline bool                                operator == (::std::complex<T> const & lhs, quaternion<T> const & rhs)
  1.1159 +        {
  1.1160 +            return    (
  1.1161 +                        (rhs.R_component_1() == lhs.real())&&
  1.1162 +                        (rhs.R_component_2() == lhs.imag())&&
  1.1163 +                        (rhs.R_component_3() == static_cast<T>(0))&&
  1.1164 +                        (rhs.R_component_4() == static_cast<T>(0))
  1.1165 +                    );
  1.1166 +        }
  1.1167 +        
  1.1168 +        
  1.1169 +        template<typename T>
  1.1170 +        inline bool                                operator == (quaternion<T> const & lhs, ::std::complex<T> const & rhs)
  1.1171 +        {
  1.1172 +            return    (
  1.1173 +                        (lhs.R_component_1() == rhs.real())&&
  1.1174 +                        (lhs.R_component_2() == rhs.imag())&&
  1.1175 +                        (lhs.R_component_3() == static_cast<T>(0))&&
  1.1176 +                        (lhs.R_component_4() == static_cast<T>(0))
  1.1177 +                    );
  1.1178 +        }
  1.1179 +        
  1.1180 +        
  1.1181 +        template<typename T>
  1.1182 +        inline bool                                operator == (quaternion<T> const & lhs, quaternion<T> const & rhs)
  1.1183 +        {
  1.1184 +            return    (
  1.1185 +                        (rhs.R_component_1() == lhs.R_component_1())&&
  1.1186 +                        (rhs.R_component_2() == lhs.R_component_2())&&
  1.1187 +                        (rhs.R_component_3() == lhs.R_component_3())&&
  1.1188 +                        (rhs.R_component_4() == lhs.R_component_4())
  1.1189 +                    );
  1.1190 +        }
  1.1191 +        
  1.1192 +        
  1.1193 +#define    BOOST_QUATERNION_NOT_EQUAL_GENERATOR  \
  1.1194 +        {                                        \
  1.1195 +            return(!(lhs == rhs));               \
  1.1196 +        }
  1.1197 +        
  1.1198 +        template<typename T>
  1.1199 +        inline bool                                operator != (T const & lhs, quaternion<T> const & rhs)
  1.1200 +        BOOST_QUATERNION_NOT_EQUAL_GENERATOR
  1.1201 +        
  1.1202 +        template<typename T>
  1.1203 +        inline bool                                operator != (quaternion<T> const & lhs, T const & rhs)
  1.1204 +        BOOST_QUATERNION_NOT_EQUAL_GENERATOR
  1.1205 +        
  1.1206 +        template<typename T>
  1.1207 +        inline bool                                operator != (::std::complex<T> const & lhs, quaternion<T> const & rhs)
  1.1208 +        BOOST_QUATERNION_NOT_EQUAL_GENERATOR
  1.1209 +        
  1.1210 +        template<typename T>
  1.1211 +        inline bool                                operator != (quaternion<T> const & lhs, ::std::complex<T> const & rhs)
  1.1212 +        BOOST_QUATERNION_NOT_EQUAL_GENERATOR
  1.1213 +        
  1.1214 +        template<typename T>
  1.1215 +        inline bool                                operator != (quaternion<T> const & lhs, quaternion<T> const & rhs)
  1.1216 +        BOOST_QUATERNION_NOT_EQUAL_GENERATOR
  1.1217 +        
  1.1218 +#undef    BOOST_QUATERNION_NOT_EQUAL_GENERATOR
  1.1219 +        
  1.1220 +        
  1.1221 +        // Note:    we allow the following formats, whith a, b, c, and d reals
  1.1222 +        //            a
  1.1223 +        //            (a), (a,b), (a,b,c), (a,b,c,d)
  1.1224 +        //            (a,(c)), (a,(c,d)), ((a)), ((a),c), ((a),(c)), ((a),(c,d)), ((a,b)), ((a,b),c), ((a,b),(c)), ((a,b),(c,d))
  1.1225 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1226 +        template<typename T>
  1.1227 +        std::istream &                            operator >> (    ::std::istream & is,
  1.1228 +                                                                quaternion<T> & q)
  1.1229 +#else
  1.1230 +        template<typename T, typename charT, class traits>
  1.1231 +        ::std::basic_istream<charT,traits> &    operator >> (    ::std::basic_istream<charT,traits> & is,
  1.1232 +                                                                quaternion<T> & q)
  1.1233 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1234 +        {
  1.1235 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1236 +            typedef char    charT;
  1.1237 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1238 +            
  1.1239 +#ifdef    BOOST_NO_STD_LOCALE
  1.1240 +#else
  1.1241 +            const ::std::ctype<charT> & ct = ::std::use_facet< ::std::ctype<charT> >(is.getloc());
  1.1242 +#endif /* BOOST_NO_STD_LOCALE */
  1.1243 +            
  1.1244 +            T    a = T();
  1.1245 +            T    b = T();
  1.1246 +            T    c = T();
  1.1247 +            T    d = T();
  1.1248 +            
  1.1249 +            ::std::complex<T>    u = ::std::complex<T>();
  1.1250 +            ::std::complex<T>    v = ::std::complex<T>();
  1.1251 +            
  1.1252 +            charT    ch = charT();
  1.1253 +            char    cc;
  1.1254 +            
  1.1255 +            is >> ch;                                        // get the first lexeme
  1.1256 +            
  1.1257 +            if    (!is.good())    goto finish;
  1.1258 +            
  1.1259 +#ifdef    BOOST_NO_STD_LOCALE
  1.1260 +            cc = ch;
  1.1261 +#else
  1.1262 +            cc = ct.narrow(ch, char());
  1.1263 +#endif /* BOOST_NO_STD_LOCALE */
  1.1264 +            
  1.1265 +            if    (cc == '(')                            // read "(", possible: (a), (a,b), (a,b,c), (a,b,c,d), (a,(c)), (a,(c,d)), ((a)), ((a),c), ((a),(c)), ((a),(c,d)), ((a,b)), ((a,b),c), ((a,b),(c)), ((a,b,),(c,d,))
  1.1266 +            {
  1.1267 +                is >> ch;                                    // get the second lexeme
  1.1268 +                
  1.1269 +                if    (!is.good())    goto finish;
  1.1270 +                
  1.1271 +#ifdef    BOOST_NO_STD_LOCALE
  1.1272 +                cc = ch;
  1.1273 +#else
  1.1274 +                cc = ct.narrow(ch, char());
  1.1275 +#endif /* BOOST_NO_STD_LOCALE */
  1.1276 +                
  1.1277 +                if    (cc == '(')                        // read "((", possible: ((a)), ((a),c), ((a),(c)), ((a),(c,d)), ((a,b)), ((a,b),c), ((a,b),(c)), ((a,b,),(c,d,))
  1.1278 +                {
  1.1279 +                    is.putback(ch);
  1.1280 +                    
  1.1281 +                    is >> u;                                // we extract the first and second components
  1.1282 +                    a = u.real();
  1.1283 +                    b = u.imag();
  1.1284 +                    
  1.1285 +                    if    (!is.good())    goto finish;
  1.1286 +                    
  1.1287 +                    is >> ch;                                // get the next lexeme
  1.1288 +                    
  1.1289 +                    if    (!is.good())    goto finish;
  1.1290 +                    
  1.1291 +#ifdef    BOOST_NO_STD_LOCALE
  1.1292 +                    cc = ch;
  1.1293 +#else
  1.1294 +                    cc = ct.narrow(ch, char());
  1.1295 +#endif /* BOOST_NO_STD_LOCALE */
  1.1296 +                    
  1.1297 +                    if        (cc == ')')                    // format: ((a)) or ((a,b))
  1.1298 +                    {
  1.1299 +                        q = quaternion<T>(a,b);
  1.1300 +                    }
  1.1301 +                    else if    (cc == ',')                // read "((a)," or "((a,b),", possible: ((a),c), ((a),(c)), ((a),(c,d)), ((a,b),c), ((a,b),(c)), ((a,b,),(c,d,))
  1.1302 +                    {
  1.1303 +                        is >> v;                            // we extract the third and fourth components
  1.1304 +                        c = v.real();
  1.1305 +                        d = v.imag();
  1.1306 +                        
  1.1307 +                        if    (!is.good())    goto finish;
  1.1308 +                        
  1.1309 +                        is >> ch;                                // get the last lexeme
  1.1310 +                        
  1.1311 +                        if    (!is.good())    goto finish;
  1.1312 +                        
  1.1313 +#ifdef    BOOST_NO_STD_LOCALE
  1.1314 +                        cc = ch;
  1.1315 +#else
  1.1316 +                        cc = ct.narrow(ch, char());
  1.1317 +#endif /* BOOST_NO_STD_LOCALE */
  1.1318 +                        
  1.1319 +                        if    (cc == ')')                    // format: ((a),c), ((a),(c)), ((a),(c,d)), ((a,b),c), ((a,b),(c)) or ((a,b,),(c,d,))
  1.1320 +                        {
  1.1321 +                            q = quaternion<T>(a,b,c,d);
  1.1322 +                        }
  1.1323 +                        else                            // error
  1.1324 +                        {
  1.1325 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1326 +                            is.setstate(::std::ios::failbit);
  1.1327 +#else
  1.1328 +                            is.setstate(::std::ios_base::failbit);
  1.1329 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1330 +                        }
  1.1331 +                    }
  1.1332 +                    else                                // error
  1.1333 +                    {
  1.1334 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1335 +                        is.setstate(::std::ios::failbit);
  1.1336 +#else
  1.1337 +                        is.setstate(::std::ios_base::failbit);
  1.1338 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1339 +                    }
  1.1340 +                }
  1.1341 +                else                                // read "(a", possible: (a), (a,b), (a,b,c), (a,b,c,d), (a,(c)), (a,(c,d))
  1.1342 +                {
  1.1343 +                    is.putback(ch);
  1.1344 +                    
  1.1345 +                    is >> a;                                // we extract the first component
  1.1346 +                    
  1.1347 +                    if    (!is.good())    goto finish;
  1.1348 +                    
  1.1349 +                    is >> ch;                                // get the third lexeme
  1.1350 +                    
  1.1351 +                    if    (!is.good())    goto finish;
  1.1352 +                    
  1.1353 +#ifdef    BOOST_NO_STD_LOCALE
  1.1354 +                    cc = ch;
  1.1355 +#else
  1.1356 +                    cc = ct.narrow(ch, char());
  1.1357 +#endif /* BOOST_NO_STD_LOCALE */
  1.1358 +                    
  1.1359 +                    if        (cc == ')')                    // format: (a)
  1.1360 +                    {
  1.1361 +                        q = quaternion<T>(a);
  1.1362 +                    }
  1.1363 +                    else if    (cc == ',')                // read "(a,", possible: (a,b), (a,b,c), (a,b,c,d), (a,(c)), (a,(c,d))
  1.1364 +                    {
  1.1365 +                        is >> ch;                            // get the fourth lexeme
  1.1366 +                        
  1.1367 +                        if    (!is.good())    goto finish;
  1.1368 +                        
  1.1369 +#ifdef    BOOST_NO_STD_LOCALE
  1.1370 +                        cc = ch;
  1.1371 +#else
  1.1372 +                        cc = ct.narrow(ch, char());
  1.1373 +#endif /* BOOST_NO_STD_LOCALE */
  1.1374 +                        
  1.1375 +                        if    (cc == '(')                // read "(a,(", possible: (a,(c)), (a,(c,d))
  1.1376 +                        {
  1.1377 +                            is.putback(ch);
  1.1378 +                            
  1.1379 +                            is >> v;                        // we extract the third and fourth component
  1.1380 +                            
  1.1381 +                            c = v.real();
  1.1382 +                            d = v.imag();
  1.1383 +                            
  1.1384 +                            if    (!is.good())    goto finish;
  1.1385 +                            
  1.1386 +                            is >> ch;                        // get the ninth lexeme
  1.1387 +                            
  1.1388 +                            if    (!is.good())    goto finish;
  1.1389 +                            
  1.1390 +#ifdef    BOOST_NO_STD_LOCALE
  1.1391 +                            cc = ch;
  1.1392 +#else
  1.1393 +                            cc = ct.narrow(ch, char());
  1.1394 +#endif /* BOOST_NO_STD_LOCALE */
  1.1395 +                            
  1.1396 +                            if    (cc == ')')                // format: (a,(c)) or (a,(c,d))
  1.1397 +                            {
  1.1398 +                                q = quaternion<T>(a,b,c,d);
  1.1399 +                            }
  1.1400 +                            else                        // error
  1.1401 +                            {
  1.1402 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1403 +                                is.setstate(::std::ios::failbit);
  1.1404 +#else
  1.1405 +                                is.setstate(::std::ios_base::failbit);
  1.1406 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1407 +                            }
  1.1408 +                        }
  1.1409 +                        else                        // read "(a,b", possible: (a,b), (a,b,c), (a,b,c,d)
  1.1410 +                        {
  1.1411 +                            is.putback(ch);
  1.1412 +                            
  1.1413 +                            is >> b;                        // we extract the second component
  1.1414 +                            
  1.1415 +                            if    (!is.good())    goto finish;
  1.1416 +                            
  1.1417 +                            is >> ch;                        // get the fifth lexeme
  1.1418 +                            
  1.1419 +                            if    (!is.good())    goto finish;
  1.1420 +                            
  1.1421 +#ifdef    BOOST_NO_STD_LOCALE
  1.1422 +                            cc = ch;
  1.1423 +#else
  1.1424 +                            cc = ct.narrow(ch, char());
  1.1425 +#endif /* BOOST_NO_STD_LOCALE */
  1.1426 +                            
  1.1427 +                            if    (cc == ')')                // format: (a,b)
  1.1428 +                            {
  1.1429 +                                q = quaternion<T>(a,b);
  1.1430 +                            }
  1.1431 +                            else if    (cc == ',')        // read "(a,b,", possible: (a,b,c), (a,b,c,d)
  1.1432 +                            {
  1.1433 +                                is >> c;                    // we extract the third component
  1.1434 +                                
  1.1435 +                                if    (!is.good())    goto finish;
  1.1436 +                                
  1.1437 +                                is >> ch;                    // get the seventh lexeme
  1.1438 +                                
  1.1439 +                                if    (!is.good())    goto finish;
  1.1440 +                                
  1.1441 +#ifdef    BOOST_NO_STD_LOCALE
  1.1442 +                                cc = ch;
  1.1443 +#else
  1.1444 +                                cc = ct.narrow(ch, char());
  1.1445 +#endif /* BOOST_NO_STD_LOCALE */
  1.1446 +                                
  1.1447 +                                if        (cc == ')')        // format: (a,b,c)
  1.1448 +                                {
  1.1449 +                                    q = quaternion<T>(a,b,c);
  1.1450 +                                }
  1.1451 +                                else if    (cc == ',')    // read "(a,b,c,", possible: (a,b,c,d)
  1.1452 +                                {
  1.1453 +                                    is >> d;                // we extract the fourth component
  1.1454 +                                    
  1.1455 +                                    if    (!is.good())    goto finish;
  1.1456 +                                    
  1.1457 +                                    is >> ch;                // get the ninth lexeme
  1.1458 +                                    
  1.1459 +                                    if    (!is.good())    goto finish;
  1.1460 +                                    
  1.1461 +#ifdef    BOOST_NO_STD_LOCALE
  1.1462 +                                    cc = ch;
  1.1463 +#else
  1.1464 +                                    cc = ct.narrow(ch, char());
  1.1465 +#endif /* BOOST_NO_STD_LOCALE */
  1.1466 +                                    
  1.1467 +                                    if    (cc == ')')        // format: (a,b,c,d)
  1.1468 +                                    {
  1.1469 +                                        q = quaternion<T>(a,b,c,d);
  1.1470 +                                    }
  1.1471 +                                    else                // error
  1.1472 +                                    {
  1.1473 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1474 +                                        is.setstate(::std::ios::failbit);
  1.1475 +#else
  1.1476 +                                        is.setstate(::std::ios_base::failbit);
  1.1477 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1478 +                                    }
  1.1479 +                                }
  1.1480 +                                else                    // error
  1.1481 +                                {
  1.1482 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1483 +                                    is.setstate(::std::ios::failbit);
  1.1484 +#else
  1.1485 +                                    is.setstate(::std::ios_base::failbit);
  1.1486 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1487 +                                }
  1.1488 +                            }
  1.1489 +                            else                        // error
  1.1490 +                            {
  1.1491 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1492 +                                is.setstate(::std::ios::failbit);
  1.1493 +#else
  1.1494 +                                is.setstate(::std::ios_base::failbit);
  1.1495 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1496 +                            }
  1.1497 +                        }
  1.1498 +                    }
  1.1499 +                    else                                // error
  1.1500 +                    {
  1.1501 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1502 +                        is.setstate(::std::ios::failbit);
  1.1503 +#else
  1.1504 +                        is.setstate(::std::ios_base::failbit);
  1.1505 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1506 +                    }
  1.1507 +                }
  1.1508 +            }
  1.1509 +            else                                        // format:    a
  1.1510 +            {
  1.1511 +                is.putback(ch);
  1.1512 +                
  1.1513 +                is >> a;                                    // we extract the first component
  1.1514 +                
  1.1515 +                if    (!is.good())    goto finish;
  1.1516 +                
  1.1517 +                q = quaternion<T>(a);
  1.1518 +            }
  1.1519 +            
  1.1520 +            finish:
  1.1521 +            return(is);
  1.1522 +        }
  1.1523 +        
  1.1524 +        
  1.1525 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1526 +        template<typename T>
  1.1527 +        ::std::ostream &                         operator << (    ::std::ostream & os,
  1.1528 +                                                                quaternion<T> const & q)
  1.1529 +#else
  1.1530 +        template<typename T, typename charT, class traits>
  1.1531 +        ::std::basic_ostream<charT,traits> &    operator << (    ::std::basic_ostream<charT,traits> & os,
  1.1532 +                                                                quaternion<T> const & q)
  1.1533 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1534 +        {
  1.1535 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1536 +            ::std::ostringstream                        s;
  1.1537 +#else
  1.1538 +            ::std::basic_ostringstream<charT,traits>    s;
  1.1539 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1540 +            
  1.1541 +            s.flags(os.flags());
  1.1542 +#ifdef    BOOST_NO_STD_LOCALE
  1.1543 +#else
  1.1544 +            s.imbue(os.getloc());
  1.1545 +#endif /* BOOST_NO_STD_LOCALE */
  1.1546 +            s.precision(os.precision());
  1.1547 +            
  1.1548 +            s << '('    << q.R_component_1() << ','
  1.1549 +                        << q.R_component_2() << ','
  1.1550 +                        << q.R_component_3() << ','
  1.1551 +                        << q.R_component_4() << ')';
  1.1552 +            
  1.1553 +            return os << s.str();
  1.1554 +        }
  1.1555 +        
  1.1556 +        
  1.1557 +        // values
  1.1558 +        
  1.1559 +        template<typename T>
  1.1560 +        inline T                                real(quaternion<T> const & q)
  1.1561 +        {
  1.1562 +            return(q.real());
  1.1563 +        }
  1.1564 +        
  1.1565 +        
  1.1566 +        template<typename T>
  1.1567 +        inline quaternion<T>                    unreal(quaternion<T> const & q)
  1.1568 +        {
  1.1569 +            return(q.unreal());
  1.1570 +        }
  1.1571 +        
  1.1572 +        
  1.1573 +#define    BOOST_QUATERNION_VALARRAY_LOADER  \
  1.1574 +            using    ::std::valarray;        \
  1.1575 +                                             \
  1.1576 +            valarray<T>    temp(4);          \
  1.1577 +                                             \
  1.1578 +            temp[0] = q.R_component_1();     \
  1.1579 +            temp[1] = q.R_component_2();     \
  1.1580 +            temp[2] = q.R_component_3();     \
  1.1581 +            temp[3] = q.R_component_4();
  1.1582 +        
  1.1583 +        
  1.1584 +        template<typename T>
  1.1585 +        inline T                                sup(quaternion<T> const & q)
  1.1586 +        {
  1.1587 +#ifdef    BOOST_NO_ARGUMENT_DEPENDENT_LOOKUP
  1.1588 +            using    ::std::abs;
  1.1589 +#endif    /* BOOST_NO_ARGUMENT_DEPENDENT_LOOKUP */
  1.1590 +            
  1.1591 +            BOOST_QUATERNION_VALARRAY_LOADER
  1.1592 +            
  1.1593 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1594 +            return((BOOST_GET_VALARRAY(T, abs(temp)).max)());
  1.1595 +#else
  1.1596 +            return((abs(temp).max)());
  1.1597 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1598 +        }
  1.1599 +        
  1.1600 +        
  1.1601 +        template<typename T>
  1.1602 +        inline T                                l1(quaternion<T> const & q)
  1.1603 +        {
  1.1604 +#ifdef    BOOST_NO_ARGUMENT_DEPENDENT_LOOKUP
  1.1605 +            using    ::std::abs;
  1.1606 +#endif    /* BOOST_NO_ARGUMENT_DEPENDENT_LOOKUP */
  1.1607 +            
  1.1608 +            BOOST_QUATERNION_VALARRAY_LOADER
  1.1609 +            
  1.1610 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1611 +            return(BOOST_GET_VALARRAY(T, abs(temp)).sum());
  1.1612 +#else
  1.1613 +            return(abs(temp).sum());
  1.1614 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1615 +        }
  1.1616 +        
  1.1617 +        
  1.1618 +        template<typename T>
  1.1619 +        inline T                                abs(quaternion<T> const & q)
  1.1620 +        {
  1.1621 +#ifdef    BOOST_NO_ARGUMENT_DEPENDENT_LOOKUP
  1.1622 +            using    ::std::abs;
  1.1623 +#endif    /* BOOST_NO_ARGUMENT_DEPENDENT_LOOKUP */
  1.1624 +            
  1.1625 +            using    ::std::sqrt;
  1.1626 +            
  1.1627 +            BOOST_QUATERNION_VALARRAY_LOADER
  1.1628 +            
  1.1629 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1630 +            T            maxim = (BOOST_GET_VALARRAY(T, abs(temp)).max)();    // overflow protection
  1.1631 +#else
  1.1632 +            T            maxim = (abs(temp).max)();    // overflow protection
  1.1633 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1634 +            
  1.1635 +            if    (maxim == static_cast<T>(0))
  1.1636 +            {
  1.1637 +                return(maxim);
  1.1638 +            }
  1.1639 +            else
  1.1640 +            {
  1.1641 +                T    mixam = static_cast<T>(1)/maxim;    // prefer multiplications over divisions
  1.1642 +                
  1.1643 +                temp *= mixam;
  1.1644 +                
  1.1645 +                temp *= temp;
  1.1646 +                
  1.1647 +                return(maxim*sqrt(temp.sum()));
  1.1648 +            }
  1.1649 +            
  1.1650 +            //return(sqrt(norm(q)));
  1.1651 +        }
  1.1652 +        
  1.1653 +        
  1.1654 +#undef    BOOST_QUATERNION_VALARRAY_LOADER
  1.1655 +        
  1.1656 +        
  1.1657 +        // Note:    This is the Cayley norm, not the Euclidian norm...
  1.1658 +        
  1.1659 +        template<typename T>
  1.1660 +        inline T                                norm(quaternion<T>const  & q)
  1.1661 +        {
  1.1662 +            return(real(q*conj(q)));
  1.1663 +        }
  1.1664 +        
  1.1665 +        
  1.1666 +        template<typename T>
  1.1667 +        inline quaternion<T>                    conj(quaternion<T> const & q)
  1.1668 +        {
  1.1669 +            return(quaternion<T>(   +q.R_component_1(),
  1.1670 +                                    -q.R_component_2(),
  1.1671 +                                    -q.R_component_3(),
  1.1672 +                                    -q.R_component_4()));
  1.1673 +        }
  1.1674 +        
  1.1675 +        
  1.1676 +        template<typename T>
  1.1677 +        inline quaternion<T>                    spherical(  T const & rho,
  1.1678 +                                                            T const & theta,
  1.1679 +                                                            T const & phi1,
  1.1680 +                                                            T const & phi2)
  1.1681 +        {
  1.1682 +            using ::std::cos;
  1.1683 +            using ::std::sin;
  1.1684 +            
  1.1685 +            //T    a = cos(theta)*cos(phi1)*cos(phi2);
  1.1686 +            //T    b = sin(theta)*cos(phi1)*cos(phi2);
  1.1687 +            //T    c = sin(phi1)*cos(phi2);
  1.1688 +            //T    d = sin(phi2);
  1.1689 +            
  1.1690 +            T    courrant = static_cast<T>(1);
  1.1691 +            
  1.1692 +            T    d = sin(phi2);
  1.1693 +            
  1.1694 +            courrant *= cos(phi2);
  1.1695 +            
  1.1696 +            T    c = sin(phi1)*courrant;
  1.1697 +            
  1.1698 +            courrant *= cos(phi1);
  1.1699 +            
  1.1700 +            T    b = sin(theta)*courrant;
  1.1701 +            T    a = cos(theta)*courrant;
  1.1702 +            
  1.1703 +            return(rho*quaternion<T>(a,b,c,d));
  1.1704 +        }
  1.1705 +        
  1.1706 +        
  1.1707 +        template<typename T>
  1.1708 +        inline quaternion<T>                    semipolar(  T const & rho,
  1.1709 +                                                            T const & alpha,
  1.1710 +                                                            T const & theta1,
  1.1711 +                                                            T const & theta2)
  1.1712 +        {
  1.1713 +            using ::std::cos;
  1.1714 +            using ::std::sin;
  1.1715 +            
  1.1716 +            T    a = cos(alpha)*cos(theta1);
  1.1717 +            T    b = cos(alpha)*sin(theta1);
  1.1718 +            T    c = sin(alpha)*cos(theta2);
  1.1719 +            T    d = sin(alpha)*sin(theta2);
  1.1720 +            
  1.1721 +            return(rho*quaternion<T>(a,b,c,d));
  1.1722 +        }
  1.1723 +        
  1.1724 +        
  1.1725 +        template<typename T>
  1.1726 +        inline quaternion<T>                    multipolar( T const & rho1,
  1.1727 +                                                            T const & theta1,
  1.1728 +                                                            T const & rho2,
  1.1729 +                                                            T const & theta2)
  1.1730 +        {
  1.1731 +            using ::std::cos;
  1.1732 +            using ::std::sin;
  1.1733 +            
  1.1734 +            T    a = rho1*cos(theta1);
  1.1735 +            T    b = rho1*sin(theta1);
  1.1736 +            T    c = rho2*cos(theta2);
  1.1737 +            T    d = rho2*sin(theta2);
  1.1738 +            
  1.1739 +            return(quaternion<T>(a,b,c,d));
  1.1740 +        }
  1.1741 +        
  1.1742 +        
  1.1743 +        template<typename T>
  1.1744 +        inline quaternion<T>                    cylindrospherical(  T const & t,
  1.1745 +                                                                    T const & radius,
  1.1746 +                                                                    T const & longitude,
  1.1747 +                                                                    T const & latitude)
  1.1748 +        {
  1.1749 +            using ::std::cos;
  1.1750 +            using ::std::sin;
  1.1751 +            
  1.1752 +            
  1.1753 +            
  1.1754 +            T    b = radius*cos(longitude)*cos(latitude);
  1.1755 +            T    c = radius*sin(longitude)*cos(latitude);
  1.1756 +            T    d = radius*sin(latitude);
  1.1757 +            
  1.1758 +            return(quaternion<T>(t,b,c,d));
  1.1759 +        }
  1.1760 +        
  1.1761 +        
  1.1762 +        template<typename T>
  1.1763 +        inline quaternion<T>                    cylindrical(T const & r,
  1.1764 +                                                            T const & angle,
  1.1765 +                                                            T const & h1,
  1.1766 +                                                            T const & h2)
  1.1767 +        {
  1.1768 +            using ::std::cos;
  1.1769 +            using ::std::sin;
  1.1770 +            
  1.1771 +            T    a = r*cos(angle);
  1.1772 +            T    b = r*sin(angle);
  1.1773 +            
  1.1774 +            return(quaternion<T>(a,b,h1,h2));
  1.1775 +        }
  1.1776 +        
  1.1777 +        
  1.1778 +        // transcendentals
  1.1779 +        // (please see the documentation)
  1.1780 +        
  1.1781 +        
  1.1782 +        template<typename T>
  1.1783 +        inline quaternion<T>                    exp(quaternion<T> const & q)
  1.1784 +        {
  1.1785 +            using    ::std::exp;
  1.1786 +            using    ::std::cos;
  1.1787 +            
  1.1788 +            using    ::boost::math::sinc_pi;
  1.1789 +            
  1.1790 +            T    u = exp(real(q));
  1.1791 +            
  1.1792 +            T    z = abs(unreal(q));
  1.1793 +            
  1.1794 +            T    w = sinc_pi(z);
  1.1795 +            
  1.1796 +            return(u*quaternion<T>(cos(z),
  1.1797 +                w*q.R_component_2(), w*q.R_component_3(),
  1.1798 +                w*q.R_component_4()));
  1.1799 +        }
  1.1800 +        
  1.1801 +        
  1.1802 +        template<typename T>
  1.1803 +        inline quaternion<T>                    cos(quaternion<T> const & q)
  1.1804 +        {
  1.1805 +            using    ::std::sin;
  1.1806 +            using    ::std::cos;
  1.1807 +            using    ::std::cosh;
  1.1808 +            
  1.1809 +            using    ::boost::math::sinhc_pi;
  1.1810 +            
  1.1811 +            T    z = abs(unreal(q));
  1.1812 +            
  1.1813 +            T    w = -sin(q.real())*sinhc_pi(z);
  1.1814 +            
  1.1815 +            return(quaternion<T>(cos(q.real())*cosh(z),
  1.1816 +                w*q.R_component_2(), w*q.R_component_3(),
  1.1817 +                w*q.R_component_4()));
  1.1818 +        }
  1.1819 +        
  1.1820 +        
  1.1821 +        template<typename T>
  1.1822 +        inline quaternion<T>                    sin(quaternion<T> const & q)
  1.1823 +        {
  1.1824 +            using    ::std::sin;
  1.1825 +            using    ::std::cos;
  1.1826 +            using    ::std::cosh;
  1.1827 +            
  1.1828 +            using    ::boost::math::sinhc_pi;
  1.1829 +            
  1.1830 +            T    z = abs(unreal(q));
  1.1831 +            
  1.1832 +            T    w = +cos(q.real())*sinhc_pi(z);
  1.1833 +            
  1.1834 +            return(quaternion<T>(sin(q.real())*cosh(z),
  1.1835 +                w*q.R_component_2(), w*q.R_component_3(),
  1.1836 +                w*q.R_component_4()));
  1.1837 +        }
  1.1838 +        
  1.1839 +        
  1.1840 +        template<typename T>
  1.1841 +        inline quaternion<T>                    tan(quaternion<T> const & q)
  1.1842 +        {
  1.1843 +            return(sin(q)/cos(q));
  1.1844 +        }
  1.1845 +        
  1.1846 +        
  1.1847 +        template<typename T>
  1.1848 +        inline quaternion<T>                    cosh(quaternion<T> const & q)
  1.1849 +        {
  1.1850 +            return((exp(+q)+exp(-q))/static_cast<T>(2));
  1.1851 +        }
  1.1852 +        
  1.1853 +        
  1.1854 +        template<typename T>
  1.1855 +        inline quaternion<T>                    sinh(quaternion<T> const & q)
  1.1856 +        {
  1.1857 +            return((exp(+q)-exp(-q))/static_cast<T>(2));
  1.1858 +        }
  1.1859 +        
  1.1860 +        
  1.1861 +        template<typename T>
  1.1862 +        inline quaternion<T>                    tanh(quaternion<T> const & q)
  1.1863 +        {
  1.1864 +            return(sinh(q)/cosh(q));
  1.1865 +        }
  1.1866 +        
  1.1867 +        
  1.1868 +        template<typename T>
  1.1869 +        quaternion<T>                            pow(quaternion<T> const & q,
  1.1870 +                                                    int n)
  1.1871 +        {
  1.1872 +            if        (n > 1)
  1.1873 +            {
  1.1874 +                int    m = n>>1;
  1.1875 +                
  1.1876 +                quaternion<T>    result = pow(q, m);
  1.1877 +                
  1.1878 +                result *= result;
  1.1879 +                
  1.1880 +                if    (n != (m<<1))
  1.1881 +                {
  1.1882 +                    result *= q; // n odd
  1.1883 +                }
  1.1884 +                
  1.1885 +                return(result);
  1.1886 +            }
  1.1887 +            else if    (n == 1)
  1.1888 +            {
  1.1889 +                return(q);
  1.1890 +            }
  1.1891 +            else if    (n == 0)
  1.1892 +            {
  1.1893 +                return(quaternion<T>(1));
  1.1894 +            }
  1.1895 +            else    /* n < 0 */
  1.1896 +            {
  1.1897 +                return(pow(quaternion<T>(1)/q,-n));
  1.1898 +            }
  1.1899 +        }
  1.1900 +        
  1.1901 +        
  1.1902 +        // helper templates for converting copy constructors (definition)
  1.1903 +        
  1.1904 +        namespace detail
  1.1905 +        {
  1.1906 +            
  1.1907 +            template<   typename T,
  1.1908 +                        typename U
  1.1909 +                    >
  1.1910 +            quaternion<T>    quaternion_type_converter(quaternion<U> const & rhs)
  1.1911 +            {
  1.1912 +                return(quaternion<T>(   static_cast<T>(rhs.R_component_1()),
  1.1913 +                                        static_cast<T>(rhs.R_component_2()),
  1.1914 +                                        static_cast<T>(rhs.R_component_3()),
  1.1915 +                                        static_cast<T>(rhs.R_component_4())));
  1.1916 +            }
  1.1917 +        }
  1.1918 +    }
  1.1919 +}
  1.1920 +
  1.1921 +
  1.1922 +#if    BOOST_WORKAROUND(__GNUC__, < 3)
  1.1923 +    #undef    BOOST_GET_VALARRAY
  1.1924 +#endif /* BOOST_WORKAROUND(__GNUC__, < 3) */
  1.1925 +
  1.1926 +
  1.1927 +#endif /* BOOST_QUATERNION_HPP */