CP090-计算物理热传导方程的差分解法
差分方法

figure
h=plot(x,uu(:,1),'linewidth',5);
set(h,'EraseMode','xor')
axis([0,1,0,0.25]);
fork=2:200
uu(2:99,2)=(1-2*c)*uu(2:99,1)+c*(uu(3:100,1)+uu(1:98,1))-b*dt/dx*(uu(3:100,1)-uu(2:99,1));
plot(u(1,:))
subplot(2,1,2)
plot(u(end,:))
差分方程所得的数值解的图形如图4所示,其中(a)是开始状态,(b)是最后状态。
(a) 初始状态
(b) 最后状态
【程序】
N=500;dx=0.01;dt=0.000001;
c=50*dt/dx/dx;
A=500;b=5;
x=linspace(0,1,100)';
方程设置是parobolic型,系数取为 。
解题的时间范围为 ,初始条件是 。
为了有足够的精度,将初始化的网格作了两次细分。而作图的选项为Contour和Animation。
作为对比,可以更改初始条件为 ,即 。
资料来源:数学物理方程与Matlab可视化.
(a) 整体图
(b) 上图:初始状态,(c)下图:最后状态
图3 解析解的图形
【程序】:
a2=50;b=5;
[x,t]=meshgrid(0:0.01:1,0:0.000001:0.0005);
Anfun=inline('2*(x-0.5).^2.*exp(5*x./2./50).*sin(n*pi*x)','x','n');
6第六讲 典型模型方程-热传导方的差分格式程

ADE Methods
同一时间步,同时左右扫描
6、Hopscotch Method
Comments
2、Richardaon’s Method: CTCS
3、Simple Implicit(Laasonen) Method
算子表示:
其中:
未知,每一个时间步长需要求解三对角方程组
放大因子:
4、Crank-Nicolson Method: famous
1 N 1 1 uN uN j 1 2u j j 1
where
G 1
3-D: ADI Methods
假设
∆t左移,右端合并即为C-N格式
其中 rx x 2 , ry y 2
at
at
3-D: ADI Methods
4、Splitting or Fractional-step Methods
5、ADE Methods
1-D:
或
Simple Explicit Method
差分格式的放大因子与精确解的放大因子比较
精度高
Simple Explicit Method
Simple Explicit Method: Example
Simple Explicit Method: Example
没有相位误差,幅值误差:1.88%
2-D: Crank-Nicolson Scheme
通常采用迭代方 法,需要比三对 角更多的计算资 源
3、2-D: ADI Methods
time
1/ 2 n 1 / 2 1/ 2 uin , uin 1, j , ui , j 1, j
1 n 1 n 1 uin, j 1 , ui , j , ui , j 1
第九章_热传导方程的差分解法_郑大昉

类似地,其偏微分用差分近似为: 类似地 其偏微分用差分近似为 近似为
∂ui, j,k ui, j,k+1 − ui, j,k = ∂t τ 2 ∂ ui, j,k ui+1, j,k − 2ui, j,k + ui−1, j,k = 2 ∂x h2 ∂2ui, j,k ui, j+1,k − 2ui, j,k + ui, j−1,k = 2 ∂y h2
∂ui,k ∂x ∂ui,k − ∂x + h
(9-18)
二阶中心差商可近似为 二阶中心差商可近似为: 可近似为
∂2ui,k ∂x
即:
2
=
−
(9-19)
ui+1,k − 2ui,k + ui−1,k ∂2u = 2 2 ∂x i,k h
(9-20)
时间的一阶差商近似为 近似为: 另, 对时间的一阶差商近似为
(9-27)
u(x, y,0) = ϕ(x, y)
(9-28)
其边界条件留待后面给出 边界条件留待后面给出. 留待后面给出
差分方法 仍设空间步长 h 仍设空间步长: 空间步长 时间步长: 时间步长 空间为: 网格. 空间为 N× M 网格
τ
则:
Nh = l,
M =s h
t = kτ , k = 0,1 2,... , x = ih, i = 0,1 N ,..., y = jh, j = 0,1,..., M
∆t
∂u ∆Q = −K(x, y, z, t)∆t∆S ∂n
(9-1)
t1
t2 t1
和
Q =∫ 1
∂u dt ∫∫ K(x, y, z, t) dS ∂n (S)
热传导方程的差分格式汇总

热传导方程的差分格式汇总1.显式差分格式:显式差分格式是最简单的一种方法,通过将导热方程时间和空间上的导数进行近似,引入差分算子,将方程转化为差分格式。
其中最常见的差分格式有:a. 前向差分法(Forward Difference Method):利用当前节点和其相邻节点的温度值进行计算。
例如,在一维离散情况下,可以使用公式:u(i,j+1)=u(i,j)+α(u(i+1,j)-2u(i,j)+u(i-1,j))b. 后向差分法(Backward Difference Method):利用当前节点和其相邻节点的温度值进行计算。
例如,在一维离散情况下,可以使用公式:u(i,j+1)=u(i,j)+α(u(i+1,j+1)-2u(i,j+1)+u(i-1,j+1))c. 中心差分法(Central Difference Method):利用当前和其相邻节点的温度值进行计算。
例如,在一维离散情况下,可以使用公式:u(i,j+1)=u(i,j)+α(u(i+1,j)-2u(i,j)+u(i-1,j))+β(u(i+1,j)-u(i-1,j))其中α和β是时间和空间步长的比例因子。
2.隐式差分格式:显式差分格式具有较大的稳定性限制。
为了克服这个问题,可以使用隐式差分格式,其中使用下一个时间步长的温度值来求解当前时间步长。
常见的隐式差分格式有:a. C-N差分法(Crank-Nicolson Method):利用前后两个时间步长的温度值进行计算。
例如,在一维离散情况下,可以使用公式:u(i,j+1)=u(i,j)+0.5α(u(i+1,j+1)-2u(i,j+1)+u(i-1,j+1))+0.5α(u(i+1,j)-2u(i,j)+u(i-1,j))b. 力学模拟法(Finite Element Method):将空间离散化后,通过引入有限元方法,将热传导问题转化为线性方程组,再通过求解线性方程组得到温度分布。
一维热传导方程的差分法

一维热传导方程的差分法一维热传导方程描述了一个物体内部热的传递规律。
这个方程可用于解决各种问题,如材料的温度分布、传热速率等。
对于一维热传导方程,可以通过差分法来求解。
差分法是一种数值求解法,通过将原方程离散化成差分形式,将导数转化为有限差分,从而得到差分方程组。
通过求解差分方程组就可以得到离散点上的数值解。
关于一维热传导方程的差分法,以下是具体步骤。
1. 确定精度和空间网格数在差分法中,需要首先确定精度和空间离散化的步长。
通常情况下,精度越高,计算量越大,但是结果也越接近真实情况。
空间网格数越多,计算量也会越大,但是离散化的结果也越接近真实情况。
因此,需要在计算效率和结果准确性之间做出权衡。
2. 离散化热传导方程将一维热传导方程离散化,得到差分方程组。
通过 Taylor 展开,将导数转化为有限差分的形式,得到如下式子:$$ \frac{T_{i+1}-2T_{i}+T_{i-1}}{\Deltax^{2}}=\frac{\partial^{2}T}{\partial x^{2}}|_{x=i\Delta x,t}=\frac{1}{\alpha}\frac{\partial T}{\partial t}|_{x=i\Delta x,t} $$其中,$T_i$ 表示在 $x=i\Delta x$ 处的温度值,$\Delta x$ 表示空间分割步长,$\frac{1}{\alpha}$ 表示材料的热扩散系数。
3. 构建差分方程组通过对差分方程组进行简单的变形,得到一个带有时间变化的差分方程组:其中,$n$ 表示时间步长,$\Delta t$ 表示时间离散化步长。
4. 初始条件和边界条件为了有效地求解差分方程组,我们需要知道初始条件和给定的边界条件。
在一维热传导方程中,初始条件是物体最初的温度分布,而边界条件通常包括物体边界的温度和热流量。
5. 使用迭代算法求解差分方程组通过使用迭代算法(如欧拉法、隐式迭代法、迭代加速法等),可以求解差分方程组的数值解。
差分方法

解题的时间范围为 ,初始条件是 。
为了有足够的精度,将初始化的网格作了两次细分。而作图的选项为Contour和Animation。
作为对比,可以更改初始条件为 ,即 。
资料来源:数学物理方程与Matlab可视化.
end
mesh(X,Y,un);
axis([-1 1 -1 1 -0.4 0]);
pause(0.1)
end
figure(2)
wn=0;
fork=1:N
wnn=2*(U0-u0)*(-1)^k.*sin(k.*pi.*RR).*exp(-k^2*pi^2*a2*TT)./(pi*k.*RR);
wn=wn+wnn;
,即
稳定且非振荡的条件为
截断误差为
另一种格式为
即
该式称为隐式格式。对任何步长都是恒稳定的。在 上取值的唯一限制是,要将截断误差保持在合理的程度上从而节约计算时间。
截断误差为
。
二、一维热传导方问题
2.1 无限长细杆的热传导
无限长细杆的热传导的定解问题是
利用Fourier变换求得问题的解是
其中取初始温度分布如下:
end
2.3输运问题(非齐次方程)
定解问题是
解析解为
其中
我们分别用解析解与差分方程的数值解画图。
计算中取 ,这个时间很短,因为这个过程的时间很短。解析解的图形如图3所示,其中图(a)是以 为变量所画的表面图,从图中可以看出变化的全貌,图(b)是初始状态,图(c)是最后的状态。解析解在初始状态所画出的图形与差分方程的解有一定的偏差。
为了作图,取 ,在级数求和中只取了前面10项。
【程序】
热传导方程三层并行差分格式初始条件的计算

1001-246X ( 2011 )04-0488-05热传导方程三层并行差分格式初始条件的计算左风丽1崔霞2袁光伟2(1.北京应用物理与计算数学研究所高性能计算中心,北京100088;2.北京应用物理与计算数学研究所计算物理重点实验室,北京100088)摘要:给出二维热传导问题的三层差分格式初始条件的一种显式计算方法,对于由此形成的内边界预估校正三层并行差分算法,证明稳定性和收敛性定理.并行数值试验表明,方法稳定,且与通常采用隐式格式计算初始条件的方法相比,易于程序实现;与已有的扰动算法相比,能大幅度减小误差.三层差分格式;初始条件;稳定显式计算方法 O241.82;O246A2010-08-122011-02-18国防基础科研项目( B1520110011)、中国工程物理研究院科学技术基金(2010A0202010)、计算物理实验室基金、国家自然科学基金(60973151,61033009)、863课题(2009AA01A134,2010AA012303)及973 课题(2011CB309702)资助项目左风丽(1971 -),女,河南,副研究员,硕士,主要从事大规模科学与工程并行计算研究.第4期490@@[ 1 ] Thomas J W. Numerical partial differential equations finite difference methods [ M]. Springer-Verlag, 1997.@@[2] Richtmyer R D,Morton K W.初值问题的差分方法[M].袁国兴,杜明笙,王汉强,译.广州:中山大学出版社,1992.@@[3]李德元,陈光南.抛物型方程差分方法引论[M].北京:科学出版社,1995.@@[4]陆金甫,关治.偏微分方程数值解法[M].第二版.北京:清华大学出版社,2004.@@[ 5 ] Yuan G W, Zuo F L. Parallel difference schemes for heat conduction equation [ J ]. International Journal of Computer Mathematics, 2003, 80 ( 8 ) : 995 - 999.@@[ 6 ] Günter S, Lackner K. A mixed implicit-explicit finite difference scheme for heat transport in magnetized plasmas [ J ]. Journal of Computational Physics, 2009, 228 : 282 - 293.@@[ 7 ] Yuan G W, Hang X D, Sheng Z Q. Parallel difference schemes with interface extrapolation terms for quasi-linear parabolic systems [J]. Science in China Series A: Mathematics, 2007, 50(2) : 253 -275.@@[ 8 ] Sheng Z Q, Yuan G W, Hang X D. Unconditional stability of parallel difference schemes with second order accuracy for parabolic equation [ J ]. Applied Mathematics and Computation, 2007, 184 : 1015 - 1031.@@[ 9 ] Yuan G W, Sheng Z Q, Hang X D. The unconditional stability of parallel difference schemes with second order convergence for nonlinear parabolic svstem [ J ]. Journal of Partial Differential Equations, 2007, 20 ( 1 ) : 1 - 20.@@[10]许士良.计算机常用算法[M].第二版.北京:清华大学出版社,1995.Initialization Method in Three-layer Parallel Difference Scheme for Heat EquationZUO FengliCUI Xia YUAN Guangwei。
2有限差分法及热传导数值计算PPT演示课件

t1
1 a11
(b1
a12t2
a13t3 )
t2
1 a 22
(b2
a 21t1 a 23t3 )
1 t3 a 33 (b3 a 31t1 a 32t2 )
•24
(2)假设一组解(迭代初场),记为: t1(0)、t2(0)并、t代3(0) 入迭代方程求得第一 次解
每次计算t1(1)均、t用2(1)、最t3(1新) 值代入。
(1) 平直边界上的节点
如图所示 边界节点 (m,n) 只能代表半个元体,若边界上有向 该元体传递的热流密度为q ,据能量守恒定律对该元体有:
tm1,n tm,n ytm,n1tm,n x
x
y 2
x
2
tm,n1 tm,n y
Φm,n
2xyyqw
0
Байду номын сангаас
xy tm ,n1 4 2 tm 1 ,n tm ,n 1 tm ,n 1 x2 Φ m ,n 2 x q w
非稳态项 的离散有三种不同的格式。如果将函数在节 点(n,i+1)对点(n,i)作泰勒展开,可有
•30
•31
由式(b)可得在点(n,i)处一阶导数的一种差分表示式 , 的向前差分:
类似地,将t在点(n,i-1)对点(n,i)作泰勒展开,可得 的向后差分的表达式:
如果将t在点(n,i+1)及(n,i-1)处的展开式相加,则可得 一阶导数的中心差分的表达式:
qw
y x
•16
(3) 内部角点
如图所示内部角点代表了 3/4 个元体,在同样的假设条 件下有
tm1,ntm,ny tm,n1tm,nx tm,n1tm,n x
x
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
function u=fai(x) u=4*x*(1-x); function u=g1(x) u=0; function u=g2(x) u=0;
9.3 二维热传导方程的差分解法
一、二维热传导方程 各向同性介质中无热源的二维热传导方程为:
u u u K ( ), , t T (9.9) t x x c x l , y s,
9.3 二维热传导方程的差分解法
例 9.2 求热传导方程混合问题:
u u u x , y , t y t x u ( x, y,) u ( x,, t ) u ( x,, t ) x , t u (, y, t ) y , t x u (, y, t ) y M , M y , t x u (, y, t ) M y M , t
9.2 一维热传导方程的差分解法
例 9.1 求热传导方程混合问题:
u u t x u ( x,) x( x) u (, t ) , u (, t )
x , t x t
的数值解,取 N=10,h=0.1,计算到 k=36 为止。
function main %热传导方程的差分解法 lda=1;l=1;h=0.05;alpha=0.5;tao=alpha* h^2/lda;T=tao*100;N=l/h;M=T/tao; for i=1:N+1 u(i,1)=fai((i-1)*h); end for k=1:M u(1,k)=g1(k*tao); u(N,k)=g2(k*tao); end for k=1:M for i=2:N u(i,k+1)=alpha*u(i+1,k)+(12*alpha)*u(i,k)+alpha*u(i-1,k); end plot([0:h:l],u(:,k+1)); hold on; pause(0.05); end
将差分格式(9.10)代入偏微分方程中得:
(9.10)
ui , j ,k ( )u i , j ,k (ui , j ,k u i , j ,k u i , j ,k ui , j ,k )
h
i ,,..., N , j ,,...,M , k ,,,...
k ,,,... i ,,,..., N j ,,,...,M
9.3 二维热传导方程的差分解法
对于节点 (i, j ) ,在 k 时刻有:
u i , j ,k u i , j ,k u i , j ,k t u i , j ,k u i , j ,k u i , j ,k u i , j ,k x h u i , j ,k u i , j ,k u i , j ,k u i , j ,k y h
第九章 热传导方程的差分解法
9.1 热传导方程概述
一、传入热量
t 时间内通过 S 横截面积传导的热量为: u Q K ( x, y, z , t )tS n u 其中,K ( x, y, z, t ) 是介质的热传导系数。 是温度沿 S 面的 n
法向微商,即温度梯度的法向分量。 通过包面 S 传入 V 的热量为:
9.1 热传导方程概述
五、三维齐次热传导方程 当介质均匀( c 、 和 K 为常数) 内无热源( F ( x, y, z, t ) )时: 、V
u c Ku, t
上式可表示为:
其中 x y z
u K u u u ( ) t c x y z
4、 一维热传导方程的差分格式
u i ,k u i ,k
整理得:
u i ,k u i ,k u i ,k h
u i ,k u i ,k ( - )u i ,k u i ,k
h
i ,,...,N , k ,,,...,M
上式为三维齐次热传导方程。
(9.2)
9.2 一维热传导方程的差分解法
一、一维热传导方程 各向同性介质中无热源的一维热传导方程为:
u u t x
二、初始、边界条件 初始条件:
u( x,) ( x)
xl
(9.4)
3、具体步骤
h
(1)给定 , h, , T , XN, YM ; (2)计算 h / , N XN / h, M YM / h, K T / ;
(3)计算初始值: ui , j , (ih, jh) ; 计算边界值: u, j ,k g (k , jh),u N , j ,k g (k , jh) ; ui ,,k g (k , ih),ui , N ,k g (k , ih) ; (4)用差分格式计算 ui , j ,k 。
(ui 1, j ,k ui 1, j ,k ui , j 1,k ui , j 1,k ) i 0,1,..., N , j 0,1,...,M i 0,1,2,..., N j 1,2,...,M 1 1, M 2 1,...,M 1 j 0,1,...,M j M 1 , M 1 1,...,M 2
function chafen2 %二维热传导方程的差分解法 lda=1;l=1;s=2;h=0.05;alpha=0.25;tao=alpha*h^2/lda; T=tao*1000;N=l/h;M=s/h;K=T/tao;M1=0.4*M;M2=0.6*M; for i=1:N+1 for j=1:M+1 u(i,j,1)=0; for k=1:K u(i,1,k)=0; u(i,M,k)=0; if M1<=j & j<=M2 u(1,j,k)=1; end end end end for k=1:K-1 for j=2:M for i=2:N u(i,j,k+1)=(1-4*alpha)*u(i,j,k)+alpha*(u(i+1,j,k)+u(i1,j,k)+u(i,j+1,k)+u(i,j-1,k)); end u(N+1,j,k+1)=u(N,j,k+1); if (2<=j&j<=M1) || (M2<=j&j<=M) u(1,j,k+1)=u(2,j,k+1); else u(1,j,k+1)=1; end end surf(u(:,:,k+1)); pause(0.05); end
的数值解,取 N=20,M=20,h=0.05,计算到 k=100 为止。
9.3 二维热传导方程的差分解法
解:将初始边界混合问题转化为差分格式
ui , j ,k 1 (1 4 )ui , j ,k ui , j , 0 0 u i , 0 , k ui , M , k 0 u0, j ,k u1, j ,k u N 1, j , k u N , j , k u0, j ,k 1
u i ,k x
和
u i ,k t
u i ,k u i ,k
u i ,k u i ,k h
u i ,k x
3、 二阶中心差商
u i ,k ui ,k x
x
h
ui ,k ui ,k u i ,k h
9.2 一维热传导方程的差分解法
Q dt F ( x, y, z, t )dV
t V
t
其中, F ( x, y, z, t ) 为内部热源密度,表示单位时间单位体积所产生 的热量。 三、消耗热量 V 内消耗热量:
Q
t
t
其中, c 为介质的比热容, 为质量密度。
u dt c dV V t
(9.11)
9.3 二维热传导方程的差分解法
初始、边界条件的差分格式:
u i , j , (ih, jh) u , j ,k g (k , jh) u N , j ,k g (k , jh) u i , ,k g (k , ih) u i , N ,k g (k , ih)
Q dt
t
t
S
u K ( x, y, z , t ) ds n
由矢量分析(高斯散度定理)可得:
Q dt [K ( x, y, z, t )u]dV
t V
t
其中, 是哈密顿算子。
9.1 热传导方程概述
二、产生热量 V 内所有热源产生的热量:
9.1 热传导方程概述
四、三维非齐次热传导方程 由能量守恒定律,即
Q Q Q
可得:
或
u t dt V [c t ( Ku) F ]dV u c ( Ku ) F ( x, y, z, t ) t
t
(9.1)
式(9.1)称为各向同性介质有热源的热传导方程,也叫三维 非齐次热传导方程。
t T
(9.6)
9.2 一维热传导方程的差分解法
边界条件: 3、 第三类边界条件:
u (, t ) x (t )u (, t ) g (t ) u (l , t ) (t )u (l , t ) g (t ) x
其中, (t )
i ,,..., N , j ,,...,M k ,,,..., j ,,...,M k ,,...,M , i ,,..., N