例6.1  求矩阵A的每行及每列的最大和最小元素，并求整个矩阵的最大和最小元素。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image002.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

A=[13,-56,78;25,63,-235;78,25,563;1,0,-1];

max(A,[],2)         %求每行最大元素

min(A,[],2)         %求每行最小元素

max(A)              %求每列最大元素

min(A)               %求每列最小元素

max(max(A))         %求整个矩阵的最大元素。也可使用命令：max(A(:))

min(min(A))         %求整个矩阵的最小元素。也可使用命令：min(A(:))

例6.2  求矩阵A的每行元素的乘积和全部元素的乘积。

A=[1,2,3,4;5,6,7,8;9,10,11,12];

S=prod(A,2)

prod(S)        %求A的全部元素的乘积。也可以使用命令prod(A(:))

例6.3  求向量X=(1！，2！，3！，…，10！)。

X=cumprod(1:10)

例6.4  对二维矩阵x，从不同维方向求出其标准方差。

x=[4,5,6;1,4,8]                % 产生一个二维矩阵x

y1=std(x,0,1)

y2=std(x,1,1)

y3=std(x,0,2)

y4=std(x,1,2)

例6.5  生成满足正态分布的10000×5随机矩阵，然后求各列元素的均值和标准方差，再求这5列随机数据的相关系数矩阵。

X=randn(10000,5);

M=mean(X)

D=std(X)

R=corrcoef(X)

例6.6  对下列矩阵做各种排序。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image004.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

A=[1,-8,5;4,12,6;13,7,-13];

sort(A)              %对A的每列按升序排序

-sort(-A,2)            %对A的每行按降序排序

 [X,I]=sort(A)        %对A按列排序，并将每个元素所在行号送矩阵I

例6.7  给出概率积分

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image006.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

的数据表如表6.1所示，用不同的插值方法计算f(0.472)。

表6.1  概率积分数据表

x=0.46:0.01:0.49;         %给出x，f(x)

f=[0.4846555,0.4937542,0.5027498,0.5116683];

format long

interp1(x,f,0.472)     %用默认方法，即线性插值方法计算f(x)

interp1(x,f,0.472,'nearest')  %用最近点插值方法计算f(x)

interp1(x,f,0.472,'spline')    %用3次样条插值方法计算f(x)

interp1(x,f,0.472,'cubic')      %用3次多项式插值方法计算f(x)

format short

例6.8  某检测参数f随时间t的采样结果如表5.1，用数据插值法计算t=2，7，12，17，22，17，32，37，42，47，52，57时的f值。

表5.1 检测参数f随时间t的采样结果

T=0:5:65;

X=2:5:57;

F=[3.2015,2.2560,879.5,1835.9,2968.8,4136.2,5237.9,6152.7,...

6725.3,6848.3,6403.5,6824.7,7328.5,7857.6];

F1=interp1(T,F,X)        %用线性插值方法插值

F1=interp1(T,F,X,'nearest')     %用最近点插值方法插值

F1=interp1(T,F,X,'spline')      %用3次样条插值方法插值

F1=interp1(T,F,X,'cubic')       %用3次多项式插值方法插值

例6.9  设z=x²+y²，对z函数在[0,1]×[0,2]区域内进行插值。

x=0:0.1:1;y=0:0.2:2;

[X,Y]=meshgrid(x,y);       %产生自变量网格坐标

Z=X.^2+Y.^2;                %求对应的函数值

interp2(x,y,Z,0.5,0.5)     %在(0.5,0.5)点插值

interp2(x,y,Z,[0.5 0.6],0.4) %在(0.5,0.4)点和(0.6，0.4)点插值

interp2(x,y,Z,[0.5 0.6],[0.4 0.5])%在(0.5,0.4)点和(0.6，0.5)点插值

%下一命令在(0.5,0.4),(0.6,0.4),(0.5,0.5)和(0.6,0.5)各点插值

interp2(x,y,Z,[0.5 0.6]',[0.4 0.5])

例6.10  某实验对一根长10米的钢轨进行热源的温度传播测试。用x表示测量点(米)，用h表示测量时间(秒)，用T表示测得各点的温度(℃)，测量结果如表6.2所示。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image007.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />表6.2  钢轨各点温度测量值

试用用3次多项式插值求出在一分钟内每隔10秒、钢轨每隔0.5米处的温度。

x=0:2.5:10;

h=[0:30:60]';

T=[95,14,0,0,0;88,48,32,12,6;67,64,54,48,41];

xi=[0:0.5:10];

hi=[0:10:60]';

temps=interp2(x,h,T,xi,hi,'cubic');

mesh(xi,hi,temps);

例6.11  用一个3次多项式在区间[0，2π]内逼近函数file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image009.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />。

X=linspace(0,2*pi,50);

Y=sin(X);

P=polyfit(X,Y,3)             %得到3次多项式的系数和误差

X=linspace(0,2*pi,20);

Y=sin(X);

Y1=polyval(P,X)

plot(X,Y,':o',X,Y1,'-*')

例6.12  设

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image011.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

(1)求f(x)+g(x)、f(x)-g(x)。

(2)求f(x)×g(x)、f(x)/g(x)。

f=[3,-5,2,-7,5,6];g=[3,5,-3];g1=[0,0,0,g];

f+g1                %求f(x)+g(x)

f-g1                 %求f(x)-g(x)

conv(f,g)             %求f(x)*g(x)

 [Q,r]=deconv(f,g)     %求f(x)/g(x)，商式送Q，余式送r。

例6.13  求有理分式的导数。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image013.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

P=[3,5,0,-8,1,-5];

Q=[10,5,0,0,6,0,0,7,-1,0,-100];

[p,q]=polyder(P,Q)

例6.14  已知多项式x⁴+8x³-10，分别取x=1.2和一个2×3矩阵为自变量计算该多项式的值。

A=[1,8,0,0,-10];             % 4次多项式系数

x=1.2;                     % 取自变量为一数值

y1=polyval(A,x)

x=[-1,1.2,-1.4;2,-1.8,1.6]      % 给出一个矩阵x

y2=polyval(A,x)     % 分别计算矩阵x中各元素为自变量的多项式之值

例6.15  仍以多项式x⁴+8x³-10为例，取一个2×2矩阵为自变量分别用polyval和polyvalm计算该多项式的值。

A=[1,8,0,0,-10];             % 多项式系数

x=[-1,1.2;2,-1.8]             % 给出一个矩阵x

y1=polyval(A,x)             % 计算代数多项式的值

y2=polyvalm(A,x)             % 计算矩阵多项式的值

例6.16  求多项式x⁴+8x³-10的根。

A=[1,8,0,0,-10];

x=roots(A)

例6.17  已知

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image015.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

(1) 计算f(x)=0 的全部根。

(2) 由方程f(x)=0的根构造一个多项式g(x)，并与f(x)进行对比。

P=[3,0,4,-5,-7.2,5];

X=roots(P)            %求方程f(x)=0的根

G=poly(X)            %求多项式g(x)

例6.18  设x由[0,2π]间均匀分布的10个点组成，求file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image009.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />的1~3阶差分。

X=linspace(0,2*pi,10);

Y=sin(X);

DY=diff(Y);       %计算Y的一阶差分

D2Y=diff(Y,2);     %计算Y的二阶差分，也可用命令diff(DY)计算

D3Y=diff(Y,3);     %计算Y的三阶差分，也可用diff(D2Y)或diff(DY,2)

例6.19  设

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image018.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

用不同的方法求函数f(x)的数值导数，并在同一个坐标系中做出f'(x)的图像。

f=inline('sqrt(x.^3+2*x.^2-x+12)+(x+5).^(1/6)+5*x+2');

g=inline('(3*x.^2+4*x-1)./sqrt(x.^3+2*x.^2-x+12)/2+1/6./(x+5).^(5/6)+5');

x=-3:0.01:3;

p=polyfit(x,f(x),5);          %用5次多项式p拟合f(x)

dp=polyder(p);                 %对拟合多项式p求导数dp

dpx=polyval(dp,x);            %求dp在假设点的函数值

dx=diff(f([x,3.01]))/0.01;   %直接对f(x)求数值导数

gx=g(x);                         %求函数f的导函数g在假设点的导数

plot(x,dpx,x,dx,'.',x,gx,'-');   %作图

例6.20  用两种不同的方法求：

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image020.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

先建立一个函数文件ex.m：

function ex=ex(x)

ex=exp(-x.^2);

然后在MATLAB命令窗口，输入命令：

format long

I=quad('ex',0,1)     %注意函数名应加字符引号

I=quadl('ex',0,1)

例6.21  用trapz函数计算：

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image020.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

X=0:0.01:1;

Y=exp(-X.^2);

trapz(X,Y)

例6.22  计算二重定积分

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image022.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

(1) 建立一个函数文件fxy.m：

function f=fxy(x,y)

global ki;

ki=ki+1;              %ki用于统计被积函数的调用次数

f=exp(-x.^2/2).*sin(x.^2+y);

(2) 调用dblquad函数求解。

global ki;ki=0;

I=dblquad('fxy',-2,2,-1,1)

ki

例6.23  给定数学函数

x(t)=12sin(2π×10t+π/4)+5cos(2π×40t)

取N=128，试对t从0~1秒采样，用FFT作快速傅立叶变换，绘制相应的振幅-频率图。

N=128;                            % 采样点数

T=1;                              % 采样时间终点

t=linspace(0,T,N);             % 给出N个采样时间ti(i=1:N)

x=12*sin(2*pi*10*t+pi/4)+5*cos(2*pi*40*t);  % 求各采样点样本值x

dt=t(2)-t(1);                   % 采样周期

f=1/dt;                          % 采样频率(Hz)

X=fft(x);                        % 计算x的快速傅立叶变换X

F=X(1:N/2+1);                   % F(k)=X(k)(k=1:N/2+1)

f=f*(0:N/2)/N;                  % 使频率轴f从零开始

plot(f,abs(F),'-*')            % 绘制振幅-频率图

xlabel('Frequency');

ylabel('|F(k)|')

例6.24  用直接解法求解下列线性方程组。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image024.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4];

b=[13,-9,6,0]';

x=A\b

例6.25  用LU分解求解例6.24中的线性方程组。

A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4];

b=[13,-9,6,0]';

[L,U]=lu(A);

x=U\(L\b)

[L,U ,P]=lu(A);

x=U\(L\P*b)

例6.26  用QR分解求解例6.24中的线性方程组。

A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4];

b=[13,-9,6,0]';

[Q,R]=qr(A);

x=R\(Q\b)

[Q,R,E]=qr(A);

x=E*(R\(Q\b))

例6.27  用Cholesky分解求解例6.24中的线性方程组。

A=[2,1,-5,1;1,-5,0,7;0,2,1,-1;1,6,-1,-4];

b=[13,-9,6,0]';

R=chol(A)

??? Error using ==> chol

Matrix must be positive definite

Jacobi迭代法的MATLAB函数文件Jacobi.m如下：

function [y,n]=jacobi(A,b,x0,eps)

if nargin==3

    eps=1.0e-6;

elseif nargin<3

    error

    return

end     

D=diag(diag(A));    %求A的对角矩阵

L=-tril(A,-1);       %求A的下三角阵

U=-triu(A,1);       %求A的上三角阵

B=D\(L+U);

f=D\b;

y=B*x0+f;

n=1;                  %迭代次数

while norm(y-x0)>=eps

    x0=y;

    y=B*x0+f;

    n=n+1;

end

例6.28  用Jacobi迭代法求解下列线性方程组。设迭代初值为0，迭代精度为10^-6。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image026.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

在命令中调用函数文件Jacobi.m，命令如下：

A=[10,-1,0;-1,10,-2;0,-2,10];

b=[9,7,6]';

[x,n]=jacobi(A,b,[0,0,0]',1.0e-6)

 

Gauss-Serdel迭代法的MATLAB函数文件gauseidel.m如下：

function [y,n]=gauseidel(A,b,x0,eps)

if nargin==3

    eps=1.0e-6;

elseif nargin<3

    error

    return

end     

D=diag(diag(A));    %求A的对角矩阵

L=-tril(A,-1);      %求A的下三角阵

U=-triu(A,1);       %求A的上三角阵

G=(D-L)\U;

f=(D-L)\b;

y=G*x0+f;

n=1;                  %迭代次数

while norm(y-x0)>=eps

    x0=y;

    y=G*x0+f;

    n=n+1;

end

例6.29  用Gauss-Serdel迭代法求解下列线性方程组。设迭代初值为0，迭代精度为10^-6。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image026.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

在命令中调用函数文件gauseidel.m，命令如下：

A=[10,-1,0;-1,10,-2;0,-2,10];

b=[9,7,6]';

[x,n]=gauseidel(A,b,[0,0,0]',1.0e-6)

例6.30  分别用Jacobi迭代和Gauss-Serdel迭代法求解下列线性方程组，看是否收敛。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image028.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

a=[1,2,-2;1,1,1;2,2,1];

b=[9;7;6];

[x,n]=jacobi(a,b,[0;0;0])

[x,n]=gauseidel(a,b,[0;0;0])

 

有了上面这些讨论，下面设计一个求解线性方程组的函数文件line_solution.m。

function [x,y]=line_solution(A,b)

[m,n]=size(A);

y=[];

if norm(b)>0                  %非齐次方程组

if rank(A)==rank([A,b])

if rank(A)==n           %有惟一解

disp('原方程组有惟一解x');

x=A\b;

else                   %方程组有无穷多个解，基础解系

disp('原方程组有无穷个解，特解为x，其齐次方程组的基础解系为y');

x=A\b;

y=null(A,'r');

end

else

disp('方程组无解');      %方程组无解

x=[];

end

else                        %齐次方程组

disp('原方程组有零解x');

x=zeros(n,1);              %0解

if rank(A)<n

disp('方程组有无穷个解，基础解系为y');    %非0解

y=null(A,'r');

end

end

例6.31  求解方程组。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image030.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

A=[1,-2,3,-1;3,-1,5,-3;2,1,2,-2];

b=[1;2;3];

[x,y]=line_solution(A,b)

例6.32  求方程组的通解。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image032.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

format  rat           %指定有理式格式输出

A=[1,1,-3,-1;3,-1,-3,4;1,5,-9,-8];

b=[1,4,0]';

[x,y]=line_solution(A,b);

x,y

format short          %恢复默认的短格式输出

例6.33 求file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image034.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />在x₀=-5和x₀=1作为迭代初值时的零点。

先建立函数文件fz.m：

function f=fz(x)

f=x-1/x+5;

然后调用fzero函数求根。：

fzero('fz',-5)          %以-5作为迭代初值

fzero('fz',1)             %以1作为迭代初值

例6.34  求下列方程组在(1，1，1)附近的解并对结果进行验证。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image036.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

	file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image038.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

 

首先建立函数文件myfun.m。

function F=myfun (X)

x=X(1);

y=X(2);

z=X(3);

F(1)=sin(x)+y+z^2*exp(x);

F(2)=x+y+z;

F(3)=x*y*z;

在给定的初值x₀=1，y₀=1，z₀=1下，调用fsolve函数求方程的根。

X=fsolve('myfun',[1,1,1],optimset('Display', 'off'))

例6.35  求圆和直线的两个交点。

圆：file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image040.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

直线：file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image042.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

先建立方程组函数文件fxyz.m：

function F=fxyz(X)

x=X(1);

y=X(2);

z=X(3);

F(1)=x^2+y^2+z^2-9;

F(2)=3*x+5*y+6*z;

F(3)=x-3*y-6*z-1;

再在MATLAB命令窗口，输入命令：

X1=fsolve('fxyz',[-1,1,-1],optimset('Display', 'off')) %求第一个交点

X2=fsolve('fxyz',[1,-1,1],optimset('Display', 'off')) %求第二个交点

例6.36  求函数

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image044.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

在区间(-10，1)和(1，10)上的最小值点。

首先建立函数文件fx.m：

function f=f(x)

f=x-1/x+5;

上述函数文件也可用一个语句函数代替：

f=inline('x-1/x+5')

再在MATLAB命令窗口，输入命令：

fminbnd('fx',-10,-1)       %求函数在(-10，-1)内的最小值点和最小值

fminbnd(f,1,10)        %求函数在(1，10)内的最小值点。注意函数名f不用加'

例6.37  设

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image046.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

求函数f在(0.5,0.5,0.5)附近的最小值。

建立函数文件fxyz.m：

function f=fxyz(u)

x=u(1);y=u(2);z=u(3);

f=x+y.^2./x/4+z.^2./y+2./z;

在MALAB命令窗口，输入命令：

[U,fmin]=fminsearch('fxyz',[0.5,0.5,0.5])   %求函数的最小值点和最小值

例6.38  求解有约束最优化问题。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image048.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

首先编写目标函数M文件fop.m。

function f=fop(x)

f=0.4*x(2)+x(1)^2+x(2)^2-x(1)*x(2)+1/30*x(1)^3;

再设定约束条件，并调用fmincon函数求解此约束最优化问题。

x0=[0.5;0.5];

A=[-1,-0.5;-0.5,-1];

b=[-0.4;-0.5];

lb=[0;0];

option=optimset; option.LargeScale='off'; option.Display ='off';

[x,f]=fmincon('fop ',x0,A,b,[],[],lb,[],[],option)

例6.39  设有初值问题：

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image050.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

y(0)=2

试求其数值解，并与精确解相比较(精确解为y(t)=file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image052.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />)。

    (1) 建立函数文件funt.m。

function y=funt(t,y)

y=(y^2-t-2)/4/(t+1);

(2) 求解微分方程。

t0=0;tf=10;

y0=2;

[t,y]=ode23('funt',[t0,tf],y0);   %求数值解

y1=sqrt(t+1)+1;              %求精确解

plot(t,y,'b.',t,y1,'r-');  %通过图形来比较

例6.40  已知一个二阶线性系统的微分方程为：

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image054.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

其中a=2，绘制系统的时间响应曲线和相平面图。

建立一个函数文件sys.m：

function xdot=sys(t,x)

xdot=[-2*x(2);x(1)];

取t0=0，tf=20，求微分方程的解：

t0=0;tf=20;

[t,x]=ode45('sys',[t0,tf],[1,0]);

[t,x]

subplot(1,2,1);plot(t,x(:,2));             %解的曲线，即t-x

subplot(1,2,2);plot(x(:,2),x(:,1))        %相平面曲线，即x-x’

axis equal

例6.41  设

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image056.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

将X转化为稀疏存储方式。

X=[2,0,0,0,0;0,0,0,0,0;0,0,0,5,0;0,1,0,0,-1;0,0,0,0,-5];

A=sparse(X)

例6.42  根据表示稀疏矩阵的矩阵A，产生一个稀疏存储方式矩阵B。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image058.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

A=[2,2,1;3,1,-1;4,3,3;5,3,8;6,6,12];

B=spconvert(A)

例6.43  求下列三对角线性代数方程组的解。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image060.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

B=[1,2,0;1,4,3;2,6,1;1,6,4;0,1,2];  %产生非0对角元素矩阵

d=[-1;0;1];     %产生非0对角元素位置向量

A=spdiags(B,d,5,5)     %产生稀疏存储的系数矩阵

f=[0;3;2;1;5];      %方程右边参数向量

x=(inv(A)*f)'       %求解

例7.1  当λ取何值时，齐次线性方程组：

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image002.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

有非零解。

syms lamda

       A=[1-lamda,-2,4;2,3-lamda,1;1,1,1-lamda];

       D=det(A);

       factor(D)

例7.2  求下列极限。

（1）file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image004.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />             （2）file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image006.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

（3）file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image008.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />        （4）file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image010.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

syms a m x;

f=(x^(1/m)-a^(1/m))/(x-a);

limit(f,x,a)                    %求极限(1)

f=(sin(a+x)-sin(a-x))/x;

limit(f)                       %求极限(2)

f=x*(sqrt(x^2+1)-x);

limit(f,x,inf,'left')               %求极限(3)

f=(sqrt(x)-sqrt(a)-sqrt(x-a))/sqrt(x*x-a*a);

limit(f,x,a,'right')                   %求极限(4)

例7.3  求下列函数的导数。

(1)file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image012.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" /> ，求y' 。                 (2)y=xcos(x)，求y''、y'''。

(3)file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image018.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />。            (4)file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image024.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />。

(5)z=f(x,y)由方程x²+y²+z²=a²定义，求file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image028.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />。

syms a b t x y z;

f=sqrt(1+exp(x));

diff(f)                  %求(1)。未指定求导变量和阶数，按默认规则处理

f=x*cos(x);

diff(f,x,2)               %求(2)。求f对x的二阶导数

diff(f,x,3)               %求(2)。求f对x的三阶导数

f1=a*cos(t);f2=b*sin(t);

diff(f2)/diff(f1)           %求(3)。按参数方程求导公式求y对x的导数

%求(3)。求y对x的二阶导数

(diff(f1)*diff(f2,2)-diff(f1,2)*diff(f2))/(diff(f1))^3

f=x*exp(y)/y^2;

diff(f,x)                 %求(4)。z对x的偏导数

diff(f,y)                 %求(4)。z对y的偏导数

f=x^2+y^2+z^2-a^2;

zx=diff(f,x)/diff(f,z)      %求(5)。

zy=diff(f,y)/diff(f,z)      %求(5)。

例7.4  在曲线y=x³+3x-2上哪一点的切线与直线y=4x-1平行。

x=sym('x');  

y=x^3+3*x-2;            %定义曲线函数

f=diff(y);                %对曲线求导数

g=f-4;

solve(g)                 %求方程f-4=0的根，即求曲线何处的导数为4

例7.5  求下列不定积分。

(1) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image030.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />                   (2) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image032.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

(3）file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image034.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />                         (4) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image036.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

x=sym('x');

f=(3-x^2)^3;

int(f)                     %求(1)

f=sin(x)^2;

-1/2*cos(x)*sin(x)+1/2*x

syms alpha t;

f=exp(alpha*t);

int(f)                     %求(3)

f=5*x*t/(1+x^2);

int(f,t)                   %求(4)

例7.6  求下列定积分。

 (1)file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image038.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />                     (2) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image040.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

(3) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image042.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />                  (4) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image044.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

x=sym('x');t=sym('t');

int(abs(1-x),1,2)                 %求(1)

f=1/(1+x^2);

int(f,-inf,inf)                    %求(2)

f=x^3/(x-1)^10;

I=int(f,2,3)                       %求(3)

double(I)                         %将上述符号结果转换为数值

int(4*x/t,t,2,sin(x))            %求(4)

例7.7  求椭球的体积。

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image046.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

syms a b c z;

f=pi*a*b*(c^2-z^2)/c^2;

V=int(f,z,-c,c)

例7.8  求空间曲线c从点(0，0，0)到点(3，3，2)的长度。设曲线c的方程是：

file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image048.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

syms t;

x=3*t;

y=3*t^2;

z=2*t^3;

f=diff([x,y,z],t);                  %求x,y,z对参数t的导数

g=sqrt(f*f');                        %计算一型积分公式中的根式部分

l=int(g,t,0,1)                       %计算曲线c的长度

例7.9  求函数y=file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image050.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />的傅立叶变换及其逆变换。

syms x t;

y=abs(x);

Ft=fourier(y,x,t)     %求y的傅立叶变换

fx=ifourier(Ft,t,x)   %求Ft的傅立叶逆变换

例7.10  计算y=x²的拉普拉斯变换及其逆变换.

x=sym('x');

y=x^2;

Ft=laplace(y,x,t)           %对函数y进行拉普拉斯变换

fx=ilaplace(Ft,t,x)         %对函数Ft进行拉普拉斯逆变换

例7.11  求数列 fn=e^-n的Z变换及其逆变换。

syms n z

fn=exp(-n);

Fz=ztrans(fn,n,z)          %求fn的Z变换

f=iztrans(Fz,z,n)          %求Fz的逆Z变换

例7.12  求下列级数之和。

(1) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image052.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

(2) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image054.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

(3) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image056.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

(4) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image058.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

n=sym('n');

s1=symsum(1/n^2,n,1,inf)            %求s1

s2=symsum((-1)^(n+1)/n,1,inf)        %求s2。未指定求和变量，默认为n

s3=symsum(n*x^n,n,1,inf)            %求s3。此处的求和变量n不能省略。

s4=symsum(n^2,1,100)              %求s4。计算有限级数的和

例7.13  求函数的泰勒级数展开式。

(1)求file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image060.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />的5阶泰勒级数展开式。

(2)将file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image062.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />在x=1处按5次多项式展开。

x=sym('x');

f1=sqrt(1-2*x+x^3)-(1-3*x+x^2)^(1/3);

f2=(1+x+x^2)/(1-x+x^2);

taylor(f1,x,5)                %求(1)。

taylor(f2,6,1)                  %求(2)。展开到x-1的5次幂时应选择n=6

例7.14   解下列方程。

(1)file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image064.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />                   (2) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image066.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

(3) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image068.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />                          (4) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image070.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

x=solve('1/(x+2)+4*x/(x^2-4)=1+2/(x-2)','x')   %解方程(1)

f=sym('x-(x^3-4*x-7)^(1/3)=1');

x=solve(f)                                            %解方程(2)

x=solve('2*sin(3*x-pi/4)=1')                      %解方程(3)

x=solve('x+x*exp(x)-10','x')                   %解方程(4)。仅标出方程的左端

例7.15  求下列方程组的解。

 

(1)file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image072.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />                     (2) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image074.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

(3) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image076.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />                     (4) file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image078.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />

命令如下：

[x y]=solve('1/x^3+1/y^3=28','1/x+1/y=4','x,y')         %解方程组(1)

[x y]=solve('x+y-98','x^(1/3)+y^(1/3)-2','x,y')        %解方程组(2)

例7.16  求下列微分方程的解。

(1)求file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image080.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />的通解。                 (2)求file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image082.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />的通解。

(3) 求file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image086.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />。       (4)求file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image088.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />的通解。

(5)求file:///C:/DOCUME~1/ADMINI~1/LOCALS~1/Temp/msohtmlclip1/01/clip_image090.gifMatlab 例题 6-7章" TITLE="刘卫国 Matlab 例题 6-7章" />的通解。

y=dsolve('Dy-(x^2+y^2)/x^2/2','x')     %解(1)。方程的右端为0时可以不写

y=dsolve('Dy*x^2+2*x*y-exp(x)','x')    %解(2)

y=dsolve('Dy-x^2/(1+y^2)','y(2)=1','x');     %解(3)

[x,y]=dsolve('Dx=4*x-2*y','Dy=2*x-y','t')    %解方程组(4)

[x,y]=dsolve('D2x-y','D2y+x','t');       %解方程组(5)