Solve 3 unknown variable system of equations

2 views (last 30 days)
Jack Dong
Jack Dong on 15 May 2021
Commented: Star Strider on 15 May 2021
syms X Z beta
deg2rad=pi/180;
theta_b=15*deg2rad;
theta_p=15*deg2rad;
R_b=1;
R_p=0.9;
B1=[R_b*cos(15*deg2rad);R_b*sin(15*deg2rad);0];
B2=[R_b*cos(105*deg2rad);R_b*sin(105*deg2rad);0];
B3=[R_b*cos(135*deg2rad);R_b*sin(135*deg2rad);0];
B4=[R_b*cos(225*deg2rad);R_b*sin(225*deg2rad);0];
B5=[R_b*cos(255*deg2rad);R_b*sin(255*deg2rad);0];
B6=[R_b*cos(345*deg2rad);R_b*sin(345*deg2rad);0];
Bx1=B1(1,:);
Bx2=B2(1,:);
Bx3=B3(1,:);
Bx4=B4(1,:);
Bx5=B5(1,:);
Bx6=B6(1,:);
By1=B1(2,:);
By2=B2(2,:);
By3=B3(2,:);
By4=B4(2,:);
By5=B5(2,:);
By6=B6(2,:);
Bz1=B1(3,:);
Bz2=B2(3,:);
Bz3=B3(3,:);
Bz4=B4(3,:);
Bz5=B5(3,:);
Bz6=B6(3,:);
P1=[R_p*cos(75*deg2rad);R_b*sin(75*deg2rad);0];
P2=[R_p*cos(165*deg2rad);R_b*sin(165*deg2rad);0];
P3=[R_p*cos(195*deg2rad);R_b*sin(195*deg2rad);0];
P4=[R_p*cos(285*deg2rad);R_b*sin(285*deg2rad);0];
P5=[R_p*cos(315*deg2rad);R_b*sin(315*deg2rad);0];
P6=[R_p*cos(405*deg2rad);R_b*sin(405*deg2rad);0];
Px1=P1(1,:);
Px2=P2(1,:);
Px3=P3(1,:);
Px4=P4(1,:);
Px5=P5(1,:);
Px6=P6(1,:);
Py1=P1(2,:);
Py2=P2(2,:);
Py3=P3(2,:);
Py4=P4(2,:);
Py5=P5(2,:);
Py6=P6(2,:);
Pz1=P1(3,:);
Pz2=P2(3,:);
Pz3=P3(3,:);
Pz4=P4(3,:);
Pz5=P5(3,:);
Pz6=P6(3,:);
Y=0;
alpha=0;
gama=0;
L1=(X+((cos(gama)*cos(beta))*Px1)+(((cos(gama)*sin(beta)*sin(alpha))-(sin(gama)*cos(alpha)))*Py1)+(((cos(gama)*sin(beta)*cos(alpha))+(sin(gama)*sin(alpha)))*Pz1)-Bx1)^2+(Y+((sin(gama)*cos(beta))*Px1)+(((sin(gama)*sin(beta)*sin(alpha))+(cos(gama)*cos(alpha)))*Py1)+(((sin(gama)*sin(beta)*cos(alpha))-cos(gama)*sin(alpha)))*Pz1)-By1)^2+(Z+((-sin(beta))*Px1)+((cos(beta)*sin(alpha))*Py1)+((cos(beta)*cos(alpha))*Pz1)-Bz1)^2;
L2=(X+((cos(gama)*cos(beta))*Px2)+(((cos(gama)*sin(beta)*sin(alpha))-(sin(gama)*cos(alpha)))*Py2)+(((cos(gama)*sin(beta)*cos(alpha))+(sin(gama)*sin(alpha)))*Pz2)-Bx2)^2+(Y+((sin(gama)*cos(beta))*Px2)+(((sin(gama)*sin(beta)*sin(alpha))+(cos(gama)*cos(alpha)))*Py2)+(((sin(gama)*sin(beta)*cos(alpha))-cos(gama)*sin(alpha)))*Pz2)-By2)^2+(Z+((-sin(beta))*Px2)+((cos(beta)*sin(alpha))*Py2)+((cos(beta)*cos(alpha))*Pz2)-Bz2)^2;
L3=(X+((cos(gama)*cos(beta))*Px3)+(((cos(gama)*sin(beta)*sin(alpha))-(sin(gama)*cos(alpha)))*Py3)+(((cos(gama)*sin(beta)*cos(alpha))+(sin(gama)*sin(alpha)))*Pz3)-Bx3)^2+(Y+((sin(gama)*cos(beta))*Px3)+(((sin(gama)*sin(beta)*sin(alpha))+(cos(gama)*cos(alpha)))*Py3)+(((sin(gama)*sin(beta)*cos(alpha))-cos(gama)*sin(alpha)))*Pz3)-By3)^2+(Z+((-sin(beta))*Px3)+((cos(beta)*sin(alpha))*Py3)+((cos(beta)*cos(alpha))*Pz3)-Bz3)^2;
A=solve(L1==1,L2==1.1,L3==1.2);

Sign in to answer this question.

Accepted Answer

Star Strider
Star Strider on 15 May 2021
Ran in:
There are likely an infinity of solutions.
Also, for whatever reason, the solve function does not like ‘L6’.
A numeric approach (with different initial parameter estimates) —
syms X Z beta
deg2rad=pi/180;
theta_b=15*deg2rad;
theta_p=15*deg2rad;
R_b=1;
R_p=0.9;
B1=[R_b*cos(15*deg2rad);R_b*sin(15*deg2rad);0];
B2=[R_b*cos(105*deg2rad);R_b*sin(105*deg2rad);0];
B3=[R_b*cos(135*deg2rad);R_b*sin(135*deg2rad);0];
B4=[R_b*cos(225*deg2rad);R_b*sin(225*deg2rad);0];
B5=[R_b*cos(255*deg2rad);R_b*sin(255*deg2rad);0];
B6=[R_b*cos(345*deg2rad);R_b*sin(345*deg2rad);0];
Bx1=B1(1,:);
Bx2=B2(1,:);
Bx3=B3(1,:);
Bx4=B4(1,:);
Bx5=B5(1,:);
Bx6=B6(1,:);
By1=B1(2,:);
By2=B2(2,:);
By3=B3(2,:);
By4=B4(2,:);
By5=B5(2,:);
By6=B6(2,:);
Bz1=B1(3,:);
Bz2=B2(3,:);
Bz3=B3(3,:);
Bz4=B4(3,:);
Bz5=B5(3,:);
Bz6=B6(3,:);
P1=[R_p*cos(75*deg2rad);R_b*sin(75*deg2rad);0];
P2=[R_p*cos(165*deg2rad);R_b*sin(165*deg2rad);0];
P3=[R_p*cos(195*deg2rad);R_b*sin(195*deg2rad);0];
P4=[R_p*cos(285*deg2rad);R_b*sin(285*deg2rad);0];
P5=[R_p*cos(315*deg2rad);R_b*sin(315*deg2rad);0];
P6=[R_p*cos(405*deg2rad);R_b*sin(405*deg2rad);0];
Px1=P1(1,:);
Px2=P2(1,:);
Px3=P3(1,:);
Px4=P4(1,:);
Px5=P5(1,:);
Px6=P6(1,:);
Py1=P1(2,:);
Py2=P2(2,:);
Py3=P3(2,:);
Py4=P4(2,:);
Py5=P5(2,:);
Py6=P6(2,:);
Pz1=P1(3,:);
Pz2=P2(3,:);
Pz3=P3(3,:);
Pz4=P4(3,:);
Pz5=P5(3,:);
Pz6=P6(3,:);
Y=0;
alpha=0;
gama=0;
L1=((X+((cos(gama)*cos(beta))*Px1)+(((cos(gama)*sin(beta)*sin(alpha))-(sin(gama)*cos(alpha)))*Py1)+(((cos(gama)*sin(beta)*cos(alpha))+(sin(gama)*sin(alpha)))*Pz1)-Bx1)^2+(Y+((sin(gama)*cos(beta))*Px1)+(((sin(gama)*sin(beta)*sin(alpha))+(cos(gama)*cos(alpha)))*Py1)+(((sin(gama)*sin(beta)*cos(alpha))-cos(gama)*sin(alpha)))*Pz1)-By1)^2+(Z+((-sin(beta))*Px1)+((cos(beta)*sin(alpha))*Py1)+((cos(beta)*cos(alpha))*Pz1)-Bz1)^2;
L2=((X+((cos(gama)*cos(beta))*Px2)+(((cos(gama)*sin(beta)*sin(alpha))-(sin(gama)*cos(alpha)))*Py2)+(((cos(gama)*sin(beta)*cos(alpha))+(sin(gama)*sin(alpha)))*Pz2)-Bx2)^2+(Y+((sin(gama)*cos(beta))*Px2)+(((sin(gama)*sin(beta)*sin(alpha))+(cos(gama)*cos(alpha)))*Py2)+(((sin(gama)*sin(beta)*cos(alpha))-cos(gama)*sin(alpha)))*Pz2)-By2)^2+(Z+((-sin(beta))*Px2)+((cos(beta)*sin(alpha))*Py2)+((cos(beta)*cos(alpha))*Pz2)-Bz2)^2;
L3=((X+((cos(gama)*cos(beta))*Px3)+(((cos(gama)*sin(beta)*sin(alpha))-(sin(gama)*cos(alpha)))*Py3)+(((cos(gama)*sin(beta)*cos(alpha))+(sin(gama)*sin(alpha)))*Pz3)-Bx3)^2+(Y+((sin(gama)*cos(beta))*Px3)+(((sin(gama)*sin(beta)*sin(alpha))+(cos(gama)*cos(alpha)))*Py3)+(((sin(gama)*sin(beta)*cos(alpha))-cos(gama)*sin(alpha)))*Pz3)-By3)^2+(Z+((-sin(beta))*Px3)+((cos(beta)*sin(alpha))*Py3)+((cos(beta)*cos(alpha))*Pz3)-Bz3)^2;
L4=((X+((cos(gama)*cos(beta))*Px4)+(((cos(gama)*sin(beta)*sin(alpha))-(sin(gama)*cos(alpha)))*Py4)+(((cos(gama)*sin(beta)*cos(alpha))+(sin(gama)*sin(alpha)))*Pz4)-Bx4)^2+(Y+((sin(gama)*cos(beta))*Px4)+(((sin(gama)*sin(beta)*sin(alpha))+(cos(gama)*cos(alpha)))*Py4)+(((sin(gama)*sin(beta)*cos(alpha))-cos(gama)*sin(alpha)))*Pz4)-By4)^2+(Z+((-sin(beta))*Px4)+((cos(beta)*sin(alpha))*Py4)+((cos(beta)*cos(alpha))*Pz4)-Bz4)^2;
L5=((X+((cos(gama)*cos(beta))*Px5)+(((cos(gama)*sin(beta)*sin(alpha))-(sin(gama)*cos(alpha)))*Py5)+(((cos(gama)*sin(beta)*cos(alpha))+(sin(gama)*sin(alpha)))*Pz5)-Bx5)^2+(Y+((sin(gama)*cos(beta))*Px5)+(((sin(gama)*sin(beta)*sin(alpha))+(cos(gama)*cos(alpha)))*Py5)+(((sin(gama)*sin(beta)*cos(alpha))-cos(gama)*sin(alpha)))*Pz5)-By5)^2+(Z+((-sin(beta))*Px5)+((cos(beta)*sin(alpha))*Py5)+((cos(beta)*cos(alpha))*Pz5)-Bz5)^2;
L6=((X+((cos(gama)*cos(beta))*Px6)+(((cos(gama)*sin(beta)*sin(alpha))-(sin(gama)*cos(alpha)))*Py6)+(((cos(gama)*sin(beta)*cos(alpha))+(sin(gama)*sin(alpha)))*Pz6)-Bx6)^2+(Y+((sin(gama)*cos(beta))*Px6)+(((sin(gama)*sin(beta)*sin(alpha))+(cos(gama)*cos(alpha)))*Py6)+(((sin(gama)*sin(beta)*cos(alpha))-cos(gama)*sin(alpha)))*Pz6)-By6)^2+(Z+((-sin(beta))*Px6)+((cos(beta)*sin(alpha))*Py6)+((cos(beta)*cos(alpha))*Pz6)-Bz6)^2;
% A=solve(L1==1,L2==1.1,L3==1.2,L4==1.3,L5==1.4,L6=1.5);
f = matlabFunction([L1-1; L2-1.1; L3-1.2; L4-1.3; L5-1.4; L6-1.5])
f = function_handle with value:
@(X,Z,beta)[(Z-sin(beta).*2.329371405922687e-1).^2+((X+cos(beta).*2.329371405922687e-1-9.659258262890683e-1).^2+7.071067811865476e-1).^2-1.0;(Z+sin(beta).*8.693332436601614e-1).^2+((X-cos(beta).*8.693332436601614e-1+2.588190451025209e-1).^2-7.071067811865473e-1).^2-1.1e+1./1.0e+1;(sqrt(2.0)./2.0-(X-cos(beta).*8.693332436601615e-1+sqrt(2.0)./2.0).^2+2.588190451025208e-1).^2+(Z+sin(beta).*8.693332436601615e-1).^2-6.0./5.0;((X+cos(beta).*2.329371405922683e-1+sqrt(2.0)./2.0).^2+sqrt(2.0)./2.0-9.659258262890684e-1).^2+(Z-sin(beta).*2.329371405922683e-1).^2-1.3e+1./1.0e+1;(sqrt(2.0).*(-1.0./2.0)+(X+sqrt(2.0).*cos(beta).*(9.0./2.0e+1)+2.588190451025206e-1).^2+9.659258262890683e-1).^2+(Z-sqrt(2.0).*sin(beta).*(9.0./2.0e+1)).^2-7.0./5.0;(sqrt(2.0)./2.0+(X+sqrt(2.0).*cos(beta).*(9.0./2.0e+1)-9.659258262890683e-1).^2+2.588190451025207e-1).^2+(Z-sqrt(2.0).*sin(beta).*(9.0./2.0e+1)).^2-3.0./2.0]
format long G
B1 = fsolve(@(b)f(b(1),b(2),b(3)), ones(3,1))
Warning: Trust-region-dogleg algorithm of FSOLVE cannot handle non-square systems; using Levenberg-Marquardt algorithm instead.
No solution found. fsolve stopped because the last step was ineffective. However, the vector of function values is not near zero, as measured by the value of the function tolerance.
B1 = 3×1
0.33876335638099 0.0804317184641243 1.05638916688373
B11 = fsolve(@(b)f(b(1),b(2),b(3)), ones(3,1)*100)
Warning: Trust-region-dogleg algorithm of FSOLVE cannot handle non-square systems; using Levenberg-Marquardt algorithm instead.
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.
B11 = 3×1
2.22171447863478 21.7903686240671 98.7112332112143
B2 = fminunc(@(b)norm(f(b(1),b(2),b(3))), ones(3,1))
Local minimum found. Optimization completed because the size of the gradient is less than the value of the optimality tolerance.
B2 = 3×1
0.338761475210389 0.0804288959290625 1.05638545964194
B21 = fminunc(@(b)norm(f(b(1),b(2),b(3))), ones(3,1)*100)
Local minimum found. Optimization completed because the size of the gradient is less than the value of the optimality tolerance.
B21 = 3×1
-0.443281068518347 -0.251219602486744 100.438173882528
format short
.
  2 Comments
Jack Dong
Jack Dong on 15 May 2021
Thanks a lot @Star Strider,
I would like to update a new script, I can run it but it returns issue that "Text exceeds maximum line length of 25,000 characters for Command Window display", is there any ways to decrease the decimals.
syms X Z beta
%#codegen
deg2rad=pi/180;
theta_b=15*deg2rad;
theta_p=15*deg2rad;
R_b=1;
R_p=0.9;
Q0=1; %initial leg length
B1=[R_b*cos(15*deg2rad);R_b*sin(15*deg2rad);0];
B2=[R_b*cos(105*deg2rad);R_b*sin(105*deg2rad);0];
B3=[R_b*cos(135*deg2rad);R_b*sin(135*deg2rad);0];
B4=[R_b*cos(225*deg2rad);R_b*sin(225*deg2rad);0];
B5=[R_b*cos(255*deg2rad);R_b*sin(255*deg2rad);0];
B6=[R_b*cos(345*deg2rad);R_b*sin(345*deg2rad);0];
Bx1=B1(1,:);
Bx2=B2(1,:);
Bx3=B3(1,:);
Bx4=B4(1,:);
Bx5=B5(1,:);
Bx6=B6(1,:);
By1=B1(2,:);
By2=B2(2,:);
By3=B3(2,:);
By4=B4(2,:);
By5=B5(2,:);
By6=B6(2,:);
Bz1=B1(3,:);
Bz2=B2(3,:);
Bz3=B3(3,:);
Bz4=B4(3,:);
Bz5=B5(3,:);
Bz6=B6(3,:);
P1=[R_p*cos(75*deg2rad);R_b*sin(75*deg2rad);0];
P2=[R_p*cos(165*deg2rad);R_b*sin(165*deg2rad);0];
P3=[R_p*cos(195*deg2rad);R_b*sin(195*deg2rad);0];
P4=[R_p*cos(285*deg2rad);R_b*sin(285*deg2rad);0];
P5=[R_p*cos(315*deg2rad);R_b*sin(315*deg2rad);0];
P6=[R_p*cos(405*deg2rad);R_b*sin(405*deg2rad);0];
Px1=P1(1,:);
Px2=P2(1,:);
Px3=P3(1,:);
Px4=P4(1,:);
Px5=P5(1,:);
Px6=P6(1,:);
Py1=P1(2,:);
Py2=P2(2,:);
Py3=P3(2,:);
Py4=P4(2,:);
Py5=P5(2,:);
Py6=P6(2,:);
Pz1=P1(3,:);
Pz2=P2(3,:);
Pz3=P3(3,:);
Pz4=P4(3,:);
Pz5=P5(3,:);
Pz6=P6(3,:);
Y=0;
alpha=0;
gama=0;
A=solve((X+((cos(gama)*cos(beta))*Px1)+(((cos(gama)*sin(beta)*sin(alpha))-(sin(gama)*cos(alpha)))*Py1)+(((cos(gama)*sin(beta)*cos(alpha))+(sin(gama)*sin(alpha)))*Pz1)-Bx1)^2+(Y+((sin(gama)*cos(beta))*Px1)+(((sin(gama)*sin(beta)*sin(alpha))+(cos(gama)*cos(alpha)))*Py1)+(((sin(gama)*sin(beta)*cos(alpha))-(cos(gama)*sin(alpha)))*Pz1)-By1)^2+(Z+((-sin(beta))*Px1)+((cos(beta)*sin(alpha))*Py1)+((cos(beta)*cos(alpha))*Pz1)-Bz1)^2==4.65,(X+((cos(gama)*cos(beta))*Px2)+(((cos(gama)*sin(beta)*sin(alpha))-(sin(gama)*cos(alpha)))*Py2)+(((cos(gama)*sin(beta)*cos(alpha))+(sin(gama)*sin(alpha)))*Pz2)-Bx2)^2+(Y+((sin(gama)*cos(beta))*Px2)+(((sin(gama)*sin(beta)*sin(alpha))+(cos(gama)*cos(alpha)))*Py2)+(((sin(gama)*sin(beta)*cos(alpha))-(cos(gama)*sin(alpha)))*Pz2)-By2)^2+(Z+((-sin(beta))*Px2)+((cos(beta)*sin(alpha))*Py2)+((cos(beta)*cos(alpha))*Pz2)-Bz2)^2==7.75,(X+((cos(gama)*cos(beta))*Px3)+(((cos(gama)*sin(beta)*sin(alpha))-(sin(gama)*cos(alpha)))*Py3)+(((cos(gama)*sin(beta)*cos(alpha))+(sin(gama)*sin(alpha)))*Pz3)-Bx3)^2+(Y+((sin(gama)*cos(beta))*Px3)+(((sin(gama)*sin(beta)*sin(alpha))+(cos(gama)*cos(alpha)))*Py3)+(((sin(gama)*sin(beta)*cos(alpha))-(cos(gama)*sin(alpha)))*Pz3)-By3)^2+(Z+((-sin(beta))*Px3)+((cos(beta)*sin(alpha))*Py3)+((cos(beta)*cos(alpha))*Pz3)-Bz3)^2==8.15);
A_x=A.X;
A_z=A.Z;
A_beta=A.beta;
Star Strider
Star Strider on 15 May 2021
As always, my pleasure!
To decrease the decimals, first use:
digits(5) % Five Digits (Choose Appropriate Number)
and suppress the output by putting a semicolon at the end of the original ‘A’, then either:
A = solve ( ... );
A_x = vpa(A.x)
A_z = vpa(A.z)
A_beta = vpa(A.beta)
or:
A_x = vpa(simplify(A.x, 500))
A_z = vpa(simplify(A.z, 500))
A_beta = vpa(simplify(A.beta, 500))
to simplify the individual results of ‘A’ (which could take a few minutes for a long expression) and then output the results. If the solutions are purely numeric (do not contain symbolic variables), you can then convert them to double values with the double function to use in other parts of the code, if necessary.
I would still go for the numeric solution, choosing as many different initial parameter estimates as necessary to provide an acceptable description of the system. Ths symbolic solution is likely to be unwieldy.
.

Sign in to comment.

More Answers (0)

Categories

Find more on Linear Algebra in Help Center and File Exchange

Tags

Products

Community Treasure Hunt

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

Start Hunting!