I will appreciate any suggestion on how I could have a solution to this.

10 ビュー (過去 30 日間)
Isaac
Isaac 2015 年 8 月 12 日
コメント済み: Walter Roberson 2015 年 8 月 14 日
SX=1000*[1 2 3];
SY=2000*[1.5 2 3];
SXY = 1258[1 2 3];
a = [0.3 0.6 0.9];
syms rb
for j=1:1:3
if pwmid(j)<=pwc(j)
SRR(j)=0.5*(SX(j)+SY(j)).*(1-(a(j).^2)/rb^2)+0.5*(SX(j)-SY(j)).*(1+(3*a(j).^4/rb^4)-(4*a(j).^2/rb^2))*cos(2*thbkso(j))...
+SXY(j).*(1+(3*a(j).^4/rb^4)-(4*a(j).^2/rb^2))*sin(2*thbkso(j))+(a(j).^2/rb^2).*pwmid(j);
STT(j)=0.5*(SX(j)+SY(j)).*(1+(a(j).^2)/rb^2)-0.5*(SX(j)-SY(j)).*(1+(3*a(j).^4/rb^4))*cos(2*thbkso(j))...
-SXY(j).*(1+(3*a(j).^4/rb^4))*sin(2*thbkso(j))-(a(j).^2/rb^2).*pwmid(j);
SRT(j)=(0.5*(SX(j)-SY(j)).*sin(2*thbkso(j))+SXY(j).*cos(2*thbkso(j))).*(1-(3*a(j).^4/rb^4)+(2*a(j).^2/rb^2));
SIGMA1A(j)=0.5*(STT(j)+SRR(j))+0.5*((STT(j)-SRR(j)).^2+4*SRT(j).^2).^0.5;
SIGMA3A(j)=0.5*(STT(j)+SRR(j))-0.5*((STT(j)-SRR(j)).^2+4*SRT(j).^2).^0.5;
C0FUN(j)=SIGMA1A(j)-SIGMA3A(j);
rbsoln{j}=double(vpasolve(C0FUN(j)==C0(j),rb));
cell(rbsoln);
rw(j)=rbsoln{j}(1);
rbkt_art(j) = rbsoln{j}(1)-a(j);
else
rw(j)=a(j);
rbkt_art(j)=rbkt_int(j);
end
end
  8 件のコメント
Isaac
Isaac 2015 年 8 月 13 日
Sorry
SXY = 125*[1 2 3]; C0 = 100*[1 2 3];
Thanks
Isaac
Isaac 2015 年 8 月 13 日
thbkso = pi/2*[1 1 1];

サインインしてコメントする。

回答 (3 件)

Isaac
Isaac 2015 年 8 月 13 日
Sorry
SXY = 125*[1 2 3]; C0 = 100*[1 2 3];
Thanks

Walter Roberson
Walter Roberson 2015 年 8 月 13 日
All solutions to those equations are strictly imaginary for the parameters you give.
For example, for j = 1, the solutions are
(3/10)*sqrt(5)*sqrt(roots([+8105,-9500,+4790,-1164,+117]))
and the negatives of those.
  2 件のコメント
Walter Roberson
Walter Roberson 2015 年 8 月 13 日
Please explain what you mean when you said you were concerned about solve or vpasolve "not giving favorable results" ?
If you want all of the results, then you may have to use solve() instead of vpasolve(), and you might have to double() the result of solve() to get numeric values. I do not have the Symbolic Toolbox so I cannot check exactly what would be returned.
Walter Roberson
Walter Roberson 2015 年 8 月 14 日

I could have made a mistake along the way, but if I got it right then:

for j = 1 : 3
  A = 32 * SXY(j) * sin(thbkso(j)) * cos(thbkso(j))^3 * SX(j) - 32 * SXY(j) * sin(thbkso(j)) * cos(thbkso(j))^3 * SY(j) - 16 * SXY(j) * sin(thbkso(j)) * cos(thbkso(j)) * SX(j) + 16 * SXY(j) * sin(thbkso(j)) * cos(thbkso(j)) * SY(j) - C0(j)^2 + SX(j)^2 - 2 * SX(j) * SY(j) + 4 * SXY(j)^2 + SY(j)^2;
    B =  - 32 * cos(thbkso(j))^4 * SX(j)^2 * a(j)^2 + 64 * cos(thbkso(j))^4 * SX(j) * SY(j) * a(j)^2 + 128 * cos(thbkso(j))^4 * SXY(j)^2 * a(j)^2 - 32 * cos(thbkso(j))^4 * SY(j)^2 * a(j)^2 - 8 * sin(thbkso(j)) * cos(thbkso(j)) * SX(j) * SXY(j) * a(j)^2 - 8 * sin(thbkso(j)) * cos(thbkso(j)) * SXY(j) * SY(j) * a(j)^2 + 16 * sin(thbkso(j)) * cos(thbkso(j)) * SXY(j) * a(j)^2 * pwmid(j) + 28 * cos(thbkso(j))^2 * SX(j)^2 * a(j)^2 - 64 * cos(thbkso(j))^2 * SX(j) * SY(j) * a(j)^2 + 8 * cos(thbkso(j))^2 * SX(j) * a(j)^2 * pwmid(j) - 128 * cos(thbkso(j))^2 * SXY(j)^2 * a(j)^2 + 36 * cos(thbkso(j))^2 * SY(j)^2 * a(j)^2 - 8 * cos(thbkso(j))^2 * SY(j) * a(j)^2 * pwmid(j) - 2 * SX(j)^2 * a(j)^2 + 8 * SX(j) * SY(j) * a(j)^2 - 4 * SX(j) * a(j)^2 * pwmid(j) + 16 * SXY(j)^2 * a(j)^2 - 6 * SY(j)^2 * a(j)^2 + 4 * SY(j) * a(j)^2 * pwmid(j);
    C = 128 * sin(thbkso(j)) * cos(thbkso(j))^3 * SX(j) * SXY(j) * a(j)^4 - 128 * sin(thbkso(j)) * cos(thbkso(j))^3 * SXY(j) * SY(j) * a(j)^4 + 48 * cos(thbkso(j))^4 * SX(j)^2 * a(j)^4 - 96 * cos(thbkso(j))^4 * SX(j) * SY(j) * a(j)^4 - 192 * cos(thbkso(j))^4 * SXY(j)^2 * a(j)^4 + 48 * cos(thbkso(j))^4 * SY(j)^2 * a(j)^4 - 48 * sin(thbkso(j)) * cos(thbkso(j)) * SX(j) * SXY(j) * a(j)^4 + 80 * sin(thbkso(j)) * cos(thbkso(j)) * SXY(j) * SY(j) * a(j)^4 - 32 * sin(thbkso(j)) * cos(thbkso(j)) * SXY(j) * a(j)^4 * pwmid(j) - 40 * cos(thbkso(j))^2 * SX(j)^2 * a(j)^4 + 96 * cos(thbkso(j))^2 * SX(j) * SY(j) * a(j)^4 - 16 * cos(thbkso(j))^2 * SX(j) * a(j)^4 * pwmid(j) + 192 * cos(thbkso(j))^2 * SXY(j)^2 * a(j)^4 - 56 * cos(thbkso(j))^2 * SY(j)^2 * a(j)^4 + 16 * cos(thbkso(j))^2 * SY(j) * a(j)^4 * pwmid(j) + 7 * SX(j)^2 * a(j)^4 - 18 * SX(j) * SY(j) * a(j)^4 + 4 * SX(j) * a(j)^4 * pwmid(j) - 8 * SXY(j)^2 * a(j)^4 + 15 * SY(j)^2 * a(j)^4 - 12 * SY(j) * a(j)^4 * pwmid(j) + 4 * a(j)^4 * pwmid(j)^2;
    E = 288 * SXY(j) * a(j)^8 * sin(thbkso(j)) * cos(thbkso(j))^3 * SX(j) - 288 * SXY(j) * a(j)^8 * sin(thbkso(j)) * cos(thbkso(j))^3 * SY(j) - 144 * SXY(j) * a(j)^8 * sin(thbkso(j)) * cos(thbkso(j)) * SX(j) + 144 * SXY(j) * a(j)^8 * sin(thbkso(j)) * cos(thbkso(j)) * SY(j) + 9 * SX(j)^2 * a(j)^8 - 18 * SY(j) * a(j)^8 * SX(j) + 36 * SXY(j)^2 * a(j)^8 + 9 * SY(j)^2 * a(j)^8;
    sols_plus = sqrt( roots([A, B, C, 0, E]) );
    sols{j} = [sols_plus; -sols_plus];
  end

I am not certain of these coefficients; I am concerned that the previous solution did not have a 0 in the x^1 position but this does.

サインインしてコメントする。


Isaac
Isaac 2015 年 8 月 13 日
Thanks Walter...yes, the solutions are all imaginary with the current inputs

カテゴリ

Help Center および File ExchangeOrdinary Differential Equations についてさらに検索

Community Treasure Hunt

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

Start Hunting!

Translated by