Discussion :

[marsupial]

Bonjour à tous,

Code demandé en 1995, la loi binomiale aussi nommée Loi de Poisson sur laquelle s'appuit la courbe de Laplace-Gauss du nom du mathématicien normand de l'ère napoléonienne qui déboucha sur l'école libre suite à la Révolution de 1789, peut s'écrire ainsi*:

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
#include <stdio.h>
#include <math.h>
 
long fact(int i)
{
       long f;
 
       if (i==0)
       {
          return 1;
       }
       else
       {
           f = fact(i - 1);
           return i*f;
       }
 
}
 
main ()
{
     int n, k;
     double e, p, q, r, t, x, y, z;
 
     printf("entrez la valeur de k entier non nul positif\n");
     scanf("%d", &k);
     printf("entrez le nombre d'essais n entier non nul positif\n");
     scanf("%d", &n);
 
     x = fact(n);
     y = fact(k);
     z = fact(n-k);
 
     r = x/(y*z);
 
     printf("la combinaison n de k est : %lf\n", r);
 
     x = n-1;
     y = k-1;
     z = x-y;
     x = fact(x);
     y = fact(y);
     z = fact(z);    
 
     e = x/(y*z);
     e = e/x;
     p = pow(e, k);
 
     printf("la probabilite de succes k de n est : %lf\n", p);
 
     x = n-k;
     y = 1-p;
     q = pow(y, x);
 
     r = r*p*q;
     printf("le resultat de la loi binomiale est : %lf\n", r);
     e = n*p;
     printf("l'esperance est de : %lf\n", e);
     t = sqrt(e*q);
     printf("l'ecart type est de : %lf\n", t);
 
     x = exp(-e);
     y = pow (e,k);
     z = x*y;
     k = fact(k);     
     x = z/k;
     printf("la loi de poisson est : %g\n", x);
}

Évidemment, même en 2006, donc a fortiori de nos jours, un compilateur actuel retourne warning et erreurs. A quelques exceptions près comme pour toute règle. Pourquoi.

1- l'architecture Intel 8,16,32,64 ( ou même 4096 ou 8192 ) bits d'une part
2- les limitations de la conception du langage C, qui date d'un demi-siècle, imprévu pour de tels calculs demandant une telle précision

Cependant, il s'agit d'une loi mathématique utilisée dans toutes les industries, y compris financière, afin de déterminer les risques en production et les délais de maintenance, ouvrant la route au big data. Il existe donc un workaround ou solution de contournement*: l'assembleur ou ASM que je ne maîtrise pas suffisamment non par méconnaissance du langage mais de l'architecture du µ.

De ces faits, je le laisse tel quel d'autant que cela produit un résultat aussi aléatoire qu'un ordinateur quantique demandé pour les PKI.

Marsu

[koala01]

Salut,

Je n'ai plus la possibilité de le faire et donc : pourrais-je te demander d'utiliser la balise [ CODE] (le bouton en forme de # qui se trouve en haut de l'éditeur) pour présenter ton code. Ce sera beaucoup plus facile à lire pour tout le monde.

PS : je cherche encore, mais... j'arrive pas à comprendre : tu as un problème avec cette loi de poisson, ou tu veux juste étaler de la confiture

[Bktero]

Je me pose aussi l'intérêt de ce message...

[Obsidian]

Bonjour à tous,

PS : je cherche encore, mais... j'arrive pas à comprendre : tu as un problème avec cette loi de poisson, ou tu veux juste étaler de la confiture

Je me pose aussi l'intérêt de ce message...

À mon avis, ni l'un ni l'autre : je pense que certaines affirmations ici étaient en fait des questions et que son programme, tant à la compilation qu'à l'exécution, doit produire des résultats qui ne sont pas ceux attendus.

Bonjour à tous,

Code demandé en 1995, la loi binomiale aussi nommée Loi de Poisson sur laquelle s'appuit la courbe de Laplace-Gauss du nom du mathématicien normand de l'ère napoléonienne qui déboucha sur l'école libre suite à la Révolution de 1789, peut s'écrire ainsi*:

[…]

Tu as l'air de faire l'hypothèse fréquente mais erronée que le langage mathématique habituel est la première et la plus parfaite façon de modéliser un phénomène, et que l'informatique n'intervient qu'en fin de processus pour obtenir le résultat. C'est une idée assez répandue dans le corps des ingénieurs mais elle n'est pourtant que partiellement vraie. Il suffit par exemple de faire intervenir l'infini dans un bête calcul de limite tel que :

… pour comprendre que c'est strictement exact au point de vue de l'analyse mais que c'est incalculable pour une machine au niveau numérique et ne faisant pas du calcul symbolique.

Évidemment, même en 2006, donc a fortiori de nos jours, un compilateur actuel retourne warning et erreurs. A quelques exceptions près comme pour toute règle. Pourquoi.

1- l'architecture Intel 8,16,32,64 ( ou même 4096 ou 8192 ) bits d'une part

Non, l'architecture Intel fonctionne comme toutes les autres au niveau de la manipulation des données en mémoire, et même la technologie des années 1960 est largement suffisante pour faire des calculs mathématiques, même très compliqués puisque jusqu'alors, on les faisait à la main et sur papier. À la fin des années 1970, les micro-processeurs encore naissants étaient pourtant une technologie déjà au point. Avant même l'apparition du produit phare d'Intel, on avait déjà envoyé Voyager 2 visiter plusieurs planètes en utilisant un ordinateur pour calculer les effets de fronde en cascade.

2- les limitations de la conception du langage C, qui date d'un demi-siècle, imprévu pour de tels calculs demandant une telle précision

Non plus. De toutes façons, la loi de Poisson est elle-même largement antérieure au langage C et aux ordinateurs d'une manière générale, quels qu'ils soient. En outre, le langage C est extrêmement proche du fonctionnel réel des micro-processeurs et c'est pour cela qu'il est encore si répandu même de nos jours. Si tu rencontres des problèmes d'imprécision, tu connaîtras les mêmes avec tous les langages qui s'appuient sur des nombres à virgule flottante normalisés, même ceux qui sont spécialisés dans les mathématiques.

Cependant, il s'agit d'une loi mathématique utilisée dans toutes les industries, y compris financière, afin de déterminer les risques en production et les délais de maintenance, ouvrant la route au big data. Il existe donc un workaround ou solution de contournement*: l'assembleur ou ASM que je ne maîtrise pas suffisamment non par méconnaissance du langage mais de l'architecture du µ.

Toujours pas. L'assembleur, c'est vraiment génial et j'encourage tout le monde à essayer. Je trouve d'ailleurs toujours bénie la période 1980-1995 durant laquelle l'assembleur était à la fois justifié et facile d'accès. Mais dans le cas précis, l'assembleur ne te donnera rien de plus que le langage C.

De ces faits, je le laisse tel quel d'autant que cela produit un résultat aussi aléatoire qu'un ordinateur quantique demandé pour les PKI.
Marsu

Si ton programme te donne des résultats aléatoires sur une application aussi simple que la loi de Poisson, c'est parce qu'il est intrinsèquement mal conçu. Il faut donc repartir de la base et comprendre pourquoi il fonctionne mal plutôt qu'espérer se transcender vers un langage plus parfait. Autrement dit, il faut déboguer ton programme. En se limitant aux premières lignes, on trouve déjà des irrégularités. Par exemple :

Code :

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

1
2
3
4
5
6
7
8
long fact(int i)
{
       long f;
 
       if (i==0)
       {
          return 1;
       }…

… et

Code :

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

1
2
3
4
5
6
7
8
     printf("entrez la valeur de k entier non nul positif\n");
     scanf("%d", &k);
     printf("entrez le nombre d'essais n entier non nul positif\n");
     scanf("%d", &n);
 
     x = fact(n);
     y = fact(k);
     z = fact(n-k);

Donc, si n et k sont tous les deux des entiers positifs mais que k est supérieur à n, alors tu vas passer un nombre négatif à ta fonction factorielle qui n'est pas prévue pour : elle aura déjà dépassé la condition d'arrêt. Mathématiquement, le paramètre de ta fonction récursive devrait plonger vers -∞ et ne jamais se terminer (ou retourner l'infini avec un signe indéterminé). Sur le plan technique, en réalité, le champ de ton nombre entier va atteindre sa valeur minimum, puis remonter par congruence à la valeur maximum et passer en revue la totalité des nombres positifs représentables avant d'atteindre le zéro et donc de se terminer. Le résultat, lui, aura explosé et fait plusieurs fois le tour du domaine des valeurs qu'il peut représenter, et sa valeur sera indéfinie, négative ou positive. Si on ajoute le fait que chaque récursion lorsque l'argument est négatif provoque un changement de signe, lancer un dé au hasard sera effectivement tout aussi fiable.

Et cela, c'est si le compilateur est suffisamment intelligent pour se rendre compte que ta fonction est récursive terminale. Si ce n'est pas le cas, tu vas très vite dépasser la taille de la pile et provoquer une segmentation fault qui fera planter tout ton programme.

Ton programme souffre probablement d'autres tares mais rien que celle-ci suffit, à elle seule, à corrompre la totalité de tes calculs. Remplace déjà « i==0 » par « i<=0 » ce qui n'est pas la panacée (car pas strictement exact), mais devrait déjà limiter les erreurs indiquées ci-dessus.

Lorsque l'on travaille en amateur, il s'agit déjà de bonnes pratiques et d'habitudes à prendre d'emblée. Lorsqu'on le fait dans le domaine professionnel, on se porte garant de la fiabilité des calculs menés par son programme pour les comptes des utilisateurs qui s'en servent et et s'y fient. Et c'est principalement pour cette raison qu'informaticien en général et programmeur en particulier sont des métiers à part entière.

[ternel]

Dis-moi, Obsidian, où as-tu trouvé ta boule de cristal télépatique, qui lit les questions dans la tête d'un posteur?

Parce que la mienne n'a pas cette capacité, et visiblement, ca me manque

[marsupial]

7 ans après, je me replonge sur ce code et je dois remercier Obsidian. Alors ça compile sur un compilateur exotique mais pas sur gcc. J'ai procédé à quelques modifications car il y avait des erreurs de typage où je passe tout en double. Mais je bloque à la compilation sur cette ligne qui demande un int pour l'exposant :

x = frexp(x, -e);

Donc la loi de poisson ne se code pas en C.

Merci à tous pour vos réponses.

edit : un petit lien pour expliquer les probabilités Loi de Poisson

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
#include <stdio.h>
#include <math.h>
 
 
 
 
double fact(double i)
{
       double f;
 
		for(i>0; (f=i*f); i--)
		{
		continue;
		}
		return f;
}
 
int LoidePoisson (double n, double k)
{
 
     double e, p, q, r, t, x, y, z;
 
     printf("entrez la valeur de k entier non nul positif\n");
     scanf("%lf", &k);
     printf("entrez le nombre d'essais n entier non nul positif\n");
     scanf("%lf", &n);
 
     x = fact(n);
     y = fact(k);
     z = fact(n-k);
 
     r = x/(y*z);
 
     printf("la combinaison n de k est : %lf\n", r);
 
     x = n-1;
     y = k-1;
     z = x-y;
     x = fact(x);
     y = fact(y);
     z = fact(z);    
 
     e = x/(y*z);
     e = e/x;
     p = pow(e, k);
 
     printf("la probabilite de succes k de n est : %lf\n", p);
 
     x = n-k;
     y = 1-p;
     q = pow(y, x);
 
     r = r*p*q;
     printf("le resultat de la loi binomiale est : %lf\n", r);
     e = n*p;
     printf("l'esperance est de : %lf\n", e);
     t = sqrt(e*q);
     printf("l'ecart type est de : %lf\n", t);
 
     x = frexp(x,int -e);
     y = pow (e,k);
     z = x*y;
     k = fact(k);     
     x = z/k;
     printf("la loi de poisson est : %lf\n", x);
 
     return 0;
}

[Guesset]

Bonjour,

Il y a confusion entre possibilités d'un langage et l'usage d'une fonction inadéquate.

Par ailleurs, calculer une factorielle en récursif est une hérésie en performances comme en occupation sur la pile (c'est hélas souvent donné comme (contre) exemple pédagogique). Une simple récurrence (boucle) suffirait.

De plus, calculer les combinaisons à partir des factorielles n'est pas le plus efficace. n!/((n-i)! (i!)) = n.(n-1)... max(n-i,i) / min(n-i, i)! On peut même éviter le calcul de la factorielle en dénominateur en précalculant k = n / min(n-i, i) pour sauter dans la boucle n.(n-1)... max(n-i, i) les k.i, k.i - k, k.i - 2.k etc.
Outre le gain d'efficacité cela retarde l'arrivée assez rapide de l'overflow.

Salutations

[Sve@r]

Bonjour

J'ai procédé à quelques modifications car il y avait des erreurs de typage où je passe tout en double.

Déjà tu aurais pensé à l'unsigned...

Mais je bloque à la compilation sur cette ligne qui demande un int pour l'exposant :

x = frexp(x, -e);

Donc la loi de poisson ne se code pas en C.

Ah carrément. Tu n'y arrives pas donc tu décrètes que c'est impossible. Remarque c'est plus facile pour son égo de le voir ainsi...

Code c :

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

1
2
3
4
int LoidePoisson(double k, double lambda) {
	double e=2.71828182846;
	return pow(e, -lambda) * pow(lambda, k) / fact(k);
}

Oups... désolé je ne savais pas que c'était pas possible alors je l'ai fait.

Donc non, ce n'est pas parce qu'on n'y arrive pas que ce n'est pas possible. D'ailleurs, de façon plus générale, tout ce qui peut se faire avec un papier et un crayon peut se faire en C. Alors on n'y arrive pas donc on se dit qu'il faut utiliser la fonction correctement. Déjà elle ne veut pas un int mais un int étoile (j'espère que tu sais que ce n'est pas la même chose). Et donc peut-être réfléchir sur le pourquoi de cet int étoile...
De là on en arrive au but de cette fonction. Elle sert à décomposer un nombre n en k*2^p (avec "k" flottant et "p" entier). Et le int étoile c'est l'adresse de la zone pour y stocker le "p". Et là on se dit alors qu'on n'est pas certain qu'elle soit utile dans la loi de Poisson.
Accessoirement ta fonction fact() ne peut pas donner le bon résultat si tu n'initialises pas "f". Quand à ce "continue" qui vient là sans rime ni raison... (en fait c'est tout le for() qu'il faut revoir car rien ne correspond !!!)

Code c :

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

1
2
3
4
5
6
double fact(double i) {
	double f;
 
	for (f=1; i > 0; f=f*i--);
	return f;
}

De plus, calculer les combinaisons à partir des factorielles n'est pas le plus efficace. n!/((n-i)! (i!)) = n.(n-1)... max(n-i,i) / min(n-i, i)!

Le plus efficace c'est monter les facteurs en parallèle. Par exemple C(3, 8) = C(5, 8) = (6*7*8)/(1*2*3)
Donc on part d'un pivot à 5 (ou 8-3 selon le cas) et d'un indice à 1 et on les monte en parallèle en multipliant et en divisant à chaque itération. Ainsi on évite effectivement les overflow prématurés.
Exemple

Code c :

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

1
2
3
4
5
6
7
unsigned long comb(n, p) {
	unsigned long res=1;
	unsigned short i, j;
	for (i=((p - n) > n ?(p - n) :n)+ 1, j=1; i <= p; i++, j++) res=res * i / j;
 
	return res ;
}

[marsupial]

Merci de m'avoir corrigé !

[Guesset]

Bonjour Sve,

...
Le plus efficace c'est monter les facteurs en parallèle. Par exemple C(3, 8) = C(5, 8) = (6*7*8)/(1*2*3)
Donc on part d'un pivot à 5 (ou 8-3 selon le cas) et d'un indice à 1 et on les monte en parallèle en multipliant et en divisant à chaque itération. Ainsi on évite effectivement les overflow prématurés.
Exemple

Code c :

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

1
2
3
4
5
6
unsigned long comb(n, p) {
	unsigned long res=1;
	unsigned short i, j;
	for (i=((p - n) > n ?(p - n) :n)+ 1, j=1; i <= p; i++, j++) res=res * i / j;
	return res ;
}

Ton code est très bien. Pour faire le malin , je dirais qu'en changeant un peu les initialisations (et en diminuant le nombre de boucles d'une unité), tu peux éviter une multiplication par 1 et une division par 1. Mais c'est du détail.

Salut

[Sve@r]

Pour faire le malin , je dirais qu'en changeant un peu les initialisations (et en diminuant le nombre de boucles d'une unité), tu peux éviter une multiplication par 1 et une division par 1. Mais c'est du détail.

Hélas... enfin la multiplication et division par 1 (multiplication et division neutre) ont au-moins l'avantage de pouvoir être implémentées dans une boucle dans un cadre général (facile de dire ainsi "je multiplie par x et si x vaut 1 tant pis, la multiplication peut quand-même se faire sans avoir besoin de tester).
Alors oui il y a une division par 1, là c'est fatal avec j qui commence à 1. Ca a l'avantage d'avoir une division toujours exacte. En effet, soit le premier terme est impair et le terme suivant pair tombera juste avec la division par 2, soit le premier terme est pair et la division par 2 fonctionne aussi ; avec quoi qu'il arrive le 3° terme multiple de 3 ce qui matchera avec la division par 3 de la 3° itération et ensuite roule.
De plus ok certes on peut faire commencer "res" à 6 et la boucle à 7 (je parle toujours de C(5, 8)) avec "j" qui commence à 2 et ça fonctionne avec C(5, 8) mais il faut faire attention aux cas extrèmes. Exemple C(8, 8) qui est égal à 1 mais pour lequel ton algo (si j'ai bien compris la façon dont tu l'imagines) ne fonctionnera pas.
D'où test supplémentaire pour détecter les cas extrèmes et là on entre dans une nouvelle discussion sur "vaut-il mieux faire le test (lui-aussi consommateur de ressource) ou pas" et bon suis pas certain que cela soit primordial dans la vie d'un dev...

[dalfab]

Bonjour,

Les applications sur lesquelles je travaille utilisent la loi de Poisson. La probabilité qu'un atome se désintègre spontanément est - heureusement pour nous - extrêmement faible, on est typiquement dans une loi de Poisson. La loi de Poisson doit être utilisée avec des valeurs comme k=10000 et lambda=3.5, et là on se rend compte que gérer la valeur 1 dans la factorielle est le détail du détail. Si on applique la formule mathématique: pow(e, -lambda) * pow(lambda, k) / fact(k);, on doit calculer fact(10000) et pow(3.5,10000) qui se retrouvent toutes deux infinies!
Dans ce cas, on passe en log:

Code :

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

double logProbaPoisson = -lambda + k * log(lambda) - log(k!);

Heureusement on a une très bonne approximation de log(factorielle) merci à Ramanujan: log(k!) = k * log(k) - k + log(8*k*k*k + 4*k*k + k + 1./30) / 6..
Ça donne pour k=10000 et lambda=3.5 => logProbaPoisson=-70000 ce qui confirme bien que c'est extrêmement rare, mais dans un gramme de terre il y a quand même un bon paquet d'atomes de ⁴⁰K, et donc que le phénomène est finalement plutôt fréquent!

@marsupial, oui on peut calculer des lois de Poisson en C, et on peut aussi faire plein d'autres choses

[Guesset]

Bonjour Sve,

...Exemple C(8, 8) qui est égal à 1 mais pour lequel ton algo (si j'ai bien compris la façon dont tu l'imagines) ne fonctionnera pas.

Je pensais à quelque chose comme cela :

Code C :

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

1
2
3
4
5
6
7
8
unsigned long comb(n, p) {
   if(n < p - n) n = p - n;
   if(n >= p) return 1;    // si n, p entiers signés, il faudrait tester également les signes
   unsigned long res = n+1;
   unsigned short i, j;
   for (i = n+2, j = 2; i <= p; i++, j++) res = res * i / j;
   return res ;
}

C'est moins élégant et un peu plus verbeux mais il n'y a pas de test en plus car il y en a un de moins dans la boucle.

Relativement, le gain en vitesse est d'autant plus faible que le nombre d'itérations est important. Mais l'overflow guette et nous incite à la modération sur les valeurs et donc le nombre d'itérations.

Salut

[Sve@r]

Code C :

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

1
2
3
4
unsigned long comb(n, p) {
   if(n < p - n) n = p - n;
   if (n >= p) return 1;    // si n, p entiers signés, il faudrait tester également les signes
}

Astucieux

C'est moins élégant et un peu plus verbeux mais il n'y a pas de test en plus car il y en a un de moins dans la boucle.

Euh non là non. Dans ma boucle je commençais par un test pour positionner le départ, et toi tu fais un test similaire en première ligne de ta fonction. N'oublions pas que la première instruction du for() n'est exécutée qu'une seule fois. Ensuite le test de fin de boucle reste le même pour tous les deux (à savoir i <= p).
Ce qui diffère c'est que tu as rajouté un test en seconde ligne pour détecter le cas C(8, 8), test qui n'existe pas dans mon cas. Et on retombe donc sur cette éternelle question "multiplication par 1 ou test n>=p" ??? (accessoirement je pense qu'écrire if (n == p) devrait être plus rapide qu'une inégalité)

Mais l'overflow guette et nous incite à la modération sur les valeurs et donc le nombre d'itérations.

Bah le nb d'itérations c'est au max p. Et là, une de plus ou de moins (p n'est pas présumé être super grand)... Quant à l'overflow lui le seul risque qu'il arrive serait que res monte trop vite durant la boucle, ce qui ne peut pas arriver justement parce qu'il est divisé à chaque itération.

Au fait, pourquoi on parle de combinatoires pour la loi de Poisson?

[Guesset]

Bonsoir Sve,

...Ce qui diffère c'est que tu as rajouté un test en seconde ligne pour détecter le cas C(8, 8), test qui n'existe pas dans mon cas.

Je faisais allusion au fait que ce test en plus est en face du test de boucle qui a lieu une fois de moins dans cette variante (on commence à 2).

Au fait, pourquoi on parle de combinatoires pour la loi de Poisson?

Je me suis aussi posé la question car je ne me rappelais pas que les processus sans mémoire étaient liés à des combinatoires. Mais j'ai eu la flemme de vérifier.

Salut