Using MATLAB (779505), страница 66
Текст из файла (страница 66)
In this example, thesolution u has only one component, but for illustrative purposes, theexample “extracts” it from the three-dimensional array. The surface plotshows 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')15-9115Differential EquationsNumerical solution computed with 20 mesh points.10.80.60.40.20211.50.810.60.40.50.2Time t00Distance x• Display a solution profile at t f , the final value of t .
In this example,t f = t = 2. See “Evaluating the Solution at Specific Points” on page 15-93for more information.figureplot(x,u(end,:))title('Solution at t = 2')xlabel('Distance x')ylabel('u(x,2)')15-92Partial Differential EquationsSolution at t = 20.140.120.1u(x,2)0.080.060.040.020−0.0200.10.20.30.40.5Distance x0.60.70.80.91Evaluating the Solution at Specific PointsAfter obtaining and plotting the solution above, you might be interested in asolution profile for a particular value of t, or the time changes of the solutionat a particular point x. The kth column u(:,k) (of the solution extracted instep 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 canevaluate the solution u and its derivative ∂u ⁄ ∂x 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 atxout = 0.
See pdeval for details.Improving PDE Solver PerformanceThe default integration properties in the MATLAB PDE solver are selected tohandle common problems. In some cases, you can improve solver performanceby overriding these defaults. You do this by supplying pdepe with one or moreproperty values in an options structure.15-9315Differential Equationssol = pdepe(m,pdefun,icfun,bcfun,xmesh,tspan,options)Use odeset to create the options structure. Only those options of theunderlying ODE solver shown in the following table are available for pdepe.The defaults obtained by leaving off the input argument options are generallysatisfactory. “Improving ODE Solver Performance” on page 15-16 tells you howto create the structure and describes the properties.PDE Property CategoriesProperties CategoryProperty NameError controlRelTol, AbsTol, NormControlStep-sizeInitialStep, MaxStepExample: Electrodynamics ProblemThis example illustrates the solution of a system of partial differentialequations.
The problem is taken from electrodynamics. It has boundary layersat both ends of the interval, and the solution changes rapidly for small t .The PDEs are2∂ u1∂u 1--------- = 0.024 ------------ – F(u 1 – u 2)2∂t∂x2∂ u2∂u 2--------- = 0.170 ------------ + F(u 1 – u 2)2∂t∂xwhere F ( y ) = exp ( 5.73y ) – exp ( – 11.46y ) . The equations hold on an interval0 ≤ x ≤ 1 for times t ≥ 0 .The solution u satisfies the initial conditionsu 1 ( x, 0 ) ≡ 1u 2 ( x, 0 ) ≡ 0and boundary conditions15-94Partial Differential Equations∂u 1--------- ( 0, t ) ≡ 0∂xu 2(0, t) ≡ 0u 1(1, t) ≡ 1∂u 2--------- ( 1, t ) ≡ 0∂xNote The demo pdex4 contains the complete code for this example.
The demouses subfunctions to place all required functions in a single M-file. To run thisexample type pdex4 at the command line. See “PDE Solver Basic Syntax” onpage 15-85 and “Representing PDE Problems” on page 15-88 for moreinformation.1 Rewrite the PDE. In the form expected by pdepe, the equations are–F ( u1 – u2 )∂ 0.024 ( ∂u 1 ⁄ ∂x )∂ u+.∗ ----- 1 = -----∂x∂t0.170 ( ∂u 2 ⁄ ∂x )F ( u1 – u2 )u211The boundary conditions on the partial derivatives of u have to be written interms of the flux.
In the form expected by pdepe, the left boundary condition is0.024 ( ∂u 1 ⁄ ∂x )010=+.∗u20.170 ( ∂u 2 ⁄ ∂x )00and the right boundary condition isu1 – 10+01.∗0.024 ( ∂u 1 ⁄ ∂x )0.170 ( ∂u 2 ⁄ ∂x )=002 Code the PDE in MATLAB. 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 theform[c,f,s] = pdefun(x,t,u,dudx)15-9515Differential Equationswhere c, f, and s correspond to the c , f , and s terms in Equation 15-3.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 theformu = icfun(x)The code below represents the initial conditions in the MATLAB functionpdex4ic.function u0 = pdex4ic(x);u0 = [1; 0];4 Code the Boundary Conditions Function.
The boundary conditions functions mustbe 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 15-5)of the boundary conditions in the MATLAB 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 .
Theprogram selects the step size in time to resolve this sharp change, but to seethis behavior in the plots, output times must be selected accordingly. There areboundary layers in the solution at both ends of [0,1], so mesh points must beplaced there to resolve these sharp changes. Often some experimentation isneeded 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];15-96Partial Differential Equations6 Apply the PDE Solver. The example calls pdepe with m = 0, the functionspdex4pde, pdex4ic, and pdex4bc, and the mesh defined by x and t at whichpdepe is to evaluate the solution. The pdepe function returns the numericalsolution in a three-dimensional array sol, where sol(i,j,k) approximates thekth component of the solution, u 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 solutioncomponents.u1 = sol(:,:,1);u2 = sol(:,:,2);figuresurf(x,t,u1)title('u1(x,t)')xlabel('Distance x')ylabel('Time t')u1(x,t)10.80.60.40.20211.50.810.60.40.50.2Time t00Distance x15-9715Differential Equationsfiguresurf(x,t,u2)title('u2(x,t)')xlabel('Distance x')ylabel('Time t')u2(x,t)0.80.70.60.50.40.30.20.10211.50.810.60.40.50.2Time t15-9800Distance xSelected BibliographySelected Bibliography[1] Ascher, U., R. Mattheij, and R. Russell, Numerical Solution of BoundaryValue Problems for Ordinary Differential Equations, SIAM, Philadelphia, PA,1995, p.
372.[2] Cebeci, T. and H. B. Keller, “Shooting and Parallel Shooting Methods forSolving the Falkner-Skan Boundary-layer Equation,” J. Comp. Phys., Vol. 7,1971, pp. 289-300.[3] Hairer, E., and G. Wanner, Solving Ordinary Differential Equations II, Stiffand Differential-Algebraic Problems, Springer-Verlag, Berlin, 1991, pp. 5-8.[4] Hindmarsh, A. C., “LSODE and LSODI, Two New Initial Value OrdinaryDifferential Equation Solvers,” SIGNUM Newsletter, Vol. 15, 1980, pp. 10-11.[5] Hindmarsh, A. C., and G. D. Byrne, “Applications of EPISODE: AnExperimental Package for the Integration of Ordinary Differential Equations,”Numerical Methods for Differential Systems, L.
Lapidus and W. E. Schiessereds., Academic Press, Orlando, FL, 1976, pp 147-166.[6] Shampine, L. F., Numerical Solution of Ordinary Differential Equations,Chapman & Hall Mathematics, 1994.[7] Shampine, L. F., and M. K. Gordon, Computer Solution of OrdinaryDifferential Equations, W.H. Freeman & Co., 1975.[8] Skeel, R. D. and M. Berzins, “A Method for the Spatial Discretization ofParabolic Equations in One Space Variable,” SIAM Journal on Scientific andStatistical Computing, Vol.
11, 1990, pp.1-32.15-9915Differential Equations15-10016Sparse MatricesFunctionSummaryIntroduction. . . . . . . . . . . . . . .Sparse Matrix Storage . . . . . . . . . . .Introduction. . . . . . . . . . . . . .Creating Sparse. MatricesSparseMatrixStorage. from. . .Outside. . . MATLAB. . . .ImportingSparseMatricesGeneral Storage Information . . . . . . . . .CreatingViewing SparseSparseMatricesMatrices. . .
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
