I want to point out that another problem with doing this in SymPy is that 
SymPy currently does not have a way to represent a derivative evaluated at a 
point, meaning that it can't for example compute the chain rule correctly for 
something like diff(f(g(x)), x).  See issue 1620.  

The reason I point this out is that it tries to compute this so that it "looks" 
like the Euler/Lagrange style, because instead of evaluating at the point g(x), 
it "takes the derivative with respect to" g(x) (even though this is both 
incorrect and will cause SymPy to fail if you try to enter it manually):

In [48]: diff(f(g(x)), x)
Out[48]: 
  d            d       
─────(f(g(x)))⋅──(g(x))
dg(x)          dx      

I am still trying to fully understand the code below, but if one of the subs is 
trying to do something like diff(f(x), x).subs(x, g(x)), then you will get an 
error in SymPy for the above reason.

Aaron Meurer

On Oct 31, 2010, at 11:19 AM, Philippe wrote:

> thanks a lot !!
> I will try to convert that to sympy.
> 
> On 31 oct, 18:14, Tim Lahey <[email protected]> wrote:
>> There's EulerLagrange code in Maple, but in general, you can't take
>> partial derivatives with respect to a function. There's a trick to handling
>> it, though.
>> 
>> My Maple code to do EulerLagrange (written when the built-in routine
>> was awful) is below. Consider it to be GPLv3. I'm not releasing it as
>> BSD since there are commercial add-ons to Maple I'd prefer not using
>> the code.
>> 
>> It should be relatively straightforward to map to sympy. Note the map
>> and map2 are because they map over different input parameters.
>> 
>> The basic trick is to convert the function to a symbol take the derivative
>> with respect to that symbol and convert back.
>> 
>> Handling general EulerLagrange is tricky, but this code does the basic
>> case.
>> 
>> Hope this helps. This function is extremely fast in Maple and I think is 
>> still
>> faster than the built-in routines. When I wrote it (around Maple 9/10/11) it
>> was two orders of magnitude faster than Maple's routine.
>> 
>> Cheers,
>> 
>> Tim.
>> 
>> # Calculate the Mass (d/dt(dL_dqv)) and Stiffness (dL_dq) contributions to 
>> the
>> # Euler-Lagrange equation. Note that the Stiffness contribution does
>> not include the negative sign.
>> EulerLagrange := proc(Lagrangian::anything, variables::list)
>>         local num_list, qv_name, vel_var, qv_subs, qv_unsubs, Lagrange_subs1,
>> Lagrange_subs2, dL_dqv1, dL_dqv2, dL_dqv, dL_dqvt, dL_dq,       dL_dq1,
>> dL_dq2, dL_dq3, q_name, q_subs, q_unsubs:
>> 
>>         # create a list of indices from 1 to the number of variables used in
>> the formulation
>>         num_list := [seq(i,i=1..nops(variables))]:
>> 
>>         # Define a list of generalized velocity and position variables
>>         qv_name := map2(cat,qv,num_list):
>>         q_name := map2(cat,q,num_list):
>> 
>>         # Equate the time derivatives of the system variable to the
>> generalized velocities
>>         # also define the reverse mapping
>>         vel_var := map(diff,variables,t):
>>         qv_subs := zip(equate,vel_var,qv_name):
>>         qv_unsubs := zip(equate,qv_name,vel_var):
>> 
>>         # Equate the generalized positions to the system variables and define
>> the reverse mapping
>>         q_subs := zip(equate,variables,q_name):
>>         q_unsubs := zip(equate,q_name,variables):
>> 
>>         # Convert the Lagrangian to the generalized position and velocity 
>> variables
>>         Lagrange_subs1 := subs(qv_subs,Lagrangian):
>>         Lagrange_subs2 := subs(q_subs,Lagrange_subs1):
>> 
>>         # Differentiate the Lagrangian with respect to the generalized
>> velocities and positions
>>         dL_dqv1 := map2(diff,Lagrange_subs2,qv_name):
>>         dL_dq1 := map2(diff,Lagrange_subs2,q_name):
>> 
>>         # Revert back to the system variables
>>         dL_dq2 := map2(subs,qv_unsubs,dL_dq1):
>>         dL_dqv2 := map2(subs,qv_unsubs,dL_dqv1):
>>         dL_dqv := map2(subs,q_unsubs,dL_dqv2):
>>         dL_dq := map2(subs,q_unsubs,dL_dq2):
>>         dL_dqvt := map(diff,dL_dqv,t):
>> 
>>         # Return the two components of the Euler-Lagrange Equation
>>         return (dL_dqvt, dL_dq):
>>     end proc:
>> 
>> --
>> Tim Lahey
>> PhD Candidate, Systems Design Engineering
>> University of Waterloohttp://www.linkedin.com/in/timlahey
> 
> -- 
> You received this message because you are subscribed to the Google Groups 
> "sympy" group.
> To post to this group, send email to [email protected].
> To unsubscribe from this group, send email to 
> [email protected].
> For more options, visit this group at 
> http://groups.google.com/group/sympy?hl=en.
> 

-- 
You received this message because you are subscribed to the Google Groups 
"sympy" group.
To post to this group, send email to [email protected].
To unsubscribe from this group, send email to 
[email protected].
For more options, visit this group at 
http://groups.google.com/group/sympy?hl=en.

Reply via email to