Rungekutta à pas adaptatif

Description

un intégrateur d'équations différentielles doit avoir un algorithme pour modifier son pas pour etre plus performant.
voici une comparaison de deux mehtodes (exemple d'un corps en gravitation avec un autre corps tres massique) :
  • Euler avec pas fixe (voyez comment cela diverge en un tour)
  • Runge Kutta avec pas adaptatif (la solution reste valable, ie une ellipse : cf 1ere loi de Kepler)


vous pouvez modifier l'epsilon avec les touches + (multiplier pas 10) et - (diviser par 10)

Source / Exemple :


//----------------------------------------------------------
BOOL _CalcRungeKuttaODE(P_ODE_RK rk)
{
double h,errMax;
int    dim;

dim   = rk->dim;
h     = rk->h;

while(TRUE)
  {
  int     i;
  double  p;

  // calcul de la methode de Runge-Kutta
  evalRK(dim,rk->y,rk->yErr,rk->x,h,rk->derivFunc,rk->param,rk->yTmp,
         rk->k1,rk->k2,rk->k3,rk->k4,rk->k5,rk->k6);
  // calcul de l'erreur
  errMax = 0.;
  for(i=0;i<dim;i++)
    {
    double err;
    err = fabs((rk->yErr->coord[i])/(rk->yScale->coord[i]));
    if(errMax < err)
      {
      errMax = err;
      }
    }
  // si le pas est correct, on a fini notre boulot
  if(errMax <= 1.)
    {
    break;
    }
  // modification du pas, qui est trop grand pour l'erreur maximale voulue
  p  = 0.90*pow(errMax,-0.25);
  h *= max(0.1,p);
  if(CmpDouble0(h))
    {
    PopErrorMacro();
    return FALSE;
    }
  }

// on avance <x>
rk->x += h;
// on ecrit le nouveau pas
rk->h  = h * (errMax > 1.89e-4 ? 0.9*pow(errMax,-0.2) : 5.0);
// on ecrit les resultats
CopyVectND(rk->y,rk->yTmp);

PopErrorMacro();
return TRUE;
} // _CalcRungeKuttaODE()

Codes Sources

A voir également

Vous n'êtes pas encore membre ?

inscrivez-vous, c'est gratuit et ça prend moins d'une minute !

Les membres obtiennent plus de réponses que les utilisateurs anonymes.

Le fait d'être membre vous permet d'avoir un suivi détaillé de vos demandes et codes sources.

Le fait d'être membre vous permet d'avoir des options supplémentaires.