1、xlsread和xlswrit函數

在MATLAB中經常會用到數據的讀取，首先是從Excel中讀取數據到MATLAB中去。下面給出原始Excel數據內容：

在MATLAB讀取結果如下：

m=xlsread('fanjufei.xls',1,'A1:C3')


m =


     1     2     3
     4     5     6
     7     8     9

其中xlsread可以直接從Excel中讀取文件，'fanjufei.xls'表示讀取文件的名稱，1表示位于sheet1；'A1：C3'表示讀取數據的范圍。

xlswrite函數:

可以從MATLAB中寫入數據到Excel中去，下面給出要寫入數據：

clear
clc
n=[1 2 3 4;5 6 7 8;9 10 11 12];
xlswrite('fanjufei.xls',n,3,'B2:E4')

2、數據擬合

2.1 多項式擬合

例如：有兩組數據為x=[1 2 3 4 5 6 7 8 9];y=[9 7 6 3 -1 2 5 7 20];寫出x與y的表達式；

clear
clc
x=[1 2 3 4 5 6 7 8 9];
y=[9 7 6 3 -1 2 5 7 20];
P=polyfit(x,y,3);
xi=010;
yi=polyval(P,xi);
plot(xi,yi,x,y,'r*');

注釋：polyfit(x,y,N)，x、y為原始數據，N為擬合最高次冪，

polyval(P,xi)，P為各項的系數，結果展示為：

P 0.148-1.403 1.85368.2698

故多項式的結果為：

2.2工具箱擬合

打開工具→基本擬合，選定擬合的方式。

x=[1 2 3 4 5 6 7 8 9];
y=[9 7 6 3 -1 2 5 7 20];
plot(x,y,'r*');

得出結果如下：

2.3自定義擬合函數

例如：要擬合數據：

clear
clc
syms t
x=[0;0.4;1.2;2;2.8;3.6;4.4;5.2;6;7.2;8;9.2;10.4;11.6;12.4;13.6;14.4;15];
y=[1;0.85;0.29;-0.27;-0.53;-0.4;-0.12;0.17;0.28;0.15;-0.03;-0.15;-0.071;0.059;0.08;0.032;-0.015;-0.02];
f=fittype('a*cos(k*t)*exp(w*t)','independent','t','coefficients',{'a','k','w'});
cfun=fit(x,y,f)
xi=020;
yi=cfun(xi);
plot(x,y,'r*',xi,yi,'b-');

結果：

cfun = 


     General model:
     cfun(t) = a*cos(k*t)*exp(w*t)
     Coefficients (with 95% confidence bounds):
       a =      0.9987  (0.9836, 1.014)
       k =      -1.001  (-1.006, -0.9958)
       w =     -0.2066  (-0.2131, -0.2002)

注釋：fittype是自定義擬合函數，cfun=fit(x,y,f)擬合數據x、y，x、y必須為列向量。

故結果為：

3、數值計算

3.1多項式

（1）多項式表示方法

（2）多項式的運算

多項式乘除運算

計算代碼：

clear 
clc


%多項式相乘conv
u=[2 4 5 6];%多項式
v=[10 20 30];%多項式
p=[1 2 3];%多項式
w=conv(u,v)  %conv為多項式相乘函數，也可以嵌套使用；
m=conv(conv(u,p),v)  


%多項式相除deconv
[q,r]=deconv(w,v)%q為商，r為余數；

結果：

w =


    20    80   190   280   270   180


m =


 20         120         410         900        1400        1560        1170         540


q =


     2     4     5     6


r =


     0     0     0     0     0     0

多項式的導函數

k=polyder(p),返回多項式p的導函數；
k=polyder(a,b),返回多項式a乘以b的導函數；

[q,d]=polyder(b,a),返回多項式b整除a的導函數，其分子返回給q，分母為d；

clear
clc
x=[1 2 3 4];
y=[1 2 3 4];
z=polyder(x,y)

z =


     6    20    40    60    50    24

多項式求值

y=polyval(p,x),代數多項式求值，若x為一數值，則求在該點的值；若為向量、矩陣，則求向量、矩陣中的每一個值；
y=polyvalm(p,x),矩陣多項式求值，要求x為方陣；

p=[1,2,3];
x=1:5;
y=polyval(p,x)

y =


     6    11    18    27    38

多項式的根

函數roots：可以求出多項式等于0的根；
函數poly：可以通過多項式等于0的根，求出多項式；

p=[1 2 1];
r=roots(p) %求p的根
v=poly(r)  %求r根的多項式

r =


    -1
    -1


v =


     1     2     1

3.2曲線擬合

曲線擬合用一個比較簡單的函數去逼近一個未知的函數，曲線擬合最優的標準采用最小二乘法原理，擬合的結果使得誤差的平方和最小。

在MATLAB上最常采用polyfit函數來求最小二乘擬合多項式的系數，再用polyval函數求出多項式在所給出點的值；

x=linspace(0,2*pi,50);
y=cos(x);
p=polyfit(x,y,6);
t=linspace(0,2*pi,50);
y1=polyval(p,t);
plot(x,y,t,y1,'r*')

從圖像上可以看出擬合比較好，紅色星號都在曲線上；

3.3 數據插值

（1）一維數據插值

插值函數：yi=interp1(x,y,xi,method)

根據在x，y處的值，計算函數在xi處的值，其中xi的值不能大于x的值；

method插值方法：linear（線性插值）、nearest（最近點插值）、cubic（3次多項式插值）、spline（3次樣條插值）；

例2：下面為1900—1990每隔10年的人口普查數據：

t=19001990;

p=[75 91 105 123 131 150 179 203 226 249]

求在1985年人口數值；

t=19001990;
p=[75 91 105 123 131 150 179 203 226 249];
yi=interp1(t,p,1985)

得出1985年的人口數為：

yi =


  237.5000

估計1900—2000年人口數值

t=19001990;
p=[75 91 105 123 131 150 179 203 226 249];
xi=1900:2000;
yi=interp1(t,p,xi,'spline');
plot(t,p,':o',xi,yi,'-r')

（2）二維數據插值

插值函數：Z1=interp2(X,Y,Z,XI,YI,method)

X,Y為原始采樣點，Z為對應的采樣值，XI,YI表示欲插值的點，method為插值方法與一維插值方法一樣；

例3：為函數peaks插入更多的線條；

[X,Y]=meshgrid(-44);
Z=peaks(X,Y);
[XI,YI]=meshgrid(-44);
ZI=interp2(X,Y,Z,XI,YI);
mesh(X,Y,Z)
hold on
mesh(XI,YI,ZI+20)

3.4 數值微積分

（1）數值微分

在MATLAB中沒有直接求數值導數的函數，只有計算向前差分的函數diff，其調用格式為：

DX=diff(X)：計算向量X的向前差分，DX(i)=X(i+1)-X(i)，i=1......n-1；

DX=diff(X,n)：計算X的n價向前差分；

DX=diff(A,n,dim)：計算矩陣A的n價差分，dim=1（默認值）

x=[3 2 1];
dx=diff(x)

dx =


    -1    -1

（2）數值積分

被積函數是解析式

MATLAB有兩種函數求解定積分，調用格式為：

quad(函數,a,b,tol,trace)

quadl(函數,a,b,tol,trace)

其中，a為下限，b為上限，tol為精度，trace是否展現積分過程；

f=inline('1./(x.^3-2*x-5)');
y=quad(f,0,2)
y1=quadl(f,0,2)

y =


-0.4605
y1 =


   -0.4605

被積函數為表格定義

用trapz(x,y)來進行計算，x為向量，y為x的函數；

x=01;
y=exp(-x.^2);
trapz(x,y)

ans =


    0.7468

二重積分數值求解

MATLAB提供的函數為：

y=dblquad(f,a,b,c,d,tol,trace)，

function f= fan(x,y)
f=x+y;
end

y=dblquad('fan',0,2,0,2)

結果：

y =


     8

3.5 線性方程組求解

（1）直接解法

對于方程Ax=b來說，可以用x=A；即x=inv(A)*b;

A=[1 2 3;3 6 7;2 6 3];
b=[8 30 25]';
x=A

x =


   17.0000
    0.0000
   -3.0000

（2）LU求解、QR求解、Cholesky求解

例8：求例7；

A=[1 2 3;3 6 7;2 6 3];
b=[8 30 25]';
[L,U]=lu(A); %LU分解
x=U(L)


[Q,R]=qr(A); %QR分解
x_val=R(Q)


R=chol(A);   %Cholesky分解
x_val_1=R(R')

3.6 常微分方程數值求解

基于龍格—庫塔法，MATLAB提供的常微分方程求解的函數為：

[t,y]=ode23('fname',tspan,y0)，二價、三價龍格—庫塔法；

[t,y]=ode45('fname',tspan,y0)，四價、五價龍格—庫塔法；

fname是定義f(t,y)的函數文件名，該函數文件必須返回一個列向量；tspan形式為[t0,tf]表示求解區間，y0是初始狀態列向量；t 給出時間向量，y為狀態向量；

function f = fan(t,x)
f=[-2*x(2);x(1)];
end

tf=25;
[t,y]=ode45('fan',[t0,tf],[1,0]);
subplot(121);
plot(t,y(:,2))
subplot(122);
plot(y(:,2),y(:,1))
axis equal

審核編輯：湯梓紅