
Boost : 
From: lums_at_[hidden]
Date: 20010525 12:53:52
Just a quick reply to the overflow question. One way to handle this
is to find the largest of ar, br, cr, dr, scale them all by that
value, compute the scaled denominator, compute the scaled at, bt, ct,
dt, and then rescale the results.
If I count correctly, there are three comparisons, four divides, and
one more multiply introduced.
Actually, to do this a little better, one would want to compute
inverse quantities and do multiplies rather than divides throughout.
Then there would be three comparisons, one divide, and six more
multiplies introduced to deal with the overflow issue. The original
method should be modified to include one divide and four multiplies
rather than four divides. This is all assuming the basic type of the
quaternion is something the compiler/cpu would deal with efficiently.
I would produce an example, but my interface to boost is via yahoo
groups and it is quite a pain to put code in. I hope the explanation
conveys the right idea.
 In boost_at_y..., Peter Schmitteckert (boost) <boost_at_s...> wrote:
> Salut,
>
> first of all I would see the library in boost, but I suggest that
problems
> with accuracy/overflows should be revised, or at least be
mentioned.
>
> Special Functions:
> ==============
>
> There are calculations of the following form
>
> if (abs(x) <= numeric_limits<T>::epsilon())
> {
> return(static_cast<T>(1));
> }
> else
> {
> return(sin(x)/x);
> }
>
> I would like to suggest to use something similar to
>
> const T eps = sqrt(sqrt(numeric_limits<T>::epsilon()));
> if (abs(x) <= eps) // or "<" ? which is better for performance?
> {
> return(static_cast<T>(1)  x*x/6);
> } else
> {
> return(sin(x)/x);
> }
>
> e.g. using t.C with g++ 2.9.3 one can see a small inaccuracy with
long
> double:
> #include <iostream.h>
> #include <math.h>
>
> using namespace std;
>
> int main()
> {
>
> cout << sizeof(double) << '\t' << sizeof(long double) << endl
<< endl;
>
> cout.precision(20);
>
> cout << "double: " << endl;
> for( double x = 1e7; x >1e8; x *=0.8 )
> {
> cout << x << '\t' << sin(x) << '\t' << sin(x) / x << '\t'
<<
> double(1.0)  x*x/6.0
> << '\t' << double(1.0)  x*x/6.0 + x*x*x*x/120.0 
x*x*x*x*x*x /
> 5040. <<endl;
> }
>
> cout << endl << "long double: " << endl;
>
> for( long double x = 1e7; x >1e8; x *=0.8 )
> {
> cout << x << '\t' << sin(x) << '\t' << sin(x) / x << '\t'
<< 1.0 
> x*x/6.0
> << '\t' << 1.0  x*x/6.0 + x*x*x*x/120.0 
x*x*x*x*x*x / 5040.
> <<endl;
> }
>
> return 0;
> }
>
> Same applies to sinhc, similar to atanh().
> By the way, why is noone using the correct name:
> arctan == arcus tangens; artanh = area tangens hyperbolicus
> Inaddition, what is the deeper meening of sinc_pi, it behaves
> like sinc. It confused me quite some time.
>
>
> Quaternions:
> ===========
>
> There seem to be problem with unneccessary overflows,
> but I didn't have the time to look into detail:
>
> template<typename X>
\
> quaternion<type> & operator /= (quaternion<X> const &
rhs) \
> {
\
> type ar = static_cast<type>(rhs.R_component_1());
\
> type br = static_cast<type>(rhs.R_component_2());
\
> type cr = static_cast<type>(rhs.R_component_3());
\
> type dr = static_cast<type>(rhs.R_component_4());
\
>
\
> type denominator = ar*ar+br*br+cr*cr+dr*dr;
\
>
\
> type at = (+a*ar+b*br+c*cr+d*dr)/denominator;
\
> type bt = (a*br+b*arc*dr+d*cr)/denominator;
\
> type ct = (a*cr+b*dr+c*ard*br)/denominator;
\
> type dt = (a*drb*cr+c*br+d*ar)/denominator;
\
>
\
> a = at;
\
> b = bt;
\
> c = ct;
\
> d = dt;
\
>
\
> return(*this);
\
> }
>
> Can this and similar expressions be reorganized to avoid overflows
if
> at,bt,ct,dt ? If this is not done for speed reasons, it should be
noted,
> that there is a danger of overflow in the denominator, allthough
the result
> itself won't overflow.
>
> Best wishes,
> Peter
Boost list run by bdawes at acm.org, gregod at cs.rpi.edu, cpdaniel at pacbell.net, john at johnmaddock.co.uk