graph for FTCS shown only one point

1 view (last 30 days)
nur afiqah
nur afiqah on 6 Jan 2024
Answered: UDAYA PEDDIRAJU on 20 Jan 2024
% Parameters for v equation
D2 = 8e-10; % Diffusion coefficient for v
L = 1; % Length of the spatial domain for v
T = 1; % Total simulation time for v
V0 = 6938; % Boundary condition at x=0 for v
dx = 0.05; % Spatial step size for v
dt = 0.01; % Temporal step size for v
% Parameters for u equation
D1 = 0.1; % Diffusion coefficient for u
chi = 1; % Coupling coefficient
u_initial = 1; % Initial condition for u
% Discretization for v equation
x_values = 0:dx:L; % Spatial grid for v
t_values = 0:dt:T; % Temporal grid for v
% Stability criterion (FTCS stability condition) for v
r_v = D2 * dt / (dx^2);
% Initial condition for v
v_values = ones(length(x_values), length(t_values)); % Matrix to store v values
v_values(:, 1) = ones(size(x_values)); % v(x, 0) = 1 for all x
% Boundary conditions for v
v_values(1, :) = V0; % v(0, t) = C0 for all t
v_values(end, :) = 0; % v(1, t) = 0 for all t
% FTCS scheme for v equation
for n = 1:length(t_values)-1
for i = 2:length(x_values)-1
v_values(i, n+1) = v_values(i, n) + r_v * (v_values(i+1, n) - 2*v_values(i, n) + v_values(i-1, n));
end
end
% Discretization for u equation
u_values = ones(length(x_values), length(t_values)) * u_initial; % Matrix to store u values
u_values(:, 1) = ones(size(x_values)) * u_initial; % Initial condition for u
% FTCS scheme for u equation
for n = 1:length(t_values)-1
for i = 2:length(x_values)-1
u_values(i, n+1) = u_values(i, n) + dt * (D1 * (u_values(i+1, n) - 2*u_values(i, n) + u_values(i-1, n))/dx^2 - chi * u_values(i, n) * (v_values(i+1, n) - v_values(i-1, n))/(2*dx));
end
end
% Plotting Numerical Solution for v
figure;
surf(t_values, x_values, v_values, 'EdgeColor', 'none');
title('Numerical Solution for v (FTCS)');
xlabel('Time');
ylabel('x');
zlabel('v(x, t)');
grid on;
% Plotting Numerical Solution for u
figure;
surf(t_values, x_values, u_values, 'EdgeColor', 'none');
title('Numerical Solution for u (FTCS)');
xlabel('Time');
ylabel('x');
zlabel('u(x, t)');
grid on;
% ... (Previous code remains unchanged)
% Time instances to plot
selected_times = [0.2, 0.4, 0.6, 0.8];
% Find the corresponding indices for the selected times
selected_indices = round(selected_times / dt) + 1;
% Plotting 2D graphs for v and u at selected times
figure;
% Plotting Numerical Solution for v at selected times
subplot(2, 1, 1);
hold on;
for i = 1:length(selected_indices)
plot(x_values, v_values(:, selected_indices(i)), 'LineWidth', 2, 'DisplayName', ['t = ', num2str(selected_times(i))]);
end
title('Numerical Solution for v at Selected Times');
xlabel('x');
ylabel('v(x, t)');
legend('Location', 'Best');
grid on;
% Plotting Numerical Solution for u at selected times
subplot(2, 1, 2);
hold on;
for i = 1:length(selected_indices)
plot(x_values, u_values(:, selected_indices(i)), 'LineWidth', 2, 'DisplayName', ['t = ', num2str(selected_times(i))]);
end
title('Numerical Solution for u at Selected Times');
xlabel('x');
ylabel('u(x, t)');
legend('Location', 'Best');
grid on;

Answers (1)

UDAYA PEDDIRAJU
UDAYA PEDDIRAJU on 20 Jan 2024
Hi Nur,
By carefully observing the code, it appears that the calculation of indices and plotting logic is correct. If you are still seeing only a single point, it's possible that the data in “v_values” and “u_values” isn't being updated as expected due to an issue within the FTCS scheme implementation.
To avoid the single point issue:
  1. Double-check the FTCS implementation for any logical errors that could cause the data to not update correctly over time.
  2. Verify that the boundary and initial conditions are applied correctly.
  3. Ensure that the stability condition for the FTCS scheme is satisfied to prevent numerical instability. I hope this gives you an idea on how to rectify the FTCS dot issue.
You can refer to the following documentation to read about how FTCS heat equation is implemented: https://www.mathworks.com/matlabcentral/fileexchange/54617-ftcs-heat-equation.

Categories

Find more on Programming in Help Center and File Exchange

Products


Release

R2023b

Community Treasure Hunt

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

Start Hunting!