Path: csiph.com!usenet.pasdenom.info!gegeweb.org!eternal-september.org!feeder.eternal-september.org!mx02.eternal-september.org!.POSTED!not-for-mail
From: Joe Riel <joer@san.rr.com>
Newsgroups: comp.soft-sys.math.maple
Subject: Re: Differentiating with respect to an expression
Date: Tue, 11 Nov 2014 17:52:51 -0800
Organization: A noiseless patient Spider
Lines: 103
Message-ID: <87h9y5kwwc.fsf@san.rr.com>
References: <m3u1gm$7ir$1@news.albasani.net> <87sihpl61w.fsf@san.rr.com> <m3u6n6$kkr$1@news.albasani.net> <87lhnhkx98.fsf@san.rr.com>
Mime-Version: 1.0
Content-Type: text/plain; charset=us-ascii
Injection-Info: mx02.eternal-september.org; posting-host="15591ad2607da309a0d1a78a1d632bc7"; logging-data="22139"; mail-complaints-to="abuse@eternal-september.org";	posting-account="U2FsdGVkX1+PD3+HLd6OE06LA9+Jh55q"
User-Agent: Gnus/5.13 (Gnus v5.13) Emacs/23.4 (gnu/linux)
Cancel-Lock: sha1:GwBuqzUMSE+OlVpm7buB2UK9Fko= sha1:g4/KuEKI3ijF4GJHw35OiInR9eM=
Xref: csiph.com comp.soft-sys.math.maple:1016

Joe Riel <joer@san.rr.com> writes:

> rouben@shadow.(none) (Rouben Rostamian) writes:
>
>> In article <87sihpl61w.fsf@san.rr.com>, Joe Riel  <joer@san.rr.com> wrote:
>>>rouben@shadow.(none) (Rouben Rostamian) writes:
>>>
>>>> The following issue comes up quite often in the context of
>>>> analytical mechanics.  I have a clunky solution for it.
>>>> I wonder if there is a clever way.
>>>>
>>>> Let's say we have L = x^2 + x'^2 * x''^2, where x is a function of t,
>>>> and in the usual mathematical notation, x' and x'' are the first and
>>>> second derivatives of x.
>>>>
>>>> We want to compute the derivative of L with respect to x'.  The
>>>> answer should be 2*x'*x''^2.  Here is the way I do it in Maple:
>>>>
>>>> restart;
>>>> L := x(t)^2 + diff(x(t),t)^2 * diff(x(t),t,t)^2;
>>>> subs(diff(x(t),t,t)=Z2, diff(x(t),t)=Z1, %);
>>>> diff(%, Z1);
>>>> subs(Z1=diff(x(t),t), Z2=diff(x(t),t,t), %);
>>>>
>>>> The L shown above is simple enough so that we don't need a
>>>> CAS to compute the derivative.  The L in a real example will
>>>> be the result of a long chain of calculations, will depend on
>>>> several functions and their derivatives, and will take up a
>>>> couple of screenfuls.
>>>>
>>>> If there is a clever way to compute that derivative, I would
>>>> like to know.
>>>
>>>A low-level way to do this is with frontend:
>>>
>>>frontend(diff, [L,diff(x(t),t)]);
>>>       2*diff(x(t),t)*diff(diff(x(t),t),t)^2
>>>
>>>An alternative method is to use the Physics package:
>>>
>>>with(Physics):
>>>
>>>diff(L, diff(x(t),t));
>>>       2*diff(x(t),t)*diff(diff(x(t),t),t)^2
>>>
>>># alternatively, with Physics
>>>
>>>diff(L, D(x)(t));
>>>     2*D(x)(t)*`@@`(D,2)(x)(t)^2
>>>
>>
>> Thanks much, Joe, this is exactly what I was hoping for.
>> I had no idea about frontend() or the Physics package.
>>
>> This brings up a somewhat related question.  Doing
>>    frontend(int, [x(t)^2, x(t)]);
>> we get x(t)^3/3, which is fine.  The following, however,
>> issues an Error message:
>>     frontend(int, [x(t)^2, x(t)=a..b]);
>>
>> Is there a way to get that to work too?
>>
>> Rouben
>
> The optional third argument of frontend can be used for that.
> See the help page.  It consists of a list of two sets.
> The first set is significant here, it contains the types 
> that frontend does *not* freeze.  By default the `+`, `*`, and `^` 
> types are not frozen.  In this case, you also do *not* want
> to freeze `=` and `..`.  So the following works:
>
>   frontend(int, [x(t)^2, x(t)=a..b], [{`+`,`^`,`*`,`..`,`=`},{}]);
>        1/3*b^3-1/3*a^3
>
> Frequently one can use the Not type to specify what you do
> want frozen.  For example,
>
>   frontend(int, [x(t)^2, x(t)=a..b], [{Not(identical(x(t)))},{}]);
>        1/3*b^3-1/3*a^3

As mentioned, frontend is a rather low-level tool.
It has its uses, but is somewhat awkward to use.
A useful trick, when you're having trouble formulating
the appropriate set of types that will not be frozen,
is to use 'print' for the operator to see what has been frozen:


  frontend(print, [x(t)^2, x(t)=a..b]);                        
               O^2, O

  # unfreeze defaults, and equations
  frontend(print, [x(t)^2, x(t)=a..b], [{`+`,`*`,`^`,`=`},{}]);
              O^2, O = O

  # unfreeze defaults, equations, and ranges
  frontend(print, [x(t)^2, x(t)=a..b], [{`+`,`*`,`^`,`=`,`..`},{}]);
              O^2, O = a .. b

That last one is what we want to pass to int, so replacing
the 'print' with 'int' gives the desired result.

-- 
Joe Riel