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

1 visualización (últimos 30 días)
Chen Yu Yuan el 24 de Feb. de 2024
Editada: John D'Errico el 24 de Feb. de 2024
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!!
0 comentariosMostrar -2 comentarios más antiguosOcultar -2 comentarios más antiguos

Iniciar sesión para comentar.

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 comentariosMostrar NingunoOcultar Ninguno
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 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.

Categorías

Más información sobre Calculus en Help Center y File Exchange.

Community Treasure Hunt

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

Start Hunting!

Translated by