Maxima: Piecewise function differentiation and simplification.

On 7/31/08, Constantine Frangos <cfrangos at telkomsa.net> wrote:

>  (%i69) vpen : penalty(u1^2*sin(x3)^2+u1^2*cos(x3)^2,0,vsqmax);
>
>  (%o69) if u1^2*sin(x3)^2+u1^2*cos(x3)^2 > vmax
>            then (u1^2*sin(x3)^2+u1^2*cos(x3)^2-vmax)^2
>            elseif u1^2*sin(x3)^2+u1^2*cos(x3)^2 < 0
>            then (u1^2*sin(x3)^2+u1^2*cos(x3)^2)^2 else 0
>
>  I want to differentiate vpen wrt u1 and then wrt x3 directly, dvpendu1 :
>  diff(vpen,u1), dvpendx3 : diff(vpen,x3), but could not do this correctly.

Assuming that F(if a then x elseif b then y elseif ... else z) is
equal to if a then F(x) elseif b then F(y) elseif ... else F(z),
we can do the following. (There are cases for which distributing
over the branches of a conditional doesn't yield the correct result.)

ifmap (f, e) := if not atom (e) and op (e) = "if" then apply ("if",
join (odd_args (e), map (f, even_args (e)))) else e;
:lisp (defun $even_args (e) `((mlist) ,@(odds (rest e) 0)))
:lisp (defun $odd_args (e) `((mlist) ,@(odds (rest e) 1)))

ifmap (foo, if a then x elseif b then y else z);
 => if a then foo(x) elseif b then foo(y) else foo(z)

ifmap (lambda ([e], diff (e, u1, 1, x3, 1)), vpen);
=> if u1^2*sin(x3)^2+u1^2*cos(x3)^2 > vmax then 0
           elseif u1^2*sin(x3)^2+u1^2*cos(x3)^2 < 0 then 0
           else 0

I think that's correct --- because of the sin and cos terms in vpen,
the derivative has terms which cancel out.

>  penalty(fx,amin,amax) := block(
>
>  [pen],
>
>  'if (fx > amax) then
>  pen : (fx - amax)^2
>  elseif (fx < amin) then
>  pen : (fx - amin)^2
>  else
>  pen : 0);

For what it's worth, this can be phrased more briefly as

penalty (fx, amin, amax) :=
  if (fx > amax)
    then (fx - amax)^2
  elseif (fx < amin)
    then (fx - amin)^2
  else 0;

I'm generally in favor of omitting unneeded stuff; it makes it easier
to see what' s going on.

best

Robert Dodier