while loop running forever during attempt at Newton Rhapson method

1 回表示 (過去 30 日間)
Linnea
Linnea 2023 年 11 月 18 日
コメント済み: Torsten 2023 年 11 月 18 日
This code keeps running seemingly forever, or at least a few minutes. I can't figure out why, could someone please explain what is going wrong? unless i wrote the actual functions of d wrong it should be expected to converge towards something until the while loop condition isn't fullfilled anymore
d=0.007;
dI=1.219;
Rfi=1.76*(10^(-4));
Rf0=Rfi;
hs=356;
ht=hs;
kw=60;
dTm=29.6;
Q=801368;
Aex=64.15;
a=((1/ht)+Rfi)/dI;
b=log(d/dI)/2*kw;
c=d*(a+b)+Rf0+(1/hs);
f=(c^(-1))-Q/(Aex*dTm);
df=-1*(c^(-2))*(a+b)*(dI/(d*2*kw))*(1/dI);
%d1-do ska va mindre än tol
%i slutet av iterationen
%så man kan sätta while
d=0.007;
d0=d;
d1=10;
abs(d1-d0);
while abs(d1-d0)>10^(-8);
d0=d;
d1=d0-f/df;
d=d1;
end
d1

採用された回答

VBBV
VBBV 2023 年 11 月 18 日
編集済み: VBBV 2023 年 11 月 18 日
d=0.007;
dI=1.219;
Rfi=1.76*(10^(-4));
Rf0=Rfi;
hs=356;
ht=hs;
kw=60;
dTm=29.6;
Q=801368;
Aex=64.15;
a=((1/ht)+Rfi)/dI;
%d1-do ska va mindre än tol
%i slutet av iterationen
%så man kan sätta while
d=0.007;
d0=d;
d1=10;
abs(d1-d0);
while abs(d1-d0)>10^(-8);
b=log(d/dI)/(2*kw);
c=d*(a+b)+Rf0+(1/hs);
f=(c^(-1))-Q/(Aex*dTm);
df=-1*(c^(-2))*(a+b)*(dI/(d*2*kw))*(1/dI);
d0=d;
d1=d0-f/df;
d=d1;
end
d1
d1 = 0.0191
  2 件のコメント
VBBV
VBBV 2023 年 11 月 18 日
Place the lines of coe inside the while loop which are function of variables d and expression for evaluating d1 inside the while loop
b=log(d/dI)/(2*kw);
c=d*(a+b)+Rf0+(1/hs);
f=(c^(-1))-Q/(Aex*dTm);
df=-1*(c^(-2))*(a+b)*(dI/(d*2*kw))*(1/dI);
Torsten
Torsten 2023 年 11 月 18 日
Are you sure it must read
b=log(d/dI)/(2*kw);
instead of
b=log(d/dI)/2*kw;
as in the original code ?

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

その他の回答 (1 件)

Torsten
Torsten 2023 年 11 月 18 日
編集済み: Torsten 2023 年 11 月 18 日
Your function doesn't seem to have a zero.
d=0.007;
dI=1.219;
Rfi=1.76*(10^(-4));
Rf0=Rfi;
hs=356;
ht=hs;
kw=60;
dTm=29.6;
Q=801368;
Aex=64.15;
a=((1/ht)+Rfi)/dI;
b=@(d)log(d/dI)/2*kw;
c=@(d)d.*(a+b(d))+Rf0+(1/hs);
f=@(d)(c(d).^(-1))-Q/(Aex*dTm);
D = 0.1:0.01:10;
plot(D,f(D))
df=@(d)-1*(c(d)^(-2))*(a+b(d))*(dI/(d*2*kw))*(1/dI);
%d1-do ska va mindre än tol
%i slutet av iterationen
%så man kan sätta while
d=0.007;
d0=d;
d1=10;
abs(d1-d0);
while abs(d1-d0)>10^(-8);
d0=d;
%df0 = (f(d0+1e-7)-f(d0))*1e7
d1=d0-f(d0)/df(d0);
% d1 = d0-f(d0)/df0;
d=d1
end
d = 2.6870
d = -2.3299e+07
d = 3.2216e+29 + 6.0366e+28i
d = -3.0746e+96 - 1.9209e+96i
d = -1.4912e+297 +1.6013e+298i
d = Inf + Infi
d = NaN + NaNi
d1
d1 = NaN + NaNi

カテゴリ

Help Center および File ExchangeData Type Identification についてさらに検索

製品


リリース

R2023b

Community Treasure Hunt

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

Start Hunting!

Translated by