Value of the array in calculation and graph is different

k=10;

Es=0.05;

t=273;

TempK=0:1500;

cp=0:0.2:1.2;

cp=0.99403+((1.671e-4)*(TempK))+((9.7215e-8)*(TempK.^2))-((9.5838e-11)*(TempK.^3))+((1.9520e-14)*(TempK.^4));

plot(TempK,cp,'r',LineWidth=2);

title('Specific Heat Capacity, C_p (kJ/kg\cdotK) VS Temperature (K)');

xlabel('Temperature, T (K)');

ylabel('Specific Heat Capacity, C_p (kJ/kg\cdotK)');

grid on

tarr=zeros(1,k);

cp1arr=zeros(1,k);

for i=1:k

T=(((0.20597+((9.5838e-11)*(t.^3))-((9.7215e-8)*(t.^2))-((1.671e-4)*t)))/(1.9520e-14)).^(1/4);

tnew=T;

Ea=abs((((tnew-t)/tnew)*100));

cp1=0.99403+((1.671e-4)*(T))+((9.7215e-8)*(T.^2))-((9.5838e-11)*(T.^3))+((1.9520e-14)*(T.^4));

if Ea<Es

disp('---------------------------------------------');

disp(['Iteration,k = ',num2str(i,k)]);

disp(['Cp = ',num2str(cp1),' kJ/kgK'])

disp(['t = ',num2str(t)]);

disp(['t new = ',num2str(tnew)]);

disp(['Ea = ',num2str(Ea), ' %']);

break;

end

t=tnew;

disp('---------------------------------------------');

disp(['Iteration,k = ',num2str(i,k)]);

disp(['Cp = ',num2str(cp1),' kJ/kgK'])

disp(['t = ',num2str(t)]);

disp(['t new = ',num2str(tnew)]);

disp(['Ea = ',num2str(Ea),' %']);

cp1arr(1,i)=cp1;

tarr(1,i)=t;

i=i+1;

end

hold on

plot(tarr,cp1arr,'g',LineWidth=2)

legend('f_1(T)','f_2(T)')

Answer 1

[VBBV]

format long
k=10;
Es=0.05;
t=273;
TempK=0:1500;
cp=0:0.2:1.2;
cp=0.99403+((1.671e-4)*(TempK))+((9.7215e-8)*(TempK.^2))-((9.5838e-11)*(TempK.^3))+((1.9520e-14)*(TempK.^4));
plot(TempK,cp,'r',LineWidth=1.5);
title('Specific Heat Capacity, C_p (kJ/kg\cdotK) VS Temperature (K)');
xlabel('Temperature, T (K)');
ylabel('Specific Heat Capacity, C_p (kJ/kg\cdotK)');
grid on
tarr=zeros(1,k);
cp1arr=zeros(1,k);
for i=1:k
    T=((((0.20597./t.^4)+((9.5838e-11)./(t))-((9.7215e-8)./(t.^2))-((1.671e-4)./t.^3)))/(1.9520e-14)).^(1/4);
    tnew=T;
    Ea=abs((((tnew-t)/tnew)*100));
   cp1=0.99403+((1.671e-4)*(T))+((9.7215e-8)*(T.^2))-((9.5838e-11)*(T.^3))+((1.9520e-14)*(T.^4));
  
     
     
%     if Ea<Es
%         disp('---------------------------------------------');
%         disp(['Iteration,k = ',num2str(i,k)]);
%         disp(['Cp = ',num2str(cp1),' kJ/kgK'])
%         disp(['t = ',num2str(t)]);
%         disp(['t new = ',num2str(tnew)]);
%         disp(['Ea = ',num2str(Ea), ' %']);
%         break;
%     end
    t=tnew;
    
%     disp('---------------------------------------------');
%     disp(['Iteration,k = ',num2str(i,k)]);
%     disp(['Cp = ',num2str(cp1),' kJ/kgK'])
%     disp(['t = ',num2str(t)]);
%     disp(['t new = ',num2str(tnew)]);
%     disp(['Ea = ',num2str(Ea),' %']);
    
    cp1arr(1,i)=cp1;
    tarr(1,i)=t;  
   
   
end
cp1arr, cp  % check using format long option

cp1arr = 1×10

0.995061235421654 1.049012085587812 0.994985384391610 1.053692028545137 0.994908170191703 1.059311490328687 0.994829372351758 1.066190044870491 0.994748692982960 1.074811200746848

cp = 1×1501

0.994030000000000 0.994197197119182 0.994364588093608 0.994532172348955 0.994699949311365 0.994867918407450 0.995036079064290 0.995204430709434 0.995372972770898 0.995541704677169 0.995710625857200 0.995879735740414 0.996049033756703 0.996218519336425 0.996388191910408 0.996558050909950 0.996728095766815 0.996898325913236 0.997068740781916 0.997239339806024 0.997410122419200 0.997581088055551 0.997752236149653 0.997923566136550 0.998095077451756 0.998266769531250 0.998438641811484 0.998610693729374 0.998782924722309 0.998955334228143

    hold on
    plot(tarr,cp1arr,'k-',LineWidth=1)
    legend('f_1(T)','f_2(T)')

[VBBV]

While obtaining new temperature profile, you may need the below changes

T=((((0.20597./t.^4)+((9.5838e-11)./(t))-((9.7215e-8)./(t.^2))-((1.671e-4)./t.^3)))/(1.9520e-14)).^(1/4);

[Amirullah Hamzah]

Thank you so much