maxima: evaluating Jacobian at a point

Question

I'm not very good with maxima, however I'm trying to learn.
I am attempting to write a Newton iterative solver for a fairly small nonlinear system. In order to do this I must evaluate the Jacobian at the current iteration. But, I can't seem to figure out a way to evaluate the Jacobian at all.

I currently have:

F1000(D5,D6,D7,D8,D9,D10,b1,b2,b3,b4,b5,b6,b7,b8,b9,d1,d2,d3,e1,e2,e3,f1,f2,f3) 
     := (-f1-3*d1-9*b1+60)*l1$
F1001(D5,D6,D7,D8,D9,D10,b1,b2,b3,b4,b5,b6,b7,b8,b9,d1,d2,d3,e1,e2,e3,f1,f2,f3) 
     := (-f2+e2-6*d2-3*b2-15)*l1$
/*further functions omitted*/

J : jacobian([ 
   F1000(D5,D6,D7,D8,D9,D10,b1,b2,b3,b4,b5,b6,b7,b8,b9,d1,d2,d3,e1,e2,e3,f1,f2,f3),
   F1001(D5,D6,D7,D8,D9,D10,b1,b2,b3,b4,b5,b6,b7,b8,b9,d1,d2,d3,e1,e2,e3,f1,f2,f3),
   /*functions omitted*/
   F0018(D5,D6,D7,D8,D9,D10,b1,b2,b3,b4,b5,b6,b7,b8,b9,d1,d2,d3,e1,e2,e3,f1,f2,f3),
   F0019(D5,D6,D7,D8,D9,D10,b1,b2,b3,b4,b5,b6,b7,b8,b9,d1,d2,d3,e1,e2,e3,f1,f2,f3)
   ],
   [D5,D6,D7,D8,D9,D10,b1,b2,b3,b4,b5,b6,b7,b8,b9,d1,d2,d3,e1,e2,e3,f1,f2,f3]
   )$

But, now i cannot find a good way to evaluate J at, for example,

u_init : [D5_0,D6_0,D7_0,D9_0,...,e3_0,f1_0,f2_0,f3_0]$

I've tried things like:

subs(J, [D5,D6,D7,D8,D9,D10,b1,b2,b3,b4,b5,b6,b7,b8,b9,d1,d2,d3,e1,e2,e3,f1,f2,f3],
     [1,1,1,1...,1,1,1,1]),

and

ev(J, [D5,D6,D7,D8,D9,D10,b1,b2,b3,b4,b5,b6,b7,b8,b9,d1,d2,d3,e1,e2,e3,f1,f2,f3],
     [1,1,1,1...,1,1,1,1]),

and I just cannot find a good way to evaluate the Jacobian, $J$, at the point $u_i$. I really appreciate any clarifications/insights.

Thanks

What language is this? What are `subs` and `ev`? I'd ask about `jacobian` as well, but your question indicates that's the entire point. — Teepeemm, Sep 14 '13 at 00:40
@Teepeemm, I apologize if i was unclear. My question is directly related to the language Maxima (or WXMaxima). `subs` is a command from matlab that i thought i should try; i believe `subst` is a similar command for maxima that i could not get to work. `ev` is a command from maxima for evaluating `ev(expr, [variables], [values])` but i've saved J as the jacobian command and `jacobian` does not work if trying to differentiate with respect to constants; nor should it. If you still have a question about jacobian matrices, I would gladly try to help. — useslessone, Sep 14 '13 at 00:55
Sorry, I thought you meant "maxima" in the sense of using Newton's Method to find a local maxima, not in the sense of the programming language. I am familiar with Jacobians. Newton's method in this case is equation x_{n+1} = x_n - (Jinv)(x_n) * F(x_n), where Jinv is the inverse of the Jacobian matrix. But now I'm not sure where your difficulty lies: are you having trouble finding the Jacobian matrix, its inverse, or its application to the current situation? — Teepeemm, Sep 14 '13 at 22:21
No worries, My trouble lies in the fact that I can find the Jacobian, J, but it is an expression with over 24 variables in it. I cannot figure out a way to evaluate J(x_n).. that is to say i have J(), i don't know how to plug in numerical values in WXMaxima. — useslessone, Sep 15 '13 at 00:09
It seems like that means your question isn't really about math or newton's method, but is instead about evaluating a function in maxima. Which takes it from my specialty to something I know nothing about. Sorry. — Teepeemm, Sep 16 '13 at 18:57

Roberto · Answer 1 · 2013-10-06T17:05:53.913

Try

subst([D5=D5_0,D6=D6_0,...,f3=f3_0], J);

An example for a smaller Jacobian:

(%i1) F1(x,y):=x^2+y;
                                           2
(%o1)                         F1(x, y) := x  + y
(%i2) F2(x,y):=-x-y^2;
                                                2
(%o2)                        F2(x, y) := - x - y
(%i3) J:jacobian([F1(x,y),F2(x,y)],[x,y]);
                                [ 2 x    1   ]
(%o3)                           [            ]
                                [ - 1  - 2 y ]
(%i4) subst([x=2,y=3],J);
                                 [  4    1  ]
(%o4)                            [          ]
                                 [ - 1  - 6 ]

maxima: evaluating Jacobian at a point

1 Answers1