|
|
aspen求助
| aspen版本: |
aspen7.2-8.8 |
这是我仿照论坛上的示例修改的子程序，请高手帮忙看一下，从语法上看看有没什么错误，多谢！编译后无法运行。
示例文件在这里http://meng.horse/forum.php?mod=viewthread&tid=76210
C User Kinetics Subroutine for Urea synthensis, 2NH3+CO2=CARB; CARB=UREA+H2O
C
SUBROUTINE USURA (SOUT, NSUBS, IDXSUB, ITYPE, NINT,
2 INT, NREAL, REAL, IDS, NPO,
3 NBOPST, NIWORK, IWORK, NWORK, WORK,
4 NC, NR, STOIC, RATES, FLUXM,
5 FLUXS, XCURR, NTCAT, RATCAT, NTSSAT,
6 RATSSA, KCALL, KFAIL, KFLASH, NCOMP,
7 IDX, Y, X, X1, X2,
8 NRALL, RATALL, NUSERV, USERV, NINTR,
9 INTR, NREALR, REALR, NIWR, IWR,
1 NWR, WR)
C
IMPLICIT NONE
C
C DECLARE VARIABLES USED IN DIMENSIONING
C
INTEGER NSUBS, NINT, NPO, NIWORK,NWORK,
+ NC, NR, NTCAT, NTSSAT,NCOMP,
+ NRALL, NUSERV,NINTR, NREALR,NIWR,
+ NWR
C
#include "pputl_ppglob.cmn"
#include "ppexec_user.cmn"
EQUIVALENCE (RMISS, USER_RUMISS)
EQUIVALENCE (IMISS, USER_IUMISS)
C
C
C
C
C.....RCSTR...
#include "rcst_rcstri.cmn"
#include "rxn_rcstrr.cmn"
C
C.....RPLUG...
#include "rplg_rplugi.cmn"
#include "rplg_rplugr.cmn"
EQUIVALENCE (XLEN, RPLUGR_UXLONG)
EQUIVALENCE (DIAM, RPLUGR_UDIAM)
C
C.....RBATCH...
#include "rbtc_rbati.cmn"
#include "rbtc_rbatr.cmn"
C
C.....PRES-RELIEF...
#include "rbtc_presrr.cmn"
C
C.....REACTOR (OR PRES-RELIEF VESSEL) PROPERTIES...
#include "rxn_rprops.cmn"
EQUIVALENCE (TEMP, RPROPS_UTEMP)
EQUIVALENCE (PRES, RPROPS_UPRES)
EQUIVALENCE (VFRAC, RPROPS_UVFRAC)
EQUIVALENCE (BETA, RPROPS_UBETA)
EQUIVALENCE (VVAP, RPROPS_UVVAP)
EQUIVALENCE (VLIQ, RPROPS_UVLIQ)
EQUIVALENCE (VLIQS, RPROPS_UVLIQS)
C
#include "shs_stwork.cmn"
EQUIVALENCE (MKBAS, STWORK_NDUM)
EQUIVALENCE (MKPHAS, STWORK_NBLM)
EQUIVALENCE (MTAPP, STWORK_NCOVAR)
EQUIVALENCE (MKBASS, STWORK_NWR)
EQUIVALENCE (MTAPPS, STWORK_NIWR)
EQUIVALENCE (SSALT, STWORK_RDUM1)
EQUIVALENCE (VSALT, STWORK_RDUM2)
EQUIVALENCE (FSALT, STWORK_FFSALT)
C
#include "dms_ncomp.cmn"
#include "dms_plex.cmn"
EQUIVALENCE (IB(1), B(1))
C
C DECLARE ARGUMENTS
C
INTEGER IDXSUB(NSUBS),ITYPE(NSUBS), INT(NINT),
+ IDS(2),NBOPST(6,NPO),IWORK(NIWORK),
+ IDX(NCOMP), INTR(NINTR), IWR(NIWR),
+ NREAL, KCALL, KFAIL, KFLASH,I
REAL*8 SOUT(1), WORK(NWORK),
+ STOIC(NC,NSUBS,NR), RATES(5),
+ FLUXM(1), FLUXS(1), RATCAT(NTCAT),
+ RATSSA(NTSSAT), Y(NCOMP),
+ X(NCOMP), X1(NCOMP), X2(NCOMP)
REAL*8 RATALL(NRALL),USERV(NUSERV),
+ REALR(NREALR),WR(NWR), XCURR
C
C
C DECLARE LOCAL VARIABLES
C
INTEGER IMISS, MKBAS, MKPHAS,MTAPP, MKBASS,
+ MTAPPS,LMW
REAL*8 REAL(NREAL), XL(7), B(1), RMISS, XLEN,
+ DIAM, TEMP, PRES, VFRAC, BETA,
+ VVAP, VLIQ, VLIQS, SSALT, VSALT,
+ FSALT, RTEMP, DENLIQ,KKA, KKB,
+ AREA, RTA, RTB, SUMR
INTEGER XMW
#include "dms_ipoff1.cmn"
C
C INITIALIZE RATES
C
C
C STATEMENT FUNCTIONS FOLLOW
C
C
C BEGIN EXECUTABLE CODE
C
XMW(I) = LMW + I
C
C SET PLEX OFFSETS
C
LMW = IPOFF1_IPOFF1(306)
DO 100 I = 1, NC
RATES(I) = 0D0
XL(I) = 0D0
100 CONTINUE
C
DO 101 I = 1, NCOMP
XL(IDX(I)) = X(I)
101 CONTINUE
C
C The structure in the array SOUT is as follows:
C
C SOUT(1) - SOUT(NCC) : Component flowrates(kg-moles/sec)
C SOUT(NCC+1) : Total flowrates(kg-moles/sec)
C SOUT(NCC+2) : Temperature(K)
C SOUT(NCC+3) : Pressure(N/SQM)
C SOUT(NCC+4) : Mass enthalpy(J/KG)
C SOUT(NCC+5) : Molar vapor fraction
C SOUT(NCC+6) : Molar liquid fraction
C SOUT(NCC+7) : Mass entropy(J/KG-K)
C SOUT(NCC+8) : Mass density(KG/CUM)
C SOUT(NCC+9) : Molecular Weight
C
C Details of this infomation can be found in the User Guide
C on page 17-10 for 8.5-6 or Reference Manual Vol. 6, User Models
C for Release 9.
C
C Set Reactor Temperature
C
RTEMP = SOUT(NCOMP_NCC+2)
C
C Get liquid mixture density
C
C DENLIQ = 1.0 / STWORK_VL
C
C Compute Concentration for reactants:
C
C CONC(i) = X(i) * RHO
C
C CONAA = XL(1) * DENLIQ
C CONACT = XL(2) * DENLIQ
C
C
C
C CALCULATE RATES OF REACTIONS
C
C The rates must be set for every single component exists in the
C system(both reactants and products) according to stoichiometry.
C
C RATES - The rates for all components
C
C Acetone + Allyl Alcohol -> N-Propyl Propionate
C XL1=UREA,XL2=CARB,XL3=NH3,XL4=CO2,XL5=H2O
C FTERM = REALR(1) * EXP( -REALR(2) / PPGLOB_RGAS / RTEMP )
C STERM = conaa * conact**0.5
C RATE = FTERM * STERM * RCSTRR_VOLRC
C RATES(1) = -rate
C RATES(2) = -rate
C RATES(3) = rate
C AREA=DFLOAT(RPLUGI_NTUBE)*3.141592654/4.*DIAM**2
KKA=XL(2)/XL(4)/XL(4)/XL(3)
KKB=XL(1)*XL(5)/XL(2)
RTA=1.0e12*EXP(-1E07/PPGLOB_RGAS/RTEMP)*(XL(4)**2*XL(3)-XL(2)/KKA)
RTB=1.5e09*EXP(-1E08/PPGLOB_RGAS/RTEMP)*(XL(2)-XL(1)*XL(5)/KKB)
RATES(1)=AREA*RTB
RATES(2)=AREA*(RTA-RTB)
RATES(3)=AREA*RTA
RATES(4)=-2*RTA*AREA
RATES(5)=RTB*AREA
C Checking Mass Balance for RATES Vector
C
SUMR = 1.0
DO 102 I=1,NCOMP_NCC
SUMR = SUMR + RATES(I)*B(XMW(I))
102 CONTINUE
IF (DABS(1.0 - SUMR) .GT. 1E-6) THEN
WRITE(USER_NHSTRY,*)
+ 'Calculated Rates Not in Mass Balance!'
WRITE(USER_NHSTRY,*)
+ 'Errors in User Kinetics Subroutines!'
ENDIF
write(USER_NHSTRY,*) 'volrc = ', RCSTRR_VOLRC
write(USER_NHSTRY,*) 'vfrrc = ', RCSTRR_VFRRC
RETURN
#undef P_NPOFF1
END
|
|
发表于 1970-1-1 08:00:00
显示全部楼层
IP卡
狗仔卡
提升卡
置顶卡
喧嚣卡
变色卡
千斤顶
显身卡
楼主