How to find angle in equcation

Question

Ali 2021-6-18

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/859655-how-to-find-angle-in-equcation

评论： Ali 2021-6-22

采纳的回答： Walter Roberson

此个问题已被 Walter Roberson 标记

显示标记

I want to find omega value in the equation.

can i anybody tell me how to solve this equation.

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Walter Roberson 2021-6-21

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/859655-how-to-find-angle-in-equcation#answer_729345

The search turns out to take over 20 seconds for each eload value, and even then it seems to fail.

Notice I ended at 0.002, which is the most I could do with this online version without timing out.

tic

syms omega

a=3.666;

k=0.467;

s=100;

r=1.4;

ef=0.08483;

N=(a*k)/(1+(a*k));

eloadvals = 0:0.001:0.002;

for eloadidx = 1 : length(eloadvals)

eload = eloadvals(eloadidx);

c1=(1-N)*(((r*a)/2)+(r/2)-1)-(r-1);

c2=(N*(r/(2*a))+(r/2)-1);

Ac2=((r*a)*(1-N))/2;

Ah1=(r*N)/(2*a);

A=-8*(1-N)*c2*sin(pi-atan(s/((1-N)*omega)))+8*N*c1*sin(2*(pi-atan(1/(N*omega))))+(32*N*ef+32*N*eload);

B=16*(1-N)*(-1*(((1-N)*omega)/(2*s))*(((Ah1*c1*sin(pi-atan(1/(N*omega)))*cos(pi-atan(1/(N*omega))))+((((N*s)/((1-N)*omega))+(N*(r/2-1)))*c2*cos((atan((((1-N)*((r/2)-1))/(s))*omega))+(pi-atan(s/((1-N)*omega))))*cos(pi-atan(s/((1-N)*omega)))))))*c2*sin(2*(pi-atan(s/((1-N)*omega))))+(16*(((N*omega)/2)*((((((1-N)/(N*omega))+((1-N)*(((r/2)-1)-(r-1))))*c1*cos((atan((((1-N)*((r/2)-1)-(r-1))/(1-N))*(N*omega)))+(pi-atan(1/(N*omega))))*cos(atan((((1-N)*((r/2)-1)-(r-1))/(1-N))*(N*omega))))+(Ac2*c2*sin(pi-atan(s/((1-N)*omega)))*cos(pi-atan(s/((1-N)*omega)))))))+2*c1)*N*c1*sin(2*(pi-atan(1/(N*omega))))+N*c1*c1*sin(4*(pi-atan(1/(N*omega))));

thissol = vpasolve(A==B);

if length(thissol) >= 1

sol(eloadidx) = thissol(1);

else

sol(eloadidx) = nan;

end

toc

Elapsed time is 54.511198 seconds.

plot(eloadvals, sol)

Nothing visibile in the plot because vpasolve() failed.

7 个评论
显示 5更早的评论隐藏 5更早的评论

Walter Roberson 2021-6-21

I cannot quite get three plots without timing out, so the plots will be split between two posts.

syms omega

Q = @(v) sym(v);

a = Q(3.666);

k = Q(0.467);

s = Q(100);

r = Q(1.4);

ef = Q(0.08483);

N = (a*k)/(1+(a*k));

syms eload

Pi = Q(pi);

c1 = (1-N)*(((r*a)/2)+(r/2)-1)-(r-1);

c2 = (N*(r/(2*a))+(r/2)-1);

Ac2 = ((r*a)*(1-N))/2;

Ah1 = (r*N)/(2*a);

A = -8*(1-N)*c2*sin(Pi-atan(s/((1-N)*omega)))+8*N*c1*sin(2*(Pi-atan(1/(N*omega))))+(32*N*ef+32*N*eload);

B = 16*(1-N)*(-1*(((1-N)*omega)/(2*s))*(((Ah1*c1*sin(Pi-atan(1/(N*omega)))*cos(Pi-atan(1/(N*omega))))+((((N*s)/((1-N)*omega))+(N*(r/2-1)))*c2*cos((atan((((1-N)*((r/2)-1))/(s))*omega))+(Pi-atan(s/((1-N)*omega))))*cos(Pi-atan(s/((1-N)*omega)))))))*c2*sin(2*(Pi-atan(s/((1-N)*omega))))+(16*(((N*omega)/2)*((((((1-N)/(N*omega))+((1-N)*(((r/2)-1)-(r-1))))*c1*cos((atan((((1-N)*((r/2)-1)-(r-1))/(1-N))*(N*omega)))+(Pi-atan(1/(N*omega))))*cos(atan((((1-N)*((r/2)-1)-(r-1))/(1-N))*(N*omega))))+(Ac2*c2*sin(Pi-atan(s/((1-N)*omega)))*cos(Pi-atan(s/((1-N)*omega)))))))+2*c1)*N*c1*sin(2*(Pi-atan(1/(N*omega))))+N*c1*c1*sin(4*(Pi-atan(1/(N*omega))));

sol = solve(A-B, eload)

sol =

double(limit(sol, omega, -inf))

ans = -0.0848

double(limit(sol, omega, inf))

ans = -0.0848

fplot(sol, [-10 10]); xlabel('omega'); ylabel('eload')

%fplot(sol, [1.35 1.5]); xlabel('omega'); ylabel('eload')

%fplot(sol, [-1.5 -1.35]); xlabel('omega'); ylabel('eload')

Walter Roberson 2021-6-21

syms omega

Q = @(v) sym(v);

a = Q(3.666);

k = Q(0.467);

s = Q(100);

r = Q(1.4);

ef = Q(0.08483);

N = (a*k)/(1+(a*k));

syms eload

Pi = Q(pi);

c1 = (1-N)*(((r*a)/2)+(r/2)-1)-(r-1);

c2 = (N*(r/(2*a))+(r/2)-1);

Ac2 = ((r*a)*(1-N))/2;

Ah1 = (r*N)/(2*a);

A = -8*(1-N)*c2*sin(Pi-atan(s/((1-N)*omega)))+8*N*c1*sin(2*(Pi-atan(1/(N*omega))))+(32*N*ef+32*N*eload);

B = 16*(1-N)*(-1*(((1-N)*omega)/(2*s))*(((Ah1*c1*sin(Pi-atan(1/(N*omega)))*cos(Pi-atan(1/(N*omega))))+((((N*s)/((1-N)*omega))+(N*(r/2-1)))*c2*cos((atan((((1-N)*((r/2)-1))/(s))*omega))+(Pi-atan(s/((1-N)*omega))))*cos(Pi-atan(s/((1-N)*omega)))))))*c2*sin(2*(Pi-atan(s/((1-N)*omega))))+(16*(((N*omega)/2)*((((((1-N)/(N*omega))+((1-N)*(((r/2)-1)-(r-1))))*c1*cos((atan((((1-N)*((r/2)-1)-(r-1))/(1-N))*(N*omega)))+(Pi-atan(1/(N*omega))))*cos(atan((((1-N)*((r/2)-1)-(r-1))/(1-N))*(N*omega))))+(Ac2*c2*sin(Pi-atan(s/((1-N)*omega)))*cos(Pi-atan(s/((1-N)*omega)))))))+2*c1)*N*c1*sin(2*(Pi-atan(1/(N*omega))))+N*c1*c1*sin(4*(Pi-atan(1/(N*omega))));

sol = solve(A-B, eload)

sol =

double(limit(sol, omega, -inf))

ans = -0.0848

double(limit(sol, omega, inf))

ans = -0.0848

%fplot(sol, [-10 10]); xlabel('omega'); ylabel('eload')

fplot(sol, [1.35 1.5]); xlabel('omega'); ylabel('eload')

fplot(sol, [-1.5 -1.35]); xlabel('omega'); ylabel('eload')

Walter Roberson 2021-6-21

So look at the top plot. The left and right halfs are symmetric around eload = -0.0848

You can see near the discontinuity at 0 that the left side close to 0 shares values with the right side close to 0, so that section is not unique.

If we look down the discontinuity on the left side close to 0, to the point where the one on the right close to 0 gives out, and then mentally follow that level near -0.1 to the left, we can see that we are not unique, there is the other side of the valley that rises up again. How high does the left valley rise? Well it rises to double(limit(sol, omega, -inf)) = -0.0848, and by following across towards the discontinuity at 0, we can see that that entire area gets us to a matching edge just to the left of 0.

What point does not have a duplicate? Only the lowest part of the valley, at approximately omega = -1.42 corresponding to eload of about -0.15342 : if you follow that one point left and right, there is no corresponding point on either side of the discontinuity.

Likewise, when we look on the right side of the discontinuity, we have a balance point at about 1.42 (same distance away from 0 as the other one), corresponding to eload of about -0.01403 .

Those are the only two omega that have unique solutions: +/- 1.42 (approximately) corresponding to eload of about -0.01403 and eload of about -0.15342 .

If there are any stable points, it would be only those two.

If omega were being controlled by some force pushing towards a balance, then if you were to disturb omega a little, it should rebalance itself. But if not... well, if there is no restoring pressure on omega, then if the system gets knocked slightly so that omega varies, then you would start getting two eload solutions.

Ali 2021-6-22