필터 지우기
필터 지우기

Why is my code not working how I want it to?

조회 수: 3 (최근 30일)
Marlon Brutus
Marlon Brutus 2017년 11월 27일
댓글: MHZ 2017년 11월 28일
Good News Everyone!
My matlab code is malfunctioning and I have no idea why... So I am trying to convert this Polymath code into a Matlab code and it is not working out. I attach the Polymath code and what is supposed to come out but I am just not getting it. Any help will be greatly appreciated. :D
This is my code:
function nonisothermal
vspan=[0 1]; %volume of the reactor dm3
x0=[0 1035]; %Conversion and Temperature(in K)
[v x]=ode45(@thermalfun,vspan,x0);
subplot(2,1,1)
plot(v,x(:,1),'m--','linewidth',1)
xlabel('Volume of the Reactor')
title('Conversion')
grid on
subplot(2,1,2)
plot(v,x(:,2),'k--','linewidth',1)
xlabel('Volume of the Reactor')
title('Temperature')
grid on
end
function dXdV=thermalfun(v,x)
Cao=18.8;
Fao=0.0376;
Cpa=163;
delCp=-9;
To=1035;
delH=80770+delCp*(x(2)-298);
ra=-Cao*3.58*exp(34222*(1/To-1/x(2)))*(1-x(1))*(To/x(2))/(1+x(1));
Ua=165000;
Ta=1150;
dXdV=[-ra/Fao; (Ua*(Ta-x(2))+ra*delH)/(Fao*(Cpa+x(1)*delCp))]
end
  댓글 수: 1
Star Strider
Star Strider 2017년 11월 28일
I do not see any significant problems with your code, with respect to the published equations. There is an insignificant problem in that ‘Ua’ should be 16500, and that the second plot is ‘-ra’ as a function of ‘v’, so you need to calculate ‘ra’ outside the ode45 call with ‘x(:,1)’, ‘x(:,2)’ and the relevant constants.
Check to see the ODE solver the authors used. The MATLAB numerical solvers are the best I’ve encountered, however some papers use much less sophisticated solvers, leading to results that are difficult to reproduce with the MATLAB solvers. (I invariably trust the MATLAB results.) Also experiment with a much shorter upper limit for ‘vspan’.

댓글을 달려면 로그인하십시오.

답변 (1개)

MHZ
MHZ 2017년 11월 27일
편집: per isakson 2017년 11월 28일
I tried out your code and ran it fine. First make sure to save 2 files in a location in matlab path. Just creat a folder to save in and then go to Matlab home> set path> add folder then hit save and close. Although you do not need to save the first one as a function and can do a normal m-script, save the following as nonisothermal
function nonisothermal
vspan=[0 1]; %volume of the reactor dm3
x0=[0 1035]; %Conversion and Temperature(in K)
[v x]=ode45(@thermalfun,vspan,x0);
subplot(2,1,1)
plot(v,x(:,1),'m--','linewidth',1)
xlabel('Volume of the Reactor')
title('Conversion')
grid on
subplot(2,1,2)
plot(v,x(:,2),'k--','linewidth',1)
xlabel('Volume of the Reactor')
title('Temperature')
grid on
end
Now save the following as thermalfun
function dXdV=thermalfun(v,x)
Cao=18.8;
Fao=0.0376;
Cpa=163;
delCp=-9;
To=1035;
delH=80770+delCp*(x(2)-298);
ra=-Cao*3.58*exp(34222*(1/To-1/x(2)))*(1-x(1))*(To/x(2))/(1+x(1));
Ua=165000;
Ta=1150;
dXdV=[-ra/Fao; (Ua*(Ta-x(2))+ra*delH)/(Fao*(Cpa+x(1)*delCp))]
end
When done, go to the command window and type nonisothermal and hit enter.
  댓글 수: 2
Marlon Brutus
Marlon Brutus 2017년 11월 28일
I don't think you understand... The code is supposed to come with results that are pretty similar to what is attached above.
MHZ
MHZ 2017년 11월 28일
Your problem is in the implementation of the integration equation. You are using x(2) instead of T. According to Matlab, x(2) will be 1035 which is your final and not the value of T. You need to separate dT/dV and actually use T.

댓글을 달려면 로그인하십시오.

카테고리

Help CenterFile Exchange에서 Programming에 대해 자세히 알아보기

태그

Community Treasure Hunt

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

Start Hunting!

Translated by