一维非稳态导热问题的数值解
传热学一维稳态和非稳态导热
11.1 通过平壁的一维稳态导热
一、第一类边界条件：表面温度
为常数
单层平壁
a. 几何条件：单层平板；s；
b. 物理条件：、cp、 已知；无内热源；
c. 时间条件：稳态导热， ∂T/∂t=0；
d. 边界条件：第一类。
微分方程式可简化：
2T x2
0
直接积分得： TC1xC2
带入边界 条件：
CC12
Tw1 (Tw2
流量，并记为qL
qL
Q Tw1 Tw2 L 1 d2 ln 传热学一维稳态和非稳态导热
2 d1
单位长度导热热阻
11.2 通过圆筒壁的一维稳态导热
多层圆筒壁
不同材料构成的多层圆筒壁，其导热 热流量可按总温差和总热阻计算
热流量
Q
Tw1 Twn+1 n 1 ln di1
i1 2i L di
单位长度的热流量
2
C1xC2
式中积分常数C1和C2可由边界条件确定，它们分别为：
C2Tw2qv s2; C10
所以，平壁内温度分布为： TTw2qv s2x2
• 可见，该条件下平壁内温度是按抛物线规律分布。令 温度分布关系式中的x=0，则得平壁中心温度为：
qv 2 T T s w
2 传热学一维稳态和非稳态导热
• 设在一管道外面包上一层绝热层(如图所示)。
• 此时单位管长的总热阻γΣ为：
d 1 11 2 1 传热1 学ln 一维d d 稳1 2 态 和非2 稳态1 导x 热lnd d 2 xd 1 x2
11.2 通过圆筒壁的一维稳态导热
• 式中：α1为管内流体与管内壁之间的给热系数，W/m2℃；α2为绝
qL
1
一维非稳态导热 圆柱体 matlab
一维非稳态导热圆柱体 matlab导热是物体中传热作用的一种。
热传导是指物质内部热量的传递及传递机制，描述的是能量（热量）在空间和时间上的传输。
而非稳态导热是指物体内部温度场和热流密度随时间和空间的发展过程。
一维非稳态导热问题是指导热物理学中，只考虑热量沿一个方向传导的问题。
圆柱体是一种常见的几何体，因其在工程领域的广泛应用，研究圆柱体的导热问题具有重要意义。
在研究一维非稳态导热圆柱体问题时，matlab是一种常用的数学软件，其强大的数学运算和可视化功能使得它成为了工程热传导问题求解的重要工具。
通过使用matlab，可以方便地求解一维非稳态导热圆柱体问题，并进行可视化展示。
下面我们将通过matlab来求解一维非稳态导热圆柱体问题，并对结果进行分析。
1. 问题建模假设圆柱体材料均匀，热传导系数为k，密度为ρ，比热容为c。
设圆柱体半径为R，长度为L。
假设圆柱体表面维持恒定的温度T0，初始时刻整个圆柱体的温度场分布为T(x,0) = f(x)，其中f(x)为已知函数。
根据热传导方程，我们可以得到一维非稳态导热圆柱体的数学模型。
2. 热传导方程根据一维热传导方程，我们可以得到圆柱体内部温度场满足的偏微分方程：ρc∂T/∂t = k∇²T3. 离散化为了利用计算机进行求解，我们需要将偏微分方程进行离散化处理。
这里我们可以使用有限差分法（finite difference method）对空间和时间进行离散化。
将圆柱体划分为若干个网格点，并采用显式差分法进行时间推进，就可以得到圆柱体温度场随时间的演化过程。
4. matlab求解在matlab中，我们可以编写程序来实现离散化求解。
首先可以定义圆柱体以及热传导材料的参数，然后通过循环计算每个时间步长内圆柱体温度场的演化，最终得到温度在空间和时间上的分布情况。
借助matlab强大的可视化功能，我们可以直观地展示圆柱体温度场的变化过程。
5. 结果分析得到圆柱体温度场的数值解之后，我们可以对结果进行分析。
传热学-学习课件-4-4 一维非稳态导热问题的数值求解
1
2
a x2
2h cx
2
a x2
t
i
N
1
2h cx
tf
③对称点
t (i)
-1
t (i)
2
传热学 Heat Transfer
2.直接用差分代替微分
①向前差分(forward difference)
i
t
t
i
n
1
t
i
n
n,i
②向后差分(backward difference)
t
t
i
n
t
i
n
1
n,i
i n,i+1
n-1,i n,i n+1,
t
i
n
1
t
n
i
a
t (i1) n 1
2
t
( n
i
1
)
x2
t (i1) n 1
(1,1)
n,i-1 i
n
x
整理成隐式格式：
传热学 Heat Transfer
传热学 Heat Transfer
主讲老师：王舫 适用专业：能源与动力工程专业
传热学 Heat Transfer
§4.4 一维非稳态导热问题的数值求解
在非稳态导热问题中，不但需要对空间区域进 行离散，还需要对时间变量进行离散，接下来以一 个一维非稳态导热问题为例，重点介绍对非稳态项 的离散方法，以及不同离散方法对计算带来的影响 等。
(完整word版)一维非稳态导热的数值计算
传热学C 程序源 二维稳态导热的数值计算2.1物理问题一矩形区域，其边长L=W=1，假设区域内无内热源，导热系数为常数，三个边温度为T1=0，一个边温度为T2=1，求该矩形区域内的温度分布。
2.2 数学描述对上述问题的微分方程及其边界条件为：2222T T0x y∂∂+=∂∂x=0，T=T 1=0x=1，T=T 1=0 y=0，T=T 1=0 y=1，T=T 2=1该问题的解析解：112121(1)sin n n n sh y T T n L x n T T n L sh W L ππππ∞=⎛⎫⋅ ⎪---⎛⎫⎝⎭=⋅ ⎪-⎛⎫⎝⎭⋅ ⎪⎝⎭∑2.3数值离散2.3.1区域离散区域离散x 方向总节点数为N ，y 方向总节点数为M ，区域内任一节点用I,j 表示。
2.3.2方程的离散对于图中所有的内部节点方程可写为：2222,,0i j i jt t x y ⎛⎫⎛⎫∂∂+= ⎪ ⎪∂∂⎝⎭⎝⎭用I,j节点的二阶中心差分代替上式中的二阶导数，得：+1,,-1,,+1,,-1222+2+0i j i j i ji j i j i j T T T T T T x y --+=上式整理成迭代形式：()()22,1,-1,,1,-12222+2()2()i j i j i j i j i j y x T T T T T x y x y ++=++++ (i=2,3……,N-1)，(j=2,3……,M-1)补充四个边界上的第一类边界条件得：1,1j T T = (j=1,2,3……,M) ,1N j T T = (j=1,2,3……,M) ,1i j T T = (i=1,2,3……,N),2i M T T (i=1,2,3……,N)#include<stdio.h> #include<math.h> #define N 10 #define K 11 main() {int i,j,l; float cha;float a,x,y,Fo,Bi; float t[N][K],b[N][K]; /*打印出题目*/printf("\t\t\t 一维非稳态导热问题\t\t"); printf("\n\t\t\t\t\t\t----何鹏举\n"); printf("\n 题目：补充材料练习题三\n");y=1;/*y 代表Δτ*/ x=0.05/(N-1);a=34.89/(7800*712); Fo=(a*y)/(x*x); Bi=233*x/34.89;printf("\n 显示格式条件:");printf("\n1、Fo=%3.1f<0.5\t",Fo);printf("\t2、1-2Fo*Bi-2Fo=%4.2f>0\n\n",1-2*Fo*Bi-2*Fo); /*时刻为零时，赋予初场温度*/ for(i=0;i<N;i++) t[i][0]=1000;/*循环开始，每次计算一个时刻*/ for(j=0;j<K-1;j++) {for(i=0;i<N;i++) b[i][j]=t[i][j];/*下面对每一个时刻进行迭代求解对应的温度分布，公式按传热学课本P178页公式*/cha=1;while(cha>0.001) {for(i=0;i<N-1;i++){if(i==0)t[i][j+1]=Fo*(t[i+1][j]+t[i+1][j])+(1-2*Fo)*t[i][j];/*当计算t[0]时，要用到t[-1],其中t[-1]=t[2]的(对称分布)*/elset[i][j+1]=Fo*(t[i+1][j]+t[i-1][j])+(1-2*Fo)*t[i][j];t[N-1][j+1]=t[N-2][j]*(1-2*Fo*Bi-2*Fo)+2*Fo*t[N-1][j]+2*Fo*Bi*20;/*边界点温度用热平衡法推导出公式*/}cha=0;for(i=0;i<N;i++)cha=cha+abs(t[i][j]-b[i][j]);cha=cha/N;}}/*输出温度分布，其中l控制输出值的排列;这个结果是横轴为x，纵轴为τ的直角坐标下从左上角开始依次的*/printf("\n经数值离散计算的温度分布为：\n");l=0;for(j=K-1;j>=0;j--)for(i=0;i<N;i++){if(t[i][j]>999.99)printf("%6.1f ",t[i][j]);elseprintf("%6.2f ",t[i][j]);l=l+1;if(l==N){printf("\n");l=0;}}getchar();/*为了是生成的exe文件结果算的后不会立即退出，方便观看*/}。
数值传热学一维非稳态导热
数值传热学一维非稳态导热
数值传热学一维非稳态导热是一个拟表达热量输运多方面考虑下的相关分析技术，例如光斑热传递，带有间断层热传导，恒定物质热传导等等。
本文将重点简要介绍一维非稳态导热模型中的理论方法，为解决该问题提供重要基础。
首先，我们讨论的一维非稳态导热模型是一维的，在这种模型中，温度的变化是由上下相邻的单元格热传导加权平均值决定的，从一个单元格到另一个单元格的变化必须满足偏微分方程的通用表达式。
其次，根据以上的假设，一维非稳态导热的数值解将以定义的步长迭代，用于求解温度在不同单元中的变化。
在数值模拟中，需要对边界条件、热导率和温度输入进行有效描述，以确定最终的解答模式。
同时，本次分析中，利用有限差分和蒙特卡罗方法来求解温度场。
这种有趣且可行的做法，不但实现了所需求解的模式，而且能够精确地给出结果。
此外，在电脑指令中，采取该方法对数值运算很有效，从而提高了计算机解的精度和实现的质量。
最后，一维非稳态导热模型是在一定物理场中进行计算的，通用性很强，其能够很好地模拟简单模型中物理场的变化。
因此，它经常被用于诸如热管道传热、滑动轴热传导、负载温度场仿真等多种领域的研究。
总而言之，一维非稳态导热的数值模拟具有良好的数学基础、使用简单的算法以及电脑指令，从而实现快速求解热传导问题的目的，是今后研究的重要课题。
第五章 导热问题的数值方法
5 热传导问题的数值方法5.1一维稳态导热一维稳态导热在直角坐标系下的控制方程可表示为：0)(=+s dxdT k dx d (5-1) 式中k 为导热系数，T 是温度，s 是单位容积的热产生率。
首先选定控制体和网格，如图5.1所示，并对方程（5-1）在所选定的控制体进行积分，即得：0)()(=+-⎰dx s dxdTk dx dT ke w w e (5-2)图5.1 控制体和网格然后进行离散化。
如果用分线段性分布来计算方程（5-2）中的微商dxdT，那么最终的方程为:0)()()()(=∆+---x s x T T k x T T k wW P w e P E e δδ (5-3)假设源项s 在任一控制体中之值可以表示为温度的线性函数，即P P c T s s s +=，则导出的离散化方程为：b T a T a T a W W E E P P ++= (5-4)式中x s b xs a a a x k a x k a c P W E P w wW ee E ∆=∆-+=δ=δ=)()( (5-5) 式(5-4)就是一维稳态导热方程的离散形式，系数a E 和a W 分别代表了节点P 与E 间及W 与P 间导热阻力的倒数，它们的大小反映了节点W 和E 处的温度对P 点的影响程度。
式中的k e 和k w 是控制容积中的e 和w 界面上的当量导热系数。
进行计算时，物理参数值存储在节点的位置上。
为了确定k e 和k w ，还需规定由节点上的物理量来计算相应界面上的量的方法。
常用的方法由两种，即算术平均法与调和平均法。
1、算术平均法假定k 与x 呈线性关系，由P 与E 点的导数系数确定e k 的公式为：eeEe e P e x x k x x k k )()()()(δδ+δδ=-+ （5-6） 2、调和平均法利用传热学的基本公式可以导出确定界面上当量导热系数的调和平均公式。
控制容积中P 和E 的导热系数不相等，但界面上热流密度应该连续，则由Fourier 定律可得：()()()()EePePE EeeE PePe e k x k x T T k x T T k x T T q +-+-δ+δ-=δ-=δ-=（5-7）而()Pe PE e k x T T q δ-=则()()()Ee Pe eek x k x k x +-+=δδδ （5-8）这就是确定界面上当量导热系数的调和平均公式，它反映了串联过程热阻的迭加原则。
一维非稳态导热CRANK-NICOLSON解法
一维非稳态导热CRANK-NICOLSON解法题目：数值计算一维非稳态导热，长度1米的不锈钢棒原来温度都是0度，一端温度突然变为300度，并保存不变，采用CRANK-NICOLSON 方法数值计算不锈钢内温度分布随时间的变化。
解法：一维导热微分方程边界条件为u（0，t）=0；u（a0，t）=300初值u（x，0）=0；主程序clcclearuX=1; %不锈钢长1米uT=2000; %时长2000秒M=10; %空间轴等分区间数N=1000; %时间轴等分区间数rou=8030; %不锈钢密度cp=502.48; %不锈钢热容kk=16.27; %不锈钢导热率D=kk/rou/cp; %扩散系数phi=inline('0'); %初值psi1=inline('0'); %左边界psi2=inline('300'); %右边界%计算步长dx=uX/M;%x的步长dt=uT/N;%t的步长x=(0:M)*dx;r=D*dt/dx/dx;%步长比Diag=zeros(1,M-1);%矩阵的对角线元素Low=zeros(1,M-2);%矩阵的下对角线元素Up=zeros(1,M-2);%矩阵的上对角线元素for i=1:M-2Diag(i)=1+r;Low(i)=-r/2;Up(i)=-r/2;endDiag(M-1)=1+r;%计算初值和边值U=zeros(M+1,N+1);for i=1:M+1U(i,1)=phi(x(i));endfor j=1:N+1U(1,j)=psi1(t(j));U(M+1,j)=psi2(t(j));endB=zeros(M-1,M-1);for i=1:M-2B(i,i)=1-r;B(i,i+1)=r/2;B(i+1,i)=r/2;endB(M-1,M-1)=1-r;%逐层求解，需要使用追赶法（调用函数EqtsForwardAndBackward）for j=1:Nb1=zeros(M-1,1);b1(1)=r*(U(1,j+1)+U(1,j))/2;b1(M-1)=r*(U(M+1,j+1)+U(M+1,j))/2;b=B*U(2:M,j)+b1;U(2:M,j+1)=zhuiganfa(Low,Diag,Up,b);endU=U';%作出图形xlabel('空间变量x')ylabel('时间变量t')shading interp程序用到了追赶法子程序，代码如下function x=zhuiganfa(L,D,U,b)%追赶法求解三对角线性方程组Ax=b%检查参数的输入是否正确n=length(D);m=length(b);n1=length(L);n2=length(U);if n-n1 ~= 1 || n-n2 ~= 1 || n ~= mdisp('输入参数有误！')x=' ';return;end%追的过程for i=2:nL(i-1)=L(i-1)/D(i-1);D(i)=D(i)-L(i-1)*U(i-1);endx=zeros(n,1);x(1)=b(1);for i=2:nx(i)=b(i)-L(i-1)*x(i-1);end%赶的过程x(n)=x(n)/D(n);for i=n-1:-1:1x(i)=(x(i)-U(i)*x(i+1))/D(i);endreturn;运行主程序，最终得到如图所示结果。
柱坐标系含内热源一维非稳态热传导问题
2.11解：对于该一维有热源的非定常热传导问题,设导热系数k 为常数有211t t r c k r u r r r R ρτ⎡⎤∂∂∂⎛⎫⎛⎫'''=+-⎢⎥ ⎪ ⎪∂∂∂⎝⎭⎝⎭⎢⎥⎣⎦初始条件(),t r o T ∞=边界条件为()00r r R r Rtr t h t T r λ=∞==⎧∂=⎪∂⎪⎨∂⎪-=-⎪∂⎩即20111r a r u r r rc R θθτρ⎡⎤∂∂∂⎛⎫⎛⎫'''=+-⎢⎥ ⎪ ⎪∂∂∂⎝⎭⎝⎭⎢⎥⎣⎦(),0r o θ=000r r R r k h r θθθ==⎧∂=⎪∂⎪⎨∂⎛⎫⎪+= ⎪⎪∂⎝⎭⎩其中ka cρ=为热扩散率或热扩散系数（thermal diffusivity ），过余温度t T θ∞=-，则 令()u f r θ=+，并使其中()f r 满足原方程及边界条件不妨设()420122416u r f r r C r C k R '''⎛⎫=--++ ⎪⎝⎭，显然其满足原方程，现用边界条件待定其中系数301112034200112222002084001684416163416r x r R r Ru f r r C C C rk R u u f r r k hf k r C h r C r C r k R k R u R u R k h C k k u ====⎛⎫'''⎛⎫∂=--+==⇒= ⎪ ⎪∂⎝⎭⎝⎭⎡⎤⎛⎫⎛⎫''''''⎛⎫⎛⎫∂⎛⎫+=--++--++⎢⎥ ⎪ ⎪ ⎪ ⎪ ⎪∂⎝⎭⎢⎥⎝⎭⎝⎭⎝⎭⎝⎭⎣⎦''''''⎛⎫⎛⎫=-+-+ ⎪ ⎪⎝⎭⎝⎭'''=-202310441344R hR C h k u R R C h k ⎛⎫++= ⎪⎝⎭'''⎛⎫⇒=+ ⎪⎝⎭故()420021341644u u R r R f r r k R h k ''''''⎛⎫⎛⎫=--++ ⎪ ⎪⎝⎭⎝⎭420021341644u u R r R u r k R h k θ''''''⎛⎫⎛⎫=--++ ⎪ ⎪⎝⎭⎝⎭ 代入原方程及边界条件中有1u u a r r r r τ∂∂∂⎛⎫= ⎪∂∂∂⎝⎭000r r R ur u k hu r==⎧∂=⎪∂⎪⎨∂⎛⎫⎪+= ⎪⎪∂⎝⎭⎩ 代入初始条件中有420020013401644u u R r R u r k R h k ττθ==''''''⎛⎫⎛⎫=--++= ⎪ ⎪⎝⎭⎝⎭即()0u r τϕ==其中()420021341644u u R r R r r k R h k ϕ''''''⎛⎫⎛⎫=--+ ⎪ ⎪⎝⎭⎝⎭ 现用分离变量法进行求解令()()(),u r R r T ττ=带入(1.4)中泛定方程中有()()()()()()1R r T a R r T rR r T r τττ''''=+⎡⎤⎣⎦ ()()()'()()R r rR r T rR r aT τλτ'''+==-即()()()'()()0,0.T aT rR r R r rR r τλτλ⎧+=⎪⎨'''++=⎪⎩由边界条件知'(0)0,'()()0.hR R R R R k=+=由初始条件有()()()0R r T r ϕ=现对k 的值进行讨论① 当0λ=时，将其带入式(1.6)中有()()'()0,0.T rR r R r τ⎧=⎪⎨'''+=⎪⎩ ()()()0ln T C R r A r BK x τ=⎧⎪⇒⎨=+⎪⎩(1.9)代入边界条件中有0A B ==显然不符合题意舍去② 当0λ<时，对于式(1.6)令x =有'2222()()0,0.T aT R R x x x R x x τλτ⎧+=⎪⎨∂∂+-=⎪∂∂⎩()()()()1200a T Ce R x AI x BK x λττ-⎧⎪=⇒⎨⎪=+⎩()()))1200a T Ce R r AI BK λττ-⎧=⎪⇒⎨⎪=+⎩注意到当τ→∞时()τ→∞T 即最终的温度→∞t ，这在实际情况下是不可能的，故该种情况应当舍去。
一维非稳态导热的数值计算
一维非稳态导热的数值计算一维非稳态导热问题是指材料的温度在时间上发生变化，且只沿一个方向进行传热的问题。
这种问题在实际工程、材料科学和热传导研究中都十分常见。
数值计算是求解这类问题的重要方法之一，接下来我们将介绍一维非稳态导热的数值计算方法。
首先，我们来定义一维非稳态导热的数学模型。
假设我们考虑的材料是一维的杆状物体，其长度为L，温度分布随时间t和空间x而变化，记作T(x,t)。
根据热传导方程，我们可以得到如下的一维非稳态导热方程：∂T/∂t=α*∂^2T/∂x^2其中，α是热扩散系数，反应了材料导热性能的指标。
我们的目标是求解在给定边界条件下的温度分布T(x,t)。
为了使用数值方法求解该方程，我们需要将其离散化。
首先，我们将时间t离散化为一系列的时间步长Δt，将空间x离散化为一系列的空间步长Δx。
然后，我们使用中心差分法来近似替代方程的二阶空间导数项和一阶时间导数项：∂T/∂t≈(T(i,j+1)-T(i,j))/Δt∂^2T/∂x^2≈(T(i+1,j)-2T(i,j)+T(i-1,j))/Δx^2其中，i和j分别表示空间和时间的离散节点索引。
将上述近似代入导热方程中，得到离散的差分方程：(T(i,j+1)-T(i,j))/Δt=α*(T(i+1,j)-2T(i,j)+T(i-1,j))/Δx^2根据上述差分方程，我们可以通过迭代计算来逐步更新温度分布。
首先，我们需要给定初始条件T(x,0)和边界条件T(0,t)和T(L,t)。
然后，我们通过迭代计算来更新温度值，直到达到所需的时间步长和空间步长。
具体来说，我们可以根据以下的更新公式进行迭代计算：T(i,j+1)=T(i,j)+α*Δt*(T(i+1,j)-2T(i,j)+T(i-1,j))/Δx^2其中，T(i,j+1)表示在第j+1个时间步长和第i个空间步长的温度值，T(i,j)表示在第j个时间步长和第i个空间步长的温度值。
总之，一维非稳态导热的数值计算方法可以使用离散化和迭代计算来求解热传导方程。
第四章_导热问题的数值方法
5 热传导问题的数值方法5.1一维稳态导热一维稳态导热在直角坐标系下的控制方程可表示为：0)(=+s dxdT k dx d (5-1) 式中k 为导热系数，T 是温度，s 是单位容积的热产生率。
首先选定控制体和网格，如图5.1所示，并对方程（5-1）在所选定的控制体进行积分，即得：0)()(=+-⎰dx s dxdTk dx dT k e w w e (5-2)图5.1 控制体和网格然后进行离散化。
如果用分线段性分布来计算方程（5-2）中的微商dxdT，那么最终的方程为:0)()()()(=∆+---x s x T T k x T T k wW P w e P E e δδ (5-3)假设源项s 在任一控制体中之值可以表示为温度的线性函数，即P P c T s s s +=，则导出的离散化方程为：b T a T a T a W W E E P P ++= (5-4)式中x s b xs a a a x k a x k a c P W E P w wW ee E ∆=∆-+=δ=δ=)()( (5-5) 式(5-4)就是一维稳态导热方程的离散形式，系数a E 和a W 分别代表了节点P 与E 间及W 与P 间导热阻力的倒数，它们的大小反映了节点W 和E 处的温度对P 点的影响程度。
式中的k e 和k w 是控制容积中的e 和w 界面上的当量导热系数。
进行计算时，物理参数值存储在节点的位置上。
为了确定k e 和k w ，还需规定由节点上的物理量来计算相应界面上的量的方法。
常用的方法由两种，即算术平均法与调和平均法。
1、算术平均法假定k 与x 呈线性关系，由P 与E 点的导数系数确定e k 的公式为：eeE e e P e x x k x x k k )()()()(δδ+δδ=-+ （5-6）2、调和平均法利用传热学的基本公式可以导出确定界面上当量导热系数的调和平均公式。
控制容积中P 和E 的导热系数不相等，但界面上热流密度应该连续，则由Fourier 定律可得：()()()()EePePE EeeE PePe e k x k x T T k x T T k x T T q +-+-δ+δ-=δ-=δ-=（5-7）而()Pe PE e k x T T q δ-=则()()()Ee Pe eek x k x k x +-+=δδδ （5-8）这就是确定界面上当量导热系数的调和平均公式，它反映了串联过程热阻的迭加原则。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
一维非稳态导热问题的数
值解
Prepared on 22 November 2020
计算传热学程序报告
题目：一维非稳态导热问题的数值解
姓名：
学号：
学院：能源与动力工程学院
专业：工程热物理
日期：2014年5月25日
一维非稳态导热问题数值解
求解下列热传导问题：
1.方程离散化
对方程进行控制体积分得到：
非稳态项：选取T 随x 阶梯式变化，有
扩散项：选取一阶导数随时间做显示变化，有
进一步取T 随x 呈分段线性变化，有
e P E e x T T x T )()(δ-=∂∂，w
W P w x T T x T )()(δ-=∂∂ 整理可以得到总的离散方程为：
2.计算空间和时间步长
取空间步长为：
h=L/N
网格Fourier 数为：
2
20x t x t
F ∆∆=∆∆=α（小于时稳定） 时间步长为：
3.建立温度矩阵与边界条件
T=ones(N+1,M+1)
T(:,1)=Ti(初始条件温度都为0）
T(1,:)=To(边界条件x=0处温度为1）
T(N+1,:)=Te(边界条件x=L 处温度为0）
4.差分法求解温度
由离散方程可得到：
转化为相应的温度矩阵形式：
5.输入界面
考虑到方程的变量，采用inputdlg函数设置5个输入变量，对这5个变量设置了默认值，如图1所示。
在计算中可以改变不同的数值，得到不同的结果，特别注意稳定条件的临界值是。
根据设置的默认值，得到的计算结果如图2所示。
图1matlab变量输入界面
图2默认值的计算结果
6.结果分析
根据上面的分析，给出了程序的输入界面，以及默认值状态下的数值解。
可以通过改变不同的输入值，得到需要的分析结果，总结出了下面4点结论：（1）取F0=，得到一维非稳态导热结果如下图所示
图2F0=时一维非稳态导热
从图中可以看出，对于长度L=1的细杆，初始时刻t=0时温度为0，边界条件x=0时，T=1,边界条件x=1时，T=0。
随着时间的增加，温度从x=0通过导热的形式传递到x=1，不同时刻不同位置杆的温度都不同，并且随着时间的增加，杆的温度也逐渐增加。
（2）取F0=，可以得到不同位置的温度响应曲线，如下图所示
图3F0=时不同x位置处的温度响应
图中红色曲线代表x=位置的温度瞬态响应，黑色曲线代表x=位置的温度瞬态响应，蓝色曲线代表x=位置的温度瞬态响应。
从图中可以看出，随着x的增加，曲线与x轴的交点值越大，温度开始传递到该位置的所需的时间越长。
随着x的增加，温度响应曲线的变化速率越慢，最终的达到的温度也越低。
（3）取F0=，得到不同位置的温度响应曲线如下图所示
图4F0=时不同x位置处的温度响应
图中三条曲线分别是x=，x=，x=位置的温度瞬态响应。
与图3的F0=进行对比，两种情况下的F0值不同，F0值越大表明热扩散系数 的值越大。
从图中可以看出热扩散系数对于导热的影响，F0=时，与F0=相比较，各位置开始响应时所需的时间较长，而且各位置响应曲线的变化速率较小，最终的达到的温度也较低，说明了热扩散系数越小，热传导越慢，传递效率越低。
（4）取F0=，得到非稳定的数值解如图所示
图5F0=时一维非稳态导热
图6F0=时不同x位置处的温度响应
从图中可以看出，对于显示格式的离散方程，并不是所有的F0值都能得到有意义的解，必须要求F0<时才能得到稳定的数值解，当F0>时，会出现物理上不真实的解。
附件：（matlab程序）
functionheat_conduction()%一维齐次热传导方程
%设置输入界面
options={'空间杆长L','空间点数N','时间点数M','扩散系数a','稳定条件的值Fo(临界值',};
topic='一维非稳态导热';%标题栏显示
lines=1;%输入行为1行
def={'1','100','1000','1',''};%默认值输入
f=inputdlg(options,topic,lines,def);%输入框设置
L=eval(f{1});%设置输入值
N=eval(f{2});
M=eval(f{3});
a=eval(f{4});
Fo=eval(f{5});%Fo的值必须小于，小于波动
%计算空间步长与时间步长
h=L/N;%空间步长
x1=0:h:L;
x=x1';
n=Fo*h^2/a;%时间步长
tm=n*M;%传导总时间
t1=0:n:tm;
t=t1';
%计算初始条件与边界条件
Ti=x.*0;%初始条件
To=1+t.*0;%x=0的边界条件
Te=t.*0;%x=L的边界条件
%建立温度矩阵T
T=ones(N+1,M+1);
T(:,1)=Ti;%第一列为初始条件
T(1,:)=To;%第一行为x=0边界条件
T(N+1,:)=Te;%最后一行为x=L边界条件
%利用差分法求解温度矩阵T
fork=1:M
m=2;
whilem<=N;
T(m,k+1)=Fo*(T(m+1,k)+T(m-1,k)-2*T(m,k))+T(m,k);
m=m+1;
end
end
%将时间空间的一维坐标转化为二维坐标[Y,X]=meshgrid(t1,x);
%根据温度矩阵T绘图
subplot(2,2,1);
mesh(X,Y,T);%三维图绘制
view([1,-1,1]);%调整视图角度
title('非稳态导热');%图像名称
xlabel('长度x');%x轴名称
ylabel('时间t');%y轴名称
zlabel('温度T');%z轴名称
subplot(2,2,2);
A=T(11,:);%取矩阵第11列的值
plot(A,'r');%二维曲线绘制
legend('A=');%显示函数名称
title('x=瞬态响应');
xlabel('时间t');
ylabel('温度T');
axis([0100001]);%坐标轴数值范围subplot(2,2,3);
B=T(21,:);%取矩阵第21列
plot(B,'k');
legend('B=');
title('x=瞬态响应');
xlabel('时间t');
ylabel('温度T');
axis([0100001]);
subplot(2,2,4);
C=T(41,:);%取矩阵第41列plot(A,'r');
holdon;%多条曲线绘制plot(B,'k');
plot(C);
holdoff;
title('瞬态响应');
xlabel('时间t');
ylabel('温度T');
axis([0100001]);
legend('A=','B=','C=');。