Skip to content

GitLab

"platforms/vscode:/vscode.git/clone" did not exist on "3d69185ab6ae9fd4c35482b1787799869ff4703a"
Unverified Commit 462a9be3 authored by Peter Eastman's avatar Peter Eastman Committed by GitHub
Browse files

Merged tests for different platforms (#4652)

parent e2f659ce
Show whitespace changes
Inline Side-by-side
Showing with 401 additions and 1068 deletions
platforms/cuda/tests/TestCudaCheckpoints.cpp
View file @ 462a9be3
......@@ -32,105 +32,5 @@
#include "CudaTests.h"
#include "TestCheckpoints.h"
void testCheckpoint() {
const int numParticles = 100;
const double boxSize = 5.0;
const double temperature = 200.0;
System system;
system.addForce(new AndersenThermostat(0.0, 100.0));
NonbondedForce* nonbonded = new NonbondedForce();
system.addForce(nonbonded);
nonbonded->setNonbondedMethod(NonbondedForce::CutoffPeriodic);
vector<Vec3> positions(numParticles);
OpenMM_SFMT::SFMT sfmt;
init_gen_rand(0, sfmt);
for (int i = 0; i < numParticles; i++) {
system.addParticle(1.0);
nonbonded->addParticle(i%2 == 0 ? 0.1 : -0.1, 0.2, 0.1);
bool clash;
do {
clash = false;
positions[i] = Vec3(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
for (int j = 0; j < i; j++) {
Vec3 delta = positions[i]-positions[j];
if (sqrt(delta.dot(delta)) < 0.1)
clash = true;
}
} while (clash);
}
VerletIntegrator integrator(0.001);
Context context(system, integrator, platform);
context.setPositions(positions);
context.setPeriodicBoxVectors(Vec3(boxSize, 0, 0), Vec3(0, boxSize, 0), Vec3(0, 0, boxSize));
context.setParameter(AndersenThermostat::Temperature(), temperature);
// Run for a little while.
integrator.step(100);
// Record the current state and make a checkpoint.
State s1 = context.getState(State::Positions | State::Velocities | State::Parameters);
stringstream stream1(ios_base::out | ios_base::in | ios_base::binary);
context.createCheckpoint(stream1);
// Continue the simulation for a few more steps and record the state again.
integrator.step(10);
State s2 = context.getState(State::Positions | State::Velocities | State::Parameters);
// Restore from the checkpoint and see if everything gets restored correctly.
context.setPeriodicBoxVectors(Vec3(2*boxSize, 0, 0), Vec3(0, 2*boxSize, 0), Vec3(0, 0, 2*boxSize));
context.setParameter(AndersenThermostat::Temperature(), temperature+10);
context.loadCheckpoint(stream1);
State s3 = context.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s1, s3);
// Now simulate from there and see if the trajectory is identical.
integrator.step(10);
State s4 = context.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s2, s4);
// Create a new Context that uses multiple devices.
string deviceIndex = platform.getPropertyValue(context, CudaPlatform::CudaDeviceIndex());
map<string, string> props;
props[CudaPlatform::CudaDeviceIndex()] = deviceIndex+","+deviceIndex;
VerletIntegrator integrator2(0.001);
Context context2(system, integrator2, platform, props);
context2.setPositions(positions);
context2.setPeriodicBoxVectors(Vec3(boxSize, 0, 0), Vec3(0, boxSize, 0), Vec3(0, 0, boxSize));
context2.setParameter(AndersenThermostat::Temperature(), temperature);
// Now repeat all of the above tests with it.
integrator2.step(100);
State s5 = context2.getState(State::Positions | State::Velocities | State::Parameters);
stringstream stream2(ios_base::out | ios_base::in | ios_base::binary);
context2.createCheckpoint(stream2);
integrator2.step(10);
State s6 = context2.getState(State::Positions | State::Velocities | State::Parameters);
context2.setPeriodicBoxVectors(Vec3(2*boxSize, 0, 0), Vec3(0, 2*boxSize, 0), Vec3(0, 0, 2*boxSize));
context2.setParameter(AndersenThermostat::Temperature(), temperature+10);
context2.loadCheckpoint(stream2);
State s7 = context2.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s5, s7);
integrator2.step(10);
State s8 = context2.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s6, s8);
// See if a checkpoint created from one Context can be loaded into a different one.
VerletIntegrator integrator3(0.001);
Context context3(system, integrator3, platform);
stream1.seekg(0, stream1.beg);
context3.loadCheckpoint(stream1);
State s9 = context3.getState(State::Positions | State::Velocities | State::Parameters | State::Energy);
compareStates(s1, s9);
}
void runPlatformTests() {
testCheckpoint();
}
platforms/cuda/tests/TestCudaCustomIntegrator.cpp
View file @ 462a9be3
......@@ -32,72 +32,5 @@
#include "CudaTests.h"
#include "TestCustomIntegrator.h"
/**
* Make sure random numbers are computed correctly when steps get merged.
*/
void testMergedRandoms() {
const int numParticles = 10;
const int numSteps = 10;
System system;
for (int i = 0; i < numParticles; i++)
system.addParticle(1.0);
CustomIntegrator integrator(0.1);
integrator.addPerDofVariable("dofUniform1", 0);
integrator.addPerDofVariable("dofUniform2", 0);
integrator.addPerDofVariable("dofGaussian1", 0);
integrator.addPerDofVariable("dofGaussian2", 0);
integrator.addGlobalVariable("globalUniform1", 0);
integrator.addGlobalVariable("globalUniform2", 0);
integrator.addGlobalVariable("globalGaussian1", 0);
integrator.addGlobalVariable("globalGaussian2", 0);
integrator.addComputePerDof("dofUniform1", "uniform");
integrator.addComputePerDof("dofUniform2", "uniform");
integrator.addComputePerDof("dofGaussian1", "gaussian");
integrator.addComputePerDof("dofGaussian2", "gaussian");
integrator.addComputeGlobal("globalUniform1", "uniform");
integrator.addComputeGlobal("globalUniform2", "uniform");
integrator.addComputeGlobal("globalGaussian1", "gaussian");
integrator.addComputeGlobal("globalGaussian2", "gaussian");
Context context(system, integrator, platform);
// See if the random numbers are computed correctly.
vector<Vec3> values1, values2;
for (int i = 0; i < numSteps; i++) {
integrator.step(1);
integrator.getPerDofVariable(0, values1);
integrator.getPerDofVariable(1, values2);
for (int i = 0; i < numParticles; i++)
for (int j = 0; j < 3; j++) {
double v1 = values1[i][j];
double v2 = values2[i][j];
ASSERT(v1 >= 0 && v1 < 1);
ASSERT(v2 >= 0 && v2 < 1);
ASSERT(v1 != v2);
}
integrator.getPerDofVariable(2, values1);
integrator.getPerDofVariable(3, values2);
for (int i = 0; i < numParticles; i++)
for (int j = 0; j < 3; j++) {
double v1 = values1[i][j];
double v2 = values2[i][j];
ASSERT(v1 >= -10 && v1 < 10);
ASSERT(v2 >= -10 && v2 < 10);
ASSERT(v1 != v2);
}
double v1 = integrator.getGlobalVariable(0);
double v2 = integrator.getGlobalVariable(1);
ASSERT(v1 >= 0 && v1 < 1);
ASSERT(v2 >= 0 && v2 < 1);
ASSERT(v1 != v2);
v1 = integrator.getGlobalVariable(2);
v2 = integrator.getGlobalVariable(3);
ASSERT(v1 >= -10 && v1 < 10);
ASSERT(v2 >= -10 && v2 < 10);
ASSERT(v1 != v2);
}
}
void runPlatformTests() {
testMergedRandoms();
}
platforms/cuda/tests/TestCudaMultipleForces.cpp
View file @ 462a9be3
......@@ -29,85 +29,8 @@
* USE OR OTHER DEALINGS IN THE SOFTWARE. *
* -------------------------------------------------------------------------- */
/**
* This tests a system with multiple forces, to make sure CudaBondedUtilities is
* processing them correctly.
*/
#include "CudaTests.h"
#include "TestMultipleForces.h"
#include "openmm/internal/AssertionUtilities.h"
#include "openmm/Context.h"
#include "CudaPlatform.h"
#include "ReferencePlatform.h"
#include "openmm/HarmonicAngleForce.h"
#include "openmm/HarmonicBondForce.h"
#include "openmm/PeriodicTorsionForce.h"
#include "openmm/RBTorsionForce.h"
#include "openmm/System.h"
#include "openmm/VerletIntegrator.h"
#include "SimTKOpenMMRealType.h"
#include "sfmt/SFMT.h"
#include <iostream>
#include <vector>
using namespace OpenMM;
using namespace std;
CudaPlatform platform;
const double TOL = 1e-4;
void testForces() {
const int numParticles = 100;
System system;
for (int i = 0; i < numParticles; i++)
system.addParticle(1.0);
HarmonicBondForce* bonds = new HarmonicBondForce();
for (int i = 0; i < numParticles-1; i++)
bonds->addBond(i, i+1, 1.0, 1.5);
system.addForce(bonds);
HarmonicAngleForce* angles = new HarmonicAngleForce();
for (int i = 0; i < numParticles-2; i++)
angles->addAngle(i, i+1, i+2, 2.0, 1.5);
system.addForce(angles);
PeriodicTorsionForce* periodic = new PeriodicTorsionForce();
for (int i = 0; i < numParticles-3; i++)
periodic->addTorsion(i, i+1, i+2, i+3, 2, PI_M/3, 1.1);
system.addForce(periodic);
RBTorsionForce* rb = new RBTorsionForce();
for (int i = 0; i < numParticles-3; i += 3)
rb->addTorsion(i, i+1, i+2, i+3, 1.0, 1.1, 1.2, 0.3, 0.4, 0.5);
system.addForce(rb);
ReferencePlatform ref;
VerletIntegrator integrator1(0.01);
Context context1(system, integrator1, ref);
VerletIntegrator integrator2(0.01);
Context context2(system, integrator2, platform);
OpenMM_SFMT::SFMT sfmt;
init_gen_rand(0, sfmt);
vector<Vec3> positions(numParticles);
for (int i = 0; i < numParticles; i++)
positions[i] = Vec3(i, genrand_real2(sfmt), genrand_real2(sfmt));
context1.setPositions(positions);
context2.setPositions(positions);
State state1 = context1.getState(State::Forces | State::Energy);
State state2 = context2.getState(State::Forces | State::Energy);
const vector<Vec3>& forces1 = state1.getForces();
const vector<Vec3>& forces2 = state2.getForces();
for (int i = 0; i < numParticles; i++)
ASSERT_EQUAL_VEC(forces1[i], forces2[i], TOL);
ASSERT_EQUAL_TOL(state1.getPotentialEnergy(), state2.getPotentialEnergy(), TOL);
}
int main(int argc, char* argv[]) {
try {
if (argc > 1)
platform.setPropertyDefaultValue("CudaPrecision", string(argv[1]));
testForces();
}
catch(const exception& e) {
cout << "exception: " << e.what() << endl;
return 1;
}
cout << "Done" << endl;
return 0;
void runPlatformTests() {
}
platforms/cuda/tests/TestCudaVirtualSites.cpp
View file @ 462a9be3
......@@ -32,112 +32,5 @@
#include "CudaTests.h"
#include "TestVirtualSites.h"
/**
* Make sure that atom reordering respects virtual sites.
*/
void testReordering() {
const double cutoff = 2.0;
const double boxSize = 20.0;
System system;
NonbondedForce* nonbonded = new NonbondedForce();
system.addForce(nonbonded);
nonbonded->setNonbondedMethod(NonbondedForce::CutoffNonPeriodic);
nonbonded->setCutoffDistance(cutoff);
vector<Vec3> positions;
OpenMM_SFMT::SFMT sfmt;
init_gen_rand(0, sfmt);
// Create linear molecules with TwoParticleAverage virtual sites.
for (int i = 0; i < 50; i++) {
int start = system.getNumParticles();
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(0.0);
system.setVirtualSite(start+2, new TwoParticleAverageSite(start, start+1, 0.4, 0.6));
system.addConstraint(start, start+1, 2.0);
for (int i = 0; i < 3; i++) {
nonbonded->addParticle(0, 0.2, 1);
for (int j = 0; j < i; j++)
nonbonded->addException(start+i, start+j, 0, 1, 0);
}
Vec3 pos(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
positions.push_back(pos);
positions.push_back(pos+Vec3(2, 0, 0));
positions.push_back(Vec3());
}
// Create planar molecules with ThreeParticleAverage virtual sites.
for (int i = 0; i < 50; i++) {
int start = system.getNumParticles();
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(0.0);
system.setVirtualSite(start+3, new ThreeParticleAverageSite(start, start+1, start+2, 0.3, 0.5, 0.2));
system.addConstraint(start, start+1, 1.0);
system.addConstraint(start, start+2, 1.0);
system.addConstraint(start+1, start+2, sqrt(2.0));
for (int i = 0; i < 4; i++) {
nonbonded->addParticle(0, 0.2, 1);
for (int j = 0; j < i; j++)
nonbonded->addException(start+i, start+j, 0, 1, 0);
}
Vec3 pos(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
positions.push_back(pos);
positions.push_back(pos+Vec3(1, 0, 0));
positions.push_back(pos+Vec3(0, 1, 0));
positions.push_back(Vec3());
}
// Create tetrahedral molecules with OutOfPlane virtual sites.
for (int i = 0; i < 50; i++) {
int start = system.getNumParticles();
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(0.0);
system.setVirtualSite(start+3, new OutOfPlaneSite(start, start+1, start+2, 0.3, 0.5, 0.2));
system.addConstraint(start, start+1, 1.0);
system.addConstraint(start, start+2, 1.0);
system.addConstraint(start+1, start+2, sqrt(2.0));
for (int i = 0; i < 4; i++) {
nonbonded->addParticle(0, 0.2, 1);
for (int j = 0; j < i; j++)
nonbonded->addException(start+i, start+j, 0, 1, 0);
}
Vec3 pos(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
positions.push_back(pos);
positions.push_back(pos+Vec3(1, 0, 0));
positions.push_back(pos+Vec3(0, 1, 0));
positions.push_back(Vec3());
}
// Simulate it and check conservation laws.
LangevinIntegrator integrator(300.0, 0.1, 0.002);
Context context(system, integrator, platform);
context.setPositions(positions);
context.applyConstraints(0.0001);
for (int i = 0; i < 1000; i++) {
State state = context.getState(State::Positions);
const vector<Vec3>& pos = state.getPositions();
for (int j = 0; j < 150; j += 3)
ASSERT_EQUAL_VEC(pos[j]*0.4+pos[j+1]*0.6, pos[j+2], 1e-5);
for (int j = 150; j < 350; j += 4)
ASSERT_EQUAL_VEC(pos[j]*0.3+pos[j+1]*0.5+pos[j+2]*0.2, pos[j+3], 1e-5);
for (int j = 350; j < 550; j += 4) {
Vec3 v12 = pos[j+1]-pos[j];
Vec3 v13 = pos[j+2]-pos[j];
Vec3 cross = v12.cross(v13);
ASSERT_EQUAL_VEC(pos[j]+v12*0.3+v13*0.5+cross*0.2, pos[j+3], 1e-5);
}
integrator.step(1);
}
}
void runPlatformTests() {
testReordering();
}
platforms/hip/tests/TestHipCheckpoints.cpp
View file @ 462a9be3
......@@ -33,105 +33,5 @@
#include "HipTests.h"
#include "TestCheckpoints.h"
void testCheckpoint() {
const int numParticles = 100;
const double boxSize = 5.0;
const double temperature = 200.0;
System system;
system.addForce(new AndersenThermostat(0.0, 100.0));
NonbondedForce* nonbonded = new NonbondedForce();
system.addForce(nonbonded);
nonbonded->setNonbondedMethod(NonbondedForce::CutoffPeriodic);
vector<Vec3> positions(numParticles);
OpenMM_SFMT::SFMT sfmt;
init_gen_rand(0, sfmt);
for (int i = 0; i < numParticles; i++) {
system.addParticle(1.0);
nonbonded->addParticle(i%2 == 0 ? 0.1 : -0.1, 0.2, 0.1);
bool clash;
do {
clash = false;
positions[i] = Vec3(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
for (int j = 0; j < i; j++) {
Vec3 delta = positions[i]-positions[j];
if (sqrt(delta.dot(delta)) < 0.1)
clash = true;
}
} while (clash);
}
VerletIntegrator integrator(0.001);
Context context(system, integrator, platform);
context.setPositions(positions);
context.setPeriodicBoxVectors(Vec3(boxSize, 0, 0), Vec3(0, boxSize, 0), Vec3(0, 0, boxSize));
context.setParameter(AndersenThermostat::Temperature(), temperature);
// Run for a little while.
integrator.step(100);
// Record the current state and make a checkpoint.
State s1 = context.getState(State::Positions | State::Velocities | State::Parameters);
stringstream stream1(ios_base::out | ios_base::in | ios_base::binary);
context.createCheckpoint(stream1);
// Continue the simulation for a few more steps and record the state again.
integrator.step(10);
State s2 = context.getState(State::Positions | State::Velocities | State::Parameters);
// Restore from the checkpoint and see if everything gets restored correctly.
context.setPeriodicBoxVectors(Vec3(2*boxSize, 0, 0), Vec3(0, 2*boxSize, 0), Vec3(0, 0, 2*boxSize));
context.setParameter(AndersenThermostat::Temperature(), temperature+10);
context.loadCheckpoint(stream1);
State s3 = context.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s1, s3);
// Now simulate from there and see if the trajectory is identical.
integrator.step(10);
State s4 = context.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s2, s4);
// Create a new Context that uses multiple devices.
string deviceIndex = platform.getPropertyValue(context, HipPlatform::HipDeviceIndex());
map<string, string> props;
props[HipPlatform::HipDeviceIndex()] = deviceIndex+","+deviceIndex;
VerletIntegrator integrator2(0.001);
Context context2(system, integrator2, platform, props);
context2.setPositions(positions);
context2.setPeriodicBoxVectors(Vec3(boxSize, 0, 0), Vec3(0, boxSize, 0), Vec3(0, 0, boxSize));
context2.setParameter(AndersenThermostat::Temperature(), temperature);
// Now repeat all of the above tests with it.
integrator2.step(100);
State s5 = context2.getState(State::Positions | State::Velocities | State::Parameters);
stringstream stream2(ios_base::out | ios_base::in | ios_base::binary);
context2.createCheckpoint(stream2);
integrator2.step(10);
State s6 = context2.getState(State::Positions | State::Velocities | State::Parameters);
context2.setPeriodicBoxVectors(Vec3(2*boxSize, 0, 0), Vec3(0, 2*boxSize, 0), Vec3(0, 0, 2*boxSize));
context2.setParameter(AndersenThermostat::Temperature(), temperature+10);
context2.loadCheckpoint(stream2);
State s7 = context2.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s5, s7);
integrator2.step(10);
State s8 = context2.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s6, s8);
// See if a checkpoint created from one Context can be loaded into a different one.
VerletIntegrator integrator3(0.001);
Context context3(system, integrator3, platform);
stream1.seekg(0, stream1.beg);
context3.loadCheckpoint(stream1);
State s9 = context3.getState(State::Positions | State::Velocities | State::Parameters | State::Energy);
compareStates(s1, s9);
}
void runPlatformTests() {
testCheckpoint();
}
platforms/hip/tests/TestHipCustomIntegrator.cpp
View file @ 462a9be3
......@@ -33,72 +33,5 @@
#include "HipTests.h"
#include "TestCustomIntegrator.h"
/**
* Make sure random numbers are computed correctly when steps get merged.
*/
void testMergedRandoms() {
const int numParticles = 10;
const int numSteps = 10;
System system;
for (int i = 0; i < numParticles; i++)
system.addParticle(1.0);
CustomIntegrator integrator(0.1);
integrator.addPerDofVariable("dofUniform1", 0);
integrator.addPerDofVariable("dofUniform2", 0);
integrator.addPerDofVariable("dofGaussian1", 0);
integrator.addPerDofVariable("dofGaussian2", 0);
integrator.addGlobalVariable("globalUniform1", 0);
integrator.addGlobalVariable("globalUniform2", 0);
integrator.addGlobalVariable("globalGaussian1", 0);
integrator.addGlobalVariable("globalGaussian2", 0);
integrator.addComputePerDof("dofUniform1", "uniform");
integrator.addComputePerDof("dofUniform2", "uniform");
integrator.addComputePerDof("dofGaussian1", "gaussian");
integrator.addComputePerDof("dofGaussian2", "gaussian");
integrator.addComputeGlobal("globalUniform1", "uniform");
integrator.addComputeGlobal("globalUniform2", "uniform");
integrator.addComputeGlobal("globalGaussian1", "gaussian");
integrator.addComputeGlobal("globalGaussian2", "gaussian");
Context context(system, integrator, platform);
// See if the random numbers are computed correctly.
vector<Vec3> values1, values2;
for (int i = 0; i < numSteps; i++) {
integrator.step(1);
integrator.getPerDofVariable(0, values1);
integrator.getPerDofVariable(1, values2);
for (int i = 0; i < numParticles; i++)
for (int j = 0; j < 3; j++) {
double v1 = values1[i][j];
double v2 = values2[i][j];
ASSERT(v1 >= 0 && v1 < 1);
ASSERT(v2 >= 0 && v2 < 1);
ASSERT(v1 != v2);
}
integrator.getPerDofVariable(2, values1);
integrator.getPerDofVariable(3, values2);
for (int i = 0; i < numParticles; i++)
for (int j = 0; j < 3; j++) {
double v1 = values1[i][j];
double v2 = values2[i][j];
ASSERT(v1 >= -10 && v1 < 10);
ASSERT(v2 >= -10 && v2 < 10);
ASSERT(v1 != v2);
}
double v1 = integrator.getGlobalVariable(0);
double v2 = integrator.getGlobalVariable(1);
ASSERT(v1 >= 0 && v1 < 1);
ASSERT(v2 >= 0 && v2 < 1);
ASSERT(v1 != v2);
v1 = integrator.getGlobalVariable(2);
v2 = integrator.getGlobalVariable(3);
ASSERT(v1 >= -10 && v1 < 10);
ASSERT(v2 >= -10 && v2 < 10);
ASSERT(v1 != v2);
}
}
void runPlatformTests() {
testMergedRandoms();
}
platforms/hip/tests/TestHipMultipleForces.cpp
View file @ 462a9be3
......@@ -30,85 +30,8 @@
* USE OR OTHER DEALINGS IN THE SOFTWARE. *
* -------------------------------------------------------------------------- */
/**
* This tests a system with multiple forces, to make sure HipBondedUtilities is
* processing them correctly.
*/
#include "HipTests.h"
#include "TestMultipleForces.h"
#include "openmm/internal/AssertionUtilities.h"
#include "openmm/Context.h"
#include "HipPlatform.h"
#include "ReferencePlatform.h"
#include "openmm/HarmonicAngleForce.h"
#include "openmm/HarmonicBondForce.h"
#include "openmm/PeriodicTorsionForce.h"
#include "openmm/RBTorsionForce.h"
#include "openmm/System.h"
#include "openmm/VerletIntegrator.h"
#include "SimTKOpenMMRealType.h"
#include "sfmt/SFMT.h"
#include <iostream>
#include <vector>
using namespace OpenMM;
using namespace std;
HipPlatform platform;
const double TOL = 1e-4;
void testForces() {
const int numParticles = 100;
System system;
for (int i = 0; i < numParticles; i++)
system.addParticle(1.0);
HarmonicBondForce* bonds = new HarmonicBondForce();
for (int i = 0; i < numParticles-1; i++)
bonds->addBond(i, i+1, 1.0, 1.5);
system.addForce(bonds);
HarmonicAngleForce* angles = new HarmonicAngleForce();
for (int i = 0; i < numParticles-2; i++)
angles->addAngle(i, i+1, i+2, 2.0, 1.5);
system.addForce(angles);
PeriodicTorsionForce* periodic = new PeriodicTorsionForce();
for (int i = 0; i < numParticles-3; i++)
periodic->addTorsion(i, i+1, i+2, i+3, 2, PI_M/3, 1.1);
system.addForce(periodic);
RBTorsionForce* rb = new RBTorsionForce();
for (int i = 0; i < numParticles-3; i += 3)
rb->addTorsion(i, i+1, i+2, i+3, 1.0, 1.1, 1.2, 0.3, 0.4, 0.5);
system.addForce(rb);
ReferencePlatform ref;
VerletIntegrator integrator1(0.01);
Context context1(system, integrator1, ref);
VerletIntegrator integrator2(0.01);
Context context2(system, integrator2, platform);
OpenMM_SFMT::SFMT sfmt;
init_gen_rand(0, sfmt);
vector<Vec3> positions(numParticles);
for (int i = 0; i < numParticles; i++)
positions[i] = Vec3(i, genrand_real2(sfmt), genrand_real2(sfmt));
context1.setPositions(positions);
context2.setPositions(positions);
State state1 = context1.getState(State::Forces | State::Energy);
State state2 = context2.getState(State::Forces | State::Energy);
const vector<Vec3>& forces1 = state1.getForces();
const vector<Vec3>& forces2 = state2.getForces();
for (int i = 0; i < numParticles; i++)
ASSERT_EQUAL_VEC(forces1[i], forces2[i], TOL);
ASSERT_EQUAL_TOL(state1.getPotentialEnergy(), state2.getPotentialEnergy(), TOL);
}
int main(int argc, char* argv[]) {
try {
if (argc > 1)
platform.setPropertyDefaultValue("HipPrecision", string(argv[1]));
testForces();
}
catch(const exception& e) {
cout << "exception: " << e.what() << endl;
return 1;
}
cout << "Done" << endl;
return 0;
void runPlatformTests() {
}
platforms/hip/tests/TestHipVirtualSites.cpp
View file @ 462a9be3
......@@ -33,112 +33,5 @@
#include "HipTests.h"
#include "TestVirtualSites.h"
/**
* Make sure that atom reordering respects virtual sites.
*/
void testReordering() {
const double cutoff = 2.0;
const double boxSize = 20.0;
System system;
NonbondedForce* nonbonded = new NonbondedForce();
system.addForce(nonbonded);
nonbonded->setNonbondedMethod(NonbondedForce::CutoffNonPeriodic);
nonbonded->setCutoffDistance(cutoff);
vector<Vec3> positions;
OpenMM_SFMT::SFMT sfmt;
init_gen_rand(0, sfmt);
// Create linear molecules with TwoParticleAverage virtual sites.
for (int i = 0; i < 50; i++) {
int start = system.getNumParticles();
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(0.0);
system.setVirtualSite(start+2, new TwoParticleAverageSite(start, start+1, 0.4, 0.6));
system.addConstraint(start, start+1, 2.0);
for (int i = 0; i < 3; i++) {
nonbonded->addParticle(0, 0.2, 1);
for (int j = 0; j < i; j++)
nonbonded->addException(start+i, start+j, 0, 1, 0);
}
Vec3 pos(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
positions.push_back(pos);
positions.push_back(pos+Vec3(2, 0, 0));
positions.push_back(Vec3());
}
// Create planar molecules with ThreeParticleAverage virtual sites.
for (int i = 0; i < 50; i++) {
int start = system.getNumParticles();
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(0.0);
system.setVirtualSite(start+3, new ThreeParticleAverageSite(start, start+1, start+2, 0.3, 0.5, 0.2));
system.addConstraint(start, start+1, 1.0);
system.addConstraint(start, start+2, 1.0);
system.addConstraint(start+1, start+2, sqrt(2.0));
for (int i = 0; i < 4; i++) {
nonbonded->addParticle(0, 0.2, 1);
for (int j = 0; j < i; j++)
nonbonded->addException(start+i, start+j, 0, 1, 0);
}
Vec3 pos(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
positions.push_back(pos);
positions.push_back(pos+Vec3(1, 0, 0));
positions.push_back(pos+Vec3(0, 1, 0));
positions.push_back(Vec3());
}
// Create tetrahedral molecules with OutOfPlane virtual sites.
for (int i = 0; i < 50; i++) {
int start = system.getNumParticles();
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(0.0);
system.setVirtualSite(start+3, new OutOfPlaneSite(start, start+1, start+2, 0.3, 0.5, 0.2));
system.addConstraint(start, start+1, 1.0);
system.addConstraint(start, start+2, 1.0);
system.addConstraint(start+1, start+2, sqrt(2.0));
for (int i = 0; i < 4; i++) {
nonbonded->addParticle(0, 0.2, 1);
for (int j = 0; j < i; j++)
nonbonded->addException(start+i, start+j, 0, 1, 0);
}
Vec3 pos(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
positions.push_back(pos);
positions.push_back(pos+Vec3(1, 0, 0));
positions.push_back(pos+Vec3(0, 1, 0));
positions.push_back(Vec3());
}
// Simulate it and check conservation laws.
LangevinIntegrator integrator(300.0, 0.1, 0.002);
Context context(system, integrator, platform);
context.setPositions(positions);
context.applyConstraints(0.0001);
for (int i = 0; i < 1000; i++) {
State state = context.getState(State::Positions);
const vector<Vec3>& pos = state.getPositions();
for (int j = 0; j < 150; j += 3)
ASSERT_EQUAL_VEC(pos[j]*0.4+pos[j+1]*0.6, pos[j+2], 1e-5);
for (int j = 150; j < 350; j += 4)
ASSERT_EQUAL_VEC(pos[j]*0.3+pos[j+1]*0.5+pos[j+2]*0.2, pos[j+3], 1e-5);
for (int j = 350; j < 550; j += 4) {
Vec3 v12 = pos[j+1]-pos[j];
Vec3 v13 = pos[j+2]-pos[j];
Vec3 cross = v12.cross(v13);
ASSERT_EQUAL_VEC(pos[j]+v12*0.3+v13*0.5+cross*0.2, pos[j+3], 1e-5);
}
integrator.step(1);
}
}
void runPlatformTests() {
testReordering();
}
platforms/opencl/tests/TestOpenCLCheckpoints.cpp
View file @ 462a9be3
......@@ -32,108 +32,5 @@
#include "OpenCLTests.h"
#include "TestCheckpoints.h"
void testCheckpoint() {
const int numParticles = 100;
const double boxSize = 5.0;
const double temperature = 200.0;
System system;
system.addForce(new AndersenThermostat(0.0, 100.0));
NonbondedForce* nonbonded = new NonbondedForce();
system.addForce(nonbonded);
nonbonded->setNonbondedMethod(NonbondedForce::CutoffPeriodic);
vector<Vec3> positions(numParticles);
OpenMM_SFMT::SFMT sfmt;
init_gen_rand(0, sfmt);
for (int i = 0; i < numParticles; i++) {
system.addParticle(1.0);
nonbonded->addParticle(i%2 == 0 ? 0.1 : -0.1, 0.2, 0.1);
bool clash;
do {
clash = false;
positions[i] = Vec3(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
for (int j = 0; j < i; j++) {
Vec3 delta = positions[i]-positions[j];
if (sqrt(delta.dot(delta)) < 0.1)
clash = true;
}
} while (clash);
}
VerletIntegrator integrator(0.001);
Context context(system, integrator, platform);
context.setPositions(positions);
context.setPeriodicBoxVectors(Vec3(boxSize, 0, 0), Vec3(0, boxSize, 0), Vec3(0, 0, boxSize));
context.setParameter(AndersenThermostat::Temperature(), temperature);
// Run for a little while.
integrator.step(100);
// Record the current state and make a checkpoint.
State s1 = context.getState(State::Positions | State::Velocities | State::Parameters);
stringstream stream1(ios_base::out | ios_base::in | ios_base::binary);
context.createCheckpoint(stream1);
// Continue the simulation for a few more steps and record the state again.
integrator.step(10);
State s2 = context.getState(State::Positions | State::Velocities | State::Parameters);
// Restore from the checkpoint and see if everything gets restored correctly.
context.setPeriodicBoxVectors(Vec3(2*boxSize, 0, 0), Vec3(0, 2*boxSize, 0), Vec3(0, 0, 2*boxSize));
context.setParameter(AndersenThermostat::Temperature(), temperature+10);
context.loadCheckpoint(stream1);
State s3 = context.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s1, s3);
// Now simulate from there and see if the trajectory is identical.
integrator.step(10);
State s4 = context.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s2, s4);
// Create a new Context that uses multiple devices.
map<string, string> props;
string deviceIndex = platform.getPropertyValue(context, OpenCLPlatform::OpenCLDeviceIndex());
props[OpenCLPlatform::OpenCLDeviceIndex()] = deviceIndex+","+deviceIndex;
string platformIndex = platform.getPropertyValue(context, OpenCLPlatform::OpenCLPlatformIndex());
props[OpenCLPlatform::OpenCLPlatformIndex()] = platformIndex;
VerletIntegrator integrator2(0.001);
Context context2(system, integrator2, platform, props);
context2.setPositions(positions);
context2.setPeriodicBoxVectors(Vec3(boxSize, 0, 0), Vec3(0, boxSize, 0), Vec3(0, 0, boxSize));
context2.setParameter(AndersenThermostat::Temperature(), temperature);
// Now repeat all of the above tests with it.
integrator2.step(100);
State s5 = context2.getState(State::Positions | State::Velocities | State::Parameters);
stringstream stream2(ios_base::out | ios_base::in | ios_base::binary);
context2.createCheckpoint(stream2);
integrator2.step(10);
State s6 = context2.getState(State::Positions | State::Velocities | State::Parameters);
context2.setPeriodicBoxVectors(Vec3(2*boxSize, 0, 0), Vec3(0, 2*boxSize, 0), Vec3(0, 0, 2*boxSize));
context2.setParameter(AndersenThermostat::Temperature(), temperature+10);
context2.loadCheckpoint(stream2);
State s7 = context2.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s5, s7);
integrator2.step(10);
State s8 = context2.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s6, s8);
// See if a checkpoint created from one Context can be loaded into a different one.
VerletIntegrator integrator3(0.001);
Context context3(system, integrator3, platform);
stream1.seekg(0, stream1.beg);
context3.loadCheckpoint(stream1);
State s9 = context3.getState(State::Positions | State::Velocities | State::Parameters | State::Energy);
compareStates(s1, s9);
}
void runPlatformTests() {
testCheckpoint();
}
platforms/opencl/tests/TestOpenCLCustomIntegrator.cpp
View file @ 462a9be3
......@@ -32,72 +32,5 @@
#include "OpenCLTests.h"
#include "TestCustomIntegrator.h"
/**
* Make sure random numbers are computed correctly when steps get merged.
*/
void testMergedRandoms() {
const int numParticles = 10;
const int numSteps = 10;
System system;
for (int i = 0; i < numParticles; i++)
system.addParticle(1.0);
CustomIntegrator integrator(0.1);
integrator.addPerDofVariable("dofUniform1", 0);
integrator.addPerDofVariable("dofUniform2", 0);
integrator.addPerDofVariable("dofGaussian1", 0);
integrator.addPerDofVariable("dofGaussian2", 0);
integrator.addGlobalVariable("globalUniform1", 0);
integrator.addGlobalVariable("globalUniform2", 0);
integrator.addGlobalVariable("globalGaussian1", 0);
integrator.addGlobalVariable("globalGaussian2", 0);
integrator.addComputePerDof("dofUniform1", "uniform");
integrator.addComputePerDof("dofUniform2", "uniform");
integrator.addComputePerDof("dofGaussian1", "gaussian");
integrator.addComputePerDof("dofGaussian2", "gaussian");
integrator.addComputeGlobal("globalUniform1", "uniform");
integrator.addComputeGlobal("globalUniform2", "uniform");
integrator.addComputeGlobal("globalGaussian1", "gaussian");
integrator.addComputeGlobal("globalGaussian2", "gaussian");
Context context(system, integrator, platform);
// See if the random numbers are computed correctly.
vector<Vec3> values1, values2;
for (int i = 0; i < numSteps; i++) {
integrator.step(1);
integrator.getPerDofVariable(0, values1);
integrator.getPerDofVariable(1, values2);
for (int i = 0; i < numParticles; i++)
for (int j = 0; j < 3; j++) {
double v1 = values1[i][j];
double v2 = values2[i][j];
ASSERT(v1 >= 0 && v1 < 1);
ASSERT(v2 >= 0 && v2 < 1);
ASSERT(v1 != v2);
}
integrator.getPerDofVariable(2, values1);
integrator.getPerDofVariable(3, values2);
for (int i = 0; i < numParticles; i++)
for (int j = 0; j < 3; j++) {
double v1 = values1[i][j];
double v2 = values2[i][j];
ASSERT(v1 >= -10 && v1 < 10);
ASSERT(v2 >= -10 && v2 < 10);
ASSERT(v1 != v2);
}
double v1 = integrator.getGlobalVariable(0);
double v2 = integrator.getGlobalVariable(1);
ASSERT(v1 >= 0 && v1 < 1);
ASSERT(v2 >= 0 && v2 < 1);
ASSERT(v1 != v2);
v1 = integrator.getGlobalVariable(2);
v2 = integrator.getGlobalVariable(3);
ASSERT(v1 >= -10 && v1 < 10);
ASSERT(v2 >= -10 && v2 < 10);
ASSERT(v1 != v2);
}
}
void runPlatformTests() {
testMergedRandoms();
}
platforms/opencl/tests/TestOpenCLMultipleForces.cpp
View file @ 462a9be3
......@@ -29,86 +29,8 @@
* USE OR OTHER DEALINGS IN THE SOFTWARE. *
* -------------------------------------------------------------------------- */
/**
* This tests a system with multiple forces, to make sure OpenCLBondedUtilities is
* processing them correctly.
*/
#include "OpenCLTests.h"
#include "TestMultipleForces.h"
#include "openmm/internal/AssertionUtilities.h"
#include "openmm/Context.h"
#include "OpenCLPlatform.h"
#include "ReferencePlatform.h"
#include "openmm/HarmonicAngleForce.h"
#include "openmm/HarmonicBondForce.h"
#include "openmm/PeriodicTorsionForce.h"
#include "openmm/RBTorsionForce.h"
#include "openmm/System.h"
#include "openmm/VerletIntegrator.h"
#include "SimTKOpenMMRealType.h"
#include "sfmt/SFMT.h"
#include <iostream>
#include <vector>
using namespace OpenMM;
using namespace std;
static OpenCLPlatform platform;
const double TOL = 1e-4;
void testForces() {
const int numParticles = 100;
System system;
for (int i = 0; i < numParticles; i++)
system.addParticle(1.0);
HarmonicBondForce* bonds = new HarmonicBondForce();
for (int i = 0; i < numParticles-1; i++)
bonds->addBond(i, i+1, 1.0, 1.5);
system.addForce(bonds);
HarmonicAngleForce* angles = new HarmonicAngleForce();
for (int i = 0; i < numParticles-2; i++)
angles->addAngle(i, i+1, i+2, 2.0, 1.5);
system.addForce(angles);
PeriodicTorsionForce* periodic = new PeriodicTorsionForce();
for (int i = 0; i < numParticles-3; i++)
periodic->addTorsion(i, i+1, i+2, i+3, 2, PI_M/3, 1.1);
system.addForce(periodic);
RBTorsionForce* rb = new RBTorsionForce();
for (int i = 0; i < numParticles-3; i += 3)
rb->addTorsion(i, i+1, i+2, i+3, 1.0, 1.1, 1.2, 0.3, 0.4, 0.5);
system.addForce(rb);
ReferencePlatform ref;
VerletIntegrator integrator1(0.01);
Context context1(system, integrator1, ref);
VerletIntegrator integrator2(0.01);
Context context2(system, integrator2, platform);
OpenMM_SFMT::SFMT sfmt;
init_gen_rand(0, sfmt);
vector<Vec3> positions(numParticles);
for (int i = 0; i < numParticles; i++)
positions[i] = Vec3(i, genrand_real2(sfmt), genrand_real2(sfmt));
context1.setPositions(positions);
context2.setPositions(positions);
State state1 = context1.getState(State::Forces | State::Energy);
State state2 = context2.getState(State::Forces | State::Energy);
const vector<Vec3>& forces1 = state1.getForces();
const vector<Vec3>& forces2 = state2.getForces();
for (int i = 0; i < numParticles; i++)
ASSERT_EQUAL_VEC(forces1[i], forces2[i], TOL);
ASSERT_EQUAL_TOL(state1.getPotentialEnergy(), state2.getPotentialEnergy(), TOL);
void runPlatformTests() {
}
int main(int argc, char* argv[]) {
try {
if (argc > 1)
platform.setPropertyDefaultValue("OpenCLPrecision", string(argv[1]));
testForces();
}
catch(const exception& e) {
cout << "exception: " << e.what() << endl;
return 1;
}
cout << "Done" << endl;
return 0;
}
platforms/opencl/tests/TestOpenCLVirtualSites.cpp
View file @ 462a9be3
......@@ -32,112 +32,5 @@
#include "OpenCLTests.h"
#include "TestVirtualSites.h"
/**
* Make sure that atom reordering respects virtual sites.
*/
void testReordering() {
const double cutoff = 2.0;
const double boxSize = 20.0;
System system;
NonbondedForce* nonbonded = new NonbondedForce();
system.addForce(nonbonded);
nonbonded->setNonbondedMethod(NonbondedForce::CutoffNonPeriodic);
nonbonded->setCutoffDistance(cutoff);
vector<Vec3> positions;
OpenMM_SFMT::SFMT sfmt;
init_gen_rand(0, sfmt);
// Create linear molecules with TwoParticleAverage virtual sites.
for (int i = 0; i < 50; i++) {
int start = system.getNumParticles();
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(0.0);
system.setVirtualSite(start+2, new TwoParticleAverageSite(start, start+1, 0.4, 0.6));
system.addConstraint(start, start+1, 2.0);
for (int i = 0; i < 3; i++) {
nonbonded->addParticle(0, 0.2, 1);
for (int j = 0; j < i; j++)
nonbonded->addException(start+i, start+j, 0, 1, 0);
}
Vec3 pos(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
positions.push_back(pos);
positions.push_back(pos+Vec3(2, 0, 0));
positions.push_back(Vec3());
}
// Create planar molecules with ThreeParticleAverage virtual sites.
for (int i = 0; i < 50; i++) {
int start = system.getNumParticles();
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(0.0);
system.setVirtualSite(start+3, new ThreeParticleAverageSite(start, start+1, start+2, 0.3, 0.5, 0.2));
system.addConstraint(start, start+1, 1.0);
system.addConstraint(start, start+2, 1.0);
system.addConstraint(start+1, start+2, sqrt(2.0));
for (int i = 0; i < 4; i++) {
nonbonded->addParticle(0, 0.2, 1);
for (int j = 0; j < i; j++)
nonbonded->addException(start+i, start+j, 0, 1, 0);
}
Vec3 pos(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
positions.push_back(pos);
positions.push_back(pos+Vec3(1, 0, 0));
positions.push_back(pos+Vec3(0, 1, 0));
positions.push_back(Vec3());
}
// Create tetrahedral molecules with OutOfPlane virtual sites.
for (int i = 0; i < 50; i++) {
int start = system.getNumParticles();
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(0.0);
system.setVirtualSite(start+3, new OutOfPlaneSite(start, start+1, start+2, 0.3, 0.5, 0.2));
system.addConstraint(start, start+1, 1.0);
system.addConstraint(start, start+2, 1.0);
system.addConstraint(start+1, start+2, sqrt(2.0));
for (int i = 0; i < 4; i++) {
nonbonded->addParticle(0, 0.2, 1);
for (int j = 0; j < i; j++)
nonbonded->addException(start+i, start+j, 0, 1, 0);
}
Vec3 pos(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
positions.push_back(pos);
positions.push_back(pos+Vec3(1, 0, 0));
positions.push_back(pos+Vec3(0, 1, 0));
positions.push_back(Vec3());
}
// Simulate it and check conservation laws.
LangevinIntegrator integrator(300.0, 0.1, 0.002);
Context context(system, integrator, platform);
context.setPositions(positions);
context.applyConstraints(0.0001);
for (int i = 0; i < 1000; i++) {
State state = context.getState(State::Positions);
const vector<Vec3>& pos = state.getPositions();
for (int j = 0; j < 150; j += 3)
ASSERT_EQUAL_VEC(pos[j]*0.4+pos[j+1]*0.6, pos[j+2], 1e-5);
for (int j = 150; j < 350; j += 4)
ASSERT_EQUAL_VEC(pos[j]*0.3+pos[j+1]*0.5+pos[j+2]*0.2, pos[j+3], 1e-5);
for (int j = 350; j < 550; j += 4) {
Vec3 v12 = pos[j+1]-pos[j];
Vec3 v13 = pos[j+2]-pos[j];
Vec3 cross = v12.cross(v13);
ASSERT_EQUAL_VEC(pos[j]+v12*0.3+v13*0.5+cross*0.2, pos[j+3], 1e-5);
}
integrator.step(1);
}
}
void runPlatformTests() {
testReordering();
}
tests/TestCheckpoints.h
View file @ 462a9be3
......@@ -6,7 +6,7 @@
* Biological Structures at Stanford, funded under the NIH Roadmap for *
* Medical Research, grant U54 GM072970. See https://simtk.org. *
* *
* Portions copyright (c) 2012-2015 Stanford University and the Authors. *
* Portions copyright (c) 2012-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -116,12 +116,117 @@ void testSetState() {
}
}
void testMultipleDevices() {
const int numParticles = 100;
const double boxSize = 5.0;
const double temperature = 200.0;
System system;
system.addForce(new AndersenThermostat(0.0, 100.0));
NonbondedForce* nonbonded = new NonbondedForce();
system.addForce(nonbonded);
nonbonded->setNonbondedMethod(NonbondedForce::CutoffPeriodic);
vector<Vec3> positions(numParticles);
OpenMM_SFMT::SFMT sfmt;
init_gen_rand(0, sfmt);
for (int i = 0; i < numParticles; i++) {
system.addParticle(1.0);
nonbonded->addParticle(i%2 == 0 ? 0.1 : -0.1, 0.2, 0.1);
bool clash;
do {
clash = false;
positions[i] = Vec3(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
for (int j = 0; j < i; j++) {
Vec3 delta = positions[i]-positions[j];
if (sqrt(delta.dot(delta)) < 0.1)
clash = true;
}
} while (clash);
}
VerletIntegrator integrator(0.001);
Context context(system, integrator, platform);
context.setPositions(positions);
context.setPeriodicBoxVectors(Vec3(boxSize, 0, 0), Vec3(0, boxSize, 0), Vec3(0, 0, boxSize));
context.setParameter(AndersenThermostat::Temperature(), temperature);
// Run for a little while.
integrator.step(100);
// Record the current state and make a checkpoint.
State s1 = context.getState(State::Positions | State::Velocities | State::Parameters);
stringstream stream1(ios_base::out | ios_base::in | ios_base::binary);
context.createCheckpoint(stream1);
// Continue the simulation for a few more steps and record the state again.
integrator.step(10);
State s2 = context.getState(State::Positions | State::Velocities | State::Parameters);
// Restore from the checkpoint and see if everything gets restored correctly.
context.setPeriodicBoxVectors(Vec3(2*boxSize, 0, 0), Vec3(0, 2*boxSize, 0), Vec3(0, 0, 2*boxSize));
context.setParameter(AndersenThermostat::Temperature(), temperature+10);
context.loadCheckpoint(stream1);
State s3 = context.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s1, s3);
// Now simulate from there and see if the trajectory is identical.
integrator.step(10);
State s4 = context.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s2, s4);
// Create a new Context that uses multiple devices.
map<string, string> props;
try {
string deviceIndex = platform.getPropertyValue(context, "DeviceIndex");
props["DeviceIndex"] = deviceIndex+","+deviceIndex;
}
catch (OpenMMException& ex) {
// This platform doesn't have a DeviceIndex property.
}
VerletIntegrator integrator2(0.001);
Context context2(system, integrator2, platform, props);
context2.setPositions(positions);
context2.setPeriodicBoxVectors(Vec3(boxSize, 0, 0), Vec3(0, boxSize, 0), Vec3(0, 0, boxSize));
context2.setParameter(AndersenThermostat::Temperature(), temperature);
// Now repeat all of the above tests with it.
integrator2.step(100);
State s5 = context2.getState(State::Positions | State::Velocities | State::Parameters);
stringstream stream2(ios_base::out | ios_base::in | ios_base::binary);
context2.createCheckpoint(stream2);
integrator2.step(10);
State s6 = context2.getState(State::Positions | State::Velocities | State::Parameters);
context2.setPeriodicBoxVectors(Vec3(2*boxSize, 0, 0), Vec3(0, 2*boxSize, 0), Vec3(0, 0, 2*boxSize));
context2.setParameter(AndersenThermostat::Temperature(), temperature+10);
context2.loadCheckpoint(stream2);
State s7 = context2.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s5, s7);
integrator2.step(10);
State s8 = context2.getState(State::Positions | State::Velocities | State::Parameters);
compareStates(s6, s8);
// See if a checkpoint created from one Context can be loaded into a different one.
VerletIntegrator integrator3(0.001);
Context context3(system, integrator3, platform);
stream1.seekg(0, stream1.beg);
context3.loadCheckpoint(stream1);
State s9 = context3.getState(State::Positions | State::Velocities | State::Parameters | State::Energy);
compareStates(s1, s9);
}
void runPlatformTests();
int main(int argc, char* argv[]) {
try {
initializeTests(argc, argv);
testSetState();
testMultipleDevices();
runPlatformTests();
}
catch(const exception& e) {
......
tests/TestCustomIntegrator.h
View file @ 462a9be3
......@@ -1283,6 +1283,72 @@ void testAnalyzeComputations() {
}
}
/**
* Make sure random numbers are computed correctly when steps get merged.
*/
void testMergedRandoms() {
const int numParticles = 10;
const int numSteps = 10;
System system;
for (int i = 0; i < numParticles; i++)
system.addParticle(1.0);
CustomIntegrator integrator(0.1);
integrator.addPerDofVariable("dofUniform1", 0);
integrator.addPerDofVariable("dofUniform2", 0);
integrator.addPerDofVariable("dofGaussian1", 0);
integrator.addPerDofVariable("dofGaussian2", 0);
integrator.addGlobalVariable("globalUniform1", 0);
integrator.addGlobalVariable("globalUniform2", 0);
integrator.addGlobalVariable("globalGaussian1", 0);
integrator.addGlobalVariable("globalGaussian2", 0);
integrator.addComputePerDof("dofUniform1", "uniform");
integrator.addComputePerDof("dofUniform2", "uniform");
integrator.addComputePerDof("dofGaussian1", "gaussian");
integrator.addComputePerDof("dofGaussian2", "gaussian");
integrator.addComputeGlobal("globalUniform1", "uniform");
integrator.addComputeGlobal("globalUniform2", "uniform");
integrator.addComputeGlobal("globalGaussian1", "gaussian");
integrator.addComputeGlobal("globalGaussian2", "gaussian");
Context context(system, integrator, platform);
// See if the random numbers are computed correctly.
vector<Vec3> values1, values2;
for (int i = 0; i < numSteps; i++) {
integrator.step(1);
integrator.getPerDofVariable(0, values1);
integrator.getPerDofVariable(1, values2);
for (int i = 0; i < numParticles; i++)
for (int j = 0; j < 3; j++) {
double v1 = values1[i][j];
double v2 = values2[i][j];
ASSERT(v1 >= 0 && v1 < 1);
ASSERT(v2 >= 0 && v2 < 1);
ASSERT(v1 != v2);
}
integrator.getPerDofVariable(2, values1);
integrator.getPerDofVariable(3, values2);
for (int i = 0; i < numParticles; i++)
for (int j = 0; j < 3; j++) {
double v1 = values1[i][j];
double v2 = values2[i][j];
ASSERT(v1 >= -10 && v1 < 10);
ASSERT(v2 >= -10 && v2 < 10);
ASSERT(v1 != v2);
}
double v1 = integrator.getGlobalVariable(0);
double v2 = integrator.getGlobalVariable(1);
ASSERT(v1 >= 0 && v1 < 1);
ASSERT(v2 >= 0 && v2 < 1);
ASSERT(v1 != v2);
v1 = integrator.getGlobalVariable(2);
v2 = integrator.getGlobalVariable(3);
ASSERT(v1 >= -10 && v1 < 10);
ASSERT(v2 >= -10 && v2 < 10);
ASSERT(v1 != v2);
}
}
void runPlatformTests();
int main(int argc, char* argv[]) {
......@@ -1314,6 +1380,7 @@ int main(int argc, char* argv[]) {
testCheckpoint();
testSaveParameters();
testAnalyzeComputations();
testMergedRandoms();
runPlatformTests();
}
catch(const exception& e) {
......
tests/TestMultipleForces.h 0 → 100644
View file @ 462a9be3
/* -------------------------------------------------------------------------- *
* OpenMM *
* -------------------------------------------------------------------------- *
* This is part of the OpenMM molecular simulation toolkit originating from *
* Simbios, the NIH National Center for Physics-Based Simulation of *
* Biological Structures at Stanford, funded under the NIH Roadmap for *
* Medical Research, grant U54 GM072970. See https://simtk.org. *
* *
* Portions copyright (c) 2008-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
* Permission is hereby granted, free of charge, to any person obtaining a *
* copy of this software and associated documentation files (the "Software"), *
* to deal in the Software without restriction, including without limitation *
* the rights to use, copy, modify, merge, publish, distribute, sublicense, *
* and/or sell copies of the Software, and to permit persons to whom the *
* Software is furnished to do so, subject to the following conditions: *
* *
* The above copyright notice and this permission notice shall be included in *
* all copies or substantial portions of the Software. *
* *
* THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR *
* IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, *
* FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL *
* THE AUTHORS, CONTRIBUTORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, *
* DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR *
* OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE *
* USE OR OTHER DEALINGS IN THE SOFTWARE. *
* -------------------------------------------------------------------------- */
#include "openmm/internal/AssertionUtilities.h"
#include "openmm/Context.h"
#include "ReferencePlatform.h"
#include "openmm/HarmonicAngleForce.h"
#include "openmm/HarmonicBondForce.h"
#include "openmm/PeriodicTorsionForce.h"
#include "openmm/RBTorsionForce.h"
#include "openmm/System.h"
#include "openmm/VerletIntegrator.h"
#include "SimTKOpenMMRealType.h"
#include "sfmt/SFMT.h"
#include <iostream>
#include <vector>
using namespace OpenMM;
using namespace std;
const double TOL = 1e-4;
/**
* This tests a system with multiple forces, to make sure BondedUtilities is
* processing them correctly.
*/
void testForces() {
const int numParticles = 100;
System system;
for (int i = 0; i < numParticles; i++)
system.addParticle(1.0);
HarmonicBondForce* bonds = new HarmonicBondForce();
for (int i = 0; i < numParticles-1; i++)
bonds->addBond(i, i+1, 1.0, 1.5);
system.addForce(bonds);
HarmonicAngleForce* angles = new HarmonicAngleForce();
for (int i = 0; i < numParticles-2; i++)
angles->addAngle(i, i+1, i+2, 2.0, 1.5);
system.addForce(angles);
PeriodicTorsionForce* periodic = new PeriodicTorsionForce();
for (int i = 0; i < numParticles-3; i++)
periodic->addTorsion(i, i+1, i+2, i+3, 2, PI_M/3, 1.1);
system.addForce(periodic);
RBTorsionForce* rb = new RBTorsionForce();
for (int i = 0; i < numParticles-3; i += 3)
rb->addTorsion(i, i+1, i+2, i+3, 1.0, 1.1, 1.2, 0.3, 0.4, 0.5);
system.addForce(rb);
ReferencePlatform ref;
VerletIntegrator integrator1(0.01);
Context context1(system, integrator1, ref);
VerletIntegrator integrator2(0.01);
Context context2(system, integrator2, platform);
OpenMM_SFMT::SFMT sfmt;
init_gen_rand(0, sfmt);
vector<Vec3> positions(numParticles);
for (int i = 0; i < numParticles; i++)
positions[i] = Vec3(i, genrand_real2(sfmt), genrand_real2(sfmt));
context1.setPositions(positions);
context2.setPositions(positions);
State state1 = context1.getState(State::Forces | State::Energy);
State state2 = context2.getState(State::Forces | State::Energy);
const vector<Vec3>& forces1 = state1.getForces();
const vector<Vec3>& forces2 = state2.getForces();
for (int i = 0; i < numParticles; i++)
ASSERT_EQUAL_VEC(forces1[i], forces2[i], TOL);
ASSERT_EQUAL_TOL(state1.getPotentialEnergy(), state2.getPotentialEnergy(), TOL);
}
void runPlatformTests();
int main(int argc, char* argv[]) {
try {
initializeTests(argc, argv);
testForces();
runPlatformTests();
}
catch(const exception& e) {
cout << "exception: " << e.what() << endl;
return 1;
}
cout << "Done" << endl;
return 0;
}
tests/TestVirtualSites.h
View file @ 462a9be3
......@@ -6,7 +6,7 @@
* Biological Structures at Stanford, funded under the NIH Roadmap for *
* Medical Research, grant U54 GM072970. See https://simtk.org. *
* *
* Portions copyright (c) 2012-2017 Stanford University and the Authors. *
* Portions copyright (c) 2012-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -512,6 +512,112 @@ void testNestedSites() {
ASSERT_EQUAL_VEC(Vec3(1*0.25 + 2*0.75, 0, 0), state.getForces()[4], 1e-6);
}
/**
* Make sure that atom reordering respects virtual sites.
*/
void testReordering() {
const double cutoff = 2.0;
const double boxSize = 20.0;
System system;
NonbondedForce* nonbonded = new NonbondedForce();
system.addForce(nonbonded);
nonbonded->setNonbondedMethod(NonbondedForce::CutoffNonPeriodic);
nonbonded->setCutoffDistance(cutoff);
vector<Vec3> positions;
OpenMM_SFMT::SFMT sfmt;
init_gen_rand(0, sfmt);
// Create linear molecules with TwoParticleAverage virtual sites.
for (int i = 0; i < 50; i++) {
int start = system.getNumParticles();
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(0.0);
system.setVirtualSite(start+2, new TwoParticleAverageSite(start, start+1, 0.4, 0.6));
system.addConstraint(start, start+1, 2.0);
for (int i = 0; i < 3; i++) {
nonbonded->addParticle(0, 0.2, 1);
for (int j = 0; j < i; j++)
nonbonded->addException(start+i, start+j, 0, 1, 0);
}
Vec3 pos(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
positions.push_back(pos);
positions.push_back(pos+Vec3(2, 0, 0));
positions.push_back(Vec3());
}
// Create planar molecules with ThreeParticleAverage virtual sites.
for (int i = 0; i < 50; i++) {
int start = system.getNumParticles();
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(0.0);
system.setVirtualSite(start+3, new ThreeParticleAverageSite(start, start+1, start+2, 0.3, 0.5, 0.2));
system.addConstraint(start, start+1, 1.0);
system.addConstraint(start, start+2, 1.0);
system.addConstraint(start+1, start+2, sqrt(2.0));
for (int i = 0; i < 4; i++) {
nonbonded->addParticle(0, 0.2, 1);
for (int j = 0; j < i; j++)
nonbonded->addException(start+i, start+j, 0, 1, 0);
}
Vec3 pos(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
positions.push_back(pos);
positions.push_back(pos+Vec3(1, 0, 0));
positions.push_back(pos+Vec3(0, 1, 0));
positions.push_back(Vec3());
}
// Create tetrahedral molecules with OutOfPlane virtual sites.
for (int i = 0; i < 50; i++) {
int start = system.getNumParticles();
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(1.0);
system.addParticle(0.0);
system.setVirtualSite(start+3, new OutOfPlaneSite(start, start+1, start+2, 0.3, 0.5, 0.2));
system.addConstraint(start, start+1, 1.0);
system.addConstraint(start, start+2, 1.0);
system.addConstraint(start+1, start+2, sqrt(2.0));
for (int i = 0; i < 4; i++) {
nonbonded->addParticle(0, 0.2, 1);
for (int j = 0; j < i; j++)
nonbonded->addException(start+i, start+j, 0, 1, 0);
}
Vec3 pos(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
positions.push_back(pos);
positions.push_back(pos+Vec3(1, 0, 0));
positions.push_back(pos+Vec3(0, 1, 0));
positions.push_back(Vec3());
}
// Simulate it and check conservation laws.
LangevinIntegrator integrator(300.0, 0.1, 0.002);
Context context(system, integrator, platform);
context.setPositions(positions);
context.applyConstraints(0.0001);
for (int i = 0; i < 1000; i++) {
State state = context.getState(State::Positions);
const vector<Vec3>& pos = state.getPositions();
for (int j = 0; j < 150; j += 3)
ASSERT_EQUAL_VEC(pos[j]*0.4+pos[j+1]*0.6, pos[j+2], 1e-5);
for (int j = 150; j < 350; j += 4)
ASSERT_EQUAL_VEC(pos[j]*0.3+pos[j+1]*0.5+pos[j+2]*0.2, pos[j+3], 1e-5);
for (int j = 350; j < 550; j += 4) {
Vec3 v12 = pos[j+1]-pos[j];
Vec3 v13 = pos[j+2]-pos[j];
Vec3 cross = v12.cross(v13);
ASSERT_EQUAL_VEC(pos[j]+v12*0.3+v13*0.5+cross*0.2, pos[j+3], 1e-5);
}
integrator.step(1);
}
}
void runPlatformTests();
int main(int argc, char* argv[]) {
......@@ -527,6 +633,7 @@ int main(int argc, char* argv[]) {
testConservationLaws();
testOverlappingSites();
testNestedSites();
testReordering();
runPlatformTests();
}
catch(const exception& e) {
......
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