help-gsl
[Top][All Lists]

## [Help-gsl] Incorrect values with Integration

 From: Philipp Basler Subject: [Help-gsl] Incorrect values with Integration Date: Thu, 23 Jul 2015 14:19:14 +0200

```Hello,

I have an Integrand which reads

static double Integrand(double k , void * params)
{
long double *p = (long double *)params;
long double MassSquared = p;
long double Temp = p;
long double sign = p;

long double EXP = std::exp(-std::sqrt(k*k+MassSquared/(Temp*Temp)));
double res = 0;
if(sign == -1)
{
res = k*k*std::log(1-exp(-std::sqrt(k*k+MassSquared/(Temp*Temp))));
}
else if(sign == +1)
{
res = k*k*std::log(1+exp(-std::sqrt(k*k+MassSquared/(Temp*Temp))));
}
return res;
};

This one is in a Class called VEffClass.

My Integration Routine would then be

long double VEffClass::NumIntInf(long double low, double sign, long double
MassSquared, long double Temp)
{
//long double  res = 0;
double up =
std::sqrt(std::log(C_threshold)*std::log(C_threshold)-MassSquared/(Temp*Temp));
long double StepSize = (up-low)/C_IntSteps;
long double par = {MassSquared, Temp, sign};

if(Temp == 0) return 0;

gsl_integration_workspace *w = gsl_integration_workspace_alloc(5e+5);
double result, error;
gsl_function F;
F.function = &VEffClass::Integrand;
F.params = &par;

gsl_integration_qags(&F, C_threshold, up, 1e-07, 0, 5e+5, w, &result,
&error);
gsl_integration_workspace_free(w);

return result;

}

C_threshold = 1e-3;

If I run NumIntInf(10^-3, +1 , 100, 500) the result is 1.83038. The Result
given by Maple is ~ 6*10^9 . If I use a simple Simpsons 1/3 Method it
differes about 1% from the Maple result.

Can anyone point me to the Error?

Best regards,
Philipp Basler.

```