![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1150145/image.png)
My code is outputting very illogical answers for the second block (my second vector loop equation).
3 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
Emily Friedman
il 30 Set 2022
Modificato: David Goodmanson
il 9 Ott 2022
The first block of code (E1) runs and comes out with a correct answer (validated in excel). It has only one input and two outputs. The second block (E2) of code runs, but comes out with an illogical answer. This block has two inputs that are dependent on the previous block's outputs. I'm expecting the result to be some sort of sine or cosine wave as this code is the position analysis of an ornithoptor.
"fsolve" in the in E2 also says "no solution found", but provides an answer. I do not understand why.
Why won't my code work properly for the E2 block?
% given values
R1 = 59.7;
R2 = 18.8;
R3 = 41.0;
R4 = 40.1;
R5 = 136.4;
R6 = 11.9;
R7 = 15.6;
R8 = 11.7;
Rgi = 13.1;
Rec = 12.7;
Rfc = 137.2;
theta1 = 0;
theta2_all = linspace(0, 2*pi, 500); % all possible input values
% E1
for i = 1:500
VLE1 = @(x0) [(R2.*cos(theta2_all(i)) + R3.*cos(x0(1)) - R4.*cos(x0(2)) - R1.*cos(theta1));
(R2.*sin(theta2_all(i)) + R3.*sin(x0(1)) - R4.*sin(x0(2)) - R1.*sin(theta1))];
x0 = [pi/2; pi/2];
x = fsolve(VLE1, x0);
theta3(i) = x(1);
theta4(i) = x(2);
end
theta3
theta4
% E2
theta_ec_off = 16.5*(pi/180);
theta_fc_off = 8.5*(pi/180);
theta_ec = theta3 + theta_ec_off;
theta_fc = theta4 + theta_fc_off;
for i = 1:500
VLE2 = @(x0) [(Rec*cos(theta_ec(i)) + Rfc*cos(theta_fc(i)) - R6*cos(x0(2)) - R5*cos(x0(1)));
(Rec*sin(theta_ec(i)) + Rfc*sin(theta_fc(i)) - R6*cos(x0(2)) - R5*cos(x0(1)))];
x0 = [pi/2; pi/2];
x = fsolve(VLE2, x0);
theta5(i) = x(1);
theta6(i) = x(2); % storing values in a vector
end
theta5
theta6
% validating
plot(theta2_all, theta5)
hold on
plot(theta2_all, theta6)
0 Commenti
Risposta accettata
David Goodmanson
il 9 Ott 2022
Modificato: David Goodmanson
il 9 Ott 2022
Hi Emily,
No one seems to have mentioned that there are two independent solutions for the angle pair angle3 and angle4, since if a side of S of the triangle S R3 R4 lies along, say, the x axis, then the triangle can lie either above or below the x axis, as determined by the angle between S and R3 being either positive or negative.
For each of those solutions there are two independent solutions for the angle pair angle5 and angle6, so there are four independent solutions in all for the set of angles 3 through 6. Of course only one of those sets will correspond to the real world geometry that you have.
Not having access to fsolve, the code below consists of algebraic vectorized solutions. The angles are stndard, positive angle being ccw from the x axis.
% VLE1 = @(x0) [(R2.*cos(theta2(i)) + R3.*cos(x0(1)) - R4.*cos(x0(2)) - R1.*cos(theta1));
% (R2.*sin(theta2(i)) + R3.*sin(x0(1)) - R4.*sin(x0(2)) - R1.*sin(theta1))];
% VLE2 = @(x0) [(Rec*cos(theta_ec(i)) + Rfc*cos(theta_fc(i)) - R6*cos(x0(2)) - R5*cos(x0(1)));
% (Rec*sin(theta_ec(i)) + Rfc*sin(theta_fc(i)) - R6*sin(x0(2)) - R5*sin(x0(1)))];
R1 = 59.7;
R2 = 18.8;
R3 = 41.0;
R4 = 40.1;
R5 = 136.4;
R6 = 11.9;
R7 = 15.6;
R8 = 11.7;
Rgi = 13.1;
Rec = 12.7;
Rfc = 137.2;
theta1 = 0;
theta2 = linspace(0, 2*pi, 500); % all possible input values
eps1 = [1 1 -1 -1];
eps2 = [1 -1 1 -1];
for n = 1:4 % loop over solutions
% E1
R7 = abs(R2*iexp(theta2)-R1*iexp(theta1));
theta7 = angle(R2*iexp(theta2)-R1*iexp(theta1));
u1 = eps1(n)*acos((R4^2-R7.^2-R3^2)./(2*R7*R3)); % pos or neg
theta3 = u1 + theta7;
theta4 = angle(R7 +R3*iexp(u1)) + theta7;
theta3 = unwrap(theta3);
theta4 = unwrap(theta4);
theta7 = unwrap(theta7);
% check complex version of VLE1, should be small
max(abs(R7.*iexp(theta7) + R3*iexp(theta3)-R4*iexp(theta4)))
% E2
theta_ec_off = 16.5*(pi/180);
theta_fc_off = 8.5*(pi/180);
theta_ec = theta3 + theta_ec_off;
theta_fc = theta4 + theta_fc_off;
R8 = abs(Rec*iexp(theta_ec) + Rfc*iexp(theta_fc));
theta8 = angle(Rec*iexp(theta_ec) + Rfc*iexp(theta_fc));
u2 = eps2(n)*acos((R8.^2+R6^2-R5^2)./(2*R8*R6)); % pos or neg
theta6 = u2 + theta8;
theta5 = angle(R8 -R6*iexp(u2)) + theta8;
theta5 = unwrap(theta5);
theta6 = unwrap(theta6);
theta8 = unwrap(theta8);
% check complex version of VLE2, should be small
max(abs(R8.*iexp(theta8) - R5*iexp(theta5) - R6*iexp(theta6)))
figure(1)
subplot(2,2,n)
plot(theta2,[theta3;theta4;theta5;theta6])
grid on
end
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1150145/image.png)
0 Commenti
Più risposte (2)
David Goodmanson
il 1 Ott 2022
Modificato: David Goodmanson
il 6 Ott 2022
Hi Emily,
If VLE1 and VLE2 are lined up typographically, you have
VLE1 = @(x0) [(R2.*cos(theta2_all(i)) + R3.*cos(x0(1)) - R4.*cos(x0(2)) - R1.*cos(theta1));
(R2.*sin(theta2_all(i)) + R3.*sin(x0(1)) - R4.*sin(x0(2)) - R1.*sin(theta1))];
VLE2 = @(x0) [(Rec*cos(theta_ec(i)) + Rfc*cos(theta_fc(i)) - R6*cos(x0(2)) - R5*cos(x0(1)));
(Rec*sin(theta_ec(i)) + Rfc*sin(theta_fc(i)) - R6*cos(x0(2)) - R5*cos(x0(1)))];
Looking at just terms involving R3 and R4, VLE1 has cosines on the first line, sines on the second line. But for R5 and R6, VLE2 has cosines on both lines, meaning that the first line and the second line have identical terms. This does not seem right.
Torsten
il 6 Ott 2022
Modificato: Torsten
il 6 Ott 2022
% given values
R1 = 59.7;
R2 = 18.8;
R3 = 41.0;
R4 = 40.1;
R5 = 136.4;
R6 = 11.9;
R7 = 15.6;
R8 = 11.7;
Rgi = 13.1;
Rec = 12.7;
Rfc = 137.2;
theta1 = 0;
theta2_all = linspace(0, 2*pi, 500); % all possible input values
options = optimset('Display','off') ;
x0 = [pi/2 pi/2];
% E1
for i = 1:500
VLE1 = @(x0) [(R2.*cos(theta2_all(i)) + R3.*cos(x0(1)) - R4.*cos(x0(2)) - R1.*cos(theta1));
(R2.*sin(theta2_all(i)) + R3.*sin(x0(1)) - R4.*sin(x0(2)) - R1.*sin(theta1))];
x = fsolve(VLE1, x0, options);
x0 = x;
theta3(i) = x(1);
theta4(i) = x(2);
end
figure(1)
hold on
plot(theta2_all, theta3)
plot(theta2_all, theta4)
hold off
% E2
theta_ec_off = 16.5*(pi/180);
theta_fc_off = 8.5*(pi/180);
theta_ec = theta3 + theta_ec_off;
theta_fc = theta4 + theta_fc_off;
x0 = [pi/2 pi/2];
for i = 1:500
VLE2 = @(x0) [(Rec*cos(theta_ec(i)) + Rfc*cos(theta_fc(i)) - R6*cos(x0(2)) - R5*cos(x0(1)));
(Rec*sin(theta_ec(i)) + Rfc*sin(theta_fc(i)) - R6*sin(x0(2)) - R5*sin(x0(1)))];
x = fsolve(VLE2, x0, options);
x0 = x;
theta5(i) = x(1);
theta6(i) = x(2); % storing values in a vector
end
figure(2)
hold on
plot(theta2_all, theta5)
plot(theta2_all, theta6)
hold off
0 Commenti
Vedere anche
Categorie
Scopri di più su Math Operations in Help Center e File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!