"vscode:/vscode.git/clone" did not exist on "309f777838bae2a07177456f7915fee0904ed1fe"
Commit c6ebf6e2 authored by Peter Eastman's avatar Peter Eastman
Browse files

Removed the word "harmonic" from lots of places it didn't belong, including class and method names

parent 58b094ce
......@@ -30,7 +30,7 @@
* -------------------------------------------------------------------------- */
/**
* This tests the Cuda implementation of CudaAmoebaHarmonicAngleForce.
* This tests the Cuda implementation of CudaAmoebaAngleForce.
*/
#include "openmm/internal/AssertionUtilities.h"
......@@ -119,22 +119,22 @@ static void getPrefactorsGivenAngleCosine( double cosine, double idealAngle, dou
return;
}
static void computeAmoebaHarmonicAngleForce(int bondIndex, std::vector<Vec3>& positions, AmoebaHarmonicAngleForce& amoebaHarmonicAngleForce,
static void computeAmoebaAngleForce(int bondIndex, std::vector<Vec3>& positions, AmoebaAngleForce& amoebaAngleForce,
std::vector<Vec3>& forces, double* energy, FILE* log ) {
int particle1, particle2, particle3;
double idealAngle;
double quadraticK;
amoebaHarmonicAngleForce.getAngleParameters(bondIndex, particle1, particle2, particle3, idealAngle, quadraticK );
amoebaAngleForce.getAngleParameters(bondIndex, particle1, particle2, particle3, idealAngle, quadraticK );
double cubicK = amoebaHarmonicAngleForce.getAmoebaGlobalHarmonicAngleCubic();
double quarticK = amoebaHarmonicAngleForce.getAmoebaGlobalHarmonicAngleQuartic();
double penticK = amoebaHarmonicAngleForce.getAmoebaGlobalHarmonicAnglePentic();
double sexticK = amoebaHarmonicAngleForce.getAmoebaGlobalHarmonicAngleSextic();
double cubicK = amoebaAngleForce.getAmoebaGlobalAngleCubic();
double quarticK = amoebaAngleForce.getAmoebaGlobalAngleQuartic();
double penticK = amoebaAngleForce.getAmoebaGlobalAnglePentic();
double sexticK = amoebaAngleForce.getAmoebaGlobalAngleSextic();
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicAngleForce: bond %d [%d %d %d] ang=%10.3f k=%10.3f [%10.3e %10.3e %10.3e %10.3e]\n",
(void) fprintf( log, "computeAmoebaAngleForce: bond %d [%d %d %d] ang=%10.3f k=%10.3f [%10.3e %10.3e %10.3e %10.3e]\n",
bondIndex, particle1, particle2, particle3, idealAngle, quadraticK, cubicK, quarticK, penticK, sexticK );
(void) fflush( log );
}
......@@ -201,7 +201,7 @@ static void computeAmoebaHarmonicAngleForce(int bondIndex, std::vector<Vec3>& p
*energy += energyTerm;
}
static void computeAmoebaHarmonicAngleForces( Context& context, AmoebaHarmonicAngleForce& amoebaHarmonicAngleForce,
static void computeAmoebaAngleForces( Context& context, AmoebaAngleForce& amoebaAngleForce,
std::vector<Vec3>& expectedForces, double* expectedEnergy, FILE* log ) {
// get positions and zero forces
......@@ -217,13 +217,13 @@ static void computeAmoebaHarmonicAngleForces( Context& context, AmoebaHarmonicAn
// calculates forces/energy
*expectedEnergy = 0.0;
for( int ii = 0; ii < amoebaHarmonicAngleForce.getNumAngles(); ii++ ){
computeAmoebaHarmonicAngleForce(ii, positions, amoebaHarmonicAngleForce, expectedForces, expectedEnergy, log );
for( int ii = 0; ii < amoebaAngleForce.getNumAngles(); ii++ ){
computeAmoebaAngleForce(ii, positions, amoebaAngleForce, expectedForces, expectedEnergy, log );
}
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicAngleForces: expected energy=%14.7e\n", *expectedEnergy );
(void) fprintf( log, "computeAmoebaAngleForces: expected energy=%14.7e\n", *expectedEnergy );
for( unsigned int ii = 0; ii < positions.size(); ii++ ){
(void) fprintf( log, "%6u [%14.7e %14.7e %14.7e]\n", ii, expectedForces[ii][0], expectedForces[ii][1], expectedForces[ii][2] );
}
......@@ -235,19 +235,19 @@ static void computeAmoebaHarmonicAngleForces( Context& context, AmoebaHarmonicAn
}
void compareWithExpectedForceAndEnergy( Context& context, AmoebaHarmonicAngleForce& amoebaHarmonicAngleForce,
void compareWithExpectedForceAndEnergy( Context& context, AmoebaAngleForce& amoebaAngleForce,
double tolerance, const std::string& idString, FILE* log) {
std::vector<Vec3> expectedForces;
double expectedEnergy;
computeAmoebaHarmonicAngleForces( context, amoebaHarmonicAngleForce, expectedForces, &expectedEnergy, log );
computeAmoebaAngleForces( context, amoebaAngleForce, expectedForces, &expectedEnergy, log );
State state = context.getState(State::Forces | State::Energy);
const std::vector<Vec3> forces = state.getForces();
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicAngleForces: expected energy=%14.7e %14.7e\n", expectedEnergy, state.getPotentialEnergy() );
(void) fprintf( log, "computeAmoebaAngleForces: expected energy=%14.7e %14.7e\n", expectedEnergy, state.getPotentialEnergy() );
for( unsigned int ii = 0; ii < forces.size(); ii++ ){
(void) fprintf( log, "%6u [%14.7e %14.7e %14.7e] [%14.7e %14.7e %14.7e]\n", ii,
expectedForces[ii][0], expectedForces[ii][1], expectedForces[ii][2], forces[ii][0], forces[ii][1], forces[ii][2] );
......@@ -272,7 +272,7 @@ void testOneAngle( FILE* log ) {
LangevinIntegrator integrator(0.0, 0.1, 0.01);
AmoebaHarmonicAngleForce* amoebaHarmonicAngleForce = new AmoebaHarmonicAngleForce();
AmoebaAngleForce* amoebaAngleForce = new AmoebaAngleForce();
double angle = 100.0;
double quadraticK = 1.0;
......@@ -280,14 +280,14 @@ void testOneAngle( FILE* log ) {
double quarticK = 1.0e-02;
double penticK = 1.0e-03;
double sexticK = 1.0e-04;
amoebaHarmonicAngleForce->addAngle(0, 1, 2, angle, quadraticK);
amoebaAngleForce->addAngle(0, 1, 2, angle, quadraticK);
amoebaHarmonicAngleForce->setAmoebaGlobalHarmonicAngleCubic(cubicK);
amoebaHarmonicAngleForce->setAmoebaGlobalHarmonicAngleQuartic(quarticK);
amoebaHarmonicAngleForce->setAmoebaGlobalHarmonicAnglePentic(penticK);
amoebaHarmonicAngleForce->setAmoebaGlobalHarmonicAngleSextic(sexticK);
amoebaAngleForce->setAmoebaGlobalAngleCubic(cubicK);
amoebaAngleForce->setAmoebaGlobalAngleQuartic(quarticK);
amoebaAngleForce->setAmoebaGlobalAnglePentic(penticK);
amoebaAngleForce->setAmoebaGlobalAngleSextic(sexticK);
system.addForce(amoebaHarmonicAngleForce);
system.addForce(amoebaAngleForce);
Context context(system, integrator, Platform::getPlatformByName( "Cuda"));
std::vector<Vec3> positions(numberOfParticles);
......@@ -297,14 +297,14 @@ void testOneAngle( FILE* log ) {
positions[2] = Vec3(0, 0, 1);
context.setPositions(positions);
compareWithExpectedForceAndEnergy( context, *amoebaHarmonicAngleForce, TOL, "testOneAngle", log );
compareWithExpectedForceAndEnergy( context, *amoebaAngleForce, TOL, "testOneAngle", log );
}
int main( int numberOfArguments, char* argv[] ) {
try {
std::cout << "TestCudaAmoebaHarmonicAngleForce running test..." << std::endl;
std::cout << "TestCudaAmoebaAngleForce running test..." << std::endl;
registerAmoebaCudaKernelFactories();
FILE* log = NULL;
testOneAngle( log );
......
......@@ -30,7 +30,7 @@
* -------------------------------------------------------------------------- */
/**
* This tests the Cuda implementation of HarmonicBondForce.
* This tests the Cuda implementation of AmoebaBondForce.
*/
#include "openmm/internal/AssertionUtilities.h"
......@@ -46,15 +46,15 @@
using namespace OpenMM;
const double TOL = 1e-5;
static void computeAmoebaHarmonicBondForce(int bondIndex, std::vector<Vec3>& positions, AmoebaHarmonicBondForce& amoebaHarmonicBondForce,
static void computeAmoebaBondForce(int bondIndex, std::vector<Vec3>& positions, AmoebaBondForce& amoebaBondForce,
std::vector<Vec3>& forces, double* energy ) {
int particle1, particle2;
double bondLength;
double quadraticK;
double cubicK = amoebaHarmonicBondForce.getAmoebaGlobalHarmonicBondCubic();
double quarticK = amoebaHarmonicBondForce.getAmoebaGlobalHarmonicBondQuartic();
amoebaHarmonicBondForce.getBondParameters(bondIndex, particle1, particle2, bondLength, quadraticK );
double cubicK = amoebaBondForce.getAmoebaGlobalBondCubic();
double quarticK = amoebaBondForce.getAmoebaGlobalBondQuartic();
amoebaBondForce.getBondParameters(bondIndex, particle1, particle2, bondLength, quadraticK );
double deltaR[3];
double r2 = 0.0;
......@@ -82,7 +82,7 @@ static void computeAmoebaHarmonicBondForce(int bondIndex, std::vector<Vec3>& po
}
static void computeAmoebaHarmonicBondForces( Context& context, AmoebaHarmonicBondForce& amoebaHarmonicBondForce,
static void computeAmoebaBondForces( Context& context, AmoebaBondForce& amoebaBondForce,
std::vector<Vec3>& expectedForces, double* expectedEnergy, FILE* log ) {
// get positions and zero forces
......@@ -98,13 +98,13 @@ static void computeAmoebaHarmonicBondForces( Context& context, AmoebaHarmonicBon
// calculates forces/energy
*expectedEnergy = 0.0;
for( int ii = 0; ii < amoebaHarmonicBondForce.getNumBonds(); ii++ ){
computeAmoebaHarmonicBondForce(ii, positions, amoebaHarmonicBondForce, expectedForces, expectedEnergy );
for( int ii = 0; ii < amoebaBondForce.getNumBonds(); ii++ ){
computeAmoebaBondForce(ii, positions, amoebaBondForce, expectedForces, expectedEnergy );
}
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicBondForces: expected energy=%15.7e\n", *expectedEnergy );
(void) fprintf( log, "computeAmoebaBondForces: expected energy=%15.7e\n", *expectedEnergy );
for( unsigned int ii = 0; ii < positions.size(); ii++ ){
(void) fprintf( log, "%6u [%15.7e %15.7e %15.7e]\n", ii, expectedForces[ii][0], expectedForces[ii][1], expectedForces[ii][2] );
}
......@@ -115,18 +115,18 @@ static void computeAmoebaHarmonicBondForces( Context& context, AmoebaHarmonicBon
}
void compareWithExpectedForceAndEnergy( Context& context, AmoebaHarmonicBondForce& amoebaHarmonicBondForce, double tolerance, const std::string& idString, FILE* log) {
void compareWithExpectedForceAndEnergy( Context& context, AmoebaBondForce& amoebaBondForce, double tolerance, const std::string& idString, FILE* log) {
std::vector<Vec3> expectedForces;
double expectedEnergy;
computeAmoebaHarmonicBondForces( context, amoebaHarmonicBondForce, expectedForces, &expectedEnergy, NULL );
computeAmoebaBondForces( context, amoebaBondForce, expectedForces, &expectedEnergy, NULL );
State state = context.getState(State::Forces | State::Energy);
const std::vector<Vec3> forces = state.getForces();
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicBondForces: expected energy=%15.7e %15.7e\n", expectedEnergy, state.getPotentialEnergy() );
(void) fprintf( log, "computeAmoebaBondForces: expected energy=%15.7e %15.7e\n", expectedEnergy, state.getPotentialEnergy() );
for( unsigned int ii = 0; ii < forces.size(); ii++ ){
(void) fprintf( log, "%6u [%15.7e %15.7e %15.7e] [%15.7e %15.7e %15.7e]\n", ii,
expectedForces[ii][0], expectedForces[ii][1], expectedForces[ii][2], forces[ii][0], forces[ii][1], forces[ii][2] );
......@@ -150,17 +150,17 @@ void testOneBond( FILE* log ) {
LangevinIntegrator integrator(0.0, 0.1, 0.01);
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* amoebaBondForce = new AmoebaBondForce();
double bondLength = 1.5;
double quadraticK = 1.0;
double cubicK = 2.0;
double quarticicK = 3.0;
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( cubicK );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( quarticicK );
amoebaHarmonicBondForce->addBond(0, 1, bondLength, quadraticK);
amoebaBondForce->setAmoebaGlobalBondCubic( cubicK );
amoebaBondForce->setAmoebaGlobalBondQuartic( quarticicK );
amoebaBondForce->addBond(0, 1, bondLength, quadraticK);
system.addForce(amoebaHarmonicBondForce);
system.addForce(amoebaBondForce);
Context context(system, integrator, Platform::getPlatformByName( "Cuda"));
std::vector<Vec3> positions(2);
......@@ -168,7 +168,7 @@ void testOneBond( FILE* log ) {
positions[1] = Vec3(0, 0, 0);
context.setPositions(positions);
compareWithExpectedForceAndEnergy( context, *amoebaHarmonicBondForce, TOL, "testOneBond", log );
compareWithExpectedForceAndEnergy( context, *amoebaBondForce, TOL, "testOneBond", log );
}
void testTwoBond( FILE* log ) {
......@@ -182,18 +182,18 @@ void testTwoBond( FILE* log ) {
LangevinIntegrator integrator(0.0, 0.1, 0.01);
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* amoebaBondForce = new AmoebaBondForce();
double bondLength = 1.5;
double quadraticK = 1.0;
double cubicK = 2.0;
double quarticicK = 3.0;
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( cubicK );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( quarticicK );
amoebaHarmonicBondForce->addBond(0, 1, bondLength, quadraticK);
amoebaHarmonicBondForce->addBond(1, 2, bondLength, quadraticK);
amoebaBondForce->setAmoebaGlobalBondCubic( cubicK );
amoebaBondForce->setAmoebaGlobalBondQuartic( quarticicK );
amoebaBondForce->addBond(0, 1, bondLength, quadraticK);
amoebaBondForce->addBond(1, 2, bondLength, quadraticK);
system.addForce(amoebaHarmonicBondForce);
system.addForce(amoebaBondForce);
Context context(system, integrator, Platform::getPlatformByName( "Cuda"));
//Context context(system, integrator, platform );
std::vector<Vec3> positions(3);
......@@ -203,13 +203,13 @@ void testTwoBond( FILE* log ) {
positions[2] = Vec3(1, 0, 1);
context.setPositions(positions);
compareWithExpectedForceAndEnergy( context, *amoebaHarmonicBondForce, TOL, "testTwoBond", log );
compareWithExpectedForceAndEnergy( context, *amoebaBondForce, TOL, "testTwoBond", log );
}
int main( int numberOfArguments, char* argv[] ) {
try {
std::cout << "TestCudaAmoebaHarmonicBondForce running test..." << std::endl;
std::cout << "TestCudaAmoebaBondForce running test..." << std::endl;
registerAmoebaCudaKernelFactories();
FILE* log = NULL;
testTwoBond( log );
......
......@@ -277,20 +277,20 @@ static void setupAndGetForcesEnergyMultipoleAmmonia( AmoebaMultipoleForce::Polar
 
// 1-2 bonds needed
/*
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* AmoebaBondForce = new AmoebaBondForce();
 
// addBond: particle1, particle2, length, quadraticK
 
amoebaHarmonicBondForce->addBond( 0, 1, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 0, 2, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 0, 3, 0.0000000e+00, 0.0000000e+00 );
AmoebaBondForce->addBond( 0, 1, 0.0000000e+00, 0.0000000e+00 );
AmoebaBondForce->addBond( 0, 2, 0.0000000e+00, 0.0000000e+00 );
AmoebaBondForce->addBond( 0, 3, 0.0000000e+00, 0.0000000e+00 );
 
amoebaHarmonicBondForce->addBond( 4, 5, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 4, 6, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 4, 7, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( -2.5500000e+01 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( 3.7931250e+02 );
system.addForce(amoebaHarmonicBondForce);
AmoebaBondForce->addBond( 4, 5, 0.0000000e+00, 0.0000000e+00 );
AmoebaBondForce->addBond( 4, 6, 0.0000000e+00, 0.0000000e+00 );
AmoebaBondForce->addBond( 4, 7, 0.0000000e+00, 0.0000000e+00 );
AmoebaBondForce->setAmoebaGlobalBondCubic( -2.5500000e+01 );
AmoebaBondForce->setAmoebaGlobalBondQuartic( 3.7931250e+02 );
system.addForce(AmoebaBondForce);
*/
std::vector<Vec3> positions(numberOfParticles);
 
......
......@@ -30,7 +30,7 @@
* -------------------------------------------------------------------------- */
/**
* This tests the Cuda implementation of CudaAmoebaHarmonicInPlaneAngleForce.
* This tests the Cuda implementation of AmoebaInPlaneAngleForce.
*/
#include "openmm/internal/AssertionUtilities.h"
......@@ -120,22 +120,22 @@ static void getPrefactorsGivenInPlaneAngleCosine( double cosine, double idealInP
return;
}
static void computeAmoebaHarmonicInPlaneAngleForce(int bondIndex, std::vector<Vec3>& positions, AmoebaHarmonicInPlaneAngleForce& amoebaHarmonicInPlaneAngleForce,
static void computeAmoebaInPlaneAngleForce(int bondIndex, std::vector<Vec3>& positions, AmoebaInPlaneAngleForce& amoebaInPlaneAngleForce,
std::vector<Vec3>& forces, double* energy, FILE* log ) {
int particle1, particle2, particle3, particle4;
double idealInPlaneAngle;
double quadraticK;
amoebaHarmonicInPlaneAngleForce.getAngleParameters(bondIndex, particle1, particle2, particle3, particle4, idealInPlaneAngle, quadraticK );
amoebaInPlaneAngleForce.getAngleParameters(bondIndex, particle1, particle2, particle3, particle4, idealInPlaneAngle, quadraticK );
double cubicK = amoebaHarmonicInPlaneAngleForce.getAmoebaGlobalHarmonicInPlaneAngleCubic();
double quarticK = amoebaHarmonicInPlaneAngleForce.getAmoebaGlobalHarmonicInPlaneAngleQuartic();
double penticK = amoebaHarmonicInPlaneAngleForce.getAmoebaGlobalHarmonicInPlaneAnglePentic();
double sexticK = amoebaHarmonicInPlaneAngleForce.getAmoebaGlobalHarmonicInPlaneAngleSextic();
double cubicK = amoebaInPlaneAngleForce.getAmoebaGlobalInPlaneAngleCubic();
double quarticK = amoebaInPlaneAngleForce.getAmoebaGlobalInPlaneAngleQuartic();
double penticK = amoebaInPlaneAngleForce.getAmoebaGlobalInPlaneAnglePentic();
double sexticK = amoebaInPlaneAngleForce.getAmoebaGlobalInPlaneAngleSextic();
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicInPlaneAngleForce: bond %d [%d %d %d %d] ang=%10.3f k=%10.3f [%10.3e %10.3e %10.3e %10.3e]\n",
(void) fprintf( log, "computeAmoebaInPlaneAngleForce: bond %d [%d %d %d %d] ang=%10.3f k=%10.3f [%10.3e %10.3e %10.3e %10.3e]\n",
bondIndex, particle1, particle2, particle3, particle4, idealInPlaneAngle, quadraticK, cubicK, quarticK, penticK, sexticK );
(void) fflush( log );
}
......@@ -182,7 +182,7 @@ static void computeAmoebaHarmonicInPlaneAngleForce(int bondIndex, std::vector<V
if( rAp2 <= 0.0 && rCp2 <= 0.0 ){
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicInPlaneAngleForce: rAp2 or rCp2 <= 0.0\n" );
(void) fprintf( log, "computeAmoebaInPlaneAngleForce: rAp2 or rCp2 <= 0.0\n" );
(void) fflush( log );
}
#endif
......@@ -278,7 +278,7 @@ static void computeAmoebaHarmonicInPlaneAngleForce(int bondIndex, std::vector<V
}
static void computeAmoebaHarmonicInPlaneAngleForces( Context& context, AmoebaHarmonicInPlaneAngleForce& amoebaHarmonicInPlaneAngleForce,
static void computeAmoebaInPlaneAngleForces( Context& context, AmoebaInPlaneAngleForce& amoebaInPlaneAngleForce,
std::vector<Vec3>& expectedForces, double* expectedEnergy, FILE* log ) {
// get positions and zero forces
......@@ -294,13 +294,13 @@ static void computeAmoebaHarmonicInPlaneAngleForces( Context& context, AmoebaHar
// calculates forces/energy
*expectedEnergy = 0.0;
for( int ii = 0; ii < amoebaHarmonicInPlaneAngleForce.getNumAngles(); ii++ ){
computeAmoebaHarmonicInPlaneAngleForce(ii, positions, amoebaHarmonicInPlaneAngleForce, expectedForces, expectedEnergy, log );
for( int ii = 0; ii < amoebaInPlaneAngleForce.getNumAngles(); ii++ ){
computeAmoebaInPlaneAngleForce(ii, positions, amoebaInPlaneAngleForce, expectedForces, expectedEnergy, log );
}
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicInPlaneAngleForces: expected energy=%14.7e\n", *expectedEnergy );
(void) fprintf( log, "computeAmoebaInPlaneAngleForces: expected energy=%14.7e\n", *expectedEnergy );
for( unsigned int ii = 0; ii < positions.size(); ii++ ){
(void) fprintf( log, "%6u [%14.7e %14.7e %14.7e]\n", ii, expectedForces[ii][0], expectedForces[ii][1], expectedForces[ii][2] );
}
......@@ -311,19 +311,19 @@ static void computeAmoebaHarmonicInPlaneAngleForces( Context& context, AmoebaHar
}
void compareWithExpectedForceAndEnergy( Context& context, AmoebaHarmonicInPlaneAngleForce& amoebaHarmonicInPlaneAngleForce,
void compareWithExpectedForceAndEnergy( Context& context, AmoebaInPlaneAngleForce& amoebaInPlaneAngleForce,
double tolerance, const std::string& idString, FILE* log) {
std::vector<Vec3> expectedForces;
double expectedEnergy;
computeAmoebaHarmonicInPlaneAngleForces( context, amoebaHarmonicInPlaneAngleForce, expectedForces, &expectedEnergy, log );
computeAmoebaInPlaneAngleForces( context, amoebaInPlaneAngleForce, expectedForces, &expectedEnergy, log );
State state = context.getState(State::Forces | State::Energy);
const std::vector<Vec3> forces = state.getForces();
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicInPlaneAngleForces: expected energy=%14.7e %14.7e\n", expectedEnergy, state.getPotentialEnergy() );
(void) fprintf( log, "computeAmoebaInPlaneAngleForces: expected energy=%14.7e %14.7e\n", expectedEnergy, state.getPotentialEnergy() );
for( unsigned int ii = 0; ii < forces.size(); ii++ ){
(void) fprintf( log, "%6u [%14.7e %14.7e %14.7e] [%14.7e %14.7e %14.7e]\n", ii,
expectedForces[ii][0], expectedForces[ii][1], expectedForces[ii][2], forces[ii][0], forces[ii][1], forces[ii][2] );
......@@ -348,7 +348,7 @@ void testOneAngle( FILE* log ) {
LangevinIntegrator integrator(0.0, 0.1, 0.01);
AmoebaHarmonicInPlaneAngleForce* amoebaHarmonicInPlaneAngleForce = new AmoebaHarmonicInPlaneAngleForce();
AmoebaInPlaneAngleForce* amoebaInPlaneAngleForce = new AmoebaInPlaneAngleForce();
double angle = 65.0;
double quadraticK = 1.0;
......@@ -356,14 +356,14 @@ void testOneAngle( FILE* log ) {
double quarticK = 0.0e-02;
double penticK = 0.0e-03;
double sexticK = 0.0e-04;
amoebaHarmonicInPlaneAngleForce->addAngle(0, 1, 2, 3, angle, quadraticK);
amoebaInPlaneAngleForce->addAngle(0, 1, 2, 3, angle, quadraticK);
amoebaHarmonicInPlaneAngleForce->setAmoebaGlobalHarmonicInPlaneAngleCubic(cubicK);
amoebaHarmonicInPlaneAngleForce->setAmoebaGlobalHarmonicInPlaneAngleQuartic(quarticK);
amoebaHarmonicInPlaneAngleForce->setAmoebaGlobalHarmonicInPlaneAnglePentic(penticK);
amoebaHarmonicInPlaneAngleForce->setAmoebaGlobalHarmonicInPlaneAngleSextic(sexticK);
amoebaInPlaneAngleForce->setAmoebaGlobalInPlaneAngleCubic(cubicK);
amoebaInPlaneAngleForce->setAmoebaGlobalInPlaneAngleQuartic(quarticK);
amoebaInPlaneAngleForce->setAmoebaGlobalInPlaneAnglePentic(penticK);
amoebaInPlaneAngleForce->setAmoebaGlobalInPlaneAngleSextic(sexticK);
system.addForce(amoebaHarmonicInPlaneAngleForce);
system.addForce(amoebaInPlaneAngleForce);
Context context(system, integrator, Platform::getPlatformByName( "Cuda"));
std::vector<Vec3> positions(numberOfParticles);
......@@ -374,14 +374,14 @@ void testOneAngle( FILE* log ) {
positions[3] = Vec3(1, 1, 1);
context.setPositions(positions);
compareWithExpectedForceAndEnergy( context, *amoebaHarmonicInPlaneAngleForce, TOL, "testOneInPlaneAngle", log );
compareWithExpectedForceAndEnergy( context, *amoebaInPlaneAngleForce, TOL, "testOneInPlaneAngle", log );
}
int main( int numberOfArguments, char* argv[] ) {
try {
std::cout << "TestCudaAmoebaHarmonicInPlaneAngleForce running test..." << std::endl;
std::cout << "TestCudaAmoebaInPlaneAngleForce running test..." << std::endl;
registerAmoebaCudaKernelFactories();
FILE* log = NULL;
testOneAngle( NULL );
......
......@@ -274,20 +274,20 @@ static void setupAndGetForcesEnergyMultipoleAmmonia( AmoebaMultipoleForce::Nonbo
// 1-2 bonds needed
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* amoebaBondForce = new AmoebaBondForce();
// addBond: particle1, particle2, length, quadraticK
amoebaHarmonicBondForce->addBond( 0, 1, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 0, 2, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 0, 3, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( 0, 1, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( 0, 2, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( 0, 3, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 4, 5, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 4, 6, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 4, 7, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( -2.5500000e+01 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( 3.7931250e+02 );
system.addForce(amoebaHarmonicBondForce);
amoebaBondForce->addBond( 4, 5, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( 4, 6, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( 4, 7, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->setAmoebaGlobalBondCubic( -2.5500000e+01 );
amoebaBondForce->setAmoebaGlobalBondQuartic( 3.7931250e+02 );
system.addForce(amoebaBondForce);
std::vector<Vec3> positions(numberOfParticles);
......@@ -611,18 +611,18 @@ static void setupAndGetForcesEnergyMultipoleWater( AmoebaMultipoleForce::Nonbond
// 1-2 bonds needed
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* amoebaBondForce = new AmoebaBondForce();
// addBond: particle1, particle2, length, quadraticK
for( unsigned int jj = 0; jj < numberOfParticles; jj += 3 ){
amoebaHarmonicBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
}
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( -2.5500000e+01 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( 3.7931250e+02 );
system.addForce(amoebaHarmonicBondForce);
amoebaBondForce->setAmoebaGlobalBondCubic( -2.5500000e+01 );
amoebaBondForce->setAmoebaGlobalBondQuartic( 3.7931250e+02 );
system.addForce(amoebaBondForce);
std::vector<Vec3> positions(numberOfParticles);
......@@ -836,18 +836,18 @@ static void testQuadrupoleValidation( FILE* log ){
}
*/
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* amoebaBondForce = new AmoebaBondForce();
// addBond: particle1, particle2, length, quadraticK
for( unsigned int jj = 0; jj < numberOfParticles; jj += 3 ){
amoebaHarmonicBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
}
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( -2.5500000e+01 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( 3.7931250e+02 );
system.addForce(amoebaHarmonicBondForce);
amoebaBondForce->setAmoebaGlobalBondCubic( -2.5500000e+01 );
amoebaBondForce->setAmoebaGlobalBondQuartic( 3.7931250e+02 );
system.addForce(amoebaBondForce);
std::vector<Vec3> positions(numberOfParticles);
......@@ -1089,18 +1089,18 @@ static void setupAndGetForcesEnergyMultipoleIonsAndWater( AmoebaMultipoleForce::
// 1-2 bonds needed
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* amoebaBondForce = new AmoebaBondForce();
// addBond: particle1, particle2, length, quadraticK
for( unsigned int jj = 2; jj < numberOfParticles; jj += 3 ){
amoebaHarmonicBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
}
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( -2.5500000e+01 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( 3.7931250e+02 );
system.addForce(amoebaHarmonicBondForce);
amoebaBondForce->setAmoebaGlobalBondCubic( -2.5500000e+01 );
amoebaBondForce->setAmoebaGlobalBondQuartic( 3.7931250e+02 );
system.addForce(amoebaBondForce);
std::vector<Vec3> positions(numberOfParticles);
......@@ -1321,18 +1321,18 @@ static void setupAndGetForcesEnergyMultipoleLargeWater( AmoebaMultipoleForce::No
// 1-2 bonds needed
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* amoebaBondForce = new AmoebaBondForce();
// addBond: particle1, particle2, length, quadraticK
for( unsigned int jj = 0; jj < numberOfParticles; jj += 3 ){
amoebaHarmonicBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
}
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( 0.0 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( 0.0 );
system.addForce(amoebaHarmonicBondForce);
amoebaBondForce->setAmoebaGlobalBondCubic( 0.0 );
amoebaBondForce->setAmoebaGlobalBondQuartic( 0.0 );
system.addForce(amoebaBondForce);
std::vector<Vec3> positions(numberOfParticles);
......
......@@ -75,9 +75,9 @@ FOREACH(subdir ${OPENMM_SOURCE_SUBDIRS})
ENDFOREACH(subdir)
INCLUDE_DIRECTORIES(BEFORE ${CMAKE_CURRENT_SOURCE_DIR}/src)
INCLUDE_DIRECTORIES(BEFORE ${CMAKE_SOURCE_DIR}/platforms/cuda2/include)
INCLUDE_DIRECTORIES(BEFORE ${CMAKE_SOURCE_DIR}/platforms/cuda2/src)
INCLUDE_DIRECTORIES(BEFORE ${CMAKE_BINARY_DIR}/platforms/cuda2/src)
INCLUDE_DIRECTORIES(BEFORE ${CMAKE_SOURCE_DIR}/platforms/cuda/include)
INCLUDE_DIRECTORIES(BEFORE ${CMAKE_SOURCE_DIR}/platforms/cuda/src)
INCLUDE_DIRECTORIES(BEFORE ${CMAKE_BINARY_DIR}/platforms/cuda/src)
# Set variables needed for encoding kernel sources into a C++ class
......@@ -95,7 +95,7 @@ INCLUDE_DIRECTORIES(${CUDA_TOOLKIT_INCLUDE})
FILE(GLOB CUDA_KERNELS ${CUDA_SOURCE_DIR}/kernels/*.cu)
ADD_CUSTOM_COMMAND(OUTPUT ${CUDA_KERNELS_CPP} ${CUDA_KERNELS_H}
COMMAND ${CMAKE_COMMAND}
ARGS -D CUDA_SOURCE_DIR=${CUDA_SOURCE_DIR} -D CUDA_KERNELS_CPP=${CUDA_KERNELS_CPP} -D CUDA_KERNELS_H=${CUDA_KERNELS_H} -D CUDA_SOURCE_CLASS=${CUDA_SOURCE_CLASS} -P ${CMAKE_SOURCE_DIR}/platforms/cuda2/EncodeCUDAFiles.cmake
ARGS -D CUDA_SOURCE_DIR=${CUDA_SOURCE_DIR} -D CUDA_KERNELS_CPP=${CUDA_KERNELS_CPP} -D CUDA_KERNELS_H=${CUDA_KERNELS_H} -D CUDA_SOURCE_CLASS=${CUDA_SOURCE_CLASS} -P ${CMAKE_SOURCE_DIR}/platforms/cuda/EncodeCUDAFiles.cmake
DEPENDS ${CUDA_KERNELS}
)
SET_SOURCE_FILES_PROPERTIES(${CUDA_KERNELS_CPP} ${CUDA_KERNELS_H} PROPERTIES GENERATED TRUE)
......
......@@ -40,9 +40,9 @@ extern "C" OPENMM_EXPORT void registerKernelFactories() {
try {
Platform& platform = Platform::getPlatformByName("CUDA");
AmoebaCudaKernelFactory* factory = new AmoebaCudaKernelFactory();
platform.registerKernelFactory(CalcAmoebaHarmonicBondForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaHarmonicAngleForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaHarmonicInPlaneAngleForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaBondForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaAngleForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaInPlaneAngleForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaPiTorsionForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaStretchBendForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaOutOfPlaneBendForceKernel::Name(), factory);
......@@ -71,14 +71,14 @@ KernelImpl* AmoebaCudaKernelFactory::createKernelImpl(std::string name, const Pl
CudaPlatform::PlatformData& data = *static_cast<CudaPlatform::PlatformData*>(context.getPlatformData());
CudaContext& cu = *data.contexts[0];
if (name == CalcAmoebaHarmonicBondForceKernel::Name())
return new CudaCalcAmoebaHarmonicBondForceKernel(name, platform, cu, context.getSystem());
if (name == CalcAmoebaBondForceKernel::Name())
return new CudaCalcAmoebaBondForceKernel(name, platform, cu, context.getSystem());
if (name == CalcAmoebaHarmonicAngleForceKernel::Name())
return new CudaCalcAmoebaHarmonicAngleForceKernel(name, platform, cu, context.getSystem());
if (name == CalcAmoebaAngleForceKernel::Name())
return new CudaCalcAmoebaAngleForceKernel(name, platform, cu, context.getSystem());
if (name == CalcAmoebaHarmonicInPlaneAngleForceKernel::Name())
return new CudaCalcAmoebaHarmonicInPlaneAngleForceKernel(name, platform, cu, context.getSystem());
if (name == CalcAmoebaInPlaneAngleForceKernel::Name())
return new CudaCalcAmoebaInPlaneAngleForceKernel(name, platform, cu, context.getSystem());
if (name == CalcAmoebaPiTorsionForceKernel::Name())
return new CudaCalcAmoebaPiTorsionForceKernel(name, platform, cu, context.getSystem());
......
......@@ -54,12 +54,12 @@ using namespace std;
}
/* -------------------------------------------------------------------------- *
* AmoebaHarmonicBond *
* AmoebaBondForce *
* -------------------------------------------------------------------------- */
class CudaCalcAmoebaHarmonicBondForceKernel::ForceInfo : public CudaForceInfo {
class CudaCalcAmoebaBondForceKernel::ForceInfo : public CudaForceInfo {
public:
ForceInfo(const AmoebaHarmonicBondForce& force) : force(force) {
ForceInfo(const AmoebaBondForce& force) : force(force) {
}
int getNumParticleGroups() {
return force.getNumBonds();
......@@ -80,20 +80,20 @@ public:
return (length1 == length2 && k1 == k2);
}
private:
const AmoebaHarmonicBondForce& force;
const AmoebaBondForce& force;
};
CudaCalcAmoebaHarmonicBondForceKernel::CudaCalcAmoebaHarmonicBondForceKernel(std::string name, const Platform& platform, CudaContext& cu, System& system) :
CalcAmoebaHarmonicBondForceKernel(name, platform), cu(cu), system(system), params(NULL) {
CudaCalcAmoebaBondForceKernel::CudaCalcAmoebaBondForceKernel(std::string name, const Platform& platform, CudaContext& cu, System& system) :
CalcAmoebaBondForceKernel(name, platform), cu(cu), system(system), params(NULL) {
}
CudaCalcAmoebaHarmonicBondForceKernel::~CudaCalcAmoebaHarmonicBondForceKernel() {
CudaCalcAmoebaBondForceKernel::~CudaCalcAmoebaBondForceKernel() {
cu.setAsCurrent();
if (params != NULL)
delete params;
}
void CudaCalcAmoebaHarmonicBondForceKernel::initialize(const System& system, const AmoebaHarmonicBondForce& force) {
void CudaCalcAmoebaBondForceKernel::initialize(const System& system, const AmoebaBondForce& force) {
cu.setAsCurrent();
int numContexts = cu.getPlatformData().contexts.size();
int startIndex = cu.getContextIndex()*force.getNumBonds()/numContexts;
......@@ -113,23 +113,23 @@ void CudaCalcAmoebaHarmonicBondForceKernel::initialize(const System& system, con
map<string, string> replacements;
replacements["COMPUTE_FORCE"] = CudaAmoebaKernelSources::amoebaBondForce;
replacements["PARAMS"] = cu.getBondedUtilities().addArgument(params->getDevicePointer(), "float2");
replacements["CUBIC_K"] = cu.doubleToString(force.getAmoebaGlobalHarmonicBondCubic());
replacements["QUARTIC_K"] = cu.doubleToString(force.getAmoebaGlobalHarmonicBondQuartic());
replacements["CUBIC_K"] = cu.doubleToString(force.getAmoebaGlobalBondCubic());
replacements["QUARTIC_K"] = cu.doubleToString(force.getAmoebaGlobalBondQuartic());
cu.getBondedUtilities().addInteraction(atoms, cu.replaceStrings(CudaKernelSources::bondForce, replacements), force.getForceGroup());
cu.addForce(new ForceInfo(force));
}
double CudaCalcAmoebaHarmonicBondForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
double CudaCalcAmoebaBondForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
return 0.0;
}
/* -------------------------------------------------------------------------- *
* AmoebaHarmonicAngle *
* AmoebaAngleForce *
* -------------------------------------------------------------------------- */
class CudaCalcAmoebaHarmonicAngleForceKernel::ForceInfo : public CudaForceInfo {
class CudaCalcAmoebaAngleForceKernel::ForceInfo : public CudaForceInfo {
public:
ForceInfo(const AmoebaHarmonicAngleForce& force) : force(force) {
ForceInfo(const AmoebaAngleForce& force) : force(force) {
}
int getNumParticleGroups() {
return force.getNumAngles();
......@@ -151,20 +151,20 @@ public:
return (angle1 == angle2 && k1 == k2);
}
private:
const AmoebaHarmonicAngleForce& force;
const AmoebaAngleForce& force;
};
CudaCalcAmoebaHarmonicAngleForceKernel::CudaCalcAmoebaHarmonicAngleForceKernel(std::string name, const Platform& platform, CudaContext& cu, System& system) :
CalcAmoebaHarmonicAngleForceKernel(name, platform), cu(cu), system(system), params(NULL) {
CudaCalcAmoebaAngleForceKernel::CudaCalcAmoebaAngleForceKernel(std::string name, const Platform& platform, CudaContext& cu, System& system) :
CalcAmoebaAngleForceKernel(name, platform), cu(cu), system(system), params(NULL) {
}
CudaCalcAmoebaHarmonicAngleForceKernel::~CudaCalcAmoebaHarmonicAngleForceKernel() {
CudaCalcAmoebaAngleForceKernel::~CudaCalcAmoebaAngleForceKernel() {
cu.setAsCurrent();
if (params != NULL)
delete params;
}
void CudaCalcAmoebaHarmonicAngleForceKernel::initialize(const System& system, const AmoebaHarmonicAngleForce& force) {
void CudaCalcAmoebaAngleForceKernel::initialize(const System& system, const AmoebaAngleForce& force) {
cu.setAsCurrent();
int numContexts = cu.getPlatformData().contexts.size();
int startIndex = cu.getContextIndex()*force.getNumAngles()/numContexts;
......@@ -184,26 +184,26 @@ void CudaCalcAmoebaHarmonicAngleForceKernel::initialize(const System& system, co
map<string, string> replacements;
replacements["COMPUTE_FORCE"] = CudaAmoebaKernelSources::amoebaAngleForce;
replacements["PARAMS"] = cu.getBondedUtilities().addArgument(params->getDevicePointer(), "float2");
replacements["CUBIC_K"] = cu.doubleToString(force.getAmoebaGlobalHarmonicAngleCubic());
replacements["QUARTIC_K"] = cu.doubleToString(force.getAmoebaGlobalHarmonicAngleQuartic());
replacements["PENTIC_K"] = cu.doubleToString(force.getAmoebaGlobalHarmonicAnglePentic());
replacements["SEXTIC_K"] = cu.doubleToString(force.getAmoebaGlobalHarmonicAngleSextic());
replacements["CUBIC_K"] = cu.doubleToString(force.getAmoebaGlobalAngleCubic());
replacements["QUARTIC_K"] = cu.doubleToString(force.getAmoebaGlobalAngleQuartic());
replacements["PENTIC_K"] = cu.doubleToString(force.getAmoebaGlobalAnglePentic());
replacements["SEXTIC_K"] = cu.doubleToString(force.getAmoebaGlobalAngleSextic());
replacements["RAD_TO_DEG"] = cu.doubleToString(180/M_PI);
cu.getBondedUtilities().addInteraction(atoms, cu.replaceStrings(CudaKernelSources::angleForce, replacements), force.getForceGroup());
cu.addForce(new ForceInfo(force));
}
double CudaCalcAmoebaHarmonicAngleForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
double CudaCalcAmoebaAngleForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
return 0.0;
}
/* -------------------------------------------------------------------------- *
* AmoebaHarmonicInPlaneAngle *
* AmoebaInPlaneAngleForce *
* -------------------------------------------------------------------------- */
class CudaCalcAmoebaHarmonicInPlaneAngleForceKernel::ForceInfo : public CudaForceInfo {
class CudaCalcAmoebaInPlaneAngleForceKernel::ForceInfo : public CudaForceInfo {
public:
ForceInfo(const AmoebaHarmonicInPlaneAngleForce& force) : force(force) {
ForceInfo(const AmoebaInPlaneAngleForce& force) : force(force) {
}
int getNumParticleGroups() {
return force.getNumAngles();
......@@ -226,20 +226,20 @@ public:
return (angle1 == angle2 && k1 == k2);
}
private:
const AmoebaHarmonicInPlaneAngleForce& force;
const AmoebaInPlaneAngleForce& force;
};
CudaCalcAmoebaHarmonicInPlaneAngleForceKernel::CudaCalcAmoebaHarmonicInPlaneAngleForceKernel(std::string name, const Platform& platform, CudaContext& cu, System& system) :
CalcAmoebaHarmonicInPlaneAngleForceKernel(name, platform), cu(cu), system(system), params(NULL) {
CudaCalcAmoebaInPlaneAngleForceKernel::CudaCalcAmoebaInPlaneAngleForceKernel(std::string name, const Platform& platform, CudaContext& cu, System& system) :
CalcAmoebaInPlaneAngleForceKernel(name, platform), cu(cu), system(system), params(NULL) {
}
CudaCalcAmoebaHarmonicInPlaneAngleForceKernel::~CudaCalcAmoebaHarmonicInPlaneAngleForceKernel() {
CudaCalcAmoebaInPlaneAngleForceKernel::~CudaCalcAmoebaInPlaneAngleForceKernel() {
cu.setAsCurrent();
if (params != NULL)
delete params;
}
void CudaCalcAmoebaHarmonicInPlaneAngleForceKernel::initialize(const System& system, const AmoebaHarmonicInPlaneAngleForce& force) {
void CudaCalcAmoebaInPlaneAngleForceKernel::initialize(const System& system, const AmoebaInPlaneAngleForce& force) {
cu.setAsCurrent();
int numContexts = cu.getPlatformData().contexts.size();
int startIndex = cu.getContextIndex()*force.getNumAngles()/numContexts;
......@@ -258,16 +258,16 @@ void CudaCalcAmoebaHarmonicInPlaneAngleForceKernel::initialize(const System& sys
params->upload(paramVector);
map<string, string> replacements;
replacements["PARAMS"] = cu.getBondedUtilities().addArgument(params->getDevicePointer(), "float2");
replacements["CUBIC_K"] = cu.doubleToString(force.getAmoebaGlobalHarmonicInPlaneAngleCubic());
replacements["QUARTIC_K"] = cu.doubleToString(force.getAmoebaGlobalHarmonicInPlaneAngleQuartic());
replacements["PENTIC_K"] = cu.doubleToString(force.getAmoebaGlobalHarmonicInPlaneAnglePentic());
replacements["SEXTIC_K"] = cu.doubleToString(force.getAmoebaGlobalHarmonicInPlaneAngleSextic());
replacements["CUBIC_K"] = cu.doubleToString(force.getAmoebaGlobalInPlaneAngleCubic());
replacements["QUARTIC_K"] = cu.doubleToString(force.getAmoebaGlobalInPlaneAngleQuartic());
replacements["PENTIC_K"] = cu.doubleToString(force.getAmoebaGlobalInPlaneAnglePentic());
replacements["SEXTIC_K"] = cu.doubleToString(force.getAmoebaGlobalInPlaneAngleSextic());
replacements["RAD_TO_DEG"] = cu.doubleToString(180/M_PI);
cu.getBondedUtilities().addInteraction(atoms, cu.replaceStrings(CudaAmoebaKernelSources::amoebaInPlaneForce, replacements), force.getForceGroup());
cu.addForce(new ForceInfo(force));
}
double CudaCalcAmoebaHarmonicInPlaneAngleForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
double CudaCalcAmoebaInPlaneAngleForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
return 0.0;
}
......
......@@ -40,22 +40,22 @@ namespace OpenMM {
class CudaCalcAmoebaGeneralizedKirkwoodForceKernel;
/**
* This kernel is invoked by AmoebaHarmonicBondForce to calculate the forces acting on the system and the energy of the system.
* This kernel is invoked by AmoebaBondForce to calculate the forces acting on the system and the energy of the system.
*/
class CudaCalcAmoebaHarmonicBondForceKernel : public CalcAmoebaHarmonicBondForceKernel {
class CudaCalcAmoebaBondForceKernel : public CalcAmoebaBondForceKernel {
public:
CudaCalcAmoebaHarmonicBondForceKernel(std::string name,
CudaCalcAmoebaBondForceKernel(std::string name,
const Platform& platform,
CudaContext& cu,
System& system);
~CudaCalcAmoebaHarmonicBondForceKernel();
~CudaCalcAmoebaBondForceKernel();
/**
* Initialize the kernel.
*
* @param system the System this kernel will be applied to
* @param force the AmoebaHarmonicBondForce this kernel will be used for
* @param force the AmoebaBondForce this kernel will be used for
*/
void initialize(const System& system, const AmoebaHarmonicBondForce& force);
void initialize(const System& system, const AmoebaBondForce& force);
/**
* Execute the kernel to calculate the forces and/or energy.
*
......@@ -74,19 +74,19 @@ private:
};
/**
* This kernel is invoked by AmoebaHarmonicAngleForce to calculate the forces acting on the system and the energy of the system.
* This kernel is invoked by AmoebaAngleForce to calculate the forces acting on the system and the energy of the system.
*/
class CudaCalcAmoebaHarmonicAngleForceKernel : public CalcAmoebaHarmonicAngleForceKernel {
class CudaCalcAmoebaAngleForceKernel : public CalcAmoebaAngleForceKernel {
public:
CudaCalcAmoebaHarmonicAngleForceKernel(std::string name, const Platform& platform, CudaContext& cu, System& system);
~CudaCalcAmoebaHarmonicAngleForceKernel();
CudaCalcAmoebaAngleForceKernel(std::string name, const Platform& platform, CudaContext& cu, System& system);
~CudaCalcAmoebaAngleForceKernel();
/**
* Initialize the kernel.
*
* @param system the System this kernel will be applied to
* @param force the AmoebaHarmonicAngleForce this kernel will be used for
* @param force the AmoebaAngleForce this kernel will be used for
*/
void initialize(const System& system, const AmoebaHarmonicAngleForce& force);
void initialize(const System& system, const AmoebaAngleForce& force);
/**
* Execute the kernel to calculate the forces and/or energy.
*
......@@ -105,19 +105,19 @@ private:
};
/**
* This kernel is invoked by AmoebaHarmonicInPlaneAngleForce to calculate the forces acting on the system and the energy of the system.
* This kernel is invoked by AmoebaInPlaneAngleForce to calculate the forces acting on the system and the energy of the system.
*/
class CudaCalcAmoebaHarmonicInPlaneAngleForceKernel : public CalcAmoebaHarmonicInPlaneAngleForceKernel {
class CudaCalcAmoebaInPlaneAngleForceKernel : public CalcAmoebaInPlaneAngleForceKernel {
public:
CudaCalcAmoebaHarmonicInPlaneAngleForceKernel(std::string name, const Platform& platform, CudaContext& cu, System& system);
~CudaCalcAmoebaHarmonicInPlaneAngleForceKernel();
CudaCalcAmoebaInPlaneAngleForceKernel(std::string name, const Platform& platform, CudaContext& cu, System& system);
~CudaCalcAmoebaInPlaneAngleForceKernel();
/**
* Initialize the kernel.
*
* @param system the System this kernel will be applied to
* @param force the AmoebaHarmonicInPlaneAngleForce this kernel will be used for
* @param force the AmoebaInPlaneAngleForce this kernel will be used for
*/
void initialize(const System& system, const AmoebaHarmonicInPlaneAngleForce& force);
void initialize(const System& system, const AmoebaInPlaneAngleForce& force);
/**
* Execute the kernel to calculate the forces and/or energy.
*
......
......@@ -30,7 +30,7 @@
* -------------------------------------------------------------------------- */
/**
* This tests the CUDA implementation of CudaAmoebaHarmonicAngleForce.
* This tests the CUDA implementation of CudaAmoebaAngleForce.
*/
#include "openmm/internal/AssertionUtilities.h"
......@@ -120,22 +120,22 @@ static void getPrefactorsGivenAngleCosine( double cosine, double idealAngle, dou
return;
}
static void computeAmoebaHarmonicAngleForce(int bondIndex, std::vector<Vec3>& positions, AmoebaHarmonicAngleForce& amoebaHarmonicAngleForce,
static void computeAmoebaAngleForce(int bondIndex, std::vector<Vec3>& positions, AmoebaAngleForce& amoebaAngleForce,
std::vector<Vec3>& forces, double* energy, FILE* log ) {
int particle1, particle2, particle3;
double idealAngle;
double quadraticK;
amoebaHarmonicAngleForce.getAngleParameters(bondIndex, particle1, particle2, particle3, idealAngle, quadraticK );
amoebaAngleForce.getAngleParameters(bondIndex, particle1, particle2, particle3, idealAngle, quadraticK );
double cubicK = amoebaHarmonicAngleForce.getAmoebaGlobalHarmonicAngleCubic();
double quarticK = amoebaHarmonicAngleForce.getAmoebaGlobalHarmonicAngleQuartic();
double penticK = amoebaHarmonicAngleForce.getAmoebaGlobalHarmonicAnglePentic();
double sexticK = amoebaHarmonicAngleForce.getAmoebaGlobalHarmonicAngleSextic();
double cubicK = amoebaAngleForce.getAmoebaGlobalAngleCubic();
double quarticK = amoebaAngleForce.getAmoebaGlobalAngleQuartic();
double penticK = amoebaAngleForce.getAmoebaGlobalAnglePentic();
double sexticK = amoebaAngleForce.getAmoebaGlobalAngleSextic();
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicAngleForce: bond %d [%d %d %d] ang=%10.3f k=%10.3f [%10.3e %10.3e %10.3e %10.3e]\n",
(void) fprintf( log, "computeAmoebaAngleForce: bond %d [%d %d %d] ang=%10.3f k=%10.3f [%10.3e %10.3e %10.3e %10.3e]\n",
bondIndex, particle1, particle2, particle3, idealAngle, quadraticK, cubicK, quarticK, penticK, sexticK );
(void) fflush( log );
}
......@@ -202,7 +202,7 @@ static void computeAmoebaHarmonicAngleForce(int bondIndex, std::vector<Vec3>& p
*energy += energyTerm;
}
static void computeAmoebaHarmonicAngleForces( Context& context, AmoebaHarmonicAngleForce& amoebaHarmonicAngleForce,
static void computeAmoebaAngleForces( Context& context, AmoebaAngleForce& amoebaAngleForce,
std::vector<Vec3>& expectedForces, double* expectedEnergy, FILE* log ) {
// get positions and zero forces
......@@ -218,13 +218,13 @@ static void computeAmoebaHarmonicAngleForces( Context& context, AmoebaHarmonicAn
// calculates forces/energy
*expectedEnergy = 0.0;
for( int ii = 0; ii < amoebaHarmonicAngleForce.getNumAngles(); ii++ ){
computeAmoebaHarmonicAngleForce(ii, positions, amoebaHarmonicAngleForce, expectedForces, expectedEnergy, log );
for( int ii = 0; ii < amoebaAngleForce.getNumAngles(); ii++ ){
computeAmoebaAngleForce(ii, positions, amoebaAngleForce, expectedForces, expectedEnergy, log );
}
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicAngleForces: expected energy=%14.7e\n", *expectedEnergy );
(void) fprintf( log, "computeAmoebaAngleForces: expected energy=%14.7e\n", *expectedEnergy );
for( unsigned int ii = 0; ii < positions.size(); ii++ ){
(void) fprintf( log, "%6u [%14.7e %14.7e %14.7e]\n", ii, expectedForces[ii][0], expectedForces[ii][1], expectedForces[ii][2] );
}
......@@ -236,19 +236,19 @@ static void computeAmoebaHarmonicAngleForces( Context& context, AmoebaHarmonicAn
}
void compareWithExpectedForceAndEnergy( Context& context, AmoebaHarmonicAngleForce& amoebaHarmonicAngleForce,
void compareWithExpectedForceAndEnergy( Context& context, AmoebaAngleForce& amoebaAngleForce,
double tolerance, const std::string& idString, FILE* log) {
std::vector<Vec3> expectedForces;
double expectedEnergy;
computeAmoebaHarmonicAngleForces( context, amoebaHarmonicAngleForce, expectedForces, &expectedEnergy, log );
computeAmoebaAngleForces( context, amoebaAngleForce, expectedForces, &expectedEnergy, log );
State state = context.getState(State::Forces | State::Energy);
const std::vector<Vec3> forces = state.getForces();
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicAngleForces: expected energy=%14.7e %14.7e\n", expectedEnergy, state.getPotentialEnergy() );
(void) fprintf( log, "computeAmoebaAngleForces: expected energy=%14.7e %14.7e\n", expectedEnergy, state.getPotentialEnergy() );
for( unsigned int ii = 0; ii < forces.size(); ii++ ){
(void) fprintf( log, "%6u [%14.7e %14.7e %14.7e] [%14.7e %14.7e %14.7e]\n", ii,
expectedForces[ii][0], expectedForces[ii][1], expectedForces[ii][2], forces[ii][0], forces[ii][1], forces[ii][2] );
......@@ -273,7 +273,7 @@ void testOneAngle( FILE* log ) {
LangevinIntegrator integrator(0.0, 0.1, 0.01);
AmoebaHarmonicAngleForce* amoebaHarmonicAngleForce = new AmoebaHarmonicAngleForce();
AmoebaAngleForce* amoebaAngleForce = new AmoebaAngleForce();
double angle = 100.0;
double quadraticK = 1.0;
......@@ -281,14 +281,14 @@ void testOneAngle( FILE* log ) {
double quarticK = 1.0e-02;
double penticK = 1.0e-03;
double sexticK = 1.0e-04;
amoebaHarmonicAngleForce->addAngle(0, 1, 2, angle, quadraticK);
amoebaAngleForce->addAngle(0, 1, 2, angle, quadraticK);
amoebaHarmonicAngleForce->setAmoebaGlobalHarmonicAngleCubic(cubicK);
amoebaHarmonicAngleForce->setAmoebaGlobalHarmonicAngleQuartic(quarticK);
amoebaHarmonicAngleForce->setAmoebaGlobalHarmonicAnglePentic(penticK);
amoebaHarmonicAngleForce->setAmoebaGlobalHarmonicAngleSextic(sexticK);
amoebaAngleForce->setAmoebaGlobalAngleCubic(cubicK);
amoebaAngleForce->setAmoebaGlobalAngleQuartic(quarticK);
amoebaAngleForce->setAmoebaGlobalAnglePentic(penticK);
amoebaAngleForce->setAmoebaGlobalAngleSextic(sexticK);
system.addForce(amoebaHarmonicAngleForce);
system.addForce(amoebaAngleForce);
Context context(system, integrator, Platform::getPlatformByName( "CUDA"));
std::vector<Vec3> positions(numberOfParticles);
......@@ -298,14 +298,14 @@ void testOneAngle( FILE* log ) {
positions[2] = Vec3(0, 0, 1);
context.setPositions(positions);
compareWithExpectedForceAndEnergy( context, *amoebaHarmonicAngleForce, TOL, "testOneAngle", log );
compareWithExpectedForceAndEnergy( context, *amoebaAngleForce, TOL, "testOneAngle", log );
}
int main( int numberOfArguments, char* argv[] ) {
try {
std::cout << "TestCudaAmoebaHarmonicAngleForce running test..." << std::endl;
std::cout << "TestCudaAmoebaAngleForce running test..." << std::endl;
registerAmoebaCudaKernelFactories();
FILE* log = NULL;
testOneAngle( log );
......
......@@ -30,7 +30,7 @@
* -------------------------------------------------------------------------- */
/**
* This tests the Cuda implementation of HarmonicBondForce.
* This tests the Cuda implementation of AmoebaBondForce.
*/
#include "openmm/internal/AssertionUtilities.h"
......@@ -47,15 +47,15 @@ using namespace OpenMM;
extern "C" void registerAmoebaCudaKernelFactories();
const double TOL = 1e-5;
static void computeAmoebaHarmonicBondForce(int bondIndex, std::vector<Vec3>& positions, AmoebaHarmonicBondForce& amoebaHarmonicBondForce,
static void computeAmoebaBondForce(int bondIndex, std::vector<Vec3>& positions, AmoebaBondForce& amoebaBondForce,
std::vector<Vec3>& forces, double* energy ) {
int particle1, particle2;
double bondLength;
double quadraticK;
double cubicK = amoebaHarmonicBondForce.getAmoebaGlobalHarmonicBondCubic();
double quarticK = amoebaHarmonicBondForce.getAmoebaGlobalHarmonicBondQuartic();
amoebaHarmonicBondForce.getBondParameters(bondIndex, particle1, particle2, bondLength, quadraticK );
double cubicK = amoebaBondForce.getAmoebaGlobalBondCubic();
double quarticK = amoebaBondForce.getAmoebaGlobalBondQuartic();
amoebaBondForce.getBondParameters(bondIndex, particle1, particle2, bondLength, quadraticK );
double deltaR[3];
double r2 = 0.0;
......@@ -83,7 +83,7 @@ static void computeAmoebaHarmonicBondForce(int bondIndex, std::vector<Vec3>& po
}
static void computeAmoebaHarmonicBondForces( Context& context, AmoebaHarmonicBondForce& amoebaHarmonicBondForce,
static void computeAmoebaBondForces( Context& context, AmoebaBondForce& amoebaBondForce,
std::vector<Vec3>& expectedForces, double* expectedEnergy, FILE* log ) {
// get positions and zero forces
......@@ -99,13 +99,13 @@ static void computeAmoebaHarmonicBondForces( Context& context, AmoebaHarmonicBon
// calculates forces/energy
*expectedEnergy = 0.0;
for( int ii = 0; ii < amoebaHarmonicBondForce.getNumBonds(); ii++ ){
computeAmoebaHarmonicBondForce(ii, positions, amoebaHarmonicBondForce, expectedForces, expectedEnergy );
for( int ii = 0; ii < amoebaBondForce.getNumBonds(); ii++ ){
computeAmoebaBondForce(ii, positions, amoebaBondForce, expectedForces, expectedEnergy );
}
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicBondForces: expected energy=%15.7e\n", *expectedEnergy );
(void) fprintf( log, "computeAmoebaBondForces: expected energy=%15.7e\n", *expectedEnergy );
for( unsigned int ii = 0; ii < positions.size(); ii++ ){
(void) fprintf( log, "%6u [%15.7e %15.7e %15.7e]\n", ii, expectedForces[ii][0], expectedForces[ii][1], expectedForces[ii][2] );
}
......@@ -116,18 +116,18 @@ static void computeAmoebaHarmonicBondForces( Context& context, AmoebaHarmonicBon
}
void compareWithExpectedForceAndEnergy( Context& context, AmoebaHarmonicBondForce& amoebaHarmonicBondForce, double tolerance, const std::string& idString, FILE* log) {
void compareWithExpectedForceAndEnergy( Context& context, AmoebaBondForce& amoebaBondForce, double tolerance, const std::string& idString, FILE* log) {
std::vector<Vec3> expectedForces;
double expectedEnergy;
computeAmoebaHarmonicBondForces( context, amoebaHarmonicBondForce, expectedForces, &expectedEnergy, NULL );
computeAmoebaBondForces( context, amoebaBondForce, expectedForces, &expectedEnergy, NULL );
State state = context.getState(State::Forces | State::Energy);
const std::vector<Vec3> forces = state.getForces();
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicBondForces: expected energy=%15.7e %15.7e\n", expectedEnergy, state.getPotentialEnergy() );
(void) fprintf( log, "computeAmoebaBondForces: expected energy=%15.7e %15.7e\n", expectedEnergy, state.getPotentialEnergy() );
for( unsigned int ii = 0; ii < forces.size(); ii++ ){
(void) fprintf( log, "%6u [%15.7e %15.7e %15.7e] [%15.7e %15.7e %15.7e]\n", ii,
expectedForces[ii][0], expectedForces[ii][1], expectedForces[ii][2], forces[ii][0], forces[ii][1], forces[ii][2] );
......@@ -151,17 +151,17 @@ void testOneBond( FILE* log ) {
LangevinIntegrator integrator(0.0, 0.1, 0.01);
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* amoebaBondForce = new AmoebaBondForce();
double bondLength = 1.5;
double quadraticK = 1.0;
double cubicK = 2.0;
double quarticicK = 3.0;
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( cubicK );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( quarticicK );
amoebaHarmonicBondForce->addBond(0, 1, bondLength, quadraticK);
amoebaBondForce->setAmoebaGlobalBondCubic( cubicK );
amoebaBondForce->setAmoebaGlobalBondQuartic( quarticicK );
amoebaBondForce->addBond(0, 1, bondLength, quadraticK);
system.addForce(amoebaHarmonicBondForce);
system.addForce(amoebaBondForce);
Context context(system, integrator, Platform::getPlatformByName( "CUDA"));
std::vector<Vec3> positions(2);
......@@ -169,7 +169,7 @@ void testOneBond( FILE* log ) {
positions[1] = Vec3(0, 0, 0);
context.setPositions(positions);
compareWithExpectedForceAndEnergy( context, *amoebaHarmonicBondForce, TOL, "testOneBond", log );
compareWithExpectedForceAndEnergy( context, *amoebaBondForce, TOL, "testOneBond", log );
}
void testTwoBond( FILE* log ) {
......@@ -182,18 +182,18 @@ void testTwoBond( FILE* log ) {
LangevinIntegrator integrator(0.0, 0.1, 0.01);
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* amoebaBondForce = new AmoebaBondForce();
double bondLength = 1.5;
double quadraticK = 1.0;
double cubicK = 2.0;
double quarticicK = 3.0;
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( cubicK );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( quarticicK );
amoebaHarmonicBondForce->addBond(0, 1, bondLength, quadraticK);
amoebaHarmonicBondForce->addBond(1, 2, bondLength, quadraticK);
amoebaBondForce->setAmoebaGlobalBondCubic( cubicK );
amoebaBondForce->setAmoebaGlobalBondQuartic( quarticicK );
amoebaBondForce->addBond(0, 1, bondLength, quadraticK);
amoebaBondForce->addBond(1, 2, bondLength, quadraticK);
system.addForce(amoebaHarmonicBondForce);
system.addForce(amoebaBondForce);
Context context(system, integrator, Platform::getPlatformByName( "CUDA"));
//Context context(system, integrator, platform );
std::vector<Vec3> positions(3);
......@@ -203,13 +203,13 @@ void testTwoBond( FILE* log ) {
positions[2] = Vec3(1, 0, 1);
context.setPositions(positions);
compareWithExpectedForceAndEnergy( context, *amoebaHarmonicBondForce, TOL, "testTwoBond", log );
compareWithExpectedForceAndEnergy( context, *amoebaBondForce, TOL, "testTwoBond", log );
}
int main( int numberOfArguments, char* argv[] ) {
try {
std::cout << "TestCudaAmoebaHarmonicBondForce running test..." << std::endl;
std::cout << "TestCudaAmoebaBondForce running test..." << std::endl;
registerAmoebaCudaKernelFactories();
FILE* log = NULL;
testTwoBond( log );
......
......@@ -278,20 +278,20 @@ static void setupAndGetForcesEnergyMultipoleAmmonia( AmoebaMultipoleForce::Polar
 
// 1-2 bonds needed
/*
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* AmoebaBondForce = new AmoebaBondForce();
 
// addBond: particle1, particle2, length, quadraticK
 
amoebaHarmonicBondForce->addBond( 0, 1, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 0, 2, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 0, 3, 0.0000000e+00, 0.0000000e+00 );
AmoebaBondForce->addBond( 0, 1, 0.0000000e+00, 0.0000000e+00 );
AmoebaBondForce->addBond( 0, 2, 0.0000000e+00, 0.0000000e+00 );
AmoebaBondForce->addBond( 0, 3, 0.0000000e+00, 0.0000000e+00 );
 
amoebaHarmonicBondForce->addBond( 4, 5, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 4, 6, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 4, 7, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( -2.5500000e+01 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( 3.7931250e+02 );
system.addForce(amoebaHarmonicBondForce);
AmoebaBondForce->addBond( 4, 5, 0.0000000e+00, 0.0000000e+00 );
AmoebaBondForce->addBond( 4, 6, 0.0000000e+00, 0.0000000e+00 );
AmoebaBondForce->addBond( 4, 7, 0.0000000e+00, 0.0000000e+00 );
AmoebaBondForce->setAmoebaGlobalBondCubic( -2.5500000e+01 );
AmoebaBondForce->setAmoebaGlobalBondQuartic( 3.7931250e+02 );
system.addForce(AmoebaBondForce);
*/
std::vector<Vec3> positions(numberOfParticles);
 
......
......@@ -30,7 +30,7 @@
* -------------------------------------------------------------------------- */
/**
* This tests the CUDA implementation of CudaAmoebaHarmonicInPlaneAngleForce.
* This tests the CUDA implementation of AmoebaInPlaneAngleForce.
*/
#include "openmm/internal/AssertionUtilities.h"
......@@ -121,22 +121,22 @@ static void getPrefactorsGivenInPlaneAngleCosine( double cosine, double idealInP
return;
}
static void computeAmoebaHarmonicInPlaneAngleForce(int bondIndex, std::vector<Vec3>& positions, AmoebaHarmonicInPlaneAngleForce& amoebaHarmonicInPlaneAngleForce,
static void computeAmoebaInPlaneAngleForce(int bondIndex, std::vector<Vec3>& positions, AmoebaInPlaneAngleForce& amoebaInPlaneAngleForce,
std::vector<Vec3>& forces, double* energy, FILE* log ) {
int particle1, particle2, particle3, particle4;
double idealInPlaneAngle;
double quadraticK;
amoebaHarmonicInPlaneAngleForce.getAngleParameters(bondIndex, particle1, particle2, particle3, particle4, idealInPlaneAngle, quadraticK );
amoebaInPlaneAngleForce.getAngleParameters(bondIndex, particle1, particle2, particle3, particle4, idealInPlaneAngle, quadraticK );
double cubicK = amoebaHarmonicInPlaneAngleForce.getAmoebaGlobalHarmonicInPlaneAngleCubic();
double quarticK = amoebaHarmonicInPlaneAngleForce.getAmoebaGlobalHarmonicInPlaneAngleQuartic();
double penticK = amoebaHarmonicInPlaneAngleForce.getAmoebaGlobalHarmonicInPlaneAnglePentic();
double sexticK = amoebaHarmonicInPlaneAngleForce.getAmoebaGlobalHarmonicInPlaneAngleSextic();
double cubicK = amoebaInPlaneAngleForce.getAmoebaGlobalInPlaneAngleCubic();
double quarticK = amoebaInPlaneAngleForce.getAmoebaGlobalInPlaneAngleQuartic();
double penticK = amoebaInPlaneAngleForce.getAmoebaGlobalInPlaneAnglePentic();
double sexticK = amoebaInPlaneAngleForce.getAmoebaGlobalInPlaneAngleSextic();
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicInPlaneAngleForce: bond %d [%d %d %d %d] ang=%10.3f k=%10.3f [%10.3e %10.3e %10.3e %10.3e]\n",
(void) fprintf( log, "computeAmoebaInPlaneAngleForce: bond %d [%d %d %d %d] ang=%10.3f k=%10.3f [%10.3e %10.3e %10.3e %10.3e]\n",
bondIndex, particle1, particle2, particle3, particle4, idealInPlaneAngle, quadraticK, cubicK, quarticK, penticK, sexticK );
(void) fflush( log );
}
......@@ -183,7 +183,7 @@ static void computeAmoebaHarmonicInPlaneAngleForce(int bondIndex, std::vector<V
if( rAp2 <= 0.0 && rCp2 <= 0.0 ){
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicInPlaneAngleForce: rAp2 or rCp2 <= 0.0\n" );
(void) fprintf( log, "computeAmoebaInPlaneAngleForce: rAp2 or rCp2 <= 0.0\n" );
(void) fflush( log );
}
#endif
......@@ -279,7 +279,7 @@ static void computeAmoebaHarmonicInPlaneAngleForce(int bondIndex, std::vector<V
}
static void computeAmoebaHarmonicInPlaneAngleForces( Context& context, AmoebaHarmonicInPlaneAngleForce& amoebaHarmonicInPlaneAngleForce,
static void computeAmoebaInPlaneAngleForces( Context& context, AmoebaInPlaneAngleForce& amoebaInPlaneAngleForce,
std::vector<Vec3>& expectedForces, double* expectedEnergy, FILE* log ) {
// get positions and zero forces
......@@ -295,13 +295,13 @@ static void computeAmoebaHarmonicInPlaneAngleForces( Context& context, AmoebaHar
// calculates forces/energy
*expectedEnergy = 0.0;
for( int ii = 0; ii < amoebaHarmonicInPlaneAngleForce.getNumAngles(); ii++ ){
computeAmoebaHarmonicInPlaneAngleForce(ii, positions, amoebaHarmonicInPlaneAngleForce, expectedForces, expectedEnergy, log );
for( int ii = 0; ii < amoebaInPlaneAngleForce.getNumAngles(); ii++ ){
computeAmoebaInPlaneAngleForce(ii, positions, amoebaInPlaneAngleForce, expectedForces, expectedEnergy, log );
}
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicInPlaneAngleForces: expected energy=%14.7e\n", *expectedEnergy );
(void) fprintf( log, "computeAmoebaInPlaneAngleForces: expected energy=%14.7e\n", *expectedEnergy );
for( unsigned int ii = 0; ii < positions.size(); ii++ ){
(void) fprintf( log, "%6u [%14.7e %14.7e %14.7e]\n", ii, expectedForces[ii][0], expectedForces[ii][1], expectedForces[ii][2] );
}
......@@ -312,19 +312,19 @@ static void computeAmoebaHarmonicInPlaneAngleForces( Context& context, AmoebaHar
}
void compareWithExpectedForceAndEnergy( Context& context, AmoebaHarmonicInPlaneAngleForce& amoebaHarmonicInPlaneAngleForce,
void compareWithExpectedForceAndEnergy( Context& context, AmoebaInPlaneAngleForce& amoebaInPlaneAngleForce,
double tolerance, const std::string& idString, FILE* log) {
std::vector<Vec3> expectedForces;
double expectedEnergy;
computeAmoebaHarmonicInPlaneAngleForces( context, amoebaHarmonicInPlaneAngleForce, expectedForces, &expectedEnergy, log );
computeAmoebaInPlaneAngleForces( context, amoebaInPlaneAngleForce, expectedForces, &expectedEnergy, log );
State state = context.getState(State::Forces | State::Energy);
const std::vector<Vec3> forces = state.getForces();
#ifdef AMOEBA_DEBUG
if( log ){
(void) fprintf( log, "computeAmoebaHarmonicInPlaneAngleForces: expected energy=%14.7e %14.7e\n", expectedEnergy, state.getPotentialEnergy() );
(void) fprintf( log, "computeAmoebaInPlaneAngleForces: expected energy=%14.7e %14.7e\n", expectedEnergy, state.getPotentialEnergy() );
for( unsigned int ii = 0; ii < forces.size(); ii++ ){
(void) fprintf( log, "%6u [%14.7e %14.7e %14.7e] [%14.7e %14.7e %14.7e]\n", ii,
expectedForces[ii][0], expectedForces[ii][1], expectedForces[ii][2], forces[ii][0], forces[ii][1], forces[ii][2] );
......@@ -349,7 +349,7 @@ void testOneAngle( FILE* log ) {
LangevinIntegrator integrator(0.0, 0.1, 0.01);
AmoebaHarmonicInPlaneAngleForce* amoebaHarmonicInPlaneAngleForce = new AmoebaHarmonicInPlaneAngleForce();
AmoebaInPlaneAngleForce* amoebaInPlaneAngleForce = new AmoebaInPlaneAngleForce();
double angle = 65.0;
double quadraticK = 1.0;
......@@ -357,14 +357,14 @@ void testOneAngle( FILE* log ) {
double quarticK = 0.0e-02;
double penticK = 0.0e-03;
double sexticK = 0.0e-04;
amoebaHarmonicInPlaneAngleForce->addAngle(0, 1, 2, 3, angle, quadraticK);
amoebaInPlaneAngleForce->addAngle(0, 1, 2, 3, angle, quadraticK);
amoebaHarmonicInPlaneAngleForce->setAmoebaGlobalHarmonicInPlaneAngleCubic(cubicK);
amoebaHarmonicInPlaneAngleForce->setAmoebaGlobalHarmonicInPlaneAngleQuartic(quarticK);
amoebaHarmonicInPlaneAngleForce->setAmoebaGlobalHarmonicInPlaneAnglePentic(penticK);
amoebaHarmonicInPlaneAngleForce->setAmoebaGlobalHarmonicInPlaneAngleSextic(sexticK);
amoebaInPlaneAngleForce->setAmoebaGlobalInPlaneAngleCubic(cubicK);
amoebaInPlaneAngleForce->setAmoebaGlobalInPlaneAngleQuartic(quarticK);
amoebaInPlaneAngleForce->setAmoebaGlobalInPlaneAnglePentic(penticK);
amoebaInPlaneAngleForce->setAmoebaGlobalInPlaneAngleSextic(sexticK);
system.addForce(amoebaHarmonicInPlaneAngleForce);
system.addForce(amoebaInPlaneAngleForce);
Context context(system, integrator, Platform::getPlatformByName( "CUDA"));
std::vector<Vec3> positions(numberOfParticles);
......@@ -375,14 +375,14 @@ void testOneAngle( FILE* log ) {
positions[3] = Vec3(1, 1, 1);
context.setPositions(positions);
compareWithExpectedForceAndEnergy( context, *amoebaHarmonicInPlaneAngleForce, TOL, "testOneInPlaneAngle", log );
compareWithExpectedForceAndEnergy( context, *amoebaInPlaneAngleForce, TOL, "testOneInPlaneAngle", log );
}
int main( int numberOfArguments, char* argv[] ) {
try {
std::cout << "TestCudaAmoebaHarmonicInPlaneAngleForce running test..." << std::endl;
std::cout << "TestCudaAmoebaInPlaneAngleForce running test..." << std::endl;
registerAmoebaCudaKernelFactories();
FILE* log = NULL;
testOneAngle( NULL );
......
......@@ -275,20 +275,20 @@ static void setupAndGetForcesEnergyMultipoleAmmonia( AmoebaMultipoleForce::Nonbo
// 1-2 bonds needed
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* amoebaBondForce = new AmoebaBondForce();
// addBond: particle1, particle2, length, quadraticK
amoebaHarmonicBondForce->addBond( 0, 1, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 0, 2, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 0, 3, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( 0, 1, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( 0, 2, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( 0, 3, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 4, 5, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 4, 6, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( 4, 7, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( -2.5500000e+01 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( 3.7931250e+02 );
system.addForce(amoebaHarmonicBondForce);
amoebaBondForce->addBond( 4, 5, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( 4, 6, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( 4, 7, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->setAmoebaGlobalBondCubic( -2.5500000e+01 );
amoebaBondForce->setAmoebaGlobalBondQuartic( 3.7931250e+02 );
system.addForce(amoebaBondForce);
std::vector<Vec3> positions(numberOfParticles);
......@@ -612,18 +612,18 @@ static void setupAndGetForcesEnergyMultipoleWater( AmoebaMultipoleForce::Nonbond
// 1-2 bonds needed
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* amoebaBondForce = new AmoebaBondForce();
// addBond: particle1, particle2, length, quadraticK
for( unsigned int jj = 0; jj < numberOfParticles; jj += 3 ){
amoebaHarmonicBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
}
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( -2.5500000e+01 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( 3.7931250e+02 );
system.addForce(amoebaHarmonicBondForce);
amoebaBondForce->setAmoebaGlobalBondCubic( -2.5500000e+01 );
amoebaBondForce->setAmoebaGlobalBondQuartic( 3.7931250e+02 );
system.addForce(amoebaBondForce);
std::vector<Vec3> positions(numberOfParticles);
......@@ -837,18 +837,18 @@ static void testQuadrupoleValidation( FILE* log ){
}
*/
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* amoebaBondForce = new AmoebaBondForce();
// addBond: particle1, particle2, length, quadraticK
for( unsigned int jj = 0; jj < numberOfParticles; jj += 3 ){
amoebaHarmonicBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
}
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( -2.5500000e+01 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( 3.7931250e+02 );
system.addForce(amoebaHarmonicBondForce);
amoebaBondForce->setAmoebaGlobalBondCubic( -2.5500000e+01 );
amoebaBondForce->setAmoebaGlobalBondQuartic( 3.7931250e+02 );
system.addForce(amoebaBondForce);
std::vector<Vec3> positions(numberOfParticles);
......@@ -1090,18 +1090,18 @@ static void setupAndGetForcesEnergyMultipoleIonsAndWater( AmoebaMultipoleForce::
// 1-2 bonds needed
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* amoebaBondForce = new AmoebaBondForce();
// addBond: particle1, particle2, length, quadraticK
for( unsigned int jj = 2; jj < numberOfParticles; jj += 3 ){
amoebaHarmonicBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
}
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( -2.5500000e+01 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( 3.7931250e+02 );
system.addForce(amoebaHarmonicBondForce);
amoebaBondForce->setAmoebaGlobalBondCubic( -2.5500000e+01 );
amoebaBondForce->setAmoebaGlobalBondQuartic( 3.7931250e+02 );
system.addForce(amoebaBondForce);
std::vector<Vec3> positions(numberOfParticles);
......@@ -1322,18 +1322,18 @@ static void setupAndGetForcesEnergyMultipoleLargeWater( AmoebaMultipoleForce::No
// 1-2 bonds needed
AmoebaHarmonicBondForce* amoebaHarmonicBondForce = new AmoebaHarmonicBondForce();
AmoebaBondForce* amoebaBondForce = new AmoebaBondForce();
// addBond: particle1, particle2, length, quadraticK
for( unsigned int jj = 0; jj < numberOfParticles; jj += 3 ){
amoebaHarmonicBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaHarmonicBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+1, 0.0000000e+00, 0.0000000e+00 );
amoebaBondForce->addBond( jj, jj+2, 0.0000000e+00, 0.0000000e+00 );
}
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondCubic( 0.0 );
amoebaHarmonicBondForce->setAmoebaGlobalHarmonicBondQuartic( 0.0 );
system.addForce(amoebaHarmonicBondForce);
amoebaBondForce->setAmoebaGlobalBondCubic( 0.0 );
amoebaBondForce->setAmoebaGlobalBondQuartic( 0.0 );
system.addForce(amoebaBondForce);
std::vector<Vec3> positions(numberOfParticles);
......
......@@ -52,9 +52,9 @@ extern "C" void initAmoebaReferenceKernels() {
AmoebaReferenceKernelFactory* factory = new AmoebaReferenceKernelFactory();
platform.registerKernelFactory(CalcAmoebaHarmonicBondForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaHarmonicAngleForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaHarmonicInPlaneAngleForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaBondForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaAngleForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaInPlaneAngleForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaPiTorsionForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaStretchBendForceKernel::Name(), factory);
platform.registerKernelFactory(CalcAmoebaOutOfPlaneBendForceKernel::Name(), factory);
......@@ -74,14 +74,14 @@ KernelImpl* AmoebaReferenceKernelFactory::createKernelImpl(std::string name, con
// create AmoebaReferenceData object if contextToAmoebaDataMap does not contain
// key equal to current context
if (name == CalcAmoebaHarmonicBondForceKernel::Name())
return new ReferenceCalcAmoebaHarmonicBondForceKernel(name, platform, context.getSystem());
if (name == CalcAmoebaBondForceKernel::Name())
return new ReferenceCalcAmoebaBondForceKernel(name, platform, context.getSystem());
if (name == CalcAmoebaHarmonicAngleForceKernel::Name())
return new ReferenceCalcAmoebaHarmonicAngleForceKernel(name, platform, context.getSystem());
if (name == CalcAmoebaAngleForceKernel::Name())
return new ReferenceCalcAmoebaAngleForceKernel(name, platform, context.getSystem());
if (name == CalcAmoebaHarmonicInPlaneAngleForceKernel::Name())
return new ReferenceCalcAmoebaHarmonicInPlaneAngleForceKernel(name, platform, context.getSystem());
if (name == CalcAmoebaInPlaneAngleForceKernel::Name())
return new ReferenceCalcAmoebaInPlaneAngleForceKernel(name, platform, context.getSystem());
if (name == CalcAmoebaPiTorsionForceKernel::Name())
return new ReferenceCalcAmoebaPiTorsionForceKernel(name, platform, context.getSystem());
......
......@@ -25,9 +25,9 @@
* -------------------------------------------------------------------------- */
#include "AmoebaReferenceKernels.h"
#include "AmoebaReferenceHarmonicBondForce.h"
#include "AmoebaReferenceHarmonicAngleForce.h"
#include "AmoebaReferenceHarmonicInPlaneAngleForce.h"
#include "AmoebaReferenceBondForce.h"
#include "AmoebaReferenceAngleForce.h"
#include "AmoebaReferenceInPlaneAngleForce.h"
#include "AmoebaReferencePiTorsionForce.h"
#include "AmoebaReferenceStretchBendForce.h"
#include "AmoebaReferenceOutOfPlaneBendForce.h"
......@@ -72,14 +72,14 @@ static RealVec& extractBoxSize(ContextImpl& context) {
// ***************************************************************************
ReferenceCalcAmoebaHarmonicBondForceKernel::ReferenceCalcAmoebaHarmonicBondForceKernel(std::string name, const Platform& platform, System& system) :
CalcAmoebaHarmonicBondForceKernel(name, platform), system(system) {
ReferenceCalcAmoebaBondForceKernel::ReferenceCalcAmoebaBondForceKernel(std::string name, const Platform& platform, System& system) :
CalcAmoebaBondForceKernel(name, platform), system(system) {
}
ReferenceCalcAmoebaHarmonicBondForceKernel::~ReferenceCalcAmoebaHarmonicBondForceKernel() {
ReferenceCalcAmoebaBondForceKernel::~ReferenceCalcAmoebaBondForceKernel() {
}
void ReferenceCalcAmoebaHarmonicBondForceKernel::initialize(const System& system, const AmoebaHarmonicBondForce& force) {
void ReferenceCalcAmoebaBondForceKernel::initialize(const System& system, const AmoebaBondForce& force) {
numBonds = force.getNumBonds();
for( int ii = 0; ii < numBonds; ii++) {
......@@ -93,30 +93,30 @@ void ReferenceCalcAmoebaHarmonicBondForceKernel::initialize(const System& system
length.push_back( static_cast<RealOpenMM>( lengthValue ) );
kQuadratic.push_back( static_cast<RealOpenMM>( kValue ) );
}
globalHarmonicBondCubic = static_cast<RealOpenMM>(force.getAmoebaGlobalHarmonicBondCubic());
globalHarmonicBondQuartic = static_cast<RealOpenMM>(force.getAmoebaGlobalHarmonicBondQuartic());
globalBondCubic = static_cast<RealOpenMM>(force.getAmoebaGlobalBondCubic());
globalBondQuartic = static_cast<RealOpenMM>(force.getAmoebaGlobalBondQuartic());
}
double ReferenceCalcAmoebaHarmonicBondForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
double ReferenceCalcAmoebaBondForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
vector<RealVec>& posData = extractPositions(context);
vector<RealVec>& forceData = extractForces(context);
AmoebaReferenceHarmonicBondForce amoebaReferenceHarmonicBondForce;
RealOpenMM energy = amoebaReferenceHarmonicBondForce.calculateForceAndEnergy( numBonds, posData, particle1, particle2, length, kQuadratic,
globalHarmonicBondCubic, globalHarmonicBondQuartic,
AmoebaReferenceBondForce amoebaReferenceBondForce;
RealOpenMM energy = amoebaReferenceBondForce.calculateForceAndEnergy( numBonds, posData, particle1, particle2, length, kQuadratic,
globalBondCubic, globalBondQuartic,
forceData );
return static_cast<double>(energy);
}
// ***************************************************************************
ReferenceCalcAmoebaHarmonicAngleForceKernel::ReferenceCalcAmoebaHarmonicAngleForceKernel(std::string name, const Platform& platform, System& system) :
CalcAmoebaHarmonicAngleForceKernel(name, platform), system(system) {
ReferenceCalcAmoebaAngleForceKernel::ReferenceCalcAmoebaAngleForceKernel(std::string name, const Platform& platform, System& system) :
CalcAmoebaAngleForceKernel(name, platform), system(system) {
}
ReferenceCalcAmoebaHarmonicAngleForceKernel::~ReferenceCalcAmoebaHarmonicAngleForceKernel() {
ReferenceCalcAmoebaAngleForceKernel::~ReferenceCalcAmoebaAngleForceKernel() {
}
void ReferenceCalcAmoebaHarmonicAngleForceKernel::initialize(const System& system, const AmoebaHarmonicAngleForce& force) {
void ReferenceCalcAmoebaAngleForceKernel::initialize(const System& system, const AmoebaAngleForce& force) {
numAngles = force.getNumAngles();
......@@ -130,29 +130,29 @@ void ReferenceCalcAmoebaHarmonicAngleForceKernel::initialize(const System& syste
angle.push_back( static_cast<RealOpenMM>( angleValue ) );
kQuadratic.push_back( static_cast<RealOpenMM>( k) );
}
globalHarmonicAngleCubic = static_cast<RealOpenMM>(force.getAmoebaGlobalHarmonicAngleCubic());
globalHarmonicAngleQuartic = static_cast<RealOpenMM>(force.getAmoebaGlobalHarmonicAngleQuartic());
globalHarmonicAnglePentic = static_cast<RealOpenMM>(force.getAmoebaGlobalHarmonicAnglePentic());
globalHarmonicAngleSextic = static_cast<RealOpenMM>(force.getAmoebaGlobalHarmonicAngleSextic());
globalAngleCubic = static_cast<RealOpenMM>(force.getAmoebaGlobalAngleCubic());
globalAngleQuartic = static_cast<RealOpenMM>(force.getAmoebaGlobalAngleQuartic());
globalAnglePentic = static_cast<RealOpenMM>(force.getAmoebaGlobalAnglePentic());
globalAngleSextic = static_cast<RealOpenMM>(force.getAmoebaGlobalAngleSextic());
}
double ReferenceCalcAmoebaHarmonicAngleForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
double ReferenceCalcAmoebaAngleForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
vector<RealVec>& posData = extractPositions(context);
vector<RealVec>& forceData = extractForces(context);
AmoebaReferenceHarmonicAngleForce amoebaReferenceHarmonicAngleForce;
RealOpenMM energy = amoebaReferenceHarmonicAngleForce.calculateForceAndEnergy( numAngles,
posData, particle1, particle2, particle3, angle, kQuadratic, globalHarmonicAngleCubic, globalHarmonicAngleQuartic, globalHarmonicAnglePentic, globalHarmonicAngleSextic, forceData );
AmoebaReferenceAngleForce amoebaReferenceAngleForce;
RealOpenMM energy = amoebaReferenceAngleForce.calculateForceAndEnergy( numAngles,
posData, particle1, particle2, particle3, angle, kQuadratic, globalAngleCubic, globalAngleQuartic, globalAnglePentic, globalAngleSextic, forceData );
return static_cast<double>(energy);
}
ReferenceCalcAmoebaHarmonicInPlaneAngleForceKernel::ReferenceCalcAmoebaHarmonicInPlaneAngleForceKernel(std::string name, const Platform& platform, System& system) :
CalcAmoebaHarmonicInPlaneAngleForceKernel(name, platform), system(system) {
ReferenceCalcAmoebaInPlaneAngleForceKernel::ReferenceCalcAmoebaInPlaneAngleForceKernel(std::string name, const Platform& platform, System& system) :
CalcAmoebaInPlaneAngleForceKernel(name, platform), system(system) {
}
ReferenceCalcAmoebaHarmonicInPlaneAngleForceKernel::~ReferenceCalcAmoebaHarmonicInPlaneAngleForceKernel() {
ReferenceCalcAmoebaInPlaneAngleForceKernel::~ReferenceCalcAmoebaInPlaneAngleForceKernel() {
}
void ReferenceCalcAmoebaHarmonicInPlaneAngleForceKernel::initialize(const System& system, const AmoebaHarmonicInPlaneAngleForce& force) {
void ReferenceCalcAmoebaInPlaneAngleForceKernel::initialize(const System& system, const AmoebaInPlaneAngleForce& force) {
numAngles = force.getNumAngles();
for (int ii = 0; ii < numAngles; ii++) {
......@@ -166,20 +166,20 @@ void ReferenceCalcAmoebaHarmonicInPlaneAngleForceKernel::initialize(const System
angle.push_back( static_cast<RealOpenMM>( angleValue ) );
kQuadratic.push_back( static_cast<RealOpenMM>( k ) );
}
globalHarmonicInPlaneAngleCubic = static_cast<RealOpenMM>(force.getAmoebaGlobalHarmonicInPlaneAngleCubic());
globalHarmonicInPlaneAngleQuartic = static_cast<RealOpenMM>(force.getAmoebaGlobalHarmonicInPlaneAngleQuartic());
globalHarmonicInPlaneAnglePentic = static_cast<RealOpenMM>(force.getAmoebaGlobalHarmonicInPlaneAnglePentic());
globalHarmonicInPlaneAngleSextic = static_cast<RealOpenMM>(force.getAmoebaGlobalHarmonicInPlaneAngleSextic());
globalInPlaneAngleCubic = static_cast<RealOpenMM>(force.getAmoebaGlobalInPlaneAngleCubic());
globalInPlaneAngleQuartic = static_cast<RealOpenMM>(force.getAmoebaGlobalInPlaneAngleQuartic());
globalInPlaneAnglePentic = static_cast<RealOpenMM>(force.getAmoebaGlobalInPlaneAnglePentic());
globalInPlaneAngleSextic = static_cast<RealOpenMM>(force.getAmoebaGlobalInPlaneAngleSextic());
}
double ReferenceCalcAmoebaHarmonicInPlaneAngleForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
double ReferenceCalcAmoebaInPlaneAngleForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
vector<RealVec>& posData = extractPositions(context);
vector<RealVec>& forceData = extractForces(context);
AmoebaReferenceHarmonicInPlaneAngleForce amoebaReferenceHarmonicInPlaneAngleForce;
RealOpenMM energy = amoebaReferenceHarmonicInPlaneAngleForce.calculateForceAndEnergy( numAngles, posData, particle1, particle2, particle3, particle4,
angle, kQuadratic, globalHarmonicInPlaneAngleCubic, globalHarmonicInPlaneAngleQuartic,
globalHarmonicInPlaneAnglePentic, globalHarmonicInPlaneAngleSextic, forceData );
AmoebaReferenceInPlaneAngleForce amoebaReferenceInPlaneAngleForce;
RealOpenMM energy = amoebaReferenceInPlaneAngleForce.calculateForceAndEnergy( numAngles, posData, particle1, particle2, particle3, particle4,
angle, kQuadratic, globalInPlaneAngleCubic, globalInPlaneAngleQuartic,
globalInPlaneAnglePentic, globalInPlaneAngleSextic, forceData );
return static_cast<double>(energy);
}
......
......@@ -35,21 +35,21 @@
namespace OpenMM {
/**
* This kernel is invoked by AmoebaHarmonicBondForce to calculate the forces acting on the system and the energy of the system.
* This kernel is invoked by AmoebaBondForce to calculate the forces acting on the system and the energy of the system.
*/
class ReferenceCalcAmoebaHarmonicBondForceKernel : public CalcAmoebaHarmonicBondForceKernel {
class ReferenceCalcAmoebaBondForceKernel : public CalcAmoebaBondForceKernel {
public:
ReferenceCalcAmoebaHarmonicBondForceKernel(std::string name,
ReferenceCalcAmoebaBondForceKernel(std::string name,
const Platform& platform,
System& system);
~ReferenceCalcAmoebaHarmonicBondForceKernel();
~ReferenceCalcAmoebaBondForceKernel();
/**
* Initialize the kernel.
*
* @param system the System this kernel will be applied to
* @param force the AmoebaHarmonicBondForce this kernel will be used for
* @param force the AmoebaBondForce this kernel will be used for
*/
void initialize(const System& system, const AmoebaHarmonicBondForce& force);
void initialize(const System& system, const AmoebaBondForce& force);
/**
* Execute the kernel to calculate the forces and/or energy.
*
......@@ -65,25 +65,25 @@ private:
std::vector<int> particle2;
std::vector<RealOpenMM> length;
std::vector<RealOpenMM> kQuadratic;
RealOpenMM globalHarmonicBondCubic;
RealOpenMM globalHarmonicBondQuartic;
RealOpenMM globalBondCubic;
RealOpenMM globalBondQuartic;
System& system;
};
/**
* This kernel is invoked by AmoebaHarmonicAngleForce to calculate the forces acting on the system and the energy of the system.
* This kernel is invoked by AmoebaAngleForce to calculate the forces acting on the system and the energy of the system.
*/
class ReferenceCalcAmoebaHarmonicAngleForceKernel : public CalcAmoebaHarmonicAngleForceKernel {
class ReferenceCalcAmoebaAngleForceKernel : public CalcAmoebaAngleForceKernel {
public:
ReferenceCalcAmoebaHarmonicAngleForceKernel(std::string name, const Platform& platform, System& system);
~ReferenceCalcAmoebaHarmonicAngleForceKernel();
ReferenceCalcAmoebaAngleForceKernel(std::string name, const Platform& platform, System& system);
~ReferenceCalcAmoebaAngleForceKernel();
/**
* Initialize the kernel.
*
* @param system the System this kernel will be applied to
* @param force the AmoebaHarmonicAngleForce this kernel will be used for
* @param force the AmoebaAngleForce this kernel will be used for
*/
void initialize(const System& system, const AmoebaHarmonicAngleForce& force);
void initialize(const System& system, const AmoebaAngleForce& force);
/**
* Execute the kernel to calculate the forces and/or energy.
*
......@@ -100,27 +100,27 @@ private:
std::vector<int> particle3;
std::vector<RealOpenMM> angle;
std::vector<RealOpenMM> kQuadratic;
RealOpenMM globalHarmonicAngleCubic;
RealOpenMM globalHarmonicAngleQuartic;
RealOpenMM globalHarmonicAnglePentic;
RealOpenMM globalHarmonicAngleSextic;
RealOpenMM globalAngleCubic;
RealOpenMM globalAngleQuartic;
RealOpenMM globalAnglePentic;
RealOpenMM globalAngleSextic;
System& system;
};
/**
* This kernel is invoked by AmoebaHarmonicInPlaneAngleForce to calculate the forces acting on the system and the energy of the system.
* This kernel is invoked by AmoebaInPlaneAngleForce to calculate the forces acting on the system and the energy of the system.
*/
class ReferenceCalcAmoebaHarmonicInPlaneAngleForceKernel : public CalcAmoebaHarmonicInPlaneAngleForceKernel {
class ReferenceCalcAmoebaInPlaneAngleForceKernel : public CalcAmoebaInPlaneAngleForceKernel {
public:
ReferenceCalcAmoebaHarmonicInPlaneAngleForceKernel(std::string name, const Platform& platform, System& system);
~ReferenceCalcAmoebaHarmonicInPlaneAngleForceKernel();
ReferenceCalcAmoebaInPlaneAngleForceKernel(std::string name, const Platform& platform, System& system);
~ReferenceCalcAmoebaInPlaneAngleForceKernel();
/**
* Initialize the kernel.
*
* @param system the System this kernel will be applied to
* @param force the AmoebaHarmonicInPlaneAngleForce this kernel will be used for
* @param force the AmoebaInPlaneAngleForce this kernel will be used for
*/
void initialize(const System& system, const AmoebaHarmonicInPlaneAngleForce& force);
void initialize(const System& system, const AmoebaInPlaneAngleForce& force);
/**
* Execute the kernel to calculate the forces and/or energy.
*
......@@ -138,10 +138,10 @@ private:
std::vector<int> particle4;
std::vector<RealOpenMM> angle;
std::vector<RealOpenMM> kQuadratic;
RealOpenMM globalHarmonicInPlaneAngleCubic;
RealOpenMM globalHarmonicInPlaneAngleQuartic;
RealOpenMM globalHarmonicInPlaneAnglePentic;
RealOpenMM globalHarmonicInPlaneAngleSextic;
RealOpenMM globalInPlaneAngleCubic;
RealOpenMM globalInPlaneAngleQuartic;
RealOpenMM globalInPlaneAnglePentic;
RealOpenMM globalInPlaneAngleSextic;
System& system;
};
......
......@@ -23,7 +23,7 @@
*/
#include "AmoebaReferenceForce.h"
#include "AmoebaReferenceHarmonicAngleForce.h"
#include "AmoebaReferenceAngleForce.h"
using std::vector;
using OpenMM::RealVec;
......@@ -46,7 +46,7 @@ using OpenMM::RealVec;
--------------------------------------------------------------------------------------- */
RealOpenMM AmoebaReferenceHarmonicAngleForce::getPrefactorsGivenAngleCosine( RealOpenMM cosine,
RealOpenMM AmoebaReferenceAngleForce::getPrefactorsGivenAngleCosine( RealOpenMM cosine,
RealOpenMM idealAngle, RealOpenMM angleK,
RealOpenMM angleCubic, RealOpenMM angleQuartic,
RealOpenMM anglePentic, RealOpenMM angleSextic,
......@@ -62,7 +62,7 @@ RealOpenMM AmoebaReferenceHarmonicAngleForce::getPrefactorsGivenAngleCosine( Rea
static const RealOpenMM five = 5.0;
static const RealOpenMM six = 6.0;
// static const std::string methodName = "AmoebaReferenceHarmonicAngleForce::getPrefactorsGivenAngleCosine";
// static const std::string methodName = "AmoebaReferenceAngleForce::getPrefactorsGivenAngleCosine";
// ---------------------------------------------------------------------------------------
......@@ -96,7 +96,7 @@ RealOpenMM AmoebaReferenceHarmonicAngleForce::getPrefactorsGivenAngleCosine( Rea
/**---------------------------------------------------------------------------------------
Calculate Amoeba harmonic angle ixn (force and energy)
Calculate Amoeba angle ixn (force and energy)
@param positionAtomA Cartesian coordinates of atom A
@param positionAtomB Cartesian coordinates of atom B
......@@ -113,7 +113,7 @@ RealOpenMM AmoebaReferenceHarmonicAngleForce::getPrefactorsGivenAngleCosine( Rea
--------------------------------------------------------------------------------------- */
RealOpenMM AmoebaReferenceHarmonicAngleForce::calculateAngleIxn( const RealVec& positionAtomA, const RealVec& positionAtomB,
RealOpenMM AmoebaReferenceAngleForce::calculateAngleIxn( const RealVec& positionAtomA, const RealVec& positionAtomB,
const RealVec& positionAtomC,
RealOpenMM angle, RealOpenMM angleK,
RealOpenMM angleCubic, RealOpenMM angleQuartic,
......@@ -122,7 +122,7 @@ RealOpenMM AmoebaReferenceHarmonicAngleForce::calculateAngleIxn( const RealVec&
// ---------------------------------------------------------------------------------------
//static const std::string methodName = "AmoebaReferenceHarmonicAngleForce::calculateHarmonicForce";
//static const std::string methodName = "AmoebaReferenceAngleForce::calculateAngleIxn";
static const RealOpenMM zero = 0.0;
static const RealOpenMM one = 1.0;
......@@ -183,7 +183,7 @@ RealOpenMM AmoebaReferenceHarmonicAngleForce::calculateAngleIxn( const RealVec&
return energy;
}
RealOpenMM AmoebaReferenceHarmonicAngleForce::calculateForceAndEnergy( int numAngles, vector<RealVec>& posData,
RealOpenMM AmoebaReferenceAngleForce::calculateForceAndEnergy( int numAngles, vector<RealVec>& posData,
const std::vector<int>& particle1,
const std::vector<int>& particle2,
const std::vector<int>& particle3,
......
......@@ -22,15 +22,15 @@
* WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#ifndef __AmoebaReferenceHarmonicAngleForce_H__
#define __AmoebaReferenceHarmonicAngleForce_H__
#ifndef __AmoebaReferenceAngleForce_H__
#define __AmoebaReferenceAngleForce_H__
#include "SimTKUtilities/RealVec.h"
#include <vector>
// ---------------------------------------------------------------------------------------
class AmoebaReferenceHarmonicAngleForce {
class AmoebaReferenceAngleForce {
public:
......@@ -40,7 +40,7 @@ public:
--------------------------------------------------------------------------------------- */
AmoebaReferenceHarmonicAngleForce( ){};
AmoebaReferenceAngleForce( ){};
/**---------------------------------------------------------------------------------------
......@@ -48,11 +48,11 @@ public:
--------------------------------------------------------------------------------------- */
~AmoebaReferenceHarmonicAngleForce( ){};
~AmoebaReferenceAngleForce( ){};
/**---------------------------------------------------------------------------------------
Calculate Amoeba harmonic angle ixns (force and energy)
Calculate Amoeba angle ixns (force and energy)
@param numBonds number of angles
@param posData particle positions
......@@ -78,10 +78,10 @@ public:
const std::vector<int>& particle3,
const std::vector<RealOpenMM>& angle,
const std::vector<RealOpenMM>& kQuadratic,
RealOpenMM globalHarmonicAngleCubic,
RealOpenMM globalHarmonicAngleQuartic,
RealOpenMM globalHarmonicAnglePentic,
RealOpenMM globalHarmonicAngleSextic,
RealOpenMM globalAngleCubic,
RealOpenMM globalAngleQuartic,
RealOpenMM globalAnglePentic,
RealOpenMM globalAngleSextic,
std::vector<OpenMM::RealVec>& forceData ) const;
private:
......@@ -111,7 +111,7 @@ private:
/**---------------------------------------------------------------------------------------
Calculate Amoeba harmonic angle ixn (force and energy)
Calculate Amoeba angle ixn (force and energy)
@param positionAtomA Cartesian coordinates of atom A
@param positionAtomB Cartesian coordinates of atom B
......@@ -139,4 +139,4 @@ private:
// ---------------------------------------------------------------------------------------
#endif // _AmoebaReferenceHarmonicAngleForce___
#endif // _AmoebaReferenceAngleForce___
......@@ -22,7 +22,7 @@
* WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
*/
#include "AmoebaReferenceHarmonicBondForce.h"
#include "AmoebaReferenceBondForce.h"
#include "AmoebaReferenceForce.h"
using std::vector;
......@@ -30,7 +30,7 @@ using OpenMM::RealVec;
/**---------------------------------------------------------------------------------------
Calculate Amoeba harmonic bond ixn (force and energy)
Calculate Amoeba bond ixn (force and energy)
@param positionAtomA Cartesian coordinates of atom A
@param positionAtomB Cartesian coordinates of atom B
......@@ -44,14 +44,14 @@ using OpenMM::RealVec;
--------------------------------------------------------------------------------------- */
RealOpenMM AmoebaReferenceHarmonicBondForce::calculateBondIxn( const RealVec& positionAtomA, const RealVec& positionAtomB,
RealOpenMM AmoebaReferenceBondForce::calculateBondIxn( const RealVec& positionAtomA, const RealVec& positionAtomB,
RealOpenMM bondLength, RealOpenMM bondK,
RealOpenMM bondCubic, RealOpenMM bondQuartic,
RealVec* forces ) const {
// ---------------------------------------------------------------------------------------
//static const std::string methodName = "AmoebaReferenceHarmonicBondForce::calculateHarmonicForce";
//static const std::string methodName = "AmoebaReferenceBondForce::calculateBondIxn";
static const RealOpenMM zero = 0.0;
static const RealOpenMM one = 1.0;
......@@ -88,14 +88,14 @@ RealOpenMM AmoebaReferenceHarmonicBondForce::calculateBondIxn( const RealVec& po
return energy;
}
RealOpenMM AmoebaReferenceHarmonicBondForce::calculateForceAndEnergy( int numBonds,
RealOpenMM AmoebaReferenceBondForce::calculateForceAndEnergy( int numBonds,
vector<RealVec>& particlePositions,
const std::vector<int>& particle1,
const std::vector<int>& particle2,
const std::vector<RealOpenMM>& length,
const std::vector<RealOpenMM>& kQuadratic,
RealOpenMM globalHarmonicBondCubic,
RealOpenMM globalHarmonicBondQuartic,
RealOpenMM globalBondCubic,
RealOpenMM globalBondQuartic,
vector<RealVec>& forceData ) const {
RealOpenMM energy = 0.0;
for( int ii = 0; ii < numBonds; ii++ ){
......@@ -106,7 +106,7 @@ RealOpenMM AmoebaReferenceHarmonicBondForce::calculateForceAndEnergy( int numBon
RealVec forces[2];
energy += calculateBondIxn( particlePositions[particle1Index], particlePositions[particle2Index],
bondLength, bondK, globalHarmonicBondCubic, globalHarmonicBondQuartic,
bondLength, bondK, globalBondCubic, globalBondQuartic,
forces );
for( int jj = 0; jj < 3; jj++ ){
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment