Problème sur une boucle en C++

Résolu/Fermé
turbochris Messages postés 4 Date d'inscription dimanche 10 mars 2013 Statut Membre Dernière intervention 10 mars 2013 - 10 mars 2013 à 10:55
turbochris Messages postés 4 Date d'inscription dimanche 10 mars 2013 Statut Membre Dernière intervention 10 mars 2013 - 10 mars 2013 à 13:49
Bonjour,

je sollicite l'aide des membres du site. en effet j'ai un programme en c++ dans lequel je fais varier un paramètre r dans une boucle. dans cette boucle, il y'a une serie d'instruction. mon pb est que lorsque je fais varier r de 0.93 à 0.90, les résultats change lorsque je commence de 0.94 à 0.90
je précise que parmi les instructions, j'appel aussi des fonctions. que dois je faire?

6 réponses

KX Messages postés 16752 Date d'inscription samedi 31 mai 2008 Statut Modérateur Dernière intervention 31 août 2024 3 019
10 mars 2013 à 11:01
Sur ce site, il n'y a pas de devin, il faut que tu nous donnes le code qui pose problème parce que tes explications ne sont pas suffisantes pour t'aider...
0
turbochris Messages postés 4 Date d'inscription dimanche 10 mars 2013 Statut Membre Dernière intervention 10 mars 2013
10 mars 2013 à 12:15
#include <cstdlib>
#include <ctime>
#include <iomanip>
#include <iostream>
#include <math.h>
#include <fstream>
#include <string.h>
#include <stdio.h>
#include <conio.h>
#include <stdlib.h>
#include <fcntl.h>

#define pi 3.1415926535897932384626433832795
#define mu 1.0
#define sig 1.0
#define mu2 0.15
#define ki 1.0
#define alfa2 0.4
#define kip 1.0
#define k0 0.01
#define lam 0.1
#define lamp 0.0

#define dlam 0.05
#define dv 0.025
#define dr 0.01

#define vfin 0.25
//#define rfin 0.25
#define lamfin 4.0

#define kmax 39000
#define ktrans 70000
#define hh 0.05
#define dim 2
#define dip 2
#define val 1

const double c41[4] = {0.0, 0.5*hh, 0.0, 0.5*hh};
const double c42[4] = {0.5*hh, 0.5*hh, hh, 0.0};
const double c43[4] = {hh/6.0, hh/3.0, hh/3.0, hh/6.0};
//const int nbv = int(ceil(vfin/dv1));
const double alfa22=alfa2*alfa2;

int lk;
double v1, v2, r1, r2;
double vx[dip][kmax], vxp[dip][kmax];
double ver[dip][kmax];

// declaration of functions
double sqr(double u);
void initu(double vx[dip][kmax], double vxp[dip][kmax], double ver[dip][kmax]);
double f2(double xt, double xx[dim], double xp[dim]);
double f4(double xt, double xx[dim], double xp[dim]);
void update_rk4_x(long int nt, double &t, double xx[dim], double xp[dim], double vx[dip][kmax], double vxp[dip][kmax], double ver[dip][kmax]);
double Abs(double u);

using namespace std;
int main ()
{
srand(time(NULL));
//declaration of variable
double t, xx[dim], xp[dim];
double er[dim];
double kxx, lam0;
double rfin;
ofstream fichx ("synchro2a.dat");

if ( fichx.is_open() ) // opening the datas file
{
r1=-0.93; r2=r1; rfin=-0.90;
do
{
v1=0.05; v2=v1;
xx[0]=0.0; xp[0]=0.0;
xx[1]=0.01; xp[1]=0.05;
er[0] = xx[1]-xx[0];
er[1] = xp[1]-xp[0];
t=0.;
update_rk4_x(ktrans,t,xx,xp,vx,vxp,ver); // transient
initu(vx,vxp,ver);
update_rk4_x(kmax,t,xx,xp,vx,vxp,ver);
fichx << ver[0][kmax-1] << " " << ver[1][kmax-1] << "\n";
cout << " r1 = " << r1 << " r2 = " << r2 << " "<< ver[0][kmax-1]<<"\n";
r1 += dr;
r2 += dr;
} while (r1<rfin);
cin >> lk;
fichx.close();
} // if no error when opening fichx
return(0);
} // end of main program


double sqr(double u)
{
return(u*u);
}

void initu(double vx[dip][kmax], double vxp[dip][kmax], double ver[dip][kmax])
{
for (int i=0; i<dip; ++i)
{
for (long int k=0; k<kmax; ++k)
{
vx[i][k] = 0.0; vxp[i][k] = 0.0;
ver[i][k] = 0.0;
}
}
}

double f2(double xt, double xx[dim], double xp[dim])
{
double result, arg1, pterm, a, b;

arg1 = 2.*pi*(xx[0]-v1*xt);
a = -mu*sig*mu2*(xp[0]-v1)-mu*ki*xx[0]*alfa22-mu*kip*k0*xx[0]*sqr(xx[0]);
a -= mu*lam*(xx[0]-xx[1])+mu*lamp*(xx[0]-xx[1])*sqr(xx[0]-xx[1]);
b = -(mu*val*(1./(2.*pi))*(sqr(1.-r1*r1))*sin(arg1))/sqr(1.+r1*r1+2.*r1*cos(arg1));
pterm = a+b;

result = pterm;
return (result);
}

double f4(double xt, double xx[dim], double xp[dim])
{
double result, pterm, arg2, a, b;

arg2 = 2.*pi*(xx[1]-v2*xt);
a = -mu2*(xp[1]-v2)-xx[1]*alfa22-k0*xx[1]*sqr(xx[1]);
a -=lam*(xx[1]-xx[0])+lamp*(xx[1]-xx[0])*sqr(xx[1]-xx[0]);
b = -1.*val/(2.*pi)*(sqr(1.-r2*r2))*sin(arg2)/sqr(1.+r2*r2+2.*r2*cos(arg2));
pterm = a+b;

result = pterm;
return (result);
}

void update_rk4_x(long int nt, double & t, double xx[dim], double xp[dim], double vx[dip][kmax], double vxp[dip][kmax], double ver[dip][kmax])
{
double kr4p[dim], kr4x[dim], xnew[dim], xrec[dim], xpnew[dim], xprec[dim], er[dim], som[dim];

som[0] = 0.0;
som[1] = 0.0;

for (long int k=0; k<nt; ++k)
{
for (int i=0; i<dim; ++i)
{
xpnew[i]=xp[i]; xprec[i]=xp[i];
xnew[i]=xx[i]; xrec[i]=xx[i];
}
for (int i4=0; i4<4; ++i4)
{
t += c41[i4];

kr4x[0]=xpnew[0];
kr4x[1]=xpnew[1];
kr4p[0]=f2(t,xnew,xpnew);
kr4p[1]=f4(t,xnew,xpnew);
xnew[0] = xrec[0]+kr4x[0]*c42[i4]; xp[0] += kr4p[0]*c43[i4]; xnew[1] = xrec[1]+kr4x[1]*c42[i4]; xx[0] += kr4x[0]*c43[i4];
xpnew[0] = xprec[0]+kr4p[0]*c42[i4]; xp[1] += kr4p[1]*c43[i4]; xpnew[1] = xprec[1]+kr4p[1]*c42[i4]; xx[1] += kr4x[1]*c43[i4];
} // end of i4 loop
er[0] = Abs(xx[1]-xx[0]);
er[1] = Abs(xp[1]-xp[0]);
if (k<kmax)
for (int i=0; i<dip; ++i)
{
som[i] += er[i];
vx[i][k] = xx[i];
vxp[i][k] = xp[i];
ver[i][k] = som[i];
}
} // end of time loop
} // ******* end of updating x **********

double Abs(double u)
{
if (u>=0)
return(u);
if (u<0)
return(-u);
}
0
KX Messages postés 16752 Date d'inscription samedi 31 mai 2008 Statut Modérateur Dernière intervention 31 août 2024 3 019
10 mars 2013 à 12:17
Donc ta boucle c'est ça :

r1=-0.93; r2=r1; rfin=-0.90; 
do 
{ 
    ...
    r1 += dr; 
    r2 += dr; 
} while (r1<rfin);

Mais quel est le problème ? Ça a l'air de bien fonctionner...
0
turbochris Messages postés 4 Date d'inscription dimanche 10 mars 2013 Statut Membre Dernière intervention 10 mars 2013
10 mars 2013 à 12:47
quand tu initialise r à -0.94; les résultats affichés par ver[0][kmax-1] changent
0

Vous n’avez pas trouvé la réponse que vous recherchez ?

Posez votre question
KX Messages postés 16752 Date d'inscription samedi 31 mai 2008 Statut Modérateur Dernière intervention 31 août 2024 3 019
10 mars 2013 à 13:31
Il faut faire attention aux approximations sur les double !

Quand tu initialises r1=-0.93, la valeur "réelle" de r1 est -0.93000000000000005
Mais quand tu initialises r1=-0.94, la valeur de r1 est -0.93999999999999995 à quoi tu rajoutes 0.01, ce qui donne -0.92999999999999994.

La valeur de r1 n'est donc jamais vraiment égale à -0.93 et son approximation peut-être différente selon la manière dont elle est calculée (-0.94 + 0.01 ≠ -0.93)

Remarque : il en est de même pour r2.

Donc comme tes calculs utilisent intensivement r1 et r2 (dans update_rk4_x qui appelle souvent f2 et f4), les erreurs d'approximations s'accumulent et cela entraîne une bifurcation du résultat final.

Il faudrait donc revoir tout tes calculs pour diminuer l'impact de ces erreurs, par exemple en manipulant des types plus précis (comme long double).

Au passage, il serait bon de revoir ta structure de code en général, car l'utilisation de variables globales comme ça, c'est vraiment très moche !
0
turbochris Messages postés 4 Date d'inscription dimanche 10 mars 2013 Statut Membre Dernière intervention 10 mars 2013
10 mars 2013 à 13:49
merci pr ton aide
0