Problems with a Newton-Raphson code

7 次查看（过去 30 天）
beau genes
beau genes 2015-4-6
回答： Jan 2015-4-6
function pnew = newtraph(p4,p1,d,p,e,f)
p1=20;
p4=500;
k1=1.4;
k4=1.4;
r1=287;
r4=287;
c1=348.92;
c4=348.92;
p=1;
dp=.0000001;
d=((k4-1)/(2*k1))*(c1/c4)
e= (k1+1)/(2*k1)
f= (2*k4)/(k4-1)
while abs(f(p)) > 1e-6
fp =(p4/p1)*(1-((d*(p-1)/(1+e*(p-1))))^f)-p
fpmdp= (p4/p1)*(1-((d*((p-dp)-1)/(1+e*((p-dp)-1)))))^f-(p-dp)
fppdp= (p4/p1)*(1-((d*((p+dp)-1)/(1+e*((p+dp)-1)))))^f-(p+dp)
dfdp= (fppdp-fpmdp)/2*dp
pnew=fp-(fp/dfdp)
p=pnew
end
I am trying to use a newton raphson method to find pnew for a shock tube in my gas dynamics class. I get the error message below.
Subscript indices must either be real positive integers or logicals.
Error in Untitled (line 34) end
With the given values from above, the iteration should produce a pnew=4.0471
Any help would be much appreciated, thanks

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

回答（2 个）

Brendan Hamm
Brendan Hamm 2015-4-6
You have a line:
while abs(f(p)) > 1e-6
which is an indexing into the variable f , with the indexing variable p. But p changes on every iteration in the while loop. It changes to the value 9.230769236418623e+13 on the first loop and you get an error when you try and use this value as an index into f at the while loop line.
Jan
Jan 2015-4-6
The debugger helps you to indentify the source of problems. Set a breakpoint in your code and step through the function line by line. Then you will find out, what Brendan has mentioned already.

类别

在 Help Center 和 File Exchange 中查找有关 Newton-Raphson Method 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by