Post Reply 
HP49-50G,VER17.5 Lambert Funktion Wk(x), k=k= 0, ±1, ±2, ±3, ±4..., x real or complex
12-18-2024, 03:15 AM (This post was last modified: 12-21-2024 10:47 PM by Gil.)
Post: #15
RE: HP49-50G,VER16.2 Lambert Funktion Wk(x), k=k= 0, ±1, ±2, ±3, ±4..., x real or complex
New version 17.5

Corrections depending on the given arguments.

Code:
\<< "Lambert (ver 17.5)

 2 inputs 
 k & x with branch k
   =0, \1771, \1772, \1773, \1774

 Or 1 input x possible
  for k = 0 by default
 when stack level 2
  .\=/0, \1771, \1772, \1773, \1774
  .or empty

 .With x=W(x)*e^W(x),
  W(x) unknown
 .x may be a complex #
   (a,b), even ext.big
   complex, i.e a or b
   =\1779.99999999999E499

 Output(s) 
 -1 output Wk(x): 
   . w
   . or '\oo \177 in\pi/2'
     for ext.big comp.
     (\=/Wolfram\-> \oo)
   . or '-\oo \177 in\pi'
     for x\->-0 or x\->0+
     (\=/Wolfram\-> -\oo)
 -Or 2 outputs,
  if -1/e \<=x \<=0
      & 'k=0 or k=-1':
   . W0 (x): w1 
   . W-1(x): w2
     
Formulae initial guess
.Albert Chan
.WIKIP: W-1 & -1/e\<=x\<=0

For branch k of Wk(x):
https://www.johndcook.com/blog/2021/11/15/lambert-w-branch-number/
" DROP
  IFERR OVER
  THEN 0
  ELSE \-> k
    \<< k TYPE 0 == k TYPE 28 == OR
      IF
      THEN k FP 0 == { SWAP } 0 IFTE
      ELSE 0
      END
    \>>
  END SWAP DUP \->NUM 0 RCLF 0 0 0 0 .367879441171 4.42321595524E-13 0 0 0 0 0 \-> k x0 x f fg xR xi t x0\175 R R2 H X Y Z h
  \<< -3 CF x IM DUP 'xi' STO 0 ==
    IF
    THEN x RE 'x' STO
    END x RE 0 \<= x RE R NEG \>= AND k -1 == AND
    IF
    THEN 0 'k' STO
    END -103 SF x RE DUP 'xR' STO ABS \->NUM \oo \->NUM == xi ABS \->NUM \oo \->NUM == OR
    IF
    THEN "W" k R\->I + "((" +
      CASE xi \oo \->NUM NEG ==
        THEN xR ABS \oo \->NUM ==
          IF
          THEN xR 0 > "\oo" "-\oo" IFTE
          ELSE xR DUP FP 0 == { R\->I } IFT
          END + ",-\oo))" + 4 k * 1 - DUP 0 > "+i\183" "-i\183" IFTE SWAP ABS R\->I + "\pi/2" +
        END xi \oo \->NUM ==
        THEN xR ABS \oo \->NUM ==
          IF
          THEN xR 0 > "\oo" "-\oo" IFTE
          ELSE xR DUP FP 0 == { R\->I } IFT
          END + ",\oo))" + 4 k * 1 + DUP 0 > "+i\183" "-i\183" IFTE SWAP ABS R\->I + "\pi/2" +
        END xR \oo \->NUM NEG == xi 0 < AND
        THEN "-\oo," + xi DUP FP 0 == { R\->I } IFT + "))" + 2 k * 1 - DUP 0 > "+i\183" "-i\183" IFTE SWAP ABS R\->I + "\pi" +
        END xR \oo \->NUM NEG ==
        THEN "-\oo," + xi DUP FP 0 == { R\->I } IFT + "))" + 2 k * 1 + DUP 0 > "+i\183" "-i\183" IFTE SWAP ABS R\->I + "\pi" +
        END "\oo," + xi DUP FP 0 == { R\->I } IFT + "))" + 2 k * DUP 0 > "+i\183" "-i\183" IFTE SWAP ABS R\->I + "\pi" +
      END "i\1830" "" SREPL 1 ==
      IF
      THEN DROP '\oo'
      ELSE "\oo " SWAP + "\pi" "\183\pi" SREPL DROP "\1831\183" "\183" SREPL DROP " " " <U>" SREPL DROP "<U>" + "<U>" 19 CHR 3 CHR OVER + + SREPL DROP
      END SWAP \->TAG
    ELSE -105 SF RAD '2*k*\pi' x ARG + i * \->NUM 't' STO
      CASE xi 0 == xR R NEG < AND xR -.38 > AND
        THEN '\v/(2*R*(x+R+R2))*(k+k+1)' \->NUM 'Y' STO 'R+Y*\v/(1+Y/(3*R))' \->NUM 'Y' STO
          DO Y 'Y+x' '(Y-R-R2)/R' \->NUM DUP TYPE 1 ==
            IF
            THEN 1 + LN
            ELSE LNP1
            END / - \->NUM 'H' STO 'Y-H' \->NUM 'Y' STO Y 'Y+H*.0001' \->NUM ==
          UNTIL
          END x Y / "W" k R\->I + "(" + x0 TYPE 9 ==
          IF
          THEN x0 + "=" +
          END x xR FP 0 == xi 0 == AND xR ABS 1.E13 < AND { R\->I } IFT + ")" + \->TAG fg STOF KILL
        END k 1 == xi 0 < AND k -1 == xi 0 \>= AND OR
        THEN x NEG LN x ABS .5 <
          IF
          THEN 0
          ELSE t 2 /
          END +
        END k 0 ==
        THEN xi 0 == xR 0 > AND
          IF
          THEN x 2 * LNP1
          ELSE x 2 * 1 + DUP 0 ==
            IF
            THEN DROP 1.6 DUP R\->C
            ELSE LN
            END
          END 2 /
        END t
      END 1 \->ARRY xi 0 == xR 0 \<= AND xR e INV NEG \>= AND k 0 == k -1 == OR AND
      IF
      THEN 'W*e^W' x - DUP 'f' STO 'W' ROT 1 GET ROOT
        CASE x e INV NEG \->NUM == x0 TYPE 9 == AND
          THEN DROP -1
          END x -.3615 \<= x e NEG INV > AND
          THEN DROP '(x+R+R2)/R' \->NUM 'H' STO '\v/(2*H)*(2*k+1)' \->NUM 'X' STO 'X*\v/(1+X/3)' \->NUM 'X' STO 'R+R*X' \->NUM 'Y' STO
            DO X DUP TYPE 1 ==
              IF
              THEN 1 + LN
              ELSE LNP1
              END 'Z' STO 'X-(X-Z+H)/Z' \->NUM 'Z' STO X Z - 'X' STO X 'X+Z*.000001' \->NUM ==
            UNTIL
            END X DUP TYPE 1 ==
            IF
            THEN 1 + LN
            ELSE LNP1
            END 1 -
          END
        END DUP FP 0 == { R\->I } IFT "W0(" x0 TYPE 9 ==
        IF
        THEN x0 + x R NEG == "~" "=" IFTE +
        END x xR FP 0 == xi 0 == AND { R\->I } IFT + ")" + \->TAG
        CASE x 0 < x e INV NEG \>= AND
          THEN x 0 < x ABS 1.E-495 < AND
            IF
            THEN x NEG LN 1 3
              START 'x0\175' STO '(1-(LN(-x0\175)-LN(-x)))/(1+1/x0\175)' \->NUM
              NEXT
            ELSE f 'W'
              IF x 4 INV NEG \<=
              THEN '-1-\v/(2*(1+e*x))'
              ELSE 'LN(-x)-LN(-LN(-x))'
              END \->NUM ROOT x e NEG INV \->NUM == x0 TYPE 9 == AND
              IF
              THEN DROP -1
              END
            END
          END x 0 ==
          THEN -105 CF "-\oo (-\oo<U>-i\pi<U>)" "<U>" 19 CHR 3 CHR OVER + + SREPL DROP
          END -1
        END x -.3636 \<= x e NEG INV > AND
        IF
        THEN -1 'k' STO DROP '(x+R+R2)/R' \->NUM 'H' STO '\v/(2*H)*(2*k+1)' \->NUM 'X' STO 'X*\v/(1+X/3)' \->NUM 'X' STO 'R+R*X' \->NUM 'Y' STO
          DO X DUP TYPE 1 ==
            IF
            THEN 1 + LN
            ELSE LNP1
            END 'Z' STO 'X-(X-Z+H)/Z' \->NUM 'Z' STO X Z - 'X' STO X 'X+Z*.000001' \->NUM ==
          UNTIL
          END X DUP TYPE 1 ==
          IF
          THEN 1 + LN
          ELSE LNP1
          END 1 -
        END "W-1(" x0 TYPE 9 ==
        IF
        THEN x0 + x R NEG == "~" "=" IFTE +
        END x xR FP 0 == xi 0 == AND { R\->I } IFT xR 0 ==
        IF
        THEN DROP "-0 (0+)"
        END + ")" + \->TAG 'W' PURGE
      ELSE x 0 \=/
        IF
        THEN 'k*2*\pi*i+LN(x)' 'W+LN(W)' - 1 \->ARRY [ 'W' ] ROT MSLV 1 GET DUP IM 0 == { RE } IFT UNROT DROP2 DUP \-> r
          \<< TYPE 1 ==
            IF
            THEN r C\->R DUP ABS DUP .0001 + \pi 2 / \>= SWAP \pi 2 / / \->NUM DUP 0 RND - ABS .0001 \<= AND ROT ABS .00000000001 < AND
              IF
              THEN 0 SWAP R\->C
              ELSE DROP r
              END
            ELSE r
            END
          \>>
        ELSE DROP "-\oo<U>" 2 k * k SIGN - DUP 0 >
          IF
          THEN "+i"
          ELSE "-i"
          END \-> kk i
          \<< i + kk ABS k 0 > 1 -1 IFTE + R\->I + "\pi<U>" + " (-\oo<U>" + i + kk ABS + "\pi<U>)" + "i1\pi" "i\pi" SREPL DROP "<U>" 19 CHR 3 CHR OVER + + SREPL DROP
          \>>
        END "W" k R\->I + "(" + x0 TYPE 9 ==
        IF
        THEN x0 + "=" +
        END x xR 0 == xi 0 == AND
        IF
        THEN DROP "\->-0 (\->0+)"
        ELSE xR FP 0 == xi 0 == AND xR ABS 1.E13 < AND { R\->I } IFT
        END + ")" + \->TAG
      END
    END fg STOF
    IFERR 0 DOERR
    THEN
    END
  \>>
\>>


Attached File(s)
.hp  lambert.17.5.hp (Size: 5.93 KB / Downloads: 0)
Find all posts by this user
Quote this message in a reply
Post Reply 


Messages In This Thread
RE: HP49-50G,VER16.2 Lambert Funktion Wk(x), k=k= 0, ±1, ±2, ±3, ±4..., x real or complex - Gil - 12-18-2024 03:15 AM



User(s) browsing this thread: 1 Guest(s)