基于malab的偏微分方程箱在求解数学物理方程中的应用

基于malab的偏微分方程箱在求解数学物理方程中的应用

nb是世界上最流行的科技应用之一。这是1984年由美国mathrs开发的一个软件包理软件。但是Matlab软件在求解数学物理方程中的应用并没有受人们的足够重视。在工程领域,一般都利用相关的其他应用软件(如ANSYS软件等)或用VB、VC++等编程语言针对具体的数学物理方程编程,这样需要消耗大量的人力、物力。基于以上的现实情况,本文主要讨论MATLAB7.0的偏微分方程(PDE)工具箱在求解数学物理方程中的一些应用。通过作者应用研究经验表明:MATLAB7.0的PDEToolbox可求解《数学物理方程》等的教材中列出来的大多数定解问题,只是这些教材中求出的是解析解,对于初边值条件的可以用一个函数符号代替,而在用(MATLAB软件求解时必须是已知的表达式或具体的数值。目前的MATLAB7.0的PDEToolbox提供了研究和求解空间二维偏微分方程问题的一个强大而又灵活实用的平台。该工具箱能求解基本方程有椭圆型方程、抛物型方程、双曲型方程、特征值方程、椭圆型方程组以及非线性椭圆型方程。1任意的工具函数Matlab工具箱就是一些M文件的集合,用户可以修改工具箱中的函数,更为重要的是用户可以通过编制M文件来任意地添加工具箱中原来没有的工具函数。许多的专业领域在MATLAB中都有自己的工具箱。下面仅就用Matlab7.0PDEToolBox求解数学物理方程需要注意的问题作简单介绍。1.1b接触条件混合边界条件椭圆型方程:-∇·(c∇u)+au=f(1)抛物型方程:d∂u∂t−∇⋅(c∇u)+au=f(2)d∂u∂t-∇⋅(c∇u)+au=f(2)双曲型方程:d∂2u∂t2−∇⋅(c∇u)+au=f(3)d∂2u∂t2-∇⋅(c∇u)+au=f(3)特征值方程:-∇·(c∇u)+au=λdu(4)其中d是定义在Ω上的复变函数,λ是待求的特征值。椭圆型方程组:⎧⎩⎨⎪⎪⎪⎪⎪⎪⎪⎪−∇⋅(c11∇u1)−∇⋅(c12∇u2)−⋯−∇⋅(c1n∇un)+a11u1+a12u2+⋯+a1nun=f1−∇⋅(c21∇u1)−∇⋅(c22∇u2)−⋯−∇⋅(c2n∇un)+a21u1+a22u2+⋯+a2nun=f2⋯⋯−∇⋅(cm1∇u1)−∇⋅(cm2∇u2)−⋯−∇⋅(cmn∇un)+am1u1+am2u2+⋯+amnun=fm(5){-∇⋅(c11∇u1)-∇⋅(c12∇u2)-⋯-∇⋅(c1n∇un)+a11u1+a12u2+⋯+a1nun=f1-∇⋅(c21∇u1)-∇⋅(c22∇u2)-⋯-∇⋅(c2n∇un)+a21u1+a22u2+⋯+a2nun=f2⋯⋯-∇⋅(cm1∇u1)-∇⋅(cm2∇u2)-⋯-∇⋅(cmn∇un)+am1u1+am2u2+⋯+amnun=fm(5)非线性椭圆方程-∇·(c(u)∇u)+a(u)u=f(6)边界条件:①Dirichlet条件:hu=r(7)②Neumann条件:n·(c∇u)+qu=g(8)要特别说明的是:对于特征值问题仅限于齐次条件;对于非线性情形,系数g,q,h和r可以依赖于u;对于抛物型方程和双曲型方程,系数可以依赖于时间t。③Robin条件:就是未知函数在边界上的法向微商与其本身的线性组合的值。④混合边界条件:指既有Dirichlet条件又有Neumann条件,甚至还有的Robin条件的情形。初值条件:一般是含有时间(t)变量的抛物型方程和双曲型方程等才涉及初值条件。在使用MATLAB7.0的PDEToolbox之前,要把原方程(组)先转化成工具箱能识别的形式。1.2求解表的一般过程根据欲求解的方程(组)及相关的初边值条件,用MATLAB7.0的PDEToolbox和使用函数命令及相关编程手段来求解时皆可总结为以下几步:①根据未知变量的取值范围画出定解区域;②设置相关的初边值条件;③设置方程类型;④网格剖分;⑤解方程。2使用示例下面仅用MATLAB7.0的PDEToolbox求解抛物型方程、双曲型方程两类典型的方程实例。2.1求解线性系统的解析解设有半径为R的薄圆盘,其上下两底绝热,若圆盘边界上的温度保持为零度,且初始温度为已知,求圆盘内的瞬时温度分布规律。这个问题可以归结为求解下述的定解问题:⎧⎩⎨⎪⎪⎪⎪a2∂u∂t=∂2u∂x2+∂u∂y2,x2+y2=R2,t>0ut=0=φ(x,y)ux2+y2=R2=0(9){a2∂u∂t=∂2u∂x2+∂u∂y2,x2+y2=R2,t>0ut=0=φ(x,y)ux2+y2=R2=0(9)这里假设R=0.5,a=1,ut=0=0.5,t>0。用MATLAB7.0的PDEToolbox经过69步,140次函数调用,18次LU分解,139次的线性系统的求解求得本问题的色彩解如下:由图1可以很直观地看出该时刻圆盘的热传导规律:从圆盘的中央向四周逐渐传热,使得四周的温度逐渐升高。由初边界条件及相关的热传递规律知其解的是符合实际情况的。也说明了其解的正确性。用分离变量法求解式(9)的解析解时可以引出著名Bessel方程和Helmholtz方程,其相应的解析解是通过分离变量法转化成求解Bessel方程的特征值问题来求解。详细过程见参考文献。2.2ca多通道波问题的解析解假设考虑如下的具体方程如下:⎧⎩⎨⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪⎪∂2u∂t2−∇2u=0,(x,y)∈{(x,y)−1<x<1,−1<y<1}ux=±0=0∂u∂ny=±1=0ut=0=arctan(cos(π2x))∂u∂tt=0=3sin(πx)esin(π2y)(10){∂2u∂t2-∇2u=0,(x,y)∈{(x,y)-1<x<1,-1<y<1}ux=±0=0∂u∂ny=±1=0ut=0=arctan(cos(π2x))∂u∂tt=0=3sin(πx)esin(π2y)(10)用MATLAB7.0的PDEToolbox经过106步,6次failedattempts,226函数调用,27次LU分解,最终求解了225次线性系统得本问题的色彩解如图2所示。由图2(b)可以直观地看出在二维平面中波的传播现象,如波峰、波谷等,如果区域选择大一些波的转播现象更能看得清楚。关于式(10)的解析解的问题,目前只能用降维法求解出相应的Cauchy问题的解析解。由于篇幅限制和数值解本身是一个很大的数值矩阵的原因,两个实例的相应的数值解未在本文中列出。3数学物理方程的新方法、新出路与其他求解方法相比,此方法有计算速度快,计算方便,求解稳定可靠、节省财力人力等优点,为广大应用科学的