Jacobian of non-linear ODE in MATLAB

Question

I'm calculating the jacobian for a two-body problem, which is defined like so

I have set up my system of equations as follows

syms y1(t) y2(t) y3(t) y4(t)

r = sqrt(y1^2 + y2^2)
y3 = diff(y1)
y4 = diff(y2)
yd = [y3; y4; -y1/r^3; -y2/r^3]

jacobian(yd, [y1 y2 y3 y4])

However, when I run the jacobian function I get the following error

The second argument must be a vector of variables.

What am I doing wrong?

EDIT:

I have also tried parametrizing y for t y(t) to no avail.

Answer 1

As the error message suggests that the second argument must be a vector of variables, whereas in your case it is: [y1, y2, 1, 1].

Also there is no need to initialize them as symfun class i.e. y1(t), y2(t), y3(t) and y4(t), you can define them as sym class instead i.e. y1, y2, y3 and y4

So, by initializing them as sym and removing the lines where you make y3 and y4 equal to 1, i.e.

syms y1 y2 y3 y4
r = sqrt(y1^2 + y2^2);
yd = [y3; y4; -y1/r^3; -y2/r^3];
jacobian(yd, [y1 y2 y3 y4])

you will get this output:

[                                                    0,                                                    0, 1, 0]
[                                                    0,                                                    0, 0, 1]
[ (3*y1^2)/(y1^2 + y2^2)^(5/2) - 1/(y1^2 + y2^2)^(3/2),                        (3*y1*y2)/(y1^2 + y2^2)^(5/2), 0, 0]
[                        (3*y1*y2)/(y1^2 + y2^2)^(5/2), (3*y2^2)/(y1^2 + y2^2)^(5/2) - 1/(y1^2 + y2^2)^(3/2), 0, 0]