help-gsl
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: [Help-gsl] problem


From: Pau Cervera Badia
Subject: Re: [Help-gsl] problem
Date: Wed, 24 Aug 2005 20:42:59 +0200
User-agent: Mozilla Thunderbird 1.0.6-1.1.fc3 (X11/20050720)

Whow! I see now why you haven't posted them before ;-)
Sorry, but I think I can't help you with that.

address@hidden wrote:

Thanks a lot. Attached is my function.

Jun Cao



Quoting Pau Cervera Badia <address@hidden>:

I think it means that gsl_odeiv_evolve_apply is trying to adjust h to obtain the prescribed presicion in the integration of the equations, but the step-size value seems too large to obtain that precision, so the routine chooses a smaller step (see [1]). It seems like if one or more of your variables or derivatives is unboundedly growing, so your solutions are diverging, but I'm not sure if it is the only reason because I don't know you're equations.

Maybe you can post your *func* function and we can take a look.

[1] http://www.gnu.org/software/gsl/manual/gsl-ref_25.html#SEC384

address@hidden wrote:

Hi Pau Cervera Badia:

I check h. I found that h is becoming smaller and smaller and finally get to
0.0. What does it mean? Does it mean my program is wrong? Thanks again.

Jun Cao


Quoting Pau Cervera Badia <address@hidden>:



Maybe is changing only a little.
Can you monitor the h values, something as

while (t < t1) {

  int status = gsl_odeiv_evolve_apply(e,c,s,&sys,&t,t1,&h,y);

  printf("%.5e %.5e %.5e\n", t, t1, h);

  if (status != GSL_SUCCESS)
     break;

}

will do. If h is becoming smaller and smaller, something wrong is
happening
(maybe your ecuations are diverging somewhere).


address@hidden wrote:

Hi Pau Cervera Badia:

Thank you for your response. I check my program again. It stuck in
t=5.790250,
t1=6.467000.  After gsl_odeiv_evolve_apply, t always is 5.790250 and does
not
change. So the program stuck in the "while" loop. I try to change h to
1e-2
and
1e-10. I get the same result. Do you have anoter idea? Thanks again.

Jun Cao


Quoting Pau Cervera Badia <address@hidden>:



gsl_odeiv_evolve_apply(e,c,s,&sys,&t,t1,&h,y) will advance from time t to the next integration time with an optimum step-size, provided that the next time is less than t1. Otherwise it will integrate the system to
t1 exactly. So after some iterations the while condition will be false.

You can check that gsl_odeiv_evolve_apply is running properly with:

while (t < t1) {

 int status = gsl_odeiv_evolve_apply(e,c,s,&sys,&t,t1,&h,y);

 printf("%.5e %.5e\n", t, t1); // this should print t less than t1

 if (status != GSL_SUCCESS)
    break;

}
printf("%.5e %.5e\n",t,t1); // this should print t equal to t1

If I'm not missing something, your program should work. Maybe the problem is that t is always less than t1 for some other reason (maybe the routine couldn't reach the prescribed precision and gsl_odeiv_evolve_apply is trying to reduce the step-size more and more).
Maybe you can check the h values.

Maybe something will be of any help.

address@hidden wrote:

Hi:

I use ordinary differential equations from GSL. The part of my program
is
like
example shown in the document.

const gsl_odeiv_step_type * T = gsl_odeiv_step_rk8pd;

gsl_odeiv_step * s = gsl_odeiv_step_alloc (T, 9); gsl_odeiv_control * c = gsl_odeiv_control_y_new (1e-6, 0.0); gsl_odeiv_evolve * e = gsl_odeiv_evolve_alloc (9);

double mu = 10;
gsl_odeiv_system sys = {func, NULL, 9, &mu};

double t = 0.0; double h = 1e-6;

for (i = 1; i <no_of_data_points; i++) {
 t1 = mytable[i];
 while (t < t1)
 {
    int status = gsl_odeiv_evolve_apply (e, c, s,
&sys, &t, t1,
                                       &h, y);

  if (status != GSL_SUCCESS)
      break;
  }
  ...
}
gsl_odeiv_evolve_free (e);
gsl_odeiv_control_free (c);
gsl_odeiv_step_free (s);

The weird thing is that at some i, the program stuck inside of "while"
loop
and
never come out. gsl_odeiv_evolve_apply will give t=t, so t is always <
t1.
As
result the program never stop. Dose anyone know how to solve this
problem?
Thanks.

Jun Cao



_______________________________________________
Help-gsl mailing list
address@hidden
http://lists.gnu.org/mailman/listinfo/help-gsl




--
Pau Cervera i Badia (e-mail address@hidden)
{
Departament de Física Fonamental                 Martí i Franqués, 1
Universitat de Barcelona                   Planta 3, despatx 346 bis
                                                     08028 Barcelona
tel: +34 934 921 155                                           Spain

"To err is human, but to really foul things up requires a computer."
}



--
Pau Cervera i Badia (e-mail address@hidden)
{
 Departament de Física Fonamental                 Martí i Franqués, 1
 Universitat de Barcelona                   Planta 3, despatx 346 bis
                                                      08028 Barcelona
 tel: +34 934 921 155                                           Spain

 "To err is human, but to really foul things up requires a computer."
}





--
Pau Cervera i Badia (e-mail address@hidden)
{
  Departament de Física Fonamental                 Martí i Franqués, 1
  Universitat de Barcelona                   Planta 3, despatx 346 bis
                                                       08028 Barcelona
  tel: +34 934 921 155                                           Spain

  "To err is human, but to really foul things up requires a computer."
}





------------------------------------------------------------------------

#include "mymodel.h"


/*******************************************************************************
/*******************************************************************************
*                                                                              *
*   function requied by imsl_f_ode_runge_kutta                                  
   *
*                                                                              *
*******************************************************************************
*******************************************************************************/
int ode_fcn (double t, const double y[], double yprime[], void *params)
{

k = (double *)params;
yprime[0] =0  ;

yprime[1] 
=(+(+(-y[37]-y[36]-y[35])*0.100000-y[3]-2.000000*y[4]-3.000000*y[5]-4.000000*y[6]-5.000000*y[7]-6.000000*y[8]-7.000000*y[9]-8.000000*y[10]-9.000000*y[11]-10.000000*y[12]-11.000000*y[13]-12.000000*y[14]-13.000000*y[15]-y[38]-y[39]-y[40]-y[41]-y[42]-y[43]-y[44]-y[45]-y[46]-y[47]-y[48]-y[49]-y[69]-y[70]-y[71]-y[72]-y[73]-y[74]-y[75]-y[76]-y[77]-y[78]-y[79]-y[80]-y[81]-y[82]-y[83])*y[1]+(+(+y[14]+y[13]+y[12]+y[11]+y[10]+y[9]+y[8]+y[7]+y[6]+y[5]+y[4]+y[3])*y[2]+(+y[12]+y[11]+y[10]+y[9]+y[8]+y[7]+y[6]+y[5]+y[3])*y[4]+(+y[13]+y[12]+y[11]+y[10]+y[9]+y[6]+y[5])*y[3]+(+y[11]+y[10]+y[9]+y[8]+y[7]+y[6])*y[5]+(+y[10]+y[9]+y[8]+y[7])*y[6]+(+y[9]+y[8])*y[7])*2.000000+(+4.000000*y[3]+y[7])*y[7]+y[8]*y[8]+y[6]*y[6]+y[5]*y[5]+y[4]*y[4]+y[3]*y[3]+y[2]*y[2])*k[3]+(-k[4]-k[0])*y[1]+k[1]*y[53]*y[53]+k[2]*y[33]*y[32]
 ;

yprime[2] 
=(+(+(+y[15]+y[14]+y[13]+y[12]+y[11]+y[10]+y[9]+y[8]+y[7]+y[6]+y[5]+y[4]+y[3])*2.000000+(+y[37]+y[36]+y[35])*0.100000+y[38]+y[39]+y[40]+y[41]+y[42]+y[43]+y[44]+y[45]+y[46]+y[47]+y[48]+y[49])*y[1]+(+(-y[14]-y[13]-y[12]-y[11]-y[10]-y[9]-y[8]-y[7]-y[6]-y[5]-y[4]-y[3]-y[2])*y[2])*2.000000)*k[3]+(+(-k[0]-k[4])*y[2])*2.000000+k[1]*y[53]*y[54]
 ;

yprime[3] 
=(+(+(+y[15]+y[14]+y[13]+y[12]+y[11]+y[10]+y[9]+y[8]+y[7]+y[6]+y[5]+y[4])*y[1]+(-y[13]-y[12]-y[11]-y[10]-y[9]-y[6]-y[5]-y[4]-y[3]-y[2])*y[3])*2.000000+(-y[1]-4.000000*y[7])*y[3]+y[2]*y[2])*k[3]+(-2.000000*k[4]-3.000000*k[0])*y[3]+(+y[54]*y[54]+y[53]*y[55])*k[1]
 ;

yprime[4] 
=(+(+(+y[15]+y[14]+y[13]+y[12]+y[11]+y[10]+y[9]+y[8]+y[7]+y[6]+y[5]-y[4])*y[1]+(-y[12]-y[11]-y[10]-y[9]-y[8]-y[7]-y[6]-y[5]-y[4]-y[3]-y[2])*y[4]+y[2]*y[3])*k[3]-k[4]*y[4])*2.000000+(+y[54]*y[55]+y[53]*y[56])*k[1]-4.000000*k[0]*y[4]
 ;

yprime[5] 
=(+(+(+y[15]+y[14]+y[13]+y[12]+y[11]+y[10]+y[9]+y[8]+y[7]+y[6])*y[1]+(-y[11]-y[10]-y[9]-y[8]-y[7]-y[6]-y[5]-y[4]-y[3]-y[2])*y[5]+y[2]*y[4])*2.000000+y[3]*y[3]-3.000000*y[1]*y[5])*k[3]+(+y[55]*y[55]+y[54]*y[56]+y[53]*y[57])*k[1]+(-2.000000*k[4]-5.000000*k[0])*y[5]
 ;

yprime[6] 
=(+(+(+y[15]+y[14]+y[13]+y[12]+y[11]+y[10]+y[9]+y[8]+y[7])*y[1]+(-y[10]-y[9]-y[8]-y[7]-y[6]-y[5]-y[4]-y[3]-y[2])*y[6]+y[2]*y[5]+y[3]*y[4])*2.000000-4.000000*y[1]*y[6])*k[3]+(+y[55]*y[56]+y[54]*y[57]+y[53]*y[58])*k[1]+(-2.000000*k[4]-6.000000*k[0])*y[6]
 ;

yprime[7] 
=(+(+(+y[15]+y[14]+y[13]+y[12]+y[11]+y[9]+y[8])*y[1]+(-y[9]-y[8]-y[7]-y[6]-y[5]-y[4]-y[2])*y[7]+y[2]*y[6]+y[3]*y[5])*2.000000+(+y[1]*y[10]-y[3]*y[7])*4.000000+y[4]*y[4]-5.000000*y[1]*y[7])*k[3]+(+y[56]*y[56]+y[55]*y[57]+y[54]*y[58]+y[53]*y[59])*k[1]+(-2.000000*k[4]-7.000000*k[0])*y[7]
 ;

yprime[8] 
=(+(+(+y[15]+y[14]+y[13]+y[12]+y[11]+y[9])*y[1]+(-y[8]-y[7]-y[6]-y[5]-y[4]-y[2])*y[8]+y[2]*y[7]+y[3]*y[6]+y[4]*y[5])*2.000000-6.000000*y[1]*y[8])*k[3]+(+y[56]*y[57]+y[55]*y[58]+y[54]*y[59]+y[53]*y[60])*k[1]+(-2.000000*k[4]-8.000000*k[0])*y[8]
 ;

yprime[9] 
=(+(+(+y[15]+y[14]+y[13]+y[12]+y[11]+y[10])*y[1]+(-y[7]-y[6]-y[5]-y[4]-y[3]-y[2])*y[9]+y[2]*y[8]+y[3]*y[7]+y[4]*y[6])*2.000000+y[5]*y[5]-7.000000*y[1]*y[9])*k[3]+(+y[57]*y[57]+y[56]*y[58]+y[55]*y[59]+y[54]*y[60]+y[53]*y[61])*k[1]+(-2.000000*k[4]-9.000000*k[0])*y[9]
 ;

yprime[10] 
=(+(+(+y[15]+y[14]+y[13]+y[12]+y[11])*y[1]+(-y[6]-y[5]-y[4]-y[3]-y[2])*y[10]+(+y[4]+y[3])*y[7]+y[2]*y[9]+y[5]*y[6])*2.000000-8.000000*y[1]*y[10])*k[3]+(+y[57]*y[58]+y[56]*y[59]+y[55]*y[60]+y[54]*y[61]+y[53]*y[62])*k[1]+(-2.000000*k[4]-10.000000*k[0])*y[10]
 ;

yprime[11] 
=(+(+(+y[15]+y[14]+y[13]+y[12])*y[1]+(-y[5]-y[4]-y[3]-y[2])*y[11]+y[2]*y[10]+y[3]*y[9]+y[4]*y[8]+y[5]*y[7])*2.000000+y[6]*y[6]-9.000000*y[1]*y[11])*k[3]+(+y[58]*y[58]+y[57]*y[59]+y[56]*y[60]+y[55]*y[61]+y[54]*y[62]+y[53]*y[63])*k[1]+(-2.000000*k[4]-11.000000*k[0])*y[11]
 ;

yprime[12] 
=(+(+(+y[15]+y[14]+y[13])*y[1]+(-y[4]-y[3]-y[2])*y[12]+y[2]*y[11]+y[3]*y[10]+y[4]*y[9]+y[5]*y[8]+y[6]*y[7])*2.000000-10.000000*y[1]*y[12])*k[3]+(+y[58]*y[59]+y[57]*y[60]+y[56]*y[61]+y[55]*y[62]+y[54]*y[63]+y[53]*y[64])*k[1]+(-2.000000*k[4]-12.000000*k[0])*y[12]
 ;

yprime[13] 
=(+(+(+y[15]+y[14])*y[1]+(-y[13]+y[12])*y[2]+(-y[13]+y[11])*y[3]+y[4]*y[10]+y[5]*y[9]+y[6]*y[8])*2.000000+y[7]*y[7]-11.000000*y[1]*y[13])*k[3]+(+y[59]*y[59]+y[58]*y[60]+y[57]*y[61]+y[56]*y[62]+y[55]*y[63]+y[54]*y[64]+y[53]*y[65])*k[1]+(-2.000000*k[4]-13.000000*k[0])*y[13]
 ;

yprime[14] 
=(+(+(-y[14]+y[13])*y[2]+y[1]*y[15]+y[3]*y[12]+y[4]*y[11]+y[5]*y[10]+y[6]*y[9]+y[7]*y[8])*2.000000-12.000000*y[1]*y[14])*k[3]+(+y[59]*y[60]+y[58]*y[61]+y[57]*y[62]+y[56]*y[63]+y[55]*y[64]+y[54]*y[65]+y[53]*y[66])*k[1]+(-2.000000*k[4]-14.000000*k[0])*y[14]
 ;

yprime[15] 
=(+y[60]*y[60]+y[59]*y[61]+y[58]*y[62]+y[57]*y[63]+y[56]*y[64]+y[55]*y[65]+y[54]*y[66]+y[53]*y[67])*k[1]+(+(+y[2]*y[14]+y[3]*y[13]+y[4]*y[12]+y[5]*y[11]+y[6]*y[10]+y[7]*y[9])*2.000000+y[8]*y[8]-13.000000*y[1]*y[15])*k[3]+(-2.000000*k[4]-15.000000*k[0])*y[15]
 ;

yprime[16] =+k[5]*y[53]+k[4]*y[1] ;

yprime[17] =-k[6]*y[17]+k[5]*y[54]+2.000000*k[4]*y[2]+k[3]*y[69]*y[1] ;

yprime[18] =(+k[4]*y[3]-k[6]*y[18])*2.000000+k[5]*y[55]+k[3]*y[70]*y[1] ;

yprime[19] =-3.000000*k[6]*y[19]+k[5]*y[56]+2.000000*k[4]*y[4]+k[3]*y[71]*y[1] ;

yprime[20] =-4.000000*k[6]*y[20]+k[5]*y[57]+2.000000*k[4]*y[5]+k[3]*y[72]*y[1] ;

yprime[21] =-5.000000*k[6]*y[21]+k[5]*y[58]+2.000000*k[4]*y[6]+k[3]*y[73]*y[1] ;

yprime[22] =-6.000000*k[6]*y[22]+k[5]*y[59]+2.000000*k[4]*y[7]+k[3]*y[74]*y[1] ;

yprime[23] =-7.000000*k[6]*y[23]+k[5]*y[60]+2.000000*k[4]*y[8]+k[3]*y[75]*y[1] ;

yprime[24] =-8.000000*k[6]*y[24]+k[5]*y[61]+2.000000*k[4]*y[9]+k[3]*y[76]*y[1] ;

yprime[25] =-9.000000*k[6]*y[25]+k[5]*y[62]+2.000000*k[4]*y[10]+k[3]*y[77]*y[1] 
;

yprime[26] 
=-10.000000*k[6]*y[26]+k[5]*y[63]+2.000000*k[4]*y[11]+k[3]*y[78]*y[1] ;

yprime[27] 
=-11.000000*k[6]*y[27]+k[5]*y[64]+2.000000*k[4]*y[12]+k[3]*y[79]*y[1] ;

yprime[28] 
=-12.000000*k[6]*y[28]+k[5]*y[65]+2.000000*k[4]*y[13]+k[3]*y[80]*y[1] ;

yprime[29] 
=-13.000000*k[6]*y[29]+k[5]*y[66]+2.000000*k[4]*y[14]+k[3]*y[81]*y[1] ;

yprime[30] 
=-14.000000*k[6]*y[30]+k[5]*y[67]+2.000000*k[4]*y[15]+k[3]*y[82]*y[1] ;

yprime[31] =-15.000000*k[6]*y[31]+k[5]*y[68]+k[3]*y[83]*y[1] ;

yprime[32] 
=(+(+y[2]+y[3]+y[4]+y[5]+y[6]+y[7]+y[8]+y[9]+y[10]+y[11]+y[12]+y[13]+y[14]+y[15])*2.000000+y[1])*k[4]+(-k[2]*y[32]+k[7])*y[33]
 ;

yprime[33] =(-k[2]*y[32]-k[7])*y[33] ;

yprime[34] =+k[5]*y[69]+0.100000*k[3]*y[35]*y[1] ;

yprime[35] 
=(-0.100000*k[3]*y[1]-k[8])*y[35]+k[5]*y[70]+0.100000*k[3]*y[36]*y[1] ;

yprime[36] 
=(-0.100000*k[3]*y[1]-k[8])*y[36]+k[5]*y[71]+0.100000*k[3]*y[37]*y[1] ;

yprime[37] =(-0.100000*k[3]*y[1]-k[8])*y[37]+k[5]*y[72]+k[3]*y[38]*y[1] ;

yprime[38] =(-k[3]*y[1]-k[8])*y[38]+k[5]*y[73]+k[3]*y[39]*y[1] ;

yprime[39] =(-k[3]*y[1]-k[8])*y[39]+k[5]*y[74]+k[3]*y[40]*y[1] ;

yprime[40] =(-k[3]*y[1]-k[8])*y[40]+k[5]*y[75]+k[3]*y[41]*y[1] ;

yprime[41] =(-k[3]*y[1]-k[8])*y[41]+k[5]*y[76]+k[3]*y[42]*y[1] ;

yprime[42] =(-k[3]*y[1]-k[8])*y[42]+k[5]*y[77]+k[3]*y[43]*y[1] ;

yprime[43] =(-k[3]*y[1]-k[8])*y[43]+k[5]*y[78]+k[3]*y[44]*y[1] ;

yprime[44] =(-k[3]*y[1]-k[8])*y[44]+k[5]*y[79]+k[3]*y[45]*y[1] ;

yprime[45] =(-k[3]*y[1]-k[8])*y[45]+k[5]*y[80]+k[3]*y[46]*y[1] ;

yprime[46] =(-k[3]*y[1]-k[8])*y[46]+k[5]*y[81]+k[3]*y[47]*y[1] ;

yprime[47] =(-k[3]*y[1]-k[8])*y[47]+k[5]*y[82]+k[3]*y[48]*y[1] ;

yprime[48] =(-k[3]*y[1]-k[8])*y[48]+k[5]*y[83]+k[3]*y[49]*y[1] ;

yprime[49] =(-k[3]*y[1]-k[8])*y[49]+k[5]*y[84] ;

yprime[50] =0  ;

yprime[51] =(+k[2]*y[32]+k[7])*y[33] ;

yprime[52] 
=(+(+y[69]+y[70]+y[71]+y[72]+y[73]+y[74]+y[75])*0.100000+(+y[76]+y[77]+y[78]+y[79]+y[80]+y[81])*0.300000+(+y[82]+y[83]+y[84])*0.500000)*k[5]+(+y[49]+y[48]+y[47]+y[46]+y[45]+y[44]+y[43]+y[42]+y[41]+y[40]+y[39]+y[38]+y[37]+y[36]+y[35])*k[8]
 ;

yprime[53] 
=(+(-2.000000*y[53]-y[54]-y[55]-y[56]-y[57]-y[58]-y[59]-y[60]-y[61]-y[62]-y[63]-y[64]-y[65]-y[66]-y[67]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*k[1]-k[5])*y[53]+(+(+y[69]+y[70]+y[71]+y[72]+y[73]+y[74]+y[75]+y[76]+y[77]+y[78]+y[79]+y[80]+y[81]+y[82]+y[83])*k[3]+2.000000*k[0])*y[1]+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24]+y[23]+y[22]+y[21]+y[20]+y[19]+y[18]+y[17])*k[6]+(+(+y[2]+y[3]+y[4]+y[5]+y[6]+y[7]+y[8]+y[9]+y[10]+y[11]+y[12]+y[13]+y[14]+y[15])*k[0])*2.000000
 ;

yprime[54] 
=(+(-y[53]-2.000000*y[54]-y[55]-y[56]-y[57]-y[58]-y[59]-y[60]-y[61]-y[62]-y[63]-y[64]-y[65]-y[66]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[54]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[53])*k[1]+(+(+y[2]+y[3]+y[4]+y[5]+y[6]+y[7]+y[8]+y[9]+y[10]+y[11]+y[12]+y[13]+y[14]+y[15])*k[0])*2.000000+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24]+y[23]+y[22]+y[21]+y[20]+y[19]+y[18])*k[6]-k[5]*y[54]
 ;

yprime[55] 
=(+(-y[53]-y[54]-2.000000*y[55]-y[56]-y[57]-y[58]-y[59]-y[60]-y[61]-y[62]-y[63]-y[64]-y[65]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[55]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[54])*k[1]+(+(+y[3]+y[4]+y[5]+y[6]+y[7]+y[8]+y[9]+y[10]+y[11]+y[12]+y[13]+y[14]+y[15])*k[0])*2.000000+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24]+y[23]+y[22]+y[21]+y[20]+y[19])*k[6]-k[5]*y[55]
 ;

yprime[56] 
=(+(-y[53]-y[54]-y[55]-2.000000*y[56]-y[57]-y[58]-y[59]-y[60]-y[61]-y[62]-y[63]-y[64]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[56]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[55])*k[1]+(+(+y[4]+y[5]+y[6]+y[7]+y[8]+y[9]+y[10]+y[11]+y[12]+y[13]+y[14]+y[15])*k[0])*2.000000+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24]+y[23]+y[22]+y[21]+y[20])*k[6]-k[5]*y[56]
 ;

yprime[57] 
=(+(-y[53]-y[54]-y[55]-y[56]-2.000000*y[57]-y[58]-y[59]-y[60]-y[61]-y[62]-y[63]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[57]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[56])*k[1]+(+(+y[5]+y[6]+y[7]+y[8]+y[9]+y[10]+y[11]+y[12]+y[13]+y[14]+y[15])*k[0])*2.000000+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24]+y[23]+y[22]+y[21])*k[6]-k[5]*y[57]
 ;

yprime[58] 
=(+(-y[53]-y[54]-y[55]-y[56]-y[57]-2.000000*y[58]-y[59]-y[60]-y[61]-y[62]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[58]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[57])*k[1]+(+(+y[6]+y[7]+y[8]+y[9]+y[10]+y[11]+y[12]+y[13]+y[14]+y[15])*k[0])*2.000000+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24]+y[23]+y[22])*k[6]-k[5]*y[58]
 ;

yprime[59] 
=(+(-y[53]-y[54]-y[55]-y[56]-y[57]-y[58]-2.000000*y[59]-y[60]-y[61]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[59]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[58])*k[1]+(+(+y[7]+y[8]+y[9]+y[10]+y[11]+y[12]+y[13]+y[14]+y[15])*k[0])*2.000000+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24]+y[23])*k[6]-k[5]*y[59]
 ;

yprime[60] 
=(+(-y[53]-y[54]-y[55]-y[56]-y[57]-y[58]-y[59]-2.000000*y[60]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[60]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[59])*k[1]+(+(+y[8]+y[9]+y[10]+y[11]+y[12]+y[13]+y[14]+y[15])*k[0])*2.000000+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24])*k[6]-k[5]*y[60]
 ;

yprime[61] 
=(+(-y[53]-y[54]-y[55]-y[56]-y[57]-y[58]-y[59]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[61]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[60])*k[1]+(+(+y[9]+y[10]+y[11]+y[12]+y[13]+y[14]+y[15])*k[0])*2.000000+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25])*k[6]-k[5]*y[61]
 ;

yprime[62] 
=(+(-y[53]-y[54]-y[55]-y[56]-y[57]-y[58]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[62]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[61])*k[1]+(+(+y[10]+y[11]+y[12]+y[13]+y[14]+y[15])*k[0])*2.000000+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26])*k[6]-k[5]*y[62]
 ;

yprime[63] 
=(+(-y[53]-y[54]-y[55]-y[56]-y[57]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[63]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[62])*k[1]+(+(+y[11]+y[12]+y[13]+y[14]+y[15])*k[0])*2.000000+(+y[31]+y[30]+y[29]+y[28]+y[27])*k[6]-k[5]*y[63]
 ;

yprime[64] 
=(+(-y[53]-y[54]-y[55]-y[56]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[64]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[63])*k[1]+(+(+y[12]+y[13]+y[14]+y[15])*k[0])*2.000000+(+y[31]+y[30]+y[29]+y[28])*k[6]-k[5]*y[64]
 ;

yprime[65] 
=(+(-y[53]-y[54]-y[55]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[65]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[64])*k[1]+(+(+y[13]+y[14]+y[15])*k[0])*2.000000+(+y[31]+y[30]+y[29])*k[6]-k[5]*y[65]
 ;

yprime[66] 
=(+(-y[53]-y[54]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[66]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[65])*k[1]+(+(+y[14]+y[15])*k[0])*2.000000+(+y[31]+y[30])*k[6]-k[5]*y[66]
 ;

yprime[67] 
=(+(-y[53]-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[67]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[66])*k[1]+2.000000*k[0]*y[15]+k[6]*y[31]-k[5]*y[67]
 ;

yprime[68] 
=(+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[67])*k[1]-k[5]*y[68] ;

yprime[69] 
=(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24]+y[23]+y[22]+y[21]+y[20]+y[19]+y[18]+y[17])*k[6]+(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*k[1]-k[3]*y[1]-1.100000*k[5])*y[69]
 ;

yprime[70] 
=(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[70]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[69])*k[1]+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24]+y[23]+y[22]+y[21]+y[20]+y[19]+y[18])*k[6]+(-k[3]*y[1]-1.100000*k[5])*y[70]
 ;

yprime[71] 
=(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[71]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[70])*k[1]+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24]+y[23]+y[22]+y[21]+y[20]+y[19])*k[6]+(-k[3]*y[1]-1.100000*k[5])*y[71]
 ;

yprime[72] 
=(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[72]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[71])*k[1]+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24]+y[23]+y[22]+y[21]+y[20])*k[6]+(-k[3]*y[1]-1.100000*k[5])*y[72]
 ;

yprime[73] 
=(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[73]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[72])*k[1]+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24]+y[23]+y[22]+y[21])*k[6]+(-k[3]*y[1]-1.100000*k[5])*y[73]
 ;

yprime[74] 
=(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[74]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[73])*k[1]+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24]+y[23]+y[22])*k[6]+(-k[3]*y[1]-1.100000*k[5])*y[74]
 ;

yprime[75] 
=(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[75]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[74])*k[1]+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24]+y[23])*k[6]+(-k[3]*y[1]-1.100000*k[5])*y[75]
 ;

yprime[76] 
=(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[76]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[75])*k[1]+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25]+y[24])*k[6]+(-k[3]*y[1]-1.300000*k[5])*y[76]
 ;

yprime[77] 
=(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[77]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[76])*k[1]+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26]+y[25])*k[6]+(-k[3]*y[1]-1.300000*k[5])*y[77]
 ;

yprime[78] 
=(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[78]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[77])*k[1]+(+y[31]+y[30]+y[29]+y[28]+y[27]+y[26])*k[6]+(-k[3]*y[1]-1.300000*k[5])*y[78]
 ;

yprime[79] 
=(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[79]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[78])*k[1]+(+y[31]+y[30]+y[29]+y[28]+y[27])*k[6]+(-k[3]*y[1]-1.300000*k[5])*y[79]
 ;

yprime[80] 
=(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[80]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[79])*k[1]+(+y[31]+y[30]+y[29]+y[28])*k[6]+(-k[3]*y[1]-1.300000*k[5])*y[80]
 ;

yprime[81] 
=(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[81]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[80])*k[1]+(+y[31]+y[30]+y[29])*k[6]+(-k[3]*y[1]-1.300000*k[5])*y[81]
 ;

yprime[82] 
=(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[82]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[81])*k[1]+(+y[31]+y[30])*k[6]+(-k[3]*y[1]-1.500000*k[5])*y[82]
 ;

yprime[83] 
=(+(-y[85]-y[86]-y[87]-y[88]-y[89]-y[90]-y[91]-y[92])*y[83]+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[82])*k[1]+(-k[3]*y[1]-1.500000*k[5])*y[83]+k[6]*y[31]
 ;

yprime[84] 
=(+(+y[85]+y[86]+y[87]+y[88]+y[89]+y[90]+y[91]+y[92])*y[83])*k[1]-1.500000*k[5]*y[84]
 ;

yprime[85] 
=(+(+y[53]+y[54]+y[55]+y[56]+y[57]+y[58]+y[59]+y[60]+y[61]+y[62]+y[63]+y[64]+y[65]+y[66]+y[67]+y[69]+y[70]+y[71]+y[72]+y[73]+y[74]+y[75]+y[76]+y[77]+y[78]+y[79]+y[80]+y[81]+y[82]+y[83])*y[86]+(-y[53]-y[54]-y[55]-y[56]-y[57]-y[58]-y[59]-y[60]-y[61]-y[62]-y[63]-y[64]-y[65]-y[66]-y[67]-y[69]-y[70]-y[71]-y[72]-y[73]-y[74]-y[75]-y[76]-y[77]-y[78]-y[79]-y[80]-y[81]-y[82]-y[83])*y[85])*k[1]
 ;

yprime[86] 
=(+(+y[53]+y[54]+y[55]+y[56]+y[57]+y[58]+y[59]+y[60]+y[61]+y[62]+y[63]+y[64]+y[65]+y[66]+y[67]+y[69]+y[70]+y[71]+y[72]+y[73]+y[74]+y[75]+y[76]+y[77]+y[78]+y[79]+y[80]+y[81]+y[82]+y[83])*y[87]+(-y[53]-y[54]-y[55]-y[56]-y[57]-y[58]-y[59]-y[60]-y[61]-y[62]-y[63]-y[64]-y[65]-y[66]-y[67]-y[69]-y[70]-y[71]-y[72]-y[73]-y[74]-y[75]-y[76]-y[77]-y[78]-y[79]-y[80]-y[81]-y[82]-y[83])*y[86])*k[1]
 ;

yprime[87] 
=(+(+y[53]+y[54]+y[55]+y[56]+y[57]+y[58]+y[59]+y[60]+y[61]+y[62]+y[63]+y[64]+y[65]+y[66]+y[67]+y[69]+y[70]+y[71]+y[72]+y[73]+y[74]+y[75]+y[76]+y[77]+y[78]+y[79]+y[80]+y[81]+y[82]+y[83])*y[88]+(-y[53]-y[54]-y[55]-y[56]-y[57]-y[58]-y[59]-y[60]-y[61]-y[62]-y[63]-y[64]-y[65]-y[66]-y[67]-y[69]-y[70]-y[71]-y[72]-y[73]-y[74]-y[75]-y[76]-y[77]-y[78]-y[79]-y[80]-y[81]-y[82]-y[83])*y[87])*k[1]
 ;

yprime[88] 
=(+(+y[53]+y[54]+y[55]+y[56]+y[57]+y[58]+y[59]+y[60]+y[61]+y[62]+y[63]+y[64]+y[65]+y[66]+y[67]+y[69]+y[70]+y[71]+y[72]+y[73]+y[74]+y[75]+y[76]+y[77]+y[78]+y[79]+y[80]+y[81]+y[82]+y[83])*y[89]+(-y[53]-y[54]-y[55]-y[56]-y[57]-y[58]-y[59]-y[60]-y[61]-y[62]-y[63]-y[64]-y[65]-y[66]-y[67]-y[69]-y[70]-y[71]-y[72]-y[73]-y[74]-y[75]-y[76]-y[77]-y[78]-y[79]-y[80]-y[81]-y[82]-y[83])*y[88])*k[1]
 ;

yprime[89] 
=(+(+y[53]+y[54]+y[55]+y[56]+y[57]+y[58]+y[59]+y[60]+y[61]+y[62]+y[63]+y[64]+y[65]+y[66]+y[67]+y[69]+y[70]+y[71]+y[72]+y[73]+y[74]+y[75]+y[76]+y[77]+y[78]+y[79]+y[80]+y[81]+y[82]+y[83])*y[90]+(-y[53]-y[54]-y[55]-y[56]-y[57]-y[58]-y[59]-y[60]-y[61]-y[62]-y[63]-y[64]-y[65]-y[66]-y[67]-y[69]-y[70]-y[71]-y[72]-y[73]-y[74]-y[75]-y[76]-y[77]-y[78]-y[79]-y[80]-y[81]-y[82]-y[83])*y[89])*k[1]
 ;

yprime[90] 
=(+(+y[53]+y[54]+y[55]+y[56]+y[57]+y[58]+y[59]+y[60]+y[61]+y[62]+y[63]+y[64]+y[65]+y[66]+y[67]+y[69]+y[70]+y[71]+y[72]+y[73]+y[74]+y[75]+y[76]+y[77]+y[78]+y[79]+y[80]+y[81]+y[82]+y[83])*y[91]+(-y[53]-y[54]-y[55]-y[56]-y[57]-y[58]-y[59]-y[60]-y[61]-y[62]-y[63]-y[64]-y[65]-y[66]-y[67]-y[69]-y[70]-y[71]-y[72]-y[73]-y[74]-y[75]-y[76]-y[77]-y[78]-y[79]-y[80]-y[81]-y[82]-y[83])*y[90])*k[1]
 ;

yprime[91] 
=(+(+y[53]+y[54]+y[55]+y[56]+y[57]+y[58]+y[59]+y[60]+y[61]+y[62]+y[63]+y[64]+y[65]+y[66]+y[67]+y[69]+y[70]+y[71]+y[72]+y[73]+y[74]+y[75]+y[76]+y[77]+y[78]+y[79]+y[80]+y[81]+y[82]+y[83])*y[92]+(-y[53]-y[54]-y[55]-y[56]-y[57]-y[58]-y[59]-y[60]-y[61]-y[62]-y[63]-y[64]-y[65]-y[66]-y[67]-y[69]-y[70]-y[71]-y[72]-y[73]-y[74]-y[75]-y[76]-y[77]-y[78]-y[79]-y[80]-y[81]-y[82]-y[83])*y[91])*k[1]
 ;

yprime[92] 
=(+(-y[53]-y[54]-y[55]-y[56]-y[57]-y[58]-y[59]-y[60]-y[61]-y[62]-y[63]-y[64]-y[65]-y[66]-y[67]-y[69]-y[70]-y[71]-y[72]-y[73]-y[74]-y[75]-y[76]-y[77]-y[78]-y[79]-y[80]-y[81]-y[82]-y[83])*y[92])*k[1]
 ;

yprime[93] 
=(+y[68]+y[67]+y[66]+y[65]+y[64]+y[63]+y[62]+y[61]+y[60]+y[59]+y[58]+y[57]+y[56]+y[55]+y[54]+y[53])*k[5]
 ;

yprime[94] =0  ;

        return GSL_SUCCESS;
}

--
Pau Cervera i Badia (e-mail address@hidden)
{
  Departament de Física Fonamental                 Martí i Franqués, 1
  Universitat de Barcelona                   Planta 3, despatx 346 bis
                                                       08028 Barcelona
  tel: +34 934 921 155                                           Spain

  "To err is human, but to really foul things up requires a computer."
}






reply via email to

[Prev in Thread] Current Thread [Next in Thread]