0

I have this code which solves the heat equation with non homogeneous boundary conditons: $$u_t=ku_{xx}$$ where k is the diffusion coefficient but now I want to add a constant which makes the equation non homogeneous, like this: $$u_t=ku_{xx}+h$$

function [u,er]=ftcs(t0,tf,nt,a,b,nx,ci,cca,ccb,dif,cada,sol)
% explicit method

%u solution
%er vector error
%t0 y tf time limits
%nt 
%nx 
%a b x value
%ci initial condition
%cca y ccb boundary conditions
%dif diffusion coefficient
%cada time lapse
%sol solution if it's known

x=linspace(a,b,nx);  x=x';  dx=x(2)-x(1); %h
t=linspace(t0,tf,nt); t=t';  dt=t(2)-t(1); %k
r=dt/dx^2*dif;  
if r>.5
    clc
    disp('stability not satisfied')
    pause
end
di=1-2*r;
cca=inline(cca,'t');
ccb=inline(ccb,'t');
A=diag(di*ones(nx-2,1))+diag(r*ones(nx-3,1),1)+diag(r*ones(nx-3,1),-1);
A=[[r;zeros(nx-3,1)] A [zeros(nx-3,1);r]];
ci=vectorize(inline(ci,'x')); 
u=ci(x);      % vectorwith a solution for t=0, initial condition
gra=plot(x,u);
axis([a b min(u) max(u)])
pause
z=u;
u(1)=cca(t(1)); 
u(end)=ccb(t(1)); 
for i=1:nt-1;  
    u=A*u;
    u=[cca(t(i+1)); u ; ccb(t(i+1))];
    if mod(i,cada)==0
        set(gra,'ydata',u');
        pause(.1)
        z(:,end+1)=u;
    end
end
if nargin==12 & nargout==2
    sol=vectorize(inline(sol,'x','t'));
    er=u-sol(x,tf);
end

pause
figure
surfc(z);shading interp; colormap(hot);set(gca,'ydir','reverse');
rotate3d
Dev-iL
  • 23,742
  • 7
  • 57
  • 99
  • So? What's your question? Where exactly are you stuck? – tryman Apr 30 '19 at 16:09
  • @tryman I don't knoe how to add that termi which is just a constant and doesn't depend on t since I don't know how to use matlab. I just want it to get a plot of one case. The only solution I see is to use this code for the homogeneous and plot it and then create another plot just for the term h and plot it. – Phys pilot Apr 30 '19 at 17:02

0 Answers0