OOPSE program

[Dan Gezelter]

Francesco Mercuri wrote:
> 
> Dear Professor Gezelter,
> I've read your paper, recently appeared on J. Comput. Chem,
> concerning the development of the OOPSE program for molecular
> dynamics simulations.
> In particular, I am interested in molecular dynamics simulations
> of rigid bodies (e.g. TIP4P water) making use of time-reversible
> algorithms, as done in OOPSE by means of the symplectic DLM
> integrator. I did already some test runs on small water clusters.
> However, I would like to ask you a few questions concerning the
> OOPSE program, not documented in the manual:
> - Would it be possible to improve the time-reversibility of the
>   integration procedure by increasing the accuracy of floating
>   point numbers (e.g. from "double" to "long double")? Indeed,
>   long simulations are affected by a numerical "noise" that can
>   compromise time reversibility. Is there any data-structure
>   that can be modified in order to improve the representation of
>   floating point numbers?
> - Are you considering the implementation of the Ewald summation?
>   It would be an important issue, expecially for solid-state systems.
> 
> Thank you in advance for your help and best regards,
> 
> Francesco Mercuri
> Max-Planck-Institut fuer Chemische Physik fester Stoffe
> Noethnitzer Str. 40, 01187 Dresden, Germany


There are a few things you can do to improve the time-reversibility
of DLM at the expense of speed.  In the DLM integrator, the small angle
formula is used to compute the rotation matrices.  If you look in
the rotate function in Integrator.cpp, you'll find the following lines:

   // use a small angle aproximation for sin and cosine

   angleSqr = angle * angle;
   angleSqrOver4 = angleSqr / 4.0;
   top = 1.0 - angleSqrOver4;
   bottom = 1.0 + angleSqrOver4;

   cosAngle = top / bottom;
   sinAngle = angle / bottom;

If you replace these lines with:

     cosAngle = cos(angle);
     sinAngle = sin(angle);

you'll get higher accuracy on the rotational propagation.  Using long 
doubles would be a difficult change to make.

I should point out that we're nearing the release of OOPSE-2.0 which
has been an almost complete re-write of the code (the input files
remain the same, however).  We haven't yet implemented an Ewald
summation because we need a fast (PME or SPME) sum which knows
about point multipoles (dipoles, quadrupoles, etc.)   It is certainly
on our to-do list.

Let me know if you have any additional questions,

  --Dan Gezelter


-- 
***********************************************
   J. Daniel Gezelter
   Assistant Professor
   Department of Chemistry and Biochemistry
   251 Nieuwland Science Hall
   University of Notre Dame
   Notre Dame, IN 46556-5670

   phone:  +1 (574) 631-7595
   fax:    +1 (574) 631-6652
   e-mail: gezelter at nd.edu
   web:    http://www.nd.edu/~gezelter
************************************************