24小时热门版块排行榜

返回列表

木木time

新虫 (初入文坛)

应助: 0 (幼儿园)
帖子: 1
在线: 8.2小时
虫号: 3197635
注册: 2014-05-11

[交流] 用1stopt拟合动力学方程参数无法运行，新手不懂如何调试，请求各位高手帮忙，谢谢！已有2人参与

程序是利用龙格库塔法模拟动力学方程的相关参数，我用低版本的软件运行的时候一直报错，一直无法运行，提示我的编写有问题，菜鸟不太懂，希望各位大神帮忙指正，谢谢啦。
下面是程序代码：
Title "cyclohexanone";
Parameters k1,E1,k2,E2,kA,EA,kH2,EH2,kCNE,ECNE,kCOL,ECOL,DEN;
Variable V0,PP,TT,yAin,yH2in,yArin,yCNEin,yCOLin;
Variable yCNE[Output],yCOL[Output];
StartProgram [VB];
Dim iter As Integer
Dim N0 As Double
Dim RK_i As Integer
Dim h As Double
Dim i As Integer
dim j As Integer
dim c As Double
dim DM As Double
dim RCNE As Double
dim RCOL As Double
dim yAr As Double
dim yA As Double
dim yCNE As Double
dim yH2 As Double
dim yCOL As Double
Dim KT1 As Double
dim KT2 As Double
dim KA As Double
dim KH2 As Double
dim KCNE As Double
dim KCOL As Double
Dim T As Double
Dim BB As Double
Dim pH2 As Double
dim pCOL As Double
dim pCNE As Double
dim pA As Double
Dim Y_Y01 As Double
Dim Y_Y02 As Double
Dim Y_Y03 As Double
dim y_y1 As Double
dim y_y2 As Double
dim y_y3 As Double
Dim x_x1 As Double
Dim x_x2 As Double
Dim x_x3 As Double
dim d_d1 As Double
dim d_d2 As Double
dim d_d3 As Double
dim e_e1 As Double
dim e_e2 As Double
dim e_e3 As Double
dim e_e4 As Double
dim e_e5 As Double
for iter=0 to Datalength-1
N0=V0(iter)
h=0.35/20
Y_Y01=yCNEin(iter)
Y_Y02=yCOLin(iter)
Y_Y03=0
For RK_i=1 to 20
      e_e1=0.5
      e_e2=0.5
      e_e3=1
      e_e4=1
      e_e5=0.5
      x_x1=Y_Y01
      x_x2=Y_Y02
      x_x3=Y_Y03
      y_y1=Y_Y01
      y_y2=Y_Y02
      y_y3=Y_Y03
      For j=1 To 4
         yCNE=x_x1
         yCOL=x_x2

         BB=(1+2*yA)/(1+2*yAin(iter))
         yA=BB*(yAin(iter)+yCNEin(iter))-yCNE
         yH2=BB*(yH2in(iter)+yCOLin(iter)+2yCNEin(iter))-2yCNE-yCOL
         yAr=BB*yArin(iter)
         T=TT(iter)+273.15

         pA=PP(iter)*yA
         pH2=PP(iter)*yH2
         pCNE=PP(iter)*yCNE
         pCOL=PP(iter)*yCOL

         KT1=k1*Exp(-E1/(8.315*T))
         KT2=k2*Exp(-E2/(8.315*T))
         KA=kA*Exp(-EA/(8.315*T))
         KH2=kH2*Exp(-EH2/(8.315*T))
         KCNE=kCNE*Exp(-ECNE/(8.315*T))
         KCOL=kCOL*Exp(-ECOL/(8.315*T))

         DEN=1+KA*PA+KH2*PH2+KCNE*PCNE+KCOL*PCOL
         RCNE=(KT1*KA*PA*(KH2*PH2)^2)/(DEN^3)
         RCOL=(KT2*KCNE*PCNE*KH2*PH2)/(DEN^2)
         DM=(1+2*yCNE)/(N0*(1+2*yCNEin(iter)))
         d_d1=DM*(1+2*yCNE)*RCNE
         d_d2=DM*(RCOL+2*yCOL*RCNE)
         d_d3=1

         c=h*d_d1
         if j=1 then
            x_x1=Y_Y01+e_e1*c
            y_y1=y_y1+e_e2*c/3
         elseif j=2 then
            x_x1=Y_Y01+e_e2*c
            y_y1=y_y1+e_e3*c/3
         elseif j=3 then
            x_x1=Y_Y01+e_e3*c
            y_y1=y_y1+e_e4*c/3
         elseif j=4 then
            x_x1=Y_Y01+e_e4*c
            y_y1=y_y1+e_e5*c/3
         end if
         c=h*d_d2
         if j=1 then
            x_x2=Y_Y02+e_e1*c
            y_y2=y_y2+e_e2*c/3
         elseif j=2 then
            x_x2=Y_Y02+e_e2*c
            y_y2=y_y2+e_e3*c/3
         elseif j=3 then
            x_x2=Y_Y02+e_e3*c
            y_y2=y_y2+e_e4*c/3
         elseif j=4 then
            x_x2=Y_Y02+e_e4*c
            y_y2=y_y2+e_e5*c/3
         end if
         c=h*d_d3
         if j=1 then
            x_x3=Y_Y03+e_e1*c
            y_y3=y_y3+e_e2*c/3
         elseif j=2 then
            x_x3=Y_Y03+e_e2*c
            y_y3=y_y3+e_e3*c/3
         elseif j=3 then
            x_x3=Y_Y03+e_e3*c
            y_y3=y_y3+e_e4*c/3
         elseif j=4 then
            x_x3=Y_Y03+e_e4*c
            y_y3=y_y3+e_e5*c/3
         end if
      Next
      Y_Y01=y_y1
      Y_Y02=y_y2
      Y_Y03=y_y3
Next
yCNE(iter)=y_y1
yCOL(iter)=y_y2
Next
EndProgram;
Data;
4.06 0.2 120 0.113 0.788 0.099 0.0001 0.0001 0.00796539 7.23E-05
2.03 0.2 120 0.114 0.787 0.099 0.0001 0.0001 0.018276944 0.000165986
1.02 0.2 120 0.117 0.784 0.099 0.0001 0.0001 0.027001782 0.000272745
1.02 0.2 120 0.115 0.786 0.099 0.0001 0.0001 0.028928857 0.000292211
0.51 0.2 120 0.114 0.787 0.099 0.0001 0.0001 0.05419921 0.000602822
1.66 0.2 141 0.133 0.771 0.096 0.0001 0.0001 0.024173136 0.000244173
4.07 0.2 140 0.115 0.787 0.098 0.0001 0.0001 0.018688673 0.000169726
2.04 0.2 140 0.117 0.785 0.098 0.0001 0.0001 0.041576064 0.00041996
1.03 0.2 140 0.127 0.775 0.098 0.0001 0.0001 0.077775178 0.001668313
1.03 0.2 140 0.12 0.781 0.099 0.0001 0.0001 0.067879227 0.000824444
0.52 0.2 141 0.125 0.777 0.098 0.0001 0.0001 0.124489331 0.002152918
4.46 0.2 162 0.122 0.776 0.102 0.0001 0.0001 0.011771229 9.49E-05
3.32 0.2 162 0.121 0.774 0.105 0.0001 0.0001 0.016481008 0.000149676
1.66 0.2 162 0.133 0.771 0.096 0.0001 0.0001 0.039389391 0.000397873
4.08 0.2 161 0.117 0.784 0.099 0.0001 0.0001 0.035950835 0.000326496
2.05 0.2 161 0.12 0.781 0.099 0.0001 0.0001 0.071190475 0.000791805
1.04 0.2 161 0.137 0.766 0.097 0.0001 0.0001 0.137103964 0.002229333
1.04 0.2 161 0.132 0.771 0.097 0.0001 0.0001 0.134541091 0.00218766
4.45 0.2 182 0.119 0.778 0.103 0.0001 0.0001 0.018492489 0.000149133
3.32 0.2 182 0.122 0.772 0.106 0.0001 0.0001 0.029740384 0.000270094
1.66 0.2 182 0.13 0.773 0.097 0.0001 0.0001 0.060335598 0.000671073
4.09 0.2 181 0.12 0.782 0.098 0.0001 0.0001 0.058478171 0.000590689
2.06 0.2 181 0.126 0.777 0.097 0.0001 0.0001 0.107592821 0.001527687
1.05 0.2 181 0.141 0.763 0.096 0.0001 0.0001 0.186241228 0.006156735
1.05 0.2 182 0.142 0.762 0.096 0.0001 0.0001 0.183768277 0.005488445
1.04 0.2 181 0.138 0.766 0.096 0.0001 0.0001 0.179658491 0.005747588
1.04 0.2 182 0.136 0.767 0.097 0.0001 0.0001 0.176490707 0.005458475

文件中是我需要模拟的动力学方程

回复此楼