热传导方程的求解

应用物理软件训练

前 言

MATLAB 是美国MathWorks公司出品的商业数学软件，用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境，主要包括MATLAB和Simulink两大部分。

MATLAB是矩阵实验室（Matrix Laboratory）的简称，和Mathematica、

- -可修编.

. -

Maple并称为三大数学软件。它在数学类科技应用软件中在数值计算方面首屈一指。MATLAB可以进行矩阵运算、绘制函数和数据、实现算法、创建用户界面、连接其他编程语言的程序等，主要应用于工程计算、控制设计、信号处理与通讯、图像处理、信号检测、金融建模设计与分析等领域。本部分主要介绍如何根据所学热传导方程的理论知识进行MATLAB数值实现可视化。本部分主要介绍如何根据所学热传导方程的理论知识进行MATLAB数值实现可视化。本部分主要介绍如何根据所学热传导方程的理论知识进行MATLAB数值实现可视化。

本部分主要介绍如何根据所学热传导方程的理论知识进行MATLAB数值实现可视化。

- -可修编.

. -

题目：热传导方程的求解

目录

一、参数说明……………………………………………………………………..1 二、基本原理……………………………………………………………………..1 三、MATLAB程序流程图………………………………………………………3 四、源程序………………………………………………………………………….3 五、程序调试情况…………………………………………………………………6 六、仿真中遇到的问题……………………………………………………………9 七、结束语…………………………………………………………………………9 八、参考文献………………………………………………………………………10

- -可修编.

. -

一、 参数说明

U=zeros(21,101) 返回一个21*101的零矩阵 x=linspace(0,1,100);将变量设成列向量 meshz(u)绘制矩阵打的三维图

axis([0 21 0 1]);横坐标从0到21，纵坐标从0到1 eps是MATLAB默认的最小浮点数精度 [X,Y]=pol2cart(R,TH);效果和上一句相同 waterfall(RR,TT,wn)瀑布图

二、 基本原理

1、 一维热传导问题

（1）无限长细杆的热传导定解问题

利用傅里叶变换求得问题的解是：

取得初始温度分布如下

这是在区间0到1之间的高度为1的一个矩形脉冲，于是得

- -可修编.

. -

（2）有限长细杆的热传导定解问题

其中0x20，即L=20，取a=10且

得的解是

（3）非齐次方程定解问题是

解析解是

其中

- -可修编.

. -

2、二维热传导问题 定解问题

Ut=k^2(Uxx+Uyy) (0xa,0yb) U(x=0,y,t)=0, u(x=a,y,t)=sin3yb U（x,y=0,t）=0, u(x,y=b,t)=sin3xxacosaU（x,y,t=0）=0

3、三维热传导问题

球体内的热传导

令u=w+Uo,则w的定解问题是

Wt=ww W（r=ro）=0 W(t=to)=uo-Uo 解为

w2(Uouo)r(1)n2a2t/r2osinnrn1nen2ro r为空间变量，并用x，y表示。

三、 MATLAB程序流程图

开始

初始化定义

-

-可修编.

. -

预设矩阵

初始条件

用for语言

绘制动态图 四、 源程序

1、一维有限长细杆的热传导

x=0:20;t=0:0.01:1;a2=10; r=a2*0.01; u=zeros(21,101);

u(10:11,1)=1; 是把上述矩阵中的第10行，11行的第一列全部设成1 for j=1:100

u(2:20,j+1)=(1-2*r)*u(2:20,j)+r*(u(1:19,j)+u(3:21,j)); plot(u(:,j));

axis([0 21 0 1]);横坐标0到21，纵坐标0到1 pause(0.1)暂停0.1秒

end meshz(u)

2、非齐次方程的定解问题

a2=50;b=5;L=1;

[x,t]=meshgrid(0:0.01:1,0:0.000001:0.0005);

- -可修编.

. -

Anfun=inline('2/L*(x-L/2).^2.*exp(-b*x/2/a2).*sin(n*pi*x/L)','x','n','L','b','a2');%定义内联函数 u=0; for n=1:30

An=quad(Anfun,0,1,[],[],n,L,b,a2);%inline函数中定义x为向量，其它为标量

un=An*exp(-(n*n*pi*pi*a2/L/L+b*b/4/a2/a2).*t).*exp(b/2/a2.*x).*sin(n*pi*x/L); u=u+un; size(u);

mesh(x,t,u);%x,t，u都为501行101列的矩阵 figure subplot(2,1,1) plot(u(1,:)) subplot(2,1,2) plot(u(end,:)) end 差分法

dx=0.01;dt=0.000001;a2=50;b=5;c=a2*dt/dx/dx; x=linspace(0,1,100);%将变量设成列向量 uu(1:100,1)=(x-0.5).^2;%初温度为零 figure

- -可修编.

. -

subplot(1,2,1)%初始状态 plot(x,uu(:,1),'linewidth',1); axis([0,1,0,0.25]); subplot(1,2,2)%演化图 h=plot(x,uu(:,1),'linewidth',1); set(h,'EraseMode','xor') for j=2:200

uu(2:99,2)=(1-2*c)*uu(2:99,1)+c*(uu(1:98,1)+ uu(3:100,1))-... b*dt/dx*(uu(3:100,1)-uu(2:99,1)); uu(1,2)==0;uu(100,2)==0;%边界条件 uu(:,1)=uu(:,2); uu(:,1)

set(h,'YData',uu(:,1)); drawnow; pause(0.01) end

三维热传导问题

U0=2; u0=0; a2=2; N=10;

r=eps:0.05:1; theta=linspace(0,2*pi,100); t=0.1:0.001:0.2; [RR,TT]=meshgrid(r,t); figure(1)

- -可修编.

. -

[R,TH]=meshgrid(theta,r); [X,Y]=pol2cart(R,TH); for tt=1:100 un=0; for k=1:N

unn=2*(U0-u0)*(-1)^k.*sin(k.*pi.*(X.^2+Y.^2).^0.5).*... exp(-k^2*pi^2*a2*t(tt))./(pi.*(X.^2+Y.^2).^0.5);

un=unn+un; end

mesh(X,Y,un);

axis([-1 1 -1 1 -0.4 0]); pause(0.1) end figure(2) wn=0; for k=1:N

wnn=2*(U0-u0)*(-1)^k.*sin(k.*pi.*RR).*... exp(-k^2*pi^2*a2*TT)./(pi*k.*RR); wn=wnn+wn; end

waterfall(RR,TT,wn) xlabel('r')

- -可修编.

. -

ylabel('t')

五、 程序调试情况

1、有限长细杆的热传导

开始时

一段时间后

2、（1）非齐次方程的解析解

- -可修编.

. -

（2）非齐次方程的数值解（差分法）

3、二维热传导问题

4、三维热传导问题

- -可修编.

. -

解析解的动画图

解析解的瀑布图

六、仿真中遇到的问题

几乎所有的工程问题都能转化成数学模型来解，而且借助MATLAB，大多数的模型的数值解的精确度均能满足要求。但是，存在的问题也不少。首先，数值解法存在许多局限性，一个解只能适用于一个或几个模型，或者一个或几个方程。而解析解的得到能使我们得出所有同类问题的通解，并且精确度高于数值解。这 是由于数学的发展程度还不足以满足自然科学的发展要求，数值解法只是一个权 宜之计。其次，MATLAB虽然能处理大量的数学问题，但其命令繁多，再加上

- -可修编.

. -

各种工具箱，要完全学会和很好的使用MATLAB不是一件容易的事情，在编辑和阅读程序时通常要借助工具书查询相关命令，这样就增加了使用难度，使得MATLAB不能广泛的普及。再者，要合理的使用MATLAB来解决数学问题，必需是建立在良好的数学基础之上的，这就势必要求MATLAB的使用者有扎实的数学功底，这又给MATLAB的普及带来了挑战。最后，由于工程中的导热问题的数学模型并不一都能很顺利的建立，这就给使用MATLAB解决导热问题增加了难度。

七、结束语

在这短短的一周内从开始的一头雾水，到自己看书学习，到同学讨论，再进行整个题目的理论分析和计算，参考课程上的代码，写出自己的代码。

我们也明白了学无止尽的道理，在我们所查的很多参考书中，很多知识是我们从没有接触过的，我们对它的了解还仅限于皮毛，对它的很多功能以及函数还不是很了解，所以在这个学习的过程中我们穿越在知识的海洋中，一点一点吸取着它的知识。在MATLAB编程中需要很多的参考书，要尽量多的熟悉matlab自带的函数及其作用，因为matlab的自带函数特别多，基本上能够满足一般的数据和矩阵的计算，所以基本上不用你自己编函数。这一点对程序非常有帮助，可以使程序简单，运行效率高，可以节省很多时间。本次课设中用了很多MATLAB自带的函数，使程序变得很简单。

把基本的知识看过之后，就需要找一个实际的程序来动手编一下，不要等所有的知识都学好之后再去编程，你要在编程的过程中学习，程序需要什么知识再去补充，编程是一点一点积累的，所以你要需做一些随手笔记什么的。

在编写程序代码时，需要什么函数，需要什么模块就应该去着重看那个知识

- -可修编.

. -

点，不要一步登天，一步一步学，如果太急于把所有东西都学到，也是不好的，更是实现不了的。所以那时一天一天积累的，慢慢地学通这个软件。

八、参考文献

《数学物理方程的MATLAB解法与可视化》 彭芳麟著 清华大学

《量子物理学中的常用算法与程序》 井孝功 赵永芳 蒿凤有 编著 XX工业大学 《计算物理基础》 彭芳麟著 高等教育

- -可修编.