0

I have the following equation system (click to see picture) , and would like to solve for X(t), Y(t), Z(t), hopefully using Octave/Matlab, which I'm familiar with, but I would not mind solving it by any other means necessary.

Now, Fsolve is useful for regular equation systems, and Ode45, Lsode are useful for differential equations. But, what about this particular system? Note that the differential equation at the bottom contains not only Y, but, also, both X and Z, and those are dependent on the two non-differential equations above.

To be honest, I'm not quite sure how to approach a basic code to solve this system, and after spending some time thinking, I decided to ask for help. I really appreciate any guidance to solve this problem somewhat efficiently, if it is even possible. Almost any reply would be of some use to me right now.

Thunderboltwarrior
  • 5
  • 4
  • What is preventing you from applying these tools? `fsolve` or `fzero` for the second equation getting `x`, then the first equation (or directly with `acos`) to get `z` and then apply to the differential equation to get the derivatives array. – Lutz Lehmann Feb 15 '20 at 10:48
  • I can't solve for "X" on the second equation because it is dependent on "Y", and to get "Y", I would need to solve the differential equation. But the differential equation is also dependent on "X", making it a 2 unknown 2 equation non-linear problem, with some differential terms. Or at least thats how I'm understanding it right now. Appreciate the comment. – Thunderboltwarrior Feb 15 '20 at 13:33

1 Answers1

0

If you know y, you can solve for x, and this even unconditionally as the second equation in monotonous in x

x = fsolve(@(x) y^2-1413.7*x-1095.2*cos(x)+2169, 0)

Then once you know x, you can solve for z using the known inverse cosine

z = acos(0.20978-cos(x))

This may actually fail to give a result if cos(x) comes close to -1. One can artificially cut out that error, introducing a possibly wrong solution

z = acos(min(1,0.20978-cos(x)))

For simplicity, assemble these operations in a helper function

function [x,z] = solve_xz(y)
    x = fsolve(@(x) y^2-1413.7*x-1095.2*cos(x)+2169, 0);
    z = acos(min(1,0.20978-cos(x)));
end

Now use that to get the ODE for y

function dy = ode_y(t,y)
    [x,z] = solve_xz(y(1));
    dy = [ y(2); y(3); 6666.6667*(z-x)-333.3333*y(1)-33.3333*y(2)-5*y(3) ];
end

and apply the ODE solver of your choice. It is quite likely that the system is stiff, so ode45 might not be the best solver.

Lutz Lehmann
  • 25,219
  • 2
  • 22
  • 51