[gmx-developers] Fluctuating charge code, Q-Chem interface, force matching.

Lee-Ping Wang leeping at stanford.edu
Mon Nov 5 08:58:23 CET 2012

Hi Gerrit and Michael,

Thanks for the replies.  Here are some more details.

- Concerning the interface with Q-Chem, we didn't have to make any changes
to the Q-Chem source code.  The two programs communicate using a file-based
interface; Gromacs writes two Q-Chem input files (qchem.in and ESPGrid) and
Q-Chem produces a file with electric fields (efield.dat).  Out of the three
features, this one is by far the easiest.  I'd be happy to provide some
example files and documentation.

I'm interested to hear about how you guys are planning to restructure the
QM/MM interface, because having a framework that works with multiple codes
is both challenging and an important problem.  

- Regarding the fluctuating charge force field, it is integrated as follows:

a) All of the computation is done in a new file, qtpie.c .  Most of the
computation is linear algebra because it involves solving for the
fluctuating charges at each time step.
b) The forces are computed in do_force_lowlevel right before do_nonbonded is
c) The extra data (electronegativities etc.) is stored in the "t_mdatoms"
structure, but I can create a separate structure if desired.
d) The fluctuating charge parameters are read using an extra section in the
topology file ([ qtpie ], which fits between [ atoms ] and [ bonds ]). 
e) One extra neighbor list is required, which is coded into ns.c .
f) Periodic boundary conditions are handled using shift functions where the
force goes smoothly to zero (i.e. the second derivative of the energy is
g) It's currently working in 4.0.7, but I'm working on porting it over to

I'm not an expert at GPU coding, but maybe there's a straightforward way to
use the GPU to execute the calls to BLAS and LAPACK.

I have a question about BLAS and LAPACK: Currently qtpie.c includes headers
from the Intel MKL library, "mkl_lapack.h" and "mkl_cblas.h".  This probably
needs to be changed.  What is the best way to call BLAS and LAPACK without
causing trouble for the user in the build process?  Should I be using the
Gromacs-provided libraries, and are there examples for how to do this?  (I
can't find many instances of linear algebra calls in the code base.)  Also,
how do I give users the option to link to external BLAS and LAPACK

- Regarding the energy and force derivatives, I've been thinking about it
some more and I should probably discard this code.  I don't think it's
possible to have them as an external tool because they require changing a
lot of low-level subroutines in bondfree.c and the nonbonded inner loops.
It facilitates force field parameterization but it isn't required;
furthermore, it requires too many changes to the code base and can be very
cumbersome to maintain.

My force field parameterization methods are a different from Greg Voth's
methods although I learned a lot from his work.  My goal is to create a
general framework where a variety of reference data can be incorporated into
the objective function, so that force field parameters can be simultaneously
optimized to fit experimental and theoretical data in a systematic way.  I'm
implementing these methods into a software package for parameterization
called ForceBalance (https://simtk.org/home/forcebalance).  It works using
file-based interfaces with Gromacs and several other codes, and I might like
to submit it as a "user contribution" at some point. 


- Lee-Ping

-----Original Message-----
From: gmx-developers-bounces at gromacs.org
[mailto:gmx-developers-bounces at gromacs.org] On Behalf Of Gerrit Groenhof
Sent: 04 November 2012 23:05
To: Discussion list for GROMACS development
Subject: Re: [gmx-developers] Fluctuating charge code, Q-Chem interface,
force matching.

Dear Lee-Ping,

On 1 and 2 we already had some exchange in the past, and I am pleased to see
that you got it done. 

Concerning the interface with Q-chem, if it is documented, and works without
too much coding on the Q-chem part (in contrast to gaussian!!), I'd be happy
to have it. There is an initiative to restructure the complete qmmm setup,
as we now are getting too many program-specific interfaces.

On 1, I think it will be interesting to have, but I think it also depends on
how well it is integrated (and integratable). I suppose it would also have
to work on GPU now?

Concerning the forcematching. Is this the same procedure as that of Voth and
co-workers. We used that, but as a seperate tool. Unless it is something
else, I'd think such code would be more suitable as a gmxtool?


On Nov 4, 2012, at 9:43 PM, Lee-Ping Wang wrote:

> Dear developers,
> I'm a postdoctoral fellow in the Stanford Chemistry department working 
> with Vijay Pande and Todd Martinez.  I earned my Ph.D. from the MIT 
> Chemistry department under the supervision of Troy Van Voorhis.  I've 
> been using MD software as a central part of my research for over five 
> years, and I am definitely a big fan of GROMACS. :)
> Over the course of my Ph.D. I did a lot of method development in Gromacs.
I'm writing to find out whether there's interest in integrating any of my
contributed code into the Gromacs code base.  
> My contributions are in three main categories:
> 1) A fluctuating charge code following the work of Steve Rick and Bruce
Berne, with further improvements following the work of Jiahao Chen and Todd
> 2) A QM/MM interface with Q-Chem, mainly following the existing QM/MM
interface with Gaussian.
> 3) Energy and force matching code for force field optimization.  This
consists mainly of adding analytic first and second parametric derivatives
for the energy and force terms into the low-level subroutines that compute
the force field contributions.
> Items (1) and (3) are incorporated into an upcoming publication which was
just accepted to JCTC.  
> These methods have been implemented into version 4.0.7.  Most of the code
(perhaps 95%) is in separate files with as little modification to the
original code as possible.  All of the functionality is controlled through
the input files and command line arguments.  There is no performance impact
for normal MD runs because the parametric derivatives are only switched on
when requested.
> If there is interest on your side, I'd be happy to do all of the work.  I
think the work involved is mainly (1) updating the code to be compatible
with the latest development branch, (2) writing unit tests, and (3) creating
> Please let me know and I'll do my best to follow the directions on the
Gromacs website and stay in touch.
> Thanks,
> - Lee-Ping--
> gmx-developers mailing list
> gmx-developers at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-developers
> Please don't post (un)subscribe requests to the list. Use the www 
> interface or send it to gmx-developers-request at gromacs.org.

gmx-developers mailing list
gmx-developers at gromacs.org
Please don't post (un)subscribe requests to the list. Use the www interface
or send it to gmx-developers-request at gromacs.org.

More information about the gromacs.org_gmx-developers mailing list