codice con ode45
1 visualización (últimos 30 días)
Mostrar comentarios más antiguos
Sto provando a risolvere un problema di carbocementazione usando ODE45, in particolr modo sto usando il metodo delle linee. Lo scopo dell'esercizio è quello di ricavare dopo un determinato istante di tempo T il grado di assorbimento del carbonio. Non riesco proprio a capire come risolvere perchè mi sembra che sia inserito tutto correttamente. qualcuno mi può aiutare?
clear
clc
%cabocemntazione
%defnizione parametri
L=1; Temp=1273; D0=20*10^-6; Q=142000; R=8.314; T=10;
nx=100; nt=100;
%controllo parametri
if any([L Temp T D0 nx-2 nt-2])<=0
error('controlla i parametri')
end
%inizializzazione
D=D0*exp(-Q/(R*Temp)); %diffusività
dx=L/nx; %ho discretizzato dx perchè mi sefve nella formula di p
p=D/dx^2;
x=linspace(0,L,nx+1);%non ho fatto il tempo non mi serve
A=[p*ones(nx,1) -2*p*ones(nx,1) p*ones(nx,1)];
AA=spdiags(A, -1:1, nx,nx);
A(end,end-1)=2*p;
a=sparse(nx,1);
%condizione iniziale
c=GetPhi(nx);
[t,y]=ode45(@system,[0 T], c(2:end),[],A,a,p);
c(2:end)=y(end,:);
c(1)=Getg1;
plot(x,c,'b*');
axis([0 L 0 Getg1]);
xlabel('x'); ylabel('%w di C');
legend('MOL');
grid
pause(0.1)
%---------subfunction----------
function phi=GetPhi(nx)
c0=0.2;
phi=c0*ones(nx+1,1);
end
function g1=Getg1
cs=1;
g1=cs;
end
function g2=Getg2
g2=0;
end
function cdot=system(t,c,A,a,p);
a(1)=p*Getg1;
a(end)=p*Getg2;
cdot=A.*c+a ;
end
clear
clc
%cabocemntazione
%defnizione parametri
L=1; Temp=1273; D0=20*10^-6; Q=142000; R=8.314; T=10;
nx=100; nt=100;
%controllo parametri
if any([L Temp T D0 nx-2 nt-2])<=0
error('controlla i parametri')
end
%inizializzazione
D=D0*exp(-Q/(R*Temp)); %diffusività
dx=L/nx; %ho discretizzato dx perchè mi sefve nella formula di p
p=D/dx^2;
x=linspace(0,L,nx+1);%non ho fatto il tempo non mi serve
A=[p*ones(nx,1) -2*p*ones(nx,1) p*ones(nx,1)];
AA=spdiags(A, -1:1, nx,nx);
A(end,end-1)=2*p;
a=sparse(nx,1);
%condizione iniziale
c=GetPhi(nx);
[t,y]=ode45(@system,[0 T], c(2:end),[],A,a,p);
c(2:end)=y(end,:);
c(1)=Getg1;
plot(x,c,'b*');
axis([0 L 0 Getg1]);
xlabel('x'); ylabel('%w di C');
legend('MOL');
grid
pause(0.1)
%---------subfunction----------
function phi=GetPhi(nx)
c0=0.2;
phi=c0*ones(nx+1,1);
end
function g1=Getg1
cs=1;
g1=cs;
end
function g2=Getg2
g2=0;
end
function cdot=system(t,c,A,a,p);
a(1)=p*Getg1;
a(end)=p*Getg2;
cdot=A.*c+[p*Getg1; a; p*Getg2] ;
end![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/877895/image.png)
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/877895/image.png)
0 comentarios
Respuestas (1)
Varun
el 23 de En. de 2024
Hey! I ran the code you shared on my end and was able to reproduce the error. While I am unsure about what exactly you are trying to achieve, I see that the reason behind this error is that the function "system" is expected to run a column vector of dimensions [n x 1]. However, in the case the "cdot" variable is actually of size 100 x 3, which results in this error.
So, refactoring your code to ensure that the "system" method returns a column vector should resolve this error. Hope this helps!
0 comentarios
Ver también
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!