|
|
【答案】应助回帖
jiakai003(金币+51): 真NB,其实我对程序一窍不通，非常感谢你，我以后也要努力学习下这方面的知识，略表心意，希望笑纳，以后还望多多指教！！！ 2011-05-14 20:27:00
function liweimin
clear
clc
format long
global Cin1 Cin2 Cin3 Cin4 Cin5 k1 k2 k3 k4 k5 k6 d1 d2 d3 d4 m Macoh b2 b4 b6 Cinh2o Cgh2o Clh2o t
global tol tog tl cpacoh cph2o cppx cpair Cinacoh Cgacoh Clacoh Fing Hvacoh Hvh2o Cinpx Hr Finl0 Fout0 L G a b c d
%parameter values
Cin1=2.85; Cin2=0; Cin3=0; Cin4=0; Cin5=0; k1=5.50e5; k2=15.88e3; k3=17.09e3; k4=3.28e3; k5=9.81e3; k6=0.54;
d1=1.4247; d2=0; d3=0; d4=4.8419; m=0.0146; Macoh=182.1; b2=0.5254; b4=0.8111; b6=0.9302;
Cinh2o=8.33; Cgh2o=3.124; Clh2o=1.98; t=469; tog=433; tol=412; tl=363; cpacoh=140; cph2o=4200; cppx=214; cpair=50.9;
Cinacoh=0.0167; Cgacoh=0.0167; Clacoh=0.0167; Fing=2.12; Hvacoh=19.9; Hvh2o=35.3; Cinpx=2.85; Hr=1332.2;
C10=0.025; C20=0.025; C30=0.025; C40=0.005; C50=0.005; C60=0.005; Finl0=1.05e5; Fout0=1.00e5; L=3000; G=1000; a=2; b=0.5; c=1.5; d=2;
Finli=Finl0:1000:Finl0+5000; %initial conditions
n=length(Finli);
for i=1:n
sol=fsolve(@Equations,[C10,C20,C30,C40,C50,C60,Finli(i),Fout0],optimset('Display','off'));
y(:,i)=real(sol);
end
y
%--------------------------------------------------------------------------
function f= Equations(y)
global Cin1 Cin2 Cin3 Cin4 Cin5 k1 k2 k3 k4 k5 k6 d1 d2 d3 d4 m Macoh b2 b4 b6 Cinh2o Cgh2o Clh2o t
global tol tog tl cpacoh cph2o cppx cpair Cinacoh Cgacoh Clacoh Fing Hvacoh Hvh2o Cinpx Hr L G a b c d
C1=y(1); C2=y(2); C3=y(3); C4=y(4); C5=y(5); C6=y(6); Finl=y(7); Fout=y(8);
f1=Fout*C1+(k1*C1/(d1*C1+d2*C2+d3*C3+d4*C4+m)+k5*C5/(d1*C1+d2*C2+d3*C3+d4*C4+m))*Macoh-Finl*Cin1;
f2=Fout*C2+(k2*C2/(d1*C1+d2*C2+d3*C3+d4*C4+m)^b2-k1*C1/(d1*C1+d2*C2+d3*C3+d4*C4+m)-k6*C6/(d1*C1+d2*C2+d3*C3+d4*C4+m)^b6)*Macoh-Finl*Cin2;
f3=Fout*C3+(k3*C3/(d1*C1+d2*C2+d3*C3+d4*C4+m)-k2*C2/(d1*C1+d2*C2+d3*C3+d4*C4+m)^b2)*Macoh-Finl*Cin3;
f4=Fout*C4+(k4*C4/(d1*C1+d2*C2+d3*C3+d4*C4+m)^b4-k3*C3/(d1*C1+d2*C2+d3*C3+d4*C4+m))*Macoh-Finl*Cin4;
f5=Fout*C5+(k6*C6/(d1*C1+d2*C2+d3*C3+d4*C4+m)^b6-k5*C5/(d1*C1+d2*C2+d3*C3+d4*C4+m))*Macoh-Finl*Cin5;
f6=Fout+G-Finl-L;
f7=Fout*Cinh2o+G*Cgh2o-Finl*Cinh2o-a*(Finl*C1-Fout*C1)+b*Fout*C4+Fout*C3+c*Fout*C2+d*Fout*C6-L*Clh2o;
f8=Finl*(t-tol)*(Cinacoh*cpacoh+Cinh2o*cph2o+Cinpx*cppx)+L*(t-tl)*(Clacoh*cpacoh+Clh2o*cph2o)+Fing*cpair*...
(t-tog)+G*(Cgacoh*Hvacoh+Cgh2o*Hvh2o)-(Finl*Cinpx-Fout*(Cinh2o+Cinacoh+Cinpx))*Hr;
f=[f1;f2;f3;f4;f5;f6;f7;f8];
y =
1.0e+003 *
Columns 1 through 4
0.000000006802276 0.000000006802276 0.000000006802276 0.000000006802276
-0.000000135512226 -0.000000135512226 -0.000000135512226 -0.000000135512226
-0.000000016422428 -0.000000016422428 -0.000000016422428 -0.000000016422428
-0.000000192485189 -0.000000192485189 -0.000000192485189 -0.000000192485189
-0.000000409983711 -0.000000409983711 -0.000000409983711 -0.000000409983711
-0.010049532361599 -0.010049532361599 -0.010049532361599 -0.010049532361599
-1.311221203993331 -1.311221203993331 -1.311221203993331 -1.311221203993331
0.688778796006669 0.688778796006670 0.688778796006670 0.688778796006669
Columns 5 through 6
0.000000006802276 0.000000006802276
-0.000000135512226 -0.000000135512226
-0.000000016422428 -0.000000016422428
-0.000000192485189 -0.000000192485189
-0.000000409983711 -0.000000409983711
-0.010049532361599 -0.010049532361599
-1.311221203993331 -1.311221203993331
0.688778796006670 0.688778796006669 |
|
回复此楼
建立一个完整的m函数，而不是在命令窗口一行一行的输入