# Analytical solution of 1-D transient heat transfer

12 ビュー (過去 30 日間)
Nour 2024 年 4 月 12 日

Hello,
I am trying to solve a 1D transient heat equation using the analytical solution for radii from 0 to 5 cm, with a convective bounday conditions as shown in the picture. I have attached my discretization method as well to give a better insight into my problem. I have already coded the solution with out including convenction boundary condition but with a constant T_inf= 298 K.
-- Governing equation used
-- Convection boundary condition to be included
clc
clear all;
%-----------------------------Input----------------------------------------
% Inputs
r0=0.05; % wall thickness(m)
k=0.1; % conductivity(W/m-K)
rho=821.42; % density(kg/m^3)
cp=1399; % specific heat capacity(J/kg-K)
T_ini=353; % initial temperature(K)
T_inf=298; % external temperature(K)
h=5; % heat transfer coefficient(w/m^2-K)
t_f=1000
lampda1=1.694 % Eigen root
C1=1.3787 % constant
R=0.05;
T_i=353;
T_inf=298;
r=linspace(0,0.05,50);
t=linspace(0,100,t_f);
%-----------------------------Main parameters -----------------------------
% Parameters
Bi=(h*r0)/k % Biot number
alpha=k/(rho*cp) % Thermal diffusivity
F0=(alpha*t)/(r0*r0) % Fourier number
r_div=50;
epsilon=1.694;
C_1=1.3787;
T=zeros(1,r_div);
u=zeros(1,r_div);
t=480;
%-----------------------------Temperature calculation ---------------------
for j=1:5
At least one END is missing. The statement beginning here does not have a matching end.
for i=1:50
u(i)=epsilon*r(i) /R;
T(i)=T_inf+(T_i-T_inf)*C_1*exp(-epsilon^2*alpha*t/R^2)*besselj(0,u(i));
end
##### 6 件のコメント4 件の古いコメントを表示4 件の古いコメントを非表示
Torsten 2024 年 4 月 13 日

I don't understand your code from above.
T is a vector - so you only compute the r-profile for one give time t (in your case t = 50000) ?
And if you plot the temperature profile for smaller times (e.g. t = 100) , the temperature at r = R_i does not equal T_w = 298 K for all times t as you claim.
Or is the code you posted an attempt to implement the convective boundary condition ?
Nour 2024 年 4 月 13 日

yes, I am now struggling to add the convection boundary condition to the same code.
Actually I am just a beginner and these are my first trials. I was changing the time every time, save the value, and then plot them as follows.
Thank you a lot for your patience!
T_0=[399.348105708207 400.220997362980 399.823899267308 400.157085055238 399.853772767351 400.138907622977 399.866656443988 400.128678977623 399.875572112861 400.120262237576 399.884092056048 400.111071455955 399.894645024123 400.098108624175 399.911899650103 400.072636544409 399.955270636652 399.981215858787 400.234856351442 397.386534761319]
T_0 = 1x20
399.3481 400.2210 399.8239 400.1571 399.8538 400.1389 399.8667 400.1287 399.8756 400.1203 399.8841 400.1111 399.8946 400.0981 399.9119 400.0726 399.9553 399.9812 400.2349 397.3865
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
T_100=[399.999999230084 400.000000268964 399.999999783514 400.000000191480 399.999999826333 400.000000157780 399.999999858242 400.000000124585 399.999999894301 400.000000083422 399.999999873225 399.999997772070 399.999945739678 399.999098778850 399.989556948428 399.914856521539 399.505742260982 397.923154547764 393.537920248508 384.648819048777]
T_100 = 1x20
400.0000 400.0000 400.0000 400.0000 400.0000 400.0000 400.0000 400.0000 400.0000 400.0000 400.0000 400.0000 399.9999 399.9991 399.9896 399.9149 399.5057 397.9232 393.5379 384.6488
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
T_600=[399.999883833492 399.999829774393 399.999604832226 399.998958507437 399.997222201791 399.992765755409 399.981860358590 399.956501011340 399.900598983991 399.783965784248 399.553876193179 399.124915245478 398.369306722075 397.111634660331 395.132978110469 392.188893420349 388.042627631852 382.509713544762 375.504394130675 367.074828552094]
T_600 = 1x20
399.9999 399.9998 399.9996 399.9990 399.9972 399.9928 399.9819 399.9565 399.9006 399.7840 399.5539 399.1249 398.3693 397.1116 395.1330 392.1889 388.0426 382.5097 375.5044 367.0748
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
T_1800=[399.190501334017 399.148212305607 399.017444481274 398.786457266990 398.435582809434 397.937169878682 397.255683223726 396.348111646569 395.164846259286 393.651172628276 391.749473812938 389.402166358795 386.555293647069 383.162591487097 379.189734947228 374.618391275875 369.449659021043 363.706481627073 357.434690829772 350.702457235095]
T_1800 = 1x20
399.1905 399.1482 399.0174 398.7865 398.4356 397.9372 397.2557 396.3481 395.1648 393.6512 391.7495 389.4022 386.5553 383.1626 379.1897 374.6184 369.4497 363.7065 357.4347 350.7025
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
T_3600=[390.805227064895 390.675779000604 390.286017493707 389.631786740209 388.706477511719 387.501491077936 386.006865145261 384.212038201660 382.106722046482 379.681846109104 376.930531814971 373.849051217311 370.437721882964 366.701690070457 362.651556927952 358.303807952486 353.681014246139 348.811784878551 343.730462361549 338.476567080242]
T_3600 = 1x20
390.8052 390.6758 390.2860 389.6318 388.7065 387.5015 386.0069 384.2120 382.1067 379.6818 376.9305 373.8491 370.4377 366.7017 362.6516 358.3038 353.6810 348.8118 343.7305 338.4766
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
T_18000=[320.784154829704 320.738258431128 320.600846434660 320.372748583693 320.055341615340 319.650540052191 319.160783399654 318.589019837082 317.938686515191 317.213686595838 316.418363192915 315.557470394871 314.636141569840 313.659855173590 312.634398298245 311.565828215802 310.460432184841 309.324685801294 308.165210184555 306.988728298620]
T_18000 = 1x20
320.7842 320.7383 320.6008 320.3727 320.0553 319.6505 319.1608 318.5890 317.9387 317.2137 316.4184 315.5575 314.6361 313.6599 312.6344 311.5658 310.4604 309.3247 308.1652 306.9887
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
T_50000=[298.890663675783 298.888869371201 298.883497300502 298.874579920011 298.862171081446 298.846345670453 298.827199104151 298.804846691191 298.779422858871 298.751080252716 298.719988714898 298.686334148649 298.650317276681 298.612152302323 298.572065482773 298.530293624494 298.487082511280 298.442685276008 298.397360727463 298.351371643918]
T_50000 = 1x20
298.8907 298.8889 298.8835 298.8746 298.8622 298.8463 298.8272 298.8048 298.7794 298.7511 298.7200 298.6863 298.6503 298.6122 298.5721 298.5303 298.4871 298.4427 298.3974 298.3514
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
figure
r = linspace(0,0.05,20)
r = 1x20
0 0.0026 0.0053 0.0079 0.0105 0.0132 0.0158 0.0184 0.0211 0.0237 0.0263 0.0289 0.0316 0.0342 0.0368 0.0395 0.0421 0.0447 0.0474 0.0500
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
plot(r,T_0,'LineWidth',2)
hold on
plot(r,T_100,'LineWidth',2)
hold on
plot(r,T_600,'LineWidth',2)
hold on
plot(r,T_1800,'LineWidth',2)
hold on
plot(r,T_3600,'LineWidth',2)
hold on
plot(r,T_50000,'LineWidth',2)
%legend('0 sec','100 sec','600 sec','1800 sec','3600 sec','50000')
xlabel('r (m)')
ylabel('T (K)')

サインインしてコメントする。

### 回答 (1 件)

Torsten 2024 年 4 月 13 日

I'd compare with the solution from "pdepe":
m = 1;
Ri = 0.05;
nd = 20;
r = linspace(0,Ri,nd);
t = [0 100 600 1800 3600 50000];
sol = pdepe(m,@heatcyl,@heatic,@heatbc,r,t);
u = sol(:,:,1);
plot(r,u(1:6,:),'LineWidth',2)
grid on
function [c,f,s] = heatcyl(x,t,u,dudx)
rho = 821.42;
cp = 1399;
k = 0.1;
c = rho*cp;
f = k*dudx;
s = 0;
end
%----------------------------------------------
function u0 = heatic(x)
u0 = 400;
end
%----------------------------------------------
function [pl,ql,pr,qr] = heatbc(xl,ul,xr,ur,t)
pl = 0; %ignored by solver since m=1
ql = 0; %ignored by solver since m=1
alpha = 5;
Ta = 298.0;
pr = alpha*(ur-Ta);
qr = 1;
end
%----------------------------------------------
But your results look ok at a quick glance.

サインインしてコメントする。

### カテゴリ

Help Center および File ExchangeHeat Transfer についてさらに検索

R2022b

### Community Treasure Hunt

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

Start Hunting!

Translated by