Fminsearch is not giving a correct answer. Why?
1 次查看(过去 30 天)
显示 更早的评论
Temp (temperature) value must be equal T value that this script answered.
Temp=550 % Temperature
P=2.2 %pressure
rho =( P * 1000)./( R.* T) ; % Density
R =8.314510; % Joule / ( mol * Kelvin)
Mass =28.9669; % 28,9586; g / mol
Tj =132.6312; %Kelvin (Maxcondentherm temperature)
rhoj =10.4477; %mol/dm3 (Maxcondentherm density)
Pj =3.78502; %MPa (Maxcondentherm pressure)
%% CONSTANTS of FUNCTION%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
N1=0.605719400e-7; N2=-0.210274769e-4; N3=-0.158860716E-3;
N4=-13.841928076; N5=17.275266575; N6=-0.195363420E-3;
N7=2.490888032; N8=0.791309509; N9=0.212236768;
N10=-0.197938904; N11=25.36365; N12=16.90741;
N13=87.31279;
x=Tj./Temp; y=rho./rhoj;
a0 = log(y) + N1.*x.^(-3) + N2.* x.^(-2) + N3.* x.^(-1) + N4 + N5.* x + ...
N6.* x.^1.5 + N7.* log(x) + N8.* log(1 - exp(-N11.* x)) + N9.* log(1 - exp(-N12.* x)) + ...
N10.* log((2/3) + exp(N13.* x));
x_da0_dx = -x.*((3*N1)./x.^4 - N5 + (2*N2)./x.^3 + N3./x.^2 - N7./x - (3*N6.*x.^(1/2))/2 + ...
(N8*N11*exp(-N11.*x))./(exp(-N11.*x) - 1) + (N9*N12*exp(-N12.*x))./(exp(-N12.*x) - 1) - ...
(N10*N13*exp(N13.*x))./(exp(N13.*x) + 2/3));
s =(R.*(-a0-x_da0_dx+1))./Mass % kj / kg Kelvin %% THIS VALUE IS ENTROPY. WE FOUND ENTROPY WITH FIRST VALUES OF TEMPERATURE AND PRESSURE
eqn =@(x,s,P,R,Mass,Tj,rhoj,N1,N2,N3,N4,N5,N6,N7,N8,N9,N10,N11,N12,N13)(R/Mass)*(-(log(( x * P * 1000)./( R* Tj*rhoj)) + N1.*x.^(-3) + N2.* x.^(-2) + N3.* x.^(-1) + N4 + N5.* x + ...
N6.* x.^1.5 + N7.* log(x) + N8.* log(1 - exp(-N11.* x)) + N9.* log(1 - exp(-N12.* x)) + ...
N10.* log((2/3) + exp(N13.* x)))-(-x.*((3*N1)./x.^4 - N5 + (2*N2)./x.^3 + N3./x.^2 - N7./x - (3*N6.*x.^(1/2))/2 + ...
(N8*N11*exp(-N11.*x))./(exp(-N11.*x) - 1) + (N9*N12*exp(-N12.*x))./(exp(-N12.*x) - 1) - ...
(N10*N13*exp(N13.*x))./(exp(N13.*x) + 2/3)))+1)-s;
% BUT WHEN WE TRY TO GET TEMPRATURE FROM THE FIRST ENTROPY AND PRESSURE VALUE , ANSWER İS NOT CORRECT. WE HAVE TO FOUR TIMES TO CLICK "RUN" BUTTON. WHY?
options = optimset('Display','iter','MaxFunEvals',2e30,'TolFun',1e-20,'TolX',1e-20,'NoStopIfFlatInfeas','on')
x0=5
x = fminsearch(@(x)norm(eqn(x,s,P,R,Mass,Tj,rhoj,N1,N2,N3,N4,N5,N6,N7,N8,N9,N10,N11,N12,N13)),x0,options)
T=Tj./x
eqn(x,s,P,R,Mass,Tj,rhoj,N1,N2,N3,N4,N5,N6,N7,N8,N9,N10,N11,N12,N13) % TO CHECK FUNCTION ERROR
0 个评论
回答(0 个)
另请参阅
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!