Commit 9f37b18b authored by peastman's avatar peastman
Browse files

Code cleanup to reference and CPU platforms

parent f9106ddb
/* Portions copyright (c) 2009-2016 Stanford University and Simbios. /* Portions copyright (c) 2009-2018 Stanford University and Simbios.
* Contributors: Peter Eastman * Contributors: Peter Eastman
* *
* Permission is hereby granted, free of charge, to any person obtaining * Permission is hereby granted, free of charge, to any person obtaining
...@@ -150,7 +150,7 @@ ReferenceCustomGBIxn::~ReferenceCustomGBIxn() { ...@@ -150,7 +150,7 @@ ReferenceCustomGBIxn::~ReferenceCustomGBIxn() {
periodicBoxVectors[2] = vectors[2]; periodicBoxVectors[2] = vectors[2];
} }
void ReferenceCustomGBIxn::calculateIxn(int numberOfAtoms, vector<Vec3>& atomCoordinates, double** atomParameters, void ReferenceCustomGBIxn::calculateIxn(int numberOfAtoms, vector<Vec3>& atomCoordinates, vector<vector<double> >& atomParameters,
const vector<set<int> >& exclusions, map<string, double>& globalParameters, vector<Vec3>& forces, const vector<set<int> >& exclusions, map<string, double>& globalParameters, vector<Vec3>& forces,
double* totalEnergy, double* energyParamDerivs) { double* totalEnergy, double* energyParamDerivs) {
for (auto& param : globalParameters) for (auto& param : globalParameters)
...@@ -193,7 +193,7 @@ void ReferenceCustomGBIxn::calculateIxn(int numberOfAtoms, vector<Vec3>& atomCoo ...@@ -193,7 +193,7 @@ void ReferenceCustomGBIxn::calculateIxn(int numberOfAtoms, vector<Vec3>& atomCoo
calculateChainRuleForces(numberOfAtoms, atomCoordinates, atomParameters, exclusions, forces, energyParamDerivs); calculateChainRuleForces(numberOfAtoms, atomCoordinates, atomParameters, exclusions, forces, energyParamDerivs);
} }
void ReferenceCustomGBIxn::calculateSingleParticleValue(int index, int numAtoms, vector<Vec3>& atomCoordinates, double** atomParameters) { void ReferenceCustomGBIxn::calculateSingleParticleValue(int index, int numAtoms, vector<Vec3>& atomCoordinates, vector<vector<double> >& atomParameters) {
values[index].resize(numAtoms); values[index].resize(numAtoms);
for (int i = 0; i < numAtoms; i++) { for (int i = 0; i < numAtoms; i++) {
expressionSet.setVariable(xIndex, atomCoordinates[i][0]); expressionSet.setVariable(xIndex, atomCoordinates[i][0]);
...@@ -217,7 +217,7 @@ void ReferenceCustomGBIxn::calculateSingleParticleValue(int index, int numAtoms, ...@@ -217,7 +217,7 @@ void ReferenceCustomGBIxn::calculateSingleParticleValue(int index, int numAtoms,
} }
} }
void ReferenceCustomGBIxn::calculateParticlePairValue(int index, int numAtoms, vector<Vec3>& atomCoordinates, double** atomParameters, void ReferenceCustomGBIxn::calculateParticlePairValue(int index, int numAtoms, vector<Vec3>& atomCoordinates, vector<vector<double> >& atomParameters,
const vector<set<int> >& exclusions, bool useExclusions) { const vector<set<int> >& exclusions, bool useExclusions) {
values[index].resize(numAtoms); values[index].resize(numAtoms);
for (int i = 0; i < numAtoms; i++) for (int i = 0; i < numAtoms; i++)
...@@ -246,7 +246,7 @@ void ReferenceCustomGBIxn::calculateParticlePairValue(int index, int numAtoms, v ...@@ -246,7 +246,7 @@ void ReferenceCustomGBIxn::calculateParticlePairValue(int index, int numAtoms, v
} }
} }
void ReferenceCustomGBIxn::calculateOnePairValue(int index, int atom1, int atom2, vector<Vec3>& atomCoordinates, double** atomParameters) { void ReferenceCustomGBIxn::calculateOnePairValue(int index, int atom1, int atom2, vector<Vec3>& atomCoordinates, vector<vector<double> >& atomParameters) {
double deltaR[ReferenceForce::LastDeltaRIndex]; double deltaR[ReferenceForce::LastDeltaRIndex];
if (periodic) if (periodic)
ReferenceForce::getDeltaRPeriodic(atomCoordinates[atom2], atomCoordinates[atom1], periodicBoxVectors, deltaR); ReferenceForce::getDeltaRPeriodic(atomCoordinates[atom2], atomCoordinates[atom1], periodicBoxVectors, deltaR);
...@@ -273,7 +273,7 @@ void ReferenceCustomGBIxn::calculateOnePairValue(int index, int atom1, int atom2 ...@@ -273,7 +273,7 @@ void ReferenceCustomGBIxn::calculateOnePairValue(int index, int atom1, int atom2
} }
void ReferenceCustomGBIxn::calculateSingleParticleEnergyTerm(int index, int numAtoms, vector<Vec3>& atomCoordinates, void ReferenceCustomGBIxn::calculateSingleParticleEnergyTerm(int index, int numAtoms, vector<Vec3>& atomCoordinates,
double** atomParameters, vector<Vec3>& forces, double* totalEnergy, double* energyParamDerivs) { vector<vector<double> >& atomParameters, vector<Vec3>& forces, double* totalEnergy, double* energyParamDerivs) {
for (int i = 0; i < numAtoms; i++) { for (int i = 0; i < numAtoms; i++) {
expressionSet.setVariable(xIndex, atomCoordinates[i][0]); expressionSet.setVariable(xIndex, atomCoordinates[i][0]);
expressionSet.setVariable(yIndex, atomCoordinates[i][1]); expressionSet.setVariable(yIndex, atomCoordinates[i][1]);
...@@ -300,7 +300,7 @@ void ReferenceCustomGBIxn::calculateSingleParticleEnergyTerm(int index, int numA ...@@ -300,7 +300,7 @@ void ReferenceCustomGBIxn::calculateSingleParticleEnergyTerm(int index, int numA
} }
} }
void ReferenceCustomGBIxn::calculateParticlePairEnergyTerm(int index, int numAtoms, vector<Vec3>& atomCoordinates, double** atomParameters, void ReferenceCustomGBIxn::calculateParticlePairEnergyTerm(int index, int numAtoms, vector<Vec3>& atomCoordinates, vector<vector<double> >& atomParameters,
const vector<set<int> >& exclusions, bool useExclusions, vector<Vec3>& forces, double* totalEnergy, double* energyParamDerivs) { const vector<set<int> >& exclusions, bool useExclusions, vector<Vec3>& forces, double* totalEnergy, double* energyParamDerivs) {
if (cutoff) { if (cutoff) {
// Loop over all pairs in the neighbor list. // Loop over all pairs in the neighbor list.
...@@ -324,7 +324,7 @@ void ReferenceCustomGBIxn::calculateParticlePairEnergyTerm(int index, int numAto ...@@ -324,7 +324,7 @@ void ReferenceCustomGBIxn::calculateParticlePairEnergyTerm(int index, int numAto
} }
} }
void ReferenceCustomGBIxn::calculateOnePairEnergyTerm(int index, int atom1, int atom2, vector<Vec3>& atomCoordinates, double** atomParameters, void ReferenceCustomGBIxn::calculateOnePairEnergyTerm(int index, int atom1, int atom2, vector<Vec3>& atomCoordinates, vector<vector<double> >& atomParameters,
vector<Vec3>& forces, double* totalEnergy, double* energyParamDerivs) { vector<Vec3>& forces, double* totalEnergy, double* energyParamDerivs) {
// Compute the displacement. // Compute the displacement.
...@@ -370,7 +370,7 @@ void ReferenceCustomGBIxn::calculateOnePairEnergyTerm(int index, int atom1, int ...@@ -370,7 +370,7 @@ void ReferenceCustomGBIxn::calculateOnePairEnergyTerm(int index, int atom1, int
energyParamDerivs[i] += energyParamDerivExpressions[index][i].evaluate(); energyParamDerivs[i] += energyParamDerivExpressions[index][i].evaluate();
} }
void ReferenceCustomGBIxn::calculateChainRuleForces(int numAtoms, vector<Vec3>& atomCoordinates, double** atomParameters, void ReferenceCustomGBIxn::calculateChainRuleForces(int numAtoms, vector<Vec3>& atomCoordinates, vector<vector<double> >& atomParameters,
const vector<set<int> >& exclusions, vector<Vec3>& forces, double* energyParamDerivs) { const vector<set<int> >& exclusions, vector<Vec3>& forces, double* energyParamDerivs) {
if (cutoff) { if (cutoff) {
// Loop over all pairs in the neighbor list. // Loop over all pairs in the neighbor list.
...@@ -429,7 +429,7 @@ void ReferenceCustomGBIxn::calculateChainRuleForces(int numAtoms, vector<Vec3>& ...@@ -429,7 +429,7 @@ void ReferenceCustomGBIxn::calculateChainRuleForces(int numAtoms, vector<Vec3>&
energyParamDerivs[k] += dEdV[j][i]*dValuedParam[j][k][i]; energyParamDerivs[k] += dEdV[j][i]*dValuedParam[j][k][i];
} }
void ReferenceCustomGBIxn::calculateOnePairChainRule(int atom1, int atom2, vector<Vec3>& atomCoordinates, double** atomParameters, void ReferenceCustomGBIxn::calculateOnePairChainRule(int atom1, int atom2, vector<Vec3>& atomCoordinates, vector<vector<double> >& atomParameters,
vector<Vec3>& forces, bool isExcluded) { vector<Vec3>& forces, bool isExcluded) {
// Compute the displacement. // Compute the displacement.
......
...@@ -116,7 +116,7 @@ void ReferenceCustomHbondIxn::setPeriodic(Vec3* vectors) { ...@@ -116,7 +116,7 @@ void ReferenceCustomHbondIxn::setPeriodic(Vec3* vectors) {
--------------------------------------------------------------------------------------- */ --------------------------------------------------------------------------------------- */
void ReferenceCustomHbondIxn::calculatePairIxn(vector<Vec3>& atomCoordinates, double** donorParameters, double** acceptorParameters, void ReferenceCustomHbondIxn::calculatePairIxn(vector<Vec3>& atomCoordinates, vector<vector<double> >& donorParameters, vector<vector<double> >& acceptorParameters,
vector<set<int> >& exclusions, const map<string, double>& globalParameters, vector<Vec3>& forces, vector<set<int> >& exclusions, const map<string, double>& globalParameters, vector<Vec3>& forces,
double* totalEnergy) const { double* totalEnergy) const {
......
/* Portions copyright (c) 2009-2014 Stanford University and Simbios. /* Portions copyright (c) 2009-2018 Stanford University and Simbios.
* Contributors: Peter Eastman * Contributors: Peter Eastman
* *
* Permission is hereby granted, free of charge, to any person obtaining * Permission is hereby granted, free of charge, to any person obtaining
...@@ -105,7 +105,7 @@ ReferenceCustomManyParticleIxn::ReferenceCustomManyParticleIxn(const CustomManyP ...@@ -105,7 +105,7 @@ ReferenceCustomManyParticleIxn::ReferenceCustomManyParticleIxn(const CustomManyP
ReferenceCustomManyParticleIxn::~ReferenceCustomManyParticleIxn() { ReferenceCustomManyParticleIxn::~ReferenceCustomManyParticleIxn() {
} }
void ReferenceCustomManyParticleIxn::calculateIxn(vector<Vec3>& atomCoordinates, double** particleParameters, void ReferenceCustomManyParticleIxn::calculateIxn(vector<Vec3>& atomCoordinates, vector<vector<double> >& particleParameters,
const map<string, double>& globalParameters, vector<Vec3>& forces, const map<string, double>& globalParameters, vector<Vec3>& forces,
double* totalEnergy) const { double* totalEnergy) const {
map<string, double> variables = globalParameters; map<string, double> variables = globalParameters;
...@@ -130,7 +130,7 @@ void ReferenceCustomManyParticleIxn::setPeriodic(Vec3* vectors) { ...@@ -130,7 +130,7 @@ void ReferenceCustomManyParticleIxn::setPeriodic(Vec3* vectors) {
} }
void ReferenceCustomManyParticleIxn::loopOverInteractions(vector<int>& particles, int loopIndex, vector<OpenMM::Vec3>& atomCoordinates, void ReferenceCustomManyParticleIxn::loopOverInteractions(vector<int>& particles, int loopIndex, vector<OpenMM::Vec3>& atomCoordinates,
double** particleParameters, map<string, double>& variables, vector<OpenMM::Vec3>& forces, vector<vector<double> >& particleParameters, map<string, double>& variables, vector<OpenMM::Vec3>& forces,
double* totalEnergy) const { double* totalEnergy) const {
int numParticles = atomCoordinates.size(); int numParticles = atomCoordinates.size();
int firstPartialLoop = (centralParticleMode ? 2 : 1); int firstPartialLoop = (centralParticleMode ? 2 : 1);
...@@ -147,7 +147,7 @@ void ReferenceCustomManyParticleIxn::loopOverInteractions(vector<int>& particles ...@@ -147,7 +147,7 @@ void ReferenceCustomManyParticleIxn::loopOverInteractions(vector<int>& particles
} }
void ReferenceCustomManyParticleIxn::calculateOneIxn(const vector<int>& particles, vector<Vec3>& atomCoordinates, void ReferenceCustomManyParticleIxn::calculateOneIxn(const vector<int>& particles, vector<Vec3>& atomCoordinates,
double** particleParameters, map<string, double>& variables, vector<Vec3>& forces, double* totalEnergy) const { vector<vector<double> >& particleParameters, map<string, double>& variables, vector<Vec3>& forces, double* totalEnergy) const {
// Select the ordering to use for the particles. // Select the ordering to use for the particles.
vector<int> permutedParticles(numParticlesPerSet); vector<int> permutedParticles(numParticlesPerSet);
......
/* Portions copyright (c) 2009-2016 Stanford University and Simbios. /* Portions copyright (c) 2009-2018 Stanford University and Simbios.
* Contributors: Peter Eastman * Contributors: Peter Eastman
* *
* Permission is hereby granted, free of charge, to any person obtaining * Permission is hereby granted, free of charge, to any person obtaining
...@@ -146,18 +146,16 @@ void ReferenceCustomNonbondedIxn::setUseSwitchingFunction(double distance) { ...@@ -146,18 +146,16 @@ void ReferenceCustomNonbondedIxn::setUseSwitchingFunction(double distance) {
@param atomParameters atom parameters atomParameters[atomIndex][paramterIndex] @param atomParameters atom parameters atomParameters[atomIndex][paramterIndex]
@param exclusions atom exclusion indices @param exclusions atom exclusion indices
exclusions[atomIndex] contains the list of exclusions for that atom exclusions[atomIndex] contains the list of exclusions for that atom
@param fixedParameters non atom parameters (not currently used)
@param globalParameters the values of global parameters @param globalParameters the values of global parameters
@param forces force array (forces added) @param forces force array (forces added)
@param energyByAtom atom energy
@param totalEnergy total energy @param totalEnergy total energy
--------------------------------------------------------------------------------------- */ --------------------------------------------------------------------------------------- */
void ReferenceCustomNonbondedIxn::calculatePairIxn(int numberOfAtoms, vector<Vec3>& atomCoordinates, void ReferenceCustomNonbondedIxn::calculatePairIxn(int numberOfAtoms, vector<Vec3>& atomCoordinates,
double** atomParameters, vector<set<int> >& exclusions, vector<vector<double> >& atomParameters, vector<set<int> >& exclusions,
double* fixedParameters, const map<string, double>& globalParameters, vector<Vec3>& forces, const map<string, double>& globalParameters, vector<Vec3>& forces,
double* energyByAtom, double* totalEnergy, double* energyParamDerivs) { double* totalEnergy, double* energyParamDerivs) {
for (auto& param : globalParameters) for (auto& param : globalParameters)
expressionSet.setVariable(expressionSet.getVariableIndex(param.first), param.second); expressionSet.setVariable(expressionSet.getVariableIndex(param.first), param.second);
...@@ -177,7 +175,7 @@ void ReferenceCustomNonbondedIxn::calculatePairIxn(int numberOfAtoms, vector<Vec ...@@ -177,7 +175,7 @@ void ReferenceCustomNonbondedIxn::calculatePairIxn(int numberOfAtoms, vector<Vec
expressionSet.setVariable(particleParamIndex[j*2], atomParameters[*atom1][j]); expressionSet.setVariable(particleParamIndex[j*2], atomParameters[*atom1][j]);
expressionSet.setVariable(particleParamIndex[j*2+1], atomParameters[*atom2][j]); expressionSet.setVariable(particleParamIndex[j*2+1], atomParameters[*atom2][j]);
} }
calculateOneIxn(*atom1, *atom2, atomCoordinates, forces, energyByAtom, totalEnergy, energyParamDerivs); calculateOneIxn(*atom1, *atom2, atomCoordinates, forces, totalEnergy, energyParamDerivs);
} }
} }
} }
...@@ -190,7 +188,7 @@ void ReferenceCustomNonbondedIxn::calculatePairIxn(int numberOfAtoms, vector<Vec ...@@ -190,7 +188,7 @@ void ReferenceCustomNonbondedIxn::calculatePairIxn(int numberOfAtoms, vector<Vec
expressionSet.setVariable(particleParamIndex[j*2], atomParameters[pair.first][j]); expressionSet.setVariable(particleParamIndex[j*2], atomParameters[pair.first][j]);
expressionSet.setVariable(particleParamIndex[j*2+1], atomParameters[pair.second][j]); expressionSet.setVariable(particleParamIndex[j*2+1], atomParameters[pair.second][j]);
} }
calculateOneIxn(pair.first, pair.second, atomCoordinates, forces, energyByAtom, totalEnergy, energyParamDerivs); calculateOneIxn(pair.first, pair.second, atomCoordinates, forces, totalEnergy, energyParamDerivs);
} }
} }
else { else {
...@@ -203,7 +201,7 @@ void ReferenceCustomNonbondedIxn::calculatePairIxn(int numberOfAtoms, vector<Vec ...@@ -203,7 +201,7 @@ void ReferenceCustomNonbondedIxn::calculatePairIxn(int numberOfAtoms, vector<Vec
expressionSet.setVariable(particleParamIndex[j*2], atomParameters[ii][j]); expressionSet.setVariable(particleParamIndex[j*2], atomParameters[ii][j]);
expressionSet.setVariable(particleParamIndex[j*2+1], atomParameters[jj][j]); expressionSet.setVariable(particleParamIndex[j*2+1], atomParameters[jj][j]);
} }
calculateOneIxn(ii, jj, atomCoordinates, forces, energyByAtom, totalEnergy, energyParamDerivs); calculateOneIxn(ii, jj, atomCoordinates, forces, totalEnergy, energyParamDerivs);
} }
} }
} }
...@@ -217,15 +215,13 @@ void ReferenceCustomNonbondedIxn::calculatePairIxn(int numberOfAtoms, vector<Vec ...@@ -217,15 +215,13 @@ void ReferenceCustomNonbondedIxn::calculatePairIxn(int numberOfAtoms, vector<Vec
@param ii the index of the first atom @param ii the index of the first atom
@param jj the index of the second atom @param jj the index of the second atom
@param atomCoordinates atom coordinates @param atomCoordinates atom coordinates
@param atomParameters atom parameters (charges, c6, c12, ...) atomParameters[atomIndex][paramterIndex]
@param forces force array (forces added) @param forces force array (forces added)
@param energyByAtom atom energy
@param totalEnergy total energy @param totalEnergy total energy
--------------------------------------------------------------------------------------- */ --------------------------------------------------------------------------------------- */
void ReferenceCustomNonbondedIxn::calculateOneIxn(int ii, int jj, vector<Vec3>& atomCoordinates, vector<Vec3>& forces, void ReferenceCustomNonbondedIxn::calculateOneIxn(int ii, int jj, vector<Vec3>& atomCoordinates, vector<Vec3>& forces,
double* energyByAtom, double* totalEnergy, double* energyParamDerivs) { double* totalEnergy, double* energyParamDerivs) {
// get deltaR, R2, and R between 2 atoms // get deltaR, R2, and R between 2 atoms
double deltaR[ReferenceForce::LastDeltaRIndex]; double deltaR[ReferenceForce::LastDeltaRIndex];
...@@ -262,14 +258,8 @@ void ReferenceCustomNonbondedIxn::calculateOneIxn(int ii, int jj, vector<Vec3>& ...@@ -262,14 +258,8 @@ void ReferenceCustomNonbondedIxn::calculateOneIxn(int ii, int jj, vector<Vec3>&
// accumulate energies // accumulate energies
if (totalEnergy || energyByAtom) {
if (totalEnergy) if (totalEnergy)
*totalEnergy += energy; *totalEnergy += energy;
if (energyByAtom) {
energyByAtom[ii] += energy;
energyByAtom[jj] += energy;
}
}
} }
...@@ -81,9 +81,9 @@ void ReferenceCustomTorsionIxn::setPeriodic(OpenMM::Vec3* vectors) { ...@@ -81,9 +81,9 @@ void ReferenceCustomTorsionIxn::setPeriodic(OpenMM::Vec3* vectors) {
--------------------------------------------------------------------------------------- */ --------------------------------------------------------------------------------------- */
void ReferenceCustomTorsionIxn::calculateBondIxn(int* atomIndices, void ReferenceCustomTorsionIxn::calculateBondIxn(vector<int>& atomIndices,
vector<Vec3>& atomCoordinates, vector<Vec3>& atomCoordinates,
double* parameters, vector<double>& parameters,
vector<Vec3>& forces, vector<Vec3>& forces,
double* totalEnergy, double* energyParamDerivs) { double* totalEnergy, double* energyParamDerivs) {
double deltaR[3][ReferenceForce::LastDeltaRIndex]; double deltaR[3][ReferenceForce::LastDeltaRIndex];
......
...@@ -70,9 +70,9 @@ void ReferenceHarmonicBondIxn::setPeriodic(OpenMM::Vec3* vectors) { ...@@ -70,9 +70,9 @@ void ReferenceHarmonicBondIxn::setPeriodic(OpenMM::Vec3* vectors) {
--------------------------------------------------------------------------------------- */ --------------------------------------------------------------------------------------- */
void ReferenceHarmonicBondIxn::calculateBondIxn(int* atomIndices, void ReferenceHarmonicBondIxn::calculateBondIxn(vector<int>& atomIndices,
vector<Vec3>& atomCoordinates, vector<Vec3>& atomCoordinates,
double* parameters, vector<double>& parameters,
vector<Vec3>& forces, vector<Vec3>& forces,
double* totalEnergy, double* energyParamDerivs) { double* totalEnergy, double* energyParamDerivs) {
double deltaR[ReferenceForce::LastDeltaRIndex]; double deltaR[ReferenceForce::LastDeltaRIndex];
......
...@@ -65,8 +65,8 @@ ReferenceLJCoulomb14::~ReferenceLJCoulomb14() { ...@@ -65,8 +65,8 @@ ReferenceLJCoulomb14::~ReferenceLJCoulomb14() {
--------------------------------------------------------------------------------------- */ --------------------------------------------------------------------------------------- */
void ReferenceLJCoulomb14::calculateBondIxn(int* atomIndices, vector<Vec3>& atomCoordinates, void ReferenceLJCoulomb14::calculateBondIxn(vector<int>& atomIndices, vector<Vec3>& atomCoordinates,
double* parameters, vector<Vec3>& forces, vector<double>& parameters, vector<Vec3>& forces,
double* totalEnergy, double* energyParamDerivs) { double* totalEnergy, double* energyParamDerivs) {
double deltaR[2][ReferenceForce::LastDeltaRIndex]; double deltaR[2][ReferenceForce::LastDeltaRIndex];
......
/* Portions copyright (c) 2006-2013 Stanford University and Simbios. /* Portions copyright (c) 2006-2018 Stanford University and Simbios.
* Contributors: Pande Group * Contributors: Pande Group
* *
* Permission is hereby granted, free of charge, to any person obtaining * Permission is hereby granted, free of charge, to any person obtaining
...@@ -176,9 +176,7 @@ void ReferenceLJCoulombIxn::setUseLJPME(double alpha, int meshSize[3]) { ...@@ -176,9 +176,7 @@ void ReferenceLJCoulombIxn::setUseLJPME(double alpha, int meshSize[3]) {
@param atomParameters atom parameters atomParameters[atomIndex][paramterIndex] @param atomParameters atom parameters atomParameters[atomIndex][paramterIndex]
@param exclusions atom exclusion indices @param exclusions atom exclusion indices
exclusions[atomIndex] contains the list of exclusions for that atom exclusions[atomIndex] contains the list of exclusions for that atom
@param fixedParameters non atom parameters (not currently used)
@param forces force array (forces added) @param forces force array (forces added)
@param energyByAtom atom energy
@param totalEnergy total energy @param totalEnergy total energy
@param includeDirect true if direct space interactions should be included @param includeDirect true if direct space interactions should be included
@param includeReciprocal true if reciprocal space interactions should be included @param includeReciprocal true if reciprocal space interactions should be included
...@@ -186,9 +184,8 @@ void ReferenceLJCoulombIxn::setUseLJPME(double alpha, int meshSize[3]) { ...@@ -186,9 +184,8 @@ void ReferenceLJCoulombIxn::setUseLJPME(double alpha, int meshSize[3]) {
--------------------------------------------------------------------------------------- */ --------------------------------------------------------------------------------------- */
void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& atomCoordinates, void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& atomCoordinates,
double** atomParameters, vector<set<int> >& exclusions, vector<vector<double> >& atomParameters, vector<set<int> >& exclusions,
double* fixedParameters, vector<Vec3>& forces, vector<Vec3>& forces, double* totalEnergy, bool includeDirect, bool includeReciprocal) const {
double* energyByAtom, double* totalEnergy, bool includeDirect, bool includeReciprocal) const {
typedef std::complex<double> d_complex; typedef std::complex<double> d_complex;
static const double epsilon = 1.0; static const double epsilon = 1.0;
...@@ -224,9 +221,6 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a ...@@ -224,9 +221,6 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a
selfEwaldEnergy -= pow(alphaDispersionEwald, 6.0) * 64.0*pow(atomParameters[atomID][SigIndex], 6.0) * pow(atomParameters[atomID][EpsIndex], 2.0) / 12.0; selfEwaldEnergy -= pow(alphaDispersionEwald, 6.0) * 64.0*pow(atomParameters[atomID][SigIndex], 6.0) * pow(atomParameters[atomID][EpsIndex], 2.0) / 12.0;
} }
totalSelfEwaldEnergy -= selfEwaldEnergy; totalSelfEwaldEnergy -= selfEwaldEnergy;
if (energyByAtom) {
energyByAtom[atomID] -= selfEwaldEnergy;
}
} }
} }
...@@ -252,10 +246,6 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a ...@@ -252,10 +246,6 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a
if (totalEnergy) if (totalEnergy)
*totalEnergy += recipEnergy; *totalEnergy += recipEnergy;
if (energyByAtom)
for (int n = 0; n < numberOfAtoms; n++)
energyByAtom[n] += recipEnergy;
pme_destroy(pmedata); pme_destroy(pmedata);
if (ljpme) { if (ljpme) {
...@@ -275,10 +265,6 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a ...@@ -275,10 +265,6 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a
} }
if (totalEnergy) if (totalEnergy)
*totalEnergy += recipDispersionEnergy; *totalEnergy += recipDispersionEnergy;
if (energyByAtom)
for (int n = 0; n < numberOfAtoms; n++)
energyByAtom[n] += recipDispersionEnergy;
pme_destroy(pmedata); pme_destroy(pmedata);
} }
} }
...@@ -374,10 +360,6 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a ...@@ -374,10 +360,6 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a
if (totalEnergy) if (totalEnergy)
*totalEnergy += recipEnergy; *totalEnergy += recipEnergy;
if (energyByAtom)
for (int n = 0; n < numberOfAtoms; n++)
energyByAtom[n] += recipEnergy;
lowrz = 1 - numRz; lowrz = 1 - numRz;
} }
lowry = 1 - numRy; lowry = 1 - numRy;
...@@ -473,11 +455,6 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a ...@@ -473,11 +455,6 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a
totalVdwEnergy += vdwEnergy; totalVdwEnergy += vdwEnergy;
totalRealSpaceEwaldEnergy += realSpaceEwaldEnergy; totalRealSpaceEwaldEnergy += realSpaceEwaldEnergy;
if (energyByAtom) {
energyByAtom[ii] += realSpaceEwaldEnergy + vdwEnergy;
energyByAtom[jj] += realSpaceEwaldEnergy + vdwEnergy;
}
} }
if (totalEnergy) if (totalEnergy)
...@@ -537,10 +514,6 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a ...@@ -537,10 +514,6 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a
} }
totalExclusionEnergy += realSpaceEwaldEnergy; totalExclusionEnergy += realSpaceEwaldEnergy;
if (energyByAtom) {
energyByAtom[ii] -= realSpaceEwaldEnergy;
energyByAtom[jj] -= realSpaceEwaldEnergy;
}
} }
} }
...@@ -558,9 +531,7 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a ...@@ -558,9 +531,7 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a
@param atomParameters atom parameters atomParameters[atomIndex][paramterIndex] @param atomParameters atom parameters atomParameters[atomIndex][paramterIndex]
@param exclusions atom exclusion indices @param exclusions atom exclusion indices
exclusions[atomIndex] contains the list of exclusions for that atom exclusions[atomIndex] contains the list of exclusions for that atom
@param fixedParameters non atom parameters (not currently used)
@param forces force array (forces added) @param forces force array (forces added)
@param energyByAtom atom energy
@param totalEnergy total energy @param totalEnergy total energy
@param includeDirect true if direct space interactions should be included @param includeDirect true if direct space interactions should be included
@param includeReciprocal true if reciprocal space interactions should be included @param includeReciprocal true if reciprocal space interactions should be included
...@@ -568,12 +539,11 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a ...@@ -568,12 +539,11 @@ void ReferenceLJCoulombIxn::calculateEwaldIxn(int numberOfAtoms, vector<Vec3>& a
--------------------------------------------------------------------------------------- */ --------------------------------------------------------------------------------------- */
void ReferenceLJCoulombIxn::calculatePairIxn(int numberOfAtoms, vector<Vec3>& atomCoordinates, void ReferenceLJCoulombIxn::calculatePairIxn(int numberOfAtoms, vector<Vec3>& atomCoordinates,
double** atomParameters, vector<set<int> >& exclusions, vector<vector<double> >& atomParameters, vector<set<int> >& exclusions,
double* fixedParameters, vector<Vec3>& forces, vector<Vec3>& forces, double* totalEnergy, bool includeDirect, bool includeReciprocal) const {
double* energyByAtom, double* totalEnergy, bool includeDirect, bool includeReciprocal) const {
if (ewald || pme || ljpme) { if (ewald || pme || ljpme) {
calculateEwaldIxn(numberOfAtoms, atomCoordinates, atomParameters, exclusions, fixedParameters, forces, energyByAtom, calculateEwaldIxn(numberOfAtoms, atomCoordinates, atomParameters, exclusions, forces,
totalEnergy, includeDirect, includeReciprocal); totalEnergy, includeDirect, includeReciprocal);
return; return;
} }
...@@ -581,7 +551,7 @@ void ReferenceLJCoulombIxn::calculatePairIxn(int numberOfAtoms, vector<Vec3>& at ...@@ -581,7 +551,7 @@ void ReferenceLJCoulombIxn::calculatePairIxn(int numberOfAtoms, vector<Vec3>& at
return; return;
if (cutoff) { if (cutoff) {
for (auto& pair : *neighborList) for (auto& pair : *neighborList)
calculateOneIxn(pair.first, pair.second, atomCoordinates, atomParameters, forces, energyByAtom, totalEnergy); calculateOneIxn(pair.first, pair.second, atomCoordinates, atomParameters, forces, totalEnergy);
} }
else { else {
for (int ii = 0; ii < numberOfAtoms; ii++) { for (int ii = 0; ii < numberOfAtoms; ii++) {
...@@ -589,7 +559,7 @@ void ReferenceLJCoulombIxn::calculatePairIxn(int numberOfAtoms, vector<Vec3>& at ...@@ -589,7 +559,7 @@ void ReferenceLJCoulombIxn::calculatePairIxn(int numberOfAtoms, vector<Vec3>& at
for (int jj = ii+1; jj < numberOfAtoms; jj++) for (int jj = ii+1; jj < numberOfAtoms; jj++)
if (exclusions[jj].find(ii) == exclusions[jj].end()) if (exclusions[jj].find(ii) == exclusions[jj].end())
calculateOneIxn(ii, jj, atomCoordinates, atomParameters, forces, energyByAtom, totalEnergy); calculateOneIxn(ii, jj, atomCoordinates, atomParameters, forces, totalEnergy);
} }
} }
} }
...@@ -603,14 +573,13 @@ void ReferenceLJCoulombIxn::calculatePairIxn(int numberOfAtoms, vector<Vec3>& at ...@@ -603,14 +573,13 @@ void ReferenceLJCoulombIxn::calculatePairIxn(int numberOfAtoms, vector<Vec3>& at
@param atomCoordinates atom coordinates @param atomCoordinates atom coordinates
@param atomParameters atom parameters (charges, c6, c12, ...) atomParameters[atomIndex][paramterIndex] @param atomParameters atom parameters (charges, c6, c12, ...) atomParameters[atomIndex][paramterIndex]
@param forces force array (forces added) @param forces force array (forces added)
@param energyByAtom atom energy
@param totalEnergy total energy @param totalEnergy total energy
--------------------------------------------------------------------------------------- */ --------------------------------------------------------------------------------------- */
void ReferenceLJCoulombIxn::calculateOneIxn(int ii, int jj, vector<Vec3>& atomCoordinates, void ReferenceLJCoulombIxn::calculateOneIxn(int ii, int jj, vector<Vec3>& atomCoordinates,
double** atomParameters, vector<Vec3>& forces, vector<vector<double> >& atomParameters, vector<Vec3>& forces,
double* energyByAtom, double* totalEnergy) const { double* totalEnergy) const {
double deltaR[2][ReferenceForce::LastDeltaRIndex]; double deltaR[2][ReferenceForce::LastDeltaRIndex];
// get deltaR, R2, and R between 2 atoms // get deltaR, R2, and R between 2 atoms
...@@ -665,9 +634,5 @@ void ReferenceLJCoulombIxn::calculateOneIxn(int ii, int jj, vector<Vec3>& atomCo ...@@ -665,9 +634,5 @@ void ReferenceLJCoulombIxn::calculateOneIxn(int ii, int jj, vector<Vec3>& atomCo
if (totalEnergy) if (totalEnergy)
*totalEnergy += energy; *totalEnergy += energy;
if (energyByAtom) {
energyByAtom[ii] += energy;
energyByAtom[jj] += energy;
}
} }
...@@ -71,9 +71,9 @@ void ReferenceProperDihedralBond::setPeriodic(OpenMM::Vec3* vectors) { ...@@ -71,9 +71,9 @@ void ReferenceProperDihedralBond::setPeriodic(OpenMM::Vec3* vectors) {
--------------------------------------------------------------------------------------- */ --------------------------------------------------------------------------------------- */
void ReferenceProperDihedralBond::calculateBondIxn(int* atomIndices, void ReferenceProperDihedralBond::calculateBondIxn(vector<int>& atomIndices,
vector<Vec3>& atomCoordinates, vector<Vec3>& atomCoordinates,
double* parameters, vector<double>& parameters,
vector<Vec3>& forces, vector<Vec3>& forces,
double* totalEnergy, double* energyParamDerivs) { double* totalEnergy, double* energyParamDerivs) {
double deltaR[3][ReferenceForce::LastDeltaRIndex]; double deltaR[3][ReferenceForce::LastDeltaRIndex];
......
...@@ -69,9 +69,9 @@ void ReferenceRbDihedralBond::setPeriodic(OpenMM::Vec3* vectors) { ...@@ -69,9 +69,9 @@ void ReferenceRbDihedralBond::setPeriodic(OpenMM::Vec3* vectors) {
--------------------------------------------------------------------------------------- */ --------------------------------------------------------------------------------------- */
void ReferenceRbDihedralBond::calculateBondIxn(int* atomIndices, void ReferenceRbDihedralBond::calculateBondIxn(vector<int>& atomIndices,
vector<Vec3>& atomCoordinates, vector<Vec3>& atomCoordinates,
double* parameters, vector<double>& parameters,
vector<Vec3>& forces, vector<Vec3>& forces,
double* totalEnergy, double* energyParamDerivs) { double* totalEnergy, double* energyParamDerivs) {
// number of parameters // number of parameters
......
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