Look carefullly at what is returned.
solution = solve(eqn_final, variable_final, "ReturnConditions", true);
solution
solution =
struct with fields:
x0: [3×1 sym]
y0: [3×1 sym]
omega0: [3×1 sym]
E0: [3×1 sym]
parameters: [1×1 sym]
conditions: [3×1 sym]
So there are THREE distinct solutions possible that it found.
solution.x0
ans =
-1
1
2*x - 1
So x0 is one of the unknowns. In the first two solutions, it is a fixed value.
Looking at the conditions return, we see why.
solution.conditions
ans =
Pm0 == 0
Pm0 == 0
in(x, 'real') & Pm0 ~= 0 & Pm0^2 + 4*x^4 == 4*x^3 & (in(Pm0, 'real') | x == 0 | x == 1)
Ah, when Pm0 is exactly zero, the problem reduces. But MATLAB does not know that Pm0 will not be zero. and when Pm0 is zero, the problem becomes considerably simpler.
So what happens when pm0 is some non-zero value? Anything but zero. Then we have solution number 3, the one you will care about.
So all seems simple. Just pick solution #2, and substiture the value of Pm0. But is it simple? Does your problem have a solution?
solution.x0(3)
ans =
2*x - 1
solution.y0(3)
ans =
(4*(- x^3 + x^2))/Pm0
solution.omega0(3)
ans =
0
solution.E0(3)
ans =
x
This tells us, that for ANY value of Pm0, there are infinitely many solutions. We can arbitrarily choose some real number for x. What is x? x just happens to be the value of E0. It seems you can arbitrarily choose the value for E0.
Could we see that by some simple manipulations of those equations? It gets complicated, because when we try this:
simplify(eqn2 - eqn1)
ans =
x0 + y0 == 0 | omega0 == 0
So that tells us the solution has EITHER omega0 == 0, or it reduces to a simple linear equation in terms of x0 and y0. If omega0 is non-zero, then what happens? Now we have a problem. No non-trivial solution will exist at all, because then we must have x0 == -y0.
That means we can assume omega0 must be zero. Then both of equations 1 and 2 effectively reduce to a unit circle in the parameters x0 and y0, thus the equation
x0^2 + y0^2 == 1
So a circle with center at the origin, and radius 1.
But then we have eqn4, which gives us x0 as a function of E0. And eqn3, which allows us to compute y0, in terms of Pm0 and E0.
subs(eqn3,omega0,0)
ans =
(333*Pm0)/100 - (333*E0*y0)/100 == 0
Effectively, this allows us to solve for y0 as
y0 == Pm0/E0
Similarly, eqn 4 allows us to solve for x0 as
x0 = 2*E0 - 1
But since we already know that x0 and y0 must lie on the unit circle, we must have
subs(x0^2 + y0^2 == 1,[x0,y0],[2*E0-1,Pm0/E0])
ans =
(2*E0 - 1)^2 + Pm0^2/E0^2 == 1
Multiply by E0^2,
simplify(subs(x0^2 + y0^2 == 1,[x0,y0],[2*E0-1,Pm0/E0])*E0^2)
ans =
4*E0^4 + Pm0^2 == 4*E0^3 | E0 == 0
So as long as E0 is non-zer0, that would seem to have 4 solutions, as effectively a 4th degree polynomial in E0. Some of those solutions may be complex of course. And that would depend on the value of Pm0. For some values of Pm0, no solution at all will exist. For example. when Pm0 is 1, we would have the polynomial:
fplot(subs(4*E0^4 + Pm0^2 - 4*E0^3,Pm0,1))
You can see this function NEVER crosses zero. Therefore E0 has no real solutions, and when Pm0 is 1, the problem never has a solution. But when Pm0 has a different value, might the polynomial have a solution?
4*E0^4 + Pm0^2 - 4*E0^3 == 0
As it turns out, since Pm0 is squared there, the constant term in that function is ALWAYS non-negative. The most extreme case is where Pm0 is 0, and then we will have a triple root at E0 == 0, and a single root at E0 == 1.
For some other cases of Pm0, when Pm0 a value between 0 and 1, then we will have exactly 2 real solutions for E0.
solve(subs(4*E0^4 + Pm0^2 - 4*E0^3,Pm0,.5),'maxdegree',4)
ans =
1/2
1/(9*((11^(1/2)*1728^(1/2))/1728 + 19/216)^(1/3)) + ((11^(1/2)*1728^(1/2))/1728 + 19/216)^(1/3) + 1/6
vpa(ans)
ans =
0.5
0.91964337760708056627592628232664
So your system of equations depends on the value of Pm0, but the number of solutions will vary.
