Для решения данной задачи воспользуемся уравнением теплопроводности:
∂T/∂t = (λ/(ρc)) * ∆²T
где T – температура, t – время, λ – коэффициент теплопроводности, ρ – плотность, c – удельная теплоемкость, ∆² – оператор Лапласа.
Учитывая сферическую симметрию задачи и условие на границе цилиндра, преобразуем уравнение в цилиндрические координаты:
∂T/∂t = (λ/(ρc)) * (1/r) * (∂/∂r) * (r * (∂T/∂r))
Также учтем начальное условие T(r, 0) = T0 и граничное условие T(R, t) = Th.
Для численного решения введем сетку по координатам r и времени t. Разобьем радиус цилиндра на N равных интервалов, и время на M интервалов.
Шаг по времени h_t = t/M
Шаг по координате h_r = R/N
Заменим производные разностными аналогами, используя центральные разности:
(T_i_j+1 – T_i_j)/h_t = (λ/(ρc)) * (1/r_i) * ((r_i+1 * (T_i+1_j – T_i_j))/(h_r * h_r) – (r_i-1 * (T_i-1_j – T_i_j))/(h_r * h_r))
где i – индекс по радиусу r, j – индекс по времени t.
Решим полученные разностные уравнения численно, начиная со времени t = 0. После каждого шага по времени t обновим значения температуры на границе цилиндра T_N_j = Th.
Повторим алгоритм для 5, 10, 30 и 60 секунд, используя полученные значения температуры на предыдущем шаге в качестве начальных условий.
В результате, мы получим температурное поле в цилиндре через указанные промежутки времени.