I wrote:
>
> If you mean support for derivatives of "jumpy" functions, this
> is reasonable easy. Namely, we need to add 'signum' and
> 'diracDelta' to standard operators and define new functions
> (this is very much like current definition of 'abs'). Then
> we need to define derivatives:
>
> D(abs(x), x) = signum(x)
> D(signum(x), x) = diracDelta(x)
>
> and derivatives of 'diracDelta' would remain in symbolic
> form. Such definitions are sound and cause no special
> trouble.
In the attachement there is a patch in this direction.
Simple things work, but I did only very light testing.
Up to now I have noticed two problem. First, I used
name 'sign' to be consistent with other uses. But
this means that 'sign' is heavily overloaded and
interpreter has problems with choosing correct
signature. In fact, there was also compilation
problem: in DoubleFloat there is implementation
of 'sign : % -> Integer'. But with the patch we
also have 'sign : % -> %' and compiler was picking
the second signature and report error.
Second problem is that after changing definition of
derivative of 'abs' several integrals containing
'abs' which used to work now fail. In a sense this
is trivial problem: our integrator can integrate
only derivatives of a few simple forms containg 'abs'
and changing derivative of 'abs' changed set of
"integrable" functions. But still, this is not
nice.
I am thinking about including (possibly improved)
version of this patch. Comments?
--
Waldek Hebisch
[email protected]
--
You received this message because you are subscribed to the Google Groups
"FriCAS - computer algebra system" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to [email protected].
To post to this group, send email to [email protected].
Visit this group at http://groups.google.com/group/fricas-devel.
For more options, visit https://groups.google.com/d/optout.
Index: src/algebra/trigcat.spad
===================================================================
--- src/algebra/trigcat.spad (revision 1791)
+++ src/algebra/trigcat.spad (working copy)
@@ -261,6 +261,10 @@
SpecialFunctionCategory() : Category == with
abs : % -> %
++ abs(x) returns the absolute value of x.
+ sign : % -> %
+ ++ sign(x) returns the sign of x.
+ diracDelta : % -> %
+ ++ diracDelta(x) is unit mass at zeros of x.
Gamma : % -> %
++ Gamma(x) is the Euler Gamma function.
Beta : (%, %)->%
Index: src/algebra/sf.spad
===================================================================
--- src/algebra/sf.spad (revision 1791)
+++ src/algebra/sf.spad (working copy)
@@ -570,7 +570,7 @@
x = ((n := wholePart x)::%) => n
"failed"
- sign(x) == retract FLOAT_-SIGN(x, 1)$Lisp
+ sign(x : %) : Integer == retract FLOAT_-SIGN(x, 1)$Lisp
abs x == abs_DF(x)$Lisp
Index: src/algebra/expr.spad
===================================================================
--- src/algebra/expr.spad (revision 1791)
+++ src/algebra/expr.spad (working copy)
@@ -302,6 +302,8 @@
acsch x == acsch(x)$EF
abs x == abs(x)$FSF
+ sign x == sign(x)$FSF
+ diracDelta x == diracDelta(x)$FSF
Gamma x == Gamma(x)$FSF
Gamma(a, x) == Gamma(a, x)$FSF
Beta(x, y) == Beta(x, y)$FSF
Index: src/algebra/combfunc.spad
===================================================================
--- src/algebra/combfunc.spad (revision 1791)
+++ src/algebra/combfunc.spad (working copy)
@@ -485,6 +485,10 @@
++ error if op is not a special function operator
abs : F -> F
++ abs(f) returns the absolute value operator applied to f
+ sign : F -> F
+ ++ sign(x) returns the sign of x.
+ diracDelta : F -> F
+ ++ diracDelta(x) is unit mass at zeros of x.
Gamma : F -> F
++ Gamma(f) returns the formal Gamma function applied to f
Gamma : (F, F) -> F
@@ -685,6 +689,8 @@
SPECIALINPUT ==> '%specialInput
iabs : F -> F
+ isign : F -> F
+ iDiracDelta : F -> F
iGamma : F -> F
iBeta : (F, F) -> F
idigamma : F -> F
@@ -708,6 +714,8 @@
opabs := operator('abs)$CommonOperators
+ opsign := operator('sign)$CommonOperators
+ opDiracDelta := operator('diracDelta)$CommonOperators
opGamma := operator('Gamma)$CommonOperators
opGamma2 := operator('Gamma2)$CommonOperators
opBeta := operator('Beta)$CommonOperators
@@ -740,6 +748,8 @@
op_erfis := operator('%erfis)$CommonOperators
abs x == opabs x
+ sign x == opsign x
+ diracDelta x == opDiracDelta x
Gamma(x) == opGamma(x)
Gamma(a, x) == opGamma2(a, x)
Beta(x, y) == opBeta(x, y)
@@ -1685,6 +1695,8 @@
operator op ==
is?(op, 'abs) => opabs
+ is?(op, 'sign) => opsign
+ is?(op, 'diracDelta) => opDiracDelta
is?(op, 'Gamma) => opGamma
is?(op, 'Gamma2) => opGamma2
is?(op, 'Beta) => opBeta
@@ -1761,6 +1773,22 @@
smaller?(x, 0) => kernel(opabs, -x)
kernel(opabs, x)
+ derivative(opabs, (x : F) : F +-> sign(x))
+
+ isign x ==
+ zero? x => 0
+ is?(x, opabs) => 1
+ is?(x, opsign) => x
+ smaller?(x, 0) => kernel(opsign, -x)
+ kernel(opsign, x)
+
+ derivative(opsign, (x : F) : F +-> (2::F)*diracDelta(x))
+
+ evaluate(opsign, isign)$BasicOperatorFunctions1(F)
+
+ iDiracDelta x == kernel(opDiracDelta, x)
+
+ evaluate(opDiracDelta, iDiracDelta)$BasicOperatorFunctions1(F)
iBeta(x, y) == kernel(opBeta, [x, y])
idigamma x == kernel(opdigamma, x)
iiipolygamma(n, x) == kernel(oppolygamma, [n, x])
@@ -2369,7 +2397,6 @@
derivative(op_erfs, d_erfs)
derivative(op_erfis, d_erfis)
- derivative(opabs, (x : F) : F +-> abs(x)*inv(x))
derivative(opGamma, (x : F) : F +-> digamma(x)*Gamma(x))
derivative(op_log_gamma, (x : F) : F +-> digamma(x))
derivative(opBeta, [iBetaGrad1, iBetaGrad2])
