Lambert W

[gege]

Bonjour,
Non il ne s'agit pas d'un acteur mais de la fonction W de Lambert, qui tente de la calculer ?
A bientôt pour quelques méthodes !
G.E.

Edit : alors la fonction W est définie par :
y*e^y = x => y = W(x)
La première idée est d'utiliser Solve ou tout équivalent sur votre machine, mais cette méthode ne sera pas simple sur nos voeux Basic ou LMS...

[dprtl]

Je me suis toujours demandé à qui pourrait bien être utiles ces fonctions intégrées dans le firmware de la WP-34S ? Extraits du manuel :

[Nikass]

Ouh... c'est pas un des modèles de la combinaison diode-résistance, ça ???

vous voulez calculer ça pour dépanner de la logique DTL/RTL... qui c'est qui a attrapé un 1401 et qui le dit pas, hein ?

a+

[gege]

Bonjour,
En fait avec cette fonction on peut facilement calculer directement les solutions de tout un tas de types d'équations par exemple :
X^n.e^(a.X) = B
X + a.log(X) = B
X + a.e^X = B
qui demanderaient l'usage de "Solve" normalement.
Bon, on ne s'en sert pas constamment non plus…
G.E.

[C.Ret]

[…] la fonction W est définie par : y*e^y = x équivaut à y = W(x)
La première idée est d'utiliser Solve ou tout équivalent sur votre machine, mais cette méthode ne sera pas simple sur nos voeux Basic ou LMS...

Comme je n'ai pas trouvé de fonction SOLVE sur mon SHARP PC-1211, j'utilise une boucle de recherche. Newton et sa méthode fonctionne bien. Je me limite aux valeurs réelles. Les réponses complexes seront envisagées avec une HP-15C

X abscisse [-1/e +oo [ , la fonction de Lambert n'étant pas définie en dessous de -1/e
Y et Z les valeurs pour W(x); la fonction est divalente pour tous les x négatifs fournit deux valeurs Y=W(x) et Z=W(x)
P précision de la détermination; la recherche des résultats se faisant par itérations, soyez modeste et patient.
U,L,E,V et W variables intermédiaires pour le calcul des approximations initiales puis pour l'affinage par la méthode de Newton


La plus grande difficulté a été que les touches W et P ne sont pas disponibles en mode DEF.
J'ai donc choisi DEF-S (en dessous de W) pour saisir l'abscisse et DEF-= (en dessous de P) pour préciser la précision.

Code : Tout sélectionner

1 "S" INPUT "W(X) X=";X : IF -X>EXP -1 PRINT X,"NO VALUE":GOTO 1
2 Y=0 , Z=0 : IF X<>0 LET L= LN ABS X , U=-L-1 , Z=-1-√ABS 2U-U , Y=.56714329 : IF L>0 LET Y=L-LN L+LN L/2L
3 "=" IF X<>0 INPUT "PRES.  ";P : W=Y : GOSUB 8 : Y=W , W=Z : IF X<0 GOSUB 8 : Z=W : PRINT Y,Z:GOTO 3
4 PRINT "W(";X;")=";Y : GOTO 3
8 V=W : IF W<>-1 LET E=EXP W , W=W-(WE-X)/(E+WE) : IF ABS (V-W)>P GOTO 8
9 RETURN
Alternativement pour un affichage plus clair remplacer la dernière ligne par:     9 W=P*INT (W/P):RETURN

En appuyant sur ENTER, la précision est redemandée afin d'avoir éventuellement plus de décimales.
Contrairement à ce qui ce passe lorsqu'une nouvelle valeur de X est saisie, les valeurs initiales de Y et Z ne sont pas réestimées.
En ne changeant que la précision les valeurs Y et Z sont réutilisées, ce qui facilite la convergence.

Quelques exemples de valeurs (à 1E-5 près)

Code : Tout sélectionner

          W(5) = 1.32672                      W(.25) = 0.20388
W(2.718281828) = 1.                             W(0) = 0.
          W(2) = 0.8526                      W(-.25) = -2.1533  et -0.35741
          W(1) = 0.56714                    W(-.367) = -1.0708  et -0.9324  
         W(.5) = 0.35173                      W(-.5)   No value
       W(0.367)= 0.27794                     W(-1/∏) = -1.63853 et -0.55383

Voilà qui est bien pratique même sur nos chers et désirés matériels vingt-ages.

[Ythunder]

Perso, je comprenais plus avec Lambert Wilson...

[gege]

Bonjour,
Bravo, avec le programme Basic que voici on obtient par exemple W(5)=1.32672466524.

10 INPUT A:Y=LOG(A)/(.6978+.02084*LOG(A))
20 IF 10000<A THEN Y=Y-5
30 IF 10000>=A THEN Y=Y-log(A)
40 B=EXP(Y):X=Y*B:Z=(X-A)/(X+B):Y=Y-Z
50 IF ABS(Z)<1E-8 THEN 40
60 PRINT Y

il ne fait que la branche principale.
C'est un Newton, qui fait 5 ou 6 boucles en général.
D'autres algorithmes sont plus rapides...
A affiner
G.E

[gege]

Bonjour,
Hop un "Halley d'ordre 3" (moi non plus je ne sais pas ce que c'est) recopié du HP Museum et traduit du RPN :

10 INPUT A:X=LN(1+A):IF X>0 THEN X=X-LN(X)
20 FOR K=1 TO 5:Y=X:B=1+1/(1+X)
30 B=(X+1)/(X-A/EXP(X))-B/2:X=X-1/B
40 IF X=Y THEN K=6
50 NEXT K:PRINT Y

Pareil branche principale seulement, il fait 4 ou 5 boucles en général.
Qu'est-ce qu'on s'amuse !!
G.E.

[C.Ret]

J'ai modifié un peu mon code pour compter le nombre de fois que la boucle se fait dans le sous-programme de la ligne 8.

Toujours pour une précision de 1.E-5, le nombre de tours varie entre 2 et 12 selon la valeur principale ou négative. Le nombre d'itération est d'autant plus important que l'on se rapproche de la valeur critique -1.

Je ne comprends pas pourquoi les programmes présentés par gégé ne fonctionnent que pour la branche principale. La courbe est très régulière, il suffit de nourrir la méthode de Newton (ou la méthode de Halley) avec une valeur initiale appropriée.
L'affinement du calcul se fait alors en quelques (4 ou 6) itérations. Il suffit juste de bien placer la première approximation sur la branche principale ou sur la branche négative (cf. courbe donnée par notre ami dprl).

W Lambert (branches réelles).gif (19.41 Kio) Consulté 5188 fois

C'est ce que fait mon code à la ligne 2:
Y est la racine approchée se trouvant sur la branche principale estimée en fonction du X donné.
Z est la racine approchée se trouvant sur la branche négative estimée en fonction du logarithme népérien du X donné.
Ces deux valeurs servent d'estimations initiales (GOSUB à la méthode de Newton. Evidemment, la recherche de la valeur sur la branche négative n'est déclenchée que sir le X donné est négatif (et supérieur à -1/e).

Pour approcher les estimations initiales, j'utilise les encadrements suivants:

Lambert W function bounds.gif (10.39 Kio) Consulté 5195 fois

Par conception, Z est calculé quelque soit le signe de X car son estimation fait intervenir des variables intermédiaires servant conjointement à l'estimation de Y. Mais la valeur Z ainsi calculée pour les X positif n'a pas de sens (et ne peut servir pour la suite si X>0).
Insidieusement, l'estimation de Y ne peut être faite (division par zéro) lorsque X=1. Dans ce cas, Y est directement initialisé avec le GRAND OMEGA

Et YThunder va regretter que je n'ai pas utilisé le GRAND MEGATRON, mais avec les frais de représentation que me demande Lambert Wilson et ce très cher Newton, je n'ai plus les moyens d'inviter les Transformers.

[Ythunder]

Ah ! Les Transformers, ça je comprends aussi !!
Merci de m'égayer sur ce post, sérieux ça me fait super plaisir !!!!!!!!

[C.Ret]

[…]Les réponses complexes seront envisagées avec une HP-15C

Chose promise, chose due !

Voici donc une version de code pour HP-15C qui tire partie de ses capacités à effectuer les calculs avec des nombres complexes.

HP15C - lambert W (real and complex).gif (43.36 Kio) Consulté 5087 fois

Ce code contient trois parties dont les deux premières permettent le calcul des valeurs réelles de la fonction W de Lambert se trouvant respectivement sur la branche principale (label A) et la branche négative (label B).
Le troisième programme (label C) permet le calcul des valeurs complexes de la fonction de Lambert : w=W(z) <=> z=w.e^w avec w = u +i.v et z=a+ i.b tous deux complexes

Pour déterminer la valeur réelle sur la branche principale de W, saisir le réel x et presser f-A.
Pour déterminer la valeur réelle sur la branche négative de W, saisir le réel x et presser f-B.
Pour déterminer la valeur complexe de la fonction W au point z, saisir le nombre complexe z et presser f-C.

Pour enter un nombre complexe z=a+i.b, procéder comme indiqué dans le manuel, saisir a, presser (ENTER), saisir b et presser f - I

L'HP-15C mouline pendant quelques dizaines de secondes puis affiche le résultat.
Guetter l'indicateur 'c' de l'affichage qui indique une valeur complexe. La partie réelle du résultat complexe est affichée, pour afficher la partie imaginaire presser f - (i) et maintenir la touche (i) enfoncée pour lire la partie imaginaire. Alternativement, utiliser f - re:im pour échanger partie réelle et complexe du registre x:

La justesse du résultat affiché peut être vérifiée en pressant (ENTER) (e^x) ( x ), la valeur ainsi calculée doit être identique à la donnée initiale.

" ERROR 0 " indique généralement une erreur de calcul due au fait que l'argument initiale n'existe pas sur la branche sélectionnée.
" ERROR 3 " indique que l'on a cherché à afficher la partie imaginaire d'un résultat réel.
" '-----.---' " un affichage clignotant indique un dépassement de capacité, la valeur demandée ne peut être déterminée.
" running " tout ce passe bien, le résultat est un cours de calcul, lors de l'itération les résultats approchés successifs s'affichent ponctuellement pour faire patienter.

Le code utilise trois registres et le mode complexe (activation/désactivation du drapeau 8 )
R0 : partie réelle de la donnée initiale
R1 : partie imaginaire de la donnée initiale
R I : compteur d'itération (initialisé à 9 à la ligne 041-)
Les trois sous-programmes sont :
Lbl 0 : calcul du résultat par itérations successives (la fonction (SOLVE) n'étant d'aucun intérêt lors de calculs avec les complexes)
Lbl 1 : calcul de l'approximation intiale et initialisation du compteur d'itération I
Lbl 2: sous-programme pour sous traction de la partie imaginaire (si flag 8 armé)

[gege]

Bonjour,
Curieux ton programme est incroyablement court et je ne vois pas où sont effectués les calculs ??
Je sors la LE...
G.E.

[ledudu]

Salut,

ton programme est incroyablement court

Oui gégé, cela s'appelle la programmation synthétique.

[gege]

Bonjour,
Ah ok le programme court c'est synthétique ?
Tu es sûr qu'on parle caltos là ?

G.E.

[bernouilli92]

Il y a un une boucle entre les pas 12 et 31 avec des calculs utilisant notamment l'exponentielle.