Re: Differentiable Forth

On Fri, 19 Jul 2024 18:01:39 +0000, minforth wrote:

I don't know if gforth has z: locals for complex numbers. If yes, they
could be used for dual fp-numbers as well to reduce code overhead and
improve readability.

Yes, z: locals exists in gforth but I haven't used it. I didn't want to
mix dual and complex numbers.
I have completed the code with the elementary functions and their
derivatives (exp, ln, trigs, inverse trig, hyperbolic, inverse
hyperbolic, ...).
The code is not optimized and can be improved. I've done it just as a
proof of concept.
I added examples for functions with multiple variables (3 variables x y
and z)
Perhaps, one must separate dual numbers, dual number valued functions of
dual numbers variables, ..., and the examples.
Here is the code:
: -frot frot frot ;
: f>dl   ( f: a -- a 0) 0e ; \ real to dual value
: f>dl_d ( f: a -- a 1) 1e ; \ real to dual value with respect to
differentiate
: dl>rl  ( f: a b -- a) fdrop ;
: dl>eps ( f: a b -- b) fnip ;
: dl. ( f: a b -- ) fswap f. ." + eps " f. ;
: dl(,). ( f: a b -- ) ." (" fswap f. ." , " f. ." )" ;
: dl+ ( f: a b c d -- a+c b+d)
    frot ( f: a c d b)
    f+   ( f: a c b+d)
    -frot ( f: b+d a c)
    f+   ( f: b+d a+c)
    fswap ( f: a+c b+d) ;
: dl- ( f: a b c d -- a-c b-d)
    frot ( f: a c d b)
    fswap f-   ( f: a c b-d)
    -frot ( f: b-d a c)
    f-   ( f: b-d a-c)
    fswap ( f: a-c b-d) ;
fvariable b*c
: dl* ( f: a b c d -- a*c  a*d+b*c)
    -frot ( f: a d b c)
    ftuck f* ( f: a d c b*c)
    b*c f! ( f: a d c)
    frot ( f: d c a)
    ftuck ( f: d a c a)
    f*    ( f: d a a*c)
    -frot  ( f: a*c d a)
    f* b*c f@ f+ ( f: a*c a*d+b*c)
;
: 1/dl ( f: a b -- 1/a -b*1/a^2)
    fswap 1/f ( f: b 1/a)
    ftuck     ( f: 1/a b 1/a)
    fdup f*   ( f: 1/a b 1/a^2)
    f* fnegate ( f: 1/a -b/a^2)
;
: dl/ ( f: a b c d -- a/c rest)
    1/dl dl*
;
: dl^f ( f: a b c -- a^c b*c)
    ftuck ( f: a c b c)
    f* -frot ( f: b*c a c)
    f** fswap ( f: a^c b*c)
;
\
: dldup ( f: a b -- a b a b) fover fover ;
: dlnip ( f: a b c d -- c d) frot fdrop frot fdrop ;
: dldrop ( f: a b -- ) fdrop fdrop ;
fvariable dlswap_temp
: dlswap ( f: a b c d -- c d a b)
    frot dlswap_temp f! ( f: a c d)
    frot ( f: c d a)
    dlswap_temp f@  ;
fvariable dlover_temp1
fvariable dlover_temp2
: dlover ( f: a b c d -- a b c d a b)
    dlswap ( f: a b c d -- c d a b)
    dlover_temp2 f! dlover_temp1 f! ( f: c d)
    dlover_temp1 f@ dlover_temp2 f@ ( f: c d a b)
    dlswap ( f: a b c d)
    dlover_temp1 f@ dlover_temp2 f@ ( f: a b c d a b)
;
: dltuck dlswap dlover ;
: dlvariable create 2 floats allot ;
: dl! ( dlvar -- ) ( f: f1 f2 -- ) dup float+ f! f! ;
: dl@ ( dlvar -- ) ( f: -- f1 f2)  dup f@ float+ f@ ;
dlvariable dlrot_temp1
dlvariable dlrot_temp2
: dlrot ( dl: d1 d2 d3 -- d2 d3 d1)
    dlrot_temp1 dl! ( dl: d1 d2)
    dlswap    ( dl: d2 d1)
    dlrot_temp2 dl! ( dl: d2)
    dlrot_temp1 dl@ ( dl: d2 d3)
    dlrot_temp2 dl@ ( dl: d2 d3 d1)
;
\ dual number funuctions of dula number variables
: dlident ( f: a b -- a b) ;
: dlexp ( f: a b -- c d) fswap fdup fexp fswap fexp frot f* ;
: dlln ( f: a b -- c d)  fswap fdup fln  fswap 1/f frot f* ;
: dl^  ( dl: a b -- a^b) \ a^b = exp(b*ln(a))
     dlswap dlln dl* dlexp
;
: dlsqrt ( dl: x --y) 0.5e f>dl dl^ ;
: dlsin ( f: a b -- c d) fswap fdup fsin fswap fcos frot f* ;
: dlcos ( f: a b -- c d) fswap fdup fcos fswap fsin fnegate frot f* ;
: dltan ( dl: x -- y) dldup dlsin dlswap dlcos dl/ ;
: dlcot ( dl: x -- y) dltan 1/dl ;
: dlsinh ( dl: x -- y) dldup dlexp dlswap -1e f>dl dl* dlexp dl- 2e f>dl
dl/ ;
: dlcosh ( dl: x -- y) dldup dlexp dlswap -1e f>dl dl* dlexp dl+ 2e f>dl
dl/ ;
: dltanh ( dl: x -- y) dldup dlsinh dlswap dlcosh dl/ ;
: dlcoth ( dl: x -- y) dltanh 1/dl ;
: dlasin ( f: a b -- c d) fswap fdup fasin fswap fdup f* 1e fswap f-
fsqrt 1/f frot f*  ;
: dlacos ( f: a b -- c d) fswap fdup facos fswap fdup f* 1e fswap f-
fsqrt 1/f fnegate frot f* ;
: dlatan ( f: a b -- c d) fswap fdup fatan fswap fdup f* 1e f+ 1/f frot
f* ;
: dlasinh ( f: a b -- c d) fswap fdup fasin fswap fdup f* 1e f+ fsqrt
1/f frot f* ;
: dlacosh ( f: a b -- c d) fswap fdup facos fswap fdup f* 1e f- fsqrt
1/f fnegate frot f* ;
: dlatanh ( f: a b -- c d) fswap fdup fatanh fswap fdup f* 1e fswap f-
1/f frot f* ;
\ ..... add others
\ derivatives
variable func
: der 1e ' func ! func @ execute dl>eps ;
\ examples
1 [if]
: dlf() 1e 0e dl+ ; \ f(x) = x + 1
: dlf2() dldup dl* ; \ f2(x) = x^2
: dlf3() dldup dlf2() dlswap dlf() dl+ ; \ f3(x) = x^2 + x + 1
: der_f3() ( f: x -- y) 2e f* 1e f+ ; \ d/dx(f3) = 2*x + 1
cr 1e der dlf3() f. \ 3. ok
cr 1e der_f3() f.   \ 3. ok
cr cr
: dlf4() dlf3() dlexp ; \ f3(x) = exp(x^2+x+1)
: der_f4() ( f: x -- y) \ d/dx(f4) = (2*x+1)*exp(x^2+x+1)
    fdup 2e f* 1e f+ fswap
    fdup fdup f* f+ 1e f+ fexp f*
;
cr 2e der dlf4() f. \ 5483.16579214229  ok, calculated at 2e
cr 2e der_f4() f.   \ 5483.16579214229  ok
cr cr
: dlf5() 1/dl ;
: der_f5() ( f: x) \ d/dx(f5) = -1/x^2
    fdup f* 1/f fnegate
;
cr 2e der dlf5() f. \ -0.25  ok calculated at 2e
cr 2e der_f5() f.   \ -0.25
cr cr
: dlf6() dldup dldup dl* dlswap dlsin dl+ 1/dl ; \ f6(x) =
1/(x^2+sin(x))
\ using the derivative calculated analytically d/dx (f6(x)) =
-(2*x+cos(x))/(x^2+sin(x))^2
: der_f6() ( f: x -- y) fdup fdup fdup f* fswap fsin f+ fdup f* 1/f
fswap fdup 2e f* fswap fcos f+ f* fnegate ;
cr 1e der dlf6() f. \ -0.749127330692909  ok calculated at x = 1
cr 1e der_f6() f.   \ -0.749127330692909  ok
cr cr
: dlf7() dldup dldup dl* dlcos dl* ; \ f7(x) = x*cos(x^2),
: der_f7() ( f: x --y)  \ its deriv: d/dx(f7) = cos(x^2)-2*x^2*sin(x^2)
    fdup f* ( f: x^2)
    fdup fsincos ( f: x^2 s c )
    -frot ( f:  c x^2 s )
    f* 2e f* ( f: c 2s*x^2)
    f-
;
cr 1e der_f7() f.   \ -1.14263966374765  ok calculated at 1e
cr 1e der dlf7() f. \ -1.14263966374765  ok
cr cr
cr 2e der_f7() f.   \ 5.40077634159981  ok calculated at 2e
cr 2e der dlf7() f. \ 5.40077634159981  ok
cr cr
: dlf8() \ f8(x) = sin(sin(sin(x)))
    dlsin dlsin dlsin
;
: der_f8() \ d/dx(f8) = (sin(sin(x)))'*cos(sin(sin(x)))
           \          = (sin(x))'*cos(sin(x)*cos(sin(sin(x)))
           \          = cos(x)*cos(sin(x))*cos(sin(sin(x)))
    ( f: x -- y)
    fsincos ( f: s c)
    fswap   ( f: c s)
    fdup    ( f: c s s)
    fcos    ( f: c s cs)
    fswap   ( f: c cs s)
    fsin    ( f: c cs ss)
    fcos    ( f: c cs css)
    f* f*
;
cr 2e der dlf8() f. \
cr 2e der_f8() f.   \
cr cr
( f9 function)
: dl_f9() ( dl: x y -- z) \ f9(x,y) = x + 5*y + x*y
     dlover dlover dl* dlswap 5e f>dl dl* dl+ dl+
;
: df9/dx ( f: x y -- z) \ df9/dx = 1 + y
     fnip 1e f+
;
: df9/dy ( f: x y -- z) \ df9/dy = 5 + x
     fdrop 5e f+
;
cr 2e f>dl_d 3e f>dl dl_f9() dl>eps f. \
cr 2e    3e    df9/dx        f. \
cr
cr 2e f>dl 3e f>dl_d dl_f9() dl>eps f. \
cr 2e    3e    df9/dy        f. \
cr cr
( f10 function)
: dl_f10() ( dl: x y -- z) \ f10(x,y) = exp(x + 5*y) * sin(x*y)
     dlover dlover dl* dlsin ( dl: x y sxy)
     dlrot dlrot             ( dl: sxy x y)
     5e f>dl dl* dl+ dlexp   ( dl: sxy e^[x+5y])
     dl*
;
: df10/dx ( f: x y -- z) \ df10/dx = exp(x+5y)*(sin(x*y)+y*cos(x*y))
     fover fover 5e f* f+ fexp ( f: x y e^[x+5y])
     frot frot ( f: e^[x+5y] x y)
     ftuck     ( f: e^[x+5y] y x y)
     f* fsincos ( f: e^[x+5y] y sxy cxy)
     frot f* f+ ( f: e^[x+5y] sxy+ycxy)
     f*
;
: df10/dy ( f: x y -- z) \ df10/dy = exp(x+5y)*(5*sin(x*y) + x*cos(x*y))
     fover fover 5e f* f+ fexp ( f: x y e^[x+5y])
     frot frot ( f: e^[x+5y] x y)
     fover f*  ( f: e^[x+5y] x yx)
     fsincos   ( f: e^[x+5y] x sxy cxy)
     frot f* fswap 5e f* f+ ( f: e^[x+5y] 5sxy+xcxy)
     f*
;
cr 2e f>dl_d 3e f>dl dl_f10() dl>eps f. \
cr 2e    3e    df10/dx        f. \
cr
cr 2e f>dl 3e f>dl_d dl_f10() dl>eps f. \
cr 2e    3e    df10/dy        f. \
cr cr
cr .( gradient)
: gradient ( f: x y -- df/dx df/dy)
     ' func !
     fover fover ( f: x y x y)
     fswap f>dl_d ( f: x y y x 1)
     frot  f>dl   ( f: x y x 1 y 0)
     func @ execute fnip ( f: x y df/dx)
     frot frot ( f: df/dx x y)
     fswap f>dl ( f: df/dx y x 0)
     frot f>dl_d ( f: df/dx x 0 y 1)
     func @ execute fnip ( f: df/dx df/dy)
;
cr 2e 3e gradient dl_f10() fswap cr f. cr f.
cr cr
cr .( gradient using values)
0e fvalue x
0e fvalue y
: grad ( f: x y -- df/dx df/dy)
    ' dup
    to y to x
    x f>dl_d y f>dl   execute fnip
    x f>dl   y f>dl_d execute fnip
;
(  defined)
cr 2e 3e grad dl_f10() fswap cr f. cr f.
cr cr
cr .( gradient using arrays)
create xy 2 floats allot
xy 2 floats erase
: grad1 ( f: x y -- df/dx df/dy)
    ' dup
    xy float+ f!
    xy f!
    xy f@ f>dl_d xy float+ f@ f>dl   execute fnip
    xy f@ f>dl   xy float+ f@ f>dl_d execute fnip
;
(  defined)
cr 2e 3e grad1 dl_f10() fswap cr f. cr f.
cr cr
cr .( gradient using arrays)
create xy2 2 floats allot
xy2 2 floats erase
: to_xy2  xy2 float+ f! xy2 f! ;
: df/dx xy2 f@ f>dl_d xy2 float+ f@ f>dl   execute fnip ;
: df/dy xy2 f@ f>dl   xy2 float+ f@ f>dl_d execute fnip ;
: grad2 ( f: x y -- df/dx df/dy)
    ' dup
    to_xy2
    df/dx
    df/dy
;
(  defined)
cr 2e 3e grad2 dl_f10() fswap cr f. cr f.
cr cr
cr .( gradient using arrays, 3 variables)
: dl_f11() ( dl : x y z -- t)  \ f11(x,y,z) = exp(-x^2)*(cos(2y)+sin(z))
     dlsin dlswap 2e f>dl dl* dlcos dl+ ( dl: x cos[2y]+sin[z])
     dlswap dldup dl* -1e f>dl dl* dlexp dl*
;
\ for verif
: df11/dx ( f: x y z -- df11/dx) \ df11/dx =
-2*x*exp(-x^2)*(cos(2y)+sin(z))
      fsin fswap 2e f* fcos f+ fover f* -2e f* ( f: x
-2*[cos[2y]+sin[z]])
      fswap fdup f* -1e f* fexp f*
;
: df11/dy ( f: x y z -- df11/dy) \ df11/dy = exp(-x^2)*(-2)*sin(2*y)
      fdrop 2e f* fsin -2e f* fswap fdup f* -1e f* fexp  f*
;
: df11/dz ( f: x y z -- df11/dz) \ df11/dz = exp(-x^2)*cos(z)
      fnip fcos fswap fdup f* -1e f* fexp f*
;
create xyz 3 floats allot
xyz 3 floats erase
: to_xyz  xyz 2 floats + f! xyz float+ f! xyz f! ;
: df/dx xyz f@ f>dl_d xyz float+ f@ f>dl   xyz 2 floats + f@ f>dl    execute fnip ;
: df/dy xyz f@ f>dl   xyz float+ f@ f>dl_d xyz 2 floats + f@ f>dl    execute fnip ;
: df/dz xyz f@ f>dl   xyz float+ f@ f>dl   xyz 2 floats + f@ f>dl_d  execute fnip ;
: grad2 ( f: x y z -- df/dx df/dy df/fz)
    ' dup dup
    to_xyz
    df/dx
    df/dy
    df/dz
;
(  defined)
cr 1e 2e 3e grad2 dl_f11() frot cr f. fswap cr f. cr f.
cr cr
cr 1e 2e 3e df11/dx f.
cr 1e 2e 3e df11/dy f.
cr 1e 2e 3e df11/dz f.
cr cr
: dl_f12() ( dl: x -- 2^x) 2e f>dl dlswap dl^ ;
: der_f12() ( f: x -- 2^x) \ df12/dx = ln(2) * 2^x
     2e fswap f** 2e fln f*
;
cr 2e der dl_f12() f.
cr 2e der_f12()    f.
cr
cr 1e der dl_f12() f.
cr 1e der_f12() f.
cr cr
: dl_f13() ( dl: x -- x^2) 2e f>dl dl^ ;
: der_f13() ( f: x -- x^2) \ df12/dx = 2*x
     2e f*
;
cr 2e der dl_f13() f.
cr 2e der_f13()    f.
cr
cr 1e der dl_f13() f.
cr 1e der_f13() f.
cr cr
cr cr cr
[then]
Here the code ends.
Ahmed