Possible QM/MM Bug

[Isaiah Sumner]

Hi,
  I have been trying to run parallel QM/MM calculations and have run
into a problem. I get different energies depending on how many
processors I use. The energy differences are significant. For example,
using 24 procs on a Cray XE6, I get,

 ENERGY| Total FORCE_EVAL ( FIST ) energy (a.u.):
6.710514008036768
 ENERGY| Total FORCE_EVAL ( QS ) energy (a.u.):
-64.534140477272018
 ENERGY| Total FORCE_EVAL ( QMMM ) energy (a.u.):
-55.496751397737185

but 48 procs gives,

 ENERGY| Total FORCE_EVAL ( FIST ) energy (a.u.):
6.710514008036785
 ENERGY| Total FORCE_EVAL ( QS ) energy (a.u.):
-64.482262744620684
 ENERGY| Total FORCE_EVAL ( QMMM ) energy (a.u.):
-55.444873665085836

The energy difference always comes in at the QS calculation. And it
seems to only occur if I use electrostatic coupling. This error
persists over different machines with different architectures and
compilers. Is this a known problem or have I made a mistake somewhere?

Thanks,
   Isaiah

 My input is the one given by Pauwels Ewald in this forum as an
example input. I've replicated the input below.

&FORCE_EVAL
  METHOD QMMM
  &DFT
    BASIS_SET_FILE_NAME GTH_BASIS_SETS
    POTENTIAL_FILE_NAME GTH_POTENTIALS
    &MGRID
      COMMENSURATE
      CUTOFF 280
    &END MGRID
    &QS
     METHOD GPW
     EXTRAPOLATION PS
     EXTRAPOLATION_ORDER 2
    &END QS
    &SCF
      EPS_SCF 5.0E-7
      MAX_SCF 20
      SCF_GUESS RESTART
      &OUTER_SCF
        EPS_SCF 5.0E-7
        MAX_SCF 20
      &END
      &OT
        MINIMIZER DIIS
        PRECONDITIONER FULL_ALL
        STEPSIZE 0.10
      &END OT
    &END SCF
    &XC
      &XC_FUNCTIONAL BLYP
      &END XC_FUNCTIONAL
    &END XC
  &END DFT
  &MM
    &FORCEFIELD
      PARMTYPE CHM
      PARM_FILE_NAME par_all27_prot_na.prm
      EI_SCALE14 1.0
      VDW_SCALE14 1.0
      &SPLINE
         EMAX_SPLINE 1.0
         RCUT_NB 10.0
      &END SPLINE
    &END FORCEFIELD
    &NEIGHBOR_LISTS
      GEO_CHECK F
    &END NEIGHBOR_LISTS
    &POISSON
      &EWALD
        EWALD_TYPE spme
        ALPHA .4
        GMAX 27
      &END EWALD
    &END POISSON
  &END MM
  &SUBSYS
    &CELL
      ABC 20.0 20.0 20.0
    &END CELL
    &TOPOLOGY
      CONN_FILE TAwater.psf
      CONNECTIVITY PSF
      COORD_FILE_NAME TAwater.pdb
      COORD_FILE_FORMAT PDB
    &END TOPOLOGY
    &KIND Cqm
      ELEMENT C
      BASIS_SET DZVP-GTH
      POTENTIAL GTH-BLYP-q4
    &END KIND
    &KIND Oqm
      ELEMENT O
      BASIS_SET DZVP-GTH
      POTENTIAL GTH-BLYP-q6
    &END KIND
    &KIND Nqm
      ELEMENT N
      BASIS_SET DZVP-GTH
      POTENTIAL GTH-BLYP-q5
    &END KIND
    &KIND Hqm
      ELEMENT H
      BASIS_SET DZVP-GTH
      POTENTIAL GTH-BLYP-q1
    &END KIND
    &KIND H_LINK
      ELEMENT H
      BASIS_SET DZVP-GTH
      POTENTIAL GTH-BLYP-q1
    &END KIND
  &END SUBSYS
  &QMMM
    &CELL
      ABC 10.0 10.0 10.0
    &END CELL
    ECOUPL GAUSS
    NOCOMPATIBILITY
    USE_GEEP_LIB 12
    NOCENTER0 FALSE
    NOCENTER TRUE
    &QM_KIND Cqm
      MM_INDEX 15 17 21
    &END QM_KIND
    &QM_KIND Oqm
      MM_INDEX 12 22
    &END QM_KIND
    &QM_KIND Nqm
      MM_INDEX 13
    &END QM_KIND
    &QM_KIND Hqm
      MM_INDEX 14 16 18 19 20
    &END QM_KIND
    &LINK
      QM_INDEX 13
      QM_KIND H
      MM_INDEX 11
      LINK_TYPE IMOMM
      ALPHA_IMOMM 1.50
    &END LINK
    &LINK
      QM_INDEX 21
      QM_KIND H
      MM_INDEX 26
      LINK_TYPE IMOMM
      ALPHA_IMOMM 1.50
    &END LINK
    &PERIODIC
      GMAX .1
      &MULTIPOLE
        RCUT 10.0
      &END MULTIPOLE
    &END PERIODIC
    &PRINT
      &QMMM_LINK_INFO
      &END QMMM_LINK_INFO
      &POTENTIAL
      &END POTENTIAL
    &END PRINT
  &END QMMM
&END FORCE_EVAL
&GLOBAL
  PREFERRED_FFT_LIBRARY FFTSG
  EXTENDED_FFT_LENGTHS T
  PROJECT TA
  RUN_TYPE GEO_OPT
  PRINT_LEVEL LOW
&END GLOBAL
&MOTION
  &GEO_OPT
    OPTIMIZER LBFGS
    MAX_ITER 4000
  &END GEO_OPT
&END MOTION