E:WHU教学2010Lecture2008Lec4Matlab基础

Download Report

Transcript E:WHU教学2010Lecture2008Lec4Matlab基础

Matlab 基础
IEEE standard
1.1*realmax;
Double precision (64-bit word)
0/0; inf/inf;
16 significant decimal digits
Realmin*eps;
| x  fl ( x) | 1
  M ,  M   1t .
| x|
2
Eps 2-52 , 2.22 e-016
Unit roundoff 1.11 e-016
Realmax;
Realmin;
Format long e/g;
Realmin*eps*0.1;
数学函数
Help elfun
Help specfun
a=1;
b=1;
while a+b~=a
b=b/2;
end
a=1.0e+308;
b=1.1e+308;
c=-1.001e+308;
a+(b+c)
a+b+c
%numerical cancellation
x=1.e-15;
((1+x)-1)/x
矩阵的基本操作
矩阵的创建
A=[1 2 3; 4 5 6];
%行内元素用_或,
%行间用;
Zeros()
Ones()
Eye()
A=rand(5);
A=magic(3);
用小矩阵作元素建立大矩阵
B=[A;[1 2 3]];
C=[A eye(size(A))];
从已存在的矩阵创建
Diag()
Tril()
Triu()
Rot90
Fliplr
Reshape
Repmat
Shiftdim
Squeeze
Cat
permute
元素的提取,赋值
向量的操作
下标
X,v是向量,x(v)也是向量
A(:,[2 3 5])=B(:,1:3);
A(:,n:-1:1);
按列拉直
B=A(:)
A(:)=11:19 与原来阶数相同
Zeros()
Ones()
X=(4:-1:1)’;
Y=exp(-x).*sin(x);
X=linspace(0,2*pi,100)
[x y]
Size()
元素与元素运算( .* 和 .^)
Whos
区别 w’*v 与w.*v
[1 2 3].^[4 5 6];
Special matrices
Basic data analysis functions
Hilb; invhilb
Max
Toeplitz
Min
Compan
Mean
Gallery
Median
Hadamard
Std
Hankel
Var
Magic
Sort
Pascal
Sum
Rosser
Prod
Vander
Cumsum
wilkinson
Cumprod
diff
稀疏矩阵
比较eye(1000)和speye(1000), (8Mb,16Kb)
Spdiags 对角稀疏矩阵
 2 1



 1 2 1


  



 1 2 

D=sparse(1:n, 1:n, -2*ones(1,n), n, n)
从外部导入
L=sparse(2:n, 1:n-1, ones(1,n-1), n, n)
Full()
A=L+D+L’
Nnz
如何生成
load dd.dat
S=spconvert(dd)
%将一个有三列的矩阵转化成一个稀疏矩阵
Sprand() %随机稀疏矩阵
Spy
典型的矩阵运算
基本的线代函数
乘法:
LU, CHOL,
A*B
QR, SVD,
左除和右除
RREF, DET,
A\b;
INV, PINV,
A/B
矩阵卷积 conv2
ORTH, NULL,
RCOND, COND,
condest, NORM,
RANK, TRACE,
张量积 kron
subspace
Dot
Cross
共轭转置 A’
仅转置 A.’
线性代数方程组求解
CG, PCG, CGS, BICG,
BICGSTAB, GMRES, QMR,
SYMMLQ, MINRES
LUINC, CHOLINC
LSQR
特征值
Eig, eigs, schur, qz, svd, hess, poly, balance
矩阵函数
EXPM
创建一个10x10矩阵,使它的第一列和最
后一列全为0;其余元素全为2.
用Matlab命令(可用左除,右除或LU)求
解方程组Ax=b,这里A是一个随机矩阵,
由A=rand(n)生成,取n=100;b是右端
项,b=A*ones(n,1).请用下面的命令给出
所求解的误差: norm(x-ones(n,1));并用
下面的命令给出残量范数: norm(bA*x,’fro’).
非线性方程的求解
fzero(),
fmin( ,0.5, 1)
对多项式 polyval() polyfit()
roots()
非线性方程组
fsolve()
多项式
x=[1 1.5 3 ];
y=[1.2 3.3 ];
p1=polyfit(x,y,1); %线性拟合
y1=polyval(p1,x);
plot(x,y1,x,y,’*’);
function y= ff(var )
y = var^3 - 3*var - 1;
return;
%%
fzero(‘ff’,1);
Newton法解非线性方程组
f=inline('[2*x-y^2+log(x); x^2-x*y-x+1]','x','y');
df=inline('[2+1/x -2*y; 2*x-y-1 -x]','x','y');
x=[1 1];
for i=1:6
b=f(x(1),x(2));
B=df(x(1),x(2));
dx=-B\b;
err=norm(dx,'inf')
x=x+dx‘
end
Matlab编程
关系运算符
控制语句
If expression
Statements
end
==
~=
&
for variable=expression
Statements
end
|
~
XOR
注释 %
续行 …
while expression
Statements
end
switch x
case 1
statement1
case
statement2
otherwise
end
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
1
11
d
 dt y1 (t )  10  y2 (t )  y1 (t ) ,

d
 y2 (t )  28y1 (t )  y2 (t )  y1 (t ) y3 (t ),
 dt
8
d
 dt y3 (t )  y1 (t ) y2 (t )  3 y3 (t ).

g=2;
for k=1:100
g=1+1/g;
end
g
function yprime = lorenzeq(t,y)
yprime = [ 10*( y(2)-y(1) )
28*y(1) -y(2) -y(1)*y(3)
y(1)*y(2) -8/3 *y(3) ];
%ODE solver
t = [0 50];y0= [0;1;0];
[t,y] = ode45('lorenzeq',t,y0);
plot( y(:,1),y(:,3) );
xlabel('y_1');ylabel('y_2');
title( 'Lorenz equations','FontSize',16)
输入输出语句
load
save
fid=fopen('Velocity.dat','w');
fprintf(fid,'varibles="x","y","u","v" \n');
fprintf(fid,['zone i=' num2str(n-1) ' j=' num2str(n-1) ' f=point \n']);
fprintf(fid,'%g %g %g %g\n',[ xs'; ys'; z(1:(n-1)^2)'; z((n-1)^2+1:2*(n-1)^2)' ]);
fclose(fid);
figure(4); surf(X,Y,P);
print -deps Oseen_Surf.eps
图形基础
fplot(‘ ‘, [ ])
ezplot()
plot(x,y,’r---+’)
surf, mesh, waterfall, contour
x=linspace(起点,终点,个数);
y=sin(2*x);
plot(x,y);
plotyy()
semilogy()
z  e10t (t 1) sin( 12t )
xlabel(‘ ‘)
title(‘ ‘)
legend(‘ ‘, ‘ ‘)
text(10,100, ‘ ‘)
gtext(‘ ‘)
axis off
clf, colse, refesh
t=0:0.005:1;
z=exp(10*t.*(t-1)).*sin(12*pi*t);
plot(t,z);
fplot('exp(10*x*(x-1))*sin(12*pi*x)',[0 1])
ezplot('exp(10*x*(x-1))*sin(12*pi*x)‘,[0,1])
% 例子: Mandelbrot Set, Z->Z.^+C
h = waitbar(0,'Computing...');
x=linspace(-2.1, 0.6, 301);
y=linspace(-1.1, 1.1, 301);
[X,Y]=meshgrid(x,y);
C=complex(X,Y);
Z_max = 1e6;
it_max =50;
Z=C;
for k=1:it_max
Z= Z.^2 +C;
waitbar(k/it_max);
End
close(h)
contourf(x,y,abs(Z)<Z_max,1);
title ('Mandelbrot Set','FontSize',16);
符号计算基础
f=‘x^2+2*x+2’;
diag
syms x
tril
diff(f)
det
int(f)
rank
taylor()
rref
null
funtool
svd
expm
syms a b c x
y=solve(a*x^2+b*x+c)
vpa(pi,50)
maple
colspace*
jordan*
A_num=gallery('frank',5);
A_sym=sym(A_num);
inv(A_sym);
eig(A_sym);
poly(A_sym);
其他
More about eigenvalue computations
Polynomials and data fitting
Numerical differentiation and integration
ODE/PDE solvers
Troubleshooting