% Define constants and parameters
kp0 = 140270;
Eap = 43900;
kd0 = 0.355;
Ead = 17400;
T = [333 334 335 336 337 338 339 340]; % Temperature range
m = 2; % Symmetry parameter for pdepe
t = linspace(0, 120, 100); % Time range
x = linspace(0, 1, 100); % Spatial domain
u0 = 5.3347; % Initial concentration
% Loop over temperatures
for i = 1:length(T)
kp = kp0 * exp(-Eap / (8.314 * T(i)));
kd = kd0 * exp(-Ead / (8.314 * T(i)));
% Define the PDE function inline
eqn11 = @(x, t, u, DuDx) deal(1, (1.34 * 10^-6) * exp(-16000 / (8.314 * T(i))) * DuDx, 0);
% Define the initial condition function inline
initial11 = @(x) u0;
% Define the boundary condition function inline
bc11 = @(xl, ul, xr, ur, t) deal(0, 1, 0, 1); % Modify as per your BCs
% Solve the PDE
sol = pdepe(m, eqn11, initial11, bc11, x, t);
% Plotting results for a specific time range (e.g., t=10 to t=20)
% Find indices of t that fall within the desired time range
t_indices = find(t >= 10 & t <= 20);
figure; % Open a new figure for each temperature
hold on;
for j = t_indices
plot(x, sol(j, :)); % Plot x vs. concentration at this time
end
hold off;
title(['Temperature = ', num2str(T(i)), ' K']);
xlabel('Radius');
ylabel('Concentration');
end
-----------------------------------------------------------------------------------------------------------------------------------------------------
If you find the solution helpful and it resolves your issue, it would be greatly appreciated if you could accept the answer. Also, leaving an upvote and a comment are also wonderful ways to provide feedback.
It's important to note that the advice and code are based on limited information and meant for educational purposes. Users should verify and adapt the code to their specific needs, ensuring compatibility and adherence to ethical standards.
Professional Interests
- Technical Services and Consulting
- Embedded Systems | Firmware Developement | Simulations
- Electrical and Electronics Engineering
Feel free to contact me.
