I do not know what's wrong with my code

1 visualización (últimos 30 días)
Nathan Bueno
Nathan Bueno el 12 de Mzo. de 2019
Respondida: Steven Lord el 12 de Mzo. de 2019
%1-D Heat equation
%example 1 at page 782
%dT/dt=c.d^2T/dx^2
%T(x,t)=temperature along the rod
%by finite difference method
%(T(x,t+dt)-T(x,t))/dt = (T(x+dx,t)-2T(x,t)+T(x-dx,t))/dx^2
%solve for T(x,t+dt) by iteration
%heat constant
close all;
%clear all;
c=1;
L=1; %length of the rod
N=5; %# of elements
%dicretize Xspace
x_vec=0:0.2:1;
dx=x_vec(2)-x_vec(1);
%discretize time
dt=0.5/50;
t_vec=0:dt:0.5;
%temperature matrix
T_mat=zeros(length(x_vec),length(t_vec));
T_mat1=zeros(length(x_vec),length(t_vec));
%boundary conditions
T_mat(1,:)=0;
T_mat(end,:)=0;
%initial conditions
T_mat(:,1)= 1-cos(pi*x_vec(idx))(400/npi)*sin*(npi.*x_vec(idx));
%T_mat(:,1)= sin*((pi.*x_vec);
%[tt,xx]=meshgrid(t_vec,x_vec);
subplot(2,1,1);
mesh(xx,tt,T_mat);
lamda=(c*dt/(dx^2));
for tdx=1:length(t_vec)-1
for idx=1:length(x_vec)
if idx==1
T_mat(idx,tdx+1)=T_mat(idx+1,tdx);
T_mat1(idx,tdx)= 1-cos(pi*x_vec(idx))(400/npi)sin(npix_vec(idx))*exp(-1*(pi^2)*t_vec(tdx)))); %modified code, 2nd line of 'if'
elseif idx==length(x_vec)
T_mat(idx,tdx+1)=T_mat(idx-1,tdx);
T_mat1(idx,tdx)= 1-cos(pi*x_vec(idx))(400/npi)sin*(npi*x_vec(idx))*exp(-1*(pi^2)*t_vec(tdx)); %modified code, 2nd line of 'elseif
else
T_mat(idx,tdx+1)=T_mat(idx,tdx)+lamda*((T_mat(idx+1,tdx)-2*T_mat(idx,tdx)+T_mat(idx-1,tdx)));
T_mat1(idx,tdx)=1-cos(pi*x_vec(idx))(400/npi)*sin(npix_vec(idx))*exp(-1*(pi^2)*t_vec(tdx)); %modified code, 2nd line of 'elseif)
end
end
end
%plot
subplot(2,1,1)
[tt,xx]=meshgrid(t_vec,x_vec);
mesh(xx,tt,T_mat);
title('Finite')
%plot
subplot(2,1,2)
[tt,xx]=meshgrid(t_vec, x_vec);
mesh(xx,tt,T_mat1);
title('Analytical');
figure
plot(x_vec,T_mat(:,1),x_vec,T_mat(:,11),x_vec,T_mat(:,21),x_vec,T_mat(:,31),x_vec,T_mat(:,41),x_vec,T_mat(:,51));
xlabel('Rod length (m)'); ylabel('Temperature (C)');
legend('Initially','At t=0.1sec','At t=0.2 secs','At t=0.3 secs',' At t=0.4 secs',' At t=0.5 secs');
figure
plot(x_vec,T_mat1(:,1),x_vec,T_mat1(:,11),x_vec,T_mat1(:,21),x_vec,T_mat1(:,31),x_vec,T_mat1(:,41),x_vec,T_mat1(:,51));
xlabel('Rod length (m)'); ylabel('Temperature (C)');
legend('Initially','At t=0.1sec','At t=0.2 secs','At t=0.3 secs',' At t=0.4 secs',' At t=0.5 secs');
title('Analytical')
xlabel('Rod Length |m'); ylabel('Temperature |C|');
legend('Initially','At t=0.1 sec','At t=0.0 sec','At t=0.3 sec','At t=0.4 sec','At t=0.5 sec');
  1 comentario
Raghunandan V
Raghunandan V el 12 de Mzo. de 2019
https://in.mathworks.com/matlabcentral/answers/6200-tutorial-how-to-ask-a-question-on-answers-and-get-a-fast-answer

Iniciar sesión para comentar.

Respuestas (1)

Steven Lord
Steven Lord el 12 de Mzo. de 2019
If you receive an error when running your code, post the full text (everything displayed in red) of the error message you receive, and indicate (if it's not obvious) the line in your code on which the error occurs.
If you don't receive an error but receive a different result than you expected, it's likely going to be difficult for us to help you since there aren't a lot of comments in the code to help us understand it. This looks like it may be the solution to a homework assignment (the "%example 1 at page 782" suggests that to me) and if that's true you may want to ask your professor and/or teaching assistant to help you. If it's not a homework assignment, or if your professor and/or TA are unavailable, step through your code using the debugging tools and identify where the code gives a different result than you expected then figure out why.

Categorías

Más información sobre Programming en Help Center y File Exchange.

Etiquetas

Productos


Versión

R2018b

Community Treasure Hunt

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

Start Hunting!

Translated by