(完整版)偏微分方程的MATLAB解法
- 格式:doc
- 大小:528.99 KB
- 文档页数:34
偏微分方程Matlab 数值解法(补充4)Matlab 可以求解一般的偏微分方程,也可以利用偏微分方程工具箱中给出的函数求解一些偏微分方程。
1 偏微分方程组求解Matlab 语言提供了pdepe()函数,可以直接求解偏微分方程(,,,)[(,,,)](,,,)m mu u u u C x t u x x f x t u s x t u x t x x x-∂∂∂∂∂=+∂∂∂∂∂ (4.1)这样,偏微分方程可以编写为以下函数的描述,其入口为[,,](,,,)x c f s pdefun x t u u =其中:pdefun 为函数名。
由给定输入变量可计算出,,c f s 这三个函数。
边界条件可以用下面的函数描述(,,)(,,).*(,,,)0up x t u q x t u f x t u x∂+=∂ (4.2) 这样的边值函数可以写为Matlab 函数[,,,](,,,)a a b b x p q p q pdebc x t u u =初始条件数学描述为00(,)u x t u =,编写一个简单的函数即可0()u pdeic x =还可以选择x 和t 的向量,再加上描述这些函数,就可以用pdepe ()函数求解次偏微分方程,需要用如下格式求解(,@,@,@,,)sol pdepe m pdefun pdeic pdebc x t =【例1】 试求下列偏微分方程2111222221220.024()0.17()u u F u u t x u u F u u tx ⎧∂∂=--⎪⎪∂∂⎨∂∂⎪=+-⎪∂∂⎩ 其中: 5.7311.46()xx F x e e -=-,且满足初始条件1(,0)1u x =,2(,1)0u x =及边界条件1221(0,)0,(0,)0,(1,)1,(1,)0u u t u t u t t x x∂∂====∂∂解:对照给出的偏微分方程和(4.1),可将原方程改写为111222120.024/()1.*10.17/()u u x F u u u u x F u u t x ∂∂--⎡⎤⎡⎤⎡⎤⎡⎤∂∂=+⎢⎥⎢⎥⎢⎥⎢⎥∂∂-∂∂⎣⎦⎣⎦⎣⎦⎣⎦可知0m =,且1122120.024/()1,,10.17/()u x F u u c f s u x F u u ∂∂--⎡⎤⎡⎤⎡⎤===⎢⎥⎢⎥⎢⎥∂∂-⎣⎦⎣⎦⎣⎦编写下面的Matlab 函数function [c,f,s]=c7mpde(x,t,u,du)c=[1;1];y=u(1)-u(2);F=exp(5.73*y)-exp(-11.46*y);s=F*[-1;1]; f=[0.024*du(1);0.17*du(2)];套用(4.2)中的边界条件,可以写出如下的边值方程左边界2010.*00f u ⎡⎤⎡⎤⎡⎤+=⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦,右边界1100.*100u f -⎡⎤⎡⎤⎡⎤+=⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦ 编写下面的Matlab 函数function [pa,qa,pb,qb]=c7mpbc(xa,ua,xb,ub,t) pa=[0;ua(2)];qa=[1;0];pb=[ub(1)-1;0];qb=[0,1]; 另外,描述初值的函数function u0=c7mpic(x) u0=[1;0];有了这三个函数,选定x 和t 向量,则可以由下面的程序直接求此微分方程,得出解1u 和2u 。
引言偏微分方程定解问题有着广泛的应用背景。
人们用偏微分方程来描述、解释或者预见各种自然现象,并用于科学和工程技术的各个领域fll。
然而,对于广大应用工作者来说,从偏微分方程模型出发,使用有限元法或有限差分法求解都要耗费很大的工作量,才能得到数值解。
现在,MATLAB PDEToolbox已实现对于空间二维问题高速、准确的求解过程。
偏微分方程如果一个微分方程中出现的未知函数只含一个自变量,这个方程叫做常微分方程,也简称微分方程;如果一个微分方程中出现多元函数的偏导数,或者说如果未知函数和几个变量有关,而且方程中出现未知函数对几个变量的导数,那么这种微分方程就是偏微分方程。
常用的方法有变分法和有限差分法。
变分法是把定解问题转化成变分问题,再求变分问题的近似解;有限差分法是把定解问题转化成代数方程,然后用计算机进行计算;还有一种更有意义的模拟法,它用另一个物理的问题实验研究来代替所研究某个物理问题的定解。
虽然物理现象本质不同,但是抽象地表示在数学上是同一个定解问题,如研究某个不规则形状的物体里的稳定温度分布问题,由于求解比较困难,可作相应的静电场或稳恒电流场实验研究,测定场中各处的电势,从而也解决了所研究的稳定温度场中的温度分布问题。
随着物理科学所研究的现象在广度和深度两方面的扩展,偏微分方程的应用范围更广泛。
从数学自身的角度看,偏微分方程的求解促使数学在函数论、变分法、级数展开、常微分方程、代数、微分几何等各方面进行发展。
从这个角度说,偏微分方程变成了数学的中心。
一、MATLAB方法简介及应用1.1 MATLAB简介MATLAB是美国MathWorks公司出品的商业数学软件,用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境,主要包括MATLAB和Simulink两大部分。
1.2 Matlab主要功能数值分析数值和符号计算工程与科学绘图控制系统的设计与仿真数字图像处理数字信号处理通讯系统设计与仿真财务与金融工程1.3 优势特点1) 高效的数值计算及符号计算功能,能使用户从繁杂的数学运算分析中解脱出来;2) 具有完备的图形处理功能,实现计算结果和编程的可视化;3) 友好的用户界面及接近数学表达式的自然化语言,使学者易于学习和掌握;4) 功能丰富的应用工具箱(如信号处理工具箱、通信工具箱等) ,为用户提供了大量方便实用的处理工具。
偏微分方程是描述自然界中动态过程的重要数学工具,在工程领域中,求解偏微分方程是很多实际问题的重要一步。
MATLAB作为一种强大的数学计算软件,提供了丰富的工具和函数来求解各种类型的偏微分方程。
本文将介绍使用MATLAB求解初边值问题的偏微分方程的方法和步骤。
一、MATLAB中的偏微分方程求解工具MATLAB提供了几种可以用来求解偏微分方程的工具和函数,主要包括:1. pdepe函数:用于求解偏微分方程初边值问题的函数,可以处理各种类型的偏微分方程,并且可以自定义边界条件和初始条件。
2. pdepeplot函数:用于绘制pdepe函数求解得到的偏微分方程的解的可视化图形,有助于直观地理解方程的解的特性。
3. pdetool工具箱:提供了一个交互式的图形用户界面,可以用来建立偏微分方程模型并进行求解,适用于一些复杂的偏微分方程求解问题。
二、使用pdepe函数求解偏微分方程初边值问题的步骤对于给定的偏微分方程初边值问题,可以按照以下步骤使用pdepe函数进行求解:1. 定义偏微分方程需要将给定的偏微分方程转化为标准形式,即将偏微分方程化为形式为c(x,t,u)∂u/∂t = x(r ∂u/∂x) + ∂(p∂u/∂x) + f(x,t,u)的形式。
2. 编写边界条件和初始条件函数根据实际问题的边界条件和初始条件,编写相应的函数来描述这些条件。
3. 设置空间网格选择合适的空间网格来离散空间变量,可以使用linspace函数来产生均匀分布的网格。
4. 调用pdepe函数求解偏微分方程将定义好的偏微分方程、边界条件和初始条件函数以及空间网格作为参数传递给pdepe函数,调用该函数求解偏微分方程。
5. 可视化结果使用pdepeplot函数绘制偏微分方程的解的可视化图形,以便对解的性质进行分析和理解。
三、实例分析考虑一维热传导方程初边值问题:∂u/∂t = ∂^2u/∂x^2, 0<x<1, 0<t<1u(0,t) = 0, u(1,t) = 0, u(x,0) = sin(πx)使用MATLAB求解该初边值问题的步骤如下:1. 定义偏微分方程将热传导方程化为标准形式,得到c(x,t,u) = 1, r = 1, p = 1, f(x,t,u) = 0。
Matlab 偏微分方程求解方法目录:§1 Function Summary on page 10-87§2 Initial Value Problems on page 10-88§3 PDE Solver on page 10-89§4 Integrator Options on page 10-92§5 Examples” on page 10-93§1 Function Summary1.1 PDE Solver” on page 10-871,2 PDE Helper Functi on” on page 10-871.3 PDE SolverThis is the MATLAB PDE solver.PDE Helper Function§2 Initial Value Problemspdepe solves systems of parabolic and elliptic PDEs in one spatial variable x and time t, of the form)xu ,u ,t ,x (s ))x u ,u ,t ,x (f x (x x t u )x u ,u ,t ,x (c m m ∂∂+∂∂∂∂=∂∂∂∂- (10-2) The PDEs hold for b x a ,t t t f 0≤≤≤≤.The interval [a, b] must be finite. mcan be 0, 1, or 2, corresponding to slab, cylindrical, or spherical symmetry,respectively. If m > 0, thena ≥0 must also hold.In Equation 10-2,)x /u ,u ,t ,x (f ∂∂ is a flux term and )x /u ,u ,t ,x (s ∂∂ is a source term. The flux term must depend on x /u ∂∂. The coupling of the partial derivatives with respect to time is restricted to multiplication by a diagonal matrix )x /u ,u ,t ,x (c ∂∂. The diagonal elements of this matrix are either identically zero or positive. An element that is identically zero corresponds to an elliptic equation and otherwise to a parabolic equation. There must be at least one parabolic equation. An element of c that corresponds to a parabolic equation can vanish at isolated values of x if they are mesh points.Discontinuities in c and/or s due to material interfaces are permitted provided that a mesh point is placed at each interface.At the initial time t = t0, for all x the solution components satisfy initial conditions of the form)x (u )t ,x (u 00= (10-3)At the boundary x = a or x = b, for all t the solution components satisfy a boundary condition of the form0)xu ,u ,t ,x (f )t ,x (q )u ,t ,x (p =∂∂+ (10-4) q(x, t) is a diagonal matrix with elements that are either identically zero or never zero. Note that the boundary conditions are expressed in terms of the f rather than partial derivative of u with respect to x-x /u ∂∂. Also, ofthe two coefficients, only p can depend on u.§3 PDE Solver3.1 The PDE SolverThe MATLAB PDE solver, pdepe, solves initial-boundary value problems for systems of parabolic and elliptic PDEs in the one space variable x and time t.There must be at least one parabolic equation in the system.The pdepe solver converts the PDEs to ODEs using a second-order accurate spatial discretization based on a fixed set of user-specified nodes. The discretization method is described in [9]. The time integration is done with ode15s. The pdepe solver exploits the capabilities of ode15s for solving the differential-algebraic equations that arise when Equation 10-2 contains elliptic equations, and for handling Jacobians with a specified sparsity pattern. ode15s changes both the time step and the formula dynamically.After discretization, elliptic equations give rise to algebraic equations. If the elements of the initial conditions vector that correspond to elliptic equations are not “consistent” with the discretization, pdepe tries to adjust them before eginning the time integration. For this reason, the solution returned for the initial time may have a discretization error comparable to that at any other time. If the mesh is sufficiently fine, pdepe can find consistent initial conditions close to the given ones. If pdepe displays amessage that it has difficulty finding consistent initial conditions, try refining the mesh. No adjustment is necessary for elements of the initial conditions vector that correspond to parabolic equations.PDE Solver SyntaxThe basic syntax of the solver is:sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan)Note Correspondences given are to terms used in “Initial Value Problems” on page 10-88.The input arguments arem: Specifies the symmetry of the problem. m can be 0 =slab, 1 = cylindrical, or 2 = spherical. It corresponds to m in Equation 10-2. pdefun: Function that defines the components of the PDE. Itcomputes the terms f,c and s in Equation 10-2, and has the form[c,f,s] = pdefun(x,t,u,dudx)where x and t are scalars, and u and dudx are vectors that approximate the solution and its partial derivative with respect to . c, f, and s are column vectors. c stores the diagonal elements of the matrix .icfun: Function that evaluates the initial conditions. It has the formu = icfun(x)When called with an argument x, icfun evaluates and returns the initial values of the solution components at x in the column vector u.bcfun:Function that evaluates the terms and of the boundary conditions. Ithas the form[pl,ql,pr,qr] = bcfun(xl,ul,xr,ur,t)where ul is the approximate solution at the left boundary xl = a and ur is the approximate solution at the right boundary xr = b. pl and ql are column vectors corresponding to p and the diagonal of q evaluated at xl. Similarly, pr and qr correspond to xr. When m>0 and a = 0, boundedness of the solution near x = 0 requires that the f vanish at a = 0. pdepe imposes this boundary condition automatically and it ignores values returned in pl and ql.xmesh:Vector [x0, x1, ..., xn] specifying the points at which a numerical solution is requested for every value in tspan. x0 and xn correspond to a and b , respectively. Second-order approximation to the solution is made on the mesh specified in xmesh. Generally, it is best to use closely spaced mesh points where the solution changes rapidly. pdepe does not select the mesh in automatically. You must provide an appropriate fixed mesh in xmesh. The cost depends strongly on the length of xmesh. When , it is not necessary to use a fine mesh near to x=0 account for the coordinate singularity.The elements of xmesh must satisfy x0 < x1 < ... < xn.The length of xmesh must be ≥3.tspan:Vector [t0, t1, ..., tf] specifying the points at which a solution is requested for every value in xmesh. t0 and tf correspond tot and f t,respectively.pdepe performs the time integration with an ODE solver that selects both the time step and formula dynamically. The solutions at the points specified in tspan are obtained using the natural continuous extension of the integration formulas. The elements of tspan merely specify where you want answers and the cost depends weakly on the length of tspan.The elements of tspan must satisfy t0 < t1 < ... < tf.The length of tspan must be ≥3.The output argument sol is a three-dimensional array, such that•sol(:,:,k) approximates component k of the solution .•sol(i,:,k) approximates component k of the solution at time tspan(i) and mesh points xmesh(:).•sol(i,j,k) approximates component k of the solution at time tspan(i) and the mesh point xmesh(j).4.2 PDE Solver OptionsFor more advanced applications, you can also specify as input arguments solver options and additional parameters that are passed to the PDE functions.options:Structure of optional parameters that change the default integration properties. This is the seventh input argument.sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan,options)See “Integrator Options” on page 10-92 for more information.Integrator OptionsThe default integration properties in the MATLAB PDE solver are selected to handle common problems. In some cases, you can improve solver performance by overriding these defaults. You do this by supplying pdepe with one or more property values in an options structure.sol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan,options)Use odeset to create the options structure. Only those options of the underlying ODE solver shown in the following table are available for pdepe.The defaults obtained by leaving off the input argument options are generally satisfactory. “Integrator Options” on page 10-9 tells you how to create the structure and describes the properties.PDE Properties§4 Examples•“Single PDE” on page 10-93•“System of PDEs” on page 10-98•“Additional Examples” on page 10-1031.Single PDE• “Solving the Equation” on page 10-93• “Evaluating the Solution” on page 10-98Solving the Equation. This example illustrates the straightforward formulation, solution, and plotting of the solution of a single PDE222x u t u ∂∂=∂∂π This equation holds on an interval 1x 0≤≤ for times t ≥ 0. At 0t = the solution satisfies the initial condition x sin )0,x (u π=.At 0x =and 1x = , the solution satisfies the boundary conditions0)t ,1(xu e ,0)t ,0(u t =∂∂+π=- Note The demo pdex1 contains the complete code for this example. The demo uses subfunctions to place all functions it requires in a single MATLAB file.To run the demo type pdex1 at the command line. See “PDE Solver Syntax” on page 10-89 for more information. 1 Rewrite the PDE. Write the PDE in the form)xu ,u ,t ,x (s ))x u ,u ,t ,x (f x (x x t u )x u ,u ,t ,x (c m m ∂∂+∂∂∂∂=∂∂∂∂- This is the form shown in Equation 10-2 and expected by pdepe. For this example, the resulting equation is0x u x x x t u 002+⎪⎭⎫ ⎝⎛∂∂∂∂=∂∂π with parameter and the terms 0m = and the term0s ,xu f ,c 2=∂∂=π=2 Code the PDE. Once you rewrite the PDE in the form shown above (Equation 10-2) and identify the terms, you can code the PDE in a function that pdepe can use. The function must be of the form[c,f,s] = pdefun(x,t,u,dudx)where c, f, and s correspond to the f ,c and s terms. The code below computes c, f, and s for the example problem.function [c,f,s] = pdex1pde(x,t,u,DuDx)c = pi^2;f = DuDx;s = 0;3 Code the initial conditions function. You must code the initial conditions in a function of the formu = icfun(x)The code below represents the initial conditions in the function pdex1ic. Partial Differential Equationsfunction u0 = pdex1ic(x)u0 = sin(pi*x);4 Code the boundary conditions function. You must also code the boundary conditions in a function of the form[pl,ql,pr,qr] = bcfun(xl,ul,xr,ur,t) The boundary conditions, written in the same form as Equation 10-4, are0x ,0)t ,0(x u .0)t ,0(u ==∂∂+and1x ,0)t ,1(xu .1e t ==∂∂+π- The code below evaluates the components )u ,t ,x (p and )u ,t ,x (q of the boundary conditions in the function pdex1bc.function [pl,ql,pr,qr] = pdex1bc(xl,ul,xr,ur,t)pl = ul;ql = 0;pr = pi * exp(-t);qr = 1;In the function pdex1bc, pl and ql correspond to the left boundary conditions (x=0 ), and pr and qr correspond to the right boundary condition (x=1).5 Select mesh points for the solution. Before you use the MATLAB PDE solver, you need to specify the mesh points at which you want pdepe to evaluate the solution. Specify the points as vectors t and x.The vectors t and x play different roles in the solver (see “PDE Solver” on page 10-89). In particular, the cost and the accuracy of the solution depend strongly on the length of the vector x. However, the computation is much less sensitive to the values in the vector t.10 CalculusThis example requests the solution on the mesh produced by 20 equally spaced points from the spatial interval [0,1] and five values of t from thetime interval [0,2].x = linspace(0,1,20);t = linspace(0,2,5);6 Apply the PDE solver. The example calls pdepe with m = 0, the functions pdex1pde, pdex1ic, and pdex1bc, and the mesh defined by x and t at which pdepe is to evaluate the solution. The pdepe function returns the numerical solution in a three-dimensional array sol, wheresol(i,j,k) approximates the kth component of the solution,u, evaluated atkt(i) and x(j).m = 0;sol = pdepe(m,@pdex1pde,@pdex1ic,@pdex1bc,x,t);This example uses @ to pass pdex1pde, pdex1ic, and pdex1bc as function handles to pdepe.Note See the function_handle (@), func2str, and str2func reference pages, and the @ section of MATLAB Programming Fundamentals for information about function handles.7 View the results. Complete the example by displaying the results:a Extract and display the first solution component. In this example, the solution has only one component, but for illustrative purposes, the example “extracts” it from the three-dimensional array. The surface plot shows the behavior of the solution.u = sol(:,:,1);surf(x,t,u)title('Numerical solution computed with 20 mesh points') xlabel('Distance x') ylabel('Time t')Distance xNumerical solution computed with 20 mesh points.Time tb Display a solution profile at f t , the final value of . In this example,2t t f ==.figure plot(x,u(end,:)) title('Solution at t = 2') xlabel('Distance x') ylabel('u(x,2)')Solutions at t = 2.Distance xu (x ,2)Evaluating the Solution. After obtaining and plotting the solution above, you might be interested in a solution profile for a particular value of t, or the time changes of the solution at a particular point x. The kth column u(:,k)(of the solution extracted in step 7) contains the time history of the solution at x(k). The jth row u(j,:) contains the solution profile at t(j). Using the vectors x and u(j,:), and the helper function pdeval, you can evaluate the solution u and its derivative at any set of points xout [uout,DuoutDx] = pdeval(m,x,u(j,:),xout)The example pdex3 uses pdeval to evaluate the derivative of the solution at xout = 0. See pdeval for details.2. System of PDEsThis example illustrates the solution of a system of partial differential equations. The problem is taken from electrodynamics. It has boundary layers at both ends of the interval, and the solution changes rapidly for small . The PDEs are)u u (F xu017.0t u )u u (F xu 024.0t u 212222212121-+∂∂=∂∂--∂∂=∂∂ where )y 46.11exp()y 73.5exp()y (F --=. The equations hold on an interval1x 0≤≤ for times 0t ≥.The solution satisfies the initial conditions0)0,x (u ,1)0,x (u 21≡≡and boundary conditions0)t ,1(xu,0)t ,1(u ,0)t ,0(u ,0)t ,0(x u 2121=∂∂===∂∂ Note The demo pdex4 contains the complete code for this example. The demo uses subfunctions to place all required functions in a single MATLAB file. To run this example type pdex4 at the command line.1 Rewrite the PDE. In the form expected by pdepe, the equations are⎥⎦⎤⎢⎣⎡---+⎥⎦⎤⎢⎣⎡∂∂∂∂∂∂=⎥⎦⎤⎢⎣⎡∂∂⎥⎦⎤⎢⎣⎡)u u (F )u u (F )x /u 170.0)x /u (024.0x u u t *.1121212121 The boundary conditions on the partial derivatives of have to be written in terms of the flux. In the form expected by pdepe, the left boundary condition is⎥⎦⎤⎢⎣⎡=⎥⎦⎤⎢⎣⎡∂∂∂∂⎥⎦⎤⎢⎣⎡+⎥⎦⎤⎢⎣⎡00)x /u (170.0)x /u (024.0*.01u 0212and the right boundary condition is⎥⎦⎤⎢⎣⎡=⎥⎦⎤⎢⎣⎡∂∂∂∂⎥⎦⎤⎢⎣⎡+⎥⎦⎤⎢⎣⎡-00)x /u (170.0)x /u (024.0*.1001u 2112 Code the PDE. After you rewrite the PDE in the form shown above, you can code it as a function that pdepe can use. The function must be of the form[c,f,s] = pdefun(x,t,u,dudx)where c, f, and s correspond to the , , and terms in Equation 10-2. function [c,f,s] = pdex4pde(x,t,u,DuDx)c = [1; 1];f = [0.024; 0.17] .* DuDx;y = u(1) - u(2);F = exp(5.73*y)-exp(-11.47*y);s = [-F; F];3 Code the initial conditions function. The initial conditions function must be of the formu = icfun(x)The code below represents the initial conditions in the function pdex4ic. function u0 = pdex4ic(x);u0 = [1; 0];4 Code the boundary conditions function. The boundary conditions functions must be of the form[pl,ql,pr,qr] = bcfun(xl,ul,xr,ur,t)The code below evaluates the components p(x,t,u) and q(x,t) (Equation 10-4) of the boundary conditions in the function pdex4bc.function [pl,ql,pr,qr] = pdex4bc(xl,ul,xr,ur,t)pl = [0; ul(2)];ql = [1; 0];pr = [ur(1)-1; 0];qr = [0; 1];5 Select mesh points for the solution. The solution changes rapidly for small t . The program selects the step size in time to resolve this sharp change, but to see this behavior in the plots, output times must be selected accordingly. There are boundary layers in the solution at both ends of [0,1], so mesh points must be placed there to resolve these sharp changes. Often some experimentation is needed to select the mesh that reveals the behavior of the solution.x = [0 0.005 0.01 0.05 0.1 0.2 0.5 0.7 0.9 0.95 0.99 0.995 1];t = [0 0.005 0.01 0.05 0.1 0.5 1 1.5 2];6 Apply the PDE solver. The example calls pdepe with m = 0, the functions pdex4pde, pdex4ic, and pdex4bc, and the mesh defined by x and t at which pdepe is to evaluate the solution. The pdepe function returns the numerical solution in a three-dimensional array sol, wheresol(i,j,k) approximates the kth component of the solution, μk, evaluated at t(i) and x(j).m = 0;sol = pdepe(m,@pdex4pde,@pdex4ic,@pdex4bc,x,t);7 View the results. The surface plots show the behavior of the solution components. u1 = sol(:,:,1); u2 = sol(:,:,2); figure surf(x,t,u1) title('u1(x,t)') xlabel('Distance x') 其输出图形为Distance xu1(x,t)Time tfigure surf(x,t,u2) title('u2(x,t)') xlabel('Distance x') ylabel('Time t')Distance xu2(x,t)Time tAdditional ExamplesThe following additional examples are available. Type edit examplename to view the code and examplename to run the example.。
偏微分方程的MATLAB求解精讲©MA TLAB求解微分/偏微分方程,一直是一个头大的问题,两个字,“难过”,由于MA TLAB对LaTeX的支持有限,所有方程必须化成MA TLAB可接受的标准形式,不支持像其他三个数学软件那样直接傻瓜式输入,这个真把人给累坏了!不抱怨了,还是言归正传,回归我们今天的主体吧!MA TLAB提供了两种方法解决PDE问题,一是pdepe()函数,它可以求解一般的PDEs,据用较大的通用性,但只支持命令行形式调用。
二是PDE工具箱,可以求解特殊PDE问题,PDEtool有较大的局限性,比如只能求解二阶PDE问题,并且不能解决偏微分方程组,但是它提供了GUI界面,从繁杂的编程中解脱出来了,同时还可以通过File->Save As直接生成M代码一、一般偏微分方程组(PDEs)的MA TLAB求解 (3)1、pdepe函数说明 (3)2、实例讲解 (4)二、PDEtool求解特殊PDE问题 (6)1、典型偏微分方程的描述 (6)(1)椭圆型 (6)(2)抛物线型 (6)(3)双曲线型 (6)(4)特征值型 (7)2、偏微分方程边界条件的描述 (8)(1)Dirichlet条件 (8)(2)Neumann条件 (8)3、求解实例 (9)一、一般偏微分方程组(PDEs)的MATLAB 求解1、pdepe 函数说明MA TLAB 语言提供了pdepe()函数,可以直接求解一般偏微分方程(组),它的调用格式为sol=pdepe(m,@pdefun,@pdeic,@pdebc,x,t)【输入参数】@pdefun :是PDE 的问题描述函数,它必须换成下面的标准形式(,,)[(,,,)](,,,)()m m u u u uc x t x x f x t u s x t u x t x x x−∂∂∂∂∂=+∂∂∂∂∂式1 这样,PDE 就可以编写下面的入口函数 [c,f,s]=pdefun(x,t,u,du)m,x,t 就是对应于(式1)中相关参数,du 是u 的一阶导数,由给定的输入变量即可表示出出c,f,s 这三个函数@pdebc :是PDE 的边界条件描述函数,必须先化为下面的形式(,,)(,,).*(,,,)0up x t u q x t u f x t u x∂+=∂ 于是边值条件可以编写下面函数描述为 [pa,qa,pb,qb]=pdebc(x,t,u,du)其中a 表示下边界,b 表示下边界@pdeic :是PDE 的初值条件,必须化为下面的形式00(,)u x t u =我们使用下面的简单的函数来描述为 u0=pdeic(x)m,x,t :就是对应于(式1)中相关参数【输出参数】sol :是一个三维数组,sol(:,:,i)表示u i 的解,换句话说u k 对应x(i)和t(j)时的解为sol(i,j,k)通过sol ,我们可以使用pdeval()直接计算某个点的函数值2、实例讲解试求解下面的偏微分2111222221220.024()0.17()u u F u u t xu u F u u tx ∂∂=−− ∂∂ ∂∂ =−− ∂∂ 其中, 5.7311.46()x x F x e e −=−,且满足初始条件12(,0)1,(,0)0u x u x ==及边界条件1221(0,)0,(0,)0,(1,)1,(1,)0u ut u t u t t x x∂∂====∂∂【解】(1)对照给出的偏微分方程,根据标注形式,则原方程可以改写为111222120.024()1.*1()0.17u u F u u x u u F u u t t x ∂−−∂∂∂=+ ∂−∂∂∂可见m=0,且1122120.024()1,,1()0.17u F u u x c f s u F u u x ∂−− ∂===∂−∂%% 目标PDE 函数function [c,f,s]=pdefun (x,t,u,du) c=[1;1];f=[0.024*du(1);0.17*du(2)]; temp=u(1)-u(2);s=[-1;1].*(exp(5.73*temp)-exp(-11.46*temp));(2)边界条件改写为12011010.*.*00000u f f u −+=+=下边界上边界%% 边界条件函数function [pa,qa,pb,qb]=pdebc(xa,ua,xb,ub,t) %a 表示下边界,b 表示上边界 pa=[0;ua(2)];qa=[1;0]; pb=[ub(1)-1;0]; qb=[0;1];(3)初值条件改写为1210u u =%% 初值条件函数function u0=pdeic(x) u0=[1;0];(4)最后编写主调函数 clcx=0:0.05:1; t=0:0.05:2; m=0;sol=pdepe(m,@pdefun,@pdeic,@pdebc,x,t);figure('numbertitle','off','name','PDE Demo ——by Matlabsky') subplot(211)surf(x,t,sol(:,:,1)) title('The Solution of u_1') xlabel('X') ylabel('T') zlabel('U') subplot(212)surf(x,t,sol(:,:,2)) title('The Solution of u_2') xlabel('X') ylabel('T') zlabel('U')二、PDEtool 求解特殊PDE 问题MATLAB 的偏微分工具箱(PDE toolbox)可以比较规范的求解各种常见的二阶偏微分方程,但是惋惜的是只能求解特殊二阶的PDE 问题,并且不支持偏微分方程组!PDE toolbox 支持命令行形式求解PDE 问题,但是要记住那些命令以及调用形式真的很累人,还好MATLAB 提供了GUI 可视交互界面pdetool ,在pdetool 中可以很方便的求解一个PDE 问题,并且可以帮我们直接生成M 代码(File->Save As)。
matlab 求解偏微分方程求解偏微分方程是数学中的一种重要问题,而MATLAB是一种功能强大的数学软件,可以用于求解各种数学问题,包括偏微分方程。
在本文中,我们将探讨如何使用MATLAB求解偏微分方程,并介绍一些常用的求解方法和技巧。
偏微分方程是描述自然界中许多现象的数学模型,例如热传导、扩散、波动等。
求解偏微分方程的目标是找到满足方程条件的未知函数。
MATLAB提供了许多内置函数和工具箱,可以方便地进行偏微分方程的求解。
我们需要定义偏微分方程的边界条件和初始条件。
边界条件是指在求解区域的边界上给定的条件,而初始条件是指在求解区域内给定的初始状态。
这些条件将帮助我们确定偏微分方程的解。
接下来,我们可以使用MATLAB中的偏微分方程求解函数来求解方程。
MATLAB提供了几种常用的求解方法,例如有限差分法、有限元法和谱方法等。
通过选择合适的求解方法,我们可以得到偏微分方程的数值解。
在使用MATLAB求解偏微分方程时,我们还可以使用一些技巧和优化方法来提高求解效率。
例如,可以使用网格剖分方法来将求解区域划分为若干小区域,从而减少计算量。
此外,还可以使用迭代方法来逐步逼近偏微分方程的解,从而提高求解精度。
除了求解偏微分方程,MATLAB还可以用于可视化偏微分方程的解。
通过使用MATLAB的绘图函数,我们可以将数值解以图形的形式展示出来,从而更直观地理解偏微分方程的解。
需要注意的是,在使用MATLAB求解偏微分方程时,我们需要考虑计算资源的限制。
由于偏微分方程的求解通常需要大量的计算和存储资源,因此我们需要合理安排计算机的内存和处理器的使用,以避免计算过程中的错误或崩溃。
总结起来,MATLAB是一种强大的数学软件,可以用于求解偏微分方程。
通过选择合适的求解方法和优化技巧,我们可以得到偏微分方程的数值解,并用图形的形式展示出来。
使用MATLAB求解偏微分方程可以帮助我们更好地理解和解决实际问题,是数学研究和工程应用中的重要工具之一。
MATLAB偏微分方程组求解介绍偏微分方程组是描述自然界中许多现象的数学模型,包括流体力学、电磁学、热传导等。
求解偏微分方程组是科学研究和工程应用中的重要问题之一。
MATLAB作为一种强大的数值计算工具,提供了丰富的函数和工具箱,可以用于求解偏微分方程组。
本文将介绍如何使用MATLAB求解偏微分方程组。
我们将从基本的概念和数学理论开始,然后介绍MATLAB中的相关函数和工具箱,最后给出一个具体的求解偏微分方程组的示例。
基本概念和数学理论偏微分方程组偏微分方程组是一个包含多个未知函数的方程组,其中每个未知函数的导数(偏导数)出现在方程中。
一般形式的偏微分方程组可以写成以下形式:F1(u1,u2,…,u n,∂u1∂x,∂u2∂x,…,∂u n∂x,∂u1∂y,∂u2∂y,…,∂u n∂y,…)=0F2(u1,u2,…,u n,∂u1∂x,∂u2∂x,…,∂u n∂x,∂u1∂y,∂u2∂y,…,∂u n∂y,…)=0⋮F m(u1,u2,…,u n,∂u1∂x,∂u2∂x,…,∂u n∂x,∂u1∂y,∂u2∂y,…,∂u n∂y,…)=0其中,u1,u2,…,u n是未知函数,∂u1∂x ,∂u2∂x,…,∂u n∂x,∂u1∂y,∂u2∂y,…,∂u n∂y,…是未知函数的偏导数。
边界条件为了求解偏微分方程组,我们需要给出适当的边界条件。
边界条件是在给定的边界上给出未知函数或其导数的值。
常见的边界条件包括:Dirichlet边界条件、Neumann边界条件和Robin边界条件。
•Dirichlet边界条件:给定未知函数在边界上的值。
•Neumann边界条件:给定未知函数的法向导数在边界上的值。
•Robin边界条件:给定未知函数和其法向导数的线性组合在边界上的值。
数值方法由于一般情况下无法找到偏微分方程组的解析解,我们需要使用数值方法来求解。
常见的数值方法包括:有限差分法、有限元法和谱方法。
•有限差分法:将偏微分方程组转化为差分方程组,通过在网格上逼近导数来近似原方程。
偏微分方程的matlab解法pdf
MATLAB提供了多种用于求解偏微分方程的方法,包括:
1.有限差分法:将空间域离散化,用有限差分代替偏导数。
2.有限元法:将空间域划分为有限元,用有限元的代数方程代替偏微分方程。
3.有限体积法:将空间域划分为有限体积,用有限体积的积分代替偏微分方程。
4.谱方法:利用空间域的正交函数来求解偏微分方程。
5.变分法:将偏微分方程转化为变分问题,然后用数值方法求解变分问题。
具体选择哪种方法,需要根据偏微分方程的类型、边界条件和初始条件等因素来决定。
文章标题:深入探讨 Matlab 中求解偏微分方程的方法和应用一、引言在现代科学和工程中,偏微分方程是一种重要的数学工具,用于描述各种自然现象和物理过程,如热传导、流体力学、电磁场等。
Matlab 是一个用于科学计算和工程应用的强大工具,提供了丰富的数值计算和数据可视化功能,其中包括求解偏微分方程的工具箱,本文将深入探讨在Matlab中求解偏微分方程的方法和应用。
二、基本概念偏微分方程(Partial Differential Equation, PDE)是关于多个变量的函数及其偏导数的方程。
在物理学和工程学中,PDE广泛应用于描述空间变量和时间变量之间的关系。
在Matlab中,求解PDE通常涉及到确定PDE类型、边界条件、初始条件和求解方法等步骤。
三、求解方法1. 有限差分法(Finite Difference Method)有限差分法是求解PDE的常用数值方法之一,它将PDE转化为差分方程组,并通过迭代求解得到数值解。
在Matlab中,可以使用pdepe 函数来求解具有一维、二维或三维空间变量的PDE,该函数可以直接处理边界条件和初始条件。
2. 有限元法(Finite Element Method)有限元法是另一种常用的数值方法,它将求解区域离散化为有限数量的单元,并通过单元之间的插值来逼近PDE的解。
Matlab提供了pdenonlin函数来求解非线性PDE,该函数支持各种复杂的几何形状和非线性材料参数。
3. 特征线法(Method of Characteristics)特征线法适用于一维双曲型PDE的求解,该方法基于特征线方程的性质来构造数值解。
在Matlab中,可以使用pdegplot函数来展示特征线,并通过构造特征线网格来求解PDE。
四、实际应用1. 热传导方程的求解假设我们需要求解一个长条形的材料中的热传导方程,可以通过在Matlab中定义边界条件和初始条件,然后使用pdepe函数来求解得到温度分布和热流线。
MATLAB学习(序列1)偏微分方程组的求解ode23 解非刚性微分方程,低精度,使用Runge-Kutta法的二三阶算法。
ode45 解非刚性微分方程,中等精度,使用Runge-Kutta法的四五阶算法。
ode113 解非刚性微分方程,变精度变阶次Adams-Bashforth-Moulton PECE算法。
ode23t 解中等刚性微分方程,使用自由内插法的梯形法则。
ode15s 解刚性微分方程,使用可变阶次的数值微分(NDFs)算法。
ode23s 解刚性微分方程,低阶方法,使用修正的Rosenbrock公式。
ode23tb 解刚性微分方程,低阶方法,使用TR-BDF2方法,即Runger-Kutta公式的第一级采用梯形法则,第二级采用Gear法。
[t,YY]=solver('F',tspan,Yo解算ODE初值问题的最简调用格式。
solver指上面的指令。
tspan=[0,30]; %时域t的范围y0=[1;0]; %y(1)y(2的初始值[tt,yy]=ode45(@DyDt,tspan,y0;plot(tt,yy(:,1,title('x(t'function ydot=DyDt(t,yydot=[y(2; 2*(1-y(1^2*y(2-y(1]刚性方程:刚性是指其Jacobian矩阵的特征值相差十分悬殊。
在解的性态上表现为,其中一些解变化缓慢,另一些变化快,且相差较悬殊,这类方程常常称为刚性方程,又称为Stiff方程。
刚性方程和非刚性方程对解法中步长选择的要求不同。
刚性方程一般不适合由ode45这类函数求解,而应该采用ode15s等。
如果不能分辨是否是刚性方程,先试用ode45,再用ode15s。
[t,YY,Te,Ye,Ie] = solver('F',tspan,Yo,options,p1,p2,…解算ODE初值问题的最完整调用格式。
为了能够解出方程,要用指令odeset确定求解的条件和要求。
Matlab 差分法解偏微分方程1.引言解偏微分方程是数学和工程领域中的一项重要课题,它在科学研究和工程实践中具有广泛的应用。
而 Matlab 差分法是一种常用的数值方法,用于求解偏微分方程。
本文将介绍 Matlab 差分法在解偏微分方程中的应用,包括原理、步骤和实例。
2. Matlab 差分法原理差分法是一种离散化求解微分方程的方法,通过近似替代微分项来求解微分方程的数值解。
在 Matlab 中,差分法可以通过有限差分法或者差分格式来实现。
有限差分法将微分方程中的导数用有限差分替代,而差分格式指的是使用不同的差分格式来近似微分方程中的各个项,通常包括前向差分、后向差分和中心差分等。
3. Matlab 差分法步骤使用 Matlab 差分法解偏微分方程一般包括以下步骤:(1)建立离散化的区域:将求解区域离散化为网格点或节点,并确定网格间距。
(2)建立离散化的时间步长:对于时间相关的偏微分方程,需要建立离散化的时间步长。
(3)建立离散化的微分方程:使用差分法将偏微分方程中的微分项转化为离散形式。
(4)建立迭代方程:根据离散化的微分方程建立迭代方程,求解数值解。
(5)编写 Matlab 代码:根据建立的迭代方程编写 Matlab 代码求解数值解。
(6)求解并分析结果:使用 Matlab 对建立的代码进行求解,并对结果进行分析和后处理。
4. Matlab 差分法解偏微分方程实例假设我们要使用 Matlab 差分法解决以下一维热传导方程:∂u/∂t = α * ∂^2u/∂x^2其中 u(x, t) 是热传导方程的温度分布,α 是热扩散系数。
4.1. 离散化区域和时间步长我们将求解区域离散化为网格点,分别为 x_i,i=1,2,...,N。
时间步长为Δt。
4.2. 离散化的微分方程使用中心差分格式将偏微分方程中的导数项离散化得到:∂u/∂t ≈ (u_i(t+Δt) - u_i(t))/Δt∂^2u/∂x^2 ≈ (u_i-1(t) - 2u_i(t) + u_i+1(t))/(Δx)^2代入原偏微分方程可得离散化的微分方程:(u_i(t+Δt) - u_i(t))/Δt = α * (u_i-1(t) - 2u_i(t) + u_i+1(t))/(Δx)^24.3. 建立迭代方程根据离散化的微分方程建立迭代方程:u_i(t+Δt) = u_i(t) + α * Δt * (u_i-1(t) - 2u_i(t) + u_i+1(t))/(Δx)^24.4. 编写 Matlab 代码使用以上建立的迭代方程编写 Matlab 代码求解热传导方程。
matlab用有限元法求解偏微分方程组使用有限元法求解偏微分方程组是一种常见的数值计算方法,它在工程领域和科学研究中广泛应用。
本文将介绍如何利用MATLAB软件进行有限元法求解偏微分方程组的基本步骤和注意事项。
我们需要了解有限元法的基本原理。
有限元法是一种将连续问题离散化为有限个小区域,通过在每个小区域内建立适当的数学模型,然后将这些小区域连接起来形成整个问题的数学模型的方法。
在有限元法中,我们通常将问题的域分割成许多小的有限元,每个有限元都具有简单的几何形状,如线段、三角形或四边形。
然后,在每个有限元上建立适当的近似函数,通过对这些函数的系数进行求解,我们可以得到问题的近似解。
在MATLAB中,有限元法的求解过程可以分为以下几个步骤:1. 离散化域:根据问题的几何形状,将问题的域进行离散化处理。
离散化可以采用三角剖分法或四边形剖分法,将域分割成许多小的有限元。
2. 建立数学模型:在每个有限元上建立适当的数学模型。
这通常涉及选择适当的近似函数,并在每个有限元上求解这些函数的系数。
3. 组装方程:将每个有限元上的数学模型组装成整个问题的数学模型。
这涉及到将有限元之间的边界条件进行匹配,并建立整个问题的刚度矩阵和载荷向量。
4. 求解方程:利用线性代数求解方法,求解得到问题的近似解。
MATLAB提供了各种求解线性方程组的函数,如“\”运算符、LU 分解和共轭梯度法等。
5. 后处理:对求解结果进行后处理,包括绘制解的图形、计算问题的误差等。
在进行有限元法求解偏微分方程组时,需要注意以下几点:1. 网格剖分的合理性:网格剖分的精细程度对结果的精确性有很大影响。
网格过于粗糙可能导致结果的不准确,而网格过于细小则会增加计算的复杂性。
因此,需要根据问题的特点和计算资源的限制选择合适的网格剖分。
2. 近似函数的选择:近似函数的选择直接影响到结果的准确性和计算的效率。
一般情况下,近似函数的阶数越高,结果的准确性越高,但计算的复杂性也越大。
使用matlab差分法解偏微分方程1. 引言差分法是一种常用的数值方法,用于求解偏微分方程(Partial Differential Equations,简称PDE)的数值解。
在工程学和科学研究中,PDE广泛应用于描述各种物理现象和过程。
本文将介绍使用MATLAB差分法来解偏微分方程的方法和步骤,并探讨其优势和局限性。
2. 差分法简介差分法是一种基于离散点的数值求解方法,它将连续的空间或时间变量离散化为有限个点,通过对这些离散点上的方程进行逼近,得到PDE的数值解。
其中,MATLAB作为一种功能强大的数值计算工具,提供了快速而高效的差分法求解PDE的功能。
3. 二阶偏微分方程的差分方法在本节中,我们将以一个简单的二阶偏微分方程为例,说明如何使用差分法来解决。
考虑一个二维的泊松方程,即:∂²u/∂x² + ∂²u/∂y² = f(x, y)其中,u是未知函数,f(x, y)是已知函数。
为了使用差分法求解该方程,我们需要将空间离散化,假设网格步长为Δx和Δy。
我们可以使用中心差分法来逼近二阶导数,从而将偏微分方程转化为一个代数方程组。
在MATLAB中,我们可以通过设置好网格步长和边界条件,构建对应的代数方程组,并使用线性代数求解方法(如直接解法或迭代解法)获得数值解。
4. 差分法的优势和局限性差分法作为一种数值方法,具有许多优势和应用范围,但也存在一些局限性。
优势:- 简单易懂:差分法的思想直观明了,易于理解和实现。
- 适应性广泛:差分法可以用于求解各种类型的偏微分方程,包括常微分方程和偏微分方程。
- 准确度可控:通过调整网格步长,可以控制数值解的精度和稳定性。
局限性:- 离散误差:当空间或时间步长过大时,差分法的数值解可能会出现较大的离散误差。
- 边界条件:合适的边界条件对于差分法的求解结果至关重要,不合理的边界条件可能导致数值解的不准确。
- 计算效率:对于复杂的偏微分方程,差分法的计算成本可能较高,需要耗费大量的计算资源和时间。
matlab 求解偏微分方程使用MATLAB求解偏微分方程摘要:偏微分方程(partial differential equation, PDE)是数学中重要的一类方程,广泛应用于物理、工程、经济、生物等领域。
MATLAB 是一种强大的数值计算软件,提供了丰富的工具箱和函数,可以用来求解各种类型的偏微分方程。
本文将介绍如何使用MATLAB来求解偏微分方程,并通过具体案例进行演示。
引言:偏微分方程是描述多变量函数的方程,其中包含了函数的偏导数。
一般来说,偏微分方程可以分为椭圆型方程、双曲型方程和抛物型方程三类。
求解偏微分方程的方法有很多,其中数值方法是最常用的一种。
MATLAB作为一种强大的数值计算软件,提供了丰富的工具箱和函数,可以用来求解各种类型的偏微分方程。
方法:MATLAB提供了多种求解偏微分方程的函数和工具箱,包括pdepe、pdetoolbox和pde模块等。
其中,pdepe函数是用来求解带有初始条件和边界条件的常微分方程组的函数,可以用来求解一维和二维的偏微分方程。
pdepe函数使用有限差分法或有限元法来离散化偏微分方程,然后通过求解离散化后的常微分方程组得到最终的解。
案例演示:考虑一维热传导方程的求解,偏微分方程为:∂u/∂t = α * ∂^2u/∂x^2其中,u(x,t)是温度分布函数,α是热扩散系数。
假设初始条件为u(x,0)=sin(pi*x),边界条件为u(0,t)=0和u(1,t)=0。
我们需要定义偏微分方程和边界条件。
在MATLAB中,可以使用匿名函数来定义偏微分方程和边界条件。
然后,我们使用pdepe函数求解偏微分方程。
```matlabfunction [c,f,s] = pde(x,t,u,DuDx)c = 1;f = DuDx;s = 0;endfunction u0 = uinitial(x)u0 = sin(pi*x);endfunction [pl,ql,pr,qr] = uboundary(xl,ul,xr,ur,t)pl = ul;ql = 0;pr = ur;qr = 0;endx = linspace(0,1,100);t = linspace(0,0.1,10);m = 0;sol = pdepe(m,@pde,@uinitial,@uboundary,x,t);u = sol(:,:,1);surf(x,t,u);xlabel('Distance x');ylabel('Time t');zlabel('Temperature u');```在上述代码中,我们首先定义了偏微分方程函数pde,其中c、f和s分别表示系数c、f和s。
matlab解偏微分方程组使用Matlab解偏微分方程组在科学与工程领域,偏微分方程组是描述自然现象和物理过程的重要数学工具。
解偏微分方程组是求解这些现象和过程的数值模拟方法之一。
Matlab作为一种高级的数值计算软件,提供了强大的功能来解决偏微分方程组。
本文将介绍如何使用Matlab来解偏微分方程组,并给出实例说明。
一、Matlab解偏微分方程组的基本原理Matlab是一种基于矩阵运算的高级数值计算软件,它提供了丰富的函数和工具箱来解决数学问题。
在解偏微分方程组时,Matlab主要采用有限差分法、有限元法和谱方法等数值方法。
这些方法将偏微分方程转化为离散的代数方程组,然后通过求解代数方程组得到数值解。
二、使用Matlab解偏微分方程组的步骤1. 定义偏微分方程组:首先需要将偏微分方程组转化为Matlab可以处理的形式。
通常将自变量和因变量离散化,并用矩阵和向量表示。
2. 离散化:将偏微分方程中的连续变量转化为离散变量,通常采用有限差分法或有限元法。
有限差分法将偏微分方程中的导数用差商表示,有限元法则将区域划分为有限个小单元。
3. 构建代数方程组:根据离散化后的方程,可以得到相应的代数方程组。
这一步需要根据边界条件和初始条件来确定代数方程的边界值和初始值。
4. 求解代数方程组:利用Matlab提供的求解函数,如\texttt{fsolve}或\texttt{ode45}等,求解代数方程组得到数值解。
5. 可视化结果:使用Matlab的绘图函数,如\texttt{plot}或\texttt{surf}等,将数值解可视化展示出来。
这可以帮助我们更好地理解解的特性和趋势。
三、一个简单的例子为了更好地理解如何使用Matlab解偏微分方程组,我们将以一个简单的热传导问题为例。
考虑一个一维热传导方程:$$\frac{{\partial u}}{{\partial t}} = \frac{{\partial^2 u}}{{\partial x^2}}$$其中$u(x,t)$是温度分布,$x$是空间变量,$t$是时间变量。
引言偏微分方程定解问题有着广泛的应用背景。
人们用偏微分方程来描述、解释或者预见各种自然现象,并用于科学和工程技术的各个领域fll。
然而,对于广大应用工作者来说,从偏微分方程模型出发,使用有限元法或有限差分法求解都要耗费很大的工作量,才能得到数值解。
现在,MATLAB PDEToolbox已实现对于空间二维问题高速、准确的求解过程。
偏微分方程如果一个微分方程中出现的未知函数只含一个自变量,这个方程叫做常微分方程,也简称微分方程;如果一个微分方程中出现多元函数的偏导数,或者说如果未知函数和几个变量有关,而且方程中出现未知函数对几个变量的导数,那么这种微分方程就是偏微分方程。
常用的方法有变分法和有限差分法。
变分法是把定解问题转化成变分问题,再求变分问题的近似解;有限差分法是把定解问题转化成代数方程,然后用计算机进行计算;还有一种更有意义的模拟法,它用另一个物理的问题实验研究来代替所研究某个物理问题的定解。
虽然物理现象本质不同,但是抽象地表示在数学上是同一个定解问题,如研究某个不规则形状的物体里的稳定温度分布问题,由于求解比较困难,可作相应的静电场或稳恒电流场实验研究,测定场中各处的电势,从而也解决了所研究的稳定温度场中的温度分布问题。
随着物理科学所研究的现象在广度和深度两方面的扩展,偏微分方程的应用范围更广泛。
从数学自身的角度看,偏微分方程的求解促使数学在函数论、变分法、级数展开、常微分方程、代数、微分几何等各方面进行发展。
从这个角度说,偏微分方程变成了数学的中心。
一、MATLAB方法简介及应用1.1 MATLAB简介MATLAB是美国MathWorks公司出品的商业数学软件,用于算法开发、数据可视化、数据分析以及数值计算的高级技术计算语言和交互式环境,主要包括MATLAB和Simulink两大部分。
1.2 Matlab主要功能数值分析数值和符号计算工程与科学绘图控制系统的设计与仿真数字图像处理数字信号处理通讯系统设计与仿真财务与金融工程1.3 优势特点1) 高效的数值计算及符号计算功能,能使用户从繁杂的数学运算分析中解脱出来;2) 具有完备的图形处理功能,实现计算结果和编程的可视化;3) 友好的用户界面及接近数学表达式的自然化语言,使学者易于学习和掌握;4) 功能丰富的应用工具箱(如信号处理工具箱、通信工具箱等) ,为用户提供了大量方便实用的处理工具。
1.4 MATLAB 产品族可以用来进行以下各种工作●数值分析 ●数值和符号计算 ●工程与科学绘图 ●控制系统的设计与仿真 ●数字图像处理技术 ●数字信号处理技术 ●通讯系统设计与仿真 ●财务与金融工程●管理与调度优化计算(运筹学)MATLAB 的应用范围非常广,包括信号和图像处理、通讯 MATLAB 在通讯系统设计与仿真的应用、控制系统设计、测试和测量、财务建模和分析以及计算生物学等众多应用领域。
附加的工具箱(单独提供的专用MATLAB 函数集)扩展了MATLAB 环境,以解决这些应用领域内特定类型的问题。
二、Laplacian 算子简介Laplacian 算子:2222yx ∂∂+∂∂=∆Poisson 方程(ellipptic):fu =∆Laplacian 算子的特征值问题:fu u =+∆λHeat equation(parabolic):ut u∆=∂∂Wave equation(hyperbolic):u t u ∆=∂∂22五点离散:hhh y x u y x u h y x u y h x u y x u y h x u y x u h 22),(),(2),(),(),(2),(),(-+-++-+-+=∆hP u S u E u W u N u p u h2)(4)()()()()(-+++=∆)(=∆P u hPoisson 方程离散:)()(P f P u h=∆特征值问题:)()(P P u u k k h λ=-∆热方程:),(),,(),,(2y x u t y x u t y x u h ∆=--+δδδ),(),,(),,(y x u t y x u t y x u h ∆+=+δδ波动方程:),(),,(),,(2),,(2y x u t y x u t y x u t y x u h ∆=-+-+δδδ),(),,(),,(y x u t y x u t y x u h ∆+=+δδ波动方程:),(),,(),,(2),,(2y x u t y x u t y x u t y x u h ∆=-+-+δδδ),(),,(),,(2),,(2y x u t y x u t y x u t y x u h ∆+--=+δδδδ椭圆方程:b Au =),(2y x f h b = 特征值方程:u Au h λ=-21热方程: )()()1(n n n Auuuσ+=+2h δσ=波动方程:)()1()()1(2n n n n Au u u u σ+-=-+2⎥⎦⎤⎢⎣⎡=h δδ)()1(n n Mu u =+热方程: )()1(n n Mu u =+ A I M σ+=波动方程:212≤⎥⎦⎤⎢⎣⎡=h δσ三、Matlab 解偏微分方程解偏微分方程不是一件轻松的事情,但是偏微分方程在自然科学和工程领域应用很广,因此,研究解偏微分方程的方法、开发解偏微分方程的工具是数学和计算机领域中的一项重要工作。
MATLAB提供了专门用于解二维偏微分方程的工具箱,使用这个工具箱,一方面解偏微分方程,另一方面,可以让我们学习如何把求解数学问题的过程与方法工程化。
应当承认,我们国家在数学软件的开发方面还比较落后,MATLAB是当今世界上最好的数学软件之一,通过对这个软件的认识,有助于研发我们自己的数学软件。
MATLAB的偏微分方程工具箱名字叫pdetool,它采用有限元法解偏微分方程。
用这个工具箱可以解如下方程。
椭圆方程抛物线方程双曲线方程特征值方程所有的方程都在二维平面Ω域上。
方程中,▽是Laplace 算子,u是待解的未知函数,c, a, f是已知的实值标量函数,d 是已知的复值函数,λ是未知的特征值。
MATLAB提供了两种方法Ⅲ解决PDE问题,一是pdepe函数,它可以求解一般的PDEs,具有较大的通用性,但只支持命令行形式调用。
二是PDE工具箱,可以求解特殊PDE问题,但有较大的局限性,比如只能求解二阶PDE问题,并且不能解决偏微分方程组。
它提供了GUI界面,可以从繁杂的编程中解脱出来,同时还可以通过FileSave As直接生成M代码。
3.1函数解法pdepe函数介绍它的调用格式为sol=pdepe(m,@pdefun,@pdeic,@pdebc,x,t) 输入参数:@pdefun:是PDE的问题描述函数@pdebc:是PDE的边界条件描述函数@pdeic:是PDE的初值条件输出参数:sol:是一个三维数组,sol(:,:,i)表示ui的解,换句话说uk 对应x(i)和t(j)时的解为sol(i,j ,k),通过sol ,我们可以使用pdeval()直接计算某个点的函数值。
实例讲解例:试求解下面的偏微分其中, ,且满足初始条件及边界条件 解:(1)对照给出的偏微分方程,根据标注形式,则原方程可以改写为可见m=0,且%目标PDE 函数function[c ,f ,s]--pdefun(x ,t ,u ,du) c=[1;1];f_[0.024*du(1);0.17*du(2)];)1()(170.0)(024.0212222212121⎪⎪⎩⎪⎪⎨⎧--∂∂=∂∂--∂∂=∂∂u u F x u t u u u F x u t ux x e e x F 46.1173.5)(--=0)0,(,1)0,(21==x u x u 1),1(,0),0(,0),0(121===∂∂t u t u t t u 0),1(,2=∂∂t tu )2()()(170.0024.0*.1121212121⎥⎦⎤⎢⎣⎡---+⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡∂∂∂∂∂∂=⎥⎦⎤⎢⎣⎡∂∂⎥⎦⎤⎢⎣⎡u u F u u F x u x u t u u t ⎥⎦⎤⎢⎣⎡---=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡∂∂∂∂=⎥⎦⎤⎢⎣⎡=)()(,170.0024.0,11212121u u F u u F S x u x u f ctemp=u(1)一u(2);s=[一1;1].*(exp(5.73*temp)一exp(-11.46*temp));(2)边界条件改写为下边界上边界%边界条件函数function[pa .qa ,pb ,qb]--pdebe(xa ,ua ,xb ,ub ,t) %a 表示下边界,b 表示上边界 pa=[0;ua(2)]; qa=[1;0 ]; pb=[ub(1)一1;O]; qb=[O ;1]; (3)初值条件改写为%初值条件函数 function u0=pdeic(x) u0=[1;0];(4)最后编写主调函数 tic m=0;x=[0 0.005 0.01 0.05 0.1 0.2 0.5 0.7 0.9 0.95⎥⎦⎤⎢⎣⎡=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡∂∂∂∂⎥⎦⎤⎢⎣⎡+⎥⎦⎤⎢⎣⎡00170.0024.0*.010212x u x u u ⎥⎦⎤⎢⎣⎡=⎥⎥⎥⎦⎤⎢⎢⎢⎣⎡∂∂∂∂⎥⎦⎤⎢⎣⎡+⎥⎦⎤⎢⎣⎡-00170.0024.0*.0101211x u x u u )3(0121⎥⎦⎤⎢⎣⎡=⎥⎦⎤⎢⎣⎡u u0.99 0.995 1];T=[0 0.005 0.01 0.05 0.1 0.5 1 1.5 2];sol=pdepe(m,@pdefun,@pdeic,@pdebc,x,t);U1=sol(:,:,1);u2=sol(:,:,2);figuresurf(x,t,u1)title(’ul(x,t)’)xlabel(’Distance x’)ylabel(’Time t’)figuresurf(x,t,u2)tiffe(’u2(x,t)’)xlabel(’Distance x’)ylabel(’Time t’)结果图3.2偏微分方程的pdetool解法3.2.1偏微分方程工具箱的功能偏微分方程工具箱(PDE Toolbox)提供了研究和求解空间二维偏微分方程问题的一个强大而又灵活实用的环境。
PDE Toolbox 的功能包括:(1)设置PDE (偏微分方程)定解问题,即设置二维定解区域、边界条件以及方程的形式和系数;(2)用有限元法 (FEM) 求解PDE数值解;(3)解的可视化。
无论是高级研究人员还是初学者,在使用PDE Too1box时都会感到非常方便。
只要PDE定解问题的提法正确,那么,启动MATLAB后,在MATLAB工作空间的命令行中键人pdetool,系统立即产生偏微分方程工具箱(PDE Toolbox)的图形用户界面(Graphical User Interface,简记为GUI),即PDE解的图形环境,这时就可以在它上面画出定解区域、设置方程和边界条件、作网格剖分、求解、作图等工作,详见1.4节中的例子。