当前位置：文档之家 › 数值积分与微分MATLAB公式

数值积分与微分MATLAB公式

  • 格式：doc
  • 大小：50.00 KB
  • 文档页数：3

数值积分与微分

实验目的：1）用matlab软件掌握梯形公式、辛普森公式和蒙特卡罗方法计算数值积分；

2）通过实例学习用数值积分和数值微分解决实际问题。

实验内容：

第一题：用梯形、辛普森和蒙特卡罗方法计算积分。改变步长（对梯形），改变精度要求（对辛普森），改变随机点数目（对蒙特卡罗），进行比较、分析。

y=21e22x ,-2x2

解：用三种方法计算积分的源程序如下：

h=4/50;

x=-2:h:2;

y=fun(x);

z1=trapz(y)*h

z2=quad('fun',-2,2)

plot(x,y,'g');

n=1000;

x1=rand(1,n);

y1=fun(x1.*2);

z3=sum(y1)*4/n

对梯形公式取h=4/50,4/100,4/10000;对辛普森分别取精度为103，710,108；对蒙特卡罗方法分别取n=1000,10000,100000.得到的结果如下：

梯形公式 辛普森公式 蒙特卡罗方法

0.95438456767789 0.95449943824154 0.93999211059586

0.95447094168964 0.95449973610735 0.95794338594866

0.95449973322412 0.95449973610373 0.95427317381756

从得到的结果可以看到对梯形公式，步长越小，计算的积分结果越准确；对于辛普森公式，在一般的103精度下结果已经很准确（小数点后前六位均为准确数字），提高精度后结果更加精确，可见辛普森具有很高的优越性，但它的局限性在于必须要有函数解析式；对于蒙特卡罗方法，虽然结果具有随机性，但随着n增大，得到的结果越来越接近准确值。

第二题：测得活塞中气体压力P和体积V的一组数据如下：

P 60 80 100 120 140 160 180

V 80.0 69.2 60.0 52.0 45.0 38.6 32.5

求V＝60，50处，V改变1时P的改变量；求从70减至40时气体作的功。

解：用中点公式计算导数k.则P＝kV。因为V＝1，所以P数值上等于k。

取h=0.1，利用三次样条计算P在V－h，V+h处的数值，从而利用中点公式计算导数。

源程序如下：

h=0.1;x=[60-h,60+h];

p=60:20:180;

v=[80.0 69.2 60.0 52.0 45.0 38.6 32.5];

y=spline(v,p,x);

p=(y(1)-y(2))/2/h

结果为 p=2.3341（2/inlbf） 同理可以算出V=50时,p=2.7891（2/inlbf）

求导的问题也可以用书后补充知识中样条求导的方法解决，计算后可以得到相同结果。

利用三次样条插值计算V在40～70之间时相应的一系列P值，然后用梯形公式计算积分即得气体作功。

源程序如下：

x=40:30/1000:70;

p=60:20:180;

v=[80.0 69.2 60.0 52.0 45.0 38.6 32.5];

y=spline(v,p,x);

w=trapz(x,y)

输出结果为：W= 3414.36(lbf.in)

第三题：冰淇淋的下部为锥体，上部为半球。设它由锥面z=22yx和球面1)1(222zyx围成，用蒙特卡罗方法计算它的体积。

解：两个曲面方程联立可以解得几何体的边界方程为单位圆：22yx＝1。

应用蒙特卡罗均值估计法计算体积的思路如下：

利用计算机每次产生两个0～1的随机数x,y，若落在单位圆内，则计算球面与锥面上在（x,y）处的z值之差，产生n次随机数，并将得到的z值累加，累加和除以n即得所求体积的四分之一。取n＝100000，matlab源程序如下：

z=0;

n=100000;

for i=1:n

x=rand(1,2);

y=0;

if x(1)^2+x(2)^2<=1

y=sqrt(1-x(1)^2-x(2)^2)+1-sqrt(x(1)^2+x(2)^2);

z=z+y;

end

end

V=4*z/n

输出结果为：V= 3.1336

即所求冰淇淋的体积为3.1336。

小结：通过本章的学习掌握了用matlab软件利用梯形公式、辛普森公式和蒙特卡罗方法计算数值积分的方法。虽然利用数值积分与微分得到的不是我们一贯接触的那种百分之百准确的“完美”的结果，但是现实当中许多问题是不能用经典的解析方法来描述的，而

只能通过数值积分和数值微分的方法得到满足我们精度要求的结果，所以学好本章的内容对我们解决许多实际问题有很大的帮助。另外在本章当中还有许多精巧的思维方法，比如蒙特卡罗方法，他的这种思想方法很朴素，但是由此推广得到的各种方法（比如均值估计法）却可以解决很多实际问题。我想我们不仅应该学好解决问题的具体方法，更应该体会、理解方法当中蕴涵的深刻的数学思想和精巧的思维方式。

合集下载

matlab公式总结

matlab公式总结

一 matlab常用函数

1、 特殊变量与常数

ans 计算结果的变量名 computer 确定运行的计算机 eps 浮点相对精度 Inf 无穷大 I 虚数单位 inputname 输入参数名

NaN 非数 nargin 输入参数个数 nargout 输出参数的数目 pi

圆周率 nargoutchk 有效的输出参数数目 realmax 最大正浮点数

realmin 最小正浮点数 varargin 实际输入 的参量 varargout 实际返回的参量 操作符与特殊字符 + 加 - 减 * 矩阵乘法 .* 数组乘（对应元素相乘） ^ 矩阵幂 .^ 数组幂（各个元素求幂） \ 左除或反斜杠 / 右除或斜面杠 ./ 数组除（对应元素除） kron

Kronecker张量积 : 冒号 () 圆括 [] 方括 . 小数点 .. 父目录 ... 继续 , 逗号（分割多条命令） ; 分号（禁止结果显示） % 注释 !

感叹号 ' 转置或引用 = 赋值 == 相等 <> 不等于 & 逻辑与

| 逻辑或 ~ 逻辑非 xor 逻辑异或

2、 基本数学函数

abs 绝对值和复数模长 acos,acodh 反余弦，反双曲余弦 acot,acoth 反余切，反双曲余切 acsc,acsch 反余割，反双曲余割 angle 相角

matlab ode45 求解带积分的常微分方程

matlab ode45 求解带积分的常微分方程

概述

在工程和科学领域中，常微分方程是一种常见的数学建模工具。其中，带积分的常微分方程更是一种需要特殊解法的方程形式。MATLAB是一种功能强大的数学工具软件，而ode45是MATLAB中用于求解常微分方程的函数之一。本文将详细介绍如何使用MATLAB中的ode45函数来求解带积分的常微分方程。

一、带积分的常微分方程简介

带积分的常微分方程是指在微分方程中出现积分形式的项，通常表现为对某个函数进行积分。这种形式的微分方程在工程和科学领域中有着广泛的应用，例如在电路分析、控制系统、生物学模型等领域中都能见到。典型的带积分的常微分方程形式如下所示：

y' = f(t,y) + ∫g(t,y)dt

其中，y'表示y对自变量t的导数，f(t,y)为已知的函数，g(t,y)为未知的函数需要求解。这种形式的微分方程要比普通的常微分方程更复杂，需要使用特定的求解方法来得到解析解或数值解。

二、MATLAB中的ode45函数介绍

MATLAB是一种被广泛应用于科学计算和工程领域的数学软件工具，其中有丰富的数值计算函数库。其中，用于求解常微分方程的ode45函数是应用较为广泛的函数之一。ode45函数可以通过数值计算的方法来求解常微分方程的数值解，其基本调用格式如下：

[t,y] = ode45(odefun,tspan,y0)

其中，odefun是定义了微分方程的函数句柄，tspan是求解的时间范围，y0是初始条件。ode45函数会返回微分方程在tspan范围内的数值解t和对应的y值。

三、使用MATLAB求解带积分的常微分方程

对于带积分的常微分方程，我们需要将其转化为标准形式，然后利用MATLAB的ode45函数进行求解。假设我们有如下形式的带积分的常微分方程：

y' = f(t,y) + ∫g(t,y)dt

我们将其转化为等价的无積分項的方程形式，例如

∂F/∂t = f(t,y) + ∫g(t,y)dt

matlab常微分方程的数值解法实验报告

matlab常微分方程的数值解法实验报告

实验四

常微分方程的数值解法

指令：

[t,y]=ode23(‘fun’,tspan,yo) 2/3阶龙格库塔方法

[t,y]=ode45(‘fun’,tspan,yo) 4/5阶龙格库塔方法

[t,y]=ode113(‘fun’,tspan,yo) 高阶微分方程数值方法

其中fun是定义函数的文件名。该函数fun必须以为dx输出量，以t,y为输入量。tspan=[t0

tfina]表示积分的起始值和终止值。yo是初始状态列向量。

考虑到初始条件有

00d, (0)0,dd, (0)0.dSSISStISIIIIt (5.24)

这就是Kermack与McKendrick的SIR仓室模型. 方程(5.24)无法求出()St和()It的解析解.我们先做数值计算。Matlab代码为：

function dy=rigid(t,y)

dy=zeros(2,1);

a=1;

b=0.3;

dy(1)=a*y(1).*y(2)-b*y(1);

dy(2)=-a*y(1).*y(2);

ts=0:.5:50;

x0=[0.02,0.98];

[T,Y]=ode45('rigid',ts,x0);

%plot(T,Y(:,1),'-',T,Y(:,2),'*')

plot(Y(:,2),Y(:,1),'b--')

xlabel('s')

ylabel('i')

任务：

1 画出i（t），

2分析各参数的影响

例57：求解两点边值问题：0)5(,0)1(,32yyxyyx 。（注意：相应的数值解法比较复杂）。

y=dsolve('x*D2y-3*Dy=x^2','y(1)=0,y(5)=0','x') ↙

使用Matlab进行公式推导、求解

使用Matlab进行公式推导、求解

使用Matlab进行公式推导、求解

MATLAB是一种强大的数值计算和科学编程软件，在MATLAB中，可以使用其丰富的数学函数和符号计算工具进行公式推导和求解。本文将以案例的形式介绍如何使用MATLAB进行公式推导和求解，包括符号计算、方程求解、微分和积分等方面的应用。

案例1：对以下公式进行“去括号展开”和“幂级数形式整理”

（1）采用符号工具箱syms进行方程写入

注意：代码建议写在实时编辑器中，这里语句都没有加分号，运行之后会直接显示结果如上图

（2）“去括号展开”采用expand，“幂级数形式整理”采用series函数

说明：

"expand" 是 MATLAB 中用于展开代数表达式的一个函数。它是 Symbolic Math Toolbox 的一部分。当应用于符号表达式时，"expand" 函数会通过展开括号和简化项来展开表达式。这在处理复杂的数学表达式或尝试简化和操作方程时非常有用。

"series" 是 MATLAB 中用于展开符号表达式的级数形式的一个函数。它也属于 Symbolic Math Toolbox，提供了处理符号表达式和符号数学的功能，可以将符号表达式展开为指定的级数形式（‘order’后加展开到的最高指数），通常是泰勒级数。它可以帮助在数学建模、分析和求解问题时处理复杂的表达式，并在需要时进行近似计算。

（3）案例完整代码

clearsymsa0a1a2a3b0b1b2b3tt_ea=[a2;a1;a0];b=[b2;b1;b0];x(t)=[t^2t1]*a+a3*t*(t-t_e)^2y(t)=[t^2t1]*b+b3*t*(t-t_e)^2x(t)=expand(x)y(t)=expand(y)x(t)=series(x,t,'Order',4)y(t)=series(y,t,'Order',4)

案例2：建立等式方程组，转换为线性矩阵表达式，然后求解

（1）利用实时编辑器建立方程组