vpasolve finds no solutions except the null solution

9 次查看(过去 30 天)
I have the following code to obtain the solutions of a system of 10 nonlinear equations (eqn1-eqn10) with 10 variables (x1-x10).
x = (0,0,...,0) and x = (1,1,1,1,1,0,0,0,0,0) are solutions of this system. I believe that this system has other solutions, that I want to obtain, but the code is not working and I don't understand why. In fact, with range1 the code gives the null solution as output, But with range2, it gives empty. Shouldn't the code give the solution x = (1,1,1,1,1,0,0,0,0,0)?
Other question I have is if I need to add the equations x1 + x2 + x3 + x4 + x5 = 0 and x6 + x7 + x8 + x9 + x10 = 0, the code gives an error. Can I not have more equations than variables? How do I deal with this situation?
syms x1 x2 x3 x4 x5 x6 x7 x8 x9 x10
A1 = (x1 - x2)^2 + (x6 - x7)^2 - (x2 - x3)^2 - (x7 - x8)^2;
A2 = (x1 - x2)^2 + (x6 - x7)^2 - (x3 - x4)^2 - (x8 - x9)^2;
A3 = (x1 - x2)^2 + (x6 - x7)^2 - (x4 - x5)^2 - (x9 - x10)^2;
A4 = (x1 - x2)^2 + (x6 - x7)^2 - (x1 - x5)^2 - (x6 - x10)^2;
D1 = (x1 - x3)^2 + (x6 - x8)^2 - (x1 - x4)^2 - (x6 - x9)^2;
D2 = (x1 - x3)^2 + (x6 - x8)^2 - (x2 - x4)^2 - (x7 - x9)^2;
D3 = (x1 - x3)^2 + (x6 - x8)^2 - (x3 - x5)^2 - (x8 - x10)^2;
D4 = (x1 - x3)^2 + (x6 - x8)^2 - (x5 - x2)^2 - (x10 - x7)^2;
eqn1 = (D1 + D2 + D3 + D4) * (x1 - x3) + (4 * D1 - D2 - D3 - ...
D4) * (x4 - x1) + (A1 + A2 + A3 + A4) * (x1 - x2) + (4 * A4 - A1 - ...
A2 - A3) * (x5 - x1) == 0;
eqn2 = (4 * D2 - D1 - D3 - D4) * (x4 - x2) + (4 * D4 - D1 - D2 - ...
D3) * (x5 - x2) + (4 * A1 - A2 - A3 - A4) * (x3 - x2) == 0;
eqn3 = (D1 + D2 + D3 + D4) * (x3 - x1) + (4 * D3 - D1 - D2 - ...
D4) * (x5 - x3) + (4 * A1 - A2 - A3 - A4) * (x2 - x3) + (4 * A2 - A1 - ...
A3 - A4) * (x4 - x3) == 0;
eqn4 = (4 * D1 - D2 - D3 - D4) * (x1 - x4) + (4 * D2 - D1 - D3 - ...
D4) * (x2 - x4) + (4 * A3 - A1 - A2 - A4) * (x5 - x4) == 0;
eqn5 = (4 * D4 - D1 - D2 - D3) * (x2 - x5) + (4 * D3 - D1 - D2 - ...
D4) * (x3 - x5) + (4 * A3 - A1 - A2 - A4) * (x4 - x5) + (4 * A4 - A1 - ...
A2 - A3) * (x1 - x5) == 0;
eqn6 = (D1 + D2 + D3 + D4) * (x6 - x8) + (4 * D1 - D2 - D3 - D4) * (x9 - x6) + ...
(A1 + A2 + A3 + A4) * (x6 - x7) + (4 * A4 - A1 -A2 - A3) * (x10 - x6) == 0;
eqn7 = (4 * D2 - D1 - D3 - D4) * (x9 - x7) + (4 * D4 - D1 - D2 - ...
D3) * (x10 - x7) + (4 * A1 - A2 - A3 - A4) * (x8 - x7) == 0;
eqn8 = (D1 + D2 + D3 + D4) * (x8 - x6) + (4 * D3 - D1 - D2 - ...
D4) * (x10 - x8) + (4 * A1 - A2 - A3 - A4) * (x7 - x8) + (4 * A2 - A1 - ...
A3 - A4) * (x9 - x8) == 0;
eqn9 = (4 * D1 - D2 - D3 - D4) * (x6 - x9) + (4 * D2 - D1 - D3 - ...
D4) * (x7 - x9) + (4 * A3 - A1 - A2 - A4) * (x10 - x9) == 0;
eqn10 = (4 * D4 - D1 - D2 - D3) * (x7 - x10) + (4 * D3 - D1 - D2 - ...
D4) * (x8 - x10) + (4 * A3 - A1 - A2 - A4) * (x9 - x10) + (4 * A4 - A1 - ...
A2 - A3) * (x6 - x10) == 0;
equations = [eqn1 eqn2 eqn3 eqn4 eqn5 eqn6 eqn7 eqn8 eqn9 eqn10];
vars = [x1 x2 x3 x4 x5 x6 x7 x8 x9 x10];
range1 = [-Inf Inf; -Inf Inf; -Inf Inf; -Inf Inf; -Inf Inf; -Inf Inf; -Inf Inf; -Inf Inf; -Inf Inf; -Inf Inf];
range2 = [0.5 1.5; 0.5 1.5; 0.5 1.5; 0.5 1.5; 0.5 1.5; -0.5 0.5; -0.5 0.5; -0.5 0.5; -0.5 0.5; -0.5 0.5];
Y = vpasolve(equations, vars, range1);
[Y.x1 Y.x2 Y.x3 Y.x4 Y.x5 Y.x6 Y.x7 Y.x8 Y.x9 Y.x10]
Y = vpasolve(equations, vars, range2);
[Y.x1 Y.x2 Y.x3 Y.x4 Y.x5 Y.x6 Y.x7 Y.x8 Y.x9 Y.x10]
  1 个评论
Matt J
Matt J 2025-3-24
I believe that this system has other solutions, that I want to obtain,
Well, you can't get them all numerically. There is an infinite continuum of solutions, i.e. too many to write down. Any x of the form,
x=(a,a,a,a,a,b,b,b,b,b)
is a solution (and there may be others).

请先登录,再进行评论。

采纳的回答

Matt J
Matt J 2025-3-24
编辑:Matt J 2025-3-24
With numerical algorithms like vpasolve, you never can tell what it will find, if anything. In any case, lsqnonlin seems to offer better performance, assuming you have the Optimizatio Toolbox. It finds a solution for range2 readily with and without the additional equations:
syms x [10,1]
A1 = (x1 - x2)^2 + (x6 - x7)^2 - (x2 - x3)^2 - (x7 - x8)^2;
A2 = (x1 - x2)^2 + (x6 - x7)^2 - (x3 - x4)^2 - (x8 - x9)^2;
A3 = (x1 - x2)^2 + (x6 - x7)^2 - (x4 - x5)^2 - (x9 - x10)^2;
A4 = (x1 - x2)^2 + (x6 - x7)^2 - (x1 - x5)^2 - (x6 - x10)^2;
D1 = (x1 - x3)^2 + (x6 - x8)^2 - (x1 - x4)^2 - (x6 - x9)^2;
D2 = (x1 - x3)^2 + (x6 - x8)^2 - (x2 - x4)^2 - (x7 - x9)^2;
D3 = (x1 - x3)^2 + (x6 - x8)^2 - (x3 - x5)^2 - (x8 - x10)^2;
D4 = (x1 - x3)^2 + (x6 - x8)^2 - (x5 - x2)^2 - (x10 - x7)^2;
eqn1 = (D1 + D2 + D3 + D4) * (x1 - x3) + (4 * D1 - D2 - D3 - ...
D4) * (x4 - x1) + (A1 + A2 + A3 + A4) * (x1 - x2) + (4 * A4 - A1 - ...
A2 - A3) * (x5 - x1) == 0;
eqn2 = (4 * D2 - D1 - D3 - D4) * (x4 - x2) + (4 * D4 - D1 - D2 - ...
D3) * (x5 - x2) + (4 * A1 - A2 - A3 - A4) * (x3 - x2) == 0;
eqn3 = (D1 + D2 + D3 + D4) * (x3 - x1) + (4 * D3 - D1 - D2 - ...
D4) * (x5 - x3) + (4 * A1 - A2 - A3 - A4) * (x2 - x3) + (4 * A2 - A1 - ...
A3 - A4) * (x4 - x3) == 0;
eqn4 = (4 * D1 - D2 - D3 - D4) * (x1 - x4) + (4 * D2 - D1 - D3 - ...
D4) * (x2 - x4) + (4 * A3 - A1 - A2 - A4) * (x5 - x4) == 0;
eqn5 = (4 * D4 - D1 - D2 - D3) * (x2 - x5) + (4 * D3 - D1 - D2 - ...
D4) * (x3 - x5) + (4 * A3 - A1 - A2 - A4) * (x4 - x5) + (4 * A4 - A1 - ...
A2 - A3) * (x1 - x5) == 0;
eqn6 = (D1 + D2 + D3 + D4) * (x6 - x8) + (4 * D1 - D2 - D3 - D4) * (x9 - x6) + ...
(A1 + A2 + A3 + A4) * (x6 - x7) + (4 * A4 - A1 -A2 - A3) * (x10 - x6) == 0;
eqn7 = (4 * D2 - D1 - D3 - D4) * (x9 - x7) + (4 * D4 - D1 - D2 - ...
D3) * (x10 - x7) + (4 * A1 - A2 - A3 - A4) * (x8 - x7) == 0;
eqn8 = (D1 + D2 + D3 + D4) * (x8 - x6) + (4 * D3 - D1 - D2 - ...
D4) * (x10 - x8) + (4 * A1 - A2 - A3 - A4) * (x7 - x8) + (4 * A2 - A1 - ...
A3 - A4) * (x9 - x8) == 0;
eqn9 = (4 * D1 - D2 - D3 - D4) * (x6 - x9) + (4 * D2 - D1 - D3 - ...
D4) * (x7 - x9) + (4 * A3 - A1 - A2 - A4) * (x10 - x9) == 0;
eqn10 = (4 * D4 - D1 - D2 - D3) * (x7 - x10) + (4 * D3 - D1 - D2 - ...
D4) * (x8 - x10) + (4 * A3 - A1 - A2 - A4) * (x9 - x10) + (4 * A4 - A1 - ...
A2 - A3) * (x6 - x10) == 0;
equations = [eqn1 eqn2 eqn3 eqn4 eqn5 eqn6 eqn7 eqn8 eqn9 eqn10];
range2 = [0.5 1.5; 0.5 1.5; 0.5 1.5; 0.5 1.5; 0.5 1.5; -0.5 0.5; -0.5 0.5; -0.5 0.5; -0.5 0.5; -0.5 0.5];
equations= lhs(equations)-rhs(equations);
F=matlabFunction(equations, File="eqnFunc.m", Vars={x});
xInitial=mean(range2,2) + randn(10,1);
opts=optimoptions('lsqnonlin', 'ConstraintTol', 1e-16,'FunctionTol',1e-16,...
'StepTol',1e-16,'OptimalityTol',1e-16);
[xnum,fval]=lsqnonlin(F,xInitial, range2(:,1), range2(:,2),opts)
Local minimum possible. lsqnonlin stopped because the final change in the sum of squares relative to its initial value is less than the value of the function tolerance.
xnum = 10×1
0.9196 0.9198 0.9197 0.9191 0.9192 0.0088 0.0088 0.0088 0.0088 0.0088
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
fval = 2.2809e-18
Can I not have more equations than variables? How do I deal with this situation?
You can, but of course you cannot be guaranteed of an exact solution,
equations=[equations, (x1 + x2 + x3 + x4 + x5) , (x6 + x7 + x8 + x9 + x10)];
F=matlabFunction(equations, File="eqnFunc.m", Vars={x});
[xnum,fval]=lsqnonlin(F, xInitial, range2(:,1), range2(:,2), opts)
Local minimum possible. lsqnonlin stopped because the size of the current step is less than the value of the step size tolerance.
xnum = 10×1
0.5000 0.5000 0.5000 0.5000 0.5000 -0.0000 0.0000 -0.0000 0.0000 -0.0000
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
fval = 6.2500
  9 个评论
Catarina Pina
Catarina Pina 2025-3-25
I realize there are an infinite number of solutions, I just wanted to find some of them. With this code I can do what I need, thanks!
Matt J
Matt J 2025-3-25
You are quite welcome, but if this answer addresses your question, please Accept-click it.

请先登录,再进行评论。

更多回答(2 个)

Paul
Paul 2025-3-25
This code exceeds the run time limit here, and is taking forever running on my computer, so don't know what it produces.
If solve produced a solution when specifying two of the variables, then try using it with ReturnConditions = true to get a parameterization of all solutions.
Might also want to assume all of the xi are real if that is, in fact, the case. Might narrow down the solution space.
syms x1 x2 x3 x4 x5 x6 x7 x8 x9 x10
A1 = (x1 - x2)^2 + (x6 - x7)^2 - (x2 - x3)^2 - (x7 - x8)^2;
A2 = (x1 - x2)^2 + (x6 - x7)^2 - (x3 - x4)^2 - (x8 - x9)^2;
A3 = (x1 - x2)^2 + (x6 - x7)^2 - (x4 - x5)^2 - (x9 - x10)^2;
A4 = (x1 - x2)^2 + (x6 - x7)^2 - (x1 - x5)^2 - (x6 - x10)^2;
D1 = (x1 - x3)^2 + (x6 - x8)^2 - (x1 - x4)^2 - (x6 - x9)^2;
D2 = (x1 - x3)^2 + (x6 - x8)^2 - (x2 - x4)^2 - (x7 - x9)^2;
D3 = (x1 - x3)^2 + (x6 - x8)^2 - (x3 - x5)^2 - (x8 - x10)^2;
D4 = (x1 - x3)^2 + (x6 - x8)^2 - (x5 - x2)^2 - (x10 - x7)^2;
eqn1 = (D1 + D2 + D3 + D4) * (x1 - x3) + (4 * D1 - D2 - D3 - ...
D4) * (x4 - x1) + (A1 + A2 + A3 + A4) * (x1 - x2) + (4 * A4 - A1 - ...
A2 - A3) * (x5 - x1) == 0;
eqn2 = (4 * D2 - D1 - D3 - D4) * (x4 - x2) + (4 * D4 - D1 - D2 - ...
D3) * (x5 - x2) + (4 * A1 - A2 - A3 - A4) * (x3 - x2) == 0;
eqn3 = (D1 + D2 + D3 + D4) * (x3 - x1) + (4 * D3 - D1 - D2 - ...
D4) * (x5 - x3) + (4 * A1 - A2 - A3 - A4) * (x2 - x3) + (4 * A2 - A1 - ...
A3 - A4) * (x4 - x3) == 0;
eqn4 = (4 * D1 - D2 - D3 - D4) * (x1 - x4) + (4 * D2 - D1 - D3 - ...
D4) * (x2 - x4) + (4 * A3 - A1 - A2 - A4) * (x5 - x4) == 0;
eqn5 = (4 * D4 - D1 - D2 - D3) * (x2 - x5) + (4 * D3 - D1 - D2 - ...
D4) * (x3 - x5) + (4 * A3 - A1 - A2 - A4) * (x4 - x5) + (4 * A4 - A1 - ...
A2 - A3) * (x1 - x5) == 0;
eqn6 = (D1 + D2 + D3 + D4) * (x6 - x8) + (4 * D1 - D2 - D3 - D4) * (x9 - x6) + ...
(A1 + A2 + A3 + A4) * (x6 - x7) + (4 * A4 - A1 -A2 - A3) * (x10 - x6) == 0;
eqn7 = (4 * D2 - D1 - D3 - D4) * (x9 - x7) + (4 * D4 - D1 - D2 - ...
D3) * (x10 - x7) + (4 * A1 - A2 - A3 - A4) * (x8 - x7) == 0;
eqn8 = (D1 + D2 + D3 + D4) * (x8 - x6) + (4 * D3 - D1 - D2 - ...
D4) * (x10 - x8) + (4 * A1 - A2 - A3 - A4) * (x7 - x8) + (4 * A2 - A1 - ...
A3 - A4) * (x9 - x8) == 0;
eqn9 = (4 * D1 - D2 - D3 - D4) * (x6 - x9) + (4 * D2 - D1 - D3 - ...
D4) * (x7 - x9) + (4 * A3 - A1 - A2 - A4) * (x10 - x9) == 0;
eqn10 = (4 * D4 - D1 - D2 - D3) * (x7 - x10) + (4 * D3 - D1 - D2 - ...
D4) * (x8 - x10) + (4 * A3 - A1 - A2 - A4) * (x9 - x10) + (4 * A4 - A1 - ...
A2 - A3) * (x6 - x10) == 0;
equations = [eqn1 eqn2 eqn3 eqn4 eqn5 eqn6 eqn7 eqn8 eqn9 eqn10];
vars = [x1 x2 x3 x4 x5 x6 x7 x8 x9 x10];
Y10 = solve(equations,vars,'ReturnConditions',true)
eqn11 = x1 + x2 + x3 + x4 + x5 == 0;
eqn12 = x6 + x7 + x8 + x9 + x10 == 0;
Y12 = solve([equations,eqn11,eqn12],vars,'ReturnConditions',true)

Alex Sha
Alex Sha 2025-6-21
For:
range2 = [0.5 1.5; 0.5 1.5; 0.5 1.5; 0.5 1.5; 0.5 1.5; -0.5 0.5; -0.5 0.5; -0.5 0.5; -0.5 0.5; -0.5 0.5];
result-1:
x1: 1.08082083137654
x2: 1.08080951998694
x3: 1.08081167383851
x4: 1.0808049613845
x5: 1.0807961790416
x6: 0.139381627943366
x7: 0.139378007206469
x8: 0.139387284204371
x9: 0.13939708413668
x10: 0.139379799060729
feval:
4.18856104328548E-14
-3.48053189694011E-15
1.9850994632769E-15
-1.31438771899855E-14
-3.90425935929704E-14
-9.46492326396661E-15
7.3761149429896E-16
-6.93623770817858E-15
3.28876546733261E-14
-1.31813703670225E-14
result-2:
x1: 1.02958234918508
x2: 1.02959012135357
x3: 1.0295940618282
x4: 1.02959069634112
x5: 1.02958674833031
x6: -0.0307085239327012
x7: -0.0306995429173505
x8: -0.0306983806107563
x9: -0.0307242623557914
x10: -0.0307004373970983
feval:
8.54065517352478E-15
5.87827032304209E-16
-6.61360139369969E-15
9.51411238197795E-15
5.9308056815423E-16
2.18297090864782E-14
4.55445556176149E-14
4.51489361938118E-14
-8.39926540201438E-14
2.78159953330793E-14
result-3:
x1: 1.11722991939913
x2: 1.11722093473141
x3: 1.11721577030959
x4: 1.11723361423116
x5: 1.11722092000214
x6: 0.350956747942537
x7: 0.350948443506383
x8: 0.350952713845785
x9: 0.350953745554248
x10: 0.350945892821742
feval:
8.43258754382745E-15
-9.218608646947E-15
-1.5533394955795E-14
6.34621251709236E-15
-3.51572661481046E-15
2.14134991818329E-16
-4.90905368161458E-16
-7.00918017460192E-15
4.90030773194566E-15
1.47945732346218E-16
...
There are infinite results.

类别

Help CenterFile Exchange 中查找有关 Programming 的更多信息

产品


版本

R2022b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by