1、一、MATLAB程序設(shè)計(jì)實(shí)踐Matlab基礎(chǔ)班級(jí)：學(xué)號(hào)：姓名：表示多晶體材料織構(gòu)的三維取向分布函數(shù)（ff（1，2）是一個(gè)非常復(fù)雜的函數(shù)，難以精確的用解析函數(shù)表達(dá)，通常采用離散空間函數(shù)值來(lái)表示取向分布函數(shù)，Data.txt是三維取向分布函數(shù)的一個(gè)實(shí)例。由于數(shù)據(jù)量非常大，不便于分析，需要借助圖形來(lái)分析。請(qǐng)你編寫(xiě)一個(gè)matlab程序畫(huà)出如下的幾種圖形來(lái)分析其取向分布特征：（1）用Slice函數(shù)給出其整體分布特征；（2）用pcolor或contour函數(shù)分別給出(20, 5, 10, 15, 20, 25, 30, 35 90)切面上f分布情況(需要用到subplot函數(shù))；(3) 用plot函數(shù)給出

2、沿取向線(1=090，45，20)的f分布情況。開(kāi)始流程圖讀取Data.txtSlice函數(shù)繪制圖形Pcolor函數(shù)繪制圖形Contour函數(shù)繪制圖形Plot函數(shù)繪制圖形結(jié)束解：（1） 將文件Data.txt內(nèi)的數(shù)據(jù)按照要求讀取到矩陣f(phi1,phi,phi2)中，代碼如下：fid=fopen('data.txt'); %讀取數(shù)據(jù)文件Data.txtfor i=1:18 tline=fgetl(fid);endphi1=1;phi=1;phi2=1;line=0;f=zeros(19,19,19);while feof(fid) tline=fgetl(fid); data

3、=str2num(tline); line=line+1; if mod(line,20)=1 phi2=(data/5)+1; phi=1; else for phi1=1:19 f(phi1,phi,phi2)=data(phi1); end phi=phi+1; endendfclose(fid);將以上代碼保存為readtext.m文件并在MATLAB中運(yùn)行，運(yùn)行結(jié)果如下圖所示：將以下代碼保存為code1_1.m文件：fopen('readtext.m');readtext;x,y,z=meshgrid(0:5:90,0:5:90,0:5:90);slice(x,y,z

4、,f,45,90,45,90,0,45) %運(yùn)用slice函數(shù)繪制圖形運(yùn)行結(jié)果如右圖所示（2） 將以下代碼保存為code1_2_1.m文件：fopen('readtext.m');readtext;for i=1:19 subplot(5,4,i) pcolor(f(:,:,i) %運(yùn)用pcolor函數(shù)繪制圖形end運(yùn)行結(jié)果如右圖所示將以下代碼保存為code1_2_2.m文件：fopen('readtext.m'); %運(yùn)用contour函數(shù)繪制圖形readtext;for i=1:19 subplot(5,4,i) contour(f(:,:,i)end運(yùn)行結(jié)

5、果如右圖所示：（3） 1=090，45，20所對(duì)應(yīng)的f(1,2)即為f(:,10,1)。將以下代碼保存為code1_3.m文件：fopen('readtext.m');readtext;plot(0:5:90,f(:,10,1),'-bo') %運(yùn)用plot函數(shù)繪制圖形text(60,6,'phi=45 phi2=0')運(yùn)行結(jié)果如下圖所示：1. 編程實(shí)現(xiàn)以下科學(xué)計(jì)算算法，并舉一例應(yīng)用之。（參考書(shū)籍精通MATLAB科學(xué)算法，王正林等著，電子工業(yè)出版社，2009年）“多項(xiàng)式擬合”。思考：多項(xiàng)式擬合是用多項(xiàng)式擬合曲線的一種方式，低次數(shù)下運(yùn)用此方法符合較

6、好，但較高次數(shù)下波動(dòng)太大，失去真實(shí)性。1.1 多項(xiàng)式曲線擬合概述對(duì)給定數(shù)據(jù)點(diǎn)（xi,yi）(i=1,2,.N），構(gòu)造m次多項(xiàng)式， P（x）=+ （m<N ）由曲線擬合定義，應(yīng)該使得下式取極小值： 通過(guò)簡(jiǎn)單的計(jì)算可得出系數(shù)是下面的線性方程組的解： =其中, ck= bk=在MATLAB中編程實(shí)現(xiàn)的多項(xiàng)式曲線擬合函數(shù)為:multifit功能：離散實(shí)驗(yàn)數(shù)據(jù)點(diǎn)的多項(xiàng)式曲線擬合。調(diào)用格式：A=multifit(X,Y,m)其中：X為實(shí)驗(yàn)數(shù)據(jù)點(diǎn)的x坐標(biāo)向量； Y為實(shí)驗(yàn)數(shù)據(jù)點(diǎn)的y坐標(biāo)向量； m為擬合多項(xiàng)式的次數(shù)； A為擬合多項(xiàng)式的系數(shù)向量。1.2 多項(xiàng)式曲線擬合編程流程圖輸入向量X，Y，多項(xiàng)式次數(shù)m.

7、M=N得出X的項(xiàng)數(shù)為N，Y的項(xiàng)數(shù)為M建立長(zhǎng)度為2m+1的零向量C和長(zhǎng)度為2m的零向量b，變量j0,k0。j<2m+1YN顯示輸入不正確k<NC(j)=c(j)+X(k)(j-1)YYYYYj<m+2b(j)=b(j)+Y(k)*X(k)(j-1)k=k+1建立矩陣C將c向量的元素依行代入將矩陣b右除C得到結(jié)果向量輸出結(jié)果YNNN1.3 多項(xiàng)式曲線擬合的MATLAB程序代碼function A=multifit(X,Y,m)%離散實(shí)驗(yàn)數(shù)據(jù)點(diǎn)的多項(xiàng)式曲線擬合%實(shí)驗(yàn)數(shù)據(jù)點(diǎn)的x坐標(biāo)向量：X%實(shí)驗(yàn)數(shù)據(jù)點(diǎn)的y坐標(biāo)向量：Y%擬合多項(xiàng)式的次數(shù)：m%擬合多項(xiàng)式的系數(shù)向量：AN=length(X

8、);M=length(Y);if(N = M) disp('數(shù)據(jù)點(diǎn)坐標(biāo)不匹配！'); return;endc(1:(2*m+1)=0;b(1:(m+1)=0;for j=1:(2*m+1) %求出c和b for k=1:N c(j)=c(j)+X(k)(j-1); if(j<(m+2) b(j)=b(j)+Y(k)*X(k)(j-1); end endendC(1,:)=c(1:(m+1);for s=2:(m+1) C(s,:)=c(s:(m+s);endA=b'C; %直接求解法求出擬合系數(shù)1.4 多項(xiàng)式曲線擬合應(yīng)用實(shí)例用二次多項(xiàng)式擬合下表所列的數(shù)據(jù)點(diǎn)。x12

9、34y41018261.4.1 操作流程圖輸入向量X，Y，多項(xiàng)式次數(shù)m.調(diào)用函數(shù)A=multifit(X,Y,m)輸出結(jié)果1.4.2 操作命令>> x=1:4;>> y=4 10 18 26;>> A=multifit(x,y,2)1.4.3 輸出結(jié)果輸出結(jié)果為：A = 0.0489 0.1612 0.5672即擬合的多項(xiàng)式為： P=0.0489+0.1612x+0.5672x2144 結(jié)果如圖、編程解決以下科學(xué)計(jì)算問(wèn)題。2.1 問(wèn)題分析解：建模：由等效電流源電路圖可知各支路導(dǎo)納為：Y1=1/R1+1/(j*XL); Y2=1/R2-1/(j*Xc1); Y

10、3=1/R3-1/(j*Xc2)均為兩并聯(lián)元件導(dǎo)納之和，按照?qǐng)D中電流方向，其電流為I1=Ua*Y1, I2=(Ub-Ua)*Y2, I3=-Ub*Y3則a，b兩點(diǎn)的電流方程為Y1Ua-Y2(Ub-Ua)=Us1/jXL+Us2/R1Y2(Ub-Ua)-Y3Ub=Us3/R3-Us4/jXc2-Us2/R2寫(xiě)成矩陣形式：即可寫(xiě)成AU=B2.2 操作流程圖 開(kāi)始 結(jié)束 輸入數(shù)據(jù)由歐姆定律得出原電路各支路電流計(jì)算Y1，Y2，Y3 求出各支路的電流I1,I2,I3列方程組，建立矩陣AI=B得到Ua,Ub 帶入數(shù)據(jù)，用向量B左除A2.3 程序代碼：function fun1R1=2;R2=3;R3=4;

11、XL=2;XC1=3;XC2=5; %給出原始數(shù)據(jù) us1=8;us2=6;us3=8;us4=15; %給出原始數(shù)據(jù) Y1=1/R1+1/( j*XL); %用復(fù)數(shù)表示各支路導(dǎo)路 Y2=1/R2-1/( j*XC1); Y3 = 1/R3-1 /( j*XC2); A= Y1+ Y2,-Y2;- Y2,Y2+Y3; %按線性方程組列ua，ub的系數(shù)矩陣 B=us1/( j*XL)+us2/R1;us3/R3+us4/(- j*XC2)-us2/R2; %列出線性方程組右端U=AB;ua=U(1),ub=U(2) %求ua,ub I1=ua*Y1,I2=( ub -ua)*Y2,I3=ub*Y

12、3, %求各支路的II1R=ua/R1 ,I1L=ua/( j*XL), I2R=(ub-ua)/R2,I2C=(ub-ua)/(-j*XC1), I3R=ub/R3,I3C=ub/(- j*XC2), W=compass(ua,ub,I1,I2,I3) %畫(huà)向量圖，設(shè)定此圖的圖柄為w set(W,'linewidth',2) %改變向量圖線寬end2.4 運(yùn)行結(jié)果如圖：運(yùn)行>> fun1ua = 4.8845 - 0.5981iub = 5.4874 + 2.5752iI1 = 2.1432 - 2.7413iI2 = -0.8568 + 1.2587iI3 =

13、0.8568 + 1.7413iI1R = 2.4422 - 0.2990iI1L = -0.2990 - 2.4422iI2R = 0.2010 + 1.0578iI2C = -1.0578 + 0.2010iI3R = 1.3718 + 0.6438iI3C = -0.5150 + 1.0975iW = 179.0024 180.0024 181.0024 182.0024 183.00242.5 運(yùn)行結(jié)果截圖2. （2）解：由題要求，可用最小二乘擬合法擬合函數(shù)開(kāi)始流程圖輸入x,y 結(jié)束 調(diào)用擬合函數(shù)a=polyfit(x,y,2) 畫(huà)圖 程序x=0.1 0.4 0.5 0.7 0.7 0.9;y=0.61 0.92 0.99 1.52 1.47 2.03;cc=polyfit(x,y,2)       %求出A與B的系數(shù) xx=x(1):0.1:x(length(x);yy=po