Page 394 - Handbook Of Multiphase Flow Assurance

Computer code (Makogon, 1994, 1997) 393

IF (Xij .LT.-xmaxh) Xij = Xij + xmax
endif
Rsq = Xij*Xij + Yij*Yij + Zij*Zij
Rij = DSQRT(Rsq)
Ecou = Co*Qi*atomdata(natom1,4)/Rij
C
IF (Rij .LE. Rcut) THEN
c set itype for SPC water hydrogen
itype=5
c if oxygen, set itype for SPC water oxygen
if(i.eq.3)itype=6
R6red = (RR(Itype,natom1)*RR(Itype,natom1)/Rsq)**3
Elj = DD(Itype,natom1)*R6red*(R6red - 2.D0)
C
C ... collect hydrogen-bond interaction, if applicable
C
IF(iatomtype(natom1,1).eq.3) then
IF (I.EQ.3) THEN
IF (Rsq .LT.Rsqoff) CALL HBOND(Xi,Yi,Rsq,Iw,natom1,Ehb)
IF(Ehb.ne.0.D0) nhbcount=nhbcount+1
ENDIF
ENDIF
C
ENDIF
Ener = Ener + Ecou + Elj + Ehb
200 CONTINUE
nhb(j)=nhb(j)+nhbcount
nhbcurr(j)=nhbcount
Es(j)=Es(j)+Ener
end do
end do
C count the total polymer-surface energy
Ener=0.D0
DO j=1,nbackbone
Ener=Ener+Es(j)
end do
RETURN
C
END
C
C-------------------------------------------------------------------
C
SUBROUTINE ENERGYP(Ener)
IMPLICIT REAL*8(A-H,O-Z)
integer a, ap