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

6 Ansichten (letzte 30 Tage)
Isaac
Isaac am 12 Aug. 2015
Kommentiert: Walter Roberson am 14 Aug. 2015
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 Kommentare
Isaac
Isaac am 13 Aug. 2015
Sorry
SXY = 125*[1 2 3]; C0 = 100*[1 2 3];
Thanks

Melden Sie sich an, um zu kommentieren.

Antworten (3)

Isaac
Isaac am 13 Aug. 2015
Sorry
SXY = 125*[1 2 3]; C0 = 100*[1 2 3];
Thanks

Walter Roberson
Walter Roberson am 13 Aug. 2015
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 Kommentare
Walter Roberson
Walter Roberson am 13 Aug. 2015
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 am 14 Aug. 2015

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.

Melden Sie sich an, um zu kommentieren.


Isaac
Isaac am 13 Aug. 2015
Thanks Walter...yes, the solutions are all imaginary with the current inputs

Community Treasure Hunt

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

Start Hunting!

Translated by