[Help-gsl] integrating a function with externally defined coefficients

[Nathan Moore]

Hello,

I'd like to use GSL to numerically integrate a function with a bunch of 
coefficients that I defined as "global" variables.  At present the 
integration doesn't work, which I imagine might have something to do 
with the "gsl_function" type - which I don't really understand.

Specifically, I define a set of coefficients as globals,

double Ax[FREQUENCY_LIMIT];

// and then inside the main function I assign them values,

main() {

    extern double Ax[FREQUENCY_LIMIT];
...

 for (i = 0; i < FREQUENCY_LIMIT; i++) {
           Ax[i] = gsl_rng_uniform(dice);
}
...

// later, I integrate a function which uses these coefficients,

    gsl_integration_workspace *w = 
gsl_integration_workspace_alloc(1000);
    double *length;
    double *abserror;

    //perform integration with GSL adaptive routine
    //specify function
    gsl_function F;
    F.function = &length_function;
    // integration bounds
    double lower_bound = 0.0;
    double upper_bound = T_parametric_max;
    // desired absolute and relative error
    double epsabs = 0.0;
    double epsrel = 1e-7;
    // maximum number of sub-intervals
    size_t limit = 1000;
    // specify integration key, GSL_INTEG_GAUSS61  (key = 6)
    int key = 6;
    gsl_integration_qag(&F, lower_bound, upper_bound, epsabs,
                        epsrel, limit, key, w, length, abserror);

...
return 1;
}

the function  I use looks like,

double length_function(double t)
{
    extern double Ax[FREQUENCY_LIMIT];
    extern double pie;

    int i = 0;
    double length = 0.0;
    double n, angle, sin_k, cos_k;
    double A_dot_A;

    // sum over frequency modes, i, up to FREQUENCY_LIMIT
    for (i = 0; i < FREQUENCY_LIMIT; i++) {
        n = (double) i;
        angle = 2.0 * pie * n * t / T_parametric_max;
        sin_k = sin(angle);
        cos_k = cos(angle);
        A_dot_A = Ax[i] * Ax[i];
        length += A_dot_A * sin_k * sin_k;    }

    length = sqrt(length);
    return length;
}

I notice in the documentation that they define the gsl_function 
function parameters with a statement like,
          F.params = &alpha;
I don't do that, is that why I get bus errors when I try to run the 
scrpt?

any suggestions would be appreciated.

regards,

Nathan Moore