基于MATLAB在数理方程的应用

《MATLAB语言》课程论文

基于MATLAB在数理方程的应用

姓名：廖威

学号：12010245212

专业：通信工程

班级：通信班

指导老师：汤全武

学院：物理电气信息学院

完成日期：2011.12.12

基于MATLAB在数理方程的应用

（廖威 12010245212 2010级通信班）

[摘要] MATLAB 是近几年传播最快、影响最大的数学类应用软件。应用MATLAB 求解《数学物理方法》中的一些题目，使原来繁琐的手工计算变得简便，而且可将数理方程的解及一些特殊函数以图形的形式显示出来，形象、直观，便于理解。数理方程当中有许多的复杂的数值及数学符号的计算《数学物理方法》是许多理工专业的必修课和重要基础课，也是一门公认的难度大的课程。因课程内容抽象，数学推导繁琐，学生学习起来感到非常枯燥。MATLAB 是高性能的数值计算型数学类科技应用软件，具有优秀的数值计算功能和强大的数据可视化能力。本文以一些典型习题为例，介绍了MATLAB 在复变函数、积分变换、数理方程和特殊函数等方面的应用。

[关键词] MATLAB 积分变换数学物理方程特殊函数图形绘制

一、问题的提出

MATLAB是近几年传播最快、影响最大的数学类应用软件。应用MATLAB求解《数学物理方法》中的一些题目，使原来繁琐的手工计算变得简便，而且可将数理方程的解及一些特殊函数以图形的形式显示出来，形象、直观，便于理解。而且MATLAB强大的科学运算、灵活的程序设计、便捷的与其他程序和语言接口的功能，显示出其很强的优越性……

二、MATLAB在解偏微分方程中的应用

应用MATLAB 求解数学物理方程，可通过编程或直接利用偏微分方程工具箱求解，直接利用偏微分方程工具箱更为简单、方便。在数理方程课上我们学习解矩形域方程的问题：

例1：在矩形域-0.5

(1)启动偏微分方程求解界面。

(2)在MATLAB 命令窗口中键入pderect ([-0.5 0.5 -0.5 0.5]), 选择Boundary→Remove All SubdomainBorders 菜单项，得出偏微分方程的求解区域。

（3）单击偏微分方程界面工具栏中的PDE 图标，选择其中的Parabolic 选项，将给定的偏微分方程的参数输入到该对话框中。

（4）边界条件由Boundary 菜单下的Specify Boundary conditions 确定，

输入边界条件h=1，r=0。

（5）单击工具栏等号按钮，得到偏微分的解，单击图形设置按钮，得出如图1所示的结果。

MATLAB程序如下：

clear all %泊松方程的解

pderect([-0.5 0.5 -0.5 0.5 ]); %微分方程的求解区域

h1=1,r=0; %初始条件

运行结果如图1所示。

图1 泊松方程的解

通过以上程序对求解方程使其为可视化，加深对问题的了解。现在我们来解决另一个问题。

三、MATLAB在解特殊方程中的应用：

例2.画出连带勒让德多项式（2 阶）的曲线分布，

MATLAB语句为：

x=-1:0.01:1 %求解范围

plot(x,LEGENDRE(2,x)) %sunlink挂函数语句，输出图形运行结果如图2所示。

图2 连带勒让德多项式(2 阶)的曲线分布

例3：画出第一类贝塞尔函数J15 （x）的曲线分布。

MATLAB 语句为：

x=0:0.01:50 %求解范围

plot(x,BESSELJ(15,x)) %函数语句，输出图形

运行结果如图3所示。

图3 贝塞尔函数（15 阶）J15 （x ）的曲线分布

以上只是我们用MATLAB 解一些经典的数理方程，MATLAB 在显示其强大的绘图和解方程上的优势。下面我们来看MATLAB 在解实际应用方程的运用，先来看在解齐次波动方程。

四、作出有界弦振动问题解析解的图形分布：

02=-xx tt u a u l x ≤≤0

0,00====l x x u u

（1）

)(),(00x u x u t t

t φ?==== 解析解：u=l

x n l at n l at n n n n πππsin )sin cos

(1B +A ∑∞=（取系数为1） （2） 例4：n=1本振解随时间变化的图形分布。

clear

a=1;l=1;x=0:0.05:1;t=0:0.001:3;u=0; n=1; %函数参数

[X,T]=meshgrid(x,t); %函数变量

u=(cos(n*pi*a*T/l)+sin(n*pi*a*T/l)).*sin(n*pi*X/l) %函数定义式

figure(1)

axis([0 1 -0.05 0.05]) %图形取点，x 分布 mesh(X,T,u)

title('n=1本振随时分布的透视图') %图像名

xlabel('x') %图像坐标轴'x'

ylabel('t') %图像坐标轴'y'

zlabel('u') %图像坐标轴'u'

运行结果如图4所示。

'

图4 n=1本振解随时间变化的图形

例5：n=2本振解随时间变化的图形分布。

clear

a=1;l=1;x=0:0.05:1;t=0:0.001:3;u=0; n=2; %函数参数

[X,T]=meshgrid(x,t); %函数变量

u=(cos(n*pi*a*T/l)+sin(n*pi*a*T/l)).*sin(n*pi*X/l) %函数定义式figure(1)

axis([0 1 -0.05 0.05]) %图形取点，x分布mesh(X,T,u)

title('n=2本振随时分布的透视图')

xlabel('x') %图像坐标轴'x' ylabel('t') %图像坐标轴'y' zlabel('u') %图像坐标轴'u'

运行结果如图5所示。

图5 n=2本振解随时间变化的图形分布

五、MATLAB在复变函数中的应用

1、复数计算

例6：求下列复数的实部、虚部、共轭复数、模与辐角。

（1）（1+i）^6 （2）i/(1+i)+(1-i)/i

定义复数后，利用简单的MATLAB 语句：real（）、imag( )、conj( )、abs( )、angle( )可直接得出复数的的实部、虚部、共轭复数、模与辐角。

解：>>a=[(1+i)^6 i/(1+i)+(1-i)/i] %函数方程

a=0-8.0000i-0.5000-0.5000i %函数化简

>>real(a) %求实部

ans=

0-0.5000

>>imag(a) %求虚部

ans=

-8.0000-0.5000

>>conj(a) %求共轭复数

ans=

0+8.0000i-0.5000+0.5000i

>>abs(a) %求复数的模

ans=

8.0000 0.7071

>>angle(a) %求辐角

ans=

-1.5708-2.3562

例7：求方程z3 ＋8＝0 的所有的根。

调用内部函数Solve 进行求解。

解：>>Solve(‘z^3+8=0’) %调用函数Solve求方程的根ans =

-2

1+i*3^(1/2)

1-i*3^(1/2)

2、复变函数的图形绘制

在复变函数的教学与研究中，对于复变函数的图形绘制往往繁冗复杂，仅凭手工很难画出来。MATLAB软件中自带了图像处理工具包，它是由一系列支持图像处理操作的函数组成的，它囊括了几乎所有主流的科学计算中所涉及的图像处理功能。借助MATLAB的绘图功能，可以快捷、准确地绘出图形，使教学变得形象、直观、生动，有利于学生观察图形的形状，掌握图形的性质。

例8:绘制函数z3的图像。

输入MATLAB 语句如下:

z=cplxgrid(20); %在复平面上面得到一个21*41的网格

w=z.^3; %函数方程

surf(real(z),imag(z),real(w),imag(w)); %绘制三维曲线图

title('z^3'); %图像名

运行结果如图6所示。

图6 z3的图像

例9:作出圆周│z│ = 3 在映射ω = z + 1./z下的像。

编写相应的m 文件yingshe． m 为:

syms x y z t %定义变量

t=-pi:0.01:pi; %求解范围

x=3*cos(t); %定义x

y=3*sin(t); %定义y

z=x+i*y; %复数z 型式

w=z+1./z; %函数方程

subplot(2,1,1); %选择2×1个区中的1号区

plot(z) ; %函数语句，输出图形

title('z=3*cos(t)+i*sin(t)'); %图像名

axis equal %纵、横坐标轴采用等长刻度 subplot(2,1,2); %选择2×1个区中的2号区

plot(w); %函数语句，输出图形

title('w=z+1./z') %图像名

axis equal %纵、横坐标轴采用等长刻度 运行结果如图7所示。

图7 z =3在映射ω=z+z

1下的像 3、Taylor 级数展开

Taylor 级数的展开可以用符号工具箱中的Taylor( )函数直接导出。 例10：将函数1/（1－z ）^2在z=0处展成幂级数。

解：>>syms z ； %定义变量

>>Taylor(1/(1-z)^2,0) %Taylor( )函数直接导出

ans =

1+2*z+3*z^2+4*z^3+5*z^4+6*z^5

4、留数计算

通过求极限的方法计算留数。

例11：确定函数f （z ）＝ z ＋2i/(z^5 ＋4z^3)的极点，并求出函数在极点处的留数。

解：通过对函数的分析，可知z=0是三阶极点，z=2i,-2i 是单极点。

>>syms z; %定义变量

>>f=(z+2*i)/(z^5+4*z^3); %函数方程

>>limit(diff(f*z^3,z,2)/prod(1:2),z,0) %求函数在极点0处的留数 ans=

-1/8*i

>>limit(f*(z-2i),z,2i) %求函数在极点2i 处的留数 ans=

1/8*i

>>limit(f*(z+2i),z,-2i) %求函数在极点-2i 处的留数

ans=

有理函数的部分分式展开。MATLAB 语言中给出了现成的数值函数residue( )求取有理函数的部分分式展开表示，该函数的调用格式为:［R ，P ，K ］= residue( B ，

A)。返回留数、极点和两个多项式比值B( s) /A( s) 的部分分式展开的直接项

()

()=A B s s ()()()()()()()s n n K P s R P s R P s R +-++-+- 2211 (3)

向量B 和A 为分子、分母以s 降幂排列的多项式系数，向量R 是返回的留数，向量P 是返回的极点，向量K 由B( s)/A( s)的商的多项式系数组成，若length( B)＜length( A)，则K 为空向量，否则，length( K) =length( B)－length( A)+1。 如果存在m 级极点则有P(j)=…=P(j+m －1)，展开式包括以下形式:

()

()()()()()m j m j j j j j P s R P s R P s R )()(121-++-+--++ (4)

例12:计算函数f(z)= (z^2－4z+3)/(z^3+2z^2+z)的留数，并求其部分分式展开。

在MATLAB 窗口输入:

>>B=[1 -4 3]; %向量B

A=[1 2 1]; %向量A

[R,P,K]=residue(B,A); %函数调用

[B,A1]=rat(R);

[B,A1,P];

ans =

－2 1 －1

－8 1 －1

3 1 0

结果表明:留数为Res ［f(z)，－1］=－2，Res ［f(z)，0］=3，部分分式展开为

()z

z z f z 3)1(8122++-++-=

5、积分计算

积分计算可以转化为留数计算，先求被积函数的留数，利用留数定理求复变函数的积分。

例13：求c 矣dz （z^2 ＋1）（z^2 ＋4）（其中C 是逆时针方向的圆周：z ＝32）

解：>>[R,P,K]=residue([1],[1,0,5,0,4]) %函数调用

R= %求被积函数的留数 0+0.0833i

0-0.0833i

0-0.1667i

0+0.1667i

P= %四个极点

0+2.0000i

0-2.0000i

0+1.0000i

0-1.0000i

K= %（展开式中无此解析项） ［］

可见在圆周z ＝32内有四个极点，所以积分值S 为：

>>S=2*pi*i*sum(R) %计算出积分的值 S=

利用MATLAB 进行复数运算、级数。

六、MATLAB 在傅里叶变换的应用

在MATLAB 语言中使用fourier 函数来实现Fourier 变换。首先定义符号变量t 描述时域表达式f,直接调用fourier 函数就可求出所需的时域函数的Fourier 变换式。

该函数的调用格式为F=fourier(f)采用默认的t 为时域变量F=fourier( f ，u ，v),用户指定时域变量u 和频域变量v 使用ifourier 函数来实现Fourier 逆变换。它的使用格式如下:

f=ifourier(F),按默认变量进行Fourier 逆变换;

f=ifourier(F ，u ，v),将u 的函数变换成v 的函数。

例14:已知f(t)=e^(－t)sin2tu(t)，求f(t)的Fourier 变换，并对结果进行Fourier 逆变换。

输入MATLAB 语句如下:

syms t w %定义变量

f=exp(－t)*sin(2*t)*sym(‘Heaviside(t)’);

%时域表达式f

F=fourier( f) %求Fourier变换

运行结果F=

2/(1+2*i+i*w)/(1－2*i+i*w)

f=ifourier(F，w，t) %求Fourier 逆变换

运行结果f=

exp(－t)*sin(2*t)*Heaviside(t)返回原函数。

上面证明了在MATLAB语言中可以使用fourier函数来实现Fourier变换。下面我们来探讨能否在MATLAB语言中使用laplace函数来实现Laplace变换。

七、MATLAB在Laplace变换的应用

在MATLAB语言中使用laplace函数来实现Laplace变换。首先定义符号变量t，描述时域表达式f，直接调用laplace函数就可求出所需的时域函数的Laplace 变换式。该函数的调用格式为:

F=laplace(f)，采用默认的t为时域变量;

F=laplace(f，u，v)，用户指定时域变量u和频域变量v。

使用ilaplace函数来实现Laplace变换。它的使用格式如下:

f=ilaplace(F)，按默认变量进行Laplace逆变换;

f=ilaplace(F，u，v)，将u的函数变换成v的函数。

例15:设函数f(t)=e^（－t）t^2.sin(t+π)，求其Laplace变换，并对结果进行Laplace逆变换。

输入MATLAB 语句如下:

syms t s %定义变量

f=exp(－t)*t^2*sin(t+pi); %时域表达式f

F=laplace(f，t，s) %求Laplace变换。

运行结果F=

－8/((s+1)^2+1)^3*(s+1)^2+2/((s+1)^2+1)^2

f=ilaplace(F，s，t) %求Laplace逆变换。

运行结果f=

－exp(－t)*t^2* sin(t)

因为sin(t+π)=－sint，所以得出的就是原函数。

八、结论

从以上利用MATLAB语言对基本数理方程（经典方程﹑有界弦振动问题等）的分析，我们不难得出以下结论：

1.应用MATLAB求解《数学物理方法》中的一些题目，使原来繁琐的手工计算变得简便，而且可将数理方程的解及一些特殊函数以图形的形式显示出来，形象、直观，便于理解。而且MATLAB强大的科学运算、灵活的程序设计、便捷的与其他程序和语言接口的功能，显示出其很强的优越性……

2.应用MATLAB求解《数学物理方法》中，在解经典方程、复变函数问题的时候解决掉了手工计算式子多，计算繁杂，求解过程复杂的问题。在解有界弦振动问题，运用MATLAB编程求解得出图形，其汇编语言、程序编写相当简单。得出的图形直观，在帮助同学们理解复杂的数理方程的意义、掌握这门与我们物电专业相关的学科有很大的帮助。

九、课程体会

经过一学期紧张而又有序的课程学习，在忙碌之余也得到了颇多的收获。下面是我的一些心得体会：

1.认识主客体

在教学过程中，应坚持以数学物理为主，MATLAB为辅，千万不能把数学课变成计算机语言课，“体脑倒挂”、主次不分。因此，无论是课堂演示、课外作业还是数学实验，都必须紧密结合教材内容，做到有的放矢。要求教师在课前必须做充分的准备，每一个实例、每一个程序都要精心设计。另外，学生不能依靠计算机语言来完全替代理论的计算过程。我们需知道MATLAB只是验证理论、计算的一种工具。

2.充分调动学习的热情

由于课程内容困难，绝大多数同学都会对教学内容有恐惧心理。因此，必须利用MATLAB的强大图形处理能力及易用性消除大家的恐惧，增加其学习的动力。每次演示的时候，向同学们展示相关内容的最优美的部分，如动画、三维图像、图像的拖曳、声音背景等，从而吸引广大同学的眼球，提高学习的主观能动性。另外，往届同学所编制的一些样本程序同样可以打破同学们的恐惧心理，提高其自信心。

3.抓住学习的本质

尽管MATLAB语言非常简单，但教授学生掌握还是需要一个过程。因此，必须掌握一个尺度，使用正确的教授方法。我们的经验是采用演示的形式，教会基本操作、基本命令，结合具体的例子进行编程。更重要的是，同学们必须学会使用MATLAB强大的联机帮助功能以及联机演示程序所给出的实例。事实上，MATLAB 里面包含众多的工具箱。任何一门课程都无法教会同学们精通所有的内容，只有教授大家熟练掌握了自学的本领，才是学习这门语言的捷径。

4.结合各种方法，提高学习效率。

除了常规的教学方法外，在利用MATLAB辅助教学的过程中，作者也试用了基于问题的学习方法等。比如，在特殊函数绘图以及数值计算方法的实验中，要求广大同学自己学习命令，绘制教材上的曲线，求解书本上的例题或作业；在这个过程中，鼓励同学们发扬团队精神，利用集体智慧解决问题；鼓励“头脑风暴法”等的使用。这些方法的应用，一方面挑战了同学们的智力，一面提升了其合作精神，是学习数学物理方法、学习MATLAB语言的有效途径。

[参考文献]

[1] 刘卫国.MATLAB程序设计与应用（第二版）[M].北京：高等教育出版社，2006 [2］求是科技编著.MATLAB7.0 从入门到精通[M].北京:人民邮电出版社,2006. [3］张志涌.精通MATLAB[M].北京:北京航天大学出版社,2000..

[4］张元林.积分变换[M].北京:高等教育出版社(第四版),2003.

[5］薛定宇.高等应用数学的MATLAB 求解[M].北京:清华大学出版社,2004. [6］彭芳麟.数学物理方程的MATLAB 解法与可视化[M].北京:清华大学出版社, 2004.

[7］梁昆淼.数学物理方法[M].北京:高等教育出版社(第三版),1998.

[8］杨华军.数学物理方法与计算机仿真[M].北京:电子工业出版社,2006.

[9] 陈育宁.数学物理方法.[M].宁夏;人民教育出版社，2008

[10]张元林.积分变换[M].北京:高等教育出版社(第四版),2003.