PLOT IS NOT SHOWING ANYTHING

1 次查看（过去 30 天）
S H
S H 2020-3-25
评论： darova 2020-3-28
Hello everyone,
I have a question regarding nothing showing up in my graph. How can I fix it?
My code is in the following:
clear all
clc
format long
%----parameters----%
global T0 Rg E1 E2 E_2 Eg k10 k20 k2_0 kgo Deo k1 k2 k_2 Ke Kg Cb Ct Ce a
a=1;
T0=800;%k Temperature
Rg=8.314;%J/k*mol;
E1=10^5;%J/mol;
E2=10^5;%J/mol;
E_2=10^2;%J/mol;
Eg=10^4;%J/mol
k10=0.1;
k20=0.1;
k2_0=0.1;
kgo=0.1;
Deo=10^-10;
k1 = k10*exp(-E1/(Rg*T0));
k2 = k20*exp(-E2/(Rg*T0));
k_2 = k2_0*exp(-E_2/(Rg*T0));
Ke=(2.01*10^(-6))*T0^2-(1.43*10^(-3))*T0+0.2544;
Kg=kgo*exp(Eg/(Rg*T0));
Cb =0;
Ct=20;
Ce=(7.67131719*10^(-47))*T0^(14.5144484);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
r=linspace(0,100e-4,5);
t=linspace(0,300,100);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
n = numel(r);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
CC0 = zeros(n,1);
CO0 = ones(n,1);
CD0 = ones(n,1);
CM0 = zeros(n,1);
CA0 = zeros(n,1);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
CO0(1)=15;
CD0(1)=5;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
y0 = [CC0;CO0;CD0;CM0;CA0];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[T, Y] = ode15s(@(t,y)fun(t,y,r,n),t,y0)%T is time
%%%%%%%%PLOT%%%%%%%%%%%%%%%%%%%
Y=length(r);
plot(r,Y)
set(gca,'YLim',[0 1e-2],'xLim',[0 5e-2]);
xlabel('Radius')
ylabel('concentration')
legend('CC','CO','CD','CM','CA')
Function:
function DyDt=fun(t,y,r,n)
global T0 Rg E1 E2 E_2 Eg k10 k20 k2_0 kgo Deo k1 k2 k_2 Ke Kg Cb Ct Ce a
CC = zeros(n,1);
CO = zeros(n,1);
CD = zeros(n,1);
CM = zeros(n,1);
CA = zeros(n,1);
DCCDt = zeros(n,1);
DCODt = zeros(n,1);
DCDDt = zeros(n,1);
DCMDt = zeros(n,1);
DCADt = zeros(n,1);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
DyDt = zeros(5*n,1);
rhalf = zeros(n-1,1);
DCCDr = zeros(n-1,1);
D2CCDr2 = zeros(n-1,1);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
CC = y(1:n);
CO = y(n+1:2*n);
CD = y(2*n+1:3*n);
CM = y(3*n+1:4*n);
CA = y(4*n+1:5*n);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
rhalf(1:n-1)=(r(1:n-1)+r(2:n))/2;
%%%%%%%B.C.Left%%%%%%%%%%%%%%%%%
DCCDr(1)=0;
D2CCDr2(1) = 1.0/(rhalf(1)-r(1))*(CC(2)-CC(1))/(r(2)-r(1));
%%%%%%%B.C.Right%%%%%%%%%%%%%%%%%
DCCDt(n)=Cb
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
DCODt(n) =-k1*(CO(n)^2 -( 1/Ke)*CD(n)^2);
DCDDt(n) =k1*(CO(n)^2 - (1/Ke)*CD(n)^2)-k2*CD(n)+k_2*CA(n)*CM(n)*CC(n);
DCMDt(n) =k2*CD(n)-k_2*CA(n)*CM(n)*CC(n);
DCADt(n) =k2*CD(n)-k_2*CA(n)*CM(n)*CC(n);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for i=2:n-1
DCCDr(i) = ((r(i)-r(i-1))/(r(i+1)-r(i))*(CC(i+1)-CC(i))+(r(i+1)-r(i))/(r(i)-r(i-1))*(r(i)-CC(i-1)))/(r(i+1)-r(i-1));
D2CCDr2(i) = (rhalf(i)*(CC(i+1)-CC(i))/(r(i+1)-r(i))-rhalf(i-1)*(CC(i)-CC(i-1))/(r(i)-r(i-1)))/(rhalf(i)-rhalf(i-1));
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for i=1:n-1
De=Deo*exp(a*CM(i)/Ct);
DCCDt(i) =2*De*DCCDr(i)+De*r(i)*D2CCDr2(i)+k2*CD(i)-k_2*CA(i)*CM(i)*CC(i);
DCODt(i) =-k1*(CO(i)^2 -( 1/Ke)*CD(i)^2);
DCDDt(i) =k1*(CO(i)^2 - (1/Ke)*CD(i)^2)-k2*CD(i)+k_2*CA(i)*CM(i)*CC(i);
DCMDt(i) =k2*CD(i)-k_2*CA(i)*CM(i)*CC(i);
DCADt(i) =k2*CD(i)-k_2*CA(i)*CM(i)*CC(i);
end
DyDt = [DCCDt;DCODt;DCDDt;DCMDt;DCADt];
end
Can anyone help? Thank you very much!
  4 个评论
显示 2更早的评论
Rik
Rik 2020-3-27
What exactly do you want to use as x and what as y?
S H
S H 2020-3-28
Hello Rik,
The picture above is what I would like to plot.

请先登录，再进行评论。

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

采纳的回答

darova
darova 2020-3-27
Use subplot to draw curvees in separate windows
subplot(2,1,1)
plot(T,Y(:,1),'r')
subplot(2,1,2)]
plot(T,Y(:,2),'r')
  3 个评论
显示 1更早的评论
Rik
Rik 2020-3-28
Did you do a simple Google search for how to plot multiple graphs in a single axes? That should turn up the hold function.
darova
darova 2020-3-28
BUt your radius is not changing. I don't understand this graph: time and conc are changing (r variable is fixed in your code)

请先登录，再进行评论。

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Spline Postprocessing 的更多信息

标签

产品

Community Treasure Hunt

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

Start Hunting!

Translated by