Boucle while infinie après passage de C++ à python

**psycopatate35** · 02/01/2022, 21h08

Bonjour et bonne année à tous !

J'ai récemment tenté de convertir en python un programme créé par un collègue en C++ (ci-après le début) afin de le modifier un peu et de le rendre plus automatique.

Code C++ :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
 
void pause();
 
 
int main(int argc, char *argv[])
 
{
   int Nmax = 4000;
   int I = 0;
   int J = 0;
   int N = 0;
   long double Psif[8100];
   long double PsiTotal[8100];
   long double Rqd[8100];
   long double Pota[8100];
   long double Dtau;
   long double Csi;
   long double Tau = 0;
   long double Tauf = 0;
   long double u;
   long double Pi = 3.14159265358979323846264338327950288419716939937510;
   long double Sum;
   long double Ct1;
   long double Ct2;
   long double Ct3;
   long double Ct4;
   long double Ct5;
   long double Lambda;
   long double Alpha = 0.5;
   long double PsiProv = 0;
   long double CsiPrime;
   long double Beta = 0; // Zero pour lest tests
   long double Gamma = 0; // Zero pour les tests
   long double Rau;
   long double CsiP;
   long double PsiP;
   long double CsiM;
   long double PsiM;
   double Ei=-1.0;
   double Ef=1.0;
   double E0=0.39;
   double S=0.0706858;
   S=S/10000;
   double C0=0.001;
   C0=C0*1000;
   double T=294.65;
   double D=0.0001;
   D=D/10000;
   double Ru=165;
   double Cd=0.000017;
   double ks=20;
   ks=ks/100;
   double v=0.05;
   long double FRT;
   long double FSC0;
   long double DFRT;
   long double F = 96485.3329;
   long double R = 8.3144621;
   long double Ep1;
   long double Ep2;
   long double Ip1;
   long double Ip2;
    FRT = F / R / T;
    FSC0 = F * S * C0;
    DFRT = D * FRT;
 
    //parameters are defined for oxidation
    u = -(Ei - E0) * FRT;
    CsiM = (Ef - E0) * FRT;
    Lambda = ks / sqrt (v * DFRT);
    Rau = FRT * Ru * FSC0 * sqrt (v *DFRT);
    Gamma = Cd * sqrt (v) / FSC0 /sqrt (DFRT);
 
    printf (" Ei = %f ; Ef = %f ; S = %e ;  C0 = %e  \n", Ei, Ef, S, C0);
    printf (" T = %f ; D = %e ; Ru = %e ;  Cd = %e  \n", T, D, Ru, Cd);
    printf (" ks = %f ; v = %f \n", ks, v);
 
 
    Psif[0] = 0;
    Psif[1] = 0;
    Dtau = (CsiM + u) / Nmax; // calculation of the time step.
    Tau = Dtau;
    Csi = -u;
    Beta = Rau * Gamma;
    Ct3 = Beta/Dtau;
 
    printf (" Lambda = %LF ; Rau = %Lf ; Gamma = %Lf ; Dtau = %Lf \n", Lambda, Rau, Gamma, Dtau);
 
    for (I = 1 ; I <= 2.0 * Nmax + 1; I++)
 
        {
            Rqd[I] = (powl ((I + 1.0) , 1.5)) - (2.0 * powl (I, 1.5)) + (powl ((I - 1) , 1.5));
        }
 
 // Forward Scan
 
 
    for (I = 1 ; I <= Nmax ; I++)
        {
            Ct1 = expl (-Alpha * Csi) / Lambda;
            Ct2 = 4.0 / 3.0 * sqrtl (Dtau / Pi) * (1.0 + expl (-Csi));
            Pota[I] = Csi;
 
            Sum = 0;
 
            for (N = 1 ; N < I ; N++)
 
                {
                    Sum = Sum + (Psif[N] * (1.0 + Ct3) - Psif[N - 1] * Ct3) * Rqd[I - N];
                }
 
            N = 0;
 
            Psif[I] = Psif[I - 1];
 
            do
 
                {
                    Psif[I] = (Psif[I] + 100.0 * PsiProv) /101.0; // Il faut calmer la convergence
                    PsiProv = Psif[I];
                    CsiPrime = Csi - Rau * Psif[I] - Beta * (1.0 - expl (-Tau/Beta));
                    Ct1 = expl (-Alpha * CsiPrime) / Lambda;
                    Ct2 = 4.0 / 3.0 * sqrtl (Dtau / Pi) * (1.0 + expl (-CsiPrime));
 
                    Psif[I] = (1.0 - Sum * Ct2 + (Ct1 * Ct3 + Ct2 * Ct3) * Psif[I - 1] ) / (Ct1 + Ct2 + Ct1 * Ct3 + Ct2 * Ct3);
 
                    N++;
 
                }
 
            while ( fabs ( Psif[I] - PsiProv ) > 0.0001 );
 
            PsiTotal[I] = Psif[I] + Gamma * (1.0 - expl (-Tau / Beta));
 
            if ((I % 200) == 0)
 
            {
 
                printf (" I= %d Csi= %Lf Psif= %Lf PsiTotal= %Lf N= %d \n", I, Csi, Psif[I], PsiTotal[I], N);
 
            }
 
            Csi = Csi + Dtau;
            Tau = Tau + Dtau;
        }
 
    Tauf = Tau - Dtau;
}

Cependant, je me heurte à un problème de boucle infinie dans mon programme python et qui n'existe pas en C++.

Code Python :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
import math
import numpy as np
import statistics
 
Nmax = int(8000)
I = int(0)
J = int(0)
N = int(0)
 
Psif = [0,0]
PsiTotal = [0,0]
Rqd = [0,0]
Pota = [0,0]
for i in range(3*Nmax):
    Psif.append(np.longdouble(0e0))
    PsiTotal.append(np.longdouble(0e0))
    Rqd.append(np.longdouble(0e0))
    Pota.append(np.longdouble(0e0))
 
Tau = float(0)
Tauf = 0
PsiProv = np.longdouble(0e0)
F = float(96485.3329)
R = float(8.3144621)
 
Ei=-1
Ef=1
E0=0.39
S=0.0706858
S=S/10000
C0=0.001
C0=C0*1000
T=294.65
D=0.0001
D=D/10000
Ru=165
Cd=0.000017
ks=20
ks=ks/100
Alpha = 0.5
v=0.05
 
FRT=F/R/T
FSC0=F*S*C0
DFRT=D*FRT
u=-(Ei-E0)*FRT
CsiM=(Ef-E0)*FRT
Lambda=ks/math.sqrt(v*DFRT)
Rau=FRT*Ru*FSC0*math.sqrt(v*DFRT)
Gamma=Cd*math.sqrt(v)/FSC0/math.sqrt(DFRT)
 
"calculation of step time"
Dtau = (CsiM+u)/Nmax
Tau=Dtau
Csi=-u
Beta=Rau*Gamma
Ct3=Beta/Dtau
 
print(f"Lambda = {Lambda} ; Rau = {Rau} ; Gamma = {Gamma} ; Dtau = {Dtau}")
 
for i in range(2*Nmax+1):
    Rqd[i]=math.pow((i+2),1.5)-(2*math.pow(i+1,1.5))+(math.pow((i),1.5))
 
"Forward Scan"
 
for i in range(Nmax+1):
    if i==0:
        pass
    else:
        Ct1=math.exp(-Alpha*Csi/Lambda)
        Ct2=np.longdouble(4/3*math.sqrt(Dtau/math.pi)*(1+math.exp(-Csi)))
        Pota[i]=Csi
        Sum=0
        for N in range(i):
            if N==0:
                pass
            else :
                Sum=Sum+(Psif[N]*(1+Ct3)-Psif[N-1]*Ct3)*Rqd[i-N]
        N=int(0)
        Psif[i]=Psif[i-1]
        while True:
            Psif
            Psif[i] = np.longdouble((Psif[i] + 100*PsiProv)/101)
            PsifProv = np.longdouble(Psif[i])
            CsiPrime = Csi - Rau*Psif[i] -Beta*(1+math.exp(-Tau/Beta))
            Ct1 = math.exp(-Alpha*CsiPrime)/Lambda
            Ct2 = np.longdouble(4/3*math.sqrt(Dtau/math.pi)*(1+math.exp(-CsiPrime)))
            Psif[i] = np.longdouble((1-Sum*Ct2+(Ct1*Ct3+Ct2*Ct3)*Psif[i-1])/(Ct1+Ct2+Ct1*Ct3+Ct2*Ct3))
            N=N+1
            if abs(Psif[i]-PsiProv)<=0.0001:
                break
        PsiTotal[i]=Psif[i]+Gamma*(1-math.exp(-Tau/Beta))
        if i%200==0 and i!=0:
            print(f"I= {i}, Csi= {Csi}, Psif= {Psif[i]}, PsiTotal= {PsiTotal[i]}, N= {N}")
        Csi=Csi+Dtau
        Tau=Tau+Dtau
Tauf=Tau-Dtau

Quelqu'un aurait la solution svp ?

**wiztricks** · 02/01/2022, 23h23

Envoyé par psycopatate35

Quelqu'un aurait la solution svp ?

Pour trouver la solution, il faut essayer de comprendre le problème.

Au pif, il y a un while True dont on ne sort qu'à la condition abs(Psif[i]-PsiProv)<=0.0001. Déjà vérifier si cette condition est réalisée sinon voir ce qui a été loupé dans la traduction du C++...

C'est juste le boulot qu'on fait pour mettre au point un programme.... et je ne vais pas le faire pour vous.

- W

**psycopatate35** · 03/01/2022, 08h58

C'est effectivement tout le problème !

La boucle Do While de C++ fonctionne parfaitement et arrive à sortir pour l'ensemble de la boucle for définie par I = [1;4000].

En python, ça fonctionne également mais seulement jusqu'à I = 2310. Jusque la, on obtient sensiblement les mêmes réponses qu'en C++. Par contre, dès que I = 2311, la boucle While n'arrive plus à atteindre la condition de sortie abs(Psif[i]-PsiProv)<=0.0001, qui pourtant ce fait bien en C++.

La programmation n'étant pas mon domaine principal d'expertise, j'ai bien tenté de comprendre avec des proches ce qui coince mais sans succès.

**wiztricks** · 03/01/2022, 09h27

C'est peut être un défaut côté précision.
Là ou il y a du boulot c'est dans la vérification de la cohérence entre types de données C++, leur traduction en types numpy et s'assurer que les opérations dans le code ne le dégradent pas.

- W

**Sve@r** · 03/01/2022, 11h13

Bonjour

Envoyé par psycopatate35

En python, ça fonctionne également mais seulement jusqu'à I = 2310. Jusque la, on obtient sensiblement les mêmes réponses qu'en C++. Par contre, dès que I = 2311, la boucle While n'arrive plus à atteindre la condition de sortie abs(Psif[i]-PsiProv)<=0.0001, qui pourtant ce fait bien en C++.

Comme le dit wiztricks, il peut y avoir un souci de précision. Souci qui est d'ailleurs parfaitement connu et qui n'est pas spécifique à Python mais dû à la façon dont un ordinateur convertit un décimal en binaire (tape par exemple 0.6 / 3 puis 0.2 * 3 dans une console Python et tu comprendras).
Bon ce qui m'étonne c'est que ce souci étant aussi présent donc en C++, ça n'aurait pas dû fonctionner en C++ mais peut-être qu'en C++ la profondeur de la précision n'est pas la même et que cette différence fait justement la différence.

Ce que je peux te proposer en retour c'est d'utiliser le module decimal, lequel travaille chiffre par chiffre et non pas en binaire. C'est plus long mais c'est garanti sans erreur. Le même exemple que tout à l'heure en utilisant ce module: decimal.Decimal("0.6") / 3 ou decimal.Decimal("0.2") * 3.

**jurassic pork** · 03/01/2022, 13h30

hello,

Envoyé par psycopatate35

Bonjour et bonne année à tous !

J'ai récemment tenté de convertir en python un programme créé par un collègue en C++ (ci-après le début) afin de le modifier un peu et de le rendre plus automatique.

Code C++ :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
 
void pause();
 
 
int main(int argc, char *argv[])
 
{
   int Nmax = 4000;
   int I = 0;
   int J = 0;
   int N = 0;
   long double Psif[8100];
   long double PsiTotal[8100];
   long double Rqd[8100];
   long double Pota[8100];
   long double Dtau;
   long double Csi;
   long double Tau = 0;
   long double Tauf = 0;
   long double u;
   long double Pi = 3.14159265358979323846264338327950288419716939937510;
   long double Sum;
   long double Ct1;
   long double Ct2;
   long double Ct3;
   long double Ct4;
   long double Ct5;
   long double Lambda;
   long double Alpha = 0.5;
   long double PsiProv = 0;
   long double CsiPrime;
   long double Beta = 0; // Zero pour lest tests
   long double Gamma = 0; // Zero pour les tests
   long double Rau;
   long double CsiP;
   long double PsiP;
   long double CsiM;
   long double PsiM;
   double Ei=-1.0;
   double Ef=1.0;
   double E0=0.39;
   double S=0.0706858;
   S=S/10000;
   double C0=0.001;
   C0=C0*1000;
   double T=294.65;
   double D=0.0001;
   D=D/10000;
   double Ru=165;
   double Cd=0.000017;
   double ks=20;
   ks=ks/100;
   double v=0.05;
   long double FRT;
   long double FSC0;
   long double DFRT;
   long double F = 96485.3329;
   long double R = 8.3144621;
   long double Ep1;
   long double Ep2;
   long double Ip1;
   long double Ip2;
    FRT = F / R / T;
    FSC0 = F * S * C0;
    DFRT = D * FRT;
 
    //parameters are defined for oxidation
    u = -(Ei - E0) * FRT;
    CsiM = (Ef - E0) * FRT;
    Lambda = ks / sqrt (v * DFRT);
    Rau = FRT * Ru * FSC0 * sqrt (v *DFRT);
    Gamma = Cd * sqrt (v) / FSC0 /sqrt (DFRT);
 
    printf (" Ei = %f ; Ef = %f ; S = %e ;  C0 = %e  \n", Ei, Ef, S, C0);
    printf (" T = %f ; D = %e ; Ru = %e ;  Cd = %e  \n", T, D, Ru, Cd);
    printf (" ks = %f ; v = %f \n", ks, v);
 
 
    Psif[0] = 0;
    Psif[1] = 0;
    Dtau = (CsiM + u) / Nmax; // calculation of the time step.
    Tau = Dtau;
    Csi = -u;
    Beta = Rau * Gamma;
    Ct3 = Beta/Dtau;
 
    printf (" Lambda = %LF ; Rau = %Lf ; Gamma = %Lf ; Dtau = %Lf \n", Lambda, Rau, Gamma, Dtau);
 
    for (I = 1 ; I <= 2.0 * Nmax + 1; I++)
 
        {
            Rqd[I] = (powl ((I + 1.0) , 1.5)) - (2.0 * powl (I, 1.5)) + (powl ((I - 1) , 1.5));
        }
 
 // Forward Scan
 
 
    for (I = 1 ; I <= Nmax ; I++)
        {
            Ct1 = expl (-Alpha * Csi) / Lambda;
            Ct2 = 4.0 / 3.0 * sqrtl (Dtau / Pi) * (1.0 + expl (-Csi));
            Pota[I] = Csi;
 
            Sum = 0;
 
            for (N = 1 ; N < I ; N++)
 
                {
                    Sum = Sum + (Psif[N] * (1.0 + Ct3) - Psif[N - 1] * Ct3) * Rqd[I - N];
                }
 
            N = 0;
 
            Psif[I] = Psif[I - 1];
 
            do
 
                {
                    Psif[I] = (Psif[I] + 100.0 * PsiProv) /101.0; // Il faut calmer la convergence
                    PsiProv = Psif[I];
                    CsiPrime = Csi - Rau * Psif[I] - Beta * (1.0 - expl (-Tau/Beta));
                    Ct1 = expl (-Alpha * CsiPrime) / Lambda;
                    Ct2 = 4.0 / 3.0 * sqrtl (Dtau / Pi) * (1.0 + expl (-CsiPrime));
 
                    Psif[I] = (1.0 - Sum * Ct2 + (Ct1 * Ct3 + Ct2 * Ct3) * Psif[I - 1] ) / (Ct1 + Ct2 + Ct1 * Ct3 + Ct2 * Ct3);
 
                    N++;
 
                }
 
            while ( fabs ( Psif[I] - PsiProv ) > 0.0001 );
 
            PsiTotal[I] = Psif[I] + Gamma * (1.0 - expl (-Tau / Beta));
 
            if ((I % 200) == 0)
 
            {
 
                printf (" I= %d Csi= %Lf Psif= %Lf PsiTotal= %Lf N= %d \n", I, Csi, Psif[I], PsiTotal[I], N);
 
            }
 
            Csi = Csi + Dtau;
            Tau = Tau + Dtau;
        }
 
    Tauf = Tau - Dtau;
}

Cependant, je me heurte à un problème de boucle infinie dans mon programme python et qui n'existe pas en C++.

Code Python :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
import math
import numpy as np
import statistics
 
Nmax = int(8000)
I = int(0)
J = int(0)
N = int(0)
 
Psif = [0,0]
PsiTotal = [0,0]
Rqd = [0,0]
Pota = [0,0]
for i in range(3*Nmax):
    Psif.append(np.longdouble(0e0))
    PsiTotal.append(np.longdouble(0e0))
    Rqd.append(np.longdouble(0e0))
    Pota.append(np.longdouble(0e0))
 
Tau = float(0)
Tauf = 0
PsiProv = np.longdouble(0e0)
F = float(96485.3329)
R = float(8.3144621)
 
Ei=-1
Ef=1
E0=0.39
S=0.0706858
S=S/10000
C0=0.001
C0=C0*1000
T=294.65
D=0.0001
D=D/10000
Ru=165
Cd=0.000017
ks=20
ks=ks/100
Alpha = 0.5
v=0.05
 
FRT=F/R/T
FSC0=F*S*C0
DFRT=D*FRT
u=-(Ei-E0)*FRT
CsiM=(Ef-E0)*FRT
Lambda=ks/math.sqrt(v*DFRT)
Rau=FRT*Ru*FSC0*math.sqrt(v*DFRT)
Gamma=Cd*math.sqrt(v)/FSC0/math.sqrt(DFRT)
 
"calculation of step time"
Dtau = (CsiM+u)/Nmax
Tau=Dtau
Csi=-u
Beta=Rau*Gamma
Ct3=Beta/Dtau
 
print(f"Lambda = {Lambda} ; Rau = {Rau} ; Gamma = {Gamma} ; Dtau = {Dtau}")
 
for i in range(2*Nmax+1):
    Rqd[i]=math.pow((i+2),1.5)-(2*math.pow(i+1,1.5))+(math.pow((i),1.5))
 
"Forward Scan"
 
for i in range(Nmax+1):
    if i==0:
        pass
    else:
        Ct1=math.exp(-Alpha*Csi/Lambda)
        Ct2=np.longdouble(4/3*math.sqrt(Dtau/math.pi)*(1+math.exp(-Csi)))
        Pota[i]=Csi
        Sum=0
        for N in range(i):
            if N==0:
                pass
            else :
                Sum=Sum+(Psif[N]*(1+Ct3)-Psif[N-1]*Ct3)*Rqd[i-N]
        N=int(0)
        Psif[i]=Psif[i-1]
        while True:
            Psif
            Psif[i] = np.longdouble((Psif[i] + 100*PsiProv)/101)
            PsifProv = np.longdouble(Psif[i])
            CsiPrime = Csi - Rau*Psif[i] -Beta*(1+math.exp(-Tau/Beta))
            Ct1 = math.exp(-Alpha*CsiPrime)/Lambda
            Ct2 = np.longdouble(4/3*math.sqrt(Dtau/math.pi)*(1+math.exp(-CsiPrime)))
            Psif[i] = np.longdouble((1-Sum*Ct2+(Ct1*Ct3+Ct2*Ct3)*Psif[i-1])/(Ct1+Ct2+Ct1*Ct3+Ct2*Ct3))
            N=N+1
            if abs(Psif[i]-PsiProv)<=0.0001:
                break
        PsiTotal[i]=Psif[i]+Gamma*(1-math.exp(-Tau/Beta))
        if i%200==0 and i!=0:
            print(f"I= {i}, Csi= {Csi}, Psif= {Psif[i]}, PsiTotal= {PsiTotal[i]}, N= {N}")
        Csi=Csi+Dtau
        Tau=Tau+Dtau
Tauf=Tau-Dtau

Quelqu'un aurait la solution svp ?

le souci c'est que math.pi de python n'a pas la précision voulue ((~14 digits). Pour avoir une meilleure précision on peut utiliser le module mpmath (présent dans pypi).
1 - On définit la précision pour les calculs ( ex:50 ):

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
from mpmath import mp
mp.prec = 50

2 - dans le code on remplace math.pi par mp.pi

Ami calmant, J.P

**marco056** · 03/01/2022, 13h58

Envoyé par jurassic pork

hello,

le souci c'est que math.pi de python n'a pas la précision voulue ((~14 digits). Pour avoir une meilleure précision on peut utiliser le module mpmath (présent dans pypi).
1 - On définit la précision pour les calculs ( ex:50 ):

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
from mpmath import mp
mp.prec = 50

2 - dans le code on remplace math.pi par mp.pi

Ami calmant, J.P

A priori, cela ne fonctionne pas mieux avec tes modifs.

**jurassic pork** · 03/01/2022, 14h00

heu c'est bizarre chez moi j'arrive à aller plus loin dans la boucle :

D:\Dev\PycharmProjects\pythonProject\venv\Scripts\python.exe D:/Dev/PycharmProjects/pythonProject/binaryread/testC++.py
Lambda = 1425.2289828701896 ; Rau = 0.6219328089590949 ; Gamma = 0.008881372761047751 ; Dtau = 0.009846019801588628
I= 200, Csi= -52.78451215631721, Psif= 1.2543724463752129e-23, PsiTotal= 0.008881372761047751, N= 1
I= 400, Csi= -50.81530819600007, Psif= 8.915267114499358e-23, PsiTotal= 0.008881372761047751, N= 1
I= 600, Csi= -48.84610423568292, Psif= 6.3847702554760075e-22, PsiTotal= 0.008881372761047751, N= 1
I= 800, Csi= -46.87690027536578, Psif= 4.574495692355078e-21, PsiTotal= 0.008881372761047751, N= 1
I= 1000, Csi= -44.90769631504863, Psif= 3.277601068629375e-20, PsiTotal= 0.008881372761047751, N= 1
I= 1200, Csi= -42.93849235473149, Psif= 2.348390683187784e-19, PsiTotal= 0.008881372761047751, N= 1
I= 1400, Csi= -40.96928839441434, Psif= 1.6826149583944407e-18, PsiTotal= 0.008881372761047753, N= 1
I= 1600, Csi= -39.0000844340972, Psif= 1.205588668328873e-17, PsiTotal= 0.008881372761047763, N= 1
I= 1800, Csi= -37.03088047378005, Psif= 8.638007379958643e-17, PsiTotal= 0.008881372761047838, N= 1
I= 2000, Csi= -35.061676513462906, Psif= 6.189106905936886e-16, PsiTotal= 0.00888137276104837, N= 1
I= 2200, Csi= -33.09247255314576, Psif= 4.434476912244294e-15, PsiTotal= 0.008881372761052185, N= 1
I= 2400, Csi= -31.123268592828296, Psif= 3.177289677331258e-14, PsiTotal= 0.008881372761079524, N= 1
I= 2600, Csi= -29.15406463251044, Psif= 2.276518717263641e-13, PsiTotal= 0.008881372761275402, N= 1
I= 2800, Csi= -27.184860672192585, Psif= 1.6311189705812705e-12, PsiTotal= 0.00888137276267887, N= 1
I= 3000, Csi= -25.21565671187473, Psif= 1.168691947401768e-11, PsiTotal= 0.00888137277273467, N= 1
I= 3200, Csi= -23.246452751556873, Psif= 8.373643419428793e-11, PsiTotal= 0.008881372844784186, N= 1
I= 3400, Csi= -21.277248791239018, Psif= 5.999690819107269e-10, PsiTotal= 0.008881373361016832, N= 1
I= 3600, Csi= -19.308044830921162, Psif= 4.298760663500684e-09, PsiTotal= 0.008881377059808414, N= 1
I= 3800, Csi= -17.338840870603306, Psif= 3.080048970131414e-08, PsiTotal= 0.008881403561537453, N= 1
I= 4000, Csi= -15.36963691028545, Psif= 2.20684488234085e-07, PsiTotal= 0.008881593445535985, N= 1
I= 4200, Csi= -13.400432949967595, Psif= 1.5811938092893269e-06, PsiTotal= 0.00888295395485704, N= 1

**marco056** · 03/01/2022, 14h21

Je n'y connais rien en C++ mais a priori :
- erreur de signe à la ligne 85
- modifier les boucles par exemple :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

for i in range(1,Nmax+1):

virer la condition si i = 0 ensuite
- es-tu sûr de cela ? Qu'est-ce que cela veut dire pour toi :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
for i in range(3*Nmax):
    Psif.append(np.longdouble(0e0))
    PsiTotal.append(np.longdouble(0e0))
    Rqd.append(np.longdouble(0e0))
    Pota.append(np.longdouble(0e0))

**marco056** · 03/01/2022, 14h22

Envoyé par jurassic pork

heu c'est bizarre chez moi j'arrive à aller plus loin dans la boucle :

Le programme se termine chez toi ?

Dans tous les cas, il faudrait connaître le résultat attendu...

**jurassic pork** · 03/01/2022, 14h29

Envoyé par marco056

Le programme se termine chez toi ?

Dans tous les cas, il faudrait connaître le résultat attendu...

ça coince à I=4600

**marco056** · 03/01/2022, 14h40

Envoyé par jurassic pork

ça coince à I=4600

Oui, pareil chez moi, donc pas trop satisfaisant sans doute.

**jurassic pork** · 03/01/2022, 15h25

heu,

finalement pareil en utilisant math.pi ça coince à I=4600
Le truc louche c'est que dans son programme c++ Nmax vaut 4000 et dans celui en python 8000
Un indice : dans le programme python avec Nmax à 4000 ça coince à 2200 , avec 8000 ça coince à 4600 alors que c'est I qui est affiché
avec Nmax à 16000 on atteint 9200 :

I= 8000, Csi= -15.36471390039043, Psif= 2.186263783958073e-07, PsiTotal= 0.008881591387426146, N= 1I= 8200, Csi= -14.380111920231503, Psif= 5.85206892410364e-07, PsiTotal= 0.008881957967940161, N= 1
I= 8400, Csi= -13.395509940072575, Psif= 1.5664475629698147e-06, PsiTotal= 0.00888293920861072, N= 1
I= 8600, Csi= -12.410907959913647, Psif= 4.192964449614133e-06, PsiTotal= 0.008885565725497365, N= 1
I= 8800, Csi= -11.42630597975472, Psif= 1.1223379580582602e-05, PsiTotal= 0.008892596140628335, N= 1
I= 9000, Csi= -10.441703999595791, Psif= 3.00412985804255e-05, PsiTotal= 0.008911414059628177, N= 1
I= 9200, Csi= -9.457102019436864, Psif= 8.040706714288705e-05, PsiTotal= 0.008961779828190639, N= 1

**wiztricks** · 03/01/2022, 16h25

Il y a un typo(*).
Ligne 84: PsifProv = np.longdouble(Psif[i]) devrait être PsiProv = np.longdouble(Psif[i]).
(*) façon gentille de dire que quelqu'un s'est endormi sur son clavier.

- W

**jurassic pork** · 03/01/2022, 16h33

Envoyé par wiztricks

Il y a un typo(*).
Ligne 84: PsifProv = np.longdouble(Psif[i]) devrait être PsiProv = np.longdouble(Psif[i]).
(*) façon gentille de dire que quelqu'un s'est endormi sur son clavier.

- W

bien vu !

avec NMax à 4000 :

I= 3400, Csi= 12.189372514369143, Psif= 0.16470310961426204, PsiTotal= 0.1735844823753098, N= 89
I= 3600, Csi= 16.12778043500451, Psif= 0.14222592690530977, PsiTotal= 0.15110729966635753, N= 65
I= 3800, Csi= 20.06618835564022, Psif= 0.12708176940682908, PsiTotal= 0.13596314216787683, N= 50
I= 4000, Csi= 24.004596276275933, Psif= 0.11596149504673671, PsiTotal= 0.12484286780778446, N= 40
Process finished with exit code 0

**marco056** · 03/01/2022, 19h20

Ah oui, bien vu.
Avec un peu de nettoyage (j'ai dit un peu) :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
import math
import numpy as np
 
Nmax = 4000
I,J,N = 0,0,0
Psif, PsiTotal, Rqd, Pota = [0,0], [0,0], [0,0], [0,0]
 
for i in range(3*Nmax):
    Psif.append(np.longdouble(0e0))
    PsiTotal.append(np.longdouble(0e0))
    Rqd.append(np.longdouble(0e0))
    Pota.append(np.longdouble(0e0))
 
Tau, Tauf = 0,0
PsiProv = np.longdouble(0e0)
F = 96485.3329
R = 8.3144621
 
Ei, Ef =-1, 1
E0=0.39
S=0.0706858
S=S/10000
C0=0.001
C0=C0*1000
T=294.65
D=1e-4
D=D/1e5
Ru=165
Cd=0.000017
ks=20
ks=ks/100
Alpha = 0.5
v=0.05
 
FRT=F/R/T
FSC0=F*S*C0
DFRT=D*FRT
u=-(Ei-E0)*FRT
CsiM=(Ef-E0)*FRT
Lambda=ks/math.sqrt(v*DFRT)
Rau=FRT*Ru*FSC0*math.sqrt(v*DFRT)
Gamma=Cd*math.sqrt(v)/FSC0/math.sqrt(DFRT)
 
"calculation of step time"
Dtau = (CsiM+u)/Nmax
Tau=Dtau
Csi=-u
Beta=Rau*Gamma
Ct3=Beta/Dtau
 
print(f"Lambda = {Lambda} ; Rau = {Rau} ; Gamma = {Gamma} ; Dtau = {Dtau}")
 
for i in range(1,2*Nmax+2):
    Rqd[i]=math.pow((i+2),1.5)-(2*math.pow(i+1,1.5))+(math.pow((i),1.5))
 
"Forward Scan"
 
for i in range(1,Nmax+1):
    Ct1=math.exp(-Alpha*Csi/Lambda)
    Ct2=np.longdouble(4/3*math.sqrt(Dtau/math.pi)*(1+math.exp(-Csi)))
    Pota[i]=Csi
    Sum=0
    for N in range(1,i):
        Sum=Sum+(Psif[N]*(1+Ct3)-Psif[N-1]*Ct3)*Rqd[i-N]
    N = 0
    Psif[i]=Psif[i-1]
    while True:
        Psif
        Psif[i] = np.longdouble((Psif[i] + 100*PsiProv)/101)
        PsiProv = np.longdouble(Psif[i])
        CsiPrime = Csi - Rau*Psif[i] -Beta*(1-math.exp(-Tau/Beta))
        Ct1 = math.exp(-Alpha*CsiPrime)/Lambda
        Ct2 = np.longdouble(4/3*math.sqrt(Dtau/math.pi)*(1+math.exp(-CsiPrime)))
        Psif[i] = np.longdouble((1-Sum*Ct2+(Ct1*Ct3+Ct2*Ct3)*Psif[i-1])/(Ct1+Ct2+Ct1*Ct3+Ct2*Ct3))
        N=N+1
        if abs(Psif[i]-PsiProv)<=0.0001:
            break
    PsiTotal[i]=Psif[i]+Gamma*(1-math.exp(-Tau/Beta))
    if i%200==0 and i!=0:
        print(f"I= {i}, Csi= {Csi}, Psif= {Psif[i]}, PsiTotal= {PsiTotal[i]}, N= {N}")
    Csi=Csi+Dtau
    Tau=Tau+Dtau
Tauf=Tau-Dtau

**psycopatate35** · 03/01/2022, 20h26

Envoyé par wiztricks

Il y a un typo(*).
Ligne 84: PsifProv = np.longdouble(Psif[i]) devrait être PsiProv = np.longdouble(Psif[i]).
(*) façon gentille de dire que quelqu'un s'est endormi sur son clavier.

- W

Oh mon dieu !! Merci !!

Merci à tous de vous être penché dessus ! C'est super cool d'avoir pu compter sur vous tous !!

**marco056** · 03/01/2022, 20h37

Envoyé par wiztricks

Il y a un typo(*).
Ligne 84: PsifProv = np.longdouble(Psif[i]) devrait être PsiProv = np.longdouble(Psif[i]).
(*) façon gentille de dire que quelqu'un s'est endormi sur son clavier.

- W

C'est marrant quand même : 2 lignes à suivre.

**lg_53** · 04/01/2022, 09h49

Envoyé par psycopatate35

La programmation n'étant pas mon domaine principal d'expertise

Dans ce cas pourquoi cherchez vous à réécrire en python ce qui fonctionne très bien en C++ ?

**psycopatate35** · 04/01/2022, 10h01

Envoyé par lg_53

Dans ce cas pourquoi cherchez vous à réécrire en python ce qui fonctionne très bien en C++ ?

Je cherche à modifier des choses dans le programme, notamment l'automatiser, et créer une petite interface graphique, et python est l'un des langages les plus simple pour cela apparemment.

Boucle while infinie après passage de C++ à python [Python 3.X]

Python

Discussions similaires

Partager

Partager