Parciální diferenciální rovnice

Rovnice vedení tepla - explicitní metoda sítí

Ilustrace závislosti stability explicitní metody na délce kroku.

Uložte si do souboru teplo.m následující program s příkladem a spusťe v Matlabu .

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%  reseni rovnice vedeni tepla explicitni metodou  
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

p = 0.3;           % soucinitel tepelne vodivosti
h = 0.25;          % prostorovy krok
tau = 0.1;         % casovy krok
sig = p*tau/(h*h); % koeficient sigma v explicitni metode
a = 0; b = 1;      % interval reseni
xh = a:h:b;        % rozdeleni intervalu s krokem h
n = length(xh);    % pocet uzlu
f0 = xh(2:n-1);    % prava strana rovnice (zdroj tepla): f0(x,t) = x
U0 = xh.^2;        % pocatecni podminka: u0 = x*x
hold off
plot(xh,U0)        % nakresleni pocatecni podminky
hold on
U1(1) = 0; U1(n) = 1;  % okrajove podminky - konstantni

   % cyklus vypoctu pro jednotlive vrstvy
disp('Pro dalsi iteraci stisknete mezernik ... ');
for i = 1:10
  U1(2:n-1) = sig*(U0(1:n-2)+U0(3:n))+(1-2*sig)*U0(2:n-1)+tau*f0;
  pause
  plot(xh,U1)
  U0 = U1;
end
disp('konec.');

   % zmena casoveho kroku a pokracovani ve vypoctu
disp('A ted pokracovani s delsim krokem: pro dalsi iteraci stisknete mezernik ... ');
tau = 0.3;
sig = p*tau/(h*h);
for i = 1:7
  U1(2:n-1) = sig*(U0(1:n-2)+U0(3:n))+(1-2*sig)*U0(2:n-1)+tau*f0;
  pause;
  plot(xh,U1,'r')
  U0 = U1;
end
disp('konec.');