Re: Differentiable Forth

I have noticed that my definition of dl^ if dual_numbers.fs is not good.
It doesn't get the good values for
a b dl^ when a = 0. ( a and b are dual numbers).
Because the definition is based on the use of dlln that uses fln and
fln(0) is -inf.
for example, with gforth I got
1e 0e 0e 0e dl^ dl. 1. + eps 0.  ok
5e 0e 0e 0e dl^ dl. 1. + eps 0.  ok
5e 1e 0e 0e dl^ dl. 1. + eps 0.  ok
0e 1e 0e 0e dl^ dl. NaN+ eps NaN ok , here must be 1e + eps 0e
0e 0e 0e 0e dl^ dl. NaN+ eps NaN ok ,      must be 1e + eps 0e
0e 0e 2e 0e dl^ dl. 0. + eps NaN ok,       must be 0e + eps 0e
I modified it and the new definition is:
: dl^  ( dl: a b -- a^b) \ a^b = exp(b*ln(a))
     dldup  dl>rl f0=  if dldrop dldrop 1e 0e   exit then
     dlover dl>rl f0<> if dlswap dlln dl* dlexp exit then
     dldrop dldrop 0e 0e
;
which corrects the results.
I also added the definition for dlabs for the absolute value:
: dlabs   ( f: a b -- c d) fswap fdup fabs fswap fdup f0<> if fdup fabs
f/ then frot f* ;
An example of usage:
: dl_f3() ( d: x -- y) \ y = |x^2-5*x+1|
dldup 2e f>dl dl^ dlswap 5e f>dl dl* dl- 1e f>dl dl+ dlabs ;
and some calculations:
0e f>dl dl_f3() dl. 1. + eps 0.  ok
0e f>dl_d dl_f3() dl. 1. + eps -5.  ok
1e f>dl_d dl_f3() dl. 3. + eps 3.  ok
2e f>dl_d dl_f3() dl. 5. + eps 1.  ok
2.5e f>dl_d dl_f3() dl. 5.25 + eps -0.000000000000000888178419700125  ok
3e f>dl_d dl_f3() dl. 5. + eps -1.  ok
4e f>dl_d dl_f3() dl. 3. + eps -3.  ok
4.5e f>dl_d dl_f3() dl. 1.25 + eps -4.  ok
4.9e f>dl_d dl_f3() dl. 0.510000000000005 + eps 4.8  ok
21e fsqrt 5e f+ 2e f/ f. 4.79128784747792  ok here, dl_f2() crosses th x
axis (dl_f3()=0)
notice how the eps part changes it sign
4.79128784747792e f>dl_d dl_f3() dl. 0.0000000000000035527136788005 +
eps -4.58257569495584  ok
4.79128784747793e f>dl_d dl_f3() dl. 0.000000000000049737991503207 + eps
4.58257569495586  ok
And I also tried to do derivation with respect to parameters.
Here an example of use:
include autodiff.fs
\ 3 variables x, y and z. The array xyz
3 dlarray xyz
2 dlarray params
\ x = 0.5, y = and z= -2.3
0.5e 0e xyz 0 dlarray!
4e 0e xyz 1 dlarray!
-2.3e 0e xyz 2 dlarray!
\ p1 = 2, y = and p2= 3
2e 0e params 0 dlarray!
3e 0e params 1 dlarray!
\ 1 func of  var
variable func_count
: addfunc >r 2dup r> -rot func_count @ funcarray! 1 func_count +! ;
: end_addfuncs 2drop 0 func_count ! ;
1 funcarray f1()
0 func_count !
dlvariable x
dlvariable y
dlvariable z
dlvariable p1
dlvariable p2
f1()
:noname ( dl: p1 p2 -- r) \ f(x,y,z; p1,p2) =
sin(p1*x^(y+z))−3*z*ln(p2*x^2*y^3)
   p2 dl! p1 dl!
   x dl@ y dl@ z dl@ dl+ dl^ p1 dl@ dl* dlsin
   x dl@ 2e f>dl dl^ y dl@ 3e f>dl dl^ dl* p2 dl@ dl* dlln z dl@ dl* 3e
f>dl dl*
   dl-
;
addfunc
end_addfuncs
xyz 0 dlarray@  x dl!
xyz 1 dlarray@  y dl!
xyz 2 dlarray@  z dl!
cr params f1() .jacobian
\ verif
( verification)
fvariable x_ 0.5e x_ f!
fvariable y_ 4e y_ f!
fvariable z_ -2.3e z_ f!
fvariable p1_ 2e p1_ f!
fvariable p2_ 3e p2_ f!
: df/dp1 x_ f@ y_ f@ z_ f@ f+ f** fdup p1_ f@ f* fcos f* ; \ df/dp1 =
x^(y+z)*cos(p1*x^(y+z))
cr df/dp1 f.
: df/dp2 z_ f@ -3e f* p2_ f@ f/   ; \ df/dp2 = -3*z/p2
cr df/dp2 f.
cr cr .( done) cr cr
Ahmed