solve and plot nonlinear system of equations with no explicit solution

I have these a nonlinear system of equations:

phi and psi are constants and will be taken as input from the user.

I tried using fsolve but i got the message:

"No solution found.

fsolve stopped because the problem appears regular as measured by the gradient,

but the vector of function values is not near zero as measured by the

value of the function tolerance."

When I try to solve using vpasolve i get the message:

"Error using mupadengine/feval_internal

More equations than variables is only supported for polynomial systems."

And when i try to use "solve" , then i get the message:

"Warning: Unable to find explicit solution" and empty sym.

Can you help me to solve this? this is the code i have when i tried fsolve

phi = input('Enter phi value: ');
psi = input('Enter psi value: ');

f = @(y) [(-1/2)*(cos(y(2))*cos(y(3))+sin(y(1))/psi);(phi*sin(y(1)))-(1/2)*cos(y(2))*sin(y(3)); phi*(cos(y(1)))+(1/2)*sin(y(2))*cos(y(3))];
fsolve(f, [1 2 3])

And the code i tried using vpa solve:

phi = input('Enter phi value: ');
psi = input('Enter psi value: ');

syms x z theta 

eqn1 = -(1/2)*(cos(x)*cos(z)+sin(theta)/psi) == 0;

Stheta = solve(eqn1,theta)

eqns = [phi*sin(Stheta)-(1/2)*cos(x)*sin(z) == 0, phi*cos(Stheta) + (1/2)*sin(x)*cos(z) == 0];
[tout,S] = vpasolve(eqns,[x z])

Answer 1

[David Goodmanson]

HI Cengizhan,

MODIFIED

The first equation above is

(-1/2)*(cos(x)*cos(z)+sin(theta)/psi) = 0    (1)

which has an extraneous factor of (-1/2) and might as well be written as

(cos(x)*cos(z)+sin(theta)/psi) = 0 

How confident are you that this is the correct equation? When I posted this answer intially, without thinking about it I assumed the equation was

(-1/2)*(cos(x)*cos(z)) +sin(theta)/psi = 0    (2)

because this form is the same as that of the remaining two equations. (It did not help that the equations were posted as an image rather than a text file and could not be directly copied).

If you start dividing equations by hand and so forth, for (2) you can arrive at

psi = 6
phi = 1/3

lam = phi*psi;
lam2 = lam^2;
sig2 = (psi^2/4)/(1+lam2);
z = atand(lam)
theta = atand(sqrt((sig2-lam2)/(1-sig2)))
x = atand(-lam/tand(theta))

% results
% z = 63.4349
% theta = 58.9091
% x = -50.3360

% check, should be small
cosd(z)*cosd(x) - 2*sind(theta)/psi
cosd(x)*sind(z) - 2*phi*sind(theta)
cosd(z)*sind(x) + 2*phi*cosd(theta)

For (1), a solution is

psi = 6
phi = 1/3
lam = phi*psi;
lam2 = lam^2;
sig2 = (psi^2)/(1+4*lam2);
z = atand(-2*lam)
theta = -atand(sqrt((sig2-4*lam2)/(1-sig2)))
x = atand(2*lam/tand(theta))

% results
% z = 75.9638
% theta = -74.1592
% x = -48.6171

% check, should be small
cosd(z)*cosd(x) + sind(theta)/psi
cosd(x)*sind(z) - 2*phi*sind(theta)
cosd(z)*sind(x) + 2*phi*cosd(theta)

With the use of atand, all angles are -90<angle<90 and the cosines of the angles are all positive.

For the x = pi/2 solution, put cos(x) = 0 into the equations, which forces sin(theta) = 0. Then sin(x) = +-1 and cos(theta) = +-1. So it appears that z = acos(+-2*phi), depending how some of the signs are chosen. Chasing around minus signs for some of the special cases can get tedious.

The algebra above uses atan and I can't really comment on the acos solutions in your comment below. Except that y = acos(1/sqrt(1+tan(y)^2)) where the square root can be of either sign, and algebraic expressions for tangent can be plugged into this.

There are only so many values of phi and psi that result in real angles. For many values of phi and psi the solution is complex. There are regions of the phi-psi plane that don't lead to real angles, and you can't "solve" that issue.

[Cengizhan Demirbas]

Hello, thanks for your answer. I divided the equations by hand and got these (slightly different than yours):

psi = 2

phi = 0.4

and plugging these into first equation:

When i use solve function to solve this, i get 3 results with a "warning: possibly spurious solutions". Solutions are:

x1 = pi/2

x2 = -asin((4*11^(1/2)*39^(1/2))/195)

x3 = asin((4*11^(1/2)*39^(1/2))/195)

What exactly is happening? I think x2 and x3 do not have real solutions for high phi and psi, and that causes a problem. If so, how can i solve this issue?

I am trying to confirm the results of a paper. In it, it finds 3 solutions correspondingly:

x = pi/2 (yay!) z = acos(-2*phi) (how?)

x = acos(gamma^(1/2)) z = atan(-2*lam) (yay again!)

x = acos(-gamma^(1/2)) z = atan(-2*lam)

How can i get these? Again, thanks a lot for your answer.

[David Goodmanson]

Hi Cengizhan, I posted a modified answer.

[Cengizhan Demirbas]

Thanks. Yes, I am certain the equations are correct. They were initially differential equations, but i am looking for a case in which d(theta,x,z)/dt = 0, so i forgot to simplify the equations. Can you tell me, how did you find sig2 and theta? I couldn't manage to write theta independent from x.

I noticed that when i use phi psi as symbols rather than give values, using solve x2 and x3 becomes

x = (+-) asin(phi*psi*(-(4*phi^2*psi^2 - 3)/(4*phi^2*psi^2 - 1))^(1/2))

Since it is sqrt(-...) i have to deal with complex numbers.

To get gamma i mentioned earlier, i used this equation and found

x = acos(((4*phi^4*psi^4+phi^2*psi^2-1)/(4*phi^2*psi^2-1))^(1/2))

Resembles the equation in paper but not equal.

Is there a possibility that solve function fails to perform when complex numbers are involved? Or am i somehow doing the math wrong? Thanks for all the help

[David Goodmanson]

Since the argument of the square root is not always negative for all phi and psi, there are some regions in the phi-psi plane that work. There are exclusion regions that don't.

I don't want to delve into acos, but I can tell you how I found the atan equations, anyway.

starting point
cosd(z)*cosd(x) = -sind(theta)/psi          (3)
cosd(x)*sind(z) =  2*phi*sind(theta)        (4)
cosd(z)*sind(x) = -2*phi*cosd(theta)        (5)

divide (4) by (3) to obtain

z = atand(-2*lam)    where lam = phi*psi

divide (5) by (3) to obtain

tand(x) = 2*lam/tand(theta)    x = atand(2*lam/tand(theta))   (6a, 6b) 

square (3), invert both sides, use 1/cos^2(y) = 1+tan^2(y) twice, do some algebra, obtain

theta = -atand(sqrt((sig2-4*lam2)/(1-sig2)))

where sig2 is shown in the answer. In the last equation there are two choices for the sign of the square root. The minus sign (which passes through atand) works for the example case I had.