matlab33.docx
- 文档编号:29371219
- 上传时间:2023-07-22
- 格式:DOCX
- 页数:18
- 大小:46.66KB
matlab33.docx
《matlab33.docx》由会员分享,可在线阅读,更多相关《matlab33.docx(18页珍藏版)》请在冰豆网上搜索。
matlab33
MATLAB具有以下几个特点:
易学、适用范围广、功能强、开放性强、网络资源丰富。
启动点击MATLAB图标,进入到MATLAB命令窗(MatlabCommandWindow)。
学会使用help命令。
学会使用demo命令。
说明其功能强大。
演示census;spinner;truss;
pend.m
plot([-0.2,0.2],[0;0],'color','y','linestyle','-','linewidth',10);
g=0.98;l=1;
theta0=pi/6;x0=l*sin(theta0);
y0=-l*cos(theta0);
axis([-0.75,0.75,-1.25,0]);
axis('off');
head=line(x0,y0,'color','r','linestyle','.','erasemode','xor','markersize',40);
body=line([0;x0],[0,y0],'color','b','linestyle','-','erasemode','xor');
t=0;
dt=0.01;
whilet<=50
t=t+dt;
theta=theta0*cos(sqrt(g/l)*t);
x=l*sin(theta);y=-l*cos(theta);
set(head,'xdata',x,'ydata',y);
set(body,'xdata',[0;x],'ydata',[0;y]);
drawnow;
end
退出
在工具栏中点击File按钮,在下拉式菜单中单击ExitMATLAB项即可。
或者,在指令窗内键入exit或quit亦可。
矩阵运算的操作(demo)
MATLAB的符号运算功能
求和
symsum(S)对通项S求和,其中k为变量。
且从0变到k-1。
symsum(S,v)对通项S求和,指定其中v为变量。
且v从0变到v-1。
symsum(S,a,b)对通项S求和,其中k为变量。
且从a变到b。
symsum(S,v,a,b)对通项S求和,指定其中v为变量。
且v从a变到b。
例:
键入k=sym('k');symsum(k)得
ans=
1/2*k^2-1/2*k
又例如:
键入symsum(k^2,0,10)得
ans=
385
又例如:
键入symsum('x'^k/sym('k!
'),k,0,inf)得
ans=
exp(x)
这最后的一个例子是无穷项求和。
ⅱ求导数
diff(S,v)求表达式S对变量v的一阶导数。
diff(S,v,n)求表达式S对变量v的n阶导数。
例如:
键入命令
A=sym('[1/(1+a),(b+x)/cos(x);1,exp(x^2)]');
diff(A,'x')得
ans=
[0,1/cos(x)+(b+x)/cos(x)^2*sin(x)]
[0,2*x*exp(x^2)]
又如求sin(x)+ex的三阶导数,键入命令
diff('sin(x)+x*exp(x)',3)得
ans=
-cos(x)+3*exp(x)+x*exp(x)
再如:
求
A=
[x*sin(y),x^n+y]
[1/x/y,exp(i*x*y)]
的先对x再对y的混合偏导数。
可键入命令:
S=sym('[x*sin(y),x^n+y;1/x/y,exp(i*x*y)]');
dsdxdy=diff(diff(S,'x'),'y')得:
dsdxdy=
[cos(y),0]
[1/x^2/y^2,i*exp(i*x*y)-y*x*exp(i*x*y)]
求y=(lnx)x的导数
可键入命令:
p='(log(x))^x';
p1=diff(p,'x')
得
p1=
log(x)^x*(log(log(x))+1/log(x))
求y=xf(x2)的导数
可键入命令:
p='x*f(x^2)';
p1=diff(p,'x')
得
p1=
f(x^2)+2*x^2*D(f)(x^2)
求xy=ex+y的导数
可键入命令:
p='x*y(x)-exp(x+y(x))';
p1=diff(p,'x')
p1=
y(x)+x*diff(y(x),x)-(1+diff(y(x),x))*exp(x+y(x))
p2='y+x*dy-(1+dy)*exp(x+y)=0';
dy=solve(p2,'dy')%把dy作为变量解方程
得
dy=
-(y-exp(x+y))/(x-exp(x+y))
ⅲ求极限
limit(P)表达式P中自变量趋于零时的极限。
limit(P,a)表达式P中自变量趋于a时的极限。
limit(P,x,a,'left')表达式P中自变量x趋于a时的左极限。
limit(P,x,a,'right')表达式P中自变量x趋于a时的右极限。
例如:
键入
P=sym('sin(x)/x');
limit(P)得
ans=
1
键入
P=sym('1/x');
limit(P,x,0,'right')得
ans=
inf
键入
P=sym('(sin(x+h)-sin(x))/h');h=sym('h');
limit(P,h,0)得
ans=
cos(x)
键入
v=sym('[(1+a/x)^x,exp(-x)]');
limit(v,x,inf,'left')得
ans=
[exp(a),0]
ⅳ求泰勒展开式
taylor(f,v)f对v的五阶Maclaurin展开。
taylor(f,v,n)f对v的n-1阶Maclaurin展开。
例如求sin(x)e-x的7阶Maclaurin展开。
可键入
f=sym('sin(x)*exp(-x)');F=taylor(f,8)得
F=
x-x^2+1/3*x^3-1/30*x^5+1/90*x^6-1/630*x^7
如果要求sin(x)e-x在x=1处的7阶Taylor展开。
可键入
f=sym('sin(x)*exp(-x)');F=taylor(f,8,1)得
F=
sin
(1)*exp(-1)+(-sin
(1)*exp(-1)+cos
(1)*exp(-1))*(x-1)
-cos
(1)*exp(-1)*(x-1)^2
+(1/3*sin
(1)*exp(-1)+1/3*cos
(1)*exp(-1))*(x-1)^3
-1/6*sin
(1)*exp(-1)*(x-1)^4
+(1/30*sin
(1)*exp(-1)-1/30*cos
(1)*exp(-1))*(x-1)^5
+1/90*cos
(1)*exp(-1)*(x-1)^6
+(-1/630*cos
(1)*exp(-1)-1/630*sin
(1)*exp(-1))*(x-1)^7
多元函数的taylor展开
MATLAB不能直接进行多元函数的taylor展开。
必须先调用MAPLE函数库中的mtaylor命令。
方法为:
在MATLAB的工作窗口中键入
maple('readlib(mtaylor)')
mtaylor的格式为
mtaylor(f,v,n)
f为欲展开的函数式。
v为变量名。
写成向量的形式:
[var1=p1,var2=p2,…,varn=pn],展开式将在(p1,p2,…,pn)处进行。
如只有变量名,将在0点处展开。
n为展开式的阶数(n-1阶)。
要完成taylor展开,只需键入maple('mtaylor(f,v,n)')即可。
例:
在(x0,y0,z0)处将F=sin(x,y,z)进行2阶taylor展开。
键入
symsx0y0z0
maple('readlib(mtaylor)');
maple('mtaylor(sin(x*y*z),[x=x0,y=y0,z=z0],2)')得:
ans=
sin(x0*y0*z0)+cos(x0*y0*z0)*y0*z0*(x-x0)+cos(x0*y0*z0)*x0*z0*(y-y0)+cos(x0*y0*z0)*x0*y0*(z-z0)
ⅴ求积分
int(P)对表达式P进行不定积分。
int(P,v)以v为积分变量对P进行不定积分。
int(P,v,a,b)以v为积分变量,以a为下限,b为上限对P进行定积分。
例如可键入int('-2*x/(1+x^2)^2')得
ans=
1/(1+x^2)
键入int('x/(1+z^2)','z')得
ans=
atan(z)*x
键入int('x*log(1+x)',0,1)得
ans=
1/4
定积分的上下限可以是(符号)函数。
例如可键入:
int('2*x','sin(t)','log(t)')
得
ans=
log(t)^2-sin(t)^2
对(符号)矩阵进行积分,例
输入int('[exp(t),exp(a*t)]'),得:
ans=
[exp(t),1/a*exp(a*t)]
⑶求符号方程的解
ⅰ线性方程组的求解
线性方程组的形式为A*X=B;其中A至少行满秩。
X=linsolve(A,B)输出方程的特解X。
例如:
键入
A=sym('[cos(t),sin(t);sin(t),cos(t)]');
B=sym('[1;1]');
c=linsolve(A,B)
c=
[1/(sin(t)+cos(t))]
[1/(sin(t)+cos(t))]
例如:
键入
a=sym('[2,7,3,1;3,5,2,2;9,4,1,7]');b=sym('[6;4;2]');
X=linsolve(a,b)
Warning:
Systemisrankdeficient.Solutionisnotunique.
X=
[0]
[0]
[2]
[0]
ⅱ代数方程的求解
solve(P,v)对方程P中的指定变量v求解。
v可省略。
solve(p1,P2,…,Pn,v1,v2,…,vn)对方程P1,P2,…Pn中的指定变量v1,
v2…vn求解。
例:
可输入
solve('p+sin(x)=r')得:
ans=
-asin(p-r)
又例:
可输入:
P1='x^2+x*y+y=3';P2='x^2-4*x+3=0';
[x,y]=solve(P1,P2)得:
x=
[1]
[3]
y=
[1]
[-3/2]
可输入:
P1='a+u^2+v^2=0';P2='u-v=1';
[u,v]=solve(P1,P2,'u','v')得:
u=
[1/2+1/2*(-1-2*a)^(1/2)]
[1/2-1/2*(-1-2*a)^(1/2)]
v=
[-1/2+1/2*(-1-2*a)^(1/2)]
[-1/2-1/2*(-1-2*a)^(1/2)]
对于有些无法求出解析解的非线性方程组,MATLAB只给出一个数值解。
这一点可以
从表示解的数字不被方括号括住而确定。
例如:
键入:
[x,y]=solve('sin(x+y)-exp(x)*y=0','x^2-y=2')得:
x=
-6.0173272500593065641097297117905
y=
34.208227234306296508646214438330
由于这两个数字没有被[]括住,所以它们是数值解。
另外,可利用solve来解线性方程组的通解。
例如:
键入
P1='2*x1+7*x2+3*x3+x4=6';
P2='3*x1+5*x2+2*x3+2*x4=4';
P3='9*x1+4*x2+x3+7*x4=2';
u=solve(P1,P2,P3,'x1','x2','x3','x4')
Warning:
3equationsin4variables.
u=
x1:
[1x1sym]
x2:
[1x1sym]
x3:
[1x1sym]
x4:
[1x1sym]
可以看到:
屏幕提示“有3个方程4个变量”,意为解不唯一。
(有时会提示解不唯一)且输出的是解的结构形式。
为进一步得到解,可输入:
u.x1,u.x2,u.x3,u.x4,得:
ans=
x1
ans=
-5*x1-4*x4
ans=
11*x1+9*x4+2
ans=
x4
这样就得到了原方程组的通解。
⑷解符号微分方程
dsolve('eq1','eq2',…)其中eq表示相互独立的常微分方程、初始条件或
指定的自变量。
默认的自变量为t。
如果输入的初
始条件少于方程的个数,则在输出结果中出现常数
c1,c2,等字符。
关于微分方程的表达式有如下的约
定:
字母y表式函数,Dy表示y对t的一阶导数;
Dny表示y对t的n阶导数。
例如:
求
的解。
可键入:
[x,y]=dsolve('Dx=y','Dy=-x')
得
x=
cos(t)*C1+sin(t)*C2
y=
-sin(t)*C1+cos(t)*C2
dsolve中的输入宗量最多只能有12个,但这并不妨碍解具有多个方程的方程组,因为
可以把多个方程或初始条件定义为一个符号变量进行输入。
例如求
f(0)=0,g(0)=1的解。
可输入指令:
P='Df=3*f+4*g,Dg=-4*f+3*g';
v='f(0)=0,g(0)=1';
[f,g]=dsolve(P,v)
f=
exp(3*t)*sin(4*t)
g=
exp(3*t)*cos(4*t)
注意:
微分方程表达式中字母D必须大写。
求解微分方程
可输入
y=dsolve('D3y=-y','y(0)=1,Dy(0)=0,D2y(0)=0','x')得:
y=
(1/3+2/3*exp(1/2*x)*cos(1/2*3^(1/2)*x)*exp(x))/exp(x)
最后看一个解非线性微分方程的例子:
dsolve('(Dy)^2+y^2=1','y(0)=0','x')
ans=
[sin(x)]
[-sin(x)]
对于无法求出解析解的非线性微分方程,屏幕将提示出错信息。
与数模有关的例
1.曲线拟合美国人口预测
1.下表是美国人口统计数据,根据这份资料预测2000年美国人口总数。
年
1790180018101820183018401850
人口(百万)
3.95.37.29.612.917.123.2
年
1860187018801890190019101920
人口(百万)
31.438.650.262.976.092.0106.5
年
1930194019501960197019801990
人口(百万)
123.2131.7150.7179.3204.0226.5251.4
Step1A=[1790,1800,1810,…;
3.9,5.3,7.2,…]’;
Step2P=polyfit(A(:
1),A(:
2),3)
Step3px=poly2str(P,'x')
Step4polyval(P,2000)
如果想了解fx与数据对x-y的拟和程度,绘出二者的图形最为直观,为此可键入:
ft=polyval(P,A(:
1));plot(A(:
1),A(:
2),'bo',A(:
1),ft,'r-')得图形。
图中蓝色小圆圈是数据对的图形;而红线是拟合多项式的图形。
最后,可与demo_sensus比较。
2.插值
“线性插值”linear
“三次样条插值”spline
“三次多项式插值”cubic
对于以上问题,也可以用这三个命令来做。
3.交通流量问题
下图给出了某城市部分单行街道的交通流量(每小时过车数)
x3100x6
300x4400200
x2x7
300x1600x8
500200400
300500
x9x10
600700
所给问题满足下列方程组
x1-x3+x4=300
x4+x5=500
x7-x6=200
x1+x2=800
x1+x5=800
x7+x8=1000
x8+x3+x6=1000(x9=400,x10=600)
Step1A=[0,1,-1,1,0,0,0,0;
0,0,0,1,1,0,0,0;
0,0,0,0,0,-1,1,0;
1,1,0,0,0,0,0,0;
1,0,0,0,1,0,0,0;
0,0,0,0,0,0,1,1;
0,0,1,0,0,1,0,1];
Step2b=[300,500,200,800,800,1000,1000]’;
Step3B1=rank(A);B2=rank([A,b]);
Step4X=linsolve(A,b)
得特解。
4.线性规划有约束极小问题
用命令x=lp(C,A,b,vlb,vub)。
例:
Findxthatminimizes
f(x)=-5x1-4x2-6x3
subjectto
x1-x2+x3≦20
3x1+2x2+4x3≦42
3x1+2x2≦30
0≦x1,0≦x2,0≦x3
First,enterthecoefficients:
f=[-5;-4;-6]
A=[1-11
324
320];
b=[20;42;30];
lb=zeros(3,1);
Next,callalinearprogrammingroutine:
x=lp(f,A,b,lb);
Enteringx
x=
0.0000
15.0000
3.0000
实际此命令改为:
x=linprog(f,A,b,Aeq,beq)
x=linprog(f,A,b,Aeq,beq,lb,ub)
对以上的问题可做如下的操作:
First,enterthecoefficients:
f=[-5;-4;-6]
A=[1-11
324
320];
b=[20;42;30];
lb=zeros(3,1);
Next,callalinearprogrammingroutine:
[x,fval,exitflag,output,lambda]=linprog(f,A,b,[],[],lb);
Enteringx,fval,lambda.ineqlin,andlambda.lowergets
x=
0.0000
15.0000
3.0000
fval=
-78.0000
和其它信息。
5.非线性规划有约束极小问题
用命令x=constr('f',x0)。
Examples
Findvaluesofxthatminimizef(x)=-x1x2x3,startingatthepointx=[10;10;10]andsubjecttotheconstraints0≤x1+2x2+2x3≤72.
-x1-2x2-2x3≤0,x1+2x2+2x3≤72,
第一步:
编写M文件
function[f,g]=myfun(x)
f=-x
(1)*x
(2)*x(3);
g
(1)=-x
(1)-2*x
(2)-2*x(3);
g
(2)=x
(1)+2*x
(2)+2*x(3)-72;
第二步:
求解
在MATLAB工作窗中键入
x0=[10,10,10];
x=constr('myfun',x0)即可
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- matlab33