MATLAB Answers

PLOT IS NOT SHOWING ANYTHING

36 views (last 30 days)
S H
S H on 25 Mar 2020 at 20:05
Commented: darova on 28 Mar 2020 at 9:14
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 Comments

Show 1 older comment
S H
S H on 27 Mar 2020 at 4:45
I think the mistake is
remove Y=length(r)
and plot(t,Y) instead of plot(r, Y)
but I would like to plot the concentration of different species vs. time in separate graphs and the concentration of different species vs radius in in separate graphs.
I have no idea how to do it? Do you have any idea?
Thanks!
Rik
Rik on 27 Mar 2020 at 7:15
What exactly do you want to use as x and what as y?
S H
S H on 28 Mar 2020 at 2:38
Hello Rik,
The picture above is what I would like to plot.

Sign in to comment.

Answers (1)

darova
darova on 27 Mar 2020 at 11:07
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 Comments

S H
S H on 28 Mar 2020 at 2:50
Hello darova,
I used the method you suggested. "plot(T, Y(:,1),'r')" is the conc. vs. time for one species at a fixed radius. I was wondering if we can plot the conc. vs. time for one species at different radius in the same graph. Also, the conc. vs. radius for one species at different times in the same graph. The picture below is what exactly I would like to plot. Thank you very much for your help!
Rik
Rik on 28 Mar 2020 at 7:11
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 on 28 Mar 2020 at 9:14
BUt your radius is not changing. I don't understand this graph: time and conc are changing (r variable is fixed in your code)

Sign in to comment.

Sign in to answer this question.

Tags

Products