数值传热_数值传热学大作业3gg

合集下载

数值传热课后学习题及答案(权威版)

数值传热课后学习题及答案(权威版)

习题 4-14
充分发展区的温度控制方程如下:
c
p
u
T x
1 r
r
(r
T r
)
对于三种无量纲定义 T Tw 、 T T 、 T Tw 进行分析如下
Tb Tw
Tw T
T Tw
1)由 T Tw 得: Tb Tw
T (Tb Tw ) Tw
由 T 可得:
T [(Tb Tw ) Tw ] Tb (1 ) Tw
x
x
x
x
T r
[(Tb
Tw ) Tw ] r
(Tb
Tw
)
r
(1 ) Tw r
由 Tb 与 r
无关、
与 x 无关以及
T x

T r
的表达式可知,除了 Tw 均匀的情况外,该无量
纲温度定义在一般情况下是不能用分离变量法的;
2)由 T T 得: Tw T
T (Tw T ) T
由 T 可得:
节点 3:
T2 4T3 75
求解结果:
T2 85 , T3 40
对整个控制容积作能量平衡,有:
qB Sx h f (T f T3 ) Sx 15 (20 40) 150 2 0
即:计算区域总体守恒要求满足
习题 4-5
在 4-2 习题中,如果 h 10 (T3 T f )0.25 ,则各节点离散方程如下:
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);

数值传热大作业

数值传热大作业

放置竖直孤立平板的二维围场内的空气流动与换热的数值分析(西安交通大学能源与动力工程学院,710049,西安)摘要:针对内部放置孤立平板的二维围场内的空气流动与换热问题,在稳态、常物性和壁面温度以及孤立平板温度恒定的条件下,采用SIMPLER算法,对围场内部的空气进行了流动与换热的数值模拟计算。

在瑞利数Ra=10000时,计算得到了二维围场内的流线、等温线以及热线。

关键词:SIMPLER算法、孤立平板、流线、等温线、热线Abstract:Inorder to investigate the fluxion and heat transfer of air in a 2D square enclosure with an isolated plate. SIMPLER algorithm was adopted based on the Reylonds conservation equations of the steady-state constant property laminar flow and a constant temperature of the isolated plate and the inner walls of the enclosure condition. Slove fluid velocity and temperature fields inthe enclosure for Ra=10000,and draw the diagrams of stream lines ,isotherms and heat lines.Key words:SIMPLER algorithm; isolated plate; stream lines; isotherms; heat lines主要符号表R瑞利数aP普朗特数rν空气运动粘度m2/sg重力加速度kg.m/s2 k空气导热系数W/(m℃) β空气体膨胀系数1/℃c空气比热容J/kg.℃pρ空气密度kg/(m3s)T金属板温度℃hT围场壁面温度℃c∆温差℃T一、引言封闭空腔内孤立物体自然对流换热是一个重要的研究课题,从某种角度讲,大空间自然对流是封闭腔内孤立物体自然对流的一个特例。

传热与流体数值计算大作业

传热与流体数值计算大作业

大作业一、假设0,1x y≤≤的方腔内充满不可压缩流体,左、右、下壁面固定,上壁面以()22161u x x=--运动。

试求腔内的定常解。

(流体的物性取20℃的水。

同时,可以使用20℃的甘油作为对比)二、求解二维圆柱坐标中的Poisson-Nernst-Plack(PNP)方程,PNP方程来描述纳米孔内带电离子在浓度梯度及电场作用下的迁移行为和离子浓度分布。

具体方程如下所示:其中i=+/-,分别代表阴阳离子。

以及连续性方程:其中Φ是局域的电动势,c i表示i种离子的浓度,左侧边界上c+=10,c-=10,右侧边界上c+=1,c-=1。

j i表示离子流,D i为离子的扩散系数2×10-9,z i为离子的带电量,z i=1,T为溶液的温度,T=300。

e是电子电量1.602×10-19,ε0×εr=80,k B为波尔兹曼常数,k B=1.38×10-23。

边界上的电势Φ由高斯定律决定:对于带电的纳米孔壁(图中红色实线所示),有σs=σ(σ为纳米孔的表面电荷密度,数值为0.05);对于其余区域有σs=0。

离子流j i在边界上的法向分量为零,即,求解φ、浓度c i以及ij的场。

(备注:求解区域为一圆柱形区域,长度为1200,直径为d=10。

建议步骤:可首先猜想浓度场c+和c-,并求解电动势场φ,通过连续性方程修正离子流场ij)大作业要求:1-3人为一组,完成以上任选一题目。

最终截止时间为12月26日。

在最终截止时间之前可以提交1次,若不满意得分可以继续修改。

大作业以报告形式提交,内容至少包括计算域的网格划分、方程的离散化、边界条件的处理、计算收敛的判据、计算的结果、结果的图形化显式、结果分析等。

源代码作为附录附在报告的最后。

数值传热_数值传热学大作业3gg

数值传热_数值传热学大作业3gg

数值传热学 2009-2010 学年第一学期大作业 3
阶 梯 形 标 量 场 ( 长 宽 均 为 1 ) 的 纯 对 流 传 递 ( θ = 340 ), 控 制 方 程 为 u ∂φ + v ∂φ = 0 ,上游的边界条件都是第一类的,即给定了φ 的分布,下游按开
∂x ∂y 口边界的方式处理。
图 1 示意图
(2)从图中可以得知:数值计算在剧烈变化区域(y=0.5 处),采用 CD、 SUD、QUICK 格式时,产生越界现象。
(3)计算过程中采用 STOIC 格式,计算效果最好。 (4)采用不同的格式时,其表达式在 CBC 线内,则会出现稳定解,否则会 出现解得越界现象,越界现象与对流稳定性不同。 (5)编程过程中,采用 SOR 低松弛迭代方法,所得结果比较理想,计算速 度远远 Gauss—Seidel 方法。在本次作业中,松弛系数取为 0.8。 (6)编程过程中,要将边值点带入循环进行计算,否则会导致计算结果出 错。进而也证明了,对流现象是具有方向性,其扰动只能沿下游传播。
要求:
1、分别利用 FUD,CD,SUD,QUICK,CLAM,EULER,MINMOD,MUSCL,OSHER,SMART, STOIC 来离散对流项,观察它们的计算结果有何不同。 2、用延迟修正进行求解。 3、写出详细的离散过程和求解方法。 4、用 TECPLOT 软件画出整个流场中φ 的分布,并画出 y = 0.5 时,φ 随 x 的 分布。 5、编程采用 C/C++或 FORTRAN 语言。 6、将源程序附于作业之后,程序中要有详细的注释,以反映出思路。 7、源程序电子版和打印版上交时间:截止 2009 年 12 月 28 日。 8、每组交一份作业,给出每位同学的贡献度。(总和为 100%)

数值传热学--作业

数值传热学--作业

数值传热学大作业—淬火过程的瞬态热分析专业:材料工程班级:研1303班学号:S2*******指导教师:孙斌煜姓名:李康一、问题描述某零件材料为45钢,按照国标GB/T6912-1999规定的45钢推荐热处理制度为840C 。

淬火.600C 。

回火,淬火介质为水,试计算零件温度随时间的我变化曲线和最后时刻的温度场云图 (1)45钢弹性模量:200GPa 泊松比:0.3质量密度:78503/m kg膨胀系数:15.5e-6m/C 。

比热:448C J kg / 导热系数:70()C m W .*/ (2)水 密度:9963/m kg 比热:4185C J kg / 导热系数:2()C m W .*/水沸腾对流换热系数:1200()C m W .*2/初始45钢温度840,水的初始温度为20C 。

,水槽宽1m,中间位零件最大截面60mm ,下图为淬火过程的零件截面。

图-1二、创建模型1.建立分析项目(1)在Windows系统下执行“开始”—“所有程序”—“ANSYS14.0”—“Mechanical APDL(ANSYS)14.0”命令,启动Mechanical APDL(ANSYS)14.0,进入主界面。

(2)选择热分析过滤菜单GUI:选择菜单Main Menu —Preprocessor,弹出分析项目对话框,选择Thermal 热分析,如图2 所示,完成后单击OK按钮结束。

2.更改分析名称和标题(1)改变工作项目标题GUI:File→Change Title,弹出对话框,输入“Thermal01”如下图,单击OK结束。

(2)更改项目名称GUI:File→Change Title,弹出对话框,输入“Thermal01”下方的复选框,如下图所示,完成单击OK完成。

3.创建材料模型要点:创建模型顺序依次为工件,水(1)添加导热系数GUI:Main Menu →Preprocessor→Material Prop→Material Models→Thermal→Conductivity→Isotropic,弹出对话框,输入导热系数70,如下图完成后单击OK 结束输入(2)添加比热容GUI:Main Menu →Preprocessor→Material Prop→Material Models→Thermal→Specific Heat.弹出比热容输入对话框,在文本框中输入工件比热容448,如下图,成后单击OK结束输入(3)添加密度GUI:Main Menu →Preprocessor→Material Prop→Material Models→Thermal→Density,弹出密度输入对话框,在文本框中输入工件比热容7785,如下图,成后单击OK结束输入(4)创建材料2依照上述步骤添加水的比热容4185,密度996,导热系数25.选择单元GUI:Main Menu →Preprocessor→Element type→add/Edit/Delete→add,弹出下图所示单元对话框,选择Thermal Solide的Quad 8node77单元,按OK键结束设置单元选项GUI:Main Menu →Preprocessor→Element type→add/Edit/Delete→add,弹出Element type对话框,单击对话框中的Option,弹出设置单元对话框,在单元形状K3文本框选择Plane Thickness,如下图,单击OK结束关闭对话框。

传热大作业-数值解法-清华-传热学

传热大作业-数值解法-清华-传热学

一维非稳态导热的数值解法一、导热问题数值解法的认识(一)背景所谓求解导热问题,就是对导热微分方程在规定的定解条件下的积分求解。

这样获得的解称为分析解。

近100年来,对大量几何形状及边界条件比较简单的问题获得了分析解。

但是,对于工程技术中遇到的许多几何形状或边界条件复杂的导热问题,由于数学上的困难目前还无法得出其分析解。

另一方面,在近几十年中,随着计算机技术的迅速发展,对物理问题进行离散求解的数值方法发展十分迅速,并得到日益广泛的应用。

这些数值方法包括有限差分法、有限元法及边界元法等。

其中,有限差分法物理概念明确,实施方法简便,本次大作业即采用有限差分法。

(二)基本思想把原来在时间、空间坐标系中连续的物理量的场,如导热物体的温度场,用有限个离散点上的值的集合来代替,将连续物理量场的求解问题转化为各离散点物理量的求解问题,将微分方程的求解问题转化为离散点被求物理量的代数方程的求解问题。

(三)基本步骤(1)建立控制方程及定解条件。

根据具体的物理模型,建立符合条件的导热微分方程和边界条件。

(2)区域离散化。

用一系列与坐标轴平行的网格线把求解区域划分成许多子区域,以网格线的交点作为需要确定温度值的空间位置,称为节点。

每一个节点都可以看成是以它为中心的一个小区域的代表,将小区域称之为元体。

(3)建立节点物理量的代数方程。

建立方法主要包括泰勒级数展开法和热平衡法。

(4)设立迭代初场。

(5)求解代数方程组。

(6)解的分析。

对于数值计算所获得的温度场及所需的一些其他物理量应作仔细分析,以获得定性或定量上的一些结论。

对于不符合实际情况的应作修正。

二、问题及求解(一)题目一厚度为0.1m 的无限大平壁,两侧均为对流换热边界条件,初始时两侧流体温度与壁内温度一致,1205f f t t t ===℃;已知两侧对流换热系数分别为h 1=11 W/m 2K 、h 2=23W/m 2K ,壁的导热系数λ=0.43W/mK ,导温系数a=0.3437×10-6 m 2/s 。

数值传热学 习题答案

数值传热学习题答案数值传热学习题答案数值传热学是热力学的一个重要分支,主要研究热量在物质中传递的机理和规律。

在实际工程中,我们经常会遇到各种与传热有关的问题,通过数值计算可以得到准确的答案。

下面我将为大家提供一些数值传热学习题的答案,希望能够帮助大家更好地理解和应用这门学科。

1. 一个铝制热交换器的表面积为10平方米,其表面温度为100摄氏度,环境温度为20摄氏度。

已知铝的导热系数为200 W/(m·K),求热交换器的传热速率。

答:根据传热定律,传热速率与传热面积、传热系数和温度差之间成正比。

传热速率 = 传热系数× 传热面积× 温度差。

将已知数据代入公式中,可得传热速率= 200 × 10 × (100 - 20) = 160,000 W。

2. 一个房间的尺寸为5米× 5米× 3米,墙壁和天花板的厚度为0.2米,墙壁和天花板的导热系数为0.5 W/(m·K),室内温度为25摄氏度,室外温度为10摄氏度。

求房间的传热损失。

答:房间的传热损失可以通过计算墙壁和天花板的传热速率来得到。

墙壁和天花板的传热速率 = 传热系数× 传热面积× 温度差。

墙壁和天花板的传热面积 = 2 × (5 × 5) + 2 × (5 × 3) = 70平方米。

将已知数据代入公式中,可得墙壁和天花板的传热速率= 0.5 × 70 × (25 - 10) = 525 W。

因此,房间的传热损失为525瓦特。

3. 一个水箱的体积为1立方米,初始温度为20摄氏度,水的密度为1000千克/立方米,比热容为4186 J/(千克·摄氏度),水箱的表面积为2平方米,表面温度为100摄氏度。

已知水的传热系数为0.6 W/(m^2·K),求水箱内水的温度随时间的变化。

数值传热学作业-第四章

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 xT T S δδλλ-=++V所以 434111. 1.36311T T T =++即 43122293T T -=采用区域离散方法B22d TT=0dx - 由控制容积法 0w edT dT T x dT dT ⎛⎫⎛⎫--∆= ⎪ ⎪⎝⎭⎝⎭ 所以代入2点4点有322121011336T T T T T ----= 即 239028T T -=544431011363T T T T T ----= 即3459902828T T T -+= 对3点采用中心差分有432322+T 013T T T --=⎛⎫⎪⎝⎭即2349901919T T T -+= 对于点5 由x=11dT dx =,得 5416T T -= (1)精确解求左端点的热流密度由 ()21x x eT e e e -=-+所以有 ()2220.64806911x xx x dT e e q e e dxe e λ-====-+=-=++ (2)由A 的一阶截差公式210.247730.743113x T T dT q dxλ=-=-==⨯= (3)由B 的一阶截差公式0.216400.649213x dTq dxλ=-=-== (4)由区域离散方法B 中的一阶截差公式:210.108460.6504()B BT T dT dx x δ-⎛⎫==⨯= ⎪⎝⎭ 通过对上述计算结果进行比较可得:区域离散B 有控制容积平衡法建立的离散方程与区域离散方程A 中具有二阶精度的格式精确度相当!4-3解: 对平板最如下处理:1 2 3 4由左向右点分别表述为1、2、3、4点,x 的正方向为由左向右; 控制方程为λd 2tdx +S =0 (1)边界条件为X=0,T=75℃;X=0.1,λdTdx +ℎ(T −T f )=0;则2、3点采用二阶截差格式,有 则有以下两式:λT3−2T2+T1∆x+S=0(2)λT4−2T3+T2∆x2+S=0(3)一阶截差公式可由λdTdx+ℎ(T−T f)=0变形得到λ(T4−T3∆x)=h(T4−T f)再变形得到T4=[T3+h×∆xλT f]/(1+h×∆xλ)(4)二阶截差公式可以联立λT5−2T4+T3∆x2+S=0和λ(T5−T32∆x)=h(T4−T f),可得以下公式T4=[T3+∆x2S2λ+h×∆xλ]/(1+h×∆xλ)(5)分别联立2、3、4式与2、3、5式,把S=50×103W/m3,λ=10W/m∙℃,h=50 W/m∙℃,T f=25℃,T1=75℃,∆x= 1/30带入到式子中,则有联立2、3、4式的解为:T2=78.58℃,T3=76.59℃,T4=69.03℃联立3、4、5式的解为:T2=80.42℃,T3=80.28℃,T4=74.58℃对控制方程进行积分,并将边界条件带入,则有关于T的方程T=−2500x2+250x+75(6)把x2=130,x3=230,x3=0.1代入上述6式则有:T2=80.56℃,T3=80.56℃,T4=75.1℃相比之下,对右端点采用二阶截差的离散更接近真实值4-4解:对平板作如下分析:1 2 3 4 5 由左向右分别对点编号为1、2、3、4、5 控制方程与4-3相同,为λd 2tdx +S =0 (1)边界条件为X=0,T=75℃;X=0.1,λdTdx +ℎ(T −T f )=0;设1点和2点的距离为∆x ,另1点对2点进行泰勒展开,有d 2t dx =(T 1−T 2+dT dx ∆x )2∆x其中dT dx=T 3−T 22∆x,则有λ2T 1−3T 2+T 3∆x 2+S =0 (2)对3点进行离散有λT 4−2T 3+T 2∆x 2+S =0 (3)对右端点有: [a p +A 1ℎ+(δx )5λ]T 4=a w T 3+[S/∆x +AT f 1ℎ+(δx )5λ]代入数据有T 3−3T 2+155.56=0 T 4−2T 3+T 2=−5.56342.85T4-300T3=1681解得:T2=78.1℃,T3=78.7℃,T4=73.8℃由导热定律有T4−T3∆x =2T5−T4∆x则有T5=71.35℃4—12编写程序:M=rand(10,3)A=M(:,1);B=M(:,2);C=M(:,3);B(10)=0;C(1)=0;T=12:21;D(1)=A(1)*T(1)-B(1)*T(2)for i=2:9;D(i)= A(i)*T(i)-B(i)*T(i+1)-C(i)*T(i-1)endD(10)= A(10)*T(10)-C(10)*T(9);P(1)=B(1)/A(1);Q(1)= D(1)/A(1);for i=2:10;P(i)=B(i)/(A(i)-C(i)*P(i-1));Q(i)=(D(i)+C(i)*Q(i-1))/(A(i)-C(i)*P(i-1)); endfor i=10:-1:2;t(10)=Q(10);t(i-1)=P(i-1)*t(i)+Q(i-1);enddisp(D(1:10))disp(T(1:10))disp(t(1:10))运行结果:由运行结果可知:无论系数怎样变化,T与t都是一致的。

西安交通大学数值传热学大作业


5
数值传热学论文
能量方程: u 边界条件:
T T 2T 2T v ( 2 2 ) x y x y
T ( x, Tp ) T ( x,0)
(4)
(5) (6) (7) (8) (9) (10)
u( x, Tp ) u( x,0)
v( x, T v( x, 0 ) p )
图 2 网格划分
百叶窗翅片通道内周期性充分发展流动与换热的控制方程如下: 连续性方程: 动量方程: u
u v 0 x y
(1) (2) (3)
u u 1 p 2 u 2u v v( 2 2 ) x y x x y u v v 1 p 2v 2v v v( 2 2 ) x y y x y
2
数值传热学论文
主要符号表
f
Nu m
摩擦因数 平均 Nusselt 数 Prandtl(普朗特)数 雷诺数 竖直平板和封闭方腔壁面间的距离,热扩散系数(定义 u r ) 表面换热系数 导热系数 温度 平均温度 内部翅片的温度 W/(m2℃) W/(m℃) ℃ ℃ ℃ m/s m m m Pa W/m2
数值传热学论文
数值传热学大作业
1
数值传热学论文
百叶窗翅片流动换热的数值模拟
(西安交通大学能源与动力工程学院,710049,西安) 摘要: 针对具有一定倾斜角度的流动和换热已经进入周期性充分发展的百叶窗换 热问题,在稳态、层流、常物性和翅片温度恒定的条件下,采用 SIMPLER 算法, 对百叶窗的一个翅片单元进行了数值模拟计算。在翅片倾角θ =25°,雷诺数 Re 在 10 到 500 范围内变化时,得到了平均 Nusselt 数与阻力系数 f 的计算结果。计 算结果表明:随着 Re 的增大,平均 Nusselt 数逐渐增大,f 却随之逐渐减少。 关键词:百叶窗;周期性发展;数值模拟;SIMPLER 算法

数值传热学习题答案(汇总版)


2-4-9
= rP rS
式(2-4-9)也可以写成 a PTP = a E TE + aW TW + b 的形式。而且两种结果是一致的。
2—6:
n n TE −TW dT P , n = 解:将 , dx 2x n n TE −2TPn + TW d 2T P , n = , dx2 x 2
dk = f (x ) 代入原方程,得: dx

2-4-4
rk rk a E = , aW = , a P = a E + aW , b x w x e
= SrP r ,
式(2-4-4)可以写成 a PTP = a E TE + aW TW + b 的形式。 2. 再用 Taylor 展开法导出 k
2 2 uE + uP u = , 2 2 e
2 2 uW + uP u = 2 2 w
t u ut N − uP y = (y ) , n n
t
t ut u p − uS y = (y ) 。 s s
t
(y ) n = (y ) s = y
n n n n TE −TW TE −2TPn + TW k + f (x ) +S=0 整理得: 2x x 2
4kT P= 2k + xf ( x)T E+2k − xf ( x)T W +2x 2 S
− 2k 时, a E 会成为负值, x 2k 当 f(x)> 时, aW 会成为负值。 x
rk dr = rk r r dr dr dr
w
e
1 d
  1. 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
  2. 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
  3. 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。

=
1 2)
(1 +
) φC
)
= φC
0

) φC

1 2
1 2

) φC
≤1
其它范围
)) φf = 2φC
MUSCL:
=
1 4
+
) φC
=1 )
= φC
0

) φC

1 4
1 4

) φC

3 4
3 4

) φC

1
其它范围
) φf
=
3 2
) φC
OSHER: = 1 )
= φC
0

) φC

2 3
五、心得体会
通过此次大作业,我们认识到大家合作的重要性,同时也锻炼了我们共同处
理问题的能力。此次作业中遇到了很困惑的一个问题,就是一开始采用规正变量 后得不到收敛解,原因是计算规正变量时分母会出现为 0 的情况,一开始的做法 是当分母为 0 时在分母上加上一个很小的数,这样做的话或导致分数值很大,导 致最终不收敛。最后我们的做法是,当分母为 0 时,所计算的边界值采用一阶迎 风形式。从这次亲身体会,我们认识到在得到离散的格式之后,不要盲目的去算, 而是应该判断其特点,然后再进行程序的编写,并且,在编程过程中要熟知语句 的用法,这样不仅会节省很多时间,而且容易发现错误。同时,我们学会了解决 问题的一些基本思想,并可以应用这些基本思想解决其它学科中所遇到的问题。
+ w
− φW
)
− ue (
f
+ e
− φP ))∆y
+
(vs (
f
+ s
− φS
)
− vn (
f
+ n
− φP ))∆x
设沿斜线的对流速度为 f,则可知:
ue = uw = f cosθ , vn = vs = f sinθ
其中对
f
+ w
,
f
+ e
,
f
+ n
,
f
+ s
引入规正变量,进行延迟修正处理,且对于不同的格
b
=
(uw (
f
+ w
− φW
)
− ue (
f
+ e
− φP ))∆y
/2+
(vs (
f
+ s
− φS
)−
vn (
f
+ n
− φP ))∆x
特殊处理:
在计算边界处的
f
+ w
,
f
+ s
时按一阶迎风格式近似处理来求解其值。
具体编程处理: 取网格数为 100×100,所以 x、y 方向空间步长分别为 0.01,进行编程计算。
,0
](
f
+ s
−φS
)
−[

(v)s
,0
](
f
− s
−φP
))∆x
若对流方向沿斜线向上,则 u 的方向为从左向右,v 的方向为从下向上。可
将以上各式进行化简,如下所示: aE = 0, aN = 0, aW = uw∆y, aS = vs ∆x,a P = ue∆y + vn∆x
b
=
(uw (
f
一.采用有限容积法对方程进行进行离散
原方程: u ∂φ + v ∂φ = 0
(1)
∂x ∂y
开始离散处理:视 u、v 为常数,对(1)式进行积分:
∫ ∫n
s
e w
u
∂φ dxdy ∂x
=
[(uφ ) e

(uφ) w
]∆y
(2)
N● n
W● w P● e E● s S●
∫ ∫e w
n s
v
∂φ ∂y
引入规正变量,采用延迟修正方法,采用 SOR 进行迭代(w=0.8),当累积误差 err<10-5 时,迭代终止,输出各网格点数据并绘制图形。
二、程序
见附页。
三、y=0.5 时的图形
见附页
四、计算结果分析
(1)计算过程中引入规正变量,采用延迟修正处理,保证编程过程中方程 对角占优,从而易于求得稳定解。
∂x ∂y 口边界的方式处理。
图 1 示意图
要求:
1、分别利用 FUD,CD,SUD,QUICK,CLAM,EULER,MINMOD,MUSCL,OSHER,SMART, STOIC 来离散对流项,观察它们的计算结果有何不同。 2、用延迟修正进行求解。 3、写出详细的离散过程和求解方法。 4、用 TECPLOT 软件画出整个流场中φ 的分布,并画出 y = 0.5 时,φ 随 x 的 分布。 5、编程采用 C/C++或 FORTRAN 语言。 6、将源程序附于作业之后,程序中要有详细的注释,以反映出思路。 7、源程序电子版和打印版上交时间:截止 2009 年 12 月 28 日。 8、每组交一份作业,给出每位同学的贡献度。(总和为 100%)
+
(
f
− w
− φP )]
(5)
(vφ ) n
= [ (v)n ,0 ][φP
+
(
f
+ n
− φP )] − [ − (v)n ,0 ][φN
+
(
f
− n
− φN )]
(6)
(vφ ) s
= [ (v)s ,0 ][φS
+
(
f
+ s
− φS )] − [ − (v)s ,0 ][φP
+
(
f
− s
− φP )]
六、每个小组成员承担工作
工作 姓名
禹国军
编程 王晏林
王岩
绘图 韩立涛
报告、ppt 刘文婷
附页
(FUD)
(SUD)
(CD) (CLAM)
(QUICK) (EULER)
(MINMOD)
(MUSCL)
(OSHER)
(SMART)
(STOIC)
− φC )
= φC
) φ=
φ −φU
φD − φU
0 ≤ ϕ)C ≤ 1 其它范围
f



UC D
) φf = EULER: = 3
4) = φC
)) ) φ(C 1− φC))3 − φC2
1− 2φC
) 0 ≤ φC ≤ 1 ) φC = 0.5 其它范围
) φf
=
3 2
) φC
MINMOD:

1 5
1 5

) φC

1 2
1 2

) φC

5 6
5 6

) φC
≤1
其它范围
由于 ∂φ = 0, ∂φ = 0 ,因此,将上边界与右边界做绝热处理,认为边界处值 ∂y ∂x
等于内点值;在左上边界处φ = 1 ,因此,将边界值可赋值为 1,紧邻边界点的内 点计算时按划分网格的一半计算系数参数来近似处理。 aE = 0, aN = 0, aW = uw∆y, aS = vs ∆x / 2,a P = ue ∆y + vn∆x / 2
将式(4)、(5)、(6)、(7)分别代入式(2)、(3),整理可得:
(([ (u)e ,0 ] + [ − (u)w ,0 ])∆y + ([ (v)n ,0 ] + [ − (v)s ,0 ])∆x)φP
= [ (u)w ,0 ]∆yφW + [ − (u)e ,0 ]∆yφE + [ − (v)n ,0 ]∆xφN + [ (v)s ,0 ]∆xφS
aS = [ (v)s ,0 ]∆x
b
=
([

(u)e ,0 ](
f
− e
−φE
)
−[
(u)e ,0
](
f
+ e
−φP
)
+[
(u)w ,0
](
f
+ w
− φW
)
−[

(u)w ,0
](
f
− w
− φP
))∆y
+
([

(v)n
,0
](
f
+ n
−φP)+ [

(v)n
,0 ](
f
− n
−φN
)
+[
(v)n
dxdy
=
[(vφ
)n
− (vφ)s ]∆x
已知:
(3)
图 1 局部网格
(uφ ) e
=
[ (u)e ,0 ][φP
+
(
f
+ e
− φP )] − [ − (u)e ,0 ][φE
+
(
f
− e
− φE )]
(4)
(uφ) w
= [ (u)w ,0 ][φW
+
(
f
+ w
− φW
)] − [ − (u)w ,0 ][φP
,0
](
f
− s

φP
))∆x
因此可知:
aP = ([ (u)e ,0 ] + [ − (u)w ,0 ])∆y + ([ (v)n ,0 ] + [ − (v) s ,0 ])∆x
相关文档
最新文档