matlab数学实验

Posted on

Matlab数学实验学习笔记

MATLAB语言的基本语法

变量名的命名规则

  • 必须以字母开头
  • 区分大小写，可由字母、数字、下划线组成（变量名尽量反映其含义)

赋值语句

基本语法 变量名=表达式

1
2
3
4
5
a=[2 5 6 7 9];

a(2)=10

a=sin(pi/3), [v1,v2,v3]=fun(m1,m2)

表达式语句

一个语句可以只有表达式
系统自动将表达式的结果赋值给MATLAB内部变量”ans”。

语句的分隔符

  • 语句分隔符: 分号或逗号；

  • 语句的末尾不使用分号时，系统会显示执行结果。

    1
    a=3; b=4, c=a^2; a(1,4)=b

常用命令、快捷键

  • who 列出当前工作空间所有变量名称
  • whos 列出当前工作空间变量更多信息（维数,占用内存字节数等）
  • clc 清除命令窗口内容
  • clear 清除工作空间中的变量
  • 快捷键：向上方向键、向下方向键,用于浏览命令窗口历史命令、语句

数组的使用

  • 使用方括号

    同一行的元素用“空格或逗号”分隔,不同行的元素用“分号或换行”分隔

    1
    2
    3
    4
    5
    6
    x=[1 3 5 6 4]
    y=[1, 3, 5, 6, 4] 
    A=[1  2  3;  4  5  6;  7  8  9];
    B=[ 1  2  3
    	4  5  6
    	7  8  9];

  • 冒号操作符

    用于创建行向量 a:step:b ;
    其中a:b等同于a:1:b;

    a起始值，step增量，b用于判断向量终点值.

    1
    2
    3
    x=1:5      %表示x=[1 2 3 4 5]，增量默认为1
    x=1:2:9    %表示x=[1 3 5 7 9]   
    x=10:-2:1  %表示x=[10 8 6 4 2]

  • linspace(a,b,n)

    n-1等分区间[a, b]的节点组成的行向量(含区间端点a, b)；

    注意：如果要产生一个区间上的均匀节点,并且指定所产生数组的元素个数，则使用linspace更为方便。

    不指定n默认100个；

    1
    x=linspace(-2, 2, 5) %表示x=[-2 -1 0 1 2]

  • 拼接

    1
    2
    C=[A B]  %横向拼接要求A,B行数相同
    D=[A; B] %纵向拼接，要求A,B列数相同.

  • 空矩阵 [ ] 产生一个空矩阵

  • 调用函数创建

    1
    2
    3
    4
    5
    6
    7
    8
    a1 = zeros(m, n) 
    b1 = zeros(m) 

    a2 = ones(m, n) 
    b2 = ones(m)

    a3 = eye(m, n)
    b3 = eye(m)

  • 提取和修改数组中的元素

    • 通过数组下标访问：

      • ≥1的整数；
      • 不能越界

    • 获取子阵:
      获取第r行 A(r, :) ；
      获取第c列 A(:, c) ；
      获取子阵 A(行下标向量,列下标向量)
      示例：

      1
      2
      3
      4
      x=[1 2 3 4 5; 6 7 8 9 10;11 12 13 14 15] %创建矩阵
      y1=x([1, 2], :)                          % 取第1,2行
      y2=x([2 3],[1 3 4])                      % 取第2,3行,第1,3,4列
      y3=x([1 3 4])               % 取第1，3，4个数，列优先编号，[1 11 2]

    • 修改元素：用赋值语句修改；

  • 删除数组元素操作
    将空矩阵赋值给相应子阵达到删除；

    1
    2
    3
    4
    5
    6
    7
    8
    x=[1 5 9; 
       2 6 10; 
       3 7 11;
       4 8 12]          %创建4行3列矩阵x
    y=x([1 3 4],[2 3])  % 取第1,3,4行,第 2,3列赋值给y
    z =x([1 3 5])       %取x第1,3,5个元素赋值给z
    x([1 2],:) = [ ]    % 删除x的第1、2行

  • end用法
    end在下标表达式中表示最后一个下标值；

    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    x=[1  5  9; 2  6  10; 3 7 11; 4 8 12];  
    x(end,2)=  0; x%将矩阵x的最后一行第2列元素赋值为0

    t = rand(1,10);
    x1 = t(1:end-1) %取第1个-倒数第2个
    x2 = t(end-2:end)%取倒数第3个-倒数第1个

    A = rand(3)
    B = A(1:end-1, : )%取A的第1行-倒数第2行
    C = A(:, [2:end]) %取A的第2列-倒数第1列

算术运算符：矩阵操作运算符

  • 矩阵转置 B.' 或B' 矩阵共轭转置 $B^{T}$

  • 矩阵加减：A+B，A-B，($A,B$同型，或有一个为标量）

  • 矩阵相乘：A * B， $A$与$B$为矩阵或有一个为标量*

  • 矩阵左除：A\B， 当$A$为方阵表示: $A^{-1} B$

  • 矩阵右除：A/B， 当$B$为方阵表示 $AB^{-1}$ ，或$B$为标量

  • 矩阵幂： A^n，$ A$为方阵\

  • inv函数：inv(A) 求$A$的逆矩阵$A^{-1}$

    1
    2
    3
    4
    5
    6
    A = rand(3,3);   
    B = rand(3,3);
    C1 = A\B;    
    C2 = A/B
    E1 = C1-inv(A)*B %inv函数求矩阵的逆
    E2 = C2-A*inv(B) %E1,E2 理论上为零矩阵

数组对应元素操作运算符

  • 数组相乘：C=A.*B

  • 数组右除：C=A./B；

  • 数组左除：C=A.\B

  • 数组幂：C=A.^B（要求：A, B同维数或其中之一为标量）

    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    A=rand(3,4);
    B=rand(3,4); 
    A.*B       %两个矩阵对应元素相乘
    A./B	   %两个矩阵对应元素相除(a(i,;,j)/b(i,:,j))
    A.\B	   %两个矩阵对应元素相除(b(i,;,j)/a(i,:,j))
    A.^B	   %两个矩阵对应元素幂运算(a(i,;,j)^b(i,:,j))
    T1=A.*2;   % 以A的每个元素与2相乘构造数组
    T2=A.^2;   % 以A的每个元素2次方构造数组
    T3=2./A;   % 以A的每个元素的倒数乘以2构造数组
    T4=2.\A;   % 以2的倒数乘以A的每个元素构造数组

算术运算符：标量之间的运算符

标量加减乘除：+ - * / \ ^

逻辑运算符：& | ~ 

1
2
3
4
5
6
7
8
9
10
11
12
x =[4     9     1     2     1];
y =[1     8     5     5     1];
t1=x>=2		% t1 =     1     1     0     1     0
u=x-y		% u  =     3     1    -4    -3     0
t2=x-y>=1	%t2 =     1     1     0     0     0

a=3; b=7;
if a>3 & b>3, 
disp('handle 1')
else 
disp('handle 2')%输出
end

数据类型

主要的数据类型：double char sym struct cell

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
a=rand(3); b='Li San';   %double型，char型
symsx, y=1 + x^2         %x,y为sym类型;对y赋值的语句含符号对象
F.name='li San', F.birth=1999, F.src=rand(3) %F为struct型
whos a b x y F 
class(a) % 'double'
class(b) % 'char'
%whos结果
Name      Size            Bytes  Class     Attributes

F         1x1               620  struct

a         3x3                72  double              

b         1x6                12  char                

x         1x1               112  sym

y         1x1               112  sym

cell数组基本用法

  • 创建数组用法：a=cell(m,n)

  • 特点：一个cell数组中的元素的类型可以互不相同

  • 存取cell数组用法示例：

    1
    2
    3
    4
    5
    a=cell(2,3)
    a{1,1}='abc';a{1,2}=rand(3);a{1,3}=cell(1,2);

    % a = 'abc'    [3x3 double]    {1x2 cell}
    %     []       []              []

基本输入与格式化输出操作函数

  • input键盘输入函数

    • input(提示信息字符数组)

      1
      2
      3
      g=input('输入您的成绩: ')
      %输入您的成绩: 95
      %g =    95

  • disp 显示数组内容函数
    特点：显示数组内容，但不输出变量名,多用于调试程序时显示数组内容

    示例：

    1
    2
    a=rand(1,3);
    disp(a)

  • sprintf 将数组内容格式化为字符串

    功能：将数据格式化输出为字符串

    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    a =[pi,sqrt(2)];
    name='Li San';
    grade=[86 95 89]; 
    s1= sprintf('%0.5f',pi)    %将该实数化为5个小数位字符串
    s2= sprintf('%10.6f', a)   %将a化为10个字符长，含6个小数位的字符串
    s3=sprintf('%8s%3d%3d%3d',name,grade)

    %运行结果：
    s1 =   3.14159 
    s2 =   3.141593    1.414214
    s3 =   Li San 86 95 89

Matlab中的常用函数

eval将字符串转化为Matlab语句执行

1
2
3
4
for i=1:20
 eval(sprintf('syms x%d',i))
end
whos

length获取向量长度函数

s=length(v) 获取向量长度s

sum求和函数

  • s=sum(v) 求向量v中元素的和
  • s=sum(A,1) s=sum(A) 求矩阵A中每列的和，返回成1个行向量
  • s=sum(A,2) 求矩阵A中每行的和，返回成1个列向量
1
2
3
4
5
6
7
8
9
10
11
12
13
a= fix(1000*rand(3,5))/1000;
s=sum(a), s2=sum(a,1) % s, s2完全相同
t =sum(a’,1)
t2= sum(a,2)          % t,t2元素值相同，差别在于行/列向量

%运行输出：
s =   1.00700   0.34800   2.46700   0.87400   1.97500
s2 =  1.00700   0.34800   2.46700   0.87400   1.97500
t =   1.8430   2.2620   2.5660
t2 = 1.8430
	2.2620
	2.5660

mean求平均值函数

  • s=mean(v) 求向量v中元素的平均值
  • s=mean(A,1) s=mean(A) 求矩阵A中每列的平均值，返回成1个行向量
  • s=mean(A,2) 求矩阵A中每行的平均值，返回成1个列向量
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
X= [1 2 3; 3 3 6; 4 6 8; 4 7 7]
v1=mean(X,1) 
v2=mean(X,2)

%运行输出：
X =   1   2   3
3   3   6
4   6   8
4   7   7
v1 =   3.0000   4.5000   6.0000
v2 =
 2
 4
 6
 6

max求最大值函数

[v, I]=max(x)

  • 如果x为向量，v为向量中的最大元素；I为最大元素在x中的下标。
  • 如果x为矩阵，v为每列的最大元素组成的行向量，I则为每列最大元素的行下标组成的向量。

min求最小值函数

[v, I]=min(x)

  • 如果x为向量，v为向量中的最小元素；I为最小元素在x中的下标。
  • 如果x为矩阵，v为每列的最小元素组成的行向量，I则为每列最小元素的行下标组成的向量。

sort排序函数

  • [B, I]=sort(v) 对向量v中元素排序;
    B为按递增排序后的元素；
    I为排序后数组B中的元素在原数组v中的位置下标.
  • 当v是矩阵时，sort(v)或sort(v,1)对列分别排序，
    sort(v,2)对行分别排序。
1
2
3
4
5
6
7
a=[1 5 6 7 4 9]
v=sort(a)		  % v =1 4 5 6 7 9
[v2,idx2]=sort(a) % v2=1 4 5 6 7 9
				 %idx2=1 5 2 3 4 6
v_byidx = a(idx2) %v_byidx = 1 4 5 6 7 9
v1_c=sort(a, 'ascend') %升序排列
v2_c=sort(a, 'descend')%降序排列

find查找函数

find函数用于查找数组中的非零元素位置,结合逻辑表达式可以返回所需要元素的位置下标.

1
2
3
4
5
6
a=[1 5 6 7 4 9];
idx=find(a>=6)  % idx = 3     4     6

x =[1 9 1 7 8];
y =[8 0 3 2 8];
s = sum((x>=y-2) & (x<=y+2) )%s= 2

Matlab函数文件

函数名命名规范

  • 必须以字母开头，区分大小写，函数名可由字母、数字、下划线组成
  • 保存函数的文件名一般与程序中定义的函数名相同；否则，调用该函数时以保存函数的文件名作为函数名来调用；
  • 函数名要有一定实际含义，便于记忆，同时避免与系统内部函数名相同。
1
2
3
4
5
6
function r=sumbetween(a,b)
%sumbetween.m任意自然数a和b之间（含a和b）所有整数的和
if a <= b,
r =(b-a+1)*(a+b)/2;
else r =(a-b+1)*(a+b)/2;
end

nargin,nargout，return和pause

  • nargin，nargout分别表示函数调用时的输入参数个数，输出参数个数。
  • return 返回调用函数
  • pause 暂停运行，按任意键执行
  • pause(n) 暂停n秒，如pause(0.5)

inline函数和匿名函数

inline函数根据expr建立内联函数，函数自变量符号根据表达式自动确定；

1
2
3
4
f1 = inline('t^2')   % 系统自动确定参数
f2 = inline('x*y*z') % 系统自动确定参数
f3 = inline('x*y*z','x','y','z') % 指定参数
f4 = inline(‘x^P1+x^P2’,2) % 指定参数x,P1,P2

匿名函数使用@符号：@(参数列表)(函数表达式)，返回与输入参数同维数的数组。

1
2
3
4
f=@(x,y)((x.^2+y.^2<=1).*(x+y)+(x.^2+y.^2>1).*(x-y));
x=-2:0.1:2;
y=-2:0.1:2;
v=f(x,y);

MATLAB 一元函数绘图方法

  • plot() 基本绘图方法，利用一元函数自变量的一系列数据和对应函数值数据绘图。具有很大灵活性;
  • ezplot()简易函数绘图方法，优点：快速方便
  • fplot() 函数绘图方法,与ezplot相似.（不同版本要求不同，2019版不能给出函数值变化范围）
  • polar() 极坐标绘图命令
  • bar() 绘条形图命令
  • bar3() 三维条形图绘制命令
1
2
3
4
5
6
7
8
9
plot(x,y) 
plot(x,y,'r-') 
plot(x1,y1,x2,y2)
plot(x1,y1,'r',x2,y2,'b')
ezplot(f,[xmin,xmax,ymin,ymax])
fplot(f,[xmin,xmax])
polar(theta,rho,'r.')
bar(x,y)
bar3(x,y)
  • legend在图形上生成一个图例框
  • title在图形顶部添加标题
  • xlabel用字符串标注x轴
  • ylabel用字符串标注y轴
  • gca返回当前图形坐标区坐标轴句柄
  • gcf返回当前图窗句柄
  • subplot另开一个绘图窗口

实例

1
2
3
4
x=0:0.1:4*pi; 
y= exp(-0.5*x) ;  y1=y .*sin(5*x);
plot(x,y1,x,y,'--r',x,-y,'--r')
legend('exp(-0.5x)sin(5x)','exp(-0.5x)')
1
2
3
4
5
6
x=linspace(0,2*pi,30);
y=sin(x);
z=cos(x);
plot(x,y,x,z)
legend('sin(x)','cos(x)')
title('Figure: Legend Example')
1
2
3
4
h=ezplot('exp(-0.5*x)*sin(5*x)',  [0,10,-1,1])
set(h, 'linewidth', 2, 'color', 'r')
set(gca,'fontsize', 14, 'fontname', '宋体')
set(gcf, 'color', 'c')
1
2
3
theta=0:pi/20:4*pi;
rho= 1 + theta*2;
polar(theta,rho,'r') %r =1+2θ
1
2
3
4
5
6
7
8
x=-2.9:0.2:2.9;
y= exp(-x.*x);
subplot(1,2,1)
bar(x,y)
title('Figure:2-D Bar Chart')
subplot(1,2,2)
bar3(x,y)
title('Figrure:3-D Bar Chart')

MATLAB 二元函数绘图方法

  • 空间曲线绘制plot3

    1
    2
    3
    4
    5
    6
    t=linspace(0,10*pi);
    plot3(sin(t),cos(t),t)
    xlabel('sin(t)')
    ylabel('cos(t)')
    zlabel('t')
    title('Figure: helix')

  • 生成平面网格点命令: [X, Y]=meshgrid(x, y)

    1
    2
    3
    x=1:6; y=1:8;       %创建两个向量
    [X,Y]=meshgrid(x,y)%将x和y分别扩充为8行6列
    [X,Y]=meshgrid(1:6,1:8)%直接创建两个矩阵X和Y

  • 绘网格命令：mesh(x,y,z) 或 mesh(z)

    1
    2
    3
    4
    [x,y]=meshgrid(-2:0.2:2);
    z=x.*exp(-x.^2-y.^2);  
    mesh(x,y,z)
    colormap([0 0 1])

    绘网格命令 meshc() meshz()

    1
    2
    3
    4
    5
    6
    7
    8
    [x,y]=meshgrid(-2:0.2:2);
    z=x.*exp(-x.^2-y.^2);  
    subplot(1,2,1)
    meshc(x,y,z)
    title('Figure 1:mesh plot with Contours')
    subplot(1,2,2)
    meshz(x,y,z)
    title('Figure 2:mesh plot with Zero Plan')

  • 等高线图contour和countour3

    1
    2
    3
    4
    5
    6
    7
    8
    [x,y]=meshgrid(-2:0.2:2);
    z=x.*exp(-x.^2-y.^2);  
    subplot(1,2,1)
    contour(x,y,z,20)
    title('Figure 1: 2-D Contour Plot')
    subplot(1,2,2)
    contour3(x,y,z,20)
    title('Figure 2: 3-D Contour Plot')

Matlab中的字符串

char(S1,S2,…) 利用给定的字符串或单元数组创建字符数组
double(S) 将字符串转化成ASC码形式
cellstr(S) 利用给定的字符串数组创建字符串单元数组
blanks(n) 生成一个由n个空格组成的字符串
deblank(S) 删除尾部的空格
eval(S),evalc(S) 使用MATLAB解释器求字符串表达式的值
ischar(S) 判断是不是字符串数组
iscellstr(C) 判断是不是字符串单元数组
isletter(S) 判断是不是字母
isspace(S) 判断是不是空格
strcat(S1,S2,…) 将多个字符串水平拼接
strvcat(S1,S2,…) 将多个字符串竖直拼接

strcmp(S1,S2) 判断字符串是否相等
strncmp(S1,S2,n) 判断前n个字符串是否相等
strcmpi(S1,S2) 判断字符串是否相等（忽略大小写）
strncmpi(S1,S2,n) 判断前n个字符串是否相等（忽略大小写）
strtrim(S1) 删除结尾的空格
findstr(S1,S2) 查找
strfind(S1,S2) 在S1查找S2
strjust(S1,type) 按照指定的type调整下一个字符串数组
strmatch(S1,S2) 查找要求的字符串下标
strrep(S1,S2,S3) 将字符串S1中出现的S2用S3代替
strtok(S1,D) 查找S1中第一个给定的分隔符之前和之后的字符串
upper(S) 将一个字符串转换成大写
lower(S) 将一个字符串转换成小写

num2str(k) 将数字转换成字符串
int2str(k) 将整数型转换成字符串
mat2str(k) 将矩阵转换为字符串，供eval使用
str2double(S) 将字符串数组转换为数值数组
sprintf(S) 创建含有指定格式的字符串
sscanf(S) 按照指定的控制格式读取字符串

字符串拼接函数：strcat和strvcat

strcat('MAT','LAB')

运行输出：ans = MATLAB
strvcat('MAT',‘MATLAB')

运行输出：ans =
MAT
MATLAB

字符串与数字数组转换函数：str2num和num2str

1
2
3
4
v=([1:0.5:2]*pi)
str=num2str(v)
v2=str2num(str)
err=v2-v

字符串查找函数：findstr和strfind

k=findstr(S1, S2)找到S1和S2两者中较短那个字符串在较长字符串中的位置起始索引号.
k=strfind(S1,pattern)返回pattern在S1中的位置索引号

Matlab文本文件操作

打开文件

在读写文件之前，必须先用fopen函数打开或创建文件，并指定对该文件进行的操
作方式。

fopen函数的调用格式为：fid = fopen(文件名，‘打开方式’)
说明：其中fid用于储存文件句柄值，如果返回的句柄值大于0，则说明文件打开成功。文
件名用字符串形式，表示待打开的数据文件。

打开文件的方式

  • ‘r’只读方式打开文件（默认方式），该文件必须已存在。
  • ‘r+’读写方式打开文件，打开后先读后写。该文件必须已存在。
  • ‘w’打开后写入数据。该文件已存在则更新；不存在则创建。
  • ‘w+’读写方式打开文件。先读后写。该文件已存在则更新；不存在则创建。
  • ‘a’在打开的文件末端添加数据。文件不存在则创建。
  • ‘a+’打开文件后，先读入数据再添加数据。文件不存在则创建。
  • 另外，在这些字符串后添加一个‘t’，如‘rt’或‘wt+’，则将该文件以文本方式打开；如果添加的是‘b’，则以二进制格式打开，这也是fopen函数默认的打开方式。

关闭文件

文件在进行完读、写等操作后，应及时关闭，以免数据丢失。

关闭文件用fclose函数，调用格式为：sta = fclose(fid)
说明：该函数关闭fid所表示的文件。sta表示关闭文件操作的返回代码，若关闭成功，返回0，否则返回-1. 如果要关闭所有已打开的文件用fclose(‘all’)

常见函数

fgets 读取一行，保留换行符
fgetl 读取一行

fprintf将数据按指定格式写入到文本文件中

1
2
3
4
5
6
fid=fopen(‘mytext.txt’,‘wt’);
for i=1:10
fprintf(fid,‘LINE%5d %10.6f\n’,i ,rand);
end
fclose(fid);
edit mytext.txt

Matlab符号变量计算

符号计算函数可以完成准确的推导、计算。如求极限、求导、求积分等。
计算中的符号变量需要先定义再调用符号计算相关函数进行处理。
符号计算部分功能

  • 复合运算、变量代换及化简
  • 线性代数：行列式，特征值等
  • 微积分：求导，求极限，定积分，微分方程的求解等。

定义符号变量syms或sym函数

1
2
3
4
5
syms arg1 arg2 ...
syms arg1 arg2 ... real
syms arg1 arg2 ... positive
syms arg1 arg2 [m, n] %arg1 arg2为m*n阶符号数组
syms f(x1,x2,x3,...) %f是x1, x2, x3, ... 的符号函数

对于函数sym(A)

  • 如果A是字符串，则产生一个符号变量；
  • 如果A是数值标量或数值矩阵，则其转为符号类型。
1
2
3
4
5
6
7
x1=sym(‘x1’) % 产生数组x=sym ('x', [m, n]) 
a=sym(sqrt(200)), v=sym([100 200]) 

%运行输出： x1 =  x1 
%a = 10*2^(1/2) 
%v=[ 100, 200]

符号表达式中的符号替换

subs函数将符号表达式s中的符号变量用其他符号表达式或数值代替，实现符号的替换，将s表达式中old变量替换为new，其使用格式为：s1=subs(s, old, new)

1
2
3
syms x y z
f=x^2+y^2+z^2
fval=subs(f,[x,y],[1,2])

simplify对表达式进行化简

1
2
3
syms x y
s1 = simplify(cos(x)^2-sin(x)^2) %s1 = cos(2*x)
s2 = simplify(x^3+3*x^2+3*x+1) %s2 = (x + 1)^3

符号计算精度及其数据类型转换

  • double(s) 将符号表达式s转化为双精度数值
  • vpa(s) 计算符号表达式s的数值结果
  • vpa(s,n) 采用n位有效数字计算精度求s的数值结果
  • digits 显示vpa计算结果的有效数字的位数
  • digits(n) 设置vpa计算结果的有效数字的位数
  • char(s) 将符号表达式s转化为字符串

compose复合函数

  • compose(f, g) 返回复合函数f(g(y))，其中f=f(x),g=g(y). x和y分别为findsym(symvar)从f、g中找到的符号变量.
  • compose(f, g, z) 返回复合函数f(g(z)), f=f(x), g=g(y), x,y含义同上一种用法. 最后用指定变量z代替变量y.
  • compose(f, g, x, y, z) 返回复合函数f(g(z)). 将x=g(y)代入f(x)中, 最后用指定的变量z代替变量y
1
2
3
4
syms x y t
f = 1/(1+x);
g = sin(y)^2
h = compose(f,g,x,y,t) %h = 1/(sin(t)^2 + 1)

limit极限函数

  • limit(f,x,a) 计算f(x)当x趋向于a的极限
  • limit(f,x,a,'right') 计算右极限
  • limit(f,x,a,'left') 计算左极限

diff求导计算函数

  • diff(s,v) 求s对自变量v的1阶导数
  • diff(s,v,n) 求s对自变量v的n阶导数 注：v 为指定的符号变量，可以缺省。
1
2
3
4
5
6
7
syms x y
f= x^2*y + 2*x*y + y*y
d1 = diff(f,x,1)
d2 = diff(f,y,1)
%运行结果：
%d1 = 2*y + 2*x*y
%d2 = x^2 + 2*x + 2*y

int符号积分函数

  • s=int(expr, var)以expr表达式中的变量var为积分变量计算不定积分
  • s=int(expr, var, a, b)以expr表达式中的变量var为积分变量计算定积分，积分上下限分别为b和a

symsum级数求和函数

使用格式为： s=symsum(expr,v,a,b)

solve求解方程(组)

直接列出符号变量的方程，注意用==连接

1
2
3
4
5
6
7
8
9
10
11
syms x y a b
s=solve(a*x+b*y==10,a*x-b*y==20,x,y)
sol_x = s.x
sol_y = s.y

%运行结果：
%s = 
%x: [1x1 sym]
%y: [1x1 sym] 
%sol_x =   15/a
%sol_y =   -5/b

dsolve求解微分方程

调用dsolve函数时，应显示指定自变量，用一个
自变量字符组成的字符串作为最后一个参数传
入；

  1. 如果不指定自变量，默认自变量为t
1
2
3
syms y(x)
d1=diff(y,1);
s=dsolve(d1==(50-0.01*y)*y,y(0)==4,x)
1
y=dsolve('Dy=(50-0.01*y)*y','y(0)=4','x')

梯形法求解积分

Q = trapz(X,Y)

对由表格形式定义的函数使用梯形法进行求解积分。

1
2
3
4
5
x=0:pi/100:pi;
y=sin(x);
trapz(x,y)
%ans =
%1.9998

基于变步长Simpso法求积分

q = quad(fun,a,b,tol)

fun 被积函数文件名或函数句柄但必须是数组对应元素操作运算符;

a, b 积分下限,积分上限;

tol 积分精度,缺省默认为10^(-6);

矩形区域二重数值积分

q = dblquad(fun,a,b,c,d,tol)

fun 被积函数文件名或函数句柄

a, b 第一个自变量积分下限,积分上限

c, d 第二个自变量积分下限，积分上限

tol 积分精度

Matlab微积分实验：泰勒展开

  • 实验任务1. 请以函数$f (x)=e^x$为例，验证随着n阶泰勒多项式$Pn(x)$的阶数增加，$Pn(x)$近似函数ex的程度越高。

    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    12
    13
    14
    syms x
    fx = exp(x); 
    x0 = 0; %展开点
    n=[2 5 7];
    for i=1:length(n) %展开阶数: n-1阶
    hx(i) = taylor(fx,x,'order',n(i),'ExpansionPoint',x0)    
    end
    xp = linspace(-3,3,40); %创建节点数组
    y  = subs(fx, x, xp);   %通过subs替换符号计算函数值
    y1 = subs(hx(1), x, xp); y4 = subs(hx(2), x, xp);
    y6 = subs(hx(3), x, xp);
    h=plot(xp,y,'r',xp,y1,'k',xp,y4,'b',xp,y6,'c','linewidth',2)
    legend(h,'source','n=1','n=4','n=6')
    set(gca,'fontsize',16)

  • 实验任务2. 利用指数函数$e^x$的泰勒多项式计算数学常数e的近似值，并绘制误差曲线.

    1
    2
    3
    4
    5
    6
    7
    8
    9
    fac=1; %初始化, 存储各项阶乘值
    n = 5; %设定n, 
    s(1) = 2; %求和变量初始化s(n)
    for i= 2:n
    fac = fac*i;%计算i!
    s(i) = s(i-1) + 1/fac;
    end
    err=abs(s-exp(1))
    semilogy(1:n,err,'-o')

Matlab微积分实验：求解非线性方程

  • solve,roots,fzero命令的应用

    1
    2
    syms x
    s=solve((x+2)^x==2, x)

    对于P为按降幂排列的多项式系数, R列向量。

    1
    2
    P=[1 0 0 1];%求方程x3+1=0的根
    R=roots(P)

    fzero的典型调用格式：

    • [x,fval,exitflag]=fzero(fun,x0)
    • [x,fval,exitflag]=fzero(fun,[a,b])
    1
    2
    3
    4
    5
    syms x
    f=2*x^3-3*x^2+4*x-5
    fun=inline((f), 'x')
    x=-5:0.1:10;  plot(x, fun(x))
    [x,val,flag]=fzero(fun,2)

  • 二分法及其实现

    1
    2
    3
    4
    5
    6
    7
    8
    9
    10
    11
    12
    13
    14
    15
    16
    17
    18
    19
    20
    21
    22
    23
    24
    25
    function [x,n]=bisection(f,a,b,d) 
    %二分法求方程的根：可以是代数方程，也可以是超越方程
    %f是所要求解的函数
    %a,b是求解区间
    %d是求解精度
    %输出x是方程的近似根，n是迭代步数 
    fa=f(a);fb=f(b);
    n=1;%计数器初始化
    if fa*fb>0,
    error('所给区间内没有实数解！');
    end
    while 1%未知循环次数，不好用for循环
    c=(a+b)/2;fc=f(c);
    if abs(fc)<d;
    x=c;
    break;
    elseif fa*fc<0,
    b=c;
    fb=fc;
    else
    a=c;
    fa=fc;
    end
    n=n+1;
    end

  • Newton法迭代公式

    image-20240105221822510

image-20240105221935381

Matlab微积分实验：多项式拟合

  • 多项式拟合函数：polyfit( )
    调用格式： p=polyfit(x,y,n)
    说明：x,y为数据点，n为多项式阶数，返回p为幂次从高到低的多项式系数向量p。
  • 多项式求值函数：polyval( )

​ 调用格式： y=polyval(p,x)
​ 说明：为返回对应自变量x在给定系数p的多项式的值。

  • 非线性最小二乘拟合：lsqcurvefit( )

    [c,resnorm,residual,exitflag]=lsqcurvefit(‘fun’, x0, xdata, ydata, lb, ub).

    fun为事先建立的函数F(c, x)的M-文件或匿名函数,;

    x0为迭代初始值;
    lb,ub：解向量的上下界限，lb <= x0 <= ub，若没有要求，就设为空lb=[ ]，ub=[];
    c:拟合的最优解;

    resnorm:残差平方和;
    residual:残差向量;

    exitflag: 条件退出标示.

    1
    2
    3
    4
    5
    6
    7
    xdata =[0.9  1.5  13.8  19.8  24.1  28.2  35.2  60.3  74.6 81.3];
    ydata =[455.2  428.6  124.1  67.3  43.2  28.1  13.1  -0.4 -1.3 -1.5];
    fun = @(c,x)c(1)*exp(c(2)*x);
    x0 = [100,-1];
    [c, resnorm, residual, exitflag]= lsqcurvefit(fun, x0, xdata, ydata)
    %要计算q=3:0.4:24的函数值:
    Y=fun(c, q)

Matlab求解规划问题

  • 线性规划问题

    [x, fval] = linprog(c,A,b,Aeq,beq,Lb,Ub)

    image-20240122213537579

    1
    2
    3
    4
    5
    6
    7
    8
    9
    C = [500,600];
    A = [-1,-1;1,0;0,1];
    b = [-5000;5000;3000];
    Aeq = [];
    beq = [];
    Lb = [0;0];
    Ub = [inf;inf];
    x = linprog(-C, A, b, Aeq, beq, Lb, Ub)
    z = C*x

  • 单变量无约束（边界约束）最优化问题

    [xmin, ymin] = fminbnd(fun, x1, x2)

    fun是目标函数

    [x1,x2]是搜索区间

    xmin,ymin分别是目标函数的极小点、极小值。

    如果想求最大值，只需要取相反数即可；

  • 求多变量无约束最优化问题的最优解

    [xmin, fmin] = fminsearch(fun, x0)
    [xmin, fmin] = fminunc(fun, x0)
    fun是目标函数;

    x0是初始迭代点，
    xmin,fmin分别是目标函数的极小点、极小值。
    fminsearch适用范围更广，但处理问题的规模较小，求解速度较慢。fminunc求解更高效。

  • 非线性规划数学模型

    [xmin, fmin] = fmincon(fun,x0,A,b,Aeq,beq,Lb,Ub,nonlcon)

    fun是目标函数

    x0是初始迭代点

    nonlcon是针对非线性约束编写的函数的文件名。

  • 随机优化算法，

    解决困难、复杂、多态最优化问题的全局极小点的有效方法，即优化问题或高维、或不可导、或具有许多局部极小点、或目标函数不具有显式解析表达式等。
    [xmin, fmin] = ga(fun, n, A,b,Aeq,beq,Lb,Ub,nonlcon)
    fun是目标函数

    n是决策变量个数

    nonlcon是针对非线性约束编写的函数的文件名。
    注：随机优化算法，往往需要多次运行程序。

Matlab概率统计应用实验

随机数的生成

rand 生成$(0,1)$区间上均匀分布的随机数
unifrnd 生成指定区间内均匀分布的随机数
randn 生成服从标准正态分布的随机数
normrnd 生成指定均值、标准差的正态分布随机数
exprnd 生成服从指数分布的随机数
注：rand是采用线性同余法得到的，具有周期性，所以上述命令常被称为伪随机数生成器.
rand(‘seed’, number)

接口调用以rand为例：

rand(m) 生成$m\times m$维的随机数
rand(m,n)生成$m\times n$维的随机数
rand([m,n,p ...]) 生成排列成$m\times n \times p…$ 多维向量的随机数

蒙特卡罗方法

利用随机数模拟随机变量，求解一系列诸如面积，体积，值域，(非)线性规划相关问题；举一个例子(估计圆周率)：

1
2
3
4
5
6
7
8
9
n=input('请输入产生点的个数：');
m=0;
for i=1:n
x=-1+2*rand; y=-1+2*rand;
if x^2+y^2<=1
m=m+1;
end
end
mypi=4*m/n