Unverified Commit 78902bed authored by Peter Eastman's avatar Peter Eastman Committed by GitHub
Browse files

Optimize updateParametersInContext() (#4610)

* Optimize CustomNonbondedForce.updateParametersInContext()

* Optimized uploading changed values to GPU

* Optimized updateParametersInContext() for lots of bonded forces

* Optimized updateParametersInContext() for CustomExternalForce

* Optimized updateParametersInContext() for NonbondedForce

* Code changes for HIP platform
parent 66064fac
......@@ -540,30 +540,33 @@ double CommonCalcHarmonicBondForceKernel::execute(ContextImpl& context, bool inc
return 0.0;
}
void CommonCalcHarmonicBondForceKernel::copyParametersToContext(ContextImpl& context, const HarmonicBondForce& force) {
void CommonCalcHarmonicBondForceKernel::copyParametersToContext(ContextImpl& context, const HarmonicBondForce& force, int firstBond, int lastBond) {
ContextSelector selector(cc);
int numContexts = cc.getNumContexts();
int startIndex = cc.getContextIndex()*force.getNumBonds()/numContexts;
int endIndex = (cc.getContextIndex()+1)*force.getNumBonds()/numContexts;
if (numBonds != endIndex-startIndex)
throw OpenMMException("updateParametersInContext: The number of bonds has changed");
if (numBonds == 0)
if (numBonds == 0 || firstBond >= endIndex || lastBond < startIndex || firstBond > lastBond)
return;
firstBond = max(firstBond, startIndex);
lastBond = min(lastBond, endIndex-1);
// Record the per-bond parameters.
vector<mm_float2> paramVector(numBonds);
for (int i = 0; i < numBonds; i++) {
int numToSet = lastBond-firstBond+1;
vector<mm_float2> paramVector(numToSet);
for (int i = 0; i < numToSet; i++) {
int atom1, atom2;
double length, k;
force.getBondParameters(startIndex+i, atom1, atom2, length, k);
force.getBondParameters(firstBond+i, atom1, atom2, length, k);
paramVector[i] = mm_float2((float) length, (float) k);
}
params.upload(paramVector);
params.uploadSubArray(paramVector.data(), firstBond-startIndex, numToSet);
// Mark that the current reordering may be invalid.
cc.invalidateMolecules(info);
cc.invalidateMolecules(info, false, true);
}
class CommonCalcCustomBondForceKernel::ForceInfo : public ComputeForceInfo {
public:
......@@ -686,27 +689,30 @@ double CommonCalcCustomBondForceKernel::execute(ContextImpl& context, bool inclu
return 0.0;
}
void CommonCalcCustomBondForceKernel::copyParametersToContext(ContextImpl& context, const CustomBondForce& force) {
void CommonCalcCustomBondForceKernel::copyParametersToContext(ContextImpl& context, const CustomBondForce& force, int firstBond, int lastBond) {
ContextSelector selector(cc);
int numContexts = cc.getNumContexts();
int startIndex = cc.getContextIndex()*force.getNumBonds()/numContexts;
int endIndex = (cc.getContextIndex()+1)*force.getNumBonds()/numContexts;
if (numBonds != endIndex-startIndex)
throw OpenMMException("updateParametersInContext: The number of bonds has changed");
if (numBonds == 0)
if (numBonds == 0 || firstBond >= endIndex || lastBond < startIndex || firstBond > lastBond)
return;
firstBond = max(firstBond, startIndex);
lastBond = min(lastBond, endIndex-1);
// Record the per-bond parameters.
vector<vector<double> > paramVector(numBonds);
int numToSet = lastBond-firstBond+1;
vector<vector<double> > paramVector(numToSet);
int atom1, atom2;
for (int i = 0; i < numBonds; i++)
force.getBondParameters(startIndex+i, atom1, atom2, paramVector[i]);
params->setParameterValues(paramVector, true);
for (int i = 0; i < numToSet; i++)
force.getBondParameters(firstBond+i, atom1, atom2, paramVector[i]);
params->setParameterValuesSubset(firstBond-startIndex, paramVector, true);
// Mark that the current reordering may be invalid.
cc.invalidateMolecules(info);
cc.invalidateMolecules(info, false, true);
}
class CommonCalcHarmonicAngleForceKernel::ForceInfo : public ComputeForceInfo {
......@@ -767,30 +773,33 @@ double CommonCalcHarmonicAngleForceKernel::execute(ContextImpl& context, bool in
return 0.0;
}
void CommonCalcHarmonicAngleForceKernel::copyParametersToContext(ContextImpl& context, const HarmonicAngleForce& force) {
void CommonCalcHarmonicAngleForceKernel::copyParametersToContext(ContextImpl& context, const HarmonicAngleForce& force, int firstAngle, int lastAngle) {
ContextSelector selector(cc);
int numContexts = cc.getNumContexts();
int startIndex = cc.getContextIndex()*force.getNumAngles()/numContexts;
int endIndex = (cc.getContextIndex()+1)*force.getNumAngles()/numContexts;
if (numAngles != endIndex-startIndex)
throw OpenMMException("updateParametersInContext: The number of angles has changed");
if (numAngles == 0)
if (numAngles == 0 || firstAngle >= endIndex || lastAngle < startIndex || firstAngle > lastAngle)
return;
firstAngle = max(firstAngle, startIndex);
lastAngle = min(lastAngle, endIndex-1);
// Record the per-angle parameters.
vector<mm_float2> paramVector(numAngles);
for (int i = 0; i < numAngles; i++) {
int numToSet = lastAngle-firstAngle+1;
vector<mm_float2> paramVector(numToSet);
for (int i = 0; i < numToSet; i++) {
int atom1, atom2, atom3;
double angle, k;
force.getAngleParameters(startIndex+i, atom1, atom2, atom3, angle, k);
force.getAngleParameters(firstAngle+i, atom1, atom2, atom3, angle, k);
paramVector[i] = mm_float2((float) angle, (float) k);
}
params.upload(paramVector);
params.uploadSubArray(paramVector.data(), firstAngle-startIndex, numToSet);
// Mark that the current reordering may be invalid.
cc.invalidateMolecules();
cc.invalidateMolecules(info, false, true);
}
class CommonCalcCustomAngleForceKernel::ForceInfo : public ComputeForceInfo {
......@@ -915,27 +924,30 @@ double CommonCalcCustomAngleForceKernel::execute(ContextImpl& context, bool incl
return 0.0;
}
void CommonCalcCustomAngleForceKernel::copyParametersToContext(ContextImpl& context, const CustomAngleForce& force) {
void CommonCalcCustomAngleForceKernel::copyParametersToContext(ContextImpl& context, const CustomAngleForce& force, int firstAngle, int lastAngle) {
ContextSelector selector(cc);
int numContexts = cc.getNumContexts();
int startIndex = cc.getContextIndex()*force.getNumAngles()/numContexts;
int endIndex = (cc.getContextIndex()+1)*force.getNumAngles()/numContexts;
if (numAngles != endIndex-startIndex)
throw OpenMMException("updateParametersInContext: The number of angles has changed");
if (numAngles == 0)
if (numAngles == 0 || firstAngle >= endIndex || lastAngle < startIndex || firstAngle > lastAngle)
return;
firstAngle = max(firstAngle, startIndex);
lastAngle = min(lastAngle, endIndex-1);
// Record the per-angle parameters.
vector<vector<double> > paramVector(numAngles);
int numToSet = lastAngle-firstAngle+1;
vector<vector<double> > paramVector(numToSet);
int atom1, atom2, atom3;
for (int i = 0; i < numAngles; i++)
force.getAngleParameters(startIndex+i, atom1, atom2, atom3, paramVector[i]);
params->setParameterValues(paramVector, true);
for (int i = 0; i < numToSet; i++)
force.getAngleParameters(firstAngle+i, atom1, atom2, atom3, paramVector[i]);
params->setParameterValuesSubset(firstAngle-startIndex, paramVector, true);
// Mark that the current reordering may be invalid.
cc.invalidateMolecules(info);
cc.invalidateMolecules(info, false, true);
}
class CommonCalcPeriodicTorsionForceKernel::ForceInfo : public ComputeForceInfo {
......@@ -997,30 +1009,33 @@ double CommonCalcPeriodicTorsionForceKernel::execute(ContextImpl& context, bool
return 0.0;
}
void CommonCalcPeriodicTorsionForceKernel::copyParametersToContext(ContextImpl& context, const PeriodicTorsionForce& force) {
void CommonCalcPeriodicTorsionForceKernel::copyParametersToContext(ContextImpl& context, const PeriodicTorsionForce& force, int firstTorsion, int lastTorsion) {
ContextSelector selector(cc);
int numContexts = cc.getNumContexts();
int startIndex = cc.getContextIndex()*force.getNumTorsions()/numContexts;
int endIndex = (cc.getContextIndex()+1)*force.getNumTorsions()/numContexts;
if (numTorsions != endIndex-startIndex)
throw OpenMMException("updateParametersInContext: The number of torsions has changed");
if (numTorsions == 0)
if (numTorsions == 0 || firstTorsion >= endIndex || lastTorsion < startIndex || firstTorsion > lastTorsion)
return;
firstTorsion = max(firstTorsion, startIndex);
lastTorsion = min(lastTorsion, endIndex-1);
// Record the per-torsion parameters.
vector<mm_float4> paramVector(numTorsions);
for (int i = 0; i < numTorsions; i++) {
int numToSet = lastTorsion-firstTorsion+1;
vector<mm_float4> paramVector(numToSet);
for (int i = 0; i < numToSet; i++) {
int atom1, atom2, atom3, atom4, periodicity;
double phase, k;
force.getTorsionParameters(startIndex+i, atom1, atom2, atom3, atom4, periodicity, phase, k);
force.getTorsionParameters(firstTorsion+i, atom1, atom2, atom3, atom4, periodicity, phase, k);
paramVector[i] = mm_float4((float) k, (float) phase, (float) periodicity, 0.0f);
}
params.upload(paramVector);
params.uploadSubArray(paramVector.data(), firstTorsion-startIndex, numToSet);
// Mark that the current reordering may be invalid.
cc.invalidateMolecules();
cc.invalidateMolecules(info, false, true);
}
class CommonCalcRBTorsionForceKernel::ForceInfo : public ComputeForceInfo {
......@@ -1113,7 +1128,7 @@ void CommonCalcRBTorsionForceKernel::copyParametersToContext(ContextImpl& contex
// Mark that the current reordering may be invalid.
cc.invalidateMolecules();
cc.invalidateMolecules(info, false, true);
}
class CommonCalcCustomTorsionForceKernel::ForceInfo : public ComputeForceInfo {
......@@ -1238,27 +1253,30 @@ double CommonCalcCustomTorsionForceKernel::execute(ContextImpl& context, bool in
return 0.0;
}
void CommonCalcCustomTorsionForceKernel::copyParametersToContext(ContextImpl& context, const CustomTorsionForce& force) {
void CommonCalcCustomTorsionForceKernel::copyParametersToContext(ContextImpl& context, const CustomTorsionForce& force, int firstTorsion, int lastTorsion) {
ContextSelector selector(cc);
int numContexts = cc.getNumContexts();
int startIndex = cc.getContextIndex()*force.getNumTorsions()/numContexts;
int endIndex = (cc.getContextIndex()+1)*force.getNumTorsions()/numContexts;
if (numTorsions != endIndex-startIndex)
throw OpenMMException("updateParametersInContext: The number of torsions has changed");
if (numTorsions == 0)
if (numTorsions == 0 || firstTorsion >= endIndex || lastTorsion < startIndex || firstTorsion > lastTorsion)
return;
firstTorsion = max(firstTorsion, startIndex);
lastTorsion = min(lastTorsion, endIndex-1);
// Record the per-torsion parameters.
vector<vector<double> > paramVector(numTorsions);
int numToSet = lastTorsion-firstTorsion+1;
vector<vector<double> > paramVector(numToSet);
int atom1, atom2, atom3, atom4;
for (int i = 0; i < numTorsions; i++)
force.getTorsionParameters(startIndex+i, atom1, atom2, atom3, atom4, paramVector[i]);
params->setParameterValues(paramVector, true);
for (int i = 0; i < numToSet; i++)
force.getTorsionParameters(firstTorsion+i, atom1, atom2, atom3, atom4, paramVector[i]);
params->setParameterValuesSubset(firstTorsion-startIndex, paramVector, true);
// Mark that the current reordering may be invalid.
cc.invalidateMolecules(info);
cc.invalidateMolecules(info, false, true);
}
class CommonCalcCMAPTorsionForceKernel::ForceInfo : public ComputeForceInfo {
......@@ -1509,27 +1527,30 @@ double CommonCalcCustomExternalForceKernel::execute(ContextImpl& context, bool i
return 0.0;
}
void CommonCalcCustomExternalForceKernel::copyParametersToContext(ContextImpl& context, const CustomExternalForce& force) {
void CommonCalcCustomExternalForceKernel::copyParametersToContext(ContextImpl& context, const CustomExternalForce& force, int firstParticle, int lastParticle) {
ContextSelector selector(cc);
int numContexts = cc.getNumContexts();
int startIndex = cc.getContextIndex()*force.getNumParticles()/numContexts;
int endIndex = (cc.getContextIndex()+1)*force.getNumParticles()/numContexts;
if (numParticles != endIndex-startIndex)
throw OpenMMException("updateParametersInContext: The number of particles has changed");
if (numParticles == 0)
if (numParticles == 0 || firstParticle >= endIndex || lastParticle < startIndex || firstParticle > lastParticle)
return;
firstParticle = max(firstParticle, startIndex);
lastParticle = min(lastParticle, endIndex-1);
// Record the per-particle parameters.
vector<vector<double> > paramVector(numParticles);
int numToSet = lastParticle-firstParticle+1;
vector<vector<double> > paramVector(numToSet);
int particle;
for (int i = 0; i < numParticles; i++)
force.getParticleParameters(startIndex+i, particle, paramVector[i]);
params->setParameterValues(paramVector, true);
for (int i = 0; i < numToSet; i++)
force.getParticleParameters(firstParticle+i, particle, paramVector[i]);
params->setParameterValuesSubset(firstParticle-startIndex, paramVector, true);
// Mark that the current reordering may be invalid.
cc.invalidateMolecules(info);
cc.invalidateMolecules(info, true, false);
}
class CommonCalcCustomCompoundBondForceKernel::ForceInfo : public ComputeForceInfo {
......@@ -1718,7 +1739,7 @@ void CommonCalcCustomCompoundBondForceKernel::copyParametersToContext(ContextImp
// Mark that the current reordering may be invalid.
cc.invalidateMolecules(info);
cc.invalidateMolecules(info, false, true);
}
class CommonCalcCustomCentroidBondForceKernel::ForceInfo : public ComputeForceInfo {
......@@ -2037,7 +2058,7 @@ void CommonCalcCustomCentroidBondForceKernel::copyParametersToContext(ContextImp
// Mark that the current reordering may be invalid.
cc.invalidateMolecules(info);
cc.invalidateMolecules(info, false, true);
}
class CommonCalcCustomNonbondedForceKernel::ForceInfo : public ComputeForceInfo {
......@@ -2707,7 +2728,7 @@ double CommonCalcCustomNonbondedForceKernel::execute(ContextImpl& context, bool
return 0;
}
void CommonCalcCustomNonbondedForceKernel::copyParametersToContext(ContextImpl& context, const CustomNonbondedForce& force) {
void CommonCalcCustomNonbondedForceKernel::copyParametersToContext(ContextImpl& context, const CustomNonbondedForce& force, int firstParticle, int lastParticle) {
ContextSelector selector(cc);
int numParticles = force.getNumParticles();
if (numParticles != cc.getNumAtoms())
......@@ -2715,17 +2736,19 @@ void CommonCalcCustomNonbondedForceKernel::copyParametersToContext(ContextImpl&
// Record the per-particle parameters.
int paddedNumParticles = cc.getPaddedNumAtoms();
int numParams = force.getNumPerParticleParameters();
vector<vector<float> > paramVector(paddedNumParticles, vector<float>(numParams, 0));
vector<double> parameters;
for (int i = 0; i < numParticles; i++) {
force.getParticleParameters(i, parameters);
paramVector[i].resize(parameters.size());
for (int j = 0; j < (int) parameters.size(); j++)
paramVector[i][j] = (float) parameters[j];
if (firstParticle <= lastParticle) {
int numToSet = lastParticle-firstParticle+1;
int numParams = force.getNumPerParticleParameters();
vector<vector<float> > paramVector(numToSet, vector<float>(numParams, 0));
vector<double> parameters;
for (int i = 0; i < numToSet; i++) {
force.getParticleParameters(firstParticle+i, parameters);
paramVector[i].resize(parameters.size());
for (int j = 0; j < (int) parameters.size(); j++)
paramVector[i][j] = (float) parameters[j];
}
params->setParameterValuesSubset(firstParticle, paramVector);
}
params->setParameterValues(paramVector);
// If necessary, recompute the long range correction.
......@@ -2750,7 +2773,7 @@ void CommonCalcCustomNonbondedForceKernel::copyParametersToContext(ContextImpl&
// Mark that the current reordering may be invalid.
cc.invalidateMolecules(info);
cc.invalidateMolecules(info, firstParticle <= lastParticle, false);
}
class CommonCalcGBSAOBCForceKernel::ForceInfo : public ComputeForceInfo {
......@@ -2952,7 +2975,7 @@ void CommonCalcGBSAOBCForceKernel::copyParametersToContext(ContextImpl& context,
// Mark that the current reordering may be invalid.
cc.invalidateMolecules(info);
cc.invalidateMolecules(info, true, false);
}
class CommonCalcCustomGBForceKernel::ForceInfo : public ComputeForceInfo {
......@@ -3971,7 +3994,7 @@ void CommonCalcCustomGBForceKernel::copyParametersToContext(ContextImpl& context
// Mark that the current reordering may be invalid.
cc.invalidateMolecules(info);
cc.invalidateMolecules(info, true, false);
}
class CommonCalcCustomHbondForceKernel::ForceInfo : public ComputeForceInfo {
......@@ -4500,7 +4523,7 @@ void CommonCalcCustomHbondForceKernel::copyParametersToContext(ContextImpl& cont
// Mark that the current reordering may be invalid.
cc.invalidateMolecules(info);
cc.invalidateMolecules(info, false, true);
}
class CommonCalcCustomManyParticleForceKernel::ForceInfo : public ComputeForceInfo {
......
......@@ -338,7 +338,7 @@ void ComputeContext::invalidateMolecules() {
return;
}
bool ComputeContext::invalidateMolecules(ComputeForceInfo* force) {
bool ComputeContext::invalidateMolecules(ComputeForceInfo* force, bool checkAtoms, bool checkGroups) {
if (numAtoms == 0 || !getNonbondedUtilities().getUseCutoff())
return false;
bool valid = true;
......@@ -362,15 +362,17 @@ bool ComputeContext::invalidateMolecules(ComputeForceInfo* force) {
// See if the atoms are identical.
Molecule& m2 = molecules[instances[j]];
int offset2 = offsets[j];
for (int i = 0; i < (int) atoms.size() && valid; i++) {
if (!force->areParticlesIdentical(atoms[i]+offset1, atoms[i]+offset2))
valid = false;
if (checkAtoms) {
int offset2 = offsets[j];
for (int i = 0; i < (int) atoms.size() && valid; i++) {
if (!force->areParticlesIdentical(atoms[i]+offset1, atoms[i]+offset2))
valid = false;
}
}
// See if the force groups are identical.
if (valid && forceIndex > -1) {
if (valid && forceIndex > -1 && checkGroups) {
for (int k = 0; k < (int) m1.groups[forceIndex].size() && valid; k++)
if (!force->areGroupsIdentical(m1.groups[forceIndex][k], m2.groups[forceIndex][k]))
valid = false;
......
......@@ -117,6 +117,11 @@ void ComputeParameterSet::getParameterValues(vector<vector<T> >& values) {
template <class T>
void ComputeParameterSet::setParameterValues(const vector<vector<T> >& values, bool convert) {
setParameterValuesSubset(0, values, convert);
}
template <class T>
void ComputeParameterSet::setParameterValuesSubset(int first, const vector<vector<T> >& values, bool convert) {
if (convert && sizeof(T) == sizeof(float) && elementSize == sizeof(double)) {
vector<vector<double> > values2(values.size());
for (int i = 0; i < values.size(); i++) {
......@@ -124,7 +129,7 @@ void ComputeParameterSet::setParameterValues(const vector<vector<T> >& values, b
for (int j = 0; j < values[i].size(); j++)
values2[i][j] = (double) values[i][j];
}
setParameterValues(values2);
setParameterValuesSubset(first, values2);
return;
}
if (convert && sizeof(T) == sizeof(double) && elementSize == sizeof(float)) {
......@@ -134,16 +139,19 @@ void ComputeParameterSet::setParameterValues(const vector<vector<T> >& values, b
for (int j = 0; j < values[i].size(); j++)
values2[i][j] = (float) values[i][j];
}
setParameterValues(values2);
setParameterValuesSubset(first, values2);
return;
}
if (sizeof(T) != elementSize)
throw OpenMMException("Called setParameterValues() with vector of wrong type");
if (first < 0 || first+values.size() > numObjects)
throw OpenMMException("Called setParameterValuesSubset() for an illegal range of objects");
int numToSet = values.size();
int base = 0;
for (int i = 0; i < (int) arrays.size(); i++) {
if (arrays[i]->getElementSize() == 4*elementSize) {
vector<T> data(4*numObjects);
for (int j = 0; j < numObjects; j++) {
vector<T> data(4*numToSet);
for (int j = 0; j < numToSet; j++) {
data[4*j] = values[j][base];
if (base+1 < numParameters)
data[4*j+1] = values[j][base+1];
......@@ -152,24 +160,24 @@ void ComputeParameterSet::setParameterValues(const vector<vector<T> >& values, b
if (base+3 < numParameters)
data[4*j+3] = values[j][base+3];
}
arrays[i]->upload(data.data());
arrays[i]->uploadSubArray(data.data(), first, numToSet);
base += 4;
}
else if (arrays[i]->getElementSize() == 2*elementSize) {
vector<T> data(2*numObjects);
for (int j = 0; j < numObjects; j++) {
vector<T> data(2*numToSet);
for (int j = 0; j < numToSet; j++) {
data[2*j] = values[j][base];
if (base+1 < numParameters)
data[2*j+1] = values[j][base+1];
}
arrays[i]->upload(data.data());
arrays[i]->uploadSubArray(data.data(), first, numToSet);
base += 2;
}
else if (arrays[i]->getElementSize() == elementSize) {
vector<T> data(numObjects);
for (int j = 0; j < numObjects; j++)
vector<T> data(numToSet);
for (int j = 0; j < numToSet; j++)
data[j] = values[j][base];
arrays[i]->upload(data.data());
arrays[i]->uploadSubArray(data.data(), first, numToSet);
base++;
}
else
......@@ -201,6 +209,8 @@ string ComputeParameterSet::getParameterSuffix(int index, const std::string& ext
namespace OpenMM {
template void ComputeParameterSet::getParameterValues<float>(vector<vector<float> >& values);
template void ComputeParameterSet::setParameterValues<float>(const vector<vector<float> >& values, bool convert);
template void ComputeParameterSet::setParameterValuesSubset<float>(int first, const vector<vector<float> >& values, bool convert);
template void ComputeParameterSet::getParameterValues<double>(vector<vector<double> >& values);
template void ComputeParameterSet::setParameterValues<double>(const vector<vector<double> >& values, bool convert);
template void ComputeParameterSet::setParameterValuesSubset<double>(int first, const vector<vector<double> >& values, bool convert);
}
\ No newline at end of file
......@@ -122,10 +122,10 @@ public:
/**
* Copy changed parameters over to a context.
*
* @param context the context to copy parameters to
* @param force the HarmonicAngleForce to copy the parameters from
* @param firstAngle the index of the first bond whose parameters might have changed
* @param lastAngle the index of the last bond whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const HarmonicAngleForce& force);
void copyParametersToContext(ContextImpl& context, const HarmonicAngleForce& force, int firstAngle, int lastAngle);
private:
CpuPlatform::PlatformData& data;
int numAngles;
......@@ -162,10 +162,12 @@ public:
/**
* Copy changed parameters over to a context.
*
* @param context the context to copy parameters to
* @param force the PeriodicTorsionForce to copy the parameters from
* @param context the context to copy parameters to
* @param force the PeriodicTorsionForce to copy the parameters from
* @param firstTorsion the index of the first torsion whose parameters might have changed
* @param lastTorsion the index of the last torsion whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const PeriodicTorsionForce& force);
void copyParametersToContext(ContextImpl& context, const PeriodicTorsionForce& force, int firstTorsion, int lastTorsion);
private:
CpuPlatform::PlatformData& data;
int numTorsions;
......@@ -243,10 +245,14 @@ public:
/**
* Copy changed parameters over to a context.
*
* @param context the context to copy parameters to
* @param force the NonbondedForce to copy the parameters from
* @param context the context to copy parameters to
* @param force the NonbondedForce to copy the parameters from
* @param firstParticle the index of the first particle whose parameters might have changed
* @param lastParticle the index of the last particle whose parameters might have changed
* @param firstException the index of the first exception whose parameters might have changed
* @param lastException the index of the last exception whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const NonbondedForce& force);
void copyParametersToContext(ContextImpl& context, const NonbondedForce& force, int firstParticle, int lastParticle, int firstException, int lastException);
/**
* Get the parameters being used for PME.
*
......@@ -317,8 +323,10 @@ public:
*
* @param context the context to copy parameters to
* @param force the CustomNonbondedForce to copy the parameters from
* @param firstParticle the index of the first particle whose parameters might have changed
* @param lastParticle the index of the last particle whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const CustomNonbondedForce& force);
void copyParametersToContext(ContextImpl& context, const CustomNonbondedForce& force, int firstParticle, int lastParticle);
private:
void createInteraction(const CustomNonbondedForce& force);
CpuPlatform::PlatformData& data;
......
......@@ -337,13 +337,13 @@ double CpuCalcHarmonicAngleForceKernel::execute(ContextImpl& context, bool inclu
return energy;
}
void CpuCalcHarmonicAngleForceKernel::copyParametersToContext(ContextImpl& context, const HarmonicAngleForce& force) {
void CpuCalcHarmonicAngleForceKernel::copyParametersToContext(ContextImpl& context, const HarmonicAngleForce& force, int firstAngle, int lastAngle) {
if (numAngles != force.getNumAngles())
throw OpenMMException("updateParametersInContext: The number of angles has changed");
// Record the values.
for (int i = 0; i < numAngles; ++i) {
for (int i = firstAngle; i <= lastAngle; ++i) {
int particle1, particle2, particle3;
double angle, k;
force.getAngleParameters(i, particle1, particle2, particle3, angle, k);
......@@ -385,13 +385,13 @@ double CpuCalcPeriodicTorsionForceKernel::execute(ContextImpl& context, bool inc
return energy;
}
void CpuCalcPeriodicTorsionForceKernel::copyParametersToContext(ContextImpl& context, const PeriodicTorsionForce& force) {
void CpuCalcPeriodicTorsionForceKernel::copyParametersToContext(ContextImpl& context, const PeriodicTorsionForce& force, int firstTorsion, int lastTorsion) {
if (numTorsions != force.getNumTorsions())
throw OpenMMException("updateParametersInContext: The number of torsions has changed");
// Record the values.
for (int i = 0; i < numTorsions; ++i) {
for (int i = firstTorsion; i <= lastTorsion; ++i) {
int particle1, particle2, particle3, particle4, periodicity;
double phase, k;
force.getTorsionParameters(i, particle1, particle2, particle3, particle4, periodicity, phase, k);
......@@ -718,7 +718,7 @@ double CpuCalcNonbondedForceKernel::execute(ContextImpl& context, bool includeFo
return energy;
}
void CpuCalcNonbondedForceKernel::copyParametersToContext(ContextImpl& context, const NonbondedForce& force) {
void CpuCalcNonbondedForceKernel::copyParametersToContext(ContextImpl& context, const NonbondedForce& force, int firstParticle, int lastParticle, int firstException, int lastException) {
if (force.getNumParticles() != numParticles)
throw OpenMMException("updateParametersInContext: The number of particles has changed");
......@@ -745,7 +745,7 @@ void CpuCalcNonbondedForceKernel::copyParametersToContext(ContextImpl& context,
// Record the values.
for (int i = 0; i < numParticles; ++i)
for (int i = firstParticle; i <= lastParticle; ++i)
force.getParticleParameters(i, baseParticleParams[i][0], baseParticleParams[i][1], baseParticleParams[i][2]);
for (int i = 0; i < num14; ++i) {
int particle1, particle2;
......@@ -1028,16 +1028,15 @@ double CpuCalcCustomNonbondedForceKernel::execute(ContextImpl& context, bool inc
return energy;
}
void CpuCalcCustomNonbondedForceKernel::copyParametersToContext(ContextImpl& context, const CustomNonbondedForce& force) {
void CpuCalcCustomNonbondedForceKernel::copyParametersToContext(ContextImpl& context, const CustomNonbondedForce& force, int firstParticle, int lastParticle) {
if (numParticles != force.getNumParticles())
throw OpenMMException("updateParametersInContext: The number of particles has changed");
// Record the values.
int numParameters = force.getNumPerParticleParameters();
vector<double> params;
for (int i = 0; i < numParticles; ++i) {
vector<double> parameters;
vector<double> parameters;
for (int i = firstParticle; i <= lastParticle; ++i) {
force.getParticleParameters(i, parameters);
for (int j = 0; j < numParameters; j++)
particleParamArray[i][j] = parameters[j];
......
......@@ -9,7 +9,7 @@
* Biological Structures at Stanford, funded under the NIH Roadmap for *
* Medical Research, grant U54 GM072970. See https://simtk.org. *
* *
* Portions copyright (c) 2008-2022 Stanford University and the Authors. *
* Portions copyright (c) 2008-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -113,10 +113,14 @@ public:
/**
* Copy changed parameters over to a context.
*
* @param context the context to copy parameters to
* @param force the NonbondedForce to copy the parameters from
* @param context the context to copy parameters to
* @param force the NonbondedForce to copy the parameters from
* @param firstParticle the index of the first particle whose parameters might have changed
* @param lastParticle the index of the last particle whose parameters might have changed
* @param firstException the index of the first exception whose parameters might have changed
* @param lastException the index of the last exception whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const NonbondedForce& force);
void copyParametersToContext(ContextImpl& context, const NonbondedForce& force, int firstParticle, int lastParticle, int firstException, int lastException);
/**
* Get the parameters being used for PME.
*
......
......@@ -9,7 +9,7 @@
* Biological Structures at Stanford, funded under the NIH Roadmap for *
* Medical Research, grant U54 GM072970. See https://simtk.org. *
* *
* Portions copyright (c) 2011-2023 Stanford University and the Authors. *
* Portions copyright (c) 2011-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -126,8 +126,10 @@ public:
*
* @param context the context to copy parameters to
* @param force the HarmonicBondForce to copy the parameters from
* @param firstBond the index of the first bond whose parameters might have changed
* @param lastBond the index of the last bond whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const HarmonicBondForce& force);
void copyParametersToContext(ContextImpl& context, const HarmonicBondForce& force, int firstBond, int lastBond);
private:
class Task;
CudaPlatform::PlatformData& data;
......@@ -164,8 +166,10 @@ public:
*
* @param context the context to copy parameters to
* @param force the CustomBondForce to copy the parameters from
* @param firstBond the index of the first bond whose parameters might have changed
* @param lastBond the index of the last bond whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const CustomBondForce& force);
void copyParametersToContext(ContextImpl& context, const CustomBondForce& force, int firstBond, int lastBond);
private:
class Task;
CudaPlatform::PlatformData& data;
......@@ -202,8 +206,10 @@ public:
*
* @param context the context to copy parameters to
* @param force the HarmonicAngleForce to copy the parameters from
* @param firstAngle the index of the first bond whose parameters might have changed
* @param lastAngle the index of the last bond whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const HarmonicAngleForce& force);
void copyParametersToContext(ContextImpl& context, const HarmonicAngleForce& force, int firstAngle, int lastAngle);
private:
class Task;
CudaPlatform::PlatformData& data;
......@@ -240,8 +246,10 @@ public:
*
* @param context the context to copy parameters to
* @param force the CustomAngleForce to copy the parameters from
* @param firstAngle the index of the first bond whose parameters might have changed
* @param lastAngle the index of the last bond whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const CustomAngleForce& force);
void copyParametersToContext(ContextImpl& context, const CustomAngleForce& force, int firstAngle, int lastAngle);
private:
class Task;
CudaPlatform::PlatformData& data;
......@@ -277,10 +285,12 @@ public:
/**
* Copy changed parameters over to a context.
*
* @param context the context to copy parameters to
* @param force the PeriodicTorsionForce to copy the parameters from
* @param context the context to copy parameters to
* @param force the PeriodicTorsionForce to copy the parameters from
* @param firstTorsion the index of the first torsion whose parameters might have changed
* @param lastTorsion the index of the last torsion whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const PeriodicTorsionForce& force);
void copyParametersToContext(ContextImpl& context, const PeriodicTorsionForce& force, int firstTorsion, int lastTorsion);
private:
CudaPlatform::PlatformData& data;
std::vector<Kernel> kernels;
......@@ -390,10 +400,12 @@ public:
/**
* Copy changed parameters over to a context.
*
* @param context the context to copy parameters to
* @param force the CustomTorsionForce to copy the parameters from
* @param context the context to copy parameters to
* @param force the CustomTorsionForce to copy the parameters from
* @param firstTorsion the index of the first torsion whose parameters might have changed
* @param lastTorsion the index of the last torsion whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const CustomTorsionForce& force);
void copyParametersToContext(ContextImpl& context, const CustomTorsionForce& force, int firstTorsion, int lastTorsion);
private:
class Task;
CudaPlatform::PlatformData& data;
......@@ -430,10 +442,14 @@ public:
/**
* Copy changed parameters over to a context.
*
* @param context the context to copy parameters to
* @param force the NonbondedForce to copy the parameters from
* @param context the context to copy parameters to
* @param force the NonbondedForce to copy the parameters from
* @param firstParticle the index of the first particle whose parameters might have changed
* @param lastParticle the index of the last particle whose parameters might have changed
* @param firstException the index of the first exception whose parameters might have changed
* @param lastException the index of the last exception whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const NonbondedForce& force);
void copyParametersToContext(ContextImpl& context, const NonbondedForce& force, int firstParticle, int lastParticle, int firstException, int lastException);
/**
* Get the parameters being used for PME.
*
......@@ -488,8 +504,10 @@ public:
*
* @param context the context to copy parameters to
* @param force the CustomNonbondedForce to copy the parameters from
* @param firstParticle the index of the first particle whose parameters might have changed
* @param lastParticle the index of the last particle whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const CustomNonbondedForce& force);
void copyParametersToContext(ContextImpl& context, const CustomNonbondedForce& force, int firstParticle, int lastParticle);
private:
class Task;
CudaPlatform::PlatformData& data;
......@@ -524,10 +542,12 @@ public:
/**
* Copy changed parameters over to a context.
*
* @param context the context to copy parameters to
* @param force the CustomExternalForce to copy the parameters from
* @param context the context to copy parameters to
* @param force the CustomExternalForce to copy the parameters from
* @param firstParticle the index of the first particle whose parameters might have changed
* @param lastParticle the index of the last particle whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const CustomExternalForce& force);
void copyParametersToContext(ContextImpl& context, const CustomExternalForce& force, int firstParticle, int lastParticle);
private:
class Task;
CudaPlatform::PlatformData& data;
......
......@@ -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) 2008-2022 Stanford University and the Authors. *
* Portions copyright (c) 2008-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -1038,7 +1038,7 @@ double CudaCalcNonbondedForceKernel::execute(ContextImpl& context, bool includeF
return energy;
}
void CudaCalcNonbondedForceKernel::copyParametersToContext(ContextImpl& context, const NonbondedForce& force) {
void CudaCalcNonbondedForceKernel::copyParametersToContext(ContextImpl& context, const NonbondedForce& force, int firstParticle, int lastParticle, int firstException, int lastException) {
// Make sure the new parameters are acceptable.
ContextSelector selector(cu);
......@@ -1079,18 +1079,32 @@ void CudaCalcNonbondedForceKernel::copyParametersToContext(ContextImpl& context,
// Record the per-particle parameters.
vector<float4> baseParticleParamVec(cu.getPaddedNumAtoms(), make_float4(0, 0, 0, 0));
const vector<int>& order = cu.getAtomIndex();
for (int i = 0; i < force.getNumParticles(); i++) {
double charge, sigma, epsilon;
force.getParticleParameters(i, charge, sigma, epsilon);
baseParticleParamVec[i] = make_float4(charge, sigma, epsilon, 0);
if (firstParticle <= lastParticle) {
vector<float4> baseParticleParamVec(cu.getPaddedNumAtoms(), make_float4(0, 0, 0, 0));
for (int i = 0; i < force.getNumParticles(); i++) {
double charge, sigma, epsilon;
force.getParticleParameters(i, charge, sigma, epsilon);
baseParticleParamVec[i] = make_float4(charge, sigma, epsilon, 0);
}
baseParticleParams.uploadSubArray(&baseParticleParamVec[firstParticle], firstParticle, lastParticle-firstParticle+1);
// Compute the self energy.
ewaldSelfEnergy = 0.0;
if (nonbondedMethod == Ewald || nonbondedMethod == PME || nonbondedMethod == LJPME) {
if (cu.getContextIndex() == 0) {
for (int i = 0; i < force.getNumParticles(); i++) {
ewaldSelfEnergy -= baseParticleParamVec[i].x*baseParticleParamVec[i].x*ONE_4PI_EPS0*alpha/sqrt(M_PI);
if (doLJPME)
ewaldSelfEnergy += baseParticleParamVec[i].z*pow(baseParticleParamVec[i].y*dispersionAlpha, 6)/3.0;
}
}
}
}
baseParticleParams.upload(baseParticleParamVec);
// Record the exceptions.
if (numExceptions > 0) {
if (firstException <= lastException) {
vector<float4> baseExceptionParamsVec(numExceptions);
for (int i = 0; i < numExceptions; i++) {
int particle1, particle2;
......@@ -1105,19 +1119,9 @@ void CudaCalcNonbondedForceKernel::copyParametersToContext(ContextImpl& context,
// Compute other values.
ewaldSelfEnergy = 0.0;
if (nonbondedMethod == Ewald || nonbondedMethod == PME || nonbondedMethod == LJPME) {
if (cu.getContextIndex() == 0) {
for (int i = 0; i < force.getNumParticles(); i++) {
ewaldSelfEnergy -= baseParticleParamVec[i].x*baseParticleParamVec[i].x*ONE_4PI_EPS0*alpha/sqrt(M_PI);
if (doLJPME)
ewaldSelfEnergy += baseParticleParamVec[i].z*pow(baseParticleParamVec[i].y*dispersionAlpha, 6)/3.0;
}
}
}
if (force.getUseDispersionCorrection() && cu.getContextIndex() == 0 && (nonbondedMethod == CutoffPeriodic || nonbondedMethod == Ewald || nonbondedMethod == PME))
dispersionCoefficient = NonbondedForceImpl::calcDispersionCorrection(context.getSystem(), force);
cu.invalidateMolecules();
cu.invalidateMolecules(info, firstParticle <= lastParticle, firstException <= lastException);
recomputeParams = true;
}
......
......@@ -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) 2011-2023 Stanford University and the Authors. *
* Portions copyright (c) 2011-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -327,9 +327,9 @@ double CudaParallelCalcHarmonicBondForceKernel::execute(ContextImpl& context, bo
return 0.0;
}
void CudaParallelCalcHarmonicBondForceKernel::copyParametersToContext(ContextImpl& context, const HarmonicBondForce& force) {
void CudaParallelCalcHarmonicBondForceKernel::copyParametersToContext(ContextImpl& context, const HarmonicBondForce& force, int firstBond, int lastBond) {
for (int i = 0; i < (int) kernels.size(); i++)
getKernel(i).copyParametersToContext(context, force);
getKernel(i).copyParametersToContext(context, force, firstBond, lastBond);
}
class CudaParallelCalcCustomBondForceKernel::Task : public CudaContext::WorkTask {
......@@ -368,9 +368,9 @@ double CudaParallelCalcCustomBondForceKernel::execute(ContextImpl& context, bool
return 0.0;
}
void CudaParallelCalcCustomBondForceKernel::copyParametersToContext(ContextImpl& context, const CustomBondForce& force) {
void CudaParallelCalcCustomBondForceKernel::copyParametersToContext(ContextImpl& context, const CustomBondForce& force, int firstBond, int lastBond) {
for (int i = 0; i < (int) kernels.size(); i++)
getKernel(i).copyParametersToContext(context, force);
getKernel(i).copyParametersToContext(context, force, firstBond, lastBond);
}
class CudaParallelCalcHarmonicAngleForceKernel::Task : public CudaContext::WorkTask {
......@@ -409,9 +409,9 @@ double CudaParallelCalcHarmonicAngleForceKernel::execute(ContextImpl& context, b
return 0.0;
}
void CudaParallelCalcHarmonicAngleForceKernel::copyParametersToContext(ContextImpl& context, const HarmonicAngleForce& force) {
void CudaParallelCalcHarmonicAngleForceKernel::copyParametersToContext(ContextImpl& context, const HarmonicAngleForce& force, int firstAngle, int lastAngle) {
for (int i = 0; i < (int) kernels.size(); i++)
getKernel(i).copyParametersToContext(context, force);
getKernel(i).copyParametersToContext(context, force, firstAngle, lastAngle);
}
class CudaParallelCalcCustomAngleForceKernel::Task : public CudaContext::WorkTask {
......@@ -450,9 +450,9 @@ double CudaParallelCalcCustomAngleForceKernel::execute(ContextImpl& context, boo
return 0.0;
}
void CudaParallelCalcCustomAngleForceKernel::copyParametersToContext(ContextImpl& context, const CustomAngleForce& force) {
void CudaParallelCalcCustomAngleForceKernel::copyParametersToContext(ContextImpl& context, const CustomAngleForce& force, int firstAngle, int lastAngle) {
for (int i = 0; i < (int) kernels.size(); i++)
getKernel(i).copyParametersToContext(context, force);
getKernel(i).copyParametersToContext(context, force, firstAngle, lastAngle);
}
class CudaParallelCalcPeriodicTorsionForceKernel::Task : public CudaContext::WorkTask {
......@@ -491,9 +491,9 @@ double CudaParallelCalcPeriodicTorsionForceKernel::execute(ContextImpl& context,
return 0.0;
}
void CudaParallelCalcPeriodicTorsionForceKernel::copyParametersToContext(ContextImpl& context, const PeriodicTorsionForce& force) {
void CudaParallelCalcPeriodicTorsionForceKernel::copyParametersToContext(ContextImpl& context, const PeriodicTorsionForce& force, int firstTorsion, int lastTorsion) {
for (int i = 0; i < (int) kernels.size(); i++)
getKernel(i).copyParametersToContext(context, force);
getKernel(i).copyParametersToContext(context, force, firstTorsion, lastTorsion);
}
class CudaParallelCalcRBTorsionForceKernel::Task : public CudaContext::WorkTask {
......@@ -614,9 +614,9 @@ double CudaParallelCalcCustomTorsionForceKernel::execute(ContextImpl& context, b
return 0.0;
}
void CudaParallelCalcCustomTorsionForceKernel::copyParametersToContext(ContextImpl& context, const CustomTorsionForce& force) {
void CudaParallelCalcCustomTorsionForceKernel::copyParametersToContext(ContextImpl& context, const CustomTorsionForce& force, int firstTorsion, int lastTorsion) {
for (int i = 0; i < (int) kernels.size(); i++)
getKernel(i).copyParametersToContext(context, force);
getKernel(i).copyParametersToContext(context, force, firstTorsion, lastTorsion);
}
class CudaParallelCalcNonbondedForceKernel::Task : public CudaContext::WorkTask {
......@@ -655,9 +655,9 @@ double CudaParallelCalcNonbondedForceKernel::execute(ContextImpl& context, bool
return 0.0;
}
void CudaParallelCalcNonbondedForceKernel::copyParametersToContext(ContextImpl& context, const NonbondedForce& force) {
void CudaParallelCalcNonbondedForceKernel::copyParametersToContext(ContextImpl& context, const NonbondedForce& force, int firstParticle, int lastParticle, int firstException, int lastException) {
for (int i = 0; i < (int) kernels.size(); i++)
getKernel(i).copyParametersToContext(context, force);
getKernel(i).copyParametersToContext(context, force, firstParticle, lastParticle, firstException, lastException);
}
void CudaParallelCalcNonbondedForceKernel::getPMEParameters(double& alpha, int& nx, int& ny, int& nz) const {
......@@ -704,9 +704,9 @@ double CudaParallelCalcCustomNonbondedForceKernel::execute(ContextImpl& context,
return 0.0;
}
void CudaParallelCalcCustomNonbondedForceKernel::copyParametersToContext(ContextImpl& context, const CustomNonbondedForce& force) {
void CudaParallelCalcCustomNonbondedForceKernel::copyParametersToContext(ContextImpl& context, const CustomNonbondedForce& force, int firstParticle, int lastParticle) {
for (int i = 0; i < (int) kernels.size(); i++)
getKernel(i).copyParametersToContext(context, force);
getKernel(i).copyParametersToContext(context, force, firstParticle, lastParticle);
}
class CudaParallelCalcCustomExternalForceKernel::Task : public CudaContext::WorkTask {
......@@ -745,9 +745,9 @@ double CudaParallelCalcCustomExternalForceKernel::execute(ContextImpl& context,
return 0.0;
}
void CudaParallelCalcCustomExternalForceKernel::copyParametersToContext(ContextImpl& context, const CustomExternalForce& force) {
void CudaParallelCalcCustomExternalForceKernel::copyParametersToContext(ContextImpl& context, const CustomExternalForce& force, int firstParticle, int lastParticle) {
for (int i = 0; i < (int) kernels.size(); i++)
getKernel(i).copyParametersToContext(context, force);
getKernel(i).copyParametersToContext(context, force, firstParticle, lastParticle);
}
class CudaParallelCalcCustomHbondForceKernel::Task : public CudaContext::WorkTask {
......
......@@ -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) 2008-2015 Stanford University and the Authors. *
* Portions copyright (c) 2008-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -37,10 +37,11 @@ void testParallelComputation() {
const int numParticles = 200;
for (int i = 0; i < numParticles; i++)
system.addParticle(1.0);
CustomAngleForce* force = new CustomAngleForce("(theta-1.1)^2");
vector<double> params;
CustomAngleForce* force = new CustomAngleForce("k*(theta-theta0)^2");
force->addPerAngleParameter("k");
force->addPerAngleParameter("theta0");
for (int i = 2; i < numParticles; i++)
force->addAngle(i-2, i-1, i, params);
force->addAngle(i-2, i-1, i, {1.0, 1.1});
system.addForce(force);
vector<Vec3> positions(numParticles);
for (int i = 0; i < numParticles; i++)
......@@ -59,6 +60,21 @@ void testParallelComputation() {
ASSERT_EQUAL_TOL(state1.getPotentialEnergy(), state2.getPotentialEnergy(), 1e-5);
for (int i = 0; i < numParticles; i++)
ASSERT_EQUAL_VEC(state1.getForces()[i], state2.getForces()[i], 1e-5);
// Try updating some parameters and see if they still match.
vector<double> params;
for (int i = 95; i < 102; i++) {
int p1, p2, p3;
force->getAngleParameters(i, p1, p2, p3, params);
force->setAngleParameters(i, p1, p2, p3, {2.0, 1.2});
}
force->updateParametersInContext(context1);
force->updateParametersInContext(context2);
State state3 = context1.getState(State::Energy);
State state4 = context2.getState(State::Energy);
ASSERT_EQUAL_TOL(state3.getPotentialEnergy(), state4.getPotentialEnergy(), 1e-5);
ASSERT(fabs(state1.getPotentialEnergy()-state3.getPotentialEnergy()) > 0.1);
}
void runPlatformTests() {
......
......@@ -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) 2008-2015 Stanford University and the Authors. *
* Portions copyright (c) 2008-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -37,10 +37,11 @@ void testParallelComputation() {
const int numParticles = 200;
for (int i = 0; i < numParticles; i++)
system.addParticle(1.0);
CustomBondForce* force = new CustomBondForce(("(r-1.1)^2"));
vector<double> params;
CustomBondForce* force = new CustomBondForce(("k*(r-r0)^2"));
force->addPerBondParameter("k");
force->addPerBondParameter("r0");
for (int i = 1; i < numParticles; i++)
force->addBond(i-1, i, params);
force->addBond(i-1, i, {1.0, 1.1});
system.addForce(force);
vector<Vec3> positions(numParticles);
for (int i = 0; i < numParticles; i++)
......@@ -59,6 +60,21 @@ void testParallelComputation() {
ASSERT_EQUAL_TOL(state1.getPotentialEnergy(), state2.getPotentialEnergy(), 1e-5);
for (int i = 0; i < numParticles; i++)
ASSERT_EQUAL_VEC(state1.getForces()[i], state2.getForces()[i], 1e-5);
// Try updating some parameters and see if they still match.
vector<double> params;
for (int i = 95; i < 102; i++) {
int p1, p2;
force->getBondParameters(i, p1, p2, params);
force->setBondParameters(i, p1, p2, {2.0, 1.2});
}
force->updateParametersInContext(context1);
force->updateParametersInContext(context2);
State state3 = context1.getState(State::Energy);
State state4 = context2.getState(State::Energy);
ASSERT_EQUAL_TOL(state3.getPotentialEnergy(), state4.getPotentialEnergy(), 1e-5);
ASSERT(fabs(state1.getPotentialEnergy()-state3.getPotentialEnergy()) > 0.1);
}
void runPlatformTests() {
......
......@@ -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) 2008-2015 Stanford University and the Authors. *
* Portions copyright (c) 2008-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -37,10 +37,11 @@ void testParallelComputation() {
const int numParticles = 200;
for (int i = 0; i < numParticles; i++)
system.addParticle(1.0);
CustomNonbondedForce* force = new CustomNonbondedForce("4*eps*((sigma/r)^12-(sigma/r)^6); sigma=0.5; eps=1");
vector<double> params;
CustomNonbondedForce* force = new CustomNonbondedForce("4*eps*((sigma/r)^12-(sigma/r)^6); sigma=0.5*(sigma1+sigma2); eps=sqrt(eps1*eps2)");
force->addPerParticleParameter("sigma");
force->addPerParticleParameter("eps");
for (int i = 0; i < numParticles; i++)
force->addParticle(params);
force->addParticle({0.5, 1.0});
system.addForce(force);
OpenMM_SFMT::SFMT sfmt;
init_gen_rand(0, sfmt);
......@@ -67,6 +68,18 @@ void testParallelComputation() {
ASSERT_EQUAL_TOL(state1.getPotentialEnergy(), state2.getPotentialEnergy(), 1e-5);
for (int i = 0; i < numParticles; i++)
ASSERT_EQUAL_VEC(state1.getForces()[i], state2.getForces()[i], 1e-5);
// Try updating some parameters and see if they still match.
vector<double> params;
for (int i = 95; i < 102; i++)
force->setParticleParameters(i, {0.6, 2.0});
force->updateParametersInContext(context1);
force->updateParametersInContext(context2);
State state3 = context1.getState(State::Energy);
State state4 = context2.getState(State::Energy);
ASSERT_EQUAL_TOL(state3.getPotentialEnergy(), state4.getPotentialEnergy(), 1e-5);
ASSERT(fabs(state1.getPotentialEnergy()-state3.getPotentialEnergy()) > 0.1);
}
void runPlatformTests() {
......
......@@ -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) 2008-2015 Stanford University and the Authors. *
* Portions copyright (c) 2008-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -37,10 +37,11 @@ void testParallelComputation() {
const int numParticles = 200;
for (int i = 0; i < numParticles; i++)
system.addParticle(1.0);
CustomTorsionForce* force = new CustomTorsionForce("sin(theta-1.1)");
vector<double> params;
CustomTorsionForce* force = new CustomTorsionForce("k*sin(theta-theta0)");
force->addPerTorsionParameter("k");
force->addPerTorsionParameter("theta0");
for (int i = 3; i < numParticles; i++)
force->addTorsion(i-3, i-2, i-1, i, params);
force->addTorsion(i-3, i-2, i-1, i, {1.0, 1.1});
system.addForce(force);
vector<Vec3> positions(numParticles);
for (int i = 0; i < numParticles; i++)
......@@ -59,6 +60,21 @@ void testParallelComputation() {
ASSERT_EQUAL_TOL(state1.getPotentialEnergy(), state2.getPotentialEnergy(), 1e-5);
for (int i = 0; i < numParticles; i++)
ASSERT_EQUAL_VEC(state1.getForces()[i], state2.getForces()[i], 1e-5);
// Try updating some parameters and see if they still match.
vector<double> params;
for (int i = 95; i < 102; i++) {
int p1, p2, p3, p4;
force->getTorsionParameters(i, p1, p2, p3, p4, params);
force->setTorsionParameters(i, p1, p2, p3, p4, {2.0, 1.2});
}
force->updateParametersInContext(context1);
force->updateParametersInContext(context2);
State state3 = context1.getState(State::Energy);
State state4 = context2.getState(State::Energy);
ASSERT_EQUAL_TOL(state3.getPotentialEnergy(), state4.getPotentialEnergy(), 1e-5);
ASSERT(fabs(state1.getPotentialEnergy()-state3.getPotentialEnergy()) > 0.1);
}
void runPlatformTests() {
......
......@@ -7,7 +7,7 @@
* Biological Structures at Stanford, funded under the NIH Roadmap for *
* Medical Research, grant U54 GM072970. See https://simtk.org. *
* *
* Portions copyright (c) 2008-2015 Stanford University and the Authors. *
* Portions copyright (c) 2008-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -60,6 +60,21 @@ void testParallelComputation() {
ASSERT_EQUAL_TOL(state1.getPotentialEnergy(), state2.getPotentialEnergy(), 1e-5);
for (int i = 0; i < numParticles; i++)
ASSERT_EQUAL_VEC(state1.getForces()[i], state2.getForces()[i], 1e-5);
// Try updating some parameters and see if they still match.
for (int i = 95; i < 102; i++) {
int p1, p2, p3;
double angle, k;
force->getAngleParameters(i, p1, p2, p3, angle, k);
force->setAngleParameters(i, p1, p2, p3, angle+0.1, 2*k);
}
force->updateParametersInContext(context1);
force->updateParametersInContext(context2);
State state3 = context1.getState(State::Energy);
State state4 = context2.getState(State::Energy);
ASSERT_EQUAL_TOL(state3.getPotentialEnergy(), state4.getPotentialEnergy(), 1e-5);
ASSERT(fabs(state1.getPotentialEnergy()-state3.getPotentialEnergy()) > 0.1);
}
void runPlatformTests() {
......
......@@ -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) 2008-2015 Stanford University and the Authors. *
* Portions copyright (c) 2008-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -59,6 +59,21 @@ void testParallelComputation() {
ASSERT_EQUAL_TOL(state1.getPotentialEnergy(), state2.getPotentialEnergy(), 1e-5);
for (int i = 0; i < numParticles; i++)
ASSERT_EQUAL_VEC(state1.getForces()[i], state2.getForces()[i], 1e-5);
// Try updating some parameters and see if they still match.
for (int i = 95; i < 102; i++) {
int p1, p2;
double length, k;
force->getBondParameters(i, p1, p2, length, k);
force->setBondParameters(i, p1, p2, length+0.1, 2*k);
}
force->updateParametersInContext(context1);
force->updateParametersInContext(context2);
State state3 = context1.getState(State::Energy);
State state4 = context2.getState(State::Energy);
ASSERT_EQUAL_TOL(state3.getPotentialEnergy(), state4.getPotentialEnergy(), 1e-5);
ASSERT(fabs(state1.getPotentialEnergy()-state3.getPotentialEnergy()) > 0.1);
}
void runPlatformTests() {
......
......@@ -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) 2008-2021 Stanford University and the Authors. *
* Portions copyright (c) 2008-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -83,13 +83,20 @@ void testParallelComputation(NonbondedForce::NonbondedMethod method) {
for (int i = 0; i < numParticles; i++)
ASSERT_EQUAL_VEC(state1.getForces()[i], state2.getForces()[i], 1e-5);
// Modify some particle parameters and see if they still agree.
// Modify some parameters and see if they still agree.
for (int i = 0; i < numParticles; i += 5) {
double charge, sigma, epsilon;
force->getParticleParameters(i, charge, sigma, epsilon);
force->setParticleParameters(i, 0.9*charge, sigma, epsilon);
}
for (int i = force->getNumExceptions()/2-10; i < force->getNumExceptions()/2+10; i++) {
int p1, p2;
double charge, sigma, epsilon;
force->getExceptionParameters(i, p1, p2, charge, sigma, epsilon);
if (epsilon != 0)
force->setExceptionParameters(i, p1, p2, charge, sigma, 2*epsilon);
}
force->updateParametersInContext(context1);
force->updateParametersInContext(context2);
state1 = context1.getState(State::Forces | State::Energy);
......
......@@ -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) 2008-2015 Stanford University and the Authors. *
* Portions copyright (c) 2008-2024 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -58,6 +58,21 @@ void testParallelComputation() {
ASSERT_EQUAL_TOL(state1.getPotentialEnergy(), state2.getPotentialEnergy(), 1e-5);
for (int i = 0; i < numParticles; i++)
ASSERT_EQUAL_VEC(state1.getForces()[i], state2.getForces()[i], 1e-5);
// Try updating some parameters and see if they still match.
for (int i = 95; i < 102; i++) {
int p1, p2, p3, p4, periodicity;
double phase, k;
force->getTorsionParameters(i, p1, p2, p3, p4, periodicity, phase, k);
force->setTorsionParameters(i, p1, p2, p3, p4, periodicity, phase+0.1, 2*k);
}
force->updateParametersInContext(context1);
force->updateParametersInContext(context2);
State state3 = context1.getState(State::Energy);
State state4 = context2.getState(State::Energy);
ASSERT_EQUAL_TOL(state3.getPotentialEnergy(), state4.getPotentialEnergy(), 1e-5);
ASSERT(fabs(state1.getPotentialEnergy()-state3.getPotentialEnergy()) > 0.1);
}
void runPlatformTests() {
......
......@@ -113,10 +113,14 @@ public:
/**
* Copy changed parameters over to a context.
*
* @param context the context to copy parameters to
* @param force the NonbondedForce to copy the parameters from
* @param context the context to copy parameters to
* @param force the NonbondedForce to copy the parameters from
* @param firstParticle the index of the first particle whose parameters might have changed
* @param lastParticle the index of the last particle whose parameters might have changed
* @param firstException the index of the first exception whose parameters might have changed
* @param lastException the index of the last exception whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const NonbondedForce& force);
void copyParametersToContext(ContextImpl& context, const NonbondedForce& force, int firstParticle, int lastParticle, int firstException, int lastException);
/**
* Get the parameters being used for PME.
*
......
......@@ -9,7 +9,7 @@
* Biological Structures at Stanford, funded under the NIH Roadmap for *
* Medical Research, grant U54 GM072970. See https://simtk.org. *
* *
* Portions copyright (c) 2011-2019 Stanford University and the Authors. *
* Portions copyright (c) 2011-2024 Stanford University and the Authors. *
* Portions copyright (c) 2020-2023 Advanced Micro Devices, Inc. *
* Authors: Peter Eastman, Nicholas Curtis *
* Contributors: *
......@@ -125,8 +125,10 @@ public:
*
* @param context the context to copy parameters to
* @param force the HarmonicBondForce to copy the parameters from
* @param firstBond the index of the first bond whose parameters might have changed
* @param lastBond the index of the last bond whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const HarmonicBondForce& force);
void copyParametersToContext(ContextImpl& context, const HarmonicBondForce& force, int firstBond, int lastBond);
private:
class Task;
HipPlatform::PlatformData& data;
......@@ -163,8 +165,10 @@ public:
*
* @param context the context to copy parameters to
* @param force the CustomBondForce to copy the parameters from
* @param firstBond the index of the first bond whose parameters might have changed
* @param lastBond the index of the last bond whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const CustomBondForce& force);
void copyParametersToContext(ContextImpl& context, const CustomBondForce& force, int firstBond, int lastBond);
private:
class Task;
HipPlatform::PlatformData& data;
......@@ -201,8 +205,10 @@ public:
*
* @param context the context to copy parameters to
* @param force the HarmonicAngleForce to copy the parameters from
* @param firstAngle the index of the first bond whose parameters might have changed
* @param lastAngle the index of the last bond whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const HarmonicAngleForce& force);
void copyParametersToContext(ContextImpl& context, const HarmonicAngleForce& force, int firstAngle, int lastAngle);
private:
class Task;
HipPlatform::PlatformData& data;
......@@ -239,8 +245,10 @@ public:
*
* @param context the context to copy parameters to
* @param force the CustomAngleForce to copy the parameters from
* @param firstAngle the index of the first bond whose parameters might have changed
* @param lastAngle the index of the last bond whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const CustomAngleForce& force);
void copyParametersToContext(ContextImpl& context, const CustomAngleForce& force, int firstAngle, int lastAngle);
private:
class Task;
HipPlatform::PlatformData& data;
......@@ -278,8 +286,10 @@ public:
*
* @param context the context to copy parameters to
* @param force the PeriodicTorsionForce to copy the parameters from
* @param firstTorsion the index of the first torsion whose parameters might have changed
* @param lastTorsion the index of the last torsion whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const PeriodicTorsionForce& force);
void copyParametersToContext(ContextImpl& context, const PeriodicTorsionForce& force, int firstTorsion, int lastTorsion);
private:
HipPlatform::PlatformData& data;
std::vector<Kernel> kernels;
......@@ -391,8 +401,10 @@ public:
*
* @param context the context to copy parameters to
* @param force the CustomTorsionForce to copy the parameters from
* @param firstTorsion the index of the first torsion whose parameters might have changed
* @param lastTorsion the index of the last torsion whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const CustomTorsionForce& force);
void copyParametersToContext(ContextImpl& context, const CustomTorsionForce& force, int firstTorsion, int lastTorsion);
private:
class Task;
HipPlatform::PlatformData& data;
......@@ -431,8 +443,12 @@ public:
*
* @param context the context to copy parameters to
* @param force the NonbondedForce to copy the parameters from
* @param firstParticle the index of the first particle whose parameters might have changed
* @param lastParticle the index of the last particle whose parameters might have changed
* @param firstException the index of the first exception whose parameters might have changed
* @param lastException the index of the last exception whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const NonbondedForce& force);
void copyParametersToContext(ContextImpl& context, const NonbondedForce& force, int firstParticle, int lastParticle, int firstException, int lastException);
/**
* Get the parameters being used for PME.
*
......@@ -487,8 +503,10 @@ public:
*
* @param context the context to copy parameters to
* @param force the CustomNonbondedForce to copy the parameters from
* @param firstParticle the index of the first particle whose parameters might have changed
* @param lastParticle the index of the last particle whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const CustomNonbondedForce& force);
void copyParametersToContext(ContextImpl& context, const CustomNonbondedForce& force, int firstParticle, int lastParticle);
private:
class Task;
HipPlatform::PlatformData& data;
......@@ -525,8 +543,10 @@ public:
*
* @param context the context to copy parameters to
* @param force the CustomExternalForce to copy the parameters from
* @param firstParticle the index of the first particle whose parameters might have changed
* @param lastParticle the index of the last particle whose parameters might have changed
*/
void copyParametersToContext(ContextImpl& context, const CustomExternalForce& force);
void copyParametersToContext(ContextImpl& context, const CustomExternalForce& force, int firstParticle, int lastParticle);
private:
class Task;
HipPlatform::PlatformData& data;
......
......@@ -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) 2008-2022 Stanford University and the Authors. *
* Portions copyright (c) 2008-2024 Stanford University and the Authors. *
* Portions copyright (c) 2020-2022 Advanced Micro Devices, Inc. *
* Authors: Peter Eastman, Nicholas Curtis *
* Contributors: *
......@@ -946,7 +946,7 @@ double HipCalcNonbondedForceKernel::execute(ContextImpl& context, bool includeFo
return energy;
}
void HipCalcNonbondedForceKernel::copyParametersToContext(ContextImpl& context, const NonbondedForce& force) {
void HipCalcNonbondedForceKernel::copyParametersToContext(ContextImpl& context, const NonbondedForce& force, int firstParticle, int lastParticle, int firstException, int lastException) {
// Make sure the new parameters are acceptable.
ContextSelector selector(cu);
......@@ -987,18 +987,32 @@ void HipCalcNonbondedForceKernel::copyParametersToContext(ContextImpl& context,
// Record the per-particle parameters.
vector<float4> baseParticleParamVec(cu.getPaddedNumAtoms(), make_float4(0, 0, 0, 0));
const vector<int>& order = cu.getAtomIndex();
for (int i = 0; i < force.getNumParticles(); i++) {
double charge, sigma, epsilon;
force.getParticleParameters(i, charge, sigma, epsilon);
baseParticleParamVec[i] = make_float4(charge, sigma, epsilon, 0);
if (firstParticle <= lastParticle) {
vector<float4> baseParticleParamVec(cu.getPaddedNumAtoms(), make_float4(0, 0, 0, 0));
for (int i = 0; i < force.getNumParticles(); i++) {
double charge, sigma, epsilon;
force.getParticleParameters(i, charge, sigma, epsilon);
baseParticleParamVec[i] = make_float4(charge, sigma, epsilon, 0);
}
baseParticleParams.uploadSubArray(&baseParticleParamVec[firstParticle], firstParticle, lastParticle-firstParticle+1);
// Compute the self energy.
ewaldSelfEnergy = 0.0;
if (nonbondedMethod == Ewald || nonbondedMethod == PME || nonbondedMethod == LJPME) {
if (cu.getContextIndex() == 0) {
for (int i = 0; i < force.getNumParticles(); i++) {
ewaldSelfEnergy -= baseParticleParamVec[i].x*baseParticleParamVec[i].x*ONE_4PI_EPS0*alpha/sqrt(M_PI);
if (doLJPME)
ewaldSelfEnergy += baseParticleParamVec[i].z*pow(baseParticleParamVec[i].y*dispersionAlpha, 6)/3.0;
}
}
}
}
baseParticleParams.upload(baseParticleParamVec);
// Record the exceptions.
if (numExceptions > 0) {
if (firstException <= lastException) {
vector<float4> baseExceptionParamsVec(numExceptions);
for (int i = 0; i < numExceptions; i++) {
int particle1, particle2;
......@@ -1013,19 +1027,9 @@ void HipCalcNonbondedForceKernel::copyParametersToContext(ContextImpl& context,
// Compute other values.
ewaldSelfEnergy = 0.0;
if (nonbondedMethod == Ewald || nonbondedMethod == PME || nonbondedMethod == LJPME) {
if (cu.getContextIndex() == 0) {
for (int i = 0; i < force.getNumParticles(); i++) {
ewaldSelfEnergy -= baseParticleParamVec[i].x*baseParticleParamVec[i].x*ONE_4PI_EPS0*alpha/sqrt(M_PI);
if (doLJPME)
ewaldSelfEnergy += baseParticleParamVec[i].z*pow(baseParticleParamVec[i].y*dispersionAlpha, 6)/3.0;
}
}
}
if (force.getUseDispersionCorrection() && cu.getContextIndex() == 0 && (nonbondedMethod == CutoffPeriodic || nonbondedMethod == Ewald || nonbondedMethod == PME))
dispersionCoefficient = NonbondedForceImpl::calcDispersionCorrection(context.getSystem(), force);
cu.invalidateMolecules();
cu.invalidateMolecules(info, firstParticle <= lastParticle, firstException <= lastException);
recomputeParams = true;
}
......
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