How to get a random possible solution from 'solve' function when getting unknown parameters z1 and conditions

1 visualización (últimos 30 días)
Here are my codes:
syms x y z;
eq1 = x*x + y*y + z*z == 0.14;
eq2 = x*y + y*z + z*x == 0.11;
equations = [eq1, eq2];
assume(x>=0);
assumeAlso(y>=0);
assumeAlso(z>=0);
assumeAlso(x<=1);
assumeAlso(y<=1);
assumeAlso(z<=1);
S=solve(equations,[x,y,z],'ReturnConditions',true);
disp(S)
S will print:
x: [4×1 sym]
y: [4×1 sym]
z: [4×1 sym]
parameters: z1
conditions: [4×1 sym]
and S.x is like:
(- 75*z1^2 - 30*z1 - 2)^(1/2)/10 - z1/2 - 3/10
- z1/2 - (- 75*z1^2 - 30*z1 - 2)^(1/2)/10 - 3/10
(- 75*z1^2 + 30*z1 - 2)^(1/2)/10 - z1/2 + 3/10
3/10 - (- 75*z1^2 + 30*z1 - 2)^(1/2)/10 - z1/2
and S.conditions is like:
5*z1 + 3 <= (- 75*z1^2 - 30*z1 - 2)^(1/2) & -13 <= 5*z1 + (- 75*z1^2 - 30*z1 - 2)^(1/2) & (- 75*z1^2 - 30*z1 - 2)^(1/2) <= 5*z1 + 13 & z1/2 + 3/10 <= -(- 75*z1^2 - 30*z1 - 2)^(1/2)/10 & 0 <= z1 & z1 <= 1
5*z1 + 3 <= (- 75*z1^2 - 30*z1 - 2)^(1/2) & -13 <= 5*z1 + (- 75*z1^2 - 30*z1 - 2)^(1/2) & (- 75*z1^2 - 30*z1 - 2)^(1/2) <= 5*z1 + 13 & z1/2 + 3/10 <= -(- 75*z1^2 - 30*z1 - 2)^(1/2)/10 & 0 <= z1 & z1 <= 1
5*z1 <= (- 75*z1^2 + 30*z1 - 2)^(1/2) + 3 & -7 <= 5*z1 + (- 75*z1^2 + 30*z1 - 2)^(1/2) & (- 75*z1^2 + 30*z1 - 2)^(1/2) <= 5*z1 + 7 & z1/2 <= 3/10 - (- 75*z1^2 + 30*z1 - 2)^(1/2)/10 & 0 <= z1 & z1 <= 1
5*z1 <= (- 75*z1^2 + 30*z1 - 2)^(1/2) + 3 & -7 <= 5*z1 + (- 75*z1^2 + 30*z1 - 2)^(1/2) & (- 75*z1^2 + 30*z1 - 2)^(1/2) <= 5*z1 + 7 & z1/2 <= 3/10 - (- 75*z1^2 + 30*z1 - 2)^(1/2)/10 & 0 <= z1 & z1 <= 1
How can i get any possible solution from parameter z1 and its conditions?
For instance, one of the possible solution would be x=0.1 y=0.2 z=0.3
But I have no idea how to get a random solution.
I'm sorry that my English is poor.
Thanks for your helping!!

Respuesta aceptada

John D'Errico
John D'Errico el 24 de Feb. de 2024
Editada: John D'Errico el 24 de Feb. de 2024
Solve cannot return a random possible solution. Wanting code to do what it is not designed to do will never be sufficient.
If you start from a specific start point though, vpasolve will generate a solution, at least some of the time.
So you can provide a random start point to vpasolve. This will work, at least if it is possible to work.
The problem is, you have only 2 equations, and 3 variables. So I'm actually surprised that solve was able to give you any solution at all.
syms x y z;
eq1 = x*x + y*y + z*z == 0.14;
eq2 = x*y + y*z + z*x == 0.11;
equations = [eq1, eq2];
assume(x>=0);
assumeAlso(y>=0);
assumeAlso(z>=0);
assumeAlso(x<=1);
assumeAlso(y<=1);
assumeAlso(z<=1);
S=solve(equations,[x,y,z],'ReturnConditions',true)
S = struct with fields:
x: [4×1 sym] y: [4×1 sym] z: [4×1 sym] parameters: z1 conditions: [4×1 sym]
It does mean that vpasolve will fail completely, as this is an under-determined problem. vpasolve generally does not handle that class of problem with more unknowns than equations.
vpasolve(equations,[x,y,z])
ans = struct with fields:
x: [0×1 sym] y: [0×1 sym] z: [0×1 sym]
If you want a random solution to be generated, you can choose a specific randome value for one of the variables. For example...
xrand = rand
xrand = 0.9908
That is, xrand will be some random value between 0 and 1.
vpasolve(subs(equations,x,xrand),[y,z])
ans = struct with fields:
y: [0×1 sym] z: [0×1 sym]
Of course this will fail most of the time. Why? Look at equation 1.
eq1 = x*x + y*y + z*z == 0.14;
So if x is greater than sqrt(0.14), NO real solution can ever exist. In fact, each of the variables x,y,z would follow the same constraint. Instead, try this:
xrand = rand*sqrt(0.14)
xrand = 0.2416
That is, xrand will be some random value between 0 and 1.
yzsol = vpasolve(subs(equations,x,xrand),[y,z])
yzsol = struct with fields:
y: [2×1 sym] z: [2×1 sym]
yzsol.y
ans = 
yzsol.z
ans = 
And we now see a pair of solutions, generated conditionally subject to the value of x.
Can you do better than that? Perhaps, if we recognize that equation 1 is the equation of the surface of a sphere with three independent variables. It is centered at the origin, at the point (0,0,0).
How about equation 2? That is not so obvious, but it will generally be a hyperbolic surface. Plotting these equations in the first octant, we see:
fimplicit3(eq1,[0 1 0 1 0 1])
hold on
fimplicit3(eq2,[0 1 0 1 0 1])
hold off
grid on
box on
view(76,13)
The solution to the system of two equations is where those two surfaces intersect. Do you see the solution locus will be roughly a circle living in the 3-d space of your variables (x,y,z)? I say roughly a circle, because I cannot trivially prove it is circular. It appears to be circular, at least by my eye.
Knowing that fact, can we derive a general solution? Sigh. Possibly. I'm not at all sure it is worth the time I would need to invest though, on what appears to be possibly just an example problem.
  2 comentarios
Chen Yu Yuan
Chen Yu Yuan el 24 de Feb. de 2024
It really helps, thanks a lot!
This has been bothering me for a long time. Thank you for helping me resolve it.
Hope you have a wonderful day!
John D'Errico
John D'Errico el 24 de Feb. de 2024
Editada: John D'Errico el 24 de Feb. de 2024
I'm happy it helps. I might conjecture from the picuture that the solution locus lies in a plane. And that it is possibly circular. Given that as an idea, it might be possible to prove that assertion, knowing what to look for.
First, a quick numerical foray would be in order. I find it is often a good idea to verify my intuitions about a problem.
syms x y z;
eq1 = x*x + y*y + z*z == 0.14;
eq2 = x*y + y*z + z*x == 0.11;
equations = [eq1, eq2];
xyz = zeros(0,3);
while size(xyz,1) < 100
xrand = rand*sqrt(0.14);
yzsol = vpasolve(subs(equations,x,xrand),[y,z]);
yz = double([yzsol.y,yzsol.z]);
yz(imag(yz(:,1)) ~= 0,:) = [];
if ~isempty(yz)
nsol = size(yz,1);
xyz = [xyz;[repmat(xrand,nsol,1),yz]];
end
end
I've generated 100 solutions there. First, do they lie in a plane? This is easy enough to test.
format long g
S = svd(xyz - mean(xyz,1))
S = 3×1
1.13871106047224 0.836680188626194 2.41834573934206e-16
So they do indeed lie exactly in a plane, to within floating point trash. What is the equation of the plane that contains those solutions?
[U,S,V] = svd(xyz - mean(xyz,1));
V(:,3)
ans = 3×1
-0.577350269189626 -0.577350269189626 -0.577350269189626
Interesting. They lie in the plane with normal vector [1 1 1]. We might also notice that all solutions had the propery that x+y+z==0.6.
Do you see that is the solutions lie in a plane, then the intersectino of a plane with the sphere of equation 1 predicts that all solutions MUST lie on a perfect circle? Anyway, this is about where one should go back, and with that idea in hand, prove my claim with some mathematical rigor. I have developed a small case of rigor mortis, so I'll stop here.

Iniciar sesión para comentar.

Más respuestas (0)

Etiquetas

Community Treasure Hunt

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

Start Hunting!

Translated by