湖南大学matlab教案02.docx
- 文档编号:6178758
- 上传时间:2023-01-04
- 格式:DOCX
- 页数:49
- 大小:1.32MB
湖南大学matlab教案02.docx
《湖南大学matlab教案02.docx》由会员分享,可在线阅读,更多相关《湖南大学matlab教案02.docx(49页珍藏版)》请在冰豆网上搜索。
湖南大学matlab教案02
第2章符号计算
符号计算是指:
解算数学表达式、方程不是在离散化的数值点上进行,而是凭借一系列恒等式,数学定理,通过推理和演绎,力求获得解析结果。
符号计算的特点:
1)符号计算定义在符号变量的基础上,符号
表达式计算前必须定义符号变量。
2)符号计算是精确计算。
3)符号计算的计算速度较慢。
4)符号计算的运算符和基本数学函数与数值计算中的运算符和基本数学函数几乎完全相同
.1符号对象和符号表达式
.1.1符号对象的创建
(1)sym(‘变量’,参数)
功能:
把变量定义为符号对象。
其中参数用来设置限定符号变量的数学特性,有三种选择:
’positive’表示为“正、实”符号变量,
’real’表示为“实”符号变量,
’unreal’表示为“非实”符号变量。
如果不限定则参数可省略。
(2)syms函数的格式为:
syms(‘arg1’,‘arg2’,…,参数)
symsarg1arg2…参数
功能:
创建多个符号变量。
10一符号数字
【例2.1-1】符号(类)数字与数值(类)数字之间的差异。
a=pi+sqrt(5)
sa=sym('pi+sqrt(5)')
Ca=class(a)
Csa=class(sa)
vpa(sa-a)
a=
5.3777
sa=
pi+sqrt(5)
Ca=
double
Csa=
sym
ans=
.138223758410852e-16
10二符号变量
1、符号表达式允许使用自由变量。
确定自由变量的原则:
1)小写字母i和j不能作为自由变量。
2)符号表达式中如果有多个字符变量,则按照以下顺序选择自由变量:
首先选择x作为自由变量;如果没有x,则选择在字母顺序中最接近x的字符变量;如果与x相同距离,则在x后面的优先。
2、确定自由变量的指令
findsym的格式为:
findsym(EXPR,n)
功能:
确定EXPR中的自由变量。
其中EXPR可以是符号表达式或符号矩阵;n为按顺序得出符号变量的个数,当n省略时,则不按顺序给出EXPR中所有的符号变量。
【例2.1-2】用符号计算研究方程
的解。
(1)
symsuvwz
Eq=u*z^2+v*z+w;
result_1=solve(Eq)%
findsym(Eq,1)
result_1=
-u*z^2-v*z
ans=
w
(2)指定变量为Z
result_2=solve(Eq,z)
result_2=
1/2/u*(-v+(v^2-4*u*w)^(1/2))
1/2/u*(-v-(v^2-4*u*w)^(1/2))
【例2.1-3】对独立自由符号变量的自动辨认。
(1)
symsabxXY
k=sym('3');
z=sym('c*sqrt(delta)+y*sin(theta)');
EXPR=a*z*X+(b*x^2+k)*Y;
(2)
findsym(EXPR)
ans=
X,Y,a,b,c,delta,theta,x,y
(3)
findsym(EXPR,1)
ans=
x
(4)
findsym(EXPR,2),findsym(EXPR,3)
ans=
x,y
ans=
x,y,theta
【例2.1-4】findsym确定自由变量是对整个矩阵进行的。
symsabtuvxy
A=[a+b*x,sin(t)+u;x*exp(-t),log(y)+v]
findsym(A,1)
A=
[a+b*x,sin(t)+u]
[x*exp(-t),log(y)+v]
ans=
x
.1.2符号计算中的算符
与数字计算的符号相同
.1.3符号计算中的函数指令
见表2.1-1,P48
注意与数值计算的函数和指令的异同
.1.4符号对象的识别
【例2.1-5】数据对象及其识别指令的使用。
(1)
clear
a=1;b=2;c=3;d=4;
Mn=[a,b;c,d]
Mc='[a,b;c,d]'
Ms=sym(Mc)
Mn=
12
34
Mc=
[a,b;c,d]
Ms=
[a,b]
[c,d]
(2)
SizeMn=size(Mn)
SizeMc=size(Mc)
SizeMs=size(Ms)
SizeMn=
22
SizeMc=
19
SizeMs=
22
(3)
CMn=class(Mn)
CMc=class(Mc)
CMs=class(Ms)
CMn=
double
CMc=
char
CMs=
sym
(4)isa判断类别,正确返回1
isa(Mn,'double')
isa(Mc,'char')
isa(Ms,'sym')
ans=
1
ans=
1
ans=
1
(5)
whosMnMcMs
NameSizeBytesClassAttributes
Mc1x918char
Mn2x232double
Ms2x2312sym
.2符号数字及表达式的操作
.2.1数值数字与符号数字之间的转换
10一数值数字向符号数字的转换
如:
a=sym(8,'r')或a=sym('8')
10二符号数字向双精度数字转换
Double(a)%a变成double型
.2.2符号数字的任意精度计算
digits(n)
功能:
设定计算精度和改变默认的有效位数函数。
其中n为所期望的有效位数,默认值为32位。
Vpa的格式为:
S=vpa(s,n)
功能:
将s表示为n位有效位数的符号对象。
例:
应用digits和vpa函数设置运算精度。
a=sym('2*sqrt(5)+pi')%创建符号对象
digits%显示默认的有效位数
vpa(a)%用默认的位数计算并显示
vpa(a,20)%按指定的精度计算并显示
digits(15)%改变默认的有效位数
vpa(a)%按digits指定的精度计算并显示
【例2.2-1】digits,vpa指令的使用。
digits%显示当前有效数字位数(十进制)
p0=sym('(1+sqrt(5))/2')
pr=sym((1+sqrt(5))/2)
%
pd=sym((1+sqrt(5))/2,'d')
%
e32r=vpa(abs(p0-pr))
e16=vpa(abs(p0-pd),16)
e32d=vpa(abs(p0-pd))
Digits=32
p0=
(1+sqrt(5))/2
pr=
7286977268806824*2^(-52)
pd=
1.6180339887498949025257388711907
e32r=
.543211520368251e-16
e16=
0.
e32d=
.543211520368251e-16
.2.3符号表达式的基本操作
指令:
simple(EXPR)
【例2.2-2】简化
。
symsx
f=(1/x^3+6/x^2+12/x+8)^(1/3);
g1=simple(f)
g2=simple(g1)
g1=
(2*x+1)/x
g2=
2+1/x
.2.4表达式中的置换操作
10一通用置换指令
subs(S)
功能:
用MATLAB工作空间中的变量替换S符号表达式中的所有变量。
subs(S,NEW)
功能:
用变量NEW替换符号表达式S中的自由变量。
subs(S,OLD,NEW)
功能:
用变量NEW替换符号表达式S中的变量OLD。
【例2.2-4】用简单算例演示subs的置换规则。
(1)
symsax;
f=a*sin(x)+5
f=
a*sin(x)+5
(2)
f1=subs(f,'sin(x)',sym('y'))
class(f1)
f1=
a*y+5
ans=
sym
(3)
f2=subs(f,{a,x},{2,sym('pi/3')})
class(f2)
f2=
3^(1/2)+5
ans=
sym
(4)
f3=subs(f,{a,x},{2,pi/3})class(f3)
f3=
6.7321
ans=
double
(5)数组置换
f4=subs(subs(f,a,2),x,0:
pi/6:
pi)class(f4)
f4=
5.00006.00006.73217.00006.73216.00005.0000
ans=
double
(6)
f5=subs(f,{a,x},{0:
6,0:
pi/6:
pi})class(f5)
f5=
5.00005.50006.73218.00008.46417.50005.0000
ans=
double
.3符号微积分
.3.1极限和导数的符号计算
1、求函数极限的函数:
limit
limit(f,x,a)
功能:
求符号函数f(x)的极限值。
即计算当自变量x趋近于常数a时,f(x)函数的极限值。
limit(f,x,a,'right')
功能:
求符号函数f的极限值。
'right'表示变量x从右边趋近于a。
limit(f,x,a,'left')
功能:
求符号函数f的极限值。
'left'表示变量x从左边趋近于a。
2、求符号表达式的微分的函数:
diff
diff(f)
功能:
求f对自由变量的一阶微分
diff(f,x)
功能:
求f对符号变量x的一阶微分
diff(f,n)
功能:
求f对自由变量的n阶微分
diff(f,x,n)
功能:
求f对符号变量t的n阶微分
3、r=taylor(f,n,v,a)
功能:
把f(v)在v=a处展开成n次幂级数
以上x缺少时由findsys自动确认,n缺少时默认为1
【例2.3-1】试求
。
symsxk
Lim_f=limit((1-1/x)^(k*x),x,inf)
Lim_f=
exp(-k)
【例2.3-2】求
求
,
,
。
symsatx;
f=[a,t^3;t*cos(x),log(x)];
df=diff(f)
dfdt2=diff(f,t,2)
dfdxdt=diff(diff(f,x),t)
df=
[0,0]
[-t*sin(x),1/x]
dfdt2=
[0,6*t]
[0,0]
dfdxdt=
[0,0]
[-sin(x),0]
【例2.3-4】
,求
,
。
(1)
f=sin(abs(x));
dfdx=diff(f,x)%
dfdx=
cos(abs(x))*abs(1,x)
(2)图形观察
xn=-3/2*pi:
pi/50:
0;
xp=0:
pi/50:
3/2*pi;
xnp=[xn,xp(2:
end)];
holdon
plot(xnp,subs(f,x,xnp),'k','LineWidth',2.5)
plot(xn,subs(df_n,x,xn),'--r','LineWidth',2.5)
plot(xp,subs(df_p,x,xp),':
r','LineWidth',2.5)
legend(char(f),char(df_n),char(df_p),'Location','NorthEast')%<16>
gridon
xlabel('x')
holdoff
图2.3-1函数及其导函数
【例2.3-5】设
,求
,隐函数求导。
(1)
clear
symsx
g=sym('cos(x+sin(y(x)))=sin(y(x))')dgdx=diff(g,x)
g=
cos(x+sin(y(x)))=sin(y(x))
dgdx=
-sin(x+sin(y(x)))*(1+cos(y(x))*diff(y(x),x))=cos(y(x))*diff(y(x),x)
(2)利用maple的isolate指令
disp(maple('isolate',dgdx,diff('y(x)',x)))diff(y(x),x)=sin(x+sin(y(x)))/(-sin(x+sin(y(x)))*cos(y(x))-cos(y(x)))
【例2.3-6】求
在
处展开的8阶Maclaurin级数。
symsx
r=taylor(x*exp(x),9,x,0)
%忽略9阶及9阶以上小量的展开
r=
x+x^2+1/2*x^3+1/6*x^4+1/24*x^5+1/120*x^6+1/720*x^7+1/5040*x^8
【例2.3-7】求
在
处的截断8阶小量的Taylor展开近似。
(二元函数展开)
TL1=maple('mtaylor(sin(x^2+y),[x=0,y=0],8)')
class(TL1)
TL1=
y+x^2-1/6*y^3-1/2*x^2*y^2+1/120*y^5-1/2*x^4*y+1/24*x^2*y^4-1/6*x^6-1/5040*y^7+1/12*x^4*y^3
ans=
char
.3.2序列/级数的符号求和
级数求和运算函数:
symsum
symsum(f,v,a,b)
功能:
计算符号表达式f的级数和。
其中f为通项式,v自变量,v省略则默认为对自由变量求和;[a,b]为参数x的取值范围,缺少时求和区间默认为[0~v-1]。
【例2.3-8】求
,
。
symskt;
f1=[tk^3];
f2=[1/(2*k-1)^2,(-1)^k/k];
s1=simple(symsum(f1))
s2=simple(symsum(f2,1,inf))
s1=
[1/2*t*(t-1),k^3*t]
s2=
[1/8*pi^2,-log
(2)]
.3.3符号积分
符号积分函数:
int
•int(f,v,a,b)
功能:
以v为自变量,以a、b分别表示定积分的下限和上限,对被积函数的符号表达式f求定积分。
•int(f,v)
功能:
以v为自变量,对被积函数的符号表达式f求不定积分。
注意:
没有指定积分变量v时,按findsym函数确定的默认变量对被积函数的符号表达式f求积分。
int函数的嵌套使用可实现二重积分的计算。
【例2.3-9】求
。
clear
symsx
f=sqrt((1+x)/x)/x
s=int(f,x)
s=simple(simple(s))
f=
((1+x)/x)^(1/2)/x
s=
((1+x)/x)^(1/2)/x*(-2*(x^2+x)^(3/2)+2*(x^2+x)^(1/2)*x^2+log(1/2+x+(x^2+x)^(1/2))*x^2)/((1+x)*x)^(1/2)
s=
log(1/2+x+((1+x)*x)^(1/2))-2*((1+x)*x)^(1/2)/x
【例2.3-10】求
。
symsabx;
f=[a*x,b*x^2;1/x,sin(x)];
disp('Theintegraloffis');
pretty(int(f))
%一行显示变两行显示
Theintegraloffis
[23]
[1/2ax1/3bx]
[]
[log(x)-cos(x)]
【例2.3-11】求积分
。
symsxyz
F2=int(int(int(x^2+y^2+z^2,z,sqrt(x*y),x^2*y),y,sqrt(x),x^2),x,1,2)
VF2=vpa(F2)
F2=
-6072064/348075*2^(1/2)+64/225*2^(3/4)+14912/4641*2^(1/4)+1610027357/6563700
VF2=
224.92153573331143159790710032805
【例2.3-12】求阿基米德(Archimedes)螺线
在
到
间的曲线长度函数,并求出
时的曲线长度。
(1)
symsarthetaphipositive
x=r*cos(theta);x=subs(x,r,a*theta);
y=r*sin(theta);y=subs(y,r,a*theta);
dLdth=sqrt(diff(x,theta)^2+diff(y,theta)^2);
L=simple(int(dLdth,theta,0,phi))
L=
1/2*a*(phi*(1+phi^2)^(1/2)+asinh(phi))
(2)
L_2pi=subs(L,[a,phi],sym('[1,2*pi]'))
L_2pi_vpa=vpa(L_2pi)
L_2pi=
pi*(1+4*pi^2)^(1/2)+1/2*asinh(2*pi)
L_2pi_vpa=
21.256294148209098800702512272566
(3)
L1=subs(L,a,sym('1'));
ezplot(L1*cos(phi),L1*sin(phi),[0,2*pi])
gridon
holdon
x1=subs(x,a,sym('1'));
y1=subs(y,a,sym('1'));
h1=ezplot(x1,y1,[0,2*pi]);
set(h1,'Color','r','LineWidth',5)
title('')
legend('螺线长度-幅角曲线','阿基米德螺线')
holdoff
图2.3-2阿基米德螺线(粗红)和螺线长度函数(细蓝)
.4微分方程的符号解法
.4.1符号解法和数值解法的互补作用
符号法适合求初值和边值问题
.4.2求微分方程符号解的一般指令
符号微分方程的求解函数:
dsolve
dsolve('eqn1','eqn2',...,'eqn',’condition1’,‘condition2’,…,‘conditionN’,'v)
功能:
求微分方程组eqn1,eqn2,...,eqn的通解。
初值条件为condition1,condition2,…,conditionN下,变量为v(缺少时默认为t)。
dsolve函数也可求解二阶微分方程。
注意:
DX表示
,DnY表示
边界或初始条件表示如:
y(a)=b,Dy(c)=d
.4.3微分方程符号解示例
【例2.4-1】求
的解。
S=dsolve('Dx=y,Dy=-x')
disp([blanks(12),'x',blanks(21),'y']),disp([S.x,S.y])
S=
x:
[1x1sym]
y:
[1x1sym]
xy
[-C1*cos(t)+C2*sin(t),C1*sin(t)+C2*cos(t)]
【例2.4-3】求解两点边值问题:
。
(1)
y=dsolve('x*D2y-3*Dy=x^2','y
(1)=0,y(5)=0','x')
y=
31/468*x^4-1/3*x^3+125/468
(2)
ezplot(y,[-1,6])
%画y函数的图像,范围[-1~6]
holdon
plot([1,5],[0,0],'.r','MarkerSize',20)
%在(1,0)和(5,0)位置画出红点
text(1,1,'y
(1)=0')
text(4,1,'y(5)=0')
title(['x*D2y-3*Dy=x^2',',y
(1)=0,y(5)=0'])
holdoff
图2.4-2两点边值问题的解曲线
.5符号变换和符号卷积
.5.1Fourier变换及其反变换
•fourier
F=fourier(f,t,w)
功能:
返回函数f(t)的fourier变换F。
其中返回结果F是符号变量w的函数,当参数w省略,默认返回结果为w的函数;当参数t省略,默认自由变量为x。
•傅里叶反变换函数:
ifourier
f=ifourier(F)
f=ifourier(F,w,t)
功能:
返回函数F(w)的fourier反变换f(t)。
参数含义同fourier函数。
【例2.5-1】求
的Fourier变换。
(1)
symstw
ut=heaviside(t);%阶跃函数
UT=fourier(ut)
UT=
pi*dirac(w)-i/w
(2)
Ut=ifourier(UT,w,t)
Ut=
heaviside(t)
【例2.5-2】根据Fourier变换定义,用积分指令求方波脉冲
的Fourier变换。
(1)
symsAtw
symstaopositive
yt=heaviside(t+tao/2)-heaviside(t-tao/2);
Yw=fourier(A*yt,t,w)
Yw=
2*A/w*sin(1/2*tao*w)
(2)
Yt=ifourier(Yw,w,t)
Yt=
A*(heaviside(t+1/2*tao)-heaviside(t-1/2*tao))
(3)设tao=3,A=1
yt3=subs(yt,tao,3)
Yw3=subs(Yw,[A,tao],[1,3])
subplot(2,1,1)
Ht=ezplot(yt3,[-3,3]);
set(Ht,'Color','r','LineWidth',3)
subplot(2,1,2),ezplot(Yw3)
yt3=
heaviside(t+3/2)-heaviside(t-3/2)
Yw3=
2/w*sin(3/2*w)
图2.5-1时域方波及其Fourier变换
【例2.5-3】
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 湖南大学 matlab 教案 02