Multivariable function differentiation

[Rafid]

Hi, I'm building a utility an application to allow symbolic differentiation of multivariable functions. Right now I'm testing CasADi with some general functions to compare them to the originals.

I'm able to define symbols and take the derivative using f.grad(), but how can I evaluate the derivative for different values for the all the variables?

This is how I'm defining my function: f = SXFunction(, )

thanks

[Joel Andersson]

So if I understand it correctly, you want to simultaneously calculate d_e/d_a, d_e/d_b, d_e/d_c, d_e/d_c, d_e/d_d?

The simplest way to do this is simply to group the inputs together:

f = SXFunction([vertcat([a,b,c,d])], [e])
f.init()
gf = f.grad()
de_da = gf[0,0]
de_db = gf[1,0]
de_dc = gf[2,0]
de_dd = gf[3,0]

There are also more low-level ways of doing the same thing, as calling the f.evalSX() with appropriate arguments in order to make one adjoint sweep.

Does this answer your question?
Joel

[Rafid]

Thanks Joel! That was actually going to be my next question. Sorry for being so ambiguous.

What I'm doing right now is trying to compare derivatives computed by CasADi to common derivatives that the utility is meant to handle. Since CasADi does not simplify the expressions(or does it and I just don't know how to?) I thought of evaluating the derivatives at different values for both equations to compare the values.

So let's say I have the following equation:

q = (1-Vu^0.5-(1-V)u)*qmax -> dq/du = qmax(-0.5Vu^(-0.5)+V-1)

Now I want to enter q in python for so to differentiate it using CasAdi and then I want to be able to give different values for variables V,u and qmax to be able to compare it to the original derivative.

Thanks

[Rafid]

Also, is it possible to just get the result and not a matrix with the result?

[Joel Andersson]

Hello!

You are right that CasADi only does limited simplifications of expressions (things like multiplications by 0 and 1 and the like) as well as constant foldings. It does not do things like common subexpression elimination or grouping linear terms together (if yoy are using the just-in-time compiler or c code generation, the compiler will likely do such things, however. Anyway, because the way symbolic expressions are represented is different from CAS-tools like Maple, Matlab Symbolic Toolbox etc. CasADi can usually avoid the need for simplifications altogether. This is especially true if your expression is constructed by differentiating expressions.

It is very important that you note that output from printing an expression to the screen does not give an accurate picture of how the expression is calculated. In particular, common subexpressions are expanded *during the printing*. If you want to see how the expression is actually calculated, you should print the function (SXFunction instance) or generate C-code.

As for your example, I guess you first create expressions for u, V and qmax:

u = ssym("u")
V = ssym("V")
qmax = ssym("qmax")

Then you use these expressions to form q:

q = (1-V*u**0.5-(1-V)*u)*qmax

Then calculate dq/du, for example through (if you need several derivatives, make sure to calculate them simultaneously as noted earlier to make use of the "cheap gradient rule" from automatic differentiation):

du_dq = gradient(q,u)

Finally form a function for calculating du/dq:

du_dq_fcn = SXFunction([u,V,qmax],[du_dq])
du_dq_fcn.init()

This function object can either be evaluated (numerically or symbolically) or you can use it to generate C-code, which is a good way of visualizing the actual calculations that take place in the virtual machine:

du_dq_fcn.generateCode("du_dq.c")

Finally, about your question "Also, is it possible to just get the result and not a matrix with the result?", it is possible, sure. You will find information about this in the documentation. For example you can simply convert it to numpy containers. But, and this is very important, CasADi is designed with an "everything is a matrix" philosophy in mind, just like Matlab. So the idea is that the user always work with SXMatrix (i.e. Matrix<SX>) in this case, also for scalars and vectors.

I hope this helps!
Joel

[Joel Andersson]

It looks like part of the code above disappeared. The second to last piece of code should have been:

du_dq_fcn = SXFunction([u,V,qmax],[du_dq])
du_dq_fcn.init()

[Joel Andersson]

One more try, with some extra spaces added:

du_dq_fcn = SXFunction (  [  u ,  V , qmax ]  , [ du_dq ] )
du_dq_fcn.init()

[Rafid]

Thanks for a very detailed and fast answer! That clears things up for me

[Joel Andersson]

I'm glad to be able to help!