Using a function and fsolve with a for loop

1 次查看(过去 30 天)
I'd like to solve the system of 4 non linear equations and 4 unknonws over a span of values for T.
I'm recieving the following error:
"This statement is not inside any function.
(It follows the END that terminates the definition of the function "eqn".)"
function F = eqn(x,T)
F(1) = 10.*(T-x(1))+0.2.*(T.^4-x(1).^4) - x(4);
F(2) = -1.*(1/x(1)-1/x(2)) - x(4);
F(3) = 2.*(x(3)-x(2))-x(4);
F(4) = 5.*(1-x(3))-x(4);
end
T = 1:1:40;
for i = 1:numel(T)
x(i,:)=fsolve(@(x)eqn(x,T(i)),[0,0,0,0]);
end

回答(2 个)

Steven Lord
Steven Lord 2023-2-10
Either remove the section of code starting with "T = 1:1:40;" from your function file and run it in the Command Window or another file, or move that section of code to before the line "function F = eqn(x,T)". This second approach may require you to change the name of your function or the name of the file that contains it.
  3 个评论
Torsten
Torsten 2023-2-10
Yes, you divide by 0 in the expressions 1/x(1) and 1/x(2) if you set [0 0 0 0] as initial guess.
Alex Sha
Alex Sha 2023-2-11
See the results below:
T x1 x2 x3 x4
1 1 1 1 -4.41868100791151E-15
2 -0.04221805467944030 -15.5355259373919 -3.72443598211183 23.6221799105584
3 -0.02152989066740660 -31.4907092033764 -8.28305977239318 46.4152988619656
4 -0.01094985711415290 -62.9166489977652 -17.2618997136447 91.309498568229
5 -0.00571215265706289 -121.539985068469 -34.0114243052774 175.057121526389
6 -0.00313248055977908 -222.461927363911 -62.8462649611353 319.231324805659
7 -0.00181745226754944 -384.152722165812 -109.043634904507 550.218174522591
8 -0.00111208392250994 -628.447784587467 -178.842224167857 899.211120839266
9 -0.000713160880548333 -980.544992126182 -279.441426321773 1402.2071316088
10 -0.000476189242038096 -1469.00333332469 -419.000952378488 2100.00476189245
11 -0.000329141847745449 -2125.74230399299 -606.640658283718 3038.20329141854
12 -0.000234345559697521 -2986.04164041889 -852.440468691099 4267.20234345554
13 -0.000171168338151485 -4088.54119817842 -1167.44034233668 5842.20171168344
14 -0.000127824907104757 -5475.2408947744 -1563.64025564984 7823.20127824912
15 -9.73235904377666E-5 -7191.50068126512 -2054.00019464718 10275.0009732358
16 -7.53738494234600E-5 -9286.04052761695 -2652.44015074771 13267.2007537385
17 -5.92620663416646E-5 -11810.9404148345 -3373.84011852413 16874.2005926207
18 -4.72250545218803E-5 -14821.6403305754 -4234.04009445011 21175.2004722506
19 -3.80891432006883E-5 -18376.940266624 -5249.84007617831 26254.2003808915
20 -3.10559002788023E-5 -22539.0002173913 -6439.00006211179 32200.000310559
21 -2.55713925803802E-5 -27373.3401789998 -7820.2400511428 39106.200255714
22 -2.12444126098758E-5 -32948.8401487108 -9413.24004248879 47071.2002124441
23 -1.77941641615093E-5 -39337.7401245591 -11238.6400355883 56198.2001779416
24 -1.50160972175629E-5 -46615.6401051127 -13318.0400300322 66595.2001501609
25 -1.27591706301686E-5 -54861.5000893142 -15674.0000255183 78375.0001275917
26 -1.09104556930710E-5 -64157.6400763732 -18330.0400218209 91655.2001091045
27 -9.38454289762297E-6 -74589.7400656918 -21310.6400187691 106558.200093845
28 -8.11614527938351E-6 -86246.8400568131 -24641.2400162323 123211.200081161
29 -7.05486284239850E-6 -99221.340049384 -28348.2400141097 141746.200070549
30 -6.16142944895955E-6 -113609.00004313 -32459.0000123228 162300.000061614
31 -5.40499053405699E-6 -129508.940037835 -37001.840010810 185014.20005405
32 -4.76110670849903E-6 -147023.640033328 -42006.0400095222 210035.200047611
33 -4.21027458483498E-6 -166258.940029472 -47501.8400084205 237514.200042103
34 -3.73682023443470E-6 -187324.040026158 -53520.4400074736 267607.200037368
35 -3.32806389840558E-6 -210331.500023296 -60094.0000066561 300475.000033281
36 -2.97368408501804E-6 -235397.240020816 -67255.6400059474 336283.200029737
37 -2.66522957466352E-6 -262640.540018657 -75039.4400053305 375202.200026652
38 -2.39574209533138E-6 -292184.04001677 -83480.4400047915 417407.200023957
39 -2.15946248375913E-6 -324153.740015116 -92614.640004319 463078.200021595
40 -1.95160031217110E-6 -358679.000013661 -102479.000003903 512400.000019516

请先登录,再进行评论。


Torsten
Torsten 2023-2-10
编辑:Torsten 2023-2-10
Look whether this might help.
syms T
x = sym('x',[1 4]);
eqn1 = 10.*(T-x(1))+0.2.*(T.^4-x(1).^4) - x(4) == 0;
eqn2 = -1.*(1/x(1)-1/x(2)) - x(4) == 0;
eqn3 = 2.*(x(3)-x(2))-x(4) == 0;
eqn4 = 5.*(1-x(3))-x(4) == 0;
solve([eqn1,eqn2,eqn3,eqn4],x)
Warning: Solutions are only valid under certain conditions. To include parameters and conditions in the solution, specify the 'ReturnConditions' value as 'true'.
ans = struct with fields:
x1: ((21*T^4)/100 + (21*T)/2 + 8)*root(z^9 - z^8*(10*T + T^4/5 + 40/7) + z^7*((400*T)/7 + (8*T^4)/7 + 390/49) - z^6*((8800*T)/49 + (176*T^4)/49 - 7200/343) + z^5*((124000*T)/343 + (2480*T^4)/343 - 273500/2401) - z^4*((1234000*T)/2401 + (24680*T^4)… x2: 1 - (7*root(z^9 - z^8*(10*T + T^4/5 + 40/7) + z^7*((400*T)/7 + (8*T^4)/7 + 390/49) - z^6*((8800*T)/49 + (176*T^4)/49 - 7200/343) + z^5*((124000*T)/343 + (2480*T^4)/343 - 273500/2401) - z^4*((1234000*T)/2401 + (24680*T^4)/2401 - 2987401/12005) … x3: 1 - root(z^9 - z^8*(10*T + T^4/5 + 40/7) + z^7*((400*T)/7 + (8*T^4)/7 + 390/49) - z^6*((8800*T)/49 + (176*T^4)/49 - 7200/343) + z^5*((124000*T)/343 + (2480*T^4)/343 - 273500/2401) - z^4*((1234000*T)/2401 + (24680*T^4)/2401 - 2987401/12005) + z… x4: root(z^9 - z^8*(10*T + T^4/5 + 40/7) + z^7*((400*T)/7 + (8*T^4)/7 + 390/49) - z^6*((8800*T)/49 + (176*T^4)/49 - 7200/343) + z^5*((124000*T)/343 + (2480*T^4)/343 - 273500/2401) - z^4*((1234000*T)/2401 + (24680*T^4)/2401 - 2987401/12005) + z^3*(…
It shows that x(1)-x(4) are related to the roots of a polynomial of degree 9 in T.

类别

Help CenterFile Exchange 中查找有关 Loops and Conditional Statements 的更多信息

Community Treasure Hunt

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

Start Hunting!

Translated by