《Matlab 第七讲：矩阵代数与数值技术.ppt》由会员分享，可在线阅读，更多相关《Matlab 第七讲：矩阵代数与数值技术.ppt（56页珍藏版）》请在金锄头文库上搜索。
1、1,矩阵代数与数值技术,数学软件 Matlab,2,矩阵代数,3,矩阵运算的一些函数,4,求代数方程的近似根(解),5,解方程（代数方程）是最常见的数学问题之一，也是众多应用领域中不可避免的问题之一,目前还没有一般的解析方法来求解非线性方程，但如果在任意给定的精度下，能够解出方程的近似解，则可以认为求解问题已基本解决，至少可以满足实际需要,问题背景和实验目的,代数方程近似求解,6,相关概念,如果 f(x) 是一次多项式，称上面的方程为线性方程；否则称之为非线性方程,线性方程 与 非线性方程,本实验主要讨论非线性方程的数值求解,7,Matlab 解方程函数,roots(p)：多项式的所有零点，p
2、 是多项式系数向量。,fzero(f,x0)：求 f=0 在 x0 附近的根，f 可以使用 inline、字符串、或 ，但不能是方程或符号表达式！,solve(f,v)：求方程关于指定自变量的解，f 可以是用字符串表示的方程、符号表达式或符号方程； solve 也可解方程组(包含非线性)； 得不到解析解时，给出数值解。,linsolve(A,b)：解线性方程组。,8,Matlab 多项式的表示,在 Matlab 中，n 次多项式是用一个长度为 n+1的向量来表示，缺少的幂次项系数为 0,注：系数中的零不能省！,将多项式转化成符号表达式：poly2sym, poly2sym(2,-1,0,3),
3、Matlab 中多项式的表示方法,9,多项式的零点,x=roots(p)：若 p 是 n 次多项式，则输出是 p=0 的 n 个根组成的 n 维向量,若已知多项式的全部零点，则可用 poly 函数给出该多项式,p=poly(x),多项式的零点,10,线性方程组求解,线性方程组求解,linsolve(A,b)：解线性方程组 Ax = b,例：解方程组, A=1 2 -1; 1 0 1; 1 3 0; b=2;3;8; x=linsolve(A,b) x=Ab,b 是列向量！,11,非线性方程的根,非线性方程的数值求解,fzero(f,x0)：求方程 f=0 在 x0 附近的根,方程可能有多个根，
4、但 fzero 只给出 x0 附近的一个,fzero 先找出一个包含 x0 的区间，使得 f 在这个区间两个端点上的函数值异号，然后再在这个区间内寻找方程 f=0 的根；如果找不到这样的区间，则返回 NaN,x0 是一个标量，为参考点，不能缺省,由于 fzero 是根据函数是否穿越横轴来决定零点，因此它无法确定函数曲线仅触及横轴但不穿越的零点，如 |sin(x)| 的所有零点,12,非线性方程的根,fzero 的另外一种调用方式,fzero(f,a,b),方程在 a,b 内可能有多个根，但 fzero 只能求一个,求方程 f=0 在 a,b 区间内的根。,参数 f 可通过以下方式给出：,字符串
5、： 内联函数：,f 不是方程，即不含等号！,f=inline(x3-3*x+1); fzero(f,2),fzero(x3-3*x+1,2),13, fzero(sin(x),10) fzero(sin,10), fzero(x3-3*x+1,1) fzero(x3-3*x+1,1,2), fzero(x3-3*x+1=0,1),X, fzero(x3-3*x+1,-2,0) f=inline(x3-3*x+1); fzero(f,-2,0),用 fzero 求零点时可以先通过作图确定零点的大致范围,例：,fzero 举例,14,定积分的近似计算,15,定积分计算的基本公式是牛顿莱布尼兹公式。
6、但当被积函数的原函数不知道时，如何计算？这时就需要利用近似计算。特别是在许多实际应用中，被积函数甚至没有解析表达式，而是一条实验记录曲线，或一组离散的采样值，此时只能用近似方法计算定积分。,问题背景和实验目的,定积分的近似计算,Matlab 求积分函数,数值积分函数：trapz、quad、dblquad 符号积分函数：int,16,定积分几何意义,17,trapz(x,y) x 为分割点（节点）组成的向量， y 为被积函数在节点上的函数值组成的向量。,trapz,trapz,18,前面的做法,例：用梯形法计算下面定积分 ( 取 n=100 ),解：,a=0, b=1, n=100, yi =
7、f (xi) = 1/( 1+xi2 ), x=0:1/100:1; y=1./(1+x.2); trapz(x, y),trapz函数,trapz(x,1./(1+x.2),trapz 举例,19,quad(f,a,b,tol) f = f(x) 为被积函数，a,b 为积分区间，tol 为计算精度,将自变量看成是向量,不用自己分割积分区间 可以指定计算精度，若不指定，缺省精度是 10-6 精度越高，函数运行的时间越长 此处的函数 f 是数值形式，应该使用数组运算，即： .* ./ . .,quad,quad,20,解：, quad(1./(1+x.2),0,1), quad(1./(1+x.
8、2),0,1,1e-10), quad(1./(1+x.2),0,1,1e-16),函数表达式一定要用 单引号 括起来！ 涉及的运算一定要用 数组运算！,例：用 quad 计算定积分：,quad 举例,21,dblquad(f,a,b,c,d,tol),tol 为计算精度，若不指定，则缺省精度为 10-6 f 可以是： 字符串；inline 定义的内联函数；函数句柄 a,b 是 第一积分变量 的积分区间， c,d 是 第二积分变量 的积分区间,按字母顺序，大写字母排在小写字母的前面,dblquad,抛物线法计算二重积分： dblquad,22, f=inline(4*x*y+3*y2); I=
9、dblquad(f,-1,1,0,2),f 中关于第一自变量的运算是数组运算，即把 x 看成是向量，y 看成是标量。也可以全部采用数组运算,例：计算二重积分, dblquad(inline(4*x*y+3*x2),-1,1,0,2), dblquad(inline(4*x*y+3*x.2),-1,1,0,2),X,例：计算二重积分,dblquad 举例,23,例：计算二重积分, dblquad(x,y)4*x*y+3*x.2 , -1, 1, 0, 2),指定 x、y 分别是第一和第二积分变量, dblquad(inline(4*x*y+3*x.2) ,-1, 1, 0, 2),被积函数 f
10、(x,y) 的另一种定义方法：匿名函数,dblquad 举例,24,int(f,a,b) 计算 f 关于默认自变量 的定积分，积分区间为a,b。,int(f) 计算 f 关于默认自变量 的不定积分。,int(f,v,a,b) 计算函数 f 关于自变量 v 的定积分，积分区间为 a, b,int(f,v) 计算函数 f 关于自变量 v 的不定积分,findsym(f,1),int,符号积分： int,25,例：用 int 函数计算定积分：,解：, syms x; f=1/(1+x2); int(f,x,0,1), f=sym(1/(1+x2); int(f,x,0,1), int(1/(1+x2
11、),x,0,1),或, int(1/(1+x2),0,1),或,或,int 举例,26,double(a) 将 a 转化为双精度型，若 a 是字符，则取对应的 ASCII 码, a=3; double(a) double(a),例：,ans = 3,ans = 97,相关函数,27, x=1:0.001:2; y=exp(x.(-2); trapz(x,y),梯形法：,抛物线法：, quad(exp(x.(-2),1,2,10e-10),符号积分法：, syms x int(exp(x(-2),x,1,2),例：用 Matlab 函数近似计算定积分,数值实验,28,抛物线法：, dblquad
12、(inline(x+y2),0,2,-1,1),符号积分法：, f=int(x+y2,y,-1,1); int(f,x,0,2),数值实验,例：用 Matlab 函数近似计算二重积分,29,求微分方程的解,30,自牛顿发明微积分以来，微分方程在描述事物运动规律上已发挥了重要的作用。实际应用问题通过数学建模所得到的方程，绝大多数是微分方程。,由于实际应用的需要，人们必须求解微分方程。然而能够求得解析解的微分方程十分有限，绝大多数微分方程需要利用数值方法来近似求解。,本实验主要研究如何用 Matlab 来计算微分方程（组）的数值解。,问题背景和实验目的,31,Matlab 解初值问题,用 Malt
13、ab自带函数 解初值问题,求解析解：dsolve,求数值解： ode45、ode23、 ode113、ode23t、ode15s、 ode23s、ode23tb,32,数值求解,T,Y = solver(odefun,tspan,y0),其中 y0 为初值条件，tspan为求解区间；Matlab在数值求解时自动对求解区间进行分割，T (列向量) 中返回的是分割点的值(自变量)，Y (数组) 中返回的是这些分割点上的近似解，其列数等于因变量的个数。 solver 为Matlab的ODE求解器（可以是 ode45、ode23、ode113、ode15s、ode23s、ode23t、ode23tb）
14、,没有一种算法可以有效地解决所有的 ODE 问题，因此MATLAB 提供了多种ODE求解器，对于不同的ODE，可以调用不同的求解器。,33,Matlab的ODE求解器,34,Runge-Kutta 方法,龙格-库塔方法,是一类求解常微分方程的数值方法,有多种不同的迭代格式,35,Runge-Kutta 方法,用得较多的是 四阶R-K方法,其中,36,四阶 R-K 方法源程序,clear; f=sym(y+2*x/y2); a=0; b=2; h=0.4; n=(b-a)/h+1; % n=(b-a)/h; x=0; y=1; szj=x,y; for i=1:n-1 % i=1:n L1=su
15、bs(f,x,y,x,y); L2=subs(f,x,y,x+h/2,y+L1*h/2); L3=subs(f,x,y,x+h/2,y+L2*h/2); L4=subs(f,x,y,x+h,y+L3*h); y=y+h*(L1+2*L2+2*L3+L4)/6; x=x+h; szj=szj;x,y; end plot(szj(:,1),szj(:,2), dg-),37,Runge-Kutta 方法,例：用 Euler 法解初值问题,解析解：,y=1/3*(-18-54*x+45*exp(3*x)(1/3),38,参数说明,odefun 为显式常微分方程，可以用命令 inline 定义，或在函数文件中定义，然后通过函数句柄调用。,fun=inline(-2*y+2*x2+2*x,x,y); x,y=ode23(fun,0,0.5,1);,注：也可以在 tspan 中指定对求解区间的分割，如：,x,y=ode23(fun,0:0.1:0.5,1); %此时 x=0:0.1:0.5,T,Y = solver(odefun,tspan,y0),39,数值求解举例,如果需求解的问题是高阶常微分方程，则需将其化为一阶常微分方程组，此时必须用函数文件来定义该常微分方程组。,令,40,数值求解举例,先编写函数
