第一篇：中南大學化工原理仿真實驗報告

化工原理計算機仿真實驗

班級：化學工程與工藝1102班

姓名：王翔

學號：1505110321

日期：2014年1月1日

本套軟件系統包括8個單元仿真實驗：

實驗一 離心泵性能的測試

實驗二 管道阻力實驗

實驗三 傳熱實驗

實驗四 吸收實驗

實驗五 流體流動形態的觀測

實驗六 柏努利方程實驗

實驗七 干燥實驗

實驗八 精餾實驗

以下是實驗模擬觀測過程和計算機生成的實驗報告。

圖1 離心泵性能的測試 觀察氣蝕現象（1）

圖2 離心泵性能的測試 觀察氣蝕現象（2）

圖3 離心泵性能的測試 離心泵特性曲線測定實驗報告（1）

圖4 離心泵性能的測試 離心泵特性曲線測定實驗報告（2）

圖5 離心泵性能的測試 離心泵特性曲線測定實驗報告（3）

圖6 離心泵性能的測試 離心泵特性曲線測定實驗報告（4）

圖7 管道阻力的測定實驗報告（1）

圖8 管道阻力的測定實驗報告（2）

圖9 管道阻力的測定實驗報告（3）

圖10 傳熱實驗

圖11 傳熱實驗報告（1）

圖12 傳熱實驗報告（2）

圖13 傳熱實驗報告（3）

圖14 傳熱實驗報告（4）

圖15 吸收實驗 觀察液泛現象

圖16 吸收實驗報告

圖17 液體流動形態的觀測 觀察滯留形態

圖18 液體流動形態的觀測實驗報告

圖19 柏努利方程實驗 觀察測壓孔與水流方向方位角與水位變化（1）

圖20 柏努利方程實驗 觀察測壓孔與水流方向方位角與水位變化（2）

圖21 干燥實驗報告（1）

圖22 干燥實驗報告（2）

圖23 干燥實驗報告（3）

圖24 干燥實驗報告（4）

圖25 精餾實驗 動態平衡調整

圖26 精餾實驗報告（1）

圖27 精餾實驗報告（2）

第二篇：中南大學系統仿真實驗報告

ans =

實驗一 MATLAB 中矩陣與多項式的基本運算 實驗任務 1.了解 MATLAB 命令窗口和程序文件的調用。?熟悉如下 MATLAB 的基本運算：

① 矩陣的產生、數據的輸入、相關元素的顯示； ② 矩陣的加法、乘法、左除、右除； ③ 特殊矩陣：單位矩陣、“ 1 ”矩陣、“ 0 ”矩陣、對角陣、隨機矩陣的產生和 運算； ④ 多項式的運算：多項式求根、多項式之間的乘除。

基本命令訓練 1、>> eye(2)ans = 1 0 0 1 >> eye(4)ans = 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1 2、>> ones(2)1 1

ans =1 >> ones(4)ans = 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 >> ones(2,2)ans = 1 1 1 1 >> ones(2,3)ans = 1 1 1 1 1 1 >> ones(4,3)ans = 1 1 1 1 1 1 1 1 1 1 1 1 3、>> zeros(2)

0 0 0 0 >> zeros(4)ans = 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 >> zeros(2,2)ans = 0 0 0 0 >> zeros(2,3)ans = 0 0 0 0 0 0 >> zeros(3,2)ans = 0 0 0 0 00 4、隨機陣 >> rand(2,3)ans = 0.2785 0.9575 0.1576 0.5469 0.9649 0.9706 >> rand(2,3)

ans = 0.9572 0.8003 0.4218 0.4854 0.1419 0.9157 5、>> diag(5)ans = 5 >> diag(5,5)ans = 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 5 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 >> diag(2,3)ans = 0 0 0 2 0 0 0 0 0 0 0 0 0 0 0 0 6、(inv(A)為求 A 的逆矩陣)>> B=[5 3 1;2 3 8;1 1 1],inv(B)3 1

1

ans =

0.6250 0.2500-2.6250-0.7500-0.5000 4.7500 0.1250 0.2500-1.1250

>> A=[2 3;4 4],B=[5 3;3 8],inv(A),inv(B);AB,A/B,inv(A)*B,B*inv(A)A = 2 3 4 4 B = 5 3 3 8 ans =-1.0000 0.7500 1.0000-0.5000 ans =-2.7500 3.0000 3.5000-1.0000 ans = 0.2258 0.2903 0.6452 0.2581 ans =-2.7500 3.0000

3.5000 -1.0000 ans =-2.0000 2.2500 5.0000-1.7500 7、>> p =[1,-6,-72,-27], roots(p)p = 1-6-72-27 ans = 12.1229-5.7345-0.3884 >> p=[2,3,6],roots(p)p = 2 3 6 ans =-0.7500 + 1.5612i-0.7500-1.5612i 8、(A 為 n*n 的方陣)>> A=[0 1 0;-4 4 0;-2 1 2],poly(A),B=sym(A),poly(B)A = 0 1 0-4 4 0-2 1 2 ans =-6 12-8

B = [ 0, 1, 0] [-4, 4, 0] [-2, 1, 2]

ans = x A 3-6*x A 2+12*x-8 9,、(conv 是多項式相乘，deconv 是多項式相除)>> u=[1 2 4 6 ],v=[5 0 0-6 7],conv(u,v)u = 1 2 4 6 v =0 0-6 7

ans =10 20 24-5-10-8 42

>> v=[1 2 4 6 ],u=[5 0 0-6 7],deconv(u,v)v = 1 2 4 6 u = 5 0 0-6 7 ans = 5-10 10、（點乘是數組的運算，沒有點的乘是矩陣運算）

>> a = [2 5;3 4], b =[3 1;4 7],a.*b,a*b a = 2 5 3 4 b = 3 1 4 7 ans = 6 5 12 28 ans = 26 37 25 31 >> a = [2 3];b = [4 7];a.*b = [8 21];a*b %錯誤 a*b“ = 29;11、（who 可以看到你用過的一些變量，來了）

>> who Your variables are: A B a ans b p u >> whos Name Size Bytes whos 是把該變量及所存儲的大小等信息都顯示出 Class Attributes 2x2 32 double

B 2x2 32 double a 1x2 16 double ans 1x2 16 double b 1x2 16 double p

1x3 24 double u 1x5 40 double v 1x4 32 double12、>> A=[2 5 3;6 5 4],disp(A),size(A),length(A)A = 2 5 3 6 5 4 2 5 3 6 5 4 ans = 2 3 ans = 3 實驗二 MATLAB 繪圖命令 實驗任務 熟悉 MATLAB 基本繪圖命令，掌握如下繪圖方法：．坐標系的選擇、圖形的繪制；?圖形注解（題目、標號、說明、分格線）的加入;3 ?圖形線型、符號、顏色的選取 基本命令訓練 1、>>t=[0:pi/360:2*pi];x=cos(t)+ cos(t*4);y=si n(t)+ sin(t*4);xlabel(”x 軸“);ylabel(”y 軸“);plot(y,x),grid;2、>> t=0:0.1:100;x=3*t;y=4*t;z=si n(2*t);plot3(x,y,z, ”g:“)

■15 i 0 5 0 05 1 1 5 2

3、>>x = linspace(-2*pi,2*pi,40);y=si n(x);stairs(x,y)4、>> t=[0:pi/360:2*pi];x=cos(t)+ cos(t*4)+ sin(t*4);y=si n(t)+ si n(t*4);plot(y,x, ”r:“);

xlabel(”x 軸“);ylabel(”y 軸“);

6、>>th=[0:pi/20:2*pi];x=exp(j*th);plot(real(x),imag(x),”r-.“);grid;text(0,0,”中心“);

5、>> th=[0:pi/1000:2*pi]”;r=cos(2*th);polar(th,r);title(“四葉草圖”)270 四葉草圖

7、>>x=-2:0.01:2;8、y=-2:0.01:2;9、[X,Y] = meshgrid(x,y);Z = Y.*exp(-X.A 2-Y.A 2);[C,h] = con tour(X, YZ);set(h,“ShowText”,“o n”,“TextStep”,get(h,“LevelStep”)*2)_ 1OS ?I ,5 2 ■n.s o o.s 8、>>x = 0:0.2:10;y = 2*x+3;subplot(411);plot(x,y);grid;title(“y 的原函數”);subplot(412);semilogy(x,y);grid;title(“對 y 取對數”);丫的原畫數 40----------1------------------------------------------1--------------------1--------------------1-----------------------------------------1--------------------1-------------------1--------

l| I p | il ■ | i| I 九 _____ 1-___ — ___ I ____ L ___ :…… ：

___ J _

_______ L ___ u i| I |l I , il _-■」 “ j I ■I __ h-_____________ I I ■ Q 」 【 I 1

F I I I I II I I IT 10 1

□ 1 2 3 4 5 6 r 6 9 10 , 對 y 取對數 對弋觀對數 subplot(413);semilogx(x,y);10 10 id 10 1lZ 10 w 10 10 40 20 0 對好對數 ?掌握循環、分支語句的編寫，學會使用 look for、help 命令

grid;title(”對 x 取對數“);subplot(414);loglog(x,y);grid;title(”對 xy 均取對數“);9、>>x =-3:0.3:3;bar(x,exp(-x.*x),”g“)實驗三 MATLAB 程序設計 實驗任務 1 ?熟悉 MATLAB 程序設計的方法和思路；

程序舉例 1、>> f=[1,1];i=1;while f(i)+f(i+1)<1000 f(i+2)=f(i)+f(i+1);i=i+1;end f,i f = Columns 1 through 14 1 1 2 3 5 8 13 21 34 55 89 144 233 377 Columns 15 through 16 610 987 i = 15 2、>> m=3;n=4;for i=1:m for j=1:n a(i,j)=1/(i+j-1);end end format rat1/2 1/3 1/2 1/3 1/4 1/3 1/4 1/5(分數格式形式。用有理數逼近顯示數據)>>m=5;n=4;for i=1:m for j=1:n a(i,j)=1/(i+j-1);end end format rat a a = 1 1/2 1/3 1/4 1/2 1/3 1/4 1/5 1/3 1/4 1/5 1/6 1/4 1/5 1/6 1/7 1/5 1/6 1/7 1/8 3、程序中沒有 format rat 命令時，如果上次運行結果沒有清除，輸出的結果就是上次運行的 結果！但是運用 clear 命令清楚之前的運行結果之后就會正常運行。

4、>> x=input（”請輸入 x 的值:’）;if x==10 y=cos(x+1)+sqrt(x*x+1);else y=x*sqrt(x+sqrt(x));1/4 1/5 1/6

end y

請輸入 x 的值 :2 y = 2391/647 x=input(“請輸入 x 的值:’);if x==10 y= fprintf(” 不在定義域內，請重新輸入：

“);return else y=1/(x-10);end y 請輸入 x 的值 :2-1/8 5、>> p=[0 0 0 1 3 0 2 0 0 9];for i=1:length(p), if p(1)==0,p=p(2:length(p));end;end;p

p = Columns 1 through 5 3

0

0 Columns 6 through 7

0 9

>>p=[0 0 0 1 3 0 2 0 0 9];p(p==0)=[];p

p =3

6、>> e2(500)

ans =1 2 3 5 8 13 21 34 55 89 144 233

377 >> lookfor ffibno e2-ffibno 計算斐波那契亞數列的函數文件 >> help e2 ffibno 計算斐波那契亞數列的函數文件 n 可取任意自然數 程序如下（用法：

lookfor 關鍵詞 在所有 M 文件中找 “關鍵詞 ”，比如：

lookfor max（即尋找關鍵詞 “ max”）

其實就和我們平時用 CTRL+F 來查找關鍵詞”是一樣的 而 help 是顯示 matlab 內置的幫助信息 用法：

help 命令，比如 help inv，作用就是調用 inv 這個命令的幫助）

程序設計題 用一個 MATLAB 語言編寫一個程序：輸入一個自然數，判斷它是否是素數，如果是，輸出“ It is one prime ”，如果不是，輸出“ It is n ot o ne prime.”。要求通 過調用子函數實現。

最好能具有如下功能：

①設計較好的人機對話界面，程序中 含有提示

性的輸入輸出語句。

②能實現循環操作，由操作者輸入相關命令來控制 是否繼續進行素數的判斷。

如果操作者希望停止這種判斷，則可以退出程序。

③ 如果所輸入的自然數是一個合數，除了給出其不是素數的結論外，還應給出至少 一種其因數分解形式。例：輸入 6，因為 6 不是素數。則程序中除了有“ It is not one prime ”的結論外，還應有：“ 6=2*3 ”的說明。

function sushu while 1 x=input（” 請輸入一個自然數 “）;if x<2 disp（” 既不是質數又不是合數 “）;Else if isprime(x)==1 disp（” 這是一個素數 “）;Else disp(” 這是一個合數，可以因式分解為：

“)for n=2:sqrt(x)if rem(x,n)==0 num3=x;num1=n;num2=x/n;disp([num2str(num3),”=“,num2str(num1),”x“,num2str(num2)])End End End

end y=input(” 是否繼續判斷？繼續請按 1 ,按任意鍵退出 ：

“)； if y~=1 break end end 實驗四 MATLAB 的符號計算與 SIMULINK 的使用 實驗任務 1.掌握 MATLAB 符號計算的特點和常用基本命令;2.掌握 SIMULINK 的使用。

程序舉例 1.求矩陣對應的行列式和特征根 >> a=sym(”[a11 a12;a21 a22]“);da=det(a)ea=eig(a)da = a11*a22-a12*a21 ea = 1/2*a11+1/2*a22+1/2*(a11 A 2-2*a11*a22+a22 A 2+4*a12*a21)A(1 /2)1/2*a11+1/2*a22-1/2*(a11 A 2-2*a11*a22+a22 A 2+4*a12*a21)A(1 /2)a=sym(”[2 3;1 5]“);

da=det(a)ea=eig(a)da = 7 ea = 7/2+1/2*21A(1/2)7/2-1/2*21A(1/2)2.求方程的解(包括精確解和一定精度的解)>> r1=solve(”xA2+x-1“)rv=vpa(r1)rv4=vpa(r1,4)

rv20=vpa(r1,20)r1 = 1/2*5 八(1/2)-1/2-1/2*5 八(1/2)-1/2 rv =.6894848260-1.689484826 rv4 =.6180-1.618 rv20 =.689484820-1.68948482 3 ． a=sym(”a“);b=sym(”b“);c=sym(w=10;x=5;y=-8;z=11;A=[a,b;c,d] B=[w,x;y,z] det(A)det(B)A = [ a, b] ”c“);d=sym(”d“);%定義 4 個符號變量 %定義 4 個數值變量 %建立符號矩陣 A %建立數值矩陣 B %計算符號矩陣 A 的行列式 %計算數值矩陣 B 的行列式

[ c, d]

X

10-8 11 ans = a*d-b*c ans = 150 4.? syms x y;s=(-7*x A 2-8*y A 2)*(-x A 2+3*y A 2);expand(s)% 對 s 展開 collect(s,x)% 對 s 按變量 x 合并同類項（無同類項）

factor(ans)% 對 ans 分解因式

ans = 7*xA4-13*xA2*yA2-24*yA4 ans = 7*xA4-13*xA2*yA2-24*yA4 ans =(8*yA2+7*xA2)*(xA2-3*yA2)5.對方程 AX=b 求解 >> A=[34,8,4;3,34,3;3,6,8];b=[4;6;2];X=linsolve(A,b)% 調用 linsolve 函數求解 Ab % 用另一種方法求解

0.0675 0.1614

diff(f)%未指定求導變量和階數，按缺省規則處理

0.1037 ans = 0.0675 0.1614 0.1037 6 ． 對方程組求解 a11*x1+a12*x2+a13*x3=b1 a21*x1+a22*x2+a23*x3=b2 a31*x1+a32*x2+a33*x3=b3 >>syms a11 a12 a13 a21 a22 a23 a31 a32 a33 b1 b2 b3;A=[a11,a12,a13;a21,a22,a23;a31,a32,a33];b=[b1;b2;b3];XX=Ab %用左除運算求解(X=linsolve(A,b)%調用 linsolve 函數求的解)XX =(a12*a23*b3-a12*b2*a33+a13*a32*b2-a13*a22*b3+b1*a22*a33-b1*a32*a23)/(a11*a22*a33-a11*a32*a23-a12*a21*a33+a32*a21*a13-a22*a31*a13+a31*a12*a 23)-(a11*a23*b3-a11*b2*a33-a21*a13*b3-a23*a31*b1+b2*a31*a13+a21*b1*a33)/(a11*a22*a33-a11*a32*a23-a12*a21*a33+a32*a21*a13-a22*a31*a13+a31*a12* a23)(a32*a21*b1-a11*a32*b2+a11*a22*b3-a22*a31*b1-a12*a21*b3+a31*a12*b2)/(a11*a22*a33-a11*a32*a23-a12*a21*a33+a32*a21*a13-a22*a31*a13+a31*a12*a 23)7 ． syms a b t x y z;f=sqrt(1+exp(x));

%求 f 對 x 的二階導數 %求 f 對 x 的三階導數 %按參數方程求導公式求 y 對 x 的導數 ans = 1/2/(1+exp(x))A(1/2)*exp(x)ans =-2*si n(x)-x*cos(x)ans =-3*cos(x)+x*si n(x)ans =-b*cos(t)/a/si n(t)三、SIMULINK 的使用

f=x*cos(x);diff(f,x,2)diff(f,x,3)f1=a*cos(t);f2=b*si n(t);diff(f2)/diff(f1)

diff(f)%未指定求導變量和階數，按缺省規則處理

仿真圖: 波形圖: 其中: G i(R

0.14

實驗五 MATLAB 在控制系統分析中的應用 實驗任務 1.掌握 MATLAB 在控制系統時間響應分析中的應用;2.掌握 MATLAB 在系統根軌跡分析中的應用； 3.掌握 MATLAB 控制系統頻率分析中的應用； 4.掌握 MATLAB 在控制系統穩定性分析中的應用 基本命令 1.step 2.impulse 3.in itial 4.Isim 5.rlocfi nd 6.bode 7.margin 8.nyquist 9.Nichols 10.cloop 程序舉例 1.求下面系統的單位階躍響應 1.5 G(s)

t >> num=[4];den=[1 , 1 , 4];step(num , den)[y , x , t]=step(num , den);tp=spli ne(y , t , max(y))% 計算峰值時間 max(y)% 計算峰值 tp = Step Response 0.5 Time(sec)p m 1.6062 ans = 1.4441 0.18 0.16

2.求如下系統的單位階躍響應 X i 0 1 X i 0 u X 2 6 5 X 2 1

>> a=[0,1;-6,-5];b=[0;1];c=[1,0];d=0;[y,x]=step(a,b,c,d);Plot(y)

1,0 X i X 2 3.求下面系統的單位脈沖響應: G(s)兒 >> num=[4];den=[1 , 1 ,4];impulse(nu m,de n)Time(sec)Response to Initial Conditions 4.已知二階系統的狀態方程為: e u p m A 0.25 0.2 0.15 0.1 0.05 0-0.05

L I 1-

r

-廠、、/ Impulse Response

c=[1 , 0] ；

d=[0];xO=[1 ,0] ； subplot(1 , 2,1);in itial(a , b , c ,d,x0)subplot(1 , 2,2);impulse(a , b , c , d):系統傳遞函數為: G(s)占 輸入正弦信號時，觀察輸出 信號的相位差。

>> num=[1];den=[1 ,1];t=0 : 0.01 : 10;u=s in(2*t);hold on plot(t,u, ”r“)lsim(nu m,de n,u,t)

Real Axis

6.有一二階系統，求出周期為 4 秒的方波的輸出響應 2s 2

5s 1 ~2 s 2s 3 >>num=[2 5 1];den=[1 2 3];t=(0:.1:10);period=4;u=(rem(t,period)>=period./2);% 看 rem 函數功能 lsim(nu m,de n,u,t);7.已知開環系統傳遞函數，繪制系統的根軌跡，并分析其穩定性

G(s)k(s 2)(s 2

4s 3)2

>>num=[1 2];den 仁 [1 4 3];den=conv(de n1,de n1);figure(1)-6 rlocus(num,den)[k,p]= rlocfi nd(nu m,de n)l L l ■ 匚一 L — C--------------------

---------

-r t r r

-10-8-6-4-2 4 6 Root Locus 0 2 2 2

0 0 2 2--CPXA vyanma n-8 2.5 G(s)Lin ear Simulati on Results1.5 e 0.5-0.5-1-1.5-2 Time(sec)p m

impulse response(k=55)

5

figure(2)k=55;num 仁 k*[1 2];den=[1 4 3];den 1=c onv(de n,den);[nu m,de n]=cloop(nu m1,de n1,-1);impulse(nu m,de n)title(”impulse resp onse(k=55)“)-1.5-1 0 200 400 600 Time(sec)800 1000 1200 5 5

0 05

0 0

o

--x 10 8!-------

impulse resp on se(k=56)figure(3)k=56;num 仁 k*[1 2];den=[1 4 3];den 1=c onv(de n,den);[nu m,de n]=cloop(nu m1,de n1,-1);impulse(nu m,de n)-2-4-6-8 L

0 500 1000 Time(sec)title(”impulse resp on se(k=56)“)1500 2000 2500

Select a point in the graphics win dow selected_po int =-2.5924-0.0248i 0.7133-3.4160-2.5918-0.9961 + 0.4306i-0.9961-0.4306i Bode Diagram 8.作如下系統的 bode 圖 G(s)>> n=[1 , 1];d=[1 , 4,11 , 7];bode(n , d),grid on Frequency(rad/sec)s 1 s 3

4s 2

11s 7 9.系統傳函如下 G(s)s 1 0.5s e(s 2)3 求有理傳函的頻率響應，然后在同一張圖上繪出以四階伯德近似表示的系統頻 率響應 >> num=[1];de n=co nv([1 2],co nv([1 2],[1 2]));

w=logspace(-1,2);t=0.5;

ylabel(”gai n“);subplot(2,1,2);semilogx(w,p1,w,p2,”g--“);grid on;xlabel(”freque ncy“);ylabel(”phase“);10.已知系統模型為

求它的幅值裕度和相角裕度 >> n=[3.5];d=[1 2 3 2];[Gm,Pm,Wcg,Wcp]=margi n(n,d)

G(s)3.5 s 3

2s 2

3s 2 bode plot [m1,p1]=bode(nu m,de n, 2);p1=p1-t*w”*180/pi;[n 2,d2]=pade(t,4);nu mt=c onv(n2,nu m);den t=(c onv(de n, d2));2 freque xlabel(“freque

Gm =

1.1433 Pm =

7.1688 Wcg =

1.7323 Wcp =

1.6541

nyq uist(n, d1);hold on nyq uist(n, d2);nyq uist(n, d3);nyq uist(n, d4);Nyquist Diagram 80 |--------------------

--------------------[--------------------[---------

11.二階系統為: G(s)n s 2n S S 令 wn=1, 分別作出 E =2,1 , 0.707 , C A 5 y

a 時的 nyquist 曲線。

m >> n=[1];d 仁 [1 , 4,1];d2=[1 , 2 , 1] ； d3=[1 , 1.414,1];d4=[1,1,1];Nyquist Diagram Real Axis 12.已知系統的開環傳遞函數為

_一-”—一 一-一_一_ ~-60

005

S S

2)

S S32

S S

v vr ra an nk ky ya a

14.一多環系統，其結構圖如下，使用 Nyquist 頻率曲線判斷系統的穩定性。

16.7s(0.85s 1)(0.25s 1)(0.0625s 1)繪制系統的 Nyqusit 圖, 并討論系統的穩定性.>> G=tf(1000,co nv([1,3,2],[1,5]));nyquist(G);axis(”square“)13.分別由 w 的自動變量和人工變量作下列系統的 nyquisF 曲線：

m

G(s)1 s(s 1)>> n=[1];d=[1 , 1 ,0];nyquist(n ,d);% 自動變量 n=[1];d=[1 , 1 ,0];w=[0.5 : 0.1 : 3];nyq uist(n , d , w);% 人工變量-2-1-1.5-0.9-0.8-0.7-0.6-0.5 Real Axis-0.4-0.3-0.2-0.1 5 5

--0 0

S SI IX XA A

u u^ ^a an nT T9 9a a

卩

G(s)

30 20 10 0

figure(2)[nu m2,de n2]=cloop(nu m,de n);impulse(nu m2,de n2);

>> k1=16.70.0125;z1=[0];p1=[-1.25-4-16];[nu m1,de n1]=zp2tf(z1,p1,k1);[nu m,de n]=cloop(nu m1,de n1);[z,p,k]=tf2zp(nu m,de n);p figure(1)nyq uist(nu m,de n)-2-1 Nyquist Diagram-0.5 0 0.5 Real Axis 1 1.5 1 0.5 s a 0 n g m-0.5-1-1.5 1.5 20 Impulse Resp onse-10.5969 +36.2148i-10.5969-36.2148i-0.0562 15.已知系統為: eanML—pm

-5-10-15 0 0.1 0.2 0.3 Time(sec)Nichols Chart 0.4 0.5 0.6-2 10 0 10 2

Open-Loop Phase(deg)Frequency(rad/sec)

G(s)s(s 1)作該系統的 nichols 曲線。

>> n=[1];d=[1 , 1 , 0];ni chols(n , d);16.已知系統的開環傳遞函數為: G(s)k s(s 1)(s 2)當 k=2 時，分別作 nichols 曲線和波特圖 >> num=1;den=conv(co nv([1 0],[1 1]),[0.5 1]);subplot(1,2,1);ni chols(nu m,de n);grid;B

% n ichols 曲線 G subplot(1,2,2);2 g=tf(nu m,de n);bode(feedback(g,1,-1));grid;Nichols Chart koflr ea M^hnaaM

msuvesaB—

50

--

00

-270 Bode Diagram

90--

% 波特圖 17.系統的開環傳遞函數為:

分別確定 k=2 和 k=10 時閉環系統的穩定性 >> d 仁 [1 , 3,2,0];n 仁 [2];[nc1 , dc1]=cloop(n1 , d1 ,-1);roots(dc1)d2=d1;n2=[10];[nc2 , dc2]=cloop(n2 , d2,-1);roots(dc2)ans =-2.5214-0.2393 + 0.8579i-0.2393-0.8579i ans =-3.3089 0.1545 + 1.7316i 0.1545-1.7316i 18.系統的狀態方程為：

X 1 4 3 0 X 1 1 X 2 1 0 0 X 2 0 u X 3 0 1 0 X 3 0

X 1

y 0 1 2 x 2

X 3

G(s)k s(s 1)(s 2)

試確定系統的穩定性。

>> a=[-4,-3,0;1,0,0;0,1,0];b=[1;0;0];c=[0,1,2];d=0;eig(a)% 求特征根 ran k(ctrb(a,b))ans = 0-1-3 ans = 3 實驗六連續系統數字仿真的基本算法 實驗任務 1.理解歐拉法和龍格-庫塔法的基本思想； 2 ?理解數值積分算法的計算精度、速度、穩定性與步長的關系;程序舉例 1.取 h=0.2 ,試分別用歐拉法、RK2 法和 RK4 法求解微分方程的數值解，并 比較計算精度。

注：解析解：y(t).1 2t >> clear t(1)=0 ； y(1)=1;y_euler(1)=1;y_rk2(1)=1;y_rk4(1)=1;h=0.001;% 步長修改為 0.001 for k=1:5 y(t)y(t)y(0)2t y(t)1

t(k+1)=t(k)+h;y(k+1)=sqrt(1+2*t(k+1));end for k=1:5 y_euler(k+1)=y_euler(k)+h*(y_euler(k)-2*t(k)/y_euler(k));end for k=1:5 k1= y_rk2(k)-2*t(k)/y_rk2(k);k2=(y_rk2(k)+h*k1)-2*(t(k)+h)/(y_rk2(k)+h*k1);y_rk2(k+1)=y_rk2(k)+h*(k1+k2)/2;end for k=1:5 k1= y_rk4(k)-2*t(k)/y_rk4(k);k2=(y_rk4(k)+h*k1/2)-2*(t(k)+h/2)/(y_rk4(k)+h*k1/2);k3=(y_rk4(k)+h*k22)-2*(t(k)+h/2)/(y_rk4(k)+h*k2/2);k4=(y_rk4(k)+h*k3)-2*(t(k)+h)/(y_rk4(k)+h*k3);

y_rk4(k+1)=y_rk4(k)+h*(k1+2*k2+2*k3+k4)/6;end disp(” 時間 解析解 yt=[t“, y”, y_euler“, y_rk2”, y_rk4“];disp(yt)歐拉法 RK2 法 RK4 法 ”)時間 解析解 歐拉法 RK2 法 RK4法0 1.0000 1.0000 1.0000 1.0000 0.0010 1.0010 1.0010 1.0010 1.0010 0.0020 1.0020 1.0020 1.0020 1.0020 0.0030 1.0030 1.0030 1.0030 1.0030 0.0040 1.0040 1.0040 1.0040 1.0040 0.0050 1.0050 1.0050 1.0050 1.0050

y(t)2Ry(t)y(t)0 在 0 t 10 上的數字仿真解(已知：

y(0)y(0)0)，并將不同步長下的仿真結果與解析解進行精度比較。

說明：

已知該微分方程的解析解分別為：

100，yt y t 100cost(當 R 0)100e 2t c°s 仝 t 10 ^ e*s in 仝 t 2 3 2(當 R 0.5)采用 RK4 法進行計算，選擇狀態變量: 2.考慮如下二階系統：

x 1 y x 2 y 則有如下狀態空間模型及初值條件 x 1 x 2 x 1(0)100 x 2 x 1 2Rx 2 x 2(0)0 y x 1 采用 RK4 法進行計算。

>> clear h=input(“ 請輸入步長 h=”);M=round(10/h);t(1)=0;y_0(1)=100;y_05(1)=100;和 y_05 分別對應于為 R=0 和 R=0.5)x1(1)=100;x2(1)=0;y_rk4_0(1)=x1(1);y_rk4_05(1)=x1(1);% 求解析解 for k=1:M t(k+1)=t(k)+h;y_0(k+1)=100*cos(t(k+1));y_05(k+1)=100*exp(-t(k+1)/2).*cos(sqrt(3)/2*t(k+1))+100*sqrt(3)/3*exp(-t(k+1)/2).*si n(sqrt(3)/2*t(k+1));end% 輸入步長 % 置總計算步數 % 置自變量初值 % 置解析解的初始值(y_0 % 置狀態向量初值 % 置數值解的初值

% 利用 RK4 法求解 % R=0 for k=1:M k11=x2(k);k12=-x1(k);k21=x2(k)+h*k12/2;k22=-(x1(k)+h*k11/2);k31=x2(k)+h*k22/2;k32=-(x1(k)+h*k21/2);k41=x2(k)+h*k32;k42=-(x1(k)+h*k31);x1(k+1)=x1(k)+h*(k11+2*k21+2*k31+k41)/6;x2(k+1)=x2(k)+h*(k12+2*k22+2*k32+k42)/6;y_rk4_0(k+1)=x1(k+1);end % R=0.5 for k=1:M k11=x2(k);k12=-x1(k)-x2(k);k21=x2(k)+h*k12/2;k22=-(x1(k)+h*k11/2)-(x2(k)+h*k12/2);k31=x2(k)+h*k22/2;k32=-(x1(k)+h*k21/2)-(x2(k)+h*k22/2);k41=x2(k)+h*k32;k42=-(x1(k)+h*k31)-(x2(k)+h*k32);x1(k+1)=x1(k)+h*(k11+2*k21+2*k31+k41)/6;x2(k+1)=x2(k)+h*(k12+2*k22+2*k32+k42)/6;y_rk4_05(k+1)=x1(k+1);end % 求出誤差最大值 err_0=max(abs(y_0-y_rk4_0));err_05=max(abs(y_05-y_rk4_05));% 輸出結果 disp（“ 最大誤差（R=0）

最大誤差（R=0.5）”）

err_max=[err_0,err_05];disp（err_max）

請輸入步長 h=0.5 最大誤差（R=0）

最大誤差（R=0.5）

0.4299 0.0460 沒運行一次程序，輸入一個步長，記得可到相應的最大誤差，將結果可列表如 下：

步長 h 0.0001 0.0005 0.001 0.005 0.01 0.05 0.1 0.5

F =0 5.4330 1.6969 1.0574 4.1107 6.6029 4.1439 6.6602 4.2988 差 最大誤 X 10-10 X 10-10 X 10-10 X 10-9 X 10-8 X 10-5 X 10-4 X 10-1

R=0.5 2.7649 6.8123 5.3753 4.0902 6.5425 4.1365 6.7152 4.5976 差 最大誤 X 10-11 X 10-12 X 10-12 X 10-10 X 10-9 X 10-6 X 10-5 X 10-2 從上表中可以看出，當步長 h=0.001 時，總誤差最小；當步長 h 小于 0.001 時，由于舍入誤差變大而使總誤差增加；當步長 h 大于 0.001 時，則由于截斷誤 差的增加也使得總誤差加大。另外，當系統的解變化激烈時（如 R=0），誤差對 步長的變化較為敏感；當系統的解變化平穩時，步長的變化對誤差的影響就要緩 和得多。數值積分算法確定以后，在選擇步長時，需要綜合考慮。

第三篇：化工原理實驗報告封面格式

課程名稱： 化工原理實驗 課題

系院：

專業：

實驗人員：

學號：

同組人員：

實驗時間：

指導教師：流體流動阻力的測定生物與化學工程系11級應用化學朱穎妍1***

47、戴軍凡

49、陳謝萍

51、徐曉麗

52、丁忱2014年2 月 25 日楊金杯

第四篇：化工原理實驗報告評分標準

化工原理實驗報告評分標準

化工原理實驗報告評分分為五部分：

預習20%，考勤10%，實驗操作10%，實驗數據處理40%，思考題20%。

1.預習評分標準

進實驗室之前，仔細閱讀實驗講義相關內容，認真寫好實驗預習報告。完整的預習報告應包括如下內容：實驗名稱、實驗班級、姓名、學號以及實驗目的、實驗內容、基本原理、設備流程(設備工藝流程圖)、實驗步驟及注意事項等。缺內容者酌情扣分。

2.考勤評分標準

按時到實驗室簽到，著裝整齊者為10分。每遲到一分鐘扣一分。

3.實驗操作評分標準

正確操作實驗裝置者得10分。

4.實驗數據處理評分標準

實驗數據及中間處理過程表格化，每一計算步驟需要有詳細的計算舉例，計算公式正確，結果合理40分。

5.思考題評分標準

結合實驗，認真完整地回答實驗講義后面的思考題，且回答正確20分。

第五篇：化工原理仿真系統研究論文

一、化工原理仿真系統的制作

化工原理實驗包括流體流動阻力測定、離心泵性能測定、傳熱、精餾、吸收與解吸、干燥、萃取等基本單元操作，分別由不同的儀器儀表和管道組合而成。在仿真軟件中，把各種設備和管道用flash畫出，再根據每一套裝置流程圖的要求，以真實、立體的效果來實現。

1.整體結構。實驗仿真系統的開發過程分為三個階段：實驗前的準備、實驗過程及數據記錄和數據處理。前兩個階段在Flash動畫制作軟件上完成，第三階段在VisualStudio2005軟件開發工具上完成，并且使用Ac-cess數據庫進行數據的存儲與交換。

2.仿真系統的實現。在計算機模擬化工原理實驗時，需要通過動態數學模型來模擬真實的實驗操作，該模型主要包括實驗指導、素材演示、仿真操作、數據處理、考題測試、幫助功能等內容。下面以離心泵性能測定為例詳細說明仿真系統的制作過程。在實驗準備階段與實驗開始階段的Flash動畫的制作過程中，考慮到實驗步驟有先后，以及更好地做到人機交互，必須使用專門為Flash腳本開發的ActionScript語言。如點擊水泵開啟按鈕必須在閥門開啟以后才能啟動，直至水灌滿后，才可以點擊關閉水泵按鈕。為了使實驗更具有真實性，需設置閥門的流量控制，分為10個級別，可以逐漸增大或減小。運用VisualStudio.Net開發環境編寫C#程序，可以通過拖動添加組件，并自動生成組件需要的代碼。在制作化工原理實驗模擬課件時，可通過VisualStudio屬性窗口設置各種開發元素屬性如外觀、名稱等，且屬性窗口中顯示的內容，隨著選擇開發元素的不同而動態改變。利用VisualStudi“o工具箱”，可以向應用程序添加標準控件。在設置好窗體和控件后，利用Vi-sualStudio的代碼編輯器編寫程序代碼。在命令窗口中，可以直接輸入并執行各種命令，調試應用程序，并通過在即時窗口的命令行中輸入表達式或變量名，可以得到它們的值。編寫程序過程中，難免會遇到一些錯誤，開發人員需要對應用程序進行調試，查找錯誤的根源，以期達到設計要求。離心泵性能測定實驗涉及到流體流動、水泵運轉、儀表變化、閥門打開或關閉等動作，在仿真系統中通過Flash動畫來實現這些動作的動態效果，使整個實驗過程表現得更加真實。用Flash中的按鈕實現動畫交互效果，控制整個實驗的操作并對數據進行采集，同時將數據傳入C#，由C#對數據庫進行讀寫操作，然后作出離心泵特性曲線圖。

3.實驗數據產生及處理。化工原理實驗過程中往往要測定溫度、壓強、濃度、流速等數據，同時必須對這些參數進行整理和分析，并運用相關的理論公式進行計算，才能達到實驗預期目的。化工原理實驗實測數據多，繪圖耗時費力，計算公式復雜，有時甚至需要進行迭代計算，借助計算機輔助程序可圓滿解決這些問題。在仿真軟件中，通過C#語言設計數據處理程序。根據各化工單元操作理論建立數學模型，使仿真數據在實際操作的數據范圍內隨機產生，以保證每個學生在進行仿真實驗時即使初始條件相同，也不會得到完全相同的實驗結果，更接近真實操作狀況。試驗完成后，點擊“記錄數據”按鈕，計算機會自動記錄數據，并在后臺進行數據傳遞，然后根據預先輸入的計算公式進行數據處理。數據處理后被保存到Access數據庫中，再通過調用，將數據輸出在DataGridView進行顯示，或據此數據繪制實驗曲線。

二、操作過程及功能概述

主界面使用VisualStudioC#中的窗體，通過添加菜單欄來控制試驗的選擇。其特點是方便、簡單易用，更重要的是為今后仿真系統的逐步完善提供了空間。首先水泵的開關按鈕是不可用的，必須在打開閥門以后，才能啟動水泵。當水泵與閥門同時啟動后，便開始灌水，在這期間禁用系統中所有的按鈕。待灌水過程結束，先關閥門再關水泵。點擊“開始實驗”按鈕，可以開啟下一個界面繼續實驗。先打開水泵，然后打開閥門，通過閥門調節流量級別，儀表數值會隨之變化。點擊“記錄數據”按鈕，將儀表的數值記錄在數據庫中。當數據記錄完畢，點擊“查看數據”按鈕，屏幕上顯示10組數據以及由公式計算得出的“揚程”、“有效功率”、“效率”數值。點擊“繪圖”按鈕，可直接繪制出H-Q、P-Q及η-Q三條特性曲線。無論是實驗結束還是中途關閉實驗窗體，都將出現一個對話框以提示實驗者“是否保存當前數據？”操作者可根據提示對實驗數據進行取舍。集合Flash動畫和C#語言優點開發的化工原理實驗仿真系統，具有界面直觀、操作簡單、支持人機交互、占用空間小等特點，能顯著提高化工原理實驗教學的效果，減少實驗設備投資和損壞，降低實驗投入成本，避免實驗事故的發生。計算機輔助教學，特別是計算機仿真實驗在化工教學過程中的應用，使學生接觸了一種全新的實驗手段，激發了學生學習的積極性和主動性，使學生創新意識得到培養，從而提高了整體教學質量。