The follow code was create to calculate the time to increase the temperature, but when I run it the result is infinite. Anyone can help me?

Question

Altibano Ortenzi 2024-6-17

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2129191-the-follow-code-was-create-to-calculate-the-time-to-increase-the-temperature-but-when-i-run-it-the

评论： Altibano Ortenzi 2024-6-18

Below are the tow functions used to calculate the time to increase the temp. according to the inputs: L1, L2, L3 and heat+power. The initial values were: 0.05m; 0.1 m; 0.05m and 500W respectively. But with any other values the result returns the same: inifinite seconds! I need the results in terms of plotting.

function time_to_increase_temperature = calculate_time_to_increase_temperature( ...
    rho_concrete, rho_air, k_concrete, k_air, ~, cp_concrete, ...
    cp_air, T_initial, T_surface, L1, L2, L3, volume_P1_P3, ~, ...
    Nx, Nt, delta_T, heat_power)
    % Calculate total thickness
    total_thickness = L1 + L2 + L3;
    % Thermal diffusivities
    alpha_concrete = k_concrete / (rho_concrete * cp_concrete);  % m²/s
    alpha_air = k_air / (rho_air * cp_air);  % m²/s
    % Discretization parameters
    dx = total_thickness / Nx;  % Spatial step size
    dt = 6 / Nt;  % Time step size
    % Preallocate temperature arrays
    T = T_initial * ones(Nx, 1);
    T_new = T;
    % Define material properties at each point
    len_P1 = round(Nx * L1 / total_thickness);  % Length of concrete layer P1
    len_P2 = round(Nx * L2 / total_thickness);  % Length of air layer P2
    len_P3 = Nx - len_P1 - len_P2;  % Length of concrete layer P3
    % Ensure len_P3 is correctly adjusted if rounding caused any issues
    if len_P1 + len_P2 + len_P3 ~= Nx
        len_P3 = Nx - len_P1 - len_P2;
    end
    % Preallocate material array
    material = zeros(Nx, 1);
    material(1:len_P1) = alpha_concrete;
    material(len_P1+1:len_P1+len_P2) = alpha_air;
    material(len_P1+len_P2+1:Nx) = alpha_concrete;
    time_to_increase_temperature = inf;  % Initialize to store time to increase temperature
    % Calculate heat input per spatial point
    heat_input_per_point = heat_power * dt / (rho_concrete * cp_concrete * volume_P1_P3 / Nx);
    % Time integration using Finite Difference Method (FDM)
    time = 0;  % Initialize time
    for n = 1:Nt
        for i = 2:Nx-1
            T_new(i) = T(i) + material(i) * dt / dx^2 * (T(i+1) - 2*T(i) + T(i-1));
        end
        % Boundary conditions
        T_new(1) = T(1) + alpha_concrete * dt / dx^2 * (T(2) - T(1)) + heat_input_per_point * (T_surface - T(1));
        T_new(Nx) = T(Nx-1);
        % Update temperature
        T = T_new;
        time = time + dt;
        % Check if the temperature at the right side of P3 reaches T_initial + delta_T
        if T(Nx) >= T_initial + delta_T
            time_to_increase_temperature = time;
            break;
        end
    end
    % Return the time to increase temperature by 1 K at the right side of P3
end

Function main_process:

function main_process()

% Parameters with fixed values

rho_concrete = 2500; % kg/m³

rho_air = 1293; % kg/m³

k_concrete = 1.35; % W/(m·K)

k_air = 0.024; % W/(m·K) assuming typical value for air

cp_concrete = 0.92 * 1000; % J/(kg·K)

cp_air = 1.005 * 1000; % J/(kg·K)

T_initial = 300; % Initial temperature in K

T_surface = 411; % Surface temperature in K

volume_P1_P3 = 0.00488; % Volume of each concrete plate in m³

Nx = 100; % Number of spatial points

Nt = 1000; % Number of time steps

delta_T = 1; % Temperature increase in K

% Prompt the user for specific input parameters

% L1 = input('Enter the thickness of the first concrete layer (m): ');

% L2 = input('Enter the thickness of the air layer (m): ');

% L3 = input('Enter the thickness of the second concrete layer (m): ');

% heat_power = input('Enter the heat power (W): ');

L1 = 1;

L2 = 1;

L3 = 1;

heat_power = 1;

% Call the function

time_to_increase_temperature = calculate_time_to_increase_temperature( ...

rho_concrete, rho_air, k_concrete, k_air, [], cp_concrete, ...

cp_air, T_initial, T_surface, L1, L2, L3, volume_P1_P3, [], ...

Nx, Nt, delta_T, heat_power);

% Display the result

fprintf('Time to increase temperature by 1 K at the right side of P3: %.2f seconds\n', time_to_increase_temperature);

% Plot the result

total_thickness = L1 + L2 + L3;

x = linspace(0, total_thickness, Nx); % Spatial points

y = linspace(0, time_to_increase_temperature, Nx); % Temporal points for plotting

plot(x, y, 'DisplayName', 'Time to Increase Temperature by 1 K');

xlabel('Position (m)');

ylabel('Time (s)');

legend('show');

title('Time to Increase Temperature by 1 K at the right side of P3');

end

main_process()

Time to increase temperature by 1 K at the right side of P3: Inf seconds

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

Altibano Ortenzi 2024-6-18

I have changed the code, with your suggestions, but didn't work. I will try other approach, since the temperature did not change with time, or the time remais unchanged.

Thanks.

请先登录，再进行评论。

请先登录，再回答此问题。

Answer 1

Animesh 2024-6-17

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2129191-the-follow-code-was-create-to-calculate-the-time-to-increase-the-temperature-but-when-i-run-it-the#answer_1472541

在 MATLAB Online 中打开

Hi @Altibano Ortenzi,

Based on my understanding, you are attempting to calculate the time required to increase the temperature by 1 K. In reviewing your code, I noticed that the variable time_to_increase_temperature is intended for this purpose. It is initially set to inf (infinity) and is later updated within an if-statement as follows:

% Check if the temperature at the right side of P3 reaches T_initial + delta_T
if T(Nx) >= T_initial + delta_T
    time_to_increase_temperature = time;
    break;
end

Given the inputs: L1(0.05), L2(0.1), L3(0.05), and for heat power(500), it appears that T(Nx) never exceeds the sum of T_initial(300) and delta_T(1).

Consequently, time_to_increase_temperature remains unchanged and retains its initial value of infinity. I suggest reviewing the code to ensure everything is configured correctly and to verify that the conditions for changing time_to_increase_temperature are met as expected.

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

Altibano Ortenzi 2024-6-17

Thanks Animesh. Actually, the temperature over the P3 surface canot be >= T(Nx), since it must be increased by the heat power source. So, I should rewrite the code to avoid this comparison.

请先登录，再进行评论。

Answer 2

surya venu 2024-6-17

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2129191-the-follow-code-was-create-to-calculate-the-time-to-increase-the-temperature-but-when-i-run-it-the#answer_1472546

在 MATLAB Online 中打开

Hi,

To address the infinite time result in your temperature increase calculation, consider these points for correction:

Heat Input Calculation: The formula for "heat_input_per_point" might not be correctly distributing the heat power across the spatial points. Simplify this calculation.
Boundary Conditions: The way heat input is added in the boundary condition might not accurately model the physical process. Review and adjust this to ensure it reflects how heat is actually being applied to the system.

Here are the things you can try:

Revise Heat Input: Directly distribute the heat power across the spatial points without the complex division. Ensure it's added correctly at the boundary or throughout the volume as intended.
Adjust Boundary Condition: Ensure the boundary condition accurately models the heat addition. This could mean directly increasing the temperature based on the heat input or correctly modeling the heat flux.

Here's a simplified approach for the heat input calculation:

% Simplified heat input calculation
heat_input_per_point = heat_power / (rho_concrete * cp_concrete * volume_P1_P3);

And for the boundary condition, directly add the heat input:

T_new(1) = T(1) + (alpha_concrete * dt / dx^2 * (T(2) - T(1))) + heat_input_per_point;

Ensure these changes reflect the actual physical setup and re-evaluate the model.

Hope it helps.

1 个评论
显示 -1更早的评论隐藏 -1更早的评论

Altibano Ortenzi 2024-6-17

Thanks @surya venu ! I'll try to modify the code with the lines you have to modify the calculi, but now I I'm seeing that the point described by @Animesh may be one of the code issues too. Anyway, I'll try your both suggestions.

请先登录，再进行评论。