Programmation du modèle de Onsager

Version imprimable

1 pièce(s) jointe(s)

Bonjour à tous,

Je souhaite programmer le modèle de Onsager. Pour cela on me donne ces infos:
voir fichier joint :)

voila ce que j'ai fais:

Code:

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
#include<stdio.h>
#include<math.h>
 
#define pie 3.1419
 
 
float k1(float,float);
float integK1(float,float,float,float);
float e1(float,float);
float integE1(float,float,float,float);
 
 
int main()
{
   float km,q,K,K1q,E1q,CK,t;
 
   float inf = 0;
   float sup = pie/2;
   float N = 0.5;
   float pas = (sup-inf)/N;
 
 
  for(km=0.5; km<=5; km+=0.1)
  {
      K = 1./km; //km = 1/k
      q=(2*sinh(2*K))/(cosh(2*K)*cosh(2*K)); //q(K)
      K1q = integK1(pas,inf,sup,q); //E1(q)
      E1q = integE1(pas,inf,sup,q); //K1(q)
      //CK = ; c'est ce que je veux au finale
 
  }
printf("%f %f\n",K1q,E1q);
 
 
  return 0;
}
 
 
float k1(float q, float t)
{
 
      return 1/(sqrt(1-(q*q)*(sin(t)*sin(t))));//ma fonction de K1(q)
 
}
 
float integK1(float pas,float inf, float sup,float q)//integrale de k1
{
     int i=0;
     float K1=0;
 
     for(i=inf; i<sup; i+=pas)
     {
          K1 = K1 + (k1(q,i)*pas);
     }
 
     //printf("%f",K1);
     return K1;
 
}
 
float e1(float q, float t)
{
      return sqrt(1-(q*q)*(sin(t)*sin(t))); //fionction de E1(q)
}
 
float integE1(float pas,float inf, float sup, float q)//integrale de ma fonction
{
     int i=0;
     float E1=0;
 
     for(i=inf; i<sup; i+=pas)
     {
          E1 = E1 + (e1(q,i)*pas);
     }
     return E1;
}

les intégrales que j'ai calculés (K1q et E1q) ne marchent pas.
Pourriez-vous m'aider svp :calim2:

10/12/2011, 13h09
herosion

personne pour m'aider? :cry:

J'aimerais juste savoir pourquoi mes intégrales E1 et K1 marchent pas.
en gros:
- j'ai écris ma fonction k1 (float k1(float q, float t) // t=phi
- puis je calcul son intégrale (integK1): j'utilise la méthode des rectangles pour cela.

Je fais la même chose pour E1.

Cordialement
Ben
10/12/2011, 15h47
jeroman

Le problème, c'est que tu nous donnes le scan d'une feuille raturée de partout. On ne sait pas ce que c'est, on ne sait pas non plus quelle partie est bonne ou mauvaise.

Donne-nous les équations en question, au propre.
10/12/2011, 18h07
herosion

1 pièce(s) jointe(s)

Voila j'ai réécris au propre les équations, j’espère que ça vous aidera.

Cordialement
Ben
10/12/2011, 20h17
diogene
Code:

1 2 float N = 0.5; float pas = (sup-inf)/N;
La valeur de N est bizarre !
10/12/2011, 20h41
herosion

effectivement :lol:, j'ai confondu N et la valeur du pas au final. lorsque j'affiche les valeurs de mon intégrale K1, elles restent les mêmes, de km=0.5 a km<5, c'est pas normale.
10/12/2011, 21h29
jeroman

N doit logiquement avoir une valeur (très) élevée. A combien l'as-tu mis ?
10/12/2011, 23h00
diogene
Code:

1 2 3 4 int i=0; float K1=0; for(i=inf; i<sup; i+=pas)
Comment veux-tu que ça marche : i est entier, alors i+= pas va arrondir à la partie entière. Si pas est <1, i ne changera pas (et pas devrait être petir devant 1)

bizarrement j'avais bien mis un float i pour E1 mais pas K1...
J'ai pris N=4,comme ça j'ai un pas relativement petit. Par contre je suis toujours pas sur que mes intégrales fonctionnent bien. enfin bon, le mieux c'est que je trace la courbe pour voir.

Je dois écrire les coordonnées des points dans un fichier puis tracer la courbe mais je ne sais pas comment faire :oops:.

Pour écrire j'ai testé de faire ça:

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
int main(int argc, char *argv[])
{
    FILE* fichier = NULL;
 
    fichier = fopen("X.txt", "r+");
 
    if (fichier != NULL)
    {
        int x= 5000;
        fprintf(fichier,"valeur de x=%d",x);
 
        fclose(fichier);
    }
    else
    {
        printf("Impossible d'ouvrir le fichier test.txt");
    }
 
    return 0;
}

Pourriez-vous m'expliquer comment maintenant écrire mes coordonnées (CK et km) s'il vous plaît.

Pour tracer je me suis renseigner sur la fonction gnuplot mais j'ai pas vraiment compris comment l'utiliser. Surtout que apparemment elle n'est pas installée dans les bibliothèques.

je remet mon code:

Code:

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
#define pi 3.1419
 
 
float k1(float,float);
float integK1(float,float,float,float);
float e1(float,float);
float integE1(float,float,float,float);
 
 
int main()
{
   float km,q,K,K1q,E1q,CK;
 
   float inf = 0.0;
   float sup = pi/2;
   float N = 4.0;
   float pas = (sup-inf)/N;
 
 
  for(km=0.5; km<5; km+=0.1)
  {
      K = 1./km; //km = 1/k
      //printf("%f\n",K);
      q=(2*sinh(2*K))/(cosh(2*K)*cosh(2*K)); //q(K)
      //printf("%f\n",q);
      K1q = integK1(pas,inf,sup,q); //E1(q)
      //printf("%f\n",K1q);
      E1q = integE1(pas,inf,sup,q); //K1(q)
      //printf("\n%f\n",K1q);
 
CK =(4/pi)*((K*(1/tanh(2*K)))*(K*(1/tanh(2*K))))*(K1q-E1q-((1-((tanh(2*K))*tanh(2*K))))*(pi/2+(2*((tanh(2*K))*(tanh(2*K)))-1)-K1q));
printf("%f\n",CK);
  }
 
  return 0;
}
 
 
float k1(float q, float t)
{
 
      return (1/(sqrt(1-(q*q)*(sin(t)*sin(t)))));//ma fonction de K1(q)
 
}
 
float integK1(float pas,float inf, float sup,float q)//integrale de k1
{
     float i=0;
     float K1=0;
 
     for(i=inf; i<sup; i+=pas) // j'utilise la methode des rectangles
     {
          K1 = K1 + (k1(i,q)*pas);
     }
 
     return K1;
 
}
 
float e1(float q, float t)
{
      return (sqrt(1-(q*q)*(sin(t)*sin(t)))); //fionction de E1(q)
}
 
float integE1(float pas,float inf, float sup, float q)//integrale de ma fonction
{
     float i=0;
     float E1=0;
 
     for(i=inf; i<sup; i+=pas)
     {
          E1 = E1 + (e1(q,i)*pas);
     }
     return E1;
}

Merci pour votre aide

11/12/2011, 21h45
diogene

Citation:

bizarrement j'avais bien mis un float i pour E1 mais pas K1...

Bizarrement NON. Les deux sont des int dans le code que tu avais posté

Citation:

J'ai pris N=4,comme ça j'ai un pas relativement petit.

Je verrais plutôt un N >100

Si je met N=100, mes valeurs de CK dans mon fichier texte sont incorrectes

Code:

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
set terminal png
 
km:           CK:
0.500000	0.021293
0.600000	0.056848
0.700000	0.111444
0.800000	0.183298
0.900000	0.272499
1.000000	0.384110
1.100000	0.536260
1.200000	0.835062
1.300000	-1.#IND00
1.400000	-1.#IND00
1.500000	-1.#IND00
1.600000	-1.#IND00
1.700000	-1.#IND00
1.800000	-1.#IND00
1.900000	-1.#IND00
2.000000	-1.#IND00
2.100000	-1.#IND00
2.200000	-1.#IND00
2.300000	-1.#IND00
2.400000	-1.#IND00
2.500000	-1.#IND00
2.600000	-1.#IND00
2.700000	-1.#IND00
2.799999	-1.#IND00
2.899999	-1.#IND00
2.999999	-1.#IND00
3.099999	-1.#IND00
3.199999	-1.#IND00
3.299999	-1.#IND00
3.399999	-1.#IND00
3.499999	-1.#IND00
3.599999	-1.#IND00
3.699999	-1.#IND00
3.799999	-1.#IND00
3.899998	-1.#IND00
3.999998	-1.#IND00
4.099998	-1.#IND00
4.199998	-1.#IND00
4.299998	-1.#IND00
4.399998	-1.#IND00
4.499998	-1.#IND00
4.599998	-1.#IND00
4.699998	-1.#IND00
4.799998	-1.#IND00
4.899998	-1.#IND00
4.999998	-1.#IND00

J'ai téléchargé gnuplot, lorsque je lance le .exe et que je tape la commande: plot "X.txt" using 1:2, on me dit:
Warning: skipping unreadable file "X.txt"
no data in plot

Auriez-vous une solution s'il vous plaît ?

Mon code finale:

Code:

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
#define pi 3.1419
 
 
float k1(float,float);
float integK1(float,float,float,float);
float e1(float,float);
float integE1(float,float,float,float);
 
 
int main()
{
   float km,q,K,K1q,E1q,CK;
 
   float inf = 0.0;
   float sup = pi/2;
   float N = 100;
   float pas = (sup-inf)/N;
 
   FILE* fichier = 0;
   fichier = fopen("X.txt", "w+");
 
   if (fichier)
     fprintf(fichier, "set terminal png\n");
 
  for(km=0.5; km<=5; km+=0.1)
  {
      K = 1./km; //km = 1/k
      //printf("%f\n",K);
      q=(2*sinh(2*K))/(cosh(2*K)*cosh(2*K)); //q(K)
      //printf("%f\n",q);
      K1q = integK1(pas,inf,sup,q); //E1(q)
      //printf("%f\n",K1q);
      E1q = integE1(pas,inf,sup,q); //K1(q)
      //printf("\n%f\n",K1q);
 
      CK =(4/pi)*((K*(1/tanh(2*K)))*(K*(1/tanh(2*K))))*(K1q-E1q-((1-((tanh(2*K))*tanh(2*K))))*(pi/2+(2*((tanh(2*K))*(tanh(2*K)))-1)-K1q));
//printf("%f\n",CK);
 
    if (fichier != 0)
    {
 
      fprintf(fichier, "%f\t",km);
      fprintf(fichier,"%f\n",CK);
    }
    else
    {
      printf("Impossible d'ouvrir le fichier X.txt");
    }
  }
  fclose(fichier);
  return 0;
}
 
 
float k1(float q, float t)
{
 
      return (1/(sqrt(1-(q*q)*(sin(t)*sin(t)))));//ma fonction de K1(q)
 
}
 
float integK1(float pas,float inf, float sup,float q)//integrale de k1
{
     float i=0;
     float K1=0;
 
     for(i=inf; i<sup; i+=pas) // j'utilise la methode des rectangles
     {
          K1 = K1 + (k1(i,q)*pas);
     }
 
     return K1;
 
}
 
float e1(float q, float t)
{
      return (sqrt(1-(q*q)*(sin(t)*sin(t)))); //fionction de E1(q)
}
 
float integE1(float pas, float inf, float sup, float q)//integrale de ma fonction
{
     float i=0;
     float E1=0;
 
     for(i=inf; i<sup; i+=pas)
     {
          E1 = E1 + (e1(q,i)*pas);
     }
     return E1;
}