Maxima - differentiating a piecewise function

问题:

Suppose you have a function defined by intervals, such as

f(x):=block(if x<0 then x^2 else x^3);

When we differentiate it with

diff(f(x),x);

we get

d/dx (if x<0 then x^2 else x^3)

whereas I'd like to get

(if x<0 then 2*x else 3*x^2)

Is there a way to obtain such result?

回答1:

This may help in a simple case:

(%i1) f(x):= charfun(x<0)*x^2 + charfun(x>=0)*x^3$

(%i2) gradef(charfun(y), 0)$

(%i3) diff(f(x),x);
                                           2
(%o3)              2 x charfun(x < 0) + 3 x  charfun(x >= 0)

charfun, gradef

You can try also Pw.mac package from Richard Hennessy.

回答2:

Here's a different approach using a simplification rule for "if" expressions. The unsolved part here is to detect discontinuities and generate delta functions for those locations. If you want to ignore those, you can define FOO to return 0. Note that I didn't attempt to implement the function discontinuities; that part is unsolved here. I can give it a try if there is interest.

(%i1) display2d : false $
(%i2) matchdeclare ([aa, bb, cc], all, xx, symbolp) $
(%i3) 'diff (if aa then bb else cc, xx) $
(%i4) tellsimpafter (''%, apply ("if", [aa, diff (bb, xx), true, diff (cc, xx)]) + FOO (aa, bb, cc, xx)) $
(%i5) FOO (a, b, c, x) := 'lsum ((ev (c, x = d) - ev (b, x = d)) * delta (d, x), d, discontinuities (a, x)) $
(%i6) diff (if x > 0 then x^2 else x^3, x);
(%o6) (if x > 0 then 2*x else 3*x^2)+'lsum((d^3-d^2)*delta(d,x),d,
                                           discontinuities(x > 0,x))

回答3:

Building on slitinov's answer I wrote this quite naive implementation for functions with more than two "pieces":

gradef(charfun(dummy),0)$

/* piecewise function definition */
itv: [[x<0],[x>=0,x<1], [x>=1]]; /* intervals */
fi:  [ 1,    x^2+1,      2*x  ]; /* local functions */

/* creation of global function f and its derivative df */
f:0;
for i: 1 thru 3 do f:f+charfun(apply("and",itv[i]))*fi[i];
df:diff(f,x);

/* display of local functions and derivatives */
for i: 1 thru 3 do (
  apply(assume,itv[i]),
  newline(),
  print(itv[i]),
  print("f = ",ev(f)),
  print("df = ",ev(df)),
  apply(forget,itv[i])
  );

plot2d([f,df],[x,-2,3],[y,-1,5],[style,[lines,4,3],[lines,2,2]]);