Loading...

help-gsl@gnu.org

[Prev] Thread [Next]  |  [Prev] Date [Next]

Re: [Help-gsl] gsl precision question Rhys Ulerich Wed Feb 22 10:01:15 2012

Hi Wally,

> For example: For Euler and complex stuff a series of i^n , n= 0,1,2,...
> is useful. Following code should do this but i guess during gsl internal
> calculation (log() ) the process loose precision and i get output as shown
> below.
>
> How can i tell gsl to use a certain precision and/or how to round ?

While this does not answer your question, perhaps it can sidestep it...

If you want to take i^n I recommend directly using the properties of
that operation rather than using gsl_complex_pow_real (which will
probably always exhibit some floating point loss).  In C99 I have used
code like the following to form alpha*i^n for integer n:

static inline
void scale_by_imaginary_power(const complex_double in,
                              complex_double * const out, int p)
{
    // Modulo-four-like operation for 2's complement p
    switch (p & 3) {
        case 0: *out =  in           ; break; // I^0 = 1
        case 1: *out =  in*_Complex_I; break; // I^1 = I = I^-3
        case 2: *out = -in           ; break; // I^2 = -1 = I^-2
        case 3: *out = -in*_Complex_I; break; // I^3 = -I = I^-1
    }
}

That same idea should be translatable into GSL.  This type of logic
will give analytically exact results.

Hope that helps,
Rhys