MATLAB函数使用

>> sym s z
>> f=0.632*z^2/(z^3-1.736*z^2+1.104*z-0.368)
f =
79/125*z^2/(z^3-217/125*z^2+138/125*z-46/125)
>> r=taylor(f)
r =
-79/46*z^2-237/46*z^3-15563/2116*z^4-5135/2116*z^5 >> r=vpa(r,3)
r =
-1.72*z^2-5.15*z^3-7.35*z^4-2.43*z^5
如果你是想展开成负指数形式：
复制内容到剪贴板
代码:
clear;clc;
syms z
f=0.632*z^2/(z^3-1.736*z^2+1.104*z-0.368);
r1=subs(f,z,1/z);
r1=taylor(r1);
r1=subs(r1,z,1/z)
disp('或者')
syms z
f=0.632*z^2/(z^3-1.736*z^2+1.104*z-0.368);
r2=subs(taylor(subs(f,z,1/z)),z,1/z)
结果：
复制内容到剪贴板
代码:
r1 =
79/125/z+17143/15625/z^2+2357281/1953125/z^3+272594477/244140625/z^4+ 30811435509/30517578125/z^5
或者
r2 =
79/125/z+17143/15625/z^2+2357281/1953125/z^3+272594477/244140625/z^4+
30811435509/30517578125/z^5
化简：
复制内容到剪贴板
代码:
r1=vpa(r1,3)
r2=vpa(r2,3)
r1 =
.632/z+1.10/z^2+1.21/z^3+1.12/z^4+1.01/z^5
r2 =
.632/z+1.10/z^2+1.21/z^3+1.12/z^4+1.01/z^5小整理：适用於向量的常用函数有：
min(x): 向量x的元素的最小值
max(x): 向量x的元素的最大值
mean(x): 向量x的元素的平均值
median(x): 向量x的元素的中位数
std(x): 向量x的元素的标准差
diff(x): 向量x的相邻元素的差
sort(x): 对向量x的元素进行排序（Sorting）length(x): 向量x的元素个数
norm(x): 向量x的欧氏（Euclidean）长度sum(x): 向量x的元素总和
prod(x): 向量x的元素总乘积
cumsum(x): 向量x的累计元素总和
cumprod(x): 向量x的累计元素总乘积
dot(x, y): 向量x和y的内积
cross(x, y): 向量x和y的外积
如何在matlab中调用maple
（不用安装maple软件就可调用）
方法1：
m aple(’m aplestatem ent’)
其中m aplestatem ent 是完整的m aple语句，由一条或几条命令组成，必须符合m aple 的语法方法2：
m aple(’function’，arg1, arg2,…)
其中function为m aple中的函数名称，arg1, arg2,…是函数function所用的参数。

注：如果方法1行不通，可尝试方法2（个人经验）。

基本代数部分
如何用matlab求阶乘
factorial（n）求n的阶乘
如何用matlab配方
没有发现m atlab有这一命令，不过我们可以调用m aple的命令，调用方法如下：
首先加载m aple中的student函数库，加载方法为：m aple（’with（student）’）
然后运行m aple中的配方命令，格式为：
m aple（’com pletesquare（f）’）把f配方，其中f为代数表达式或代数方程
m aple(’com pletesquare（f，x）’)把f按指定的变量x配方，其中f同上
m aple(’com pletesquare（f，{x，y，．．．}）’)把f按指定的变量x，y，．．．配方
m aple(’com pletesquare（f，[x，y，．．．]）’)把f按指定的变量x，y，．．．配方，
如何用matlab进行多项式运算
（1）合并同类项
sym s 表达式中包含的变量 collect(表达式，指定的变量)
（2）因式分解
sym s 表达式中包含的变量 factor(表达式)
（3）展开
sym s 表达式中包含的变量 expand(表达式)
我们也可在m atlab中调用m aple的命令进行多项式的运算，调用格式如下：
m aple（’m aple中多项式的运算命令’）
如何用matlab进行分式运算
发现m atlab只有一条处理分式问题的命令，其使用格式如下：
[n，d]=num den（f）把符号表达式f化简为有理形式，其中分子和分母的系数为整数且分子分母不含公约项，返回结果n为分子，d为分母。

注意：f必须为符号表达式
不过我们可以调用m aple的命令，调用方法如下：
maple（’denom（f）’）提取分式f的分母
maple(’numer（f）’)提取分式f的分子
maple(’normal（f）’ )把分式f的分子与分母约分成最简形式
maple(’expand（f）’)把分式f的分子展开,分母不变且被看成单项。

maple(’factor（f）’)把分式f的分母和分子因式分解，并进行约分。

如何用Matlab进行因式分解
syms 表达式中包含的变量 factor(表达式)
如何用Matlab展开
syms 表达式中包含的变量 expand(表达式)
如何用Matlab进行化简
syms 表达式中包含的变量 simplify(表达式)
如何用Matlab合并同类项
syms 表达式中包含的变量 collect(表达式，指定的变量)
如何用Matlab进行数学式的转换
调用Maple中数学式的转换命令，调用格式如下：
maple(‘Maple的数学式转换命令’)
即：maple(‘convert(表达式，form)’)将表达式转换成form的表示方式
maple(‘convert(表达式，form, x)’)指定变量为x，将依赖于变量x的函数转换成form的表示方式（此指令仅对form为exp与sincos的转换式有用）
如何用Matlab进行变量替换
syms 表达式和代换式中包含的所有变量 subs(表达式，要替换的变量或式子，代换式)
如何用matlab进行复数运算
a+b*i 或 a +b*j表示复数a+bi 或 a+bj
real（z）求复数z的实部
imag（z）求复数z的虚部
abs（z）求复数z的模
angle（z）求复数z的辐角，
conj（z）求复数z的共轭复数
exp（z）复数的指数函数，表示e^z
如何在matlab中表示集合
[a, b, c,…]表示由a, b, c,…组成的集合（注意：元素之间也可用空格隔开）
unique(A) 表示集合A的最小等效集合（每个元素只出现一次）
也可调用m aple的命令，格式如下：
m aple('{a, b, c,…}')表示由a, b, c,…组成的集合
下列命令可以生成特殊的集合：
m aple(‘{seq(f(i),i=n..m)}’)生成集合{f(n), f(n+1), f(n+2), … , f(m)}
如何用Matlab求集合的交集、并集、差集和补集
我们也可以调用Maple的相应功能，调用方法如下：
m aple('用Maple求集合的交集、并集、差集和补集的命令 ')
具体地说，共有下列几个调用命令：
m aple( 'A union B' ) 求集合A和B的并集
m aple( 'A intersect B ' ) 求集合A和B的交集
m aple('A minus B ' ) 求差集A-B
m aple( '全集I m inus A ' ) 求集合A关于全集I的补集
如何用matlab排序
sort（v）将向量v的元素从小到大排列（升序排列）
MATLAB做拉式变换
2009-08-14 12:09
Poly(vector1)
Conv(vector1,vector2)
Roots(vector1)
Residue(num,den)
Syms: symbolic objects
Laplace()
Interp1(y,t,0.1)
s=tf('s')
G=1/s^2
G_z=c2d(G,0.1,'zoh')
[num,den]=tfdata(G_z,'v')
Tf2ss(num,den)
利用Matlab实现多项式互质分解的Euclid算法
摘要：利用多项式工具箱实现Euclid算法，完成了基于Matlab语言的系统特征多项式互质分解。

目前的Matlab鲁棒控制工具箱主要是针对状态空间矩阵描述方法而设计的，而对于单输入单输出(SISO)系统，在经典控制理论和工程中都更习惯采用传递函数形式来分析与综合系统。

因此，开发针对多项式的鲁棒控制工具，可更好地将鲁棒控制思想应用于工程实践。

在SISO系统鲁棒控制器设计过程中，首先一步就是对系统的特征多项式进行互质分解，利用Euclid算法[1]可实现。

在计算过程中需要进行多项式的四则运算及求解余数多项式，计算量较大，手算时易出错。

而利用其他语言实现多项式的四则运算，首先要自己编写多项式的运算规则，前期准备工作多。

本文利用[2]处所提供的多项式运算工具箱中给出的多项式四则运算函数，在Matlab中通过Euclid算法实现了多项式互质分解，提高了控制器设计的效率，值得推广。

1 互质分解
对于如图1所示的单位反馈系统，这里的P是严格正则的（传函分母阶次大于分子阶次），C是正则的（传函分母阶次大于等于分子阶次）。

图1：单位反馈系统
现在假定P是不稳定的，我们想要找到一个可镇定的C。

现作如下尝试，将P写成互质多项式比的形式，有
P=N/M
根据Euclid算法可得其他两个多项式X、Y，满足方程
NX+MY=1
Euclid算法计算两个给定的多项式n(λ)，m(λ)的最大公约数。

当n和m是互质的，该算法可以计算x(λ)和y(λ)满足
nx+my=1
步骤A：Euclid算法
输入：多项式n，m
初始化：如果阶次n≥阶次m不成立，交换n和m。

第1步用m除n得商q1和余式r1：
n=mq1+r1，
阶次r1<阶次m。

第2步用r1除m得商q2和余式r2：
m=r1q2+r2，
阶次r2<阶次r1。

第3步用r2除r1：
r1=r2q3+r3，
阶次r3<阶次r2。

继续。

当rk为非零常数，在第k步结束。

下面的例子说明如何求得x和y 。

这里k=3。

方程如下：
n=mq1+r1，
m=r1q2+r2，
r1=r2q3+r3，
即
解出r3（如用高斯消去法）：
r3=(1+q2q3)n+[-q3-q1(1+q2q3)]m
取
x=(1+q2q3)/r3, y=[-q3-q1(1+q2q3)]/r3
例题对于n(λ)= λ2，m(λ)=6λ2－5λ+1，算法如下：
q1(λ)=1/6, r1(λ)=5/6λ-1/6, q2(λ)=36λ/5-114/25, r2(λ)=6/25
由于r2是一非零常数，故在第2步结束，方程为
n=mq1+r1
m=r1q2+r2
解得
r2=(1+q1q2)m－q2n
因此取
x=-q2/r2, y=(1+q1q2)/r2
即
x(λ)=－30λ+19, y(λ)=5λ+1
下面是做G的互质分解的步骤，主要思路是做变量代换，即s→λ，以使λ的多项式对应于稳定的、正则的、实有理的函数。

步骤B
输入：G
第1步如果G是稳定的，取N=G，M=1，X=0，Y=1，停止。

否则继续。

第2步利用映射s=(1－λ)/ λ将G(s)变换成Ĝ(λ)。

将Ĝ写成互质多项式的比：
G(λ)=n(λ)/m(λ)
第3步利用Euclid算法，求得多项式x(λ)，y(λ)使得
nx+my=1
第4步利用映射λ=1/(s+1)，将n(λ)，m(λ)，x(λ)，y(λ)变换成N(s)，M(s)，X(s)，Y(s)。

此步骤中利用的映射并不唯一。

但必须保证n，m，x，y映射到N，M，X，Y后，它们属于稳定的、正则的、实有理的函数。

2 利用Matlab实现互质分解
利用[2]提供的多项式运算工具箱，可进行加、减、乘、除等多项式运算。

通过在matlab 中编程实现了上述步骤A中的多项式互质分解的Euclid算法。

其中polyadd为多项式加；polysub为多项式减；polymul为多项式乘；polydiv为多项式除。

对于上面的例题用如下程序可解
n=[-1 1 0];
m=[18 -29 10];
[q1,r1]=deconv(n,m);
r1=polysub(n,polymul(q1,m));
[q2,r2]=deconv(m,r1);
r2=polysub(m,polymul(q2,r1));
[x,r3]=deconv(q2,r2);
[y,r4]=deconv((polyadd(1,polymul(q1,q2))),r2);
-x
y
现在再来看一个k=3的例子。

设
p(s)=s/[(s-1)(10s+1)]
代入映射s=(1－λ)/ λ后
n(λ)=－λ2+λ，m(λ)=18λ2－29λ+10
利用如下程序
n=[-1 1 0];
m=[18 -29 10];
q1=polydiv(n,m)
r1=polysub(n,polymul(q1,m))
q2=polydiv(m,r1)
r2=polysub(m,polymul(q2,r1))
x=-polydiv(q2,r2)
y=polydiv((polyadd(1,polymul(q1,q2))),r2)
解得
x(λ)=－19.8λ+13.9，y(λ)=－1.1λ+0.1
再代入映射λ=1/(s+1)，得
X (s)=(13.9-5.9)/(s+1), Y(s)=(0.1s-1)/(s+1)
再通过如下程序得到可镇定系统的控制器C。

X=tf([13.9 -5.9],[1 1]);
Y=tf([0.1 -1],[1 1]);
C=minreal(X/Y)
C(s)=(139s-59)/(s-10)
第3讲 MATLAB的符号计算
2008-12-03 08:57
MATLAB是一种强大的数学计算软件，它能进行数值计算和符号计算。

符号计算是对还没有赋值的符号，如常数、变量、表达式等，进行运算处理。

它的优点是：运算准确，能得到封闭解。

符号运算可以实现微积分运算、表达式的化简、积分变换、求解代数方程和微分方程等多种功能。

3.1 符号表达式的建立
数值计算中的变量必须被赋值，符号计算中的变量必须被说明。

（1）说明符号常量的方法：
定义符号常量的命令是：sym('w1')，例如：
sym('tan(pi/3)')
把符号转成对应的数字的命令是：sym(w1)，例如：
sym(tan(pi/3))，sym(tan(pi/3),'d')
（2）创建符号变量和表达式的方法
定义符号变量和符号表达式的命令是：sym('x')和syms x y z 例如：
x=sym('x'),%定义单个符号是符号变量
y=sin(x)*cos(x)+cos(x)*sin(x),%计算符号表达式
有多个符号变量的时候：
syms a b c t;%定义多个符号变量
r1=solve('a*t^2+b*t+c=0'),%解单引号中的一元二次方程的根还可以用下面的方法求解方程的根：
syms x a b c;
y=a*x^2+b*x+c;
r=solve(y)
（3）创建符号矩阵的方法：
用命令sym或syms可以创建符号矩阵，
A=sym('[a,b;c,d]'),
syms x y z;
B=[x+3*x,5*z+6*z;y-y,z/z],
3.2 符号表达式的代数运算
例题1. 已知二元一次方程组
请求出它的解。

解：（1）先写出它的数学模型，用矩阵的形式表示为
或者
（2）然后用MATLAB语言写出它的程序表达式，
A=sym('[a11,a12;a21,a22]')，
B=sym('[b1;b2]')，
X=A\B
（3）程序的运行结果是
X =
[ (a22*b1-b2*a12)/(-a21*a12+a11*a22)]
[-(a21*b1-a11*b2)/(-a21*a12+a11*a22)]
3.3 符号表达式的操作和表达式的转换
表达式的操作是指确定表达式中的自变量，展开、因式分解或化简表达式。

例题2. 请用MATLAB将表达式
化简成通俗的形式。

解程序如下：
syms a b t
y=sin((a+b)*t)+sin((a-b)*t)
pretty(y)
expand(y)
例题3. 请用MATLAB因式分解表达式
解程序如下：
syms x;
f=x^3-6*x^2+11*x-6,
y=factor(f)
例题4. 拟合是指用一个多项式来表示一组数据的变化规律，其准则是什么？有一组温室温度的变化数据，当时间t=0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10
小时时，温室的温度f=14, 15, 16, 15, 15, 17, 19, 20, 22, 24, 25℃。

请用一个合适的多项式来表示它的变化规律。

解：先写出MATLAB的程序
t=[0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10];%设置时间变量
f=[14, 15, 16, 15, 15, 17, 19, 20, 22, 24, 25];%设置温度变量
plot(t,f);grid;hold on;%画出温度随时间变化的曲线并保持它
n=1;%设置用来模拟温度变化的多项式的阶数
order={n,'阶'}%用元胞数组写出多项式的阶数
p=polyfit(t,f,n)%计算模拟温度变化的多项式的系数矩阵
y=polyval(p,t);%计算这个多项式的函数值
plot(t,y,'r');%画出模拟多项式的红色曲线
y=poly2sym(p)%将多项式的按降幂排列的行向量转换成符号表达式
pretty(y)%把多项是写成常见的格式
3.4 符号极限、微积分和级数求和
（1）求极限时，用函数limit。

采用help limit来查看符号极限的基本函数和使用方法，
syms x a t h;
limit(sin(x)/x) returns 1
limit((x-2)/(x^2-4),2) returns 1/4
limit((1+2*t/x)^(3*x),x,inf) returns exp(6*t)
limit(1/x,x,0,'right') returns inf
limit(1/x,x,0,'left') returns -inf
limit((sin(x+h)-sin(x))/h,h,0) returns cos(x)
v = [(1 + a/x)^x, exp(-x)];
limit(v,x,inf,'left') returns [exp(a), 0]（2）求微分时，用函数diff。

可以参考书上的例题3.15：
f=sym('a*x^2+b*x+c')
diff(f)
diff(f,'a')
diff(f,'x',2)
diff(f,3)
（3）求积分时，用函数int。

可以参考书上的例题3.16续：
syms t x;
int(g)
g=[2*x t^2;t*sin(x) exp(x)]
int(g,'a','b')
3.5 符号方程的求解
代数方程的求解用函数solve。

例题5. 已知二元一次方程组
请求出它的解。

解：该题和例题1.是相同的内容，其MATLAB程序如下，
eq1=sym('a11*x1+a12*x2=b1');
eq2=sym('a21*x1+a22*x2=b2');
[x1,x2]=solve(eq1,eq2)
结果是一样的。

3.6 符号函数的可视化
1. 绘制符号函数的曲线的命令
在一般画图命令的前面加上前缀“ez”即可，“ez”表示容易使用“easy”的含义。

例题5. 已知一元二次方程
请画出它在[-1，1]的曲线。

解：该题的MATLAB程序是
ezplot(sym('2*x^2'),[-1,1])
书上表3.3的绘图命令演示：
显示的结果很好看。

2. 图形化的符号函数计算器
在MATLAB的命令窗中输入命令“funtool“就会出现窗口形式的符号计算器，在这个窗口上，可以直接对符号函数进行计算和画图。

一、在 MATLAB 中的窗函数,十分简单：
（1）矩形窗（Rectangle Window）调用格式：w=boxcar(n)，根据长度 n 产生一个矩形窗 w。

（2）三角窗（Triangular Window）调用格式：w=triang(n)，根据长度 n 产生一个三角窗 w。

（3）汉宁窗（Hanning Window）调用格式：w=hanning(n)，根据长度 n 产生一个汉宁窗 w。

（4）海明窗（Hamming Window）调用格式：w=hamming(n)，根据长度n 产生一个海明窗 w。

（5）布拉克曼窗（Blackman Window）调用格式：w=blackman(n)，根据长度 n 产生一个布拉克曼窗 w。

（6）恺撒窗（Kaiser Window）调用格式：w=kaiser(n,beta)，根据长度 n 和影响窗函数旁瓣的β参数产生一个恺撒窗w。

二、基于窗函数的 FIR 滤波器设计
利用 MATLAB 提供的函数 fir1来实现
调用格式：fir1(n,Wn,’ftype’,Window)，n 为阶数、Wn 是截止频率（如果输入是形如[W1 W2]的矢量时，本函数将设计带通滤波器，其通带为W1<ω<>
[例]设计一个长度为 8 的线性相位 FIR 滤波器。

其理想幅频特性满足
用矩形窗:
Window=boxcar(8);
b=fir1(7,0.4,Window);
freqz(b,1)
用blackman窗:
Window=blackman(8);
b=fir1(7,0.4,Window);
freqz(b,1)
[例] 设计线性相位带通滤波器，其长度N=15，上下边带截止频率分别为W1= 0.3π，w2=0.5π
Window=blackman(16);
b=fir1(15,[0.3 0.5],Window);
freqz(b,1)。