Merge branch 'master' of https://github.com/SimTk/openmm

plugins/rpmd/platforms/opencl/src/kernels/rpmdContraction.cl

@@ -113,6 +113,6 @@ __kernel void contractForces(__global real4* force, __global real4* contracted)
         
         // Store results.
 
-        force[index] = FORCE_SCALE*freal[indexInBlock];
+        force[index] = convert_real4(FORCE_SCALE*freal[indexInBlock]);
     }
 }

wrappers/python/simtk/openmm/app/modeller.py

@@ -39,7 +39,6 @@ from simtk.openmm.app.topology import Residue
 from simtk.openmm.vec3 import Vec3
 from simtk.openmm import System, Context, NonbondedForce, VerletIntegrator, LocalEnergyMinimizer
 from simtk.unit import nanometer, molar, elementary_charge, amu, gram, liter, degree, sqrt, acos, is_quantity, dot, norm
-import simtk.unit as unit
 import element as elem
 import os
 import random
@@ -809,7 +808,6 @@ class Modeller(object):
         bondedToAtomNoEP = []
         for atom in self.topology.atoms():
             bondedToAtom.append(set())
-            if atom.element is not None:
             bondedToAtomNoEP.append(set())
         for atom1, atom2 in self.topology.bonds():
             bondedToAtom[atom1.index].add(atom2.index)
@@ -889,7 +887,7 @@ class Modeller(object):
                     templateAtomPositions = len(template.atoms)*[None]
                     for index, atom in enumerate(template.atoms):
                         if atom in matchingAtoms:
-                            templateAtomPositions[index] = self.positions[matchingAtoms[atom].index]
+                            templateAtomPositions[index] = self.positions[matchingAtoms[atom].index].value_in_unit(nanometer)
                     for index, atom in enumerate(template.atoms):
                         if atom.element is None:
                             newTopology.addAtom(atom.name, None, newResidue)
@@ -903,7 +901,10 @@ class Modeller(object):
                                     elif site.type == 'average3':
                                         position = site.weights[0]*templateAtomPositions[index+site.atoms[0]] + site.weights[1]*templateAtomPositions[index+site.atoms[1]] + site.weights[2]*templateAtomPositions[index+site.atoms[2]]
                                     elif site.type == 'outOfPlane':
-                                        position = site.weights[0]*templateAtomPositions[index+site.atoms[0]] + site.weights[1]*templateAtomPositions[index+site.atoms[1]] + site.weights[2]*templateAtomPositions[index+site.atoms[2]]
+                                        v1 = templateAtomPositions[index+site.atoms[1]] - templateAtomPositions[index+site.atoms[0]]
+                                        v2 = templateAtomPositions[index+site.atoms[2]] - templateAtomPositions[index+site.atoms[0]]
+                                        cross = Vec3(v1[1]*v2[2]-v1[2]*v2[1], v1[2]*v2[0]-v1[0]*v2[2], v1[0]*v2[1]-v1[1]*v2[0])
+                                        position = templateAtomPositions[index+site.atoms[0]] + site.weights[0]*v1 + site.weights[1]*v2 + site.weights[2]*cross
                             if position is None and atom.type in drudeTypeMap:
                                 # This is a Drude particle.  Put it on top of its parent atom.
                                 
@@ -915,8 +916,8 @@ class Modeller(object):
                                 # and hope that energy minimization will fix it.
                                 
                                 knownPositions = [x for x in templateAtomPositions if x is not None]
-                                position = unit.sum(knownPositions)/len(knownPositions)
-                            newPositions.append(position)
+                                position = sum(knownPositions)/len(knownPositions)
+                            newPositions.append(position*nanometer)
         for bond in self.topology.bonds():
             if bond[0] in newAtoms and bond[1] in newAtoms:
                 newTopology.addBond(newAtoms[bond[0]], newAtoms[bond[1]])

wrappers/python/simtk/openmm/app/pdbfile.py

@@ -278,7 +278,12 @@ class PDBFile(object):
                     else:
                         atomName = atom.name
                     coords = positions[posIndex]
-                    print >>file, "ATOM  %5d %-4s %3s %s%4d    %8.3f%8.3f%8.3f  1.00  0.00" % (atomIndex%100000, atomName, resName, chainName, (resIndex+1)%10000, coords[0], coords[1], coords[2])
+                    line = "ATOM  %5d %-4s %3s %s%4d    %s%s%s  1.00  0.00" % (
+                        atomIndex%100000, atomName, resName, chainName,
+                        (resIndex+1)%10000, _format_83(coords[0]),
+                        _format_83(coords[1]), _format_83(coords[2]))
+                    assert len(line) == 66, 'Fixed width overflow detected'
+                    print >>file, line
                     posIndex += 1
                     atomIndex += 1
                 if resIndex == len(residues)-1:
@@ -297,3 +302,15 @@ class PDBFile(object):
         """
         print >>file, "END"
 
+
+def _format_83(f):
+    """Format a single float into a string of width 8, with ideally 3 decimal
+    places of precision. If the number is a little too large, we can
+    gracefully degrade the precision by lopping off some of the decimal
+    places. If it's much too large, we throw a ValueError"""
+    if -999.999 < f < 9999.999:
+        return '%8.3f' % f
+    if -9999999 < f < 99999999:
+        return ('%8.3f' % f)[:8]
+    raise ValueError('coordinate "%s" could not be represnted '
+                     'in a width-8 field' % f)

wrappers/python/src/swig_doxygen/swigInputBuilder.py

@@ -218,6 +218,7 @@ class SwigInputBuilder:
         for name in sorted(integratorSubclassList):
             self.fOut.write(",\n         OpenMM::%s" % name)
         self.fOut.write(");\n\n")
+        self.fOut.write("%factory(OpenMM::VirtualSite& OpenMM::System::getVirtualSite, OpenMM::TwoParticleAverageSite, OpenMM::ThreeParticleAverageSite, OpenMM::OutOfPlaneSite);\n\n")
         self.fOut.write("\n")
 
     def writeGlobalConstants(self):