GNU Scientific library in C. Error when solving a system of ODEs. evolve.c:317: ERROR: step dimension must match evolution size

Hi everyone

I'm trying to solve a system of differential equations using C with GSL version 2.5. For it, I used the gsl_odeiv2_evolve_apply_fixed_step function.

The function that generates this message is the following, but the problematic part of the code is inside the while loop:

int cftes(double *params, double *vec, double j, int n){

double g = params[0]; double xo = params[1];
double A = params[2]; double Omega = params[3];
double w = params[4]; double phioM = params[5];
double phioL = params[6]; double lambda = params[7];
double params2[9] = {g,xo,A,Omega,w,phioM,phioL,lambda,j};

size_t dim = size_t(2*j+1), dimtot = 2*dim*dim;
gsl_odeiv2_system sys = {edo1, NULL, 2*dim*dim, params2};

int i=0;
double t = 0.0;
double tf = 2*M_PI/Omega;
double y[dimtot];
double h = tf/(double(n));

const gsl_odeiv2_step_type *T
= gsl_odeiv2_step_rk8pd;

gsl_odeiv2_step *s = gsl_odeiv2_step_alloc(T,dimtot);
gsl_odeiv2_control * c = gsl_odeiv2_control_y_new(1e-12,0);
gsl_odeiv2_evolve * e = gsl_odeiv2_evolve_alloc(dimtot);

for(int l=0;l<dim;l++){
for(int k=0;k<dim;k++){
if(k==l){
y[2*k+2*l*dim]=1; y[2*l*dim+2*k+1]=0;
}

else{
y[2*l*dim+2*k]=0;
y[2*l*dim+2*k+1]=0;
}
}
}

while(t<tf){
int status = gsl_odeiv2_evolve_apply_fixed_step(e,c,s,
&sys,&t,h,y);

printf("%d %zu\n",i,s->dimension);
if(status != GSL_SUCCESS){
printf("error, return value = %d",status);
break;
}


for(int l=0;l<dimtot;l++){
vec[l+dimtot*i] = y[l];
} i++;

}


gsl_odeiv2_evolve_free (e);
gsl_odeiv2_control_free (c);
gsl_odeiv2_step_free (s);

return 0;
}

The dimtot stands for the size of the system, and is equal in the evolve and step functions.

After compiling and running the code, 380 steps of the while loop execute properly and the the following message appears:

gsl: evolve.c:317: ERROR: step dimension must match evolution size
Default GSL error handler invoked.


Checking the ode-initval/evolve.c on github this message appears whenever evolve and step functions have not the same dimension. Nevertheless, in all the succesfull steps in the while loop the dimension of the step and evolve functions remain constant and equal to eight (I checked it on distinct runs). I don't understand why this message appears and why would the dimension of e or s changesduring the while loop.

I would really appreciate any suggestions on how to solve this problem.

Thanks
Diego Morachis Galindo