Matlab coding issue with solving PDEs

In summary: if( rem(istep,plot_step) < 1 )ttplot(:,iplot) = tt(:);ttplot(iplot) = istep*tau;ttplot(:,iplot) = istep*tau;endend
  • #1
jkthejetplane
29
4
Homework Statement
I was tasked with putting a few methods for solving PDEs into matlab. I started with the given code that solved using FTCS method. Below is was i have modified it to. I understand i need to have it run FTCS for the first step. I am not sure if the way i have written the equations into the code makes sense but i have gotten the Richardson method to produce results. There is an error with the Dufort method though
"Array indices must be positive integers or
logical values.

Error in problem3 (line 45)
tt(2:(N-1)) = tt(2:(N-1)) +
2*coeff(tt(3:N) + tt(1:(N-2)) -
tt(2:(N-1)))/(1+ 2*coeff);
"
Relevant Equations
FTCS equation is the same form is copied below from my text
1605925246416.png


FTCS scheme
1605925514213.png

My Matlab code:
%Problem 3
%Solve diffusion problem using Richardson scheme or DuFort-Frankel scheme
clear all;

scheme = menu('Choose method of solving diffusion equation:', 'Richardson', 'DuFort-Frankel');
tau = input('Enter time step: ');
N = input('Enter the number of grid points: ');
L = 1.; % The system extends from x=-L/2 to x=L/2
h = L/(N-1); % Grid size
k = 1.; % Diffusion coefficient
coeff = k*tau/h^2;

if( coeff < 0.5 )
disp('Solution is expected to be stable');
else
disp('WARNING: Solution is expected to be unstable');
end
%* Set initial and boundary conditions.
tt = zeros(N,1); % Initialize temperature to zero at all points
tt(round(N/2)) = 1/h; % Initial cond. is delta function in center
% The boundary conditions are tt(1) = tt(N) = 0
%* Set up loop and plot variables.
xplot = (0:N-1)*h - L/2; % Record the x scale for plots
iplot = 1; % Counter used to count plots
nstep = 300; % Maximum number of iterations
nplots = 50; % Number of snapshots (plots) to take
plot_step = nstep/nplots; % Number of time steps between plots

%Run First steps using FTCS
tt(2:(N-1)) = tt(2:(N-1)) + coeff*(tt(3:N) + tt(1:(N-2)) - 2*tt(2:(N-1)))
%Run selected scheme
if (scheme==1) %richardson
for istep=2:nstep
tt(2:(N-1)) = tt(2:(N-1)) + 2*coeff*(tt(3:N) + tt(1:(N-2)) - 2*tt(2:(N-1)));

if( rem(istep,plot_step) < 1 ) % Every plot_step steps
ttplot(:,iplot) = tt(:); % record tt(i) for plotting
tplot(iplot) = istep*tau; % Record time for plots
iplot = iplot+1;
end
end
else (scheme==2)%dufort
for istep=2:nstep
tt(2:(N-1)) = tt(2:(N-1)) + 2*coeff(tt(3:N) + tt(1:(N-2)) - tt(2:(N-1)))/(1+ 2*coeff);

if( rem(istep,plot_step) < 1 ) % Every plot_step steps
ttplot(:,iplot) = tt(:); % record tt(i) for plotting
tplot(iplot) = istep*tau; % Record time for plots
iplot = iplot+1;
end
end
end

%* Plot temperature versus x and t as wire-mesh and contour plots.
figure(1); clf;
mesh(tplot,xplot,ttplot); % Wire-mesh surface plot
xlabel('Time'); ylabel('x'); zlabel('T(x,t)');
title('Diffusion of a delta spike');
pause(1);
figure(2); clf;
contourLevels = 0:0.5:10; contourLabels = 0:5;
cs = contour(tplot,xplot,ttplot,contourLevels); % Contour plot
clabel(cs,contourLabels); % Add labels to selected contour levels
xlabel('Time'); ylabel('x');
title('Temperature contour plot');
 
Physics news on Phys.org
  • #2
"Array indices must be positive integers...

Error in problem3 (line 45)
tt(2:(N-1)) = tt(2:(N-1)) +
2*coeff(tt(3:N) + tt(1:(N-2)) -
tt(2:(N-1)))/(1+ 2*coeff);
"for istep=2:nstep
tt(2:(N-1)) = tt(2:(N-1)) + 2*coeff*(tt(3:N) + tt(1:(N-2)) - 2*tt(2:(N-1)));
 

FAQ: Matlab coding issue with solving PDEs

1. What is Matlab coding and how is it used in solving PDEs?

Matlab coding is a programming language that is commonly used in scientific and engineering applications. It is particularly useful in solving partial differential equations (PDEs) due to its powerful numerical computation capabilities.

2. What are the common issues faced when using Matlab for solving PDEs?

Some common issues that may arise when using Matlab for solving PDEs include incorrect input parameters, convergence issues, and errors in the coding logic. It is important to carefully check and debug the code to ensure accurate results.

3. How can I improve the efficiency of my Matlab code for solving PDEs?

There are several ways to improve the efficiency of Matlab code for solving PDEs. These include vectorizing the code, preallocating arrays, and using built-in functions instead of loops. It is also helpful to optimize the code for the specific problem at hand.

4. Can I use Matlab to solve any type of PDE?

Yes, Matlab can be used to solve a wide range of PDEs, including linear and nonlinear equations, initial value problems, and boundary value problems. However, the complexity of the problem and the available computational resources may impact the feasibility of using Matlab for certain cases.

5. Are there any resources available to help with coding issues in Matlab for solving PDEs?

Yes, there are various online resources, tutorials, and forums available to assist with coding issues in Matlab for solving PDEs. The official Matlab documentation and MathWorks community are great places to start. Additionally, seeking help from experienced Matlab users or attending workshops can also be beneficial.

Back
Top