% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% % heatSln % ____________________________________________________________________________ % 18.385 (Nonlinear Dynamical Systems). MIT, Fall 1999. R. R. Rosales. % % Demo illustrating the solution of the heat equation T_t = T_xx with periodic % boundary conditions on 0 <= x <= 2*pi and initial conditions u(x, 0) = F(x) % as given by the same options in the fourierSC script (see help fourierSC for % more information). A spectral method is used to solve the equation. % % At the end of the run: % x (row array) contains the points in space. % fct (""""""""") contains the initial values. % tt (column array) contains the computed times. % uu contains the solution, each row corresponding to an entry in tt. % dt contains the separation between computed times. % t_f contains the final time. .... Thus: tt = (0:dt:t_f)' % %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% %========== Vector sizes and useful text variables: ========================== M = 512; N = 2*M; ind = [0:M,(1-M):(-1)]; FN = 'FontName'; HNB = 'HelveticaNarrowBold'; FS = 'FontSize'; LW = 'LineWidth'; MS = 'MarkerSize'; Space = ' '; % %========== Enter inputs needed and print on-screen info: ==================== % fprintf('\n Solution to the heat equation T_t = T_xx.') fprintf('\n The solution will be plotted for t=0:dt:t_f \n') % dt = input([Space, 'Enter dt ......................... dt = ']); z = exp(-dt*ind.^2); t_f = input([Space, 'Enter t_f ....................... t_f = ']); % %========== Construct array of x's =========================================== % dx = 2*pi/N; x = 0:dx:2*pi-dx; % %========== Construct function =============================================== % fprintf('\n NOW CHOOSE AN INITIAL CONDITION. \n') FSoption; % %========= Plotting Ranges =================================================== % Max = max(fct); Min = min(fct); range = Max - Min; Max = Max + 0.1*range; Min = Min - 0.1*range; % %========= READY TO CALCULATE AND PLOT ======================================= % fprintf('\n PRESS ANY KEY TO MOVE FORWARD IN TIME.\n') u = fct; uu = u; t=0; tt = 0; figure whitebg('k') plot(x, u, '-r', LW, 3) title('Heat equation: initial condition.', FN, HNB, FS, 22) xlabel('x', FN, HNB, FS, 18) ylabel('Solution', FN, HNB, FS, 18) grid on axis([0 2*pi Min Max]) pause % while t <= t_f - 0.1*dt u = real(ifft(fft(u).*z)); uu = [uu; u]; t = t+dt; tt = [tt; t]; TEXT = num2str(t); plot(x, u, '-r', LW, 3) title(['Heat equation solution, time t = ', TEXT], FN, HNB, FS, 22) xlabel('x', FN, HNB, FS, 18) ylabel('Solution', FN, HNB, FS, 18) grid on axis([0 2*pi Min Max]) pause end figure whitebg('k') mesh(x, tt, uu) shading interp colormap(hsv) axis([0 2*pi 0 t_f Min Max]) grid on title('Solution of the heat equation.', FN, HNB, FS, 22) xlabel('Space', FN, HNB, FS, 18) ylabel('Time', FN, HNB, FS, 18) zlabel('Solution', FN, HNB, FS, 18) fprintf('\n DONE. \n') % clear MM OPTION dyt smooth z t MS clear Space FN Max TEXT h FS Min clear ans h_fig u HNB N ind LW clear NN dx jn yb M OPT dyb clear range yt BGC CLR % % EOF