Trying to compile the gromos96 package ...

Hi

I am trying to compile the Gromos96 (molecular dynamics simulation package) code on an OS X v. 10.4.11 dual core intel Xeon.
I am using the next part of the Makefile supplied with the package

#CC=cc

#this for SGI machines
#FC = f77
#FFLAGS = -r8 -g
#FFLAGS = -O3
#FFLAGS = -OPT:roundoff=3:IEEE_arithmetic=3:fast_sqrt=OFF -lfastm
#for compiling a parallel version on the Power Challenge
#MPFLAGS = -mp
#PCFLAGS = -mp -DMP -O3 -SWP:body_ins_count_max=0 -TENV:X=3

and substituted by

CC = pgcc
FC = pgf77
FFLAGS = -fast -mp
FFLAGS = -Munroll -O3 -r8 -Mnoframe -Mdclchk -byteswapio -Mcache_align
#FFLAGS = -Munroll -O3 -r8 -Mnoframe -Mdclchk -byteswapio -Mcache_align
#FFLAGS = -OPT:roundoff=3:IEEE_arithmetic=3:fast_sqrt=OFF -lfastm

for compiling a parallel version on the Max OS X Server intel or

Mac Pro dual core intel xeon

MPFLAGS = -mp
PCFLAGS = -mp -DMP

for compiling a parallel version on the Max OS X Server intel or Mac Pro dual core intel xeon

When I tried to obtain the executable, I got next errors.

make promd_pc
pgf77 -Munroll -O3 -r8 -Mnoframe -Mdclchk -byteswapio -Mcache_align -c -o promd.o promd.f

.
.
.
pgf77 -Munroll -O3 -r8 -Mnoframe -Mdclchk -byteswapio -Mcache_align -c -o cobond.o cobond.f
pgf77 -Munroll -O3 -r8 -Mnoframe -Mdclchk -byteswapio -Mcache_align -c -o le.o le.f
pgf77 -Munroll -O3 -r8 -Mnoframe -Mdclchk -byteswapio -Mcache_align -mp -DMP -c -o nbpmlpc.o nbpmlpc.F
pgf77 -Munroll -O3 -r8 -Mnoframe -Mdclchk -byteswapio -Mcache_align -mp -DMP -c -o nonbmlpc.o nonbmlpc.F
pgf77 -Munroll -O3 -r8 -Mnoframe -Mdclchk -byteswapio -Mcache_align -mp -o …//promd_pc.64 ./promd.o ./mdutils.o ./rdmd.o ./fileio.o ./rdtopo.o ./posio.o ./blockio.o ./runmd.o ./runem.o ./shake.o ./traco.o ./cenmas.o ./concg.o ./conat.o ./force.o ./string.o ./rdpert.o ./gauss.o ./random.o ./restj.o ./restx.o ./disre.o ./fric.o ./nonbpl.o ./fourha.o ./box.o ./cobond.o ./le.o ./nbpmlpc.o ./nonbmlpc.o
Undefined symbols:
omp_num_threads, referenced from:
nbpml in nbpmlpc.o
nonbml in nonbmlpc.o
ld64-62.1 failed: symbol(s) not found
make: *** […//promd_pc.64] Error 2

Could give me any clues to compile with and without mp directives?

Many thanks in advance,

J. Cabrera

Hi,

Can you post a source code in npbmlpc.? where it has reference to omp_num_threads? Normally OMP_NUM_THREADS is an environment variable for openmp program.

Perhaps it is just a configuration issue but we don’t have source here.

Hongyon

Hi,

The source code is
::::::::::::
nbpmlpc.F:
::::::::::::
C $Id: nbpmlpc_to_pp.F,v 1.2 1996/11/22 12:53:08 wscott Exp $ --fortran--
C
COMMSUBR NBPMLPC
CCCCCCC P.HUNENBERGER, ZUERICH, MAY. 1996 CCCCCCCCCCCCCCCCCCCCC
C SGI Power-Challenge version - parallelized, optimized for software
C pipelining
C
C SUBROUTINE NBPML(NATTOT,NPM,NSM,NRAGTO,NCALCD,
C $ XCOORD,XCG,XR,FLR,
C $ NSPM,NSP,VIRLR,
C $ INB,NSZPL,JNB,RCUTP,RCUTL,
C $ MAXNRE,NRELKP,IAGRP,
C $ RCRF,EPSRF,APPAK,
C $ ELGLEL,ELGLRF,ELGLRC,ELGLLJ,
C $ EL34EL,EL34RF,EL34RC,EL34LJ,
C $ NUNRE2,ELREL,ELRRF,ELRRC,ELRLJ,
C $ LDOPER,L3D4D,RLAM,RMU,ALPHLJ,ALPHC,NLAM,MMU,
C $ LPION,LEVERY,LFREE)
C
C NBPML makes a list of non-bonded pairs using a cut-off criterion
C for the centres of the geometry of solute charge groups
C and the first atoms of solvent molecules, by scanning all
C possible pairs (using 3D distances only). The cut-off radius
C for mating the list is L.
C The pair list contains ordered sequence numbers of solute charge atom
C groups and solvent molecules (NI<NJ).
C Periodic images can be included in the list. In this case, the
C atoms of a charge group (solute) must lie within BOX/2 or
C BOXSQRT(3)/4 of each other.
C
C If L > L, a potential is evaluated for the twin-range
C region in subroutine L during the pair list construction.
C Otherwise the pair list is constructed without long range forces
C in subroutine L. The term long-range in contrast to GROMOS87
C includes Lennard-Jones interactions.
C
C The potential in L reads:
C V(r_ij) = - C6/(r_ij6) + C12/(r_ij12)
C + q_i
q_j/r_ij - (0.5C1/RCRF**3)q_iq_j(r_ij2)
C - ((1-0.5C1)/RCRF)q_iq_j
C
C where RCRF is the radius of the reaction field boundary (in principle
C RCUTL) and
C
C C1 = ( (2-2
EPSRF)(1+APPAKRCRF) - EPSRF*(APPAKRCRF)**2 ) /
C ( (1+2
EPSRF)(1+APPAKRCRF) + EPSRF*(APPAKRCRF)**2 )
C
C is the coefficient for the reaction-field correction beyond RCRF
C (with EPSRF.GE.1.0 and APPAK.GE.0.0).
C If EPSRF.LE.0.0, C1 is set to -1, corresponding to an infinite
C continuum permittivity ( a warning will be printed at first call )
C
C If LDOPER.OR.L3D4D (NTG.NE.L), a soft-core type,
C L and L dependent perturbation potential is applied
C ( see manual )
C It is evaluated for atom pairs involving the atoms II for which
C IPERT(II) .NE. L .
C The alpha coefficients for the soft-core perturbation are set to zero
C except for atom pairs involving the atoms II for which
C ISCLJ(II) .NE. 0 in which case alpha(LJ) = L
C and/or
C ISCCB(II) .NE. 0 in which case alpha(Coulomb/rx field) = L
C
C When LDOPER.AND…NOT.L3D4D (NTG.EQ.L), solely the lambda
C derivative of the potential will be returned.
C When .NOT.LDOPER.AND.L3D4D (NTG.EQ.L), solely the mu
C derivative of the potential will be returned.
C When LDOPER.AND.L3D4D (NTG.EQ.L), both will be returned.
C Derivatives not explicitely required will be returned as zero.
C
C NCALCD LDOPER/L3D4D calculated
C
C 3 .F./.F. Normal interactions in 3D
C 3 .T./.F. Normal interactions in 3D
C For pairs involving perturbed atoms:
C Soft-core lambda dependent potential in 3D,
C mu has to be zero (if not, set locally to zero)
C 3 .F./.T. - forbidden -
C 3 .T./.T. - forbidden -
C 4 .F./.F. Normal interactions in 4D
C 4 .T./.F. Normal interactions in 4D
C For pairs involving perturbed atoms:
C Soft-core lambda dependent potential,
C 3D-4D mixed at given rmu, lam deriv. returned.
C 4 .F./.T. Normal interactions in 4D
C Soft-core lambda dependent potential,
C 3D-4D mixed at given rmu, mu deriv. returned.
C 4 .T./.T. Normal interactions in 4D
C Soft-core lambda dependent potential,
C 3D-4D mixed at given rmu, both deriv. returned.
C >4 .F./.F. In principle allowed
C
C IMPORTANT NOTE:
C
C With NCALCD.EQ.4, no check is done whether an
C atom is in 4D or not according to 4D specification list (C4D array).
C This means that all atoms (solute/solvent) which are to be handled
C in 3D MUST HAVE ZERO 4thD COORDINATES !
C
C This, however, won’t prevent them from having a 4thD component
C in the force acting on them, due to interaction with 4D atoms.
C
C If ALPHLJ.NE.0.0 and for an atom II with ISCLJ(II).NE.0, the algorithm
C assumes that (either in state A or B), if C6.EQ.0.0, then C12.EQ.0.0.
C If this is not the case, a non soft-core perturbation of C12 is done.
C A warning will be printed at first call in this case.
C
C The forces are calculated in NCALCD dimensions (L.LE.L).
C The reaction-field correction is, however, always calculated in 3D.
C If L.EQ.3, whatever L, a local value of zero will be used.
C A warning will be printed at first call in this case. Additionally,
C L3D4D (NTG.EQ.L.OR.NTG.EQ.L) will not be
C The cut-off criterion is always applied in 3D.
C
C Periodic boundary conditions can be applied, depending on L,
C where the periodic box can be a truncated octahedron, rectangular
C or monoclinic.
C
C The periodic box can be a truncated octahedron, rectangular or
C monoclinic, depending on L.
C If the long range virial is to be calculated
C (ABS(L .EQ. L), it must be specified whether the
C solute consists of separate submolecules or not.
C This can be done by setting the NSPM > 1 and setting the corresponding
C values into the NSP array. The virial is only calculated
C in 3D.
C
C
C NATTOT total number of atoms
C NPM number of (identical) solute molecules
C NRAGT total number of charge groups
C NRAGTO total number of charge groups
C NCALCD the number of dimensions to calculate distance and
C forces in. Must be 3 or 4.
C Note that the 3D distance between charge groups
C I and J is always used for deciding whether to put
C the pair into the list or not, independently of NCALCD.
C The distance used in the force calculation, however,
C does depend on NCALCD.
C XCOORD(NDIM
NATTOT)
C the atom cartesian coordinates
C XCG(NDIMNRAGTO)
C the geometric centre of solute charge groups and the
C first atoms of solvent molecules.
C XR(NDIM
NATTOT)
C atom cartesian coordinates relative to submolecular
C centres of mass. This array is only accessed if
C the virial is calculated, i.e. ABS(NTB) = L.
C FLR(NDIMNATTOT)
C delivered with the long range coulomb forces if
C RCUTL > RCUTP.
C NSPM If the virial is to calculated, then L and L are
C accessed in order to calculate the molecular virial.
C NSPM is the number of separate (sub)molecules forming
C one solute molecule
C NSP(NSPM+1)
C atom sequence numbers of the last atoms of the
C submolecules.
C only elements 1…NSPM are considered to be filled
C with data. We need NSPM+1 elements in order to avoid
C a double IF statement when checking for submolecules
C VIRLR(NDIM)
C delivered with the X-, Y-, and Z- components of the
C long range virial if L = L.
C INB(NRAGTO)
C delivered with the pointer list;
C INB(I) specifies the number of charge groups or
C solvent molecules J with J>I forming a pair with charge
C group or solvent molecule I
C NSZPL delivered with the number of elements in JNB
C JNB(MAXJNB)
C delivered with the pointer list;
C it specifies the solute charge groups or solvent
C molecules J forming a pair with solute charge group
C or solvent molecule I in ascending order; for each I the
C pointer in JNB is incremented by INB(I)
C MAXNRE maximum number of energy groups
C NRELKP(MAXNRE,MAXNRE)
C index of pairs of energy groups in energy arrays
C IAGRP(MAXNAT)
C for each atom, the energy type code its
C contribution is added to.
C RCRF radius of the reaction field boundary (in principle
C RCUTL) and
C EPSRF dielectric permittivity of the reaction field
C continuum outside the boundary
C APPAK inverse Debye screening length
C NUNRE2 number of energy froups
C ELREL(MXNRE2)
C delivered with long range coulomb energies.
C ELRRF(MXNRE2)
C delivered with the long range distance dependent
C Poisson-Boltzmann reaction field energies
C ELRRC(MXNRE2)
C delivered with the long range distance independent
C Poisson-Boltzmann reaction field energies
C ELGLEL delivered with dV/d(L) at fix L
C for the long range Coulomb term if LDOPER = .TRUE.
C ELGLRF delivered with with dV/d(L) at fix L
C for the distance dependent long range
C Poisson-Boltzmann term if LDOPER = .TRUE.
C ELGLRC delivered with with dV/d(L) at fix L
C for the distance independent long range
C Poisson-Boltzmann term if LDOPER = .TRUE.
C EL34EL delivered with the 3D to 4D long range free energy contribution
C for the Coulomb term if LDOPER = .TRUE. .AND. NCALCD = 4
C EL34RF delivered with the 3D to 4D long range free energy contribution
C for the long range distance dependent
C Poisson-Boltzmann reaction field
C term if LDOPER = .TRUE. .AND. NCALCD = 4
C EL34RC delivered with the 3D to 4D long range free energy contribution
C for the long range distance independent
C Poisson-Boltzmann reaction field
C term if LDOPER = .TRUE. .AND. NCALCD = 4
C LDOPER perturbation is applied, lambda derivative computed
C LDO34D perturbation is applied, mu derivative computed
C RLAM actual value of the lambda parameter
C RMU actual value of the mu parameter
C ALPHLJ soft-core alpha parameter for the Lennard-Jones component
C ALPHC soft-core alpha parameter for the electrostatics component
C NLAM exponent of lambda in the Hamiltonian coupling scheme
C MMU exponent of mu in the Hamiltonian coupling scheme
C LPION a path-integral calculation is performed
C LEVERY
C .FALSE. a number of quantities depending on L, L,
C L, L, L, L, L, L,
C L and L are only calculated at the
C first call
C .TRUE. they are calculated at every subr. call.
C LFREE
C .TRUE. the pair list and interaction is calculated solely
C for perturbed atoms. This is used in free energy
C extrapolation.
C
C NOT passed as arguments, but accessed in common blocks are
C the following:
C
C RCUTP cut-off radius for the pair list; since no check on
C neighbour exclusions is performed, and since only charge
C group pairs I<J are scanned, L should be chosen
C large enough when L > L.
C RCUTL cut-off radius for the long range coulomb force;
C The following conditions must hold, depending on the
C boundary conditions specified.
C If L <0, L, L < BOX(M)SQRT(3)/4
C If L >0, L, L < BOX(M)/2
C If L is larger, only the nearest image is calculated.
C
C If L .EQ. .TRUE., then L, L, L, L
C and L
C are accessed in order to calculate a perturbation dependent
C potential.
C IPERT(MAXNRP)
C If IPERT(II) =0, the charge of atom II is not perturbed.
C otherwise it denotes where the new charge can be found
C in L.
C CGB(MAXPAT)
C charges determining the perturbation potential for
C state B.
C
COMMEND
SUBROUTINE NBPML(NATTOT,NPM,NSM,NRAGTO,NCALCD,
$ XCOORD,XCG,XR,FLR,
$ NSPM,NSP,VIRLR,
$ INB,NSZPL,JNB,RCUTP,RCUTL,
$ MAXNRE,NRELKP,IAGRP,
$ RCRF,EPSRF,APPAK,
$ ELGLEL,ELGLRF,ELGLRC,ELGLLJ,
$ EL34EL,EL34RF,EL34RC,EL34LJ,
$ NUNRE2,ELREL,ELRRF,ELRRC,ELRLJ,
$ LDOPER,L3D4D,RLAM,RMU,ALPHLJ,ALPHC,NLAM,MMU,
$ LPION,LEVERY,LFREE)
C
C INCLUDES
C
INCLUDE ‘toposz.h’
INCLUDE ‘topoar.h’
INCLUDE ‘coordsz.h’
INCLUDE ‘box.h’
INCLUDE ‘forcesz.h’
INCLUDE ‘pertsz.h’
INCLUDE ‘pertar.h’
C
C ARGUMENTS
C
LOGICAL LDOPER,L3D4D,LPION,LEVERY,LFREE
INTEGER NATTOT,NRAGTO,NPM,NSM,NCALCD,NSZPL
INTEGER INB(NRAGTO),JNB(MAXJNB)
INTEGER NSPM,NSP(NSPM+1)
INTEGER MAXNRE,NRELKP(MAXNRE,MAXNRE)
INTEGER NLAM,MMU
INTEGER NUNRE2
REAL RCUTP,RCUTL,RCRF,EPSRF,APPAK
REAL XCOORD(NDIM
NATTOT),XCG(NDIM
NRAGTO)
REAL XR(NDIMNATTOT),FLR(NDIMNATTOT),VIRLR(NDIM)
REAL ELREL(NUNRE2),ELGLEL,EL34EL
REAL ELRRF(NUNRE2),ELGLRF,EL34RF
REAL ELRRC(NUNRE2),ELGLRC,EL34RC
REAL ELRLJ(NUNRE2),ELGLLJ,EL34LJ
REAL RLAM,RMU,ALPHLJ,ALPHC
CHARACTER IAGRP(NATTOT)
C
C LOCAL VARIABLES
C
LOGICAL LSTOP
INTEGER IFTOPA,ILTOPA,ITOPCG,ISUBM
INTEGER I,II,M
INTEGER ITMPA,ITMPB,JTYPE,INTIJA,INTIJB
INTEGER IP_T,ICPU,NI,IAUX,ITMP,NPJ_T
INTEGER NRAGT
REAL RMULOC
REAL RCUTP2,RCUTL2
CHARACTER PRGSTR*(5)
C
C
C+++ ALL CPU VARIABLES
C The following variables are common to all CPU’s, in NONBML and
C NBPML
C
INTEGER MAXCPU
PARAMETER (MAXCPU=32)
C
REAL ONE
PARAMETER (ONE = 1.0)
REAL ENEPS
PARAMETER (ENEPS = 1.E-20)
REAL ZEROC6
PARAMETER (ZEROC6 = 1.0E-20)
REAL ZEROAL
PARAMETER (ZEROAL = 1.0E-20)
C
LOGICAL LFULL4,LPART4,LFULL3,LPART3,LPERTL
LOGICAL LMONO,LRECT,LOCTO,LVAC,LDOVIR,LDOTRA,L4D,LNOSVT
LOGICAL LEMIS,L14MIS,LMXAT,LMXCG,LMXPL
DIMENSION LEMIS(MAXCPU),L14MIS(MAXCPU),
$ LMXAT(MAXCPU),LMXCG(MAXCPU),LMXPL(MAXCPU)
INTEGER ITSUBM(MAXCAG)
INTEGER NRAGP,NRPT,IOFF
INTEGER NCPU
REAL RLMA, RL2A, RLMB, RL2B, RLLA, RLLB, RLB, RLA, RDLB, RDLA
REAL RLARMU, RLBRMU, RM4D, RM3D, RDM3D, RDM4D
REAL RMCST,RMDER
REAL RFF, RFE, RFC, C1, RCRF2
REAL PININV
REAL COSB,COSB2,BOXOH,BOXOQ,PYE,CONV,BETAR
C
COMMON /NBRREL/ RLMA,RL2A,RLMB,RL2B,RLA,RLB,
$ RLLB,RLLA,RDLA,RDLB,RM3D,RM4D,RDM3D,RDM4D,RLARMU,RLBRMU,
$ RMCST,RMDER,
$ RFF, RFE, RFC, C1, RCRF2,
$ COSB,COSB2,BOXOH,BOXOQ,PYE,CONV,BETAR,
$ PININV
C
COMMON /NBRLOG/ LFULL3,LFULL4,LPART3,LPART4,LPERTL,
$ LMONO,LRECT,LOCTO,LVAC,LDOVIR,LDOTRA,L4D,LNOSVT,
$ LEMIS,L14MIS,LMXAT,LMXCG,LMXPL
C
COMMON /NBRINT/ NCPU,NRAGP,NRPT,IOFF,ITSUBM
C
C+++ DONE: ALL CPU VARIABLES
C
C
INTEGER NPCPU(MAXCPU)
C
C+++ We need temporary storage per processor, dynamically allocated
C For the force…
C
C 4 REALS FOR ELGLEL,ELGLRF,ELGLRC,ELGLLJ
C 4 REALS FOR EL34EL,EL34RF,EL34RC,EL34LJ
C 4 * NUNRE2 REALS FOR ELREL(),ELRRF(),ELRRC(),ELRLJ()
C 1 * NDIM REALS FOR VIR()
C 1 * NDIMNATTOT REALS FOR F()
C
POINTER (P_F_TMP,F_TMP)
REAL F_TMP
DIMENSION F_TMP (1)
C
C+++ … and for the temporary pair lists
C
POINTER (P_JNB_TMP,JNB_TMP)
INTEGER JNB_TMP
DIMENSION JNB_TMP(1)
C
COMMON /SGI_COM/MPROC,MATOM,MX_JNB,P_F_TMP,P_JNB_TMP
INTEGER MPROC,MATOM,MX_JNB
C
C SAVED AND DATA
C
LOGICAL LFIRST,LDOLNG,LMPROC
SAVE LFIRST,LDOLNG,LMPROC
DATA LFIRST/.TRUE./
DATA PRGSTR /‘NBPML’/
#ifdef MP
CJMCT-begin
INTEGER OMP_NUM_THREADS
INTEGER MALLOC
EXTERNAL MALLOC
CJMCT-end
EXTERNAL OMP_NUM_THREADS
DATA LMPROC /.TRUE./
#else
DATA LMPROC /.FALSE./
#endif
C
C------------BEGIN NBPML
C
#ifdef MP
C
IF (LMPROC.AND.(P_F_TMP.EQ.0 .OR. P_JNB_TMP.EQ.0)) THEN
MPROC = OMP_NUM_THREADS()
IF ( MPROC.GT.MAXCPU ) THEN
PRINT ,PRGSTR,’ INCREASE MAXCPU TO ‘,MPROC
STOP
ENDIF
WRITE(0,’(2A,I4,A)’) PRGSTR,’ RUNNING ‘,MPROC,’ PROCESSORS’
IF (MPROC.EQ.1) THEN
C
C+++ The unique CPU gets normal arrays
C
WRITE (0,’(3A)’) PRGSTR,’ NO LOCAL MEMORY ',
$ ‘ALLOCATION REQUIRED’
LMPROC = .FALSE.
ELSE
IF (P_F_TMP.EQ.0) THEN
C
C+++ allocate space for per processor arrays if required
C
MATOM = 4 + 4 + 4
NUNRE2 + NDIM + NDIM
NATTOT
P_F_TMP = MALLOC(8*(MPROC-1)MATOM)
IF (P_F_TMP.EQ.0) THEN
WRITE(0,
) PRGSTR,'F_TMP, ALLOCATION FAILED FOR ‘,
$ 8*(MPROC-1)MATOM,’ BYTES’
STOP
ENDIF
ENDIF
C
IF (P_JNB_TMP.EQ.0) THEN
C
C+++ allocate space for per processor pair list if required
C give every thread a secure (10%)SQRT(MPROC) more space
C
MX_JNB = MAXJNB/MPROC +
$ (MAXJNB/(10.0
MPROC))SQRT(1.0MPROC)
WRITE (0,’(A,I7,A,I7,A)’) 'MAXJNB: ‘,MAXJNB,
$ ’ ; ALLOCATING PER PROCESSOR: ‘,
$ MX_JNB,’ REALS’
P_JNB_TMP = MALLOC(4
MPROCMX_JNB)
IF (P_JNB_TMP.EQ.0) THEN
WRITE(0,’(3A,I8,A)’) PRGSTR,'JNB_TMP, ALLOCATION ',
$ 'FAILED FOR ',
$ 4
MPROCMX_JNB,’ BYTES’
STOP
ENDIF
WRITE (0,’(A,I8,A)’) 'ALLOCATED: ',4
MPROC*MX_JNB,
$ ’ BYTES FOR LOCAL PAIR LISTS’
ENDIF
C
ENDIF
C
ENDIF
C
#endif
C
LDOLNG = (RCUTL .GT. RCUTP)
NRAGT = NRAGTO
RCUTP2 = RCUTP
2
RCUTL2 = RCUTL2
C
C
C+++ PERFORM FIRST CALL SETTINGS AND CHECKING
C
IF (LEVERY .OR. LFIRST) THEN
C
C+++ calculate NTB/BOX/BETA dependent quantities
C
LVAC = (NTB .EQ. 0)
LMONO = (NTB .GT. 0)
LOCTO = (NTB. LT. 0)
LDOVIR= (ABS(NTB) .EQ. NTBVIR)
IF (LMONO) THEN
PYE = 4.E0ATAN(ONE)
CONV = 1.8E2/PYE
BETAR = BETA/CONV
COSB = COS(BETAR)
LDOTRA = (ABS(COSB).GE.1.E-4)
COSB2 = 2.0
COSB
ELSE
LDOTRA = .FALSE.
ENDIF
LRECT = (LMONO.AND…NOT.LDOTRA)
IF (LOCTO) THEN
BOXOH = BOXH(1)
BOXOQ = BOX(1)0.75E0
ELSE
BOXOH = 0.0
BOXOQ = 0.0
ENDIF
C
ENDIF
C
IF (LFIRST) THEN
C
C+++ calculate other quantities
C
LFIRST = .FALSE.
NCPU = MPROC
NRAGP = NPM
NCAG
NRPT = NPM*NRP
IOFF = 1 - NDIM
IF (NPID .NE. 0) THEN
PININV = ONE / NPID
ENDIF
LNOSVT = (NRAGP.EQ.NRAGT)
LPERTL = (LDOPER.OR.L3D4D)
L4D = (NCALCD .EQ. 4)
C
C+++ a few first call consistency checks
C IF LPERTL AND ALPHLJ.NE.0.0 CHECK THAT C6.EQ.0 => C12.EQ.0 FOR ATOMS II
C WITH ISCLJ(II).NE.0 AND ANY OTHER ATOM ( STATE A AND B )
C
IF (LPERTL .AND. ABS(ALPHLJ) .GT. ZEROAL) THEN
DO 26 II=1,NRP
IF (IPERT(II) .NE. 0) THEN
IF (ISCLJ(IPERT(II)) .NE. 0) THEN
ITMPA = IAC(II)
ITMPB = IACB(IPERT(II))
DO 28 JTYPE=1,NRATT
INTIJA = MPAC(ITMPA,JTYPE)
INTIJB = MPAC(ITMPB,JTYPE)
IF (C6(INTIJA) .LT. ZEROC6.AND.
$ C12(INTIJA) .GT. ZEROC6) THEN
PRINT 900,PRGSTR,II
PRINT 910,‘A’,JTYPE
ENDIF
IF (C6(INTIJB) .LT. ZEROC6.AND.
$ C12(INTIJB) .GT. ZEROC6) THEN
PRINT 900,PRGSTR,II
PRINT 910,‘B’,JTYPE
ENDIF
IF (CS6(INTIJA) .LT. ZEROC6.AND.
$ CS12(INTIJA) .GT. ZEROC6) THEN
PRINT 900,PRGSTR,II
PRINT 915,‘A’,JTYPE
ENDIF
IF (CS6(INTIJB) .LT. ZEROC6.AND.
$ CS12(INTIJB) .GT. ZEROC6) THEN
PRINT 900,PRGSTR,II
PRINT 915,‘B’,JTYPE
ENDIF
28 CONTINUE
ENDIF
ENDIF
26 CONTINUE
ENDIF
C
C+++ check that NCALCD.LE.NDIM
C
IF ( NCALCD.GT.NDIM ) THEN
PRINT *,PRGSTR,’ ERROR: NCALCD.GT.NDIM’
STOP
ENDIF
C
C+++ no mu derivative possible in pure 3D calculation
C
IF (.NOT. L4D .AND. L3D4D) THEN
PRINT *,PRGSTR,’: NCALCD.EQ.3.AND.’,
$ ‘(NTG.EQ.NTGMU.OR.NTG.EQ.NTGBOT) NOT ALLOWED’
ENDIF
C
C+++ for this PC routine, NCALCD>4 is forbidden
C
IF ( NCALCD.GT.4 ) THEN
PRINT *,PRGSTR,’ ERROR: NCALCD.GT.4 NOT ALLOWED IN THE ',
$ ‘PC VERSION’
STOP
ENDIF
C
C+++ for this PC routine, LPERTL is incompatible with LPION
C
IF ( LPERTL.AND.LPION ) THEN
PRINT ,PRGSTR,’ ERROR: PATH INTEGRAL AND PERTURBATION ',
$ ‘ARE INCOMPATIBLE’
STOP
ENDIF
C
C+++ make a local array with the submolecules of solute charge groups
C
IF (LDOVIR) THEN
IFTOPA = 1
DO 200 ITOPCG=1,NCAG
ILTOPA = INC(ITOPCG)
ISUBM = 1
IF ( NSPM.GT.1 ) THEN
250 IF (ISUBM .LE. NSPM .AND.
$ IFTOPA .GT. NSP(ISUBM)) THEN
ISUBM = ISUBM+1
GOTO 250
ELSEIF (ISUBM .GT. NSPM) THEN
PRINT ,PRGSTR,’ ERROR IN NSP() ARRAY’
STOP
ENDIF
ENDIF
ITSUBM(ITOPCG) = ISUBM
IFTOPA = ILTOPA + 1
200 CONTINUE
ENDIF
C
C+++ set reaction field variables
C
IF (EPSRF .LE. 0.0) THEN
C1 = -1.0
PRINT ,PRGSTR, ’ WARNING, EPSRF.LE.0 , PERMITTIVITY SET ',
$ ‘TO INFINITY’
ELSE
C1 = ((2.0 - 2.0 * EPSRF) * (1.0 + APPAK
RCRF)
$ - EPSRF * (APPAK * RCRF)**2)
$ /((1.0 + 2.0 * EPSRF) * (1.0 + APPAK
RCRF)
$ + EPSRF * (APPAK * RCRF)2)
ENDIF
RFF = C1 / RCRF
3
RFE = -0.5 * RFF
RFC = -(1.0-0.5
C1)/RCRF
RCRF2 = RCRF
2
C
900 FORMAT (1X,A6,’: ATOM’,I5,’ SOFT-CORE PERTURBED’)
910 FORMAT (1X,‘AND C6 .EQ. 0 BUT C12 .NE. 0 WITH STATE ‘,
$ A1,’ ATOM TYPE’,I5)
915 FORMAT (1X,‘AND CS6 .EQ. 0 BUT CS12 .NE. 0 WITH STATE ‘,
$ A1,’ ATOM TYPE’,I5)
C
C+++ DONE: PERFORM FIRST CALL SETTINGS AND CHECKING
C
CJMCT ENDIF
C-----I—
C
IF ((LDOLNG.OR.LFREE).AND.LPERTL) THEN
C
C+++ SET GENERAL PERTURBATION VARIABLES/LOGICALS
C this code is gone through if LPERTL
C 3D calculation should imply mu=0, otherwise, set local value to zero.
C
IF ( .NOT.L4D.AND.RMU.NE.0.0 ) THEN
PRINT ,PRGSTR,’ WARNING, NCALCD .EQ. 3 .AND. RMU .NE. 0.0’
PRINT ,’ LOCAL RMU WILL BE SET TO ZERO’
RMULOC = 0.0
ELSE
RMULOC = RMU
ENDIF
CJMCT
ENDIF
C
RLMA = RLAM2 + RMULOC2
RL2A = RLAM2
RLMB = (1.0 - RLAM)2 + RMULOC2
RL2B = (1.0 - RLAM)2
C
IF (RLAM.EQ.0.0 .AND. NLAM.EQ.1) THEN
RLB = 0.0
RDLB = NLAM
ELSE
RLB = RLAM
NLAM
RDLB = RLAM
(NLAM - 1) * NLAM
ENDIF
IF (RLAM.EQ.1.0 .AND. NLAM.EQ.1) THEN
RLA = 0.0
RDLA = -NLAM
ELSE
RLA = (1.0 - RLAM)NLAM
RDLA = - (1.0 - RLAM)
(NLAM - 1) * NLAM
ENDIF
C
RLLB = -(1.0 - RLAM) * RLB
RLLA = RLAM * RLA
RLARMU = RLA * RMULOC
RLBRMU = RLB * RMULOC
C
IF (RMULOC.EQ.0.0 .AND. MMU.EQ.1) THEN
RM4D = 0.0
RDM4D = MMU
ELSE
RM4D = RMULOCMMU
RDM4D = RMULOC
(MMU - 1) * MMU
ENDIF
IF (RMULOC.EQ.1.0 .AND. MMU.EQ.1) THEN
RM3D = 0.0
RDM3D = -MMU
ELSE
RM3D = (1.0 - RMULOC)MMU
RDM3D = - (1.0 - RMULOC)
(MMU - 1) * MMU
ENDIF
C
RMCST = RM3D + RM4D
RMDER = RDM3D + RDM4D
C
IF (.NOT. L4D) THEN
RDM3D = 0.0
RMCST = 1.0
RMDER = 0.0
ENDIF
C
LFULL4 = (L4D .AND. (RMULOC .NE. 0.0))
LPART4 = (RMULOC .EQ. 0.0 .AND. L3D4D .AND. MMU.EQ.1)
LFULL3 = ((L4D .AND. (RMULOC .NE. 1.0)).OR. (.NOT. L4D))
LPART3 = (RMULOC .EQ. 1.0 .AND. L3D4D .AND. MMU.EQ.1)
C
C+++ DONE: SET GENERAL PERTURBATION VARIABLES/LOGICALS
C
C-------I—
ENDIF
C
IF (LDOTRA) THEN
CALL TRACO(NATTOT,0,XCOORD,BETA,1,LEVERY)
ENDIF
C
C+++ set arrays to zero
C
ELGLEL = 0.0
ELGLRF = 0.0
ELGLRC = 0.0
ELGLLJ = 0.0
EL34EL = 0.0
EL34RF = 0.0
EL34RC = 0.0
EL34LJ = 0.0
C
DO 10 I=1,NDIM * NATTOT
FLR(I) = 0.0
10 CONTINUE
C
DO 20 I=1,NDIM
VIRLR(I) = 0.0
20 CONTINUE
C
DO 30 I=1,NUNRE2
ELREL(I) = 0.0
ELRRF(I) = 0.0
ELRRC(I) = 0.0
ELRLJ(I) = 0.0
30 CONTINUE
C
C+++ calculate centers of geometry of charge groups
C
CALL CLCGEO(NATTOT,NRAGTO,NPM,NSM,NDIM,XCOORD,XCG)
C
C+++ now make the pair list and possibly calculate interaction
C
IF (LDOLNG .OR. LFREE) THEN
C
C+++ Now the MP stuff (long-range)
C
IF (.NOT.LMPROC) THEN
C
C+++ the unique processor gets the normal pair list and interaction
C arrays
C
CALL NBWITH_MP(NATTOT,NPM,NSM,NRAGTO,NCALCD,
$ XCG,XCOORD,XR,FLR,NSPM,NSP,VIRLR,
$ RCUTP2,RCUTL2,INB,JNB,NSZPL,MAXNRE,NRELKP,IAGRP,
$ ELGLEL,ELGLRF,ELGLRC,ELGLLJ,
$ EL34EL,EL34RF,EL34RC,EL34LJ,
$ NUNRE2,ELREL,ELRRF,ELRRC,ELRLJ,
$ ALPHLJ,ALPHC,LPION,LEVERY,LFREE,MAXJNB,1)
C
ELSE
C
CJMCT CSCO REMOVED JNB_TMP,F_TMP,
C$DOACROSS
C$& LOCAL(ICPU,NPJ_T,IP_T,I),
C$& SHARE(NATTOT,NRAGTO,NPM,NSM,NCALCD,
C$& XCG,XCOORD,XR,FLR,NSPM,NSP,VIRLR,
C$& RCUTP2,RCUTL2,INB,JNB,NSZPL,MX_JNB,
C$& MAXNRE,NRELKP,IAGRP,
C$& ELGLEL,ELGLRF,ELGLRC,ELGLLJ,
C$& EL34EL,EL34RF,EL34RC,EL34LJ,
C$& NUNRE2,ELREL,ELRRF,ELRRC,ELRLJ,
C$& ALPHLJ,ALPHC,LPION,LEVERY,LFREE,
C$& MATOM,NDIM,NCPU,NPCPU)
C
DO 1500 ICPU = 1 , NCPU
C
NPJ_T = MX_JNB
(ICPU-1)
C
IF (ICPU.EQ.1) THEN
C
C+++ the first (not unique) processor gets the normal arrays for forces,
C energies and virial, and a temporary array for the pair list
C
CALL NBWITH_MP(NATTOT,NPM,NSM,NRAGTO,NCALCD,
$ XCG,XCOORD,XR,FLR,
$ NSPM,NSP,VIRLR,
$ RCUTP2,RCUTL2,INB,JNB_TMP(NPJ_T+1),NPCPU(ICPU),
$ MAXNRE,NRELKP,IAGRP,
$ ELGLEL,ELGLRF,ELGLRC,ELGLLJ,
$ EL34EL,EL34RF,EL34RC,EL34LJ,
$ NUNRE2,ELREL,ELRRF,ELRRC,ELRLJ,
$ ALPHLJ,ALPHC,LPION,LEVERY,LFREE,MX_JNB,1)
C
ELSE
C
C+++ all the the others get temporary arrays for the pair list
C
IP_T = MATOM
(ICPU-2)
DO 890 I=1,MATOM
F_TMP(I+IP_T) = 0.0
890 CONTINUE
CALL NBWITH_MP(NATTOT,NPM,NSM,NRAGTO,NCALCD,
$ XCG,XCOORD,XR,
$ F_TMP(IP_T + 4NUNRE2 + NDIM + 9 ),NSPM,NSP,
$ F_TMP(IP_T + 4
NUNRE2 + 9),
$ RCUTP2,RCUTL2,INB,JNB_TMP(NPJ_T+1),NPCPU(ICPU),
$ MAXNRE,NRELKP,IAGRP,
$ F_TMP(IP_T + 1),F_TMP(IP_T + 2),
$ F_TMP(IP_T + 3),F_TMP(IP_T + 4),
$ F_TMP(IP_T + 5),F_TMP(IP_T + 6),
$ F_TMP(IP_T + 7),F_TMP(IP_T + 8),
$ NUNRE2,F_TMP(IP_T + 9),F_TMP(IP_T + NUNRE2 + 9),
$ F_TMP(IP_T + 2NUNRE2 + 9),
$ F_TMP(IP_T + 3
NUNRE2 + 9),
$ ALPHLJ,ALPHC,LPION,LEVERY,LFREE,MX_JNB,ICPU)
C
ENDIF
C
1500 CONTINUE
C
C+++ now collect energies and virial from all the processors
C
DO 805 ICPU = 2, MPROC
IP_T = MATOM*(ICPU-2)

ELGLEL = ELGLEL + F_TMP(IP_T + 1)
ELGLRF = ELGLRF + F_TMP(IP_T + 2)
ELGLRC = ELGLRC + F_TMP(IP_T + 3)
ELGLLJ = ELGLLJ + F_TMP(IP_T + 4)
EL34EL = EL34EL + F_TMP(IP_T + 5)
EL34RF = EL34RF + F_TMP(IP_T + 6)
EL34RC = EL34RC + F_TMP(IP_T + 7)
EL34LJ = EL34LJ + F_TMP(IP_T + 8)
C
DO 800 I=0,NUNRE2-1
ELREL(I+1) = ELREL(I+1) + F_TMP(IP_T + 9 + I)
ELRRF(I+1) = ELRRF(I+1) + F_TMP(IP_T + NUNRE2 + 9 + I)
ELRRC(I+1) = ELRRC(I+1) + F_TMP(IP_T + 2NUNRE2 + 9 + I)
ELRLJ(I+1) = ELRLJ(I+1) + F_TMP(IP_T + 3
NUNRE2 + 9 + I)
800 CONTINUE
C
DO 810 I=0,NDIM-1
VIRLR(I+1) = VIRLR(I+1) +
$ F_TMP(IP_T + 4*NUNRE2 + 9 + I)
810 CONTINUE

805 CONTINUE
C
C+++ now collect forces from all the processors
C
CJMCT CSCO removed F_TMP,
C$DOACROSS
C$& LOCAL(ICPU,I,IP_T),
C$& SHARE(MATOM,FLR,NATTOT,NDIM)
C
DO 950 I = 0 , NATTOTNDIM - 1
DO 970 ICPU = 2, MPROC
IP_T = MATOM
(ICPU-2) + 4NUNRE2 + NDIM + 9
FLR(I+1) = FLR(I+1) + F_TMP(IP_T + I)
970 CONTINUE
950 CONTINUE
C
C+++ now sort the interleaved pair lists into one single big list
C
DO 610 ICPU = 1,NCPU
NPCPU(ICPU) = 0
610 CONTINUE
C
NSZPL = 0
DO 700 NI = 1,NRAGTO
ICPU = MOD(NI-1,NCPU)+1
NPJ_T = MX_JNB
(ICPU-1)
IAUX = NPCPU(ICPU)
ITMP = INB(NI)
DO 600 I = IAUX+1,IAUX+ITMP
NSZPL = NSZPL + 1
JNB(NSZPL) = JNB_TMP(NPJ_T+I)
600 CONTINUE
NPCPU(ICPU) = IAUX+ITMP
700 CONTINUE
C
ENDIF
C
C+++ done MP stuff
C
ELSE
C
C+++ Now the MP stuff (no long-range)
C
IF (.NOT.LMPROC) THEN
C
C+++ the unique processor gets the normal pair list arrays
C
CALL NBNONE_MP(NRAGTO,3,XCG,RCUTP2,
$ INB,JNB,NPCPU(1),MAXJNB,1)
NSZPL = NPCPU(1)
C
ELSE
C
CJMCT CSCO REMOVED JNB_TMP,
C$DOACROSS
C$& LOCAL(ICPU,NPJ_T),
C$& SHARE(NATTOT,NRAGTO,
C$& XCG,RCUTP2,INB,JNB,NSZPL,MX_JNB,
C$& NCPU,NPCPU)
C
DO 1000 ICPU = 1 , NCPU
C
C+++ all the processors get temporary arrays for the pair list
C
NPJ_T = MX_JNB*(ICPU-1)
CALL NBNONE_MP(NRAGTO,3,XCG,RCUTP2,INB,
$ JNB_TMP(NPJ_T+1),NPCPU(ICPU),MX_JNB,
$ ICPU)
C
1000 CONTINUE
C
C+++ now sort the interleaved pair lists into one single big list
C
DO 310 ICPU = 1,NCPU
NPCPU(ICPU) = 0
310 CONTINUE
C
NSZPL = 0
DO 400 NI = 1,NRAGTO
ICPU = MOD(NI-1,NCPU)+1
NPJ_T = MX_JNB*(ICPU-1)
IAUX = NPCPU(ICPU)
ITMP = INB(NI)
DO 300 I = IAUX+1,IAUX+ITMP
NSZPL = NSZPL + 1
JNB(NSZPL) = JNB_TMP(NPJ_T+I)
300 CONTINUE
NPCPU(ICPU) = IAUX+ITMP
400 CONTINUE
C
ENDIF
C
C+++ done MP stuff
C
C
ENDIF
C
C+++ flag errors and weird cases
C
LSTOP = .FALSE.
DO 301 I=1,MPROC
C
IF (LMXPL(I)) THEN
PRINT ,PRGSTR,’ CPU ',I
PRINT ,‘PAIRLIST OVERFLOW:’
PRINT ,‘MAXIMUM SIZE PER PROCESSOR: MX_JNB =’,MX_JNB
LSTOP = .TRUE.
ENDIF
C
IF (LDOLNG .OR. LFREE) THEN
IF (LMXAT(I)) THEN
PRINT ,PRGSTR,’ CPU ',I
PRINT ,PRGSTR,’ TOO MANY ATOMS INTERACTING WITH A ',
$ ’ PRIMARY CHARGE GROUP, INCREASE MXCGAT IN ',PRGSTR
LSTOP = .TRUE.
ENDIF
C
IF (LMXCG(I)) THEN
PRINT ,PRGSTR,’ CPU ',I
PRINT ,PRGSTR,’ TOO MANY CHARGE GROUPS INTERACTING ',
$ ’ WITH A ',
$ ’ PRIMARY CHARGE GROUP, INCREASE MXCGCG IN ',PRGSTR
LSTOP = .TRUE.
ENDIF
ENDIF
C
301 CONTINUE
C
IF (LSTOP) THEN
STOP
ENDIF
C
C+++ correct virial if required
C
IF (LDOVIR) THEN
DO 911 M=1,NDIM
VIRLR(M) = VIRLR(M)0.5
911 CONTINUE
ENDIF
C
C+++ return derivative variables according to request
C
IF (.NOT.L3D4D) THEN
EL34EL = 0.0
EL34RF = 0.0
EL34RC = 0.0
EL34LJ = 0.0
ENDIF
IF (.NOT.LDOPER) THEN
ELGLEL = 0.0
ELGLRF = 0.0
ELGLRC = 0.0
ELGLLJ = 0.0
ENDIF
C
PRINT ‘(1X,’‘NBPML:’’,I10,’’ ELEMENTS IN THE PAIRLIST’’)’,
$ NSZPL
C
IF (LDOTRA) THEN
CALL TRACO(NATTOT,0,XCOORD,BETA,-1,LEVERY)
ENDIF
C
RETURN
C
C------------END NBPML
C
END
C
C
COMMSUBR CLCGEO
C SUBROUTINE CLCGEO(NATTOT,NRAGTO,NPM,NSM,NDIM,XSRC,XDST)
C
C CLCGEO calulates the geometric centres of
C all CG of the solute and the first atom coordinate of
C solute molecules.
C This routine is called by L
C
C+++ Unchanged, not much to optimize
C
COMMEND
SUBROUTINE CLCGEO(NATTOT,NRAGTO,NPM,NSM,NDIM,XSRC,XDST)
C
C INCLUDES
C
INCLUDE ‘coordsz.h’
INCLUDE ‘toposz.h’
INCLUDE ‘topoar.h’
C
C ARGUMENTS
C
INTEGER NATTOT,NRAGTO,NPM,NSM,NDIM
REAL XSRC(NDIM
NATTOT),XDST(NDIM
NRAGTO)
C
C LOCAL VARIABLES
C
INTEGER M,J,NSRC3,NDST3,NCG,NPP
INTEGER NINC,ILACG,IFICG,NUMCG
REAL XTEMP(MAXDIM)
C
C------------ BEGIN CLCGEO
C
NSRC3 = 0
NDST3 = 0
C
C------------ FIRST DO THE SOLUTE ( CENTRE OF GEOMETRY )
C
DO 10 NPP=1,NPM
IFICG = 1
DO 20 NCG=1,NCAG
ILACG = INC(NCG)
NUMCG = ILACG - IFICG + 1
IF (NUMCG .GT. 1) THEN
DO 30 M=1,NDIM
XTEMP(M) = 0.0
30 CONTINUE
C
DO 40 J=1,NUMCG
DO 50 M=1,NDIM
XTEMP(M) = XTEMP(M) + XSRC(NSRC3+M)
50 CONTINUE
NSRC3 = NSRC3 + NDIM
40 CONTINUE
C
DO 60 M=1,NDIM
XDST(NDST3+M) = XTEMP(M)/NUMCG
60 CONTINUE
ELSE
DO 70 M=1,NDIM
XDST(NDST3+M) = XSRC(NSRC3+M)
70 CONTINUE
NSRC3 = NSRC3 + NDIM
ENDIF
NDST3 = NDST3 + NDIM
IFICG = ILACG + 1
20 CONTINUE
10 CONTINUE
C
C------------ NOW DO THE SOLVENT (COORDINATE OF THE FIRST ATOM)
C
NINC = NDIM
NRAM
DO 100 NPP=1,NSM
DO 110 M=1,NDIM
XDST(NDST3+M) = XSRC(NSRC3+M)
110 CONTINUE
NDST3 = NDST3 + NDIM
NSRC3 = NSRC3 + NINC
100 CONTINUE
C
C------------ END CLCGEO
C
END
C
C
COMMSUBR CCCCCC P.HUNENBERGER, ZUERICH, MAY. 1996 CCCCCCCCCCCCCCCCCCCCCC
C
C SUBROUTINE NBNONE_MP(NRAGTO,NCALCD,XCG,RCUTP2,
C $ INB,JNB,NJPTR,MAXPLO,ICPU)
C
C NBNONE_MP will calculate the part of the pair-list
C allocated to processor ICPU (.LE.NCPU), that is, charge groups NI
C with MOD(NI,NCPU).EQ.(ICPU-1).
C No long range force is calculated.
C
COMMEND
C
SUBROUTINE NBNONE_MP(NRAGTO,NCALCD,XCG,RCUTP2,
$ INB,JNB,NJPTR,MAXPLO,ICPU)
C
C INCLUDES
C
INCLUDE ‘toposz.h’
INCLUDE ‘coordsz.h’
INCLUDE ‘forcesz.h’
INCLUDE ‘box.h’
C
C ARGUMENTS
C
INTEGER NRAGTO,NCALCD,INB(NRAGTO),JNB(MAXJNB),NJPTR
INTEGER ICPU
INTEGER MAXPLO
REAL XCG(NRAGTO
NDIM)
REAL RCUTP2
C
C LOCAL VARIABLES
C
INTEGER IX3,JX3,NI,NJ,JSAVE
REAL XIJ1,XIJ2,XIJ3,XCGI1,XCGI2,XCGI3,RIJCG2,DSTTMP
CHARACTER PRGSTR
(6)
C
C
C+++ ALL CPU VARIABLES
C The following variables are common to all CPU’s, in NONBML and
C NBPML
C
INTEGER MAXCPU
PARAMETER (MAXCPU=32)
C
REAL ONE
PARAMETER (ONE = 1.0)
REAL ENEPS
PARAMETER (ENEPS = 1.E-20)
REAL ZEROC6
PARAMETER (ZEROC6 = 1.0E-20)
REAL ZEROAL
PARAMETER (ZEROAL = 1.0E-20)
C
LOGICAL LFULL4,LPART4,LFULL3,LPART3,LPERTL
LOGICAL LMONO,LRECT,LOCTO,LVAC,LDOVIR,LDOTRA,L4D,LNOSVT
LOGICAL LEMIS,L14MIS,LMXAT,LMXCG,LMXPL
DIMENSION LEMIS(MAXCPU),L14MIS(MAXCPU),
$ LMXAT(MAXCPU),LMXCG(MAXCPU),LMXPL(MAXCPU)
INTEGER ITSUBM(MAXCAG)
INTEGER NRAGP,NRPT,IOFF
INTEGER NCPU
REAL RLMA, RL2A, RLMB, RL2B, RLLA, RLLB, RLB, RLA, RDLB, RDLA
REAL RLARMU, RLBRMU, RM4D, RM3D, RDM3D, RDM4D
REAL RMCST,RMDER
REAL RFF, RFE, RFC, C1, RCRF2
REAL PININV
REAL COSB,COSB2,BOXOH,BOXOQ,PYE,CONV,BETAR
C
COMMON /NBRREL/ RLMA,RL2A,RLMB,RL2B,RLA,RLB,
$ RLLB,RLLA,RDLA,RDLB,RM3D,RM4D,RDM3D,RDM4D,RLARMU,RLBRMU,
$ RMCST,RMDER,
$ RFF, RFE, RFC, C1, RCRF2,
$ COSB,COSB2,BOXOH,BOXOQ,PYE,CONV,BETAR,
$ PININV
C
COMMON /NBRLOG/ LFULL3,LFULL4,LPART3,LPART4,LPERTL,
$ LMONO,LRECT,LOCTO,LVAC,LDOVIR,LDOTRA,L4D,LNOSVT,
$ LEMIS,L14MIS,LMXAT,LMXCG,LMXPL
C
COMMON /NBRINT/ NCPU,NRAGP,NRPT,IOFF,ITSUBM
C
C+++ DONE: ALL CPU VARIABLES
C
C
C DATA
C
DATA PRGSTR /‘NBNONE’/
C
C------------BEGIN NBNONE_MP
C
LMXPL(ICPU) = .FALSE.
NJPTR = 0
C
C+++ only the primary charge-groups for this proc
C
DO 100 NI=ICPU,NRAGTO,NCPU
C
C+++ CG NI can have at most NRAGTO-NI neighbours in the pair list:
C check here that there is no risk of overflow.
C
C
C+++ SCAN CHARGE GROUPS FOR MAKING PAIR LISTS (SOLUTE+SOLVENT)
C Solute charge groups are scanned with all other groups NJ>NI.
C If the center of geometry to center of geometry square distance
C is less than RCUTP2, the pair is sticked into the pair list
C JNB(). If MAXPLO is in danger to be exceeded, LMXPL(ICPU) is
C set and a RETURN is performed.
C
IF (NJPTR+NRAGTO-NI .GE. MAXPLO) THEN
LMXPL(ICPU) = .TRUE.
RETURN
ENDIF
JSAVE = NJPTR
IX3 = IOFF + NDIM
NI
XCGI1 = XCG(IX3 )
XCGI2 = XCG(IX3+1)
XCGI3 = XCG(IX3+2)
C
C+++ The pair list is always made in 3D - but the coordinates may be 4D
C
IF (LVAC) THEN
C
CDIR$ IVDEP
DO 200 NJ=NI+1, NRAGTO
JX3 = IOFF + NDIM
NJ
XIJ1 = XCGI1 - XCG(JX3 )
XIJ2 = XCGI2 - XCG(JX3+1)
XIJ3 = XCGI3 - XCG(JX3+2)
IF (XIJ1
XIJ1 + XIJ2XIJ2 + XIJ3XIJ3
$ .LE. RCUTP2) THEN
NJPTR = NJPTR + 1
JNB(NJPTR) = NJ
ENDIF
200 CONTINUE
C
ELSEIF (LOCTO) THEN
C
CDIR$ IVDEP
DO 1200 NJ=NI+1, NRAGTO
JX3 = IOFF + NDIMNJ
XIJ1 = XCGI1 - XCG(JX3 )
XIJ2 = XCGI2 - XCG(JX3+1)
XIJ3 = XCGI3 - XCG(JX3+2)
IF (XIJ1 .GE. BOXH(1)) THEN
XIJ1 = XIJ1 - BOX(1)
ELSEIF (XIJ1 .LT. -BOXH(1)) THEN
XIJ1 = XIJ1 + BOX(1)
ENDIF
IF (XIJ2 .GE. BOXH(2)) THEN
XIJ2 = XIJ2 - BOX(2)
ELSEIF (XIJ2 .LT. -BOXH(2)) THEN
XIJ2 = XIJ2 + BOX(2)
ENDIF
IF (XIJ3 .GE. BOXH(3)) THEN
XIJ3 = XIJ3 - BOX(3)
ELSEIF (XIJ3 .LT. -BOXH(3)) THEN
XIJ3 = XIJ3 + BOX(3)
ENDIF
RIJCG2 = XIJ1
XIJ1 + XIJ2XIJ2 + XIJ3XIJ3
DSTTMP = BOXOQ - ABS(XIJ1)-ABS(XIJ2)-ABS(XIJ3)
IF (DSTTMP .LT. 0.0) THEN
RIJCG2 = RIJCG2 + DSTTMPBOX(1)
ENDIF
IF (RIJCG2 .LE. RCUTP2) THEN
NJPTR = NJPTR + 1
JNB(NJPTR) = NJ
ENDIF
1200 CONTINUE
C
ELSEIF (LRECT) THEN
C
CDIR$ IVDEP
DO 2200 NJ=NI+1, NRAGTO
JX3 = IOFF + NDIM
NJ
XIJ1 = XCGI1 - XCG(JX3 )
XIJ2 = XCGI2 - XCG(JX3+1)
XIJ3 = XCGI3 - XCG(JX3+2)
IF (XIJ1 .GE. BOXH(1)) THEN
XIJ1 = XIJ1 - BOX(1)
ELSEIF (XIJ1 .LT. -BOXH(1)) THEN
XIJ1 = XIJ1 + BOX(1)
ENDIF
IF (XIJ2 .GE. BOXH(2)) THEN
XIJ2 = XIJ2 - BOX(2)
ELSEIF (XIJ2 .LT. -BOXH(2)) THEN
XIJ2 = XIJ2 + BOX(2)
ENDIF
IF (XIJ3 .GE. BOXH(3)) THEN
XIJ3 = XIJ3 - BOX(3)
ELSEIF (XIJ3 .LT. -BOXH(3)) THEN
XIJ3 = XIJ3 + BOX(3)
ENDIF
IF (XIJ1XIJ1 + XIJ2XIJ2 + XIJ3XIJ3
$ .LE. RCUTP2) THEN
NJPTR = NJPTR + 1
JNB(NJPTR) = NJ
ENDIF
2200 CONTINUE
C
ELSE
C
CDIR$ IVDEP
DO 3200 NJ=NI+1, NRAGTO
JX3 = IOFF + NDIM
NJ
XIJ1 = XCGI1 - XCG(JX3 )
XIJ2 = XCGI2 - XCG(JX3+1)
XIJ3 = XCGI3 - XCG(JX3+2)
IF (XIJ1 .GE. BOXH(1)) THEN
XIJ1 = XIJ1 - BOX(1)
ELSEIF (XIJ1 .LT. -BOXH(1)) THEN
XIJ1 = XIJ1 + BOX(1)
ENDIF
IF (XIJ2 .GE. BOXH(2)) THEN
XIJ2 = XIJ2 - BOX(2)
ELSEIF (XIJ2 .LT. -BOXH(2)) THEN
XIJ2 = XIJ2 + BOX(2)
ENDIF
IF (XIJ3 .GE. BOXH(3)) THEN
XIJ3 = XIJ3 - BOX(3)
ELSEIF (XIJ3 .LT. -BOXH(3)) THEN
XIJ3 = XIJ3 + BOX(3)
ENDIF
IF (XIJ1XIJ1 + XIJ2XIJ2 + XIJ3XIJ3
$ + COSB2
XIJ1XIJ3 .LE. RCUTP2) THEN
NJPTR = NJPTR + 1
JNB(NJPTR) = NJ
ENDIF
3200 CONTINUE
C
ENDIF
C
INB(NI) = NJPTR - JSAVE
C
C+++ SCAN CHARGE GROUPS FOR MAKING PAIR LISTS (SOLUTE+SOLVENT)
C
C-------I—
C
100 CONTINUE
INB(NRAGTO) = 0
C
C------------END NBNONE_MP
C
END
C
C
COMMSUBR CCCCCC P.HUNENBERGER, ZUERICH, MAY. 1996 CCCCCCCCCCCCCCCCCCCCCC
C
C SUBROUTINE NBWITH_MP(NATTOT,NPM,NSM,NRAGTO,NCALCD,
C $ XCG,XCOORD,XR,F,
C $ NSPM,NSP,VIR,
C $ RCUTP2,RCUTL2,
C $ INB,JNB,NJPTR,
C $ MAXNRE,NRELKP,IAGRP,
C $ EGLEL,EGLRF,EGLRC,EGLLJ,
C $ EG34EL,EG34RF,EG34RC,EG34LJ,
C $ NUNRE2,EEL,ERF,ERC,ELJ,
C $ ALPHLJ,ALPHC,
C $ LPION,LEVERY,LFREE,
C $ MAXPLO,ICPU)
C
C
C NBWITH_MP will calculate the part of the pair-list
C allocated to processor ICPU (.LE.NCPU), that is, charge groups NI
C with MOD(NI,NCPU).EQ.(ICPU-1). A twin-range force is also
C calculated. A local pair-list is created for each charge-group
C NI listing charge-groups NJ in the twin-range. This list (array
C KCGLOC()) is then processed in much the same way as in NONBML_MP.
C
C Speed-up is achieved by building for each primary charge group
C a J-atom vector containing
C - for each NI information about all atoms J interacting with
C atoms of primary charge group NI
C - for each atom I in NI, successively, information about pairs IJ
C The integers are stored in KLOC(), the reals in DLOC(),
C perturbation indicator in LPLOC() and path-integral exclusions in
C LEXCL(). No topology exclusions and third neighbours to handle.
C
C The vector is processed step by step by small loops which allows
C - a good software pipelining
C - moving many IF statements from the J atom loop to the I atom loop
C
C General scheme:
C
C Loop over NI solute
C Loop over NJ from pair-list
C Loop over J from NJ
C Store topology info for J
C Loop over I from NI
C Loop over J-vector to get IJ interaction params (and PI exclusions)
C Loop over J-vector to get IJ distances and vector
C Loop over J-vector to get IJ interaction (perturbed/unperturbed)
C Loop over J-vector to get IJ force and virial
C
C Loop over NI solvent
C Loop over NJ from pair-list
C Loop over J from NJ
C Store topology info for J
C Loop over I from NI
C Loop over J-vector to get IJ interaction params
C Loop over J-vector to get IJ distances and vector
C Loop over J-vector to get IJ interaction (unperturbed)
C Loop over J-vector to get IJ force and virial
C
COMMEND
C
SUBROUTINE NBWITH_MP(NATTOT,NPM,NSM,NRAGTO,NCALCD,
$ XCG,XCOORD,XR,F,
$ NSPM,NSP,VIR,
$ RCUTP2,RCUTL2,
$ INB,JNB,NJPTR,
$ MAXNRE,NRELKP,IAGRP,
$ EGLEL,EGLRF,EGLRC,EGLLJ,
$ EG34EL,EG34RF,EG34RC,EG34LJ,
$ NUNRE2,EEL,ERF,ERC,ELJ,
$ ALPHLJ,ALPHC,
$ LPION,LEVERY,LFREE,
$ MAXPLO,ICPU)
C
C INCLUDES
C
INCLUDE ‘coordsz.h’
INCLUDE ‘box.h’
INCLUDE ‘toposz.h’
INCLUDE ‘topoar.h’
INCLUDE ‘forcesz.h’
INCLUDE ‘pertsz.h’
INCLUDE ‘pertar.h’
C
C ARGUMENTS
C
LOGICAL LEVERY,LPION,LFREE
INTEGER NRAGTO,NCALCD,INB(NRAGTO)
INTEGER JNB(MAXJNB),NJPTR,NATTOT,NPM,NSM
INTEGER NSPM,NSP(NSPM+1),NUNRE2
INTEGER MAXNRE,NRELKP(MAXNRE,MAXNRE)
INTEGER MAXPLO,ICPU
REAL ALPHLJ,ALPHC
REAL XCG(NDIM
NRAGTO)
REAL XCOORD(NDIMNATTOT),F(NDIMNATTOT),XR(NDIMNATTOT)
REAL VIR(NDIM)
REAL EEL(NUNRE2),EGLEL,EG34EL
REAL ERF(NUNRE2),EGLRF,EG34RF
REAL ERC(NUNRE2),EGLRC,EG34RC
REAL ELJ(NUNRE2),EGLLJ,EG34LJ
REAL RCUTP2,RCUTL2
CHARACTER IAGRP(NATTOT)
C
C LOCAL VARIABLES
C
LOGICAL LJSOLV,LIZERO,LJZERO,LICGP
INTEGER NRAGT,NIPROT,IAUX,JSAVE,JTOPCG,NJPROT
INTEGER ITOPCG,NI,II,JJ,I,IPROT,IX3,JX3,NJ,NNJ,NUMIJ
INTEGER IM,IILOC3
INTEGER ILATYP,ILBTYP,INTIJA,INTIJB,ILGRP,JLGRP,NDXGRP
INTEGER IFTOPA,ILTOPA,IFCORA,ILCORA,JFTOPA,JLTOPA
INTEGER IATOFF,ICOPTR,ICGPTR,INUMAT,INUMA3
INTEGER ILPICI,ISCELI,ISUBM,JSUBM
INTEGER JJLOC3,IAUX2,ISCLJI
INTEGER ITMP,JLPICI,NICPU
REAL RIJ2,RIJCG2,RIJ4, RIJ6,DSTTMP,PIFACT
REAL SOFTLJ,SFLJ,DF3
REAL XCGI1,XCGI2,XCGI3
REAL XI1,XI2,XI3,XI4,XRI1,XRI2,XRI3,XIJ1,XIJ2,XIJ3,XIJ4
REAL FI1,FI2,FI3,FI4
REAL CGILA,CGILB
REAL XH
REAL RIJ3D2, VELA4D, VRFA4D, DFA4D, DFRFA4
REAL VLJA4D
REAL DLLJA4, DLELA4, DLRFA4, VELB4D
REAL VLJB4D, VRFB4D, DFB4D, DFRFB4, DLLJB4, DLELB4
REAL DLRFB4
REAL VELA3D, VLJA3D, VRFA3D, DFA3D, DFRFA3, DLLJA3
REAL DLELA3, DLRFA3, VELB3D, VLJB3D, VRFB3D
REAL DFB3D, DFRFB3, DLLJB3, DLELB3, DLRFB3
REAL VEL4D, VRF4D, VLJ4D, DF4D4D, DF4D, DLLJ4D, DLEL4D
REAL DLRF4D
REAL VEL3D, VRF3D, VLJ3D, DF3D, DLLJ3D, DLEL3D, DLRF3D
REAL DDM3EL, DDM3RF, DDM3LJ
REAL VRCA,VRCB,VRC
REAL QIQJA,QIQJB,QRFR
REAL SOFTEL, SFEL
REAL RSOFT, RSOFT2, RRF, RRF2, RFSQRT
REAL RIJINV, RIJIN2, RIJIN6
REAL CA12,CA6,CB12,CB6
REAL CRA1,CRA12
REAL FACRF
CHARACTER PRGSTR
(6)
C
LOGICAL LLIPER(MXATCG),LIPERT
INTEGER ITYPEA(MXATCG),ITYPEB(MXATCG)
INTEGER IGRP(MXATCG),ISCLJL(MXATCG),ISCELL(MXATCG)
INTEGER ILPIC(MXATCG)
REAL CGIA(MXATCG),CGIB(MXATCG)
REAL XI(MAXLAT)
REAL FI(MAXLAT)
REAL XRI(MAXLAT)
C
C
C
C+++ ALL CPU VARIABLES
C The following variables are common to all CPU’s, in NONBML and
C NBPML
C
INTEGER MAXCPU
PARAMETER (MAXCPU=32)
C
REAL ONE
PARAMETER (ONE = 1.0)
REAL ENEPS
PARAMETER (ENEPS = 1.E-20)
REAL ZEROC6
PARAMETER (ZEROC6 = 1.0E-20)
REAL ZEROAL
PARAMETER (ZEROAL = 1.0E-20)
C
LOGICAL LFULL4,LPART4,LFULL3,LPART3,LPERTL
LOGICAL LMONO,LRECT,LOCTO,LVAC,LDOVIR,LDOTRA,L4D,LNOSVT
LOGICAL LEMIS,L14MIS,LMXAT,LMXCG,LMXPL
DIMENSION LEMIS(MAXCPU),L14MIS(MAXCPU),
$ LMXAT(MAXCPU),LMXCG(MAXCPU),LMXPL(MAXCPU)
INTEGER ITSUBM(MAXCAG)
INTEGER NRAGP,NRPT,IOFF
INTEGER NCPU
REAL RLMA, RL2A, RLMB, RL2B, RLLA, RLLB, RLB, RLA, RDLB, RDLA
REAL RLARMU, RLBRMU, RM4D, RM3D, RDM3D, RDM4D
REAL RMCST,RMDER
REAL RFF, RFE, RFC, C1, RCRF2
REAL PININV
REAL COSB,COSB2,BOXOH,BOXOQ,PYE,CONV,BETAR
C
COMMON /NBRREL/ RLMA,RL2A,RLMB,RL2B,RLA,RLB,
$ RLLB,RLLA,RDLA,RDLB,RM3D,RM4D,RDM3D,RDM4D,RLARMU,RLBRMU,
$ RMCST,RMDER,
$ RFF, RFE, RFC, C1, RCRF2,
$ COSB,COSB2,BOXOH,BOXOQ,PYE,CONV,BETAR,
$ PININV
C
COMMON /NBRLOG/ LFULL3,LFULL4,LPART3,LPART4,LPERTL,
$ LMONO,LRECT,LOCTO,LVAC,LDOVIR,LDOTRA,L4D,LNOSVT,
$ LEMIS,L14MIS,LMXAT,LMXCG,LMXPL
C
COMMON /NBRINT/ NCPU,NRAGP,NRPT,IOFF,ITSUBM
C
C+++ DONE: ALL CPU VARIABLES
C
C
C $Id: vectordef.inc,v 1.2 1997/10/17 08:42:40 wscott Exp $
C
C+++ SPECIFIC VARIABLES FOR THE OPTIMIZED VERSION (J-ATOM VECTOR)
C
C+++ maximal number of CG (solute/solvent) a CG may interact with,
C within the twin range
C
INTEGER MXCGCG
CJMCT
C PARAMETER (MXCGCG=100)
PARAMETER (MXCGCG=200)
C
C
C+++ maximal number of atoms a CG may interact with within the
C twin range - roughly, MXCGCG * average number of atoms in a CG
C
INTEGER MXCGAT
PARAMETER (MXCGAT=2000)
C
C+++ a few arrays
C
C KLOC contains info for a series of MXKLOC atoms J interacting
C with the atoms I of a primary charge group. Indices are defined
C by parameters for ease of modification.
C For each atom J, this is a set of NKLOC reals, i.e.
C
C SOLUTE loop:
C 1 pointer to first coordinate of atom j in the X,F arrays
C 2 j energy group j belongs to
C 3 j integer atom code in state A
C 4 j integer atom code in state B
C 5 j soft-core Lennard-Jones indicator
C 6 j soft-core Electrostatics indicator
C 7 j path-integral bead group number
C 8 ij energy group of pair IJ
C 9 ij pair Lennard-Jones soft-core indicator
C 10 ij pair Electrostatics soft-core indicator
C LPLOC is the corresponding array for perturbation indicator
C
C SOLVENT loop:
C 1 pointer to first coordinate of atom j in the X,F arrays
C 2 pointer to location of atom j in the topology arrays
C
INTEGER ICRD3J,IGRPJ,IACJA,IACJB,
$ ISCJLJ,ISCJEL,IPIBD,INDXIJ,
$ ISCLJP,ISCELP
INTEGER ISCRDJ,ISTOPJ,ISGRPJ
INTEGER NKLOC
PARAMETER (ICRD3J=1, IGRPJ=2, IACJA=3, IACJB=4,
$ ISCJLJ=5, ISCJEL=6, IPIBD=7, INDXIJ=8,
$ ISCLJP=9, ISCELP=10)
PARAMETER (ISCRDJ=1, ISTOPJ=2, ISGRPJ=3)
PARAMETER (NKLOC = ISCELP)
C
INTEGER KLOC(MXCGATNKLOC)
LOGICAL LPLOC(MXCGAT)
C
INTEGER MXKLOC,IKLOC
C
C
C DLOC contains info for a series of MXKLOC atom pairs ij.
C Indices are defined by parameters for ease of modification.
C For each pair ij, NDLOC reals, i.e.
C
C SOLUTE loop only:
C 0-3 ij vector (periodicity applied), 3 unused in 3D
C 4 ij square distance (4-dim)
C 5 ij square distance (3-dim)
C 6 i charge in state A
C 7 i charge in state B
C 8 ij charge product in state A
C 9 ij charge product in state B ( pert only )
C 10 ij C12 in state A
C 11 ij C12 in state B ( pert only )
C 12 ij C6 in state A
C 13 ij C6 in state B ( pert only )
C 14 ij C126 in state A ( pert only )
C 15 ij C126 in state B ( pert only )
C 16 ij force to vector ratio, 3D
C 17 ij force to vector ratio, 4D
C
C SOLVENT loop:
C 0-4 IDEM
C
INTEGER IDF3,IDF4TH,ICGJA,ICGJB,
$ ICA12,ICB12,ICA6,ICB6,
$ ICA126,ICB126,IQIQJA,IQIQJB
INTEGER NDLOC
PARAMETER (ICGJA=6, ICGJB=7, IQIQJA=8, IQIQJB=9,
$ ICA12=10, ICB12=11, ICA6=12, ICB6=13,
$ ICA126=14, ICB126=15, IDF3=16, IDF4TH=17)
PARAMETER (NDLOC = IDF4TH+1)
C
REAL DLOC(NDLOC
MXCGAT)
C
INTEGER IDLOC
C
INTEGER JVILOC,JEXLOC
INTEGER KCGLOC(MAXCAG)
LOGICAL LGENPA(MXCGAT)
LOGICAL LEXCL(MXCGAT)
INTEGER IILOC
PARAMETER (IILOC=0)
C
C+++ DONE: SPECIFIC VARIABLES FOR THE OPTIMIZED VERSION (J-ATOM VECTOR)
C
C-----I—
C
C DATA
C
DATA PRGSTR /‘NBWITH’/
C
C------------BEGIN NBWITH
C
NRAGT = NRAGP + NSM
LMXPL(ICPU) = .FALSE.
LMXCG(ICPU) = .FALSE.
LMXAT(ICPU) = .FALSE.
C
NJPTR = 0
C
DO 1 I=1,MXCGAT
LGENPA(I) = .TRUE.
1 CONTINUE
C
DO 100 NI=ICPU,NRAGP,NCPU
C
IPROT = (NI-1) / NCAG
NIPROT = IPROT + 1
ITOPCG = NI - IPROTNCAG
ICGPTR = IPROT
NCAG + ITOPCG
ISUBM = ITSUBM(ICGPTR)
IF (ITOPCG.EQ.1) THEN
IFTOPA = 1
ELSE
IFTOPA = INC(ITOPCG-1)+1
ENDIF
ILTOPA = INC(ITOPCG)
IATOFF = IPROT * NRP
ICOPTR = NDIM * (IFTOPA + IATOFF - 1)
IFCORA = IFTOPA + IATOFF
ILCORA = ILTOPA + IATOFF
INUMAT = ILTOPA - IFTOPA + 1
INUMA3 = INUMATNDIM
C
IM = 1
LICGP = .FALSE.
DO 11 II = IFTOPA,ILTOPA
CGIA(IM) = CG(II)
ITYPEA(IM) = IAC(II)
ILPIC(IM) = IPIC(II)
IAUX = IPERT(II)
IF (LPERTL .AND. IAUX .NE. NOPERT) THEN
LLIPER(IM) = .TRUE.
LICGP = .TRUE.
CGIB(IM) = CGB(IAUX)
ITYPEB(IM) = IACB(IAUX)
ISCLJL(IM) = 1 - ISCLJ(IAUX)
ISCELL(IM) = 1 - ISCC(IAUX)
ELSE
LLIPER(IM) = .FALSE.
CGIB(IM) = CGIA(IM)
ITYPEB(IM) = ITYPEA(IM)
ISCLJL(IM) = 1
ISCELL(IM) = 1
ENDIF
IM = IM+1
11 CONTINUE
C
IM = 1
DO 12 II=IFCORA,ILCORA
IGRP(IM) = ICHAR(IAGRP(II))
IM = IM + 1
12 CONTINUE
DO 15 II=1,INUMA3
XI(II) = XCOORD(ICOPTR+II)
XRI(II) = XR(ICOPTR+II)
FI(II)= 0.0
15 CONTINUE
C
C
C
C+++ SCAN CHARGE GROUPS FOR MAKING PAIR LISTS (SOLUTE)
C Solute charge groups are scanned with all other groups NJ>NI.
C If the center of geometry to center of geometry square distance
C is less than RCUTP2, the pair is sticked into the pair list
C JNB(). If in the twin-range, RCUTP2<…<RCUTL2, it is sticked
C into the local pair list KCGLOC(). If MAXPLO is in danger to be
C exceeded, LMXPL(ICPU) is set and a RETURN is performed.
C
IF (NJPTR+NRAGTO-NI .GE. MAXPLO) THEN
LMXPL(ICPU) = .TRUE.
RETURN
ENDIF
IX3 = IOFF + NDIM
NI
XCGI1 = XCG(IX3 )
XCGI2 = XCG(IX3+1)
XCGI3 = XCG(IX3+2)
JSAVE = NJPTR
NUMIJ = 0
C
C+++ The pair list is always made in 3D - but the coordinates may be 4D
C
IF(LFREE) THEN
C
C+++ In this case, we limit the pair list/interactions to pairs including
C perturbed atoms
C
CDIR$ IVDEP
DO 902 NJ=NI+1, NRAGTO
LGENPA(NJ) = .FALSE.
LJSOLV = (NJ .GT. NRAGP)
IF (LICGP) THEN
LGENPA(NJ) = .TRUE.
ELSEIF (.NOT.LJSOLV) THEN
JTOPCG = MOD(NJ-1,NCAG) + 1
IF (JTOPCG .EQ. 1) THEN
JFTOPA = 1
ELSE
JFTOPA = INC(JTOPCG-1) + 1
ENDIF
JLTOPA = INC(JTOPCG)
DO 901 JJ=JFTOPA,JLTOPA
IAUX = IPERT(JJ)
IF (IAUX .NE. NOPERT) THEN
LGENPA(NJ) = .TRUE.
ENDIF
901 CONTINUE
ENDIF
902 CONTINUE
C
IF (LVAC) THEN
C
CDIR$ IVDEP
DO 5202 NJ=NI+1, NRAGTO
IF (LGENPA(NJ)) THEN
JX3 = IOFF + NDIMNJ
XIJ1 = XCGI1 - XCG(JX3 )
XIJ2 = XCGI2 - XCG(JX3+1)
XIJ3 = XCGI3 - XCG(JX3+2)
RIJCG2 = XIJ1
XIJ1 + XIJ2XIJ2 + XIJ3XIJ3
IF (RIJCG2 .LE. RCUTP2) THEN
NJPTR = NJPTR + 1
JNB(NJPTR) = NJ
ELSEIF (RIJCG2 .LE. RCUTL2) THEN
NUMIJ = NUMIJ + 1
KCGLOC(NUMIJ) = NJ
ENDIF
ENDIF
5202 CONTINUE
C
ELSEIF (LOCTO) THEN
C
CDIR$ IVDEP
DO 6200 NJ=NI+1, NRAGTO
IF (LGENPA(NJ)) THEN
JX3 = IOFF + NDIMNJ
XIJ1 = XCGI1 - XCG(JX3 )
XIJ2 = XCGI2 - XCG(JX3+1)
XIJ3 = XCGI3 - XCG(JX3+2)
IF (XIJ1 .GE. BOXH(1)) THEN
XIJ1 = XIJ1 - BOX(1)
ELSEIF (XIJ1 .LT. -BOXH(1)) THEN
XIJ1 = XIJ1 + BOX(1)
ENDIF
IF (XIJ2 .GE. BOXH(2)) THEN
XIJ2 = XIJ2 - BOX(2)
ELSEIF (XIJ2 .LT. -BOXH(2)) THEN
XIJ2 = XIJ2 + BOX(2)
ENDIF
IF (XIJ3 .GE. BOXH(3)) THEN
XIJ3 = XIJ3 - BOX(3)
ELSEIF (XIJ3 .LT. -BOXH(3)) THEN
XIJ3 = XIJ3 + BOX(3)
ENDIF
RIJCG2 = XIJ1
XIJ1 + XIJ2XIJ2 + XIJ3XIJ3
DSTTMP = BOXOQ - ABS(XIJ1)-ABS(XIJ2)-ABS(XIJ3)
IF (DSTTMP .LT. 0.0) THEN
RIJCG2 = RIJCG2 + DSTTMPBOX(1)
ENDIF
IF (RIJCG2 .LE. RCUTP2) THEN
NJPTR = NJPTR + 1
JNB(NJPTR) = NJ
ELSEIF (RIJCG2 .LE. RCUTL2) THEN
NUMIJ = NUMIJ + 1
KCGLOC(NUMIJ) = NJ
ENDIF
ENDIF
6200 CONTINUE
C
ELSEIF (LRECT) THEN
C
CDIR$ IVDEP
DO 7200 NJ=NI+1, NRAGTO
IF (LGENPA(NJ)) THEN
JX3 = IOFF + NDIM
NJ
XIJ1 = XCGI1 - XCG(JX3 )
XIJ2 = XCGI2 - XCG(JX3+1)
XIJ3 = XCGI3 - XCG(JX3+2)
IF (XIJ1 .GE. BOXH(1)) THEN
XIJ1 = XIJ1 - BOX(1)
ELSEIF (XIJ1 .LT. -BOXH(1)) THEN
XIJ1 = XIJ1 + BOX(1)
ENDIF
IF (XIJ2 .GE. BOXH(2)) THEN
XIJ2 = XIJ2 - BOX(2)
ELSEIF (XIJ2 .LT. -BOXH(2)) THEN
XIJ2 = XIJ2 + BOX(2)
ENDIF
IF (XIJ3 .GE. BOXH(3)) THEN
XIJ3 = XIJ3 - BOX(3)
ELSEIF (XIJ3 .LT. -BOXH(3)) THEN
XIJ3 = XIJ3 + BOX(3)
ENDIF
RIJCG2 = XIJ1XIJ1 + XIJ2XIJ2 + XIJ3XIJ3
IF (RIJCG2 .LE. RCUTP2) THEN
NJPTR = NJPTR + 1
JNB(NJPTR) = NJ
ELSEIF (RIJCG2 .LE. RCUTL2) THEN
NUMIJ = NUMIJ + 1
KCGLOC(NUMIJ) = NJ
ENDIF
ENDIF
7200 CONTINUE
C
ELSE
C
CDIR$ IVDEP
DO 8200 NJ=NI+1, NRAGTO
IF (LGENPA(NJ)) THEN
JX3 = IOFF + NDIM
NJ
XIJ1 = XCGI1 - XCG(JX3 )
XIJ2 = XCGI2 - XCG(JX3+1)
XIJ3 = XCGI3 - XCG(JX3+2)
IF (XIJ1 .GE. BOXH(1)) THEN
XIJ1 = XIJ1 - BOX(1)
ELSEIF (XIJ1 .LT. -BOXH(1)) THEN
XIJ1 = XIJ1 + BOX(1)
ENDIF
IF (XIJ2 .GE. BOXH(2)) THEN
XIJ2 = XIJ2 - BOX(2)
ELSEIF (XIJ2 .LT. -BOXH(2)) THEN
XIJ2 = XIJ2 + BOX(2)
ENDIF
IF (XIJ3 .GE. BOXH(3)) THEN
XIJ3 = XIJ3 - BOX(3)
ELSEIF (XIJ3 .LT. -BOXH(3)) THEN
XIJ3 = XIJ3 + BOX(3)
ENDIF
RIJCG2 = XIJ1XIJ1 + XIJ2XIJ2 + XIJ3XIJ3
$ + COSB2
XIJ1XIJ3
IF (RIJCG2 .LE. RCUTP2) THEN
NJPTR = NJPTR + 1
JNB(NJPTR) = NJ
ELSEIF (RIJCG2 .LE. RCUTL2) THEN
NUMIJ = NUMIJ + 1
KCGLOC(NUMIJ) = NJ
ENDIF
ENDIF
8200 CONTINUE
C
ENDIF
C
ELSE
C
C+++ In this case, we make full pair list/interactions
C
IF (LVAC) THEN
C
CDIR$ IVDEP
DO 202 NJ=NI+1, NRAGTO
JX3 = IOFF + NDIM
NJ
XIJ1 = XCGI1 - XCG(JX3 )
XIJ2 = XCGI2 - XCG(JX3+1)
XIJ3 = XCGI3 - XCG(JX3+2)
RIJCG2 = XIJ1XIJ1 + XIJ2XIJ2 + XIJ3XIJ3
IF (RIJCG2 .LE. RCUTP2) THEN
NJPTR = NJPTR + 1
JNB(NJPTR) = NJ
ELSEIF (RIJCG2 .LE. RCUTL2) THEN
NUMIJ = NUMIJ + 1
KCGLOC(NUMIJ) = NJ
ENDIF
202 CONTINUE
C
ELSEIF (LOCTO) THEN
C
CDIR$ IVDEP
DO 1200 NJ=NI+1, NRAGTO
JX3 = IOFF + NDIM
NJ
XIJ1 = XCGI1 - XCG(JX3 )
XIJ2 = XCGI2 - XCG(JX3+1)
XIJ3 = XCGI3 - XCG(JX3+2)
IF (XIJ1 .GE. BOXH(1)) THEN
XIJ1 = XIJ1 - BOX(1)
ELSEIF (XIJ1 .LT. -BOXH(1)) T

Hi,

It looks like omp_num_threads is a user defined function. It should be somewhere in source code. You can grep for a function in source to see which file it is defined. And perhaps see which configuration flag would set it.

Hongyon

Sorry for including all file content, I have used grep to extract information about “thread” from nbpmlpc.F. Next lines appear

#ifdef MP
INTEGER OMP_NUM_THREADS
INTEGER MALLOC
EXTERNAL MALLOC
EXTERNAL OMP_NUM_THREADS
DATA LMPROC /.TRUE./
#else
DATA LMPROC /.FALSE./
#endif
C
C------------BEGIN NBPML
C
#ifdef MP
C
IF (LMPROC.AND.(P_F_TMP.EQ.0 .OR. P_JNB_TMP.EQ.0)) THEN
MPROC = OMP_NUM_THREADS()

Thank you for your support

Hi,

You will need to look at other files if the author defined omp_num_threads anywhere. If not, you will need to ask the author what to do.

Hongyon