Commit b2e9b5ac authored by Peter Eastman's avatar Peter Eastman
Browse files

Fixed bug in CustomGBForce

parent 7cdd6d16
...@@ -373,7 +373,7 @@ private: ...@@ -373,7 +373,7 @@ private:
class CpuCalcCustomGBForceKernel : public CalcCustomGBForceKernel { class CpuCalcCustomGBForceKernel : public CalcCustomGBForceKernel {
public: public:
CpuCalcCustomGBForceKernel(std::string name, const Platform& platform, CpuPlatform::PlatformData& data) : CpuCalcCustomGBForceKernel(std::string name, const Platform& platform, CpuPlatform::PlatformData& data) :
CalcCustomGBForceKernel(name, platform), data(data) { CalcCustomGBForceKernel(name, platform), data(data), ixn(NULL), neighborList(NULL) {
} }
~CpuCalcCustomGBForceKernel(); ~CpuCalcCustomGBForceKernel();
/** /**
...@@ -406,6 +406,7 @@ private: ...@@ -406,6 +406,7 @@ private:
double **particleParamArray; double **particleParamArray;
double nonbondedCutoff; double nonbondedCutoff;
CpuCustomGBForce* ixn; CpuCustomGBForce* ixn;
CpuNeighborList* neighborList;
std::vector<std::set<int> > exclusions; std::vector<std::set<int> > exclusions;
std::vector<std::string> particleParameterNames, globalParameterNames, energyParamDerivNames, valueNames; std::vector<std::string> particleParameterNames, globalParameterNames, energyParamDerivNames, valueNames;
std::vector<OpenMM::CustomGBForce::ComputationType> valueTypes; std::vector<OpenMM::CustomGBForce::ComputationType> valueTypes;
......
...@@ -1018,6 +1018,8 @@ CpuCalcCustomGBForceKernel::~CpuCalcCustomGBForceKernel() { ...@@ -1018,6 +1018,8 @@ CpuCalcCustomGBForceKernel::~CpuCalcCustomGBForceKernel() {
} }
if (ixn != NULL) if (ixn != NULL)
delete ixn; delete ixn;
if (neighborList != NULL)
delete neighborList;
} }
void CpuCalcCustomGBForceKernel::initialize(const System& system, const CustomGBForce& force) { void CpuCalcCustomGBForceKernel::initialize(const System& system, const CustomGBForce& force) {
...@@ -1064,7 +1066,7 @@ void CpuCalcCustomGBForceKernel::initialize(const System& system, const CustomGB ...@@ -1064,7 +1066,7 @@ void CpuCalcCustomGBForceKernel::initialize(const System& system, const CustomGB
nonbondedMethod = CalcCustomGBForceKernel::NonbondedMethod(force.getNonbondedMethod()); nonbondedMethod = CalcCustomGBForceKernel::NonbondedMethod(force.getNonbondedMethod());
nonbondedCutoff = force.getCutoffDistance(); nonbondedCutoff = force.getCutoffDistance();
if (nonbondedMethod != NoCutoff) if (nonbondedMethod != NoCutoff)
data.requestNeighborList(nonbondedCutoff, 0.25*nonbondedCutoff, force.getNumExclusions() > 0, exclusions); neighborList = new CpuNeighborList(4);
// Create custom functions for the tabulated functions. // Create custom functions for the tabulated functions.
...@@ -1171,7 +1173,8 @@ double CpuCalcCustomGBForceKernel::execute(ContextImpl& context, bool includeFor ...@@ -1171,7 +1173,8 @@ double CpuCalcCustomGBForceKernel::execute(ContextImpl& context, bool includeFor
ixn->setPeriodic(extractBoxSize(context)); ixn->setPeriodic(extractBoxSize(context));
if (nonbondedMethod != NoCutoff) { if (nonbondedMethod != NoCutoff) {
vector<set<int> > noExclusions(numParticles); vector<set<int> > noExclusions(numParticles);
ixn->setUseCutoff(nonbondedCutoff, *data.neighborList); neighborList->computeNeighborList(numParticles, data.posq, noExclusions, boxVectors, data.isPeriodic, nonbondedCutoff, data.threads);
ixn->setUseCutoff(nonbondedCutoff, *neighborList);
} }
map<string, double> globalParameters; map<string, double> globalParameters;
for (auto& name : globalParameterNames) for (auto& name : globalParameterNames)
......
...@@ -1474,7 +1474,8 @@ double ReferenceCalcCustomGBForceKernel::execute(ContextImpl& context, bool incl ...@@ -1474,7 +1474,8 @@ double ReferenceCalcCustomGBForceKernel::execute(ContextImpl& context, bool incl
if (periodic) if (periodic)
ixn.setPeriodic(extractBoxVectors(context)); ixn.setPeriodic(extractBoxVectors(context));
if (nonbondedMethod != NoCutoff) { if (nonbondedMethod != NoCutoff) {
computeNeighborListVoxelHash(*neighborList, numParticles, posData, exclusions, extractBoxVectors(context), periodic, nonbondedCutoff, 0.0); vector<set<int> > empty(context.getSystem().getNumParticles()); // Don't omit exclusions from the neighbor list
computeNeighborListVoxelHash(*neighborList, numParticles, posData, empty, extractBoxVectors(context), periodic, nonbondedCutoff, 0.0);
ixn.setUseCutoff(nonbondedCutoff, *neighborList); ixn.setUseCutoff(nonbondedCutoff, *neighborList);
} }
map<string, double> globalParameters; map<string, double> globalParameters;
......
...@@ -118,6 +118,7 @@ void testOBC(GBSAOBCForce::NonbondedMethod obcMethod, CustomGBForce::NonbondedMe ...@@ -118,6 +118,7 @@ void testOBC(GBSAOBCForce::NonbondedMethod obcMethod, CustomGBForce::NonbondedMe
params[1] = 0.1; params[1] = 0.1;
custom->addParticle(params); custom->addParticle(params);
} }
custom->addExclusion(2*i, 2*i+1);
positions[2*i] = Vec3(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt)); positions[2*i] = Vec3(boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt), boxSize*genrand_real2(sfmt));
positions[2*i+1] = Vec3(positions[2*i][0]+1.0, positions[2*i][1], positions[2*i][2]); positions[2*i+1] = Vec3(positions[2*i][0]+1.0, positions[2*i][1], positions[2*i][2]);
velocities[2*i] = Vec3(genrand_real2(sfmt), genrand_real2(sfmt), genrand_real2(sfmt)); velocities[2*i] = Vec3(genrand_real2(sfmt), genrand_real2(sfmt), genrand_real2(sfmt));
......
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