数值传热学陶文铨第四章作业(完整资料).doc

合集下载

计算传热学_高等教育-实验设计

坐标系不同，控制方程的形式不尽相同

必要的简化与化简
2.1 控制方程

传热的三种模式（Modes of heat transfer)

热传导（Thermal conduction) 热对流（Thermal advection)

对流换热（Convection heat transfer) 辐射换热（Radiation heat transfer)
数值方法
分析解法与实验研究

分析解法

成本最低结果最理想影响因素表达清楚缺点：局限与非常简单的问题成本较低：数值实验适用范围宽缺点：可靠性差，表达困难可靠成本高

数值方法

实验研究

将三种方法有机结合，互为补充，必然会取得相得益彰的效果
第2讲
传热问题的数学描述
1) 2)
将上面的数学模型无量纲化，并给出其分析解；取β=1，就 PeL=(ρuL)/Γ=1、10、100 三种情况分别用三点中心差分格式、迎风格式、幂律格式和 QUICK 格式进行计算，并与分析解比较（计算时节点数目可取为 10 ~ 20）； 3) 改变参数β，譬如取β=10，重复 2）中的计算；分析 2）和 3）中得到的结果，对各种格式进行比较。
Tf h A B
Tf h
δ
δ
计算传热学习题之三
考虑下述一维稳态对流-扩散问题,
d d dU ( ρuU ) = (Γ )+s dx dx dx U x=0 = U 0 U
x=L
= UL
其中 u 是流速，Γ和ρ均为常数，而 s 是 x 的单值函数，
s = 0.5 β U 0 −U L L

传热学第四版课后作业答案(杨世铭-陶文铨)]

1-9 一砖墙的表面积为122m ，厚为260mm ，平均导热系数为1.5W/（m.K ）。

设面向室内的表面温度为25℃，而外表面温度为-5℃，试确定次砖墙向外界散失的热量。

解：根据傅立叶定律有：WtA9.207626.05)(25125.1=--⨯⨯=∆=Φδλ1-12 在一次测定空气横向流过单根圆管的对流换热实验中，得到下列数据：管壁平均温度t w =69℃，空气温度t f =20℃，管子外径 d=14mm ，加热段长 80mm ，输入加热段的功率8.5w ，如果全部热量通过对流换热传给空气，试问此时的对流换热表面传热系数多大？解：根据牛顿冷却公式()f w t t rlh q -=π2所以()f w t t d qh -=π＝49.33W/(m 2.k)1-20 半径为0.5 m 的球状航天器在太空中飞行，其表面发射率为0.8。

航天器内电子元件的散热总共为175W 。

假设航天器没有从宇宙空间接受任何辐射能量，试估算其表面的平均温度。

解：电子原件的发热量＝航天器的辐射散热量即：4T Q εσ=4A QT εσ=∴ =187K 热阻分析 ;;2-4 一烘箱的炉门由两种保温材料A 及B 组成，且B A δδ2=(见附图)。

已知)./(1.0K m W A =λ,)./(06.0K m W B =λ,烘箱内空气温度4001=f t ℃，内壁面的总表面传热系数)./(501K m W h =。

为安全起见，希望烘箱炉门的外表面温度不得高于50℃。

设可把炉门导热作为一维问题处理，试决定所需保温材料的厚度。

环境温度=2f t 25℃，外表面总传热系数)./(5.922K m W h =。

解：热损失为()()22111f f BBA A fwf t t h t t h t t q -+-=+-=λδλδ又50=fw t ℃；B A δδ=联立得m m B A 039.0;078.0==δδ2-9 双层玻璃窗系由两层厚为6mm 的玻璃及其间的空气隙所组成，空气隙厚度为8mm 。

《传热学》第四章导热数值解法基础

——以右边界为例 1.第一类边界条件：
边界
2.第二类边界条件：
Байду номын сангаас
Δx＝Δy时简化为：
绝热边界：
3.第三类边界条件：
Δx＝Δy时简化为：
其他情况的节点方程 ——见教材表4－1
外拐角与内拐角节点
对流边界内部拐角节点热平衡：
节点方程式推导实例 ——对流边界外部拐角节点
Δx＝Δy时简化为：
数值导热离散方程组＝内节点离散方程＋边界节点离散方程
二、常用计算软件
1.MATLAB——矩阵计算软件
matlab软件主界面
2.FLUENT——流体流动通用数值计算软件
3. FLUENT AIRPAK ——人工环境系统分析软件，暖通空调专业和传热学领域必备软件
AIRPAK模拟温度场
第四章重点： 1.有限差分方程的建立 2.高斯－赛德尔迭代方法
谢谢观看
《传热学》
第四章导热数值解法基础
本章研究的目的 ——利用计算机求解难以用分析解求解的导热问题基本思想 ——把原来在时间、空间坐标系中连续的物理量的场，用有限个离散点的值的集合来代替，通过求解按一定方法建立起来的关于这些值的代数方程，来获得离散点上被求物理量的值。研究手段——有限差分法
数值导热离散方程组内节点离散方程边界节点离散方程三节点离散方程组的求解迭代法迭代法的原理离散方程组的求解方法消元法方程过多时计算机内存不足迭代法假定初值根据假定的初值求新值并重复此步骤若干次两次计算值足够接近认为达到真实值简单迭代法每次迭代时使用上次迭代的结果允许误差简单迭代法的缺点由于每次迭代中使用与真实值偏差较大的上次迭代的旧值使运算过程接近真实值的时间增加高斯赛德尔迭代法将本次迭代的最新结果立刻代入本次迭代过程计算其他未知值高斯赛德尔迭代法的优点由于每次迭代中使用与真实值偏差较小的本次迭代的新值使运算过程接近真实值的时间缩短第三节非稳态导热的数值计算一显式差分格式研究对象一维非稳态导热问题一维非稳态导热内节点差分方程

数值传热学陶文铨主编第二版习题答案

1
x1=x; x=t(3,1); end tcal=t
习题 4-12 的 Matlab 程序
%代数方程形式 AiTi=CiTi+1+BiTi-1+Di mdim=10;%计算的节点数 x=linspace(1,3,mdim);%生成 A、C、B、T 数据的基数； A=cos(x);%TDMA 的主对角元素 B=sin(x);%TDMA 的下对角线元素 C=cos(x)+exp(x); %TDMA 的上对角线元素 T=exp(x).*cos(x); %温度数据 %由 A、B、C 构成 TDMA coematrix=eye(mdim,mdim); for n=1:mdim coematrix(n,n)=A(1,n); if n>=2 coematrix(n,n-1)=-1*B(1,n); end if n<mdim coematrix(n,n+1)=-1*C(1,n); end end %计算 D 矢量 D=(coematrix*T')'; %由已知的 A、B、C、D 用 TDMA 方法求解 T %消元 P(1,1)=C(1,1)/A(1,1); Q(1,1)=D(1,1)/A(1,1); for n=2:mdim P(1,n)=C(1,n)/(A(1,n)-B(1,n)*P(1,n-1)); Q(1,n)=(D(1,n)+B(1,n)*Q(1,n-1))/(A(1,n)-B(1,n)*P(1,n-1)); end %回迭 Tcal(1,mdim)=Q(1,mdim); for n=(mdim-1):-1:1 Tcal(1,n)=P(1,n)*Tcal(1,n+1)+Q(1,n); end Tcom=[T;Tcal]; %绘图比较给定 T 值和计算 T 值 plot(Tcal,'r*') hold on plot(T)

数值传热学第4章作业

习题4－2图4-22 习题4-2插图[解]一维稳态导热问题的控制方程为：0=+⎪⎭⎫ ⎝⎛S dx dT dx d λ 4-2-1 该问题的边界条件为：()⎪⎩⎪⎨⎧=-=-==2,0,1001x T T h dx dT x T f λ 4-2-2分别对节点2，3进行离散，将已知数据代入离散格式中，得到方程组：130232=-T T 4-2-375432=+-T T 4-2-4 联立式（4-2-3）、式（4-2-4），可以解出2T ，3T ： 852=T ，403=T 。

下面验证总体守恒性：4-2-5右端3放出的热量为：()()30020401533=-⨯=-=f T T h Q 4-2-6在总体容积内部产生的热量为：2.0150 2.0300S Q S x =⨯∆=⨯=还需要证明左端是绝热条件：节点2的热平衡为：21851000.550.5150757501T T xS x λ--+∆=+⨯=-+=∆ 左端绝热，所以计算结果符合总体能量守恒。

习题 4-5[解] 根据习题4-2的分析，可以得到节点2的离散方程：130232+=T T 4-5-1对于节点3，应用边界条件：()()1324330.510f f T T S T T T T xλδ--+⨯=-- 4-5-2式（4-5-2）可以整理成：()5432355751020T T T =+-- 4-5-3采用局部线性化方法，可以得到：()()()()515***444333331020102012.520T T T T T -=-+-- 4-5-4节点3的离散方程表示成：()()()51**44323335575 2.52012.52020T T T T T =++---- 4-5-5迭代求解得出：2382.82;35.64T T == 检验热平衡：内热源生成热1300φ=；右端散热5/4210(35.6420)311.0h T φ=∆=-=左端散热382.821000.5150510.91φ-=⨯+⨯=-所以123()30031110.90φφφ-+=-+≅不作热平衡扣0.5 分。

数值传热学陶文铨第四章作业

4-1解：采用区域离散方法A 时；网格划分如右图。

内点采用中心差分123278.87769.9T T T === 22d T T=0dx - 有 i+1i 122+T 0i i T T T x ---=∆ 将2点,3点带入321222+T 0T T T x --=∆ 即321209T T -+= 432322+T 0T T T x --=∆4321322+T 0T T T x --=∆ 即4321209T T T -+-= 边界点4（1）一阶截差由x=1 1dT dx =，得 4313T T -= （2）二阶截差 11B M M q x x x T T S δδλλ-=++V 所以 434111. 1.36311T T T =++ 即 43122293T T -=采用区域离散方法B 22d T T=0dx - 由控制容积法 0w e dT dT T x dT dT ⎛⎫⎛⎫--∆= ⎪ ⎪⎝⎭⎝⎭ 所以代入2点4点有322121011336T T T T T ----= 即 239028T T -= 544431011363T T T T T ----= 即 34599 02828T T T -+=对3点采用中心差分有 432322+T 013T T T --=⎛⎫ ⎪⎝⎭ 即 2349901919T T T -+= 对于点5 由x=1 1dT dx =，得 5416T T -= （1）精确解求左端点的热流密度由 ()21x x e T e e e -=-+ 所以有 ()220020.64806911x x x x dTe e q e e dx e e λ-====-+=-=++ （2）由A 的一阶截差公式（3）由B 的一阶截差公式（4）由区域离散方法B 中的一阶截差公式：通过对上述计算结果进行比较可得：区域离散B 有控制容积平衡法建立的离散方程与区域离散方程A 中具有二阶精度的格式精确度相当！ 4-3解：将平板沿厚度方向3等分，如图由题可知该导热过程可看作无限大平板的一维稳态有源导热问题，则控制方程为x=0， T 0=75℃x=0.1 dT =h(T-T )dx f λ-1点，2点采用中心差分有21022+T 0T T S xλ-+=∆ （1） 32122+T 0T T S x λ-+=∆ （2）右端点采用一阶截差的离散231f hx T T T x h λλ⎡⎤+⎢⎥⎣⎦=⎛⎫+ ⎪⎝⎭V （3）右端点采用二阶截差的离散代入（1）（2）（3）得1223132280.62 5.67625T T T T T T T -=--=-= 解得123278.87769.9T T T ===代入（4）得解得 12380.6380.6675.1T T T ===精确解 22d T +S=0dxλ （4） x=0， T 0=75℃ （5） x=0.1 dT =h(T-T )dxf λ- （6）代入数据积分的将 x 1=10.13⨯，x 2=20.13⨯, x 3=0.1T 1=80.56 T 2=80.56 T 3=75.1通过比较可得右端点采用二阶截差的离散更接近真实值。

【免费下载】传热学补充教材

目录
一、教学进度计划表……………………………………...…..………….………..3 二、思考题和练习题………………………………………………...…………….4 三、补充习题…………………………………………………….……...…..……..4 四、复习提纲………………………………………………………………………7 五、计算机实习指导书…………………………………………….…….....……..11
1. 练习题一：一维稳态导热的数值计算………………………………………..11 2. 练习题二：二维稳态导热的数值计算……………………………………….14 3. 练习题三：一维非稳态导热的数值计算…………………………………….15
六、补充内容：绕流平板的流动和换热,无量纲参数,比拟和相似………….19
1. 物理问题……………………………………………………..……………19 2. 数学模型…………………………………………………………………..19 3. 数学模型的无量纲化………………………………………………………19 4. 局部对流换热系数和平均对流换热系数……………………………………20 5. 局部阻力系数和平均阻力系数……………………………………………..22 6. 无量纲参数及其关联……………………………………………………....22 7. 比拟和相似………………………………………………………………..23
1
对全部高中资料试卷电气设备，在安装过程中以及安装结束后进行高中资料试卷调整试验；通电检查所有设备高中资料电试力卷保相护互装作置用调与试相技互术关，系电，力根通保据过护生管高产线中工敷资艺设料高技试中术卷资，配料不置试仅技卷可术要以是求解指，决机对吊组电顶在气层进设配行备置继进不电行规保空范护载高与中带资负料荷试下卷高问总中题体资，配料而置试且时卷可，调保需控障要试各在验类最；管大对路限设习度备题内进到来行位确调。保整在机使管组其路高在敷中正设资常过料工程试况中卷下，安与要全过加，度强并工看且作护尽下关可都于能可管地以路缩正高小常中故工资障作料高；试中对卷资于连料继接试电管卷保口破护处坏进理范行高围整中，核资或对料者定试对值卷某，弯些审扁异核度常与固高校定中对盒资图位料纸置试，.卷保编工护写况层复进防杂行腐设自跨备动接与处地装理线置，弯高尤曲中其半资要径料避标试免高卷错等调误，试高要方中求案资技，料术编试交写5、卷底重电保。要气护管设设装线备备置敷4高、调动设中电试作技资气高，术料课中并3中试、件资且包卷管中料拒含试路调试绝线验敷试卷动槽方设技作、案技术，管以术来架及避等系免多统不项启必方动要式方高，案中为；资解对料决整试高套卷中启突语动然文过停电程机气中。课高因件中此中资，管料电壁试力薄卷高、电中接气资口设料不备试严进卷等行保问调护题试装，工置合作调理并试利且技用进术管行，线过要敷关求设运电技行力术高保。中护线资装缆料置敷试做设卷到原技准则术确：指灵在导活分。。线对对盒于于处调差，试动当过保不程护同中装电高置压中高回资中路料资交试料叉卷试时技卷，术调应问试采题技用，术金作是属为指隔调发板试电进人机行员一隔，变开需压处要器理在组；事在同前发一掌生线握内槽图部内纸故，资障强料时电、，回设需路备要须制进同造行时厂外切家部断出电习具源题高高电中中源资资，料料线试试缆卷卷敷试切设验除完报从毕告而，与采要相用进关高行技中检术资查资料和料试检，卷测并主处且要理了保。解护现装场置设。备高中资料试卷布置情况与有关高中资料试卷电气系统接线等情况，然后根据规范与规程规定，制定设备调试高中资料试卷方案。

传热学思考题参考答案(陶文铨第四版)

传热学思考题参考答案第一章：1、用铝制水壶烧开水时，尽管炉火很旺，但水壶仍安然无恙。

而一旦壶内的水烧干后水壶很快就被烧坏。

试从传热学的观点分析这一现象。

答：当壶内有水时，可以对壶底进行很好的冷却（水对壶底的对流换热系数大），壶底的热量被很快传走而不至于温度升得很高；当没有水时，和壶底发生对流换热的是气体，因为气体发生对流换热的表面换热系数小，壶底的热量不能很快被传走，故此壶底升温很快，容易被烧坏。

2、什么是串联热阻叠加原则，它在什么前提下成立以固体中的导热为例，试讨论有哪些情况可能使热量传递方向上不同截面的热流量不相等。

答：在一个串联的热量传递过程中，如果通过每个环节的热流量都相同，则各串联环节的总热阻等于各串联环节热阻的和。

例如：三块无限大平板叠加构成的平壁。

例如通过圆筒壁，对于各个传热环节的传热面积不相等，可能造成热量传递方向上不同截面的热流量不相等。

第二章：1、扩展表面中的导热问题可以按一维问题处理的条件是什么有人认为,只要扩展表面细长,就可按一维问题处理,你同意这种观点吗|答：条件：（1）材料的导热系数，表面传热系数以及沿肋高方向的横截面积均各自为常数（2）肋片温度在垂直纸面方向（即长度方向）不发生变化，因此可取一个截面（即单位长度）来分析（3）表面上的换热热阻远远大于肋片中的导热热阻，因而在任一截面上肋片温度可认为是均匀的（4）肋片顶端可视为绝热。

并不是扩展表面细长就可以按一维问题处理，必须满足上述四个假设才可视为一维问题。

2、肋片高度增加引起两种效果:肋效率下降及散热表面积增加。

因而有人认为随着肋片高度的增加会出现一个临界高度，超过这个高度后，肋片导热热流量会下降，试分析该观点的正确性。

答：的确肋片高度增加会导致肋效率下降及散热表面积增加，但是总的导热量是增加的，只是增加的部分的效率有所减低，所以我们要选择经济的肋片高度。

第三章：1、由导热微分方程可知，非稳态导热只与热扩散率有关,而与导热系数无关。

(完整版)数值传热学陶文铨主编第二版习题答案

数值传热学4-9章习题答案习题4－2一维稳态导热问题的控制方程：022=+∂∂S xTλ依据本题给定条件，对节点2节点3采用第三类边界条件具有二阶精度的差分格式，最后得到各节点的离散方程：节点1：1001=T 节点2：1505105321-=+-T T T 节点3：75432=+-T T 求解结果：，852=T 403=T 对整个控制容积作能量平衡，有：2150)4020(15)(3=⨯+-⨯=∆+-=∆+x S T T h x S q f f B 即：计算区域总体守恒要求满足习题4－5在4－2习题中，如果，则各节点离散方程如下：25.03)(10f T T h -⨯=节点1：1001=T 节点2：1505105321-=+-T T T 节点3：25.03325.032)20(4015])20(21[-⨯+=-⨯++-T T T T 对于节点3中的相关项作局部线性化处理，然后迭代计算；求解结果：，（迭代精度为10-4）818.822=T 635.353=T 迭代计算的Matlab 程序如下：x=30;x1=20;while abs(x1-x)>0.0001a=[1 0 0;5 -10 5;0 -1 1+2*(x-20)^(0.25)]; b=[100;-150; 15+40*(x-20)^(0.25)]; t=a^(-1)*b;x1=x;x=t(3,1);endtcal=t习题4-12的Matlab程序%代数方程形式A i T i=C i T i+1+B i T i-1+D imdim=10;%计算的节点数x=linspace(1,3,mdim);%生成A、C、B、T数据的基数；A=cos(x);%TDMA的主对角元素B=sin(x);%TDMA的下对角线元素C=cos(x)+exp(x); %TDMA的上对角线元素T=exp(x).*cos(x); %温度数据%由A、B、C构成TDMAcoematrix=eye(mdim,mdim);for n=1:mdimcoematrix(n,n)=A(1,n);if n>=2coematrix(n,n-1)=-1*B(1,n);endif n<mdimcoematrix(n,n+1)=-1*C(1,n);endend%计算D矢量D=(coematrix*T')';%由已知的A、B、C、D用TDMA方法求解T%消元P(1,1)=C(1,1)/A(1,1);Q(1,1)=D(1,1)/A(1,1);for n=2:mdimP(1,n)=C(1,n)/(A(1,n)-B(1,n)*P(1,n-1));Q(1,n)=(D(1,n)+B(1,n)*Q(1,n-1))/(A(1,n)-B(1,n)*P(1,n-1)); end%回迭Tcal(1,mdim)=Q(1,mdim);for n=(mdim-1):-1:1Tcal(1,n)=P(1,n)*Tcal(1,n+1)+Q(1,n);endTcom=[T;Tcal];%绘图比较给定T值和计算T值plot(Tcal,'r*')hold onplot(T)n gin th a r e 结果比较如下，由比较可知两者值非常切合（在小数点后8位之后才有区别）：习题4-14充分发展区的温度控制方程如下：)(1rTr r r x T uc p ∂∂∂∂=∂∂λρ对于三种无量纲定义、、进行分析如下w b w T T T T --=Θ∞∞--=ΘT T T T w ww T T T T --=Θ∞1）由得：wb wT T T T --=Θww b T T T T +Θ-=)(由可得：T x T x T x T T T x T w b w w b ∂∂Θ-+∂∂Θ=∂+Θ-∂=∂∂)1(])[(rT r T T r T T T r T w w b w w b ∂∂Θ-+∂Θ∂-=∂+Θ-∂=∂∂)1()(])[(由与无关、与无关以及、的表达式可知，除了均匀的情况外，该无量b T r Θx x T ∂∂rT∂∂w T 纲温度定义在一般情况下是不能用分离变量法的；2）由得：∞∞--=ΘT T T T w ∞∞+Θ-=T T T T w )(由可得：T xT x T T T x T w w ∂∂Θ=∂+Θ-∂=∂∂∞∞])[(rT r T T r T T T r T w w w ∂∂Θ+∂Θ∂-=∂+Θ-∂=∂∂∞∞∞)(])[(由与无关、与无关以及、的表达式可知，在常见的四种边界条件中除了b T r Θx x T ∂∂rT ∂∂轴向及周向均匀热流的情况外，有，则该无量纲温度定义是可以用分const q w =0=∂∂rT w离变量法的；3）由得：wwT T T T --=Θ∞ww T T T T +Θ-=∞)(由可得：T xT x T T T x T w w w ∂∂Θ-=∂+Θ-∂=∂∂∞)1(])[(r T T r T T T r T w w w -+∂Θ∂-=∂+Θ-∂=∂∂∞∞1()(])[(同2）分析可知，除了轴向及周向均匀热流const q w =温度定义是可以用分离变量法的；习题4－181）采用柱坐标分析，写出统一的稳态柱坐标形式动量方程：S r r r r r r x x w r v r r r u x +∂∂∂∂+∂∂∂∂+∂∂∂∂=∂∂+∂∂+∂∂(1)(1)()(1)(1)(θφλθφλφλφρθφρφρ、和分别是圆柱坐标的3个坐标轴，、和分别是其对应的速度分量，其中x r θu v w 是管内的流动方向；x 对于管内的层流充分发展有：、，；0=v 0=w 0=∂∂xu并且方向的源项：x x pS ∂∂-=方向的源项：r r pS ∂∂-=方向的源项：θθ∂∂-=pr S 1由以上分析可得到圆柱坐标下的动量方程：方向：x 0)(1)(1=∂∂-∂∂∂∂+∂∂∂∂x pu r r r u r r r θλθλ方向：r 0=∂∂r p 方向：θ0=∂∂θp 边界条件：，R r =0=u ，；对称线上，0=r 0=∂∂r u 0=∂∂θu 不考虑液体的轴向导热，并简化分析可以得到充分发展的能量方程为：)(1(1θλθλρ∂∂∂∂+∂∂∂∂=∂∂Tr r r T r r r x T uc p 边界条件：，；，R r =w q r T =∂∂λ0=r 0=∂∂rT，πθ/0=0=∂∂-θλT2）定义无量纲流速：dxdp R uU 2-=λ并定义无量纲半径：；将无量纲流速和无量纲半径代入方向的动量方程得：R r /=ηx 0))1((1)1((122=∂∂-∂-∂∂∂+∂-∂∂∂xp U dx dp R R R R U dx dp R RR R θληλθηηλληηη上式化简得：011(1(1=+∂∂∂∂+∂∂∂∂θηθηηηηηU U 边界条件：，1=η0=U ，；对称线上，0=η0=∂∂ηU 0=∂∂θU定义无量纲温度：λ/0R q T T b-=Θ其中，是折算到管壁表面上的平均热流密度，即：；0q Rq q wπ=0由无量纲温度定义可得：bT Rq T +Θ=λ0将表达式和无量纲半径代入能量方程得：T η(1)(100θληλθηηλληηηρ∂Θ∂∂∂+∂Θ∂∂∂=∂∂R q R R R R q R R R x T uc b p 化简得：（1）)1(1)(10θηθηηηηηρ∂Θ∂∂∂+∂Θ∂∂∂=∂∂x T u c q R b p 由热平衡条件关系可以得：mm m b m p b p p RU U q R u u R q A u u dx dT A u c x T u c x T uc 020221221)(===∂∂=∂∂ππρρρ将上式代入式（1）可得：)1(1)(12θηθηηηηη∂Θ∂∂∂+∂Θ∂∂∂=m U U 边界条件：，；，0=η0=∂Θ∂η1=ηR q q w πη10==∂Θ∂，；，0=θ0=∂Θ∂θπθ=0=∂Θ∂θ单值条件：由定义可知：且： 0/0=-=ΘλR q T T b b b ⎰⎰Θ=ΘAAb UdAUdA 即得单值性条件：=Θ⎰⎰AA UdAUdA 3）由阻力系数及定义有：f Re 228)(21/Re ⎪⎭⎫ ⎝⎛=⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎣⎡-=D D U D u u dx dp D f e m e m me νρ且：m W b m W b m W R q T T D T T q Nu ,0,,0~2)/(2Θ=-=-=λλ5－21．一维稳态无源项的对流－扩散方程如下所示：（取常物性）xx u 22∂∂Γ=∂∂φφρ边界条件如下：LL x x φφφφ====,;,00上述方程的精确解如下： 11)/(00--=--⋅Pe L x Pe L e e φφφφΓ=/uL Pe ρ2．将分成20等份，所以有：L ∆=P Pe 20 1 2 3 4 5 6……………………… 17 18 19 20 21对于中心差分、一阶迎风、混合格式和QUICK 格式分别分析如下：1）中心差分中间节点： 2)5.01()5.01(11-∆+∆++-=i i i P P φφφ20,2 =i 2）一阶迎风中间节点： ∆-∆++++=P P i i i 2)1(11φφφ20,2 =i 3）混合格式当时，中间节点： 1=∆P 2)5.01()5.01(11-∆+∆++-=i i i P P φφφ 20,2 =i 当时，中间节点： 10,5=∆P 1-=i i φφ20,2 =i 4）QUICK 格式*12111)35(8122121⎥⎦⎤⎢⎣⎡---++++++=+--∆∆-∆∆+∆i i i i i i i P P P P P φφφφφφφ2≠i*1111)336(8122121⎥⎦⎤⎢⎣⎡--++++++=+-∆∆-∆∆+∆i i i i i i P P P P P φφφφφφ2=i 数值计算结果与精确解的计算程序如下：%except for HS, any other scheme doesnt take Pe<0 into consideration %expression of exact solutiony=dsolve('a*b*Dy=c*D2y','y(0)=y0,y(L)=yL','x')y=subs(y,'L*a*b/c','t')y=simple(subs(y,'a*b/c*x','t*X'));ysim=simple(sym(strcat('(',char(y),'-y0)','/(yL-y0)')))y=sym(strcat('(',char(ysim),')*(yL-y0)','+y0'))% in the case of Pe=0y1=dsolve('D2y=0','y(0)=y0,y(L)=yL','x')y1=subs(y1,'-(y0-yL)/L*x','(-y0+yL)*X')%grid Pe number tt=[1 5 10];%dimensionless length m=20;%mdim is the number of inner node mdim=m-1;X=linspace(0,1,m+1);%initial value of variable during calculation y0=1;yL=2;%cal exact solution for n=1:size(tt,2) t=m*tt(1,n); if t==0 yval1(n,:)=eval(y1); else yval1(n,:)=eval(y); end end%extra treatment because max number in MATLAB is 10^308if max(isnan(yval1(:))) yval1=yval1'; yval1=yval1(:);indexf=find(isnan(yval1)); for n=1:size(indexf,1) if rem(indexf(n,1),size(X,2))==0 yval1(indexf(n),1)=yL; else yval1(indexf(n),1)=y0; endendyval1=reshape(yval1,size(X,2),size(yval1,1)/size(X,2));yval1=yval1';end%CD solutiond=zeros(size(tt,2),mdim);a=repmat([1],size(tt,2),mdim);for n=1:size(tt,2)t=tt(1,n);b(n,:)=repmat([0.5*(1-0.5*t)],1,mdim);c(n,:)=repmat([0.5*(1+0.5*t)],1,mdim);d(n,1)=0.5*(1+0.5*tt(1,n))*y0;d(n,mdim)=0.5*(1-0.5*tt(1,n))*yL;endc(:,1)=0;b(:,mdim)=0;%numerical cal by using TDMA subfuctionyval2=TDMA(a,b,c,d,mdim);yval2=[repmat([1],size(tt,2),1),yval2,repmat([2],size(tt,2),1)]; Fig(1,X,yval1,yval2,tt);title('CD Vs. Exact Solution')% FUS solutiond=zeros(size(tt,2),mdim);a=repmat([1],size(tt,2),mdim);for n=1:size(tt,2)t=tt(1,n);b(n,:)=repmat([1/(2+t)],1,mdim);c(n,:)=repmat([(1+t)/(2+t)],1,mdim);d(n,1)=(1+tt(1,n))/(2+tt(1,n))*y0;d(n,mdim)=1/(2+tt(1,n))*yL;endc(:,1)=0;b(:,mdim)=0;%numerical cal by using TDMA subfuctionyval3=TDMA(a,b,c,d,mdim);yval3=[repmat([1],size(tt,2),1),yval3,repmat([2],size(tt,2),1)]; Fig(2,X,yval1,yval3,tt);title('FUS Vs. Exact Solution')% HS solutiond=zeros(size(tt,2),mdim);a=repmat([1],size(tt,2),mdim);for n=1:size(tt,2)t=tt(1,n);if t>2b(n,:)=repmat([0],1,mdim);c(n,:)=repmat([1],1,mdim);d(n,1)=y0;elseif t<-2b(n,:)=repmat([1],1,mdim);c(n,:)=repmat([0],1,mdim);d(n,mdim)=yL;elseb(n,:)=repmat([0.5*(1-0.5*t)],1,mdim);c(n,:)=repmat([0.5*(1+0.5*t)],1,mdim);d(n,1)=0.5*(1+0.5*t)*y0;d(n,mdim)=0.5*(1-0.5*t)*yL;endendc(:,1)=0;b(:,mdim)=0;% numerical cal by using TDMA subfuctionyval4=TDMA(a,b,c,d,mdim);yval4=[repmat([1],size(tt,2),1),yval4,repmat([2],size(tt,2),1)]; Fig(3,X,yval1,yval4,tt);title('HS Vs. Exact Solution')%QUICK Solutiond=zeros(size(tt,2),mdim);a=repmat([1],size(tt,2),mdim);for n=1:size(tt,2)t=tt(1,n);b(n,:)=repmat([1/(2+t)],1,mdim);c(n,:)=repmat([(1+t)/(2+t)],1,mdim);d(n,1)=(1+tt(1,n))/(2+tt(1,n))*y0;d(n,mdim)=1/(2+tt(1,n))*yL;endc(:,1)=0;b(:,mdim)=0;%numerical cal by using TDMA subfuctionyval5=zeros(size(tt,2),mdim);yval5com=yval5+1;counter=1;%iterativewhile max(max(abs(yval5-yval5com)))>10^-10if counter==1yval5com=TDMA(a,b,c,d,mdim);endfor nn=1:size(tt,2)for nnn=1:mdimif nnn==1d(nn,nnn)=((6*yval5com(nn,nnn)-3*y0-3*yval5com(nn,nnn+1))*tt(1,nn))/(8*(2+tt(1,nn)))+((1+tt(1,nn))/(2+tt(1,nn))*y0);elseif nnn==2d(nn,nnn)=((5*yval5com(nn,nnn)-3*yval5com(nn,nnn+1)-yval5com(nn,nnn-1)-y0)*tt(1,nn))/(8*(2+tt(1,nn)));elseif nnn==mdimd(nn,nnn)=((5*yval5com(nn,nnn)-3*yL-yval5com(nn,nnn-1)-yval5com(nn,nnn-2))*tt(1,nn))/(8*(2+tt(1,nn)))+(1/(2+tt(1,nn))*yL);elsed(nn,nnn)=((5*yval5com(nn,nnn)-3*yval5com(nn,nnn+1)-yval5com(nn,nnn-1)-yval5com(nn,nnn-2))*tt(1,nn))/(8*(2+tt(1,nn)));endendendyval5=TDMA(a,b,c,d,mdim);temp=yval5;yval5=yval5com;yval5com=temp;counter=counter+1;endyval5=yval5com;yval5=[repmat([1],size(tt,2),1),yval5,repmat([2],size(tt,2),1)];Fig(4,X,yval1,yval5,tt);title('QUICK Vs. Exact Solution')%-------------TDMA SubFunction------------------function y=TDMA(a,b,c,d,mdim)%form a b c d resolve yval2 by using TDMA%eliminationp(:,1)=b(:,1)./a(:,1);q(:,1)=d(:,1)./a(:,1);for n=2:mdimp(:,n)=b(:,n)./(a(:,n)-c(:,n).*p(:,n-1));q(:,n)=(d(:,n)+c(:,n).*q(:,n-1))./(a(:,n)-c(:,n).*p(:,n-1));end%iterativey(:,mdim)=q(:,mdim);for n=(mdim-1):-1:1y(:,n)=p(:,n).*y(:,n+1)+q(:,n);end%-------------ResultCom SubFunction------------------function y=ResultCom (a,b,c)for n=1:max(size(c,2))y(2*n-1,:)=a(n,:);y(2*n,:)=b(n,:);end%-------------Fig SubFunction------------------function y=Fig(n,a,b,c,d)figure(n);plot(a,b);hold onplot(a,c,'*');str='''legend(';for n=1:size(d,2)if n==size(d,2)str=strcat(str,'''''Pe=',num2str(d(1,n)),''''')''');elsestr=strcat(str,'''''Pe=',num2str(d(1,n)),''''',');endendeval(eval(str));a n d A l l t h i n g s i n t h ei r b e i n g a r e g 13精确解与数值解的对比图，其中边界条件给定，。

第四版传热学第四章习题解答

第四章复习题1、试简要说明对导热问题进行有限差分数值计算的基本思想与步骤。

2、试说明用热平衡法建立节点温度离散方程的基本思想。

3、推导导热微分方程的步骤和过程与用热平衡法建立节点温度离散方程的过程十分相似，为什么前者得到的是精确描述，而后者解出的确实近似解。

4、第三类边界条件边界节点的离散那方程，也可用将第三类边界条件表达式中的一阶导数用差分公式表示来建立。

试比较这样建立起来的离散方程与用热平衡建立起来的离散方程的异同与优劣。

5．对绝热边界条件的数值处理本章采用了哪些方法？试分析比较之．6．什么是非稳态导热问题的显示格式？什么是显示格式计算中的稳定性问题？7．用高斯－塞德尔迭代法求解代数方程时是否一定可以得到收敛德解？不能得出收敛的解时是否因为初场的假设不合适而造成？8．有人对一阶导数()()()221,253x t t t x t i n i n i n in ∆-+-≈∂∂++你能否判断这一表达式是否正确，为什么？一般性数值计算4-1、采用计算机进行数值计算不仅是求解偏微分方程的有力工具，而且对一些复杂的经验公式及用无穷级数表示的分析解，也常用计算机来获得数值结果。

试用数值方法对Bi=0.1,1,10的三种情况计算下列特征方程的根:)6,2,1( =n n μ3,2,1,tan ==n Binn μμ并用计算机查明，当2.02≥=δτa Fo 时用式（3-19）表示的级数的第一项代替整个级数（计算中用前六项之和来替代）可能引起的误差。

Bi n n =μμtanFo=0.2及0.24时计算结果的对比列于下表：第一项的值 0.94513 0.61108 0.10935 前六项的值 0.94688 0.6198 0.11117 比值 0.998140.986940.98364 Bi=0.1 Bi=1 Bi=10 第一项的值 0.99277 0.93698 0.77311 前六项和的值0.99101 0.92791 0.76851 比值1.001771.009781.005984-2、试用数值计算证实，对方程组⎪⎭⎪⎬⎫⎪⎩⎪⎨⎧=++=++=-+5223122321321321x x x x x x x x x用高斯-赛德尔迭代法求解，其结果是发散的，并分析其原因。

传热学杨世铭陶文铨第四章热传导问题数值解法

(m+1,n)
y
(m,n-1)
y
x
o
第四章导热问题的数值解法
x
x
15
以二维、稳态、有内热源的导热问题为例此时：
Φ上 Φ下 Φ左＋Φ右 Φv 0
dt dt 左 A y dx dx
可见：当温度场还没有求出来之前，我们并不知道 dt dx 所以，必须假设相邻节点间的温度分布形式，这里我们
无内热源时变为：
x 2

Φ
4tm,n tm1,n tm1,n tm,n1 tm,n1
重要说明：所求节点的温度前的系数一定等于其他所有相邻节点温度前的系数之和。这一结论也适用于边界节点。但这里不包括热流(或热流密度)前的
系数。
第四章导热问题的数值解法
19
4-2 边界节点离散方程的建立及代数方程的求解
第四章导热问题的数值解法
4
§4-1 导热问题数值求解的基本思想及内部节点离散方程的建立
1 物理问题的数值求解过程
建立控制方程及定解条件确定节点（区域离散化）
设立温度场的迭代初值
建立节点物理量的代数方程
求解代数方程
改进初场
是否收敛是解的分析
否
第四章导热问题的数值解法
5
针对二维矩形区域内的稳态无内热源导热问题外部与温度为t的流体对流换热换热系数为h请建立外部角点的温度离散方程并化简到最后的形式inclassproblems第四章导热问题的数值解法26quickreviewquickreview1导热数值解法的重要意义2导热数值解法的基本思想3网格划分区域离散的过程及涉及的基本概念4代数方程离散方程的建立方法和过程第四章导热问题的数值解法2742代数方程组的求解代数方程组的求解一维无限大平板稳态常物性无内热源左侧第一类边条右侧第三类如右图所示将其均匀分成三个控制体试建立离散方程内部节点内部节点边界节点第四章导热问题的数值解法28形成如下代数方程组

数值传热学第四章课件陶文铨

主讲陶文铨西安交通大学能源与动力工程学院热流中心CFD-NHT-EHT CENTER 2010年9月27日, 西安数值传热学第四章扩散方程的数值解及其应用(1)4.1 一维导热问题4.1.1一维稳态导热的通用控制方程4.1.3界面导热系数的确定方法4.1.4 一维非稳态导热控制方程的离散化4.1.2通用控制方程控制容积积分法的离散4.1.5 数学上的稳定未必导致物理上有意义的解一维稳态导热问题不同坐标系通用控制方程0 P P()0P x x Δ=i调和平均已经广泛为国内外学术界所接受。

≤1数学上的稳定未必导致物理上有意义的解无内热源一维非稳态导热，初场均匀，两表面0]T +代入下式：P（全隐格式）才能满足。

结论：数学上的稳定未必导致物理上有意义的解；推=xΔa TP P极坐标均可以表示成为：2.解决通用化的一种方案为写出适合于三种坐标系中系数的通用表达式，特引进两个辅助变量：（1）x –方向标尺因子，scaling factor ，x-方向的距离表示成为sx x δi 。

对直角、圆柱坐标规定1;sx ≡（2）y-方向引入一个名义半径，R 。

对直角坐标R ＝1，据此，东西导热距离为：sx xδi 东西导热面积为：R /y sxΔ对极坐标取;sx r =对圆柱与极坐标R ＝r三种二维正交坐标系中离散方程的统一表达式按这种方式编制程序时，只要设置一个变量MODE，4.3 源项与边界条件的处理4.3.1非常数源项的线性化处理1. 线性化方法4.3.2第二、三类边界条件使方程组封闭的处理2. 线性化方法讨论3. 线性化方法应用实例1. 补充以边界节点代数方程的方法2. 附加源项法S= P2. 线性化方法讨论（1）对与被求解变量有关的非常数源项，线性化比假定为常数更合理：用*()PS f T =来表示P 的源项比落后一个迭代步；P C P T S S S =+（2）任何复杂的函数总可以用线性函数来近似逼近；线性又是建立线性代数方程所必须的；（3）是为保证代数方程迭代求解收敛所必须；0P S ≤P P nb nb a a b φφ=+∑P nb a a ≥∑P nb P a a S V =−Δ∑代数方程迭代求解收敛的充分条件是，因为可以确保代数方程迭代求解收敛。

传热学第四章

k及k+1表示迭代次数；
t
(k) max
—第k次迭代得到的最大值
华北电力大学
刘彦丰
传热学 Heat Transfer
4-3 一维非稳态导热问题的数值求解
在非稳态导热问题中，不但需要对空间区域进行离散，还需要对时间变量进行离散，接下来以一个一维非稳态导热问题为例，重点介绍对非稳态项的离散方法，以及不同离散方法对计算带来的影响等。
n
y
m
x
华北电力大学
刘彦丰
传热学 Heat Transfer 2.建立节点物理量的代数方程
每一个节点都与它周围相邻的节点存在一定的关系，通过相应的物理定律，可建立它们之间的关系式（属于代数方程式），此关系式又称作节点的离散方程。
(m,n+1)
(m,n)
(m-1,n)
(m+1,n)
华北电力大学
(m,n-1)
h,tf
0 q=0 H x
华北电力大学
刘彦丰
传热学 Heat Transfer
二、数学描述
华北电力大学
∂2t ∂x 2
+
∂2t ∂y 2
+
Φ& λ
=
0
x=0 x=H y=0
−
∂t ∂x
=
0
−
λ
∂t ∂x
=
h(t
−tf
)
−
∂t ∂y
=
0
y =W
−
λ
∂t ∂y
=
h(t
−tf
)
刘彦丰
传热学 Heat Transfer
tm,n+1 −tm,n ∆y
(m,n+1)

Chapter_4(2)_

主讲陶文铨西安交通大学能源与动力工程学院热流中心CFD-NHT-EHT CENTER 2010年9月27日, 西安数值传热学第四章扩散方程的数值解及其应用(2)4.4 求解代数方程的TDMA 及ADI 方法4.4.1求解一维导热问题代数方程的三对角阵算法4.4.2求解多维非稳态导热全隐格式的ADI 方法1.求解方法概述1.一维导热问题代数方程通用形式2.Peaceman-Rachford 的ADI 迭代2.Thomas 算法3.第一类边界条件的处理4.4 方法4.4.1 求解一维导热问题代数方程的三对角阵算法Thomas算法的一般形式将上式改写为：所谓消元就是要找出系数间的关系。

(b)乘以CiP=11;i−−1,i=B D3.第一类边界条件下Thomas 算法的实施第一类边界条件下，求解区域为i=2,….M1-1=M2。

将消元公式用于i=1, 注意T 1是给定的：1121T PT Q =+10;P =11Q T =因T M1已知，消元从M 2开始：2212M M M T P T Q =+注意：采用附加源项法来处理第二类，第三类边界条件时，均将第二类，第三类边界条件问题视为第一类边界条件问题，数学上的处理与此相同。

4.4.2 求解多维导热问题代数方程的方法求解二维非稳态导热全隐格式代数方程的方法(1) 五对角阵算法(Penta-diagonal ,PDMA)(2) 交替方向隐式方法(Alternative-directionImplicit, ADI)2. 3-D Peaceman-Rachford方法将tΔ三等分：tΔX方向为隐式，第一个/3y,z方向为显式式；第二，三个/3tΔ分别在y,z方向实施隐式;2－D交替方向隐式u i,j,k , v 表示方向二阶导数的中心差分；nT用von Neumann分析方法可以证明稳定性条件为：4.5 管道内充分发展对流换热概说4.5.1管道内充分发展对流换热的定义1. 简单的充分发展对流换热2. 复杂的充分发展对流换热4.5.2能实现充分发展对流换热的边界条件4.5.3部分算例汇总4.5 管道内充分发展对流换热概说4.5.1管道内充分发展对流换热的定义平直通道中的充分发展对流换热属于这一类。