Skip to content

GitLab

"plugins/vscode:/vscode.git/clone" did not exist on "a99a2d0e57a4930686ba8c737969c00bf90183d6"
Commit 6ed5bc4e authored by Rafal Wiewiora's avatar Rafal Wiewiora Committed by GitHub
Browse files

Merge branch 'master' into master

parents 656d0e3b fbf193fe
Expand all Show whitespace changes
Inline Side-by-side
Showing with 446 additions and 477 deletions
platforms/cuda/src/kernels/customHbondForce.cu
View file @ 6ed5bc4e
......@@ -105,8 +105,8 @@ extern "C" __global__ void computeDonorForces(unsigned long long* __restrict__ f
__syncthreads();
if (donorIndex < NUM_DONORS) {
for (int index = 0; index < blockSize; index++) {
#ifdef USE_EXCLUSIONS
int acceptorIndex = acceptorStart+index;
#ifdef USE_EXCLUSIONS
if (acceptorIndex == exclusionIndices.x || acceptorIndex == exclusionIndices.y || acceptorIndex == exclusionIndices.z || acceptorIndex == exclusionIndices.w)
continue;
#endif
......@@ -193,8 +193,8 @@ extern "C" __global__ void computeAcceptorForces(unsigned long long* __restrict_
__syncthreads();
if (acceptorIndex < NUM_ACCEPTORS) {
for (int index = 0; index < blockSize; index++) {
#ifdef USE_EXCLUSIONS
int donorIndex = donorStart+index;
#ifdef USE_EXCLUSIONS
if (donorIndex == exclusionIndices.x || donorIndex == exclusionIndices.y || donorIndex == exclusionIndices.z || donorIndex == exclusionIndices.w)
continue;
#endif
......
platforms/opencl/src/OpenCLBondedUtilities.cpp
View file @ 6ed5bc4e
......@@ -184,8 +184,7 @@ void OpenCLBondedUtilities::initialize(const System& system) {
// Create the kernels.
for (vector<vector<int> >::const_iterator iter = forceSets.begin(); iter != forceSets.end(); ++iter) {
const vector<int>& set = *iter;
for (auto& set : forceSets) {
int setSize = set.size();
stringstream s;
s<<"#ifdef SUPPORTS_64_BIT_ATOMICS\n";
......
platforms/opencl/src/OpenCLContext.cpp
View file @ 6ed5bc4e
......@@ -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) 2009-2016 Stanford University and the Authors. *
* Portions copyright (c) 2009-2017 Stanford University and the Authors. *
* Authors: Peter Eastman *
* Contributors: *
* *
......@@ -318,8 +318,8 @@ OpenCLContext::OpenCLContext(const System& system, int platformIndex, int device
OpenCLArray valuesArray(*this, 20, sizeof(mm_float8), "values");
vector<mm_float8> values(valuesArray.getSize());
float nextValue = 1e-4f;
for (int i = 0; i < (int) values.size(); ++i) {
values[i].s0 = nextValue;
for (auto& val : values) {
val.s0 = nextValue;
nextValue *= (float) M_PI;
}
valuesArray.upload(values);
......@@ -328,14 +328,14 @@ OpenCLContext::OpenCLContext(const System& system, int platformIndex, int device
executeKernel(accuracyKernel, values.size());
valuesArray.download(values);
double maxSqrtError = 0.0, maxRsqrtError = 0.0, maxRecipError = 0.0, maxExpError = 0.0, maxLogError = 0.0;
for (int i = 0; i < (int) values.size(); ++i) {
double v = values[i].s0;
for (auto& val : values) {
double v = val.s0;
double correctSqrt = sqrt(v);
maxSqrtError = max(maxSqrtError, fabs(correctSqrt-values[i].s1)/correctSqrt);
maxRsqrtError = max(maxRsqrtError, fabs(1.0/correctSqrt-values[i].s2)*correctSqrt);
maxRecipError = max(maxRecipError, fabs(1.0/v-values[i].s3)/values[i].s3);
maxExpError = max(maxExpError, fabs(exp(v)-values[i].s4)/values[i].s4);
maxLogError = max(maxLogError, fabs(log(v)-values[i].s5)/values[i].s5);
maxSqrtError = max(maxSqrtError, fabs(correctSqrt-val.s1)/correctSqrt);
maxRsqrtError = max(maxRsqrtError, fabs(1.0/correctSqrt-val.s2)*correctSqrt);
maxRecipError = max(maxRecipError, fabs(1.0/v-val.s3)/val.s3);
maxExpError = max(maxExpError, fabs(exp(v)-val.s4)/val.s4);
maxLogError = max(maxLogError, fabs(log(v)-val.s5)/val.s5);
}
compilationDefines["SQRT"] = (maxSqrtError < 1e-6) ? "native_sqrt" : "sqrt";
compilationDefines["RSQRT"] = (maxRsqrtError < 1e-6) ? "native_rsqrt" : "rsqrt";
......@@ -412,14 +412,14 @@ OpenCLContext::OpenCLContext(const System& system, int platformIndex, int device
}
OpenCLContext::~OpenCLContext() {
for (int i = 0; i < (int) forces.size(); i++)
delete forces[i];
for (int i = 0; i < (int) reorderListeners.size(); i++)
delete reorderListeners[i];
for (int i = 0; i < (int) preComputations.size(); i++)
delete preComputations[i];
for (int i = 0; i < (int) postComputations.size(); i++)
delete postComputations[i];
for (auto force : forces)
delete force;
for (auto listener : reorderListeners)
delete listener;
for (auto computation : preComputations)
delete computation;
for (auto computation : postComputations)
delete computation;
if (pinnedBuffer != NULL)
delete pinnedBuffer;
if (posq != NULL)
......@@ -458,8 +458,8 @@ void OpenCLContext::initialize() {
bonded->initialize(system);
numForceBuffers = platformData.contexts.size();
numForceBuffers = std::max(numForceBuffers, bonded->getNumForceBuffers());
for (int i = 0; i < (int) forces.size(); i++)
numForceBuffers = std::max(numForceBuffers, forces[i]->getRequiredForceBuffers());
for (auto force : forces)
numForceBuffers = std::max(numForceBuffers, force->getRequiredForceBuffers());
int energyBufferSize = max(numThreadBlocks*ThreadBlockSize, nonbonded->getNumEnergyBuffers());
if (useDoublePrecision) {
forceBuffers = OpenCLArray::create<mm_double4>(*this, paddedNumAtoms*numForceBuffers, "forceBuffers");
......@@ -525,17 +525,17 @@ string OpenCLContext::replaceStrings(const string& input, const std::map<std::st
symbolChars.insert(c);
}
string result = input;
for (map<string, string>::const_iterator iter = replacements.begin(); iter != replacements.end(); iter++) {
for (auto& pair : replacements) {
int index = 0;
int size = iter->first.size();
int size = pair.first.size();
do {
index = result.find(iter->first, index);
index = result.find(pair.first, index);
if (index != result.npos) {
if ((index == 0 || symbolChars.find(result[index-1]) == symbolChars.end()) && (index == result.size()-size || symbolChars.find(result[index+size]) == symbolChars.end())) {
// We have found a complete symbol, not part of a longer symbol.
result.replace(index, size, iter->second);
index += iter->second.size();
result.replace(index, size, pair.second);
index += pair.second.size();
}
else
index++;
......@@ -554,10 +554,10 @@ cl::Program OpenCLContext::createProgram(const string source, const map<string,
stringstream src;
if (!options.empty())
src << "// Compilation Options: " << options << endl << endl;
for (map<string, string>::const_iterator iter = compilationDefines.begin(); iter != compilationDefines.end(); ++iter) {
src << "#define " << iter->first;
if (!iter->second.empty())
src << " " << iter->second;
for (auto& pair : compilationDefines) {
src << "#define " << pair.first;
if (!pair.second.empty())
src << " " << pair.second;
src << endl;
}
if (!compilationDefines.empty())
......@@ -588,10 +588,10 @@ cl::Program OpenCLContext::createProgram(const string source, const map<string,
src << "typedef float3 mixed3;\n";
src << "typedef float4 mixed4;\n";
}
for (map<string, string>::const_iterator iter = defines.begin(); iter != defines.end(); ++iter) {
src << "#define " << iter->first;
if (!iter->second.empty())
src << " " << iter->second;
for (auto& pair : defines) {
src << "#define " << pair.first;
if (!pair.second.empty())
src << " " << pair.second;
src << endl;
}
if (!defines.empty())
......@@ -856,10 +856,10 @@ void OpenCLContext::findMoleculeGroups() {
atomBonds[particle1].push_back(particle2);
atomBonds[particle2].push_back(particle1);
}
for (int i = 0; i < (int) forces.size(); i++) {
for (int j = 0; j < forces[i]->getNumParticleGroups(); j++) {
for (auto force : forces) {
for (int j = 0; j < force->getNumParticleGroups(); j++) {
vector<int> particles;
forces[i]->getParticlesInGroup(j, particles);
force->getParticlesInGroup(j, particles);
for (int k = 0; k < (int) particles.size(); k++)
for (int m = 0; m < (int) particles.size(); m++)
if (k != m)
......@@ -1076,8 +1076,8 @@ bool OpenCLContext::invalidateMolecules(OpenCLForceInfo* force) {
}
atomIndexDevice->upload(atomIndex);
findMoleculeGroups();
for (int i = 0; i < (int) reorderListeners.size(); i++)
reorderListeners[i]->execute();
for (auto listener : reorderListeners)
listener->execute();
reorderAtoms();
return true;
}
......@@ -1138,10 +1138,9 @@ void OpenCLContext::reorderAtomsImpl() {
vector<Real4> newPosqCorrection(paddedNumAtoms, Real4(0,0,0,0));
vector<Mixed4> newVelm(paddedNumAtoms, Mixed4(0,0,0,0));
vector<mm_int4> newCellOffsets(numAtoms);
for (int group = 0; group < (int) moleculeGroups.size(); group++) {
for (auto& mol : moleculeGroups) {
// Find the center of each molecule.
MoleculeGroup& mol = moleculeGroups[group];
int numMolecules = mol.offsets.size();
vector<int>& atoms = mol.atoms;
vector<Real4> molPos(numMolecules);
......@@ -1235,9 +1234,9 @@ void OpenCLContext::reorderAtomsImpl() {
// Reorder the atoms.
for (int i = 0; i < numMolecules; i++) {
for (int j = 0; j < (int)atoms.size(); j++) {
int oldIndex = mol.offsets[molBins[i].second]+atoms[j];
int newIndex = mol.offsets[i]+atoms[j];
for (int atom : atoms) {
int oldIndex = mol.offsets[molBins[i].second]+atom;
int newIndex = mol.offsets[i]+atom;
originalIndex[newIndex] = atomIndex[oldIndex];
newPosq[newIndex] = oldPosq[oldIndex];
if (useMixedPrecision)
......@@ -1259,8 +1258,8 @@ void OpenCLContext::reorderAtomsImpl() {
posqCorrection->upload(newPosqCorrection);
velm->upload(newVelm);
atomIndexDevice->upload(atomIndex);
for (int i = 0; i < (int) reorderListeners.size(); i++)
reorderListeners[i]->execute();
for (auto listener : reorderListeners)
listener->execute();
}
void OpenCLContext::addReorderListener(ReorderListener* listener) {
......
platforms/opencl/src/OpenCLExpressionUtilities.cpp
View file @ 6ed5bc4e
......@@ -448,12 +448,12 @@ void OpenCLExpressionUtilities::processExpression(stringstream& out, const Expre
vector<bool> hasAssigned(powers.size(), false);
exponents.push_back((int) fabs(exponent));
names.push_back(name);
for (map<int, const ExpressionTreeNode*>::const_iterator iter = powers.begin(); iter != powers.end(); ++iter) {
if (iter->first != exponent) {
exponents.push_back(iter->first >= 0 ? iter->first : -iter->first);
for (auto& power : powers) {
if (power.first != exponent) {
exponents.push_back(power.first >= 0 ? power.first : -power.first);
string name2 = prefix+context.intToString(temps.size());
names.push_back(name2);
temps.push_back(make_pair(*iter->second, name2));
temps.push_back(make_pair(*power.second, name2));
out << tempType << " " << name2 << " = 0.0f;\n";
}
}
......
platforms/opencl/src/OpenCLKernels.cpp
View file @ 6ed5bc4e
This diff is collapsed.
platforms/opencl/src/kernels/customHbondForce.cl
View file @ 6ed5bc4e
......@@ -91,8 +91,8 @@ __kernel void computeDonorForces(__global real4* restrict forceBuffers, __global
barrier(CLK_LOCAL_MEM_FENCE);
if (donorIndex < NUM_DONORS) {
for (int index = 0; index < blockSize; index++) {
#ifdef USE_EXCLUSIONS
int acceptorIndex = acceptorStart+index;
#ifdef USE_EXCLUSIONS
if (acceptorIndex == exclusionIndices.x || acceptorIndex == exclusionIndices.y || acceptorIndex == exclusionIndices.z || acceptorIndex == exclusionIndices.w)
continue;
#endif
......@@ -179,8 +179,8 @@ __kernel void computeAcceptorForces(__global real4* restrict forceBuffers, __glo
barrier(CLK_LOCAL_MEM_FENCE);
if (acceptorIndex < NUM_ACCEPTORS) {
for (int index = 0; index < blockSize; index++) {
#ifdef USE_EXCLUSIONS
int donorIndex = donorStart+index;
#ifdef USE_EXCLUSIONS
if (donorIndex == exclusionIndices.x || donorIndex == exclusionIndices.y || donorIndex == exclusionIndices.z || donorIndex == exclusionIndices.w)
continue;
#endif
......
platforms/reference/src/ReferenceKernels.cpp
View file @ 6ed5bc4e
......@@ -159,8 +159,8 @@ static void validateVariables(const Lepton::ExpressionTreeNode& node, const set<
const Lepton::Operation& op = node.getOperation();
if (op.getId() == Lepton::Operation::VARIABLE && variables.find(op.getName()) == variables.end())
throw OpenMMException("Unknown variable in expression: "+op.getName());
for (int i = 0; i < (int) node.getChildren().size(); i++)
validateVariables(node.getChildren()[i], variables);
for (auto& child : node.getChildren())
validateVariables(child, variables);
}
/**
......@@ -214,8 +214,8 @@ void ReferenceCalcForcesAndEnergyKernel::beginComputation(ContextImpl& context,
}
else
savedForces = forceData;
for (map<string, double>::const_iterator iter = context.getParameters().begin(); iter != context.getParameters().end(); ++iter)
extractEnergyParameterDerivatives(context)[iter->first] = 0;
for (auto& param : context.getParameters())
extractEnergyParameterDerivatives(context)[param.first] = 0;
}
double ReferenceCalcForcesAndEnergyKernel::finishComputation(ContextImpl& context, bool includeForces, bool includeEnergy, int groups, bool& valid) {
......@@ -466,8 +466,8 @@ double ReferenceCalcCustomBondForceKernel::execute(ContextImpl& context, bool in
vector<Vec3>& forceData = extractForces(context);
double energy = 0;
map<string, double> globalParameters;
for (int i = 0; i < (int) globalParameterNames.size(); i++)
globalParameters[globalParameterNames[i]] = context.getParameter(globalParameterNames[i]);
for (auto& name : globalParameterNames)
globalParameters[name] = context.getParameter(name);
ReferenceCustomBondIxn bond(energyExpression, forceExpression, parameterNames, globalParameters, energyParamDerivExpressions);
if (usePeriodic)
bond.setPeriodic(extractBoxVectors(context));
......@@ -600,8 +600,8 @@ double ReferenceCalcCustomAngleForceKernel::execute(ContextImpl& context, bool i
vector<Vec3>& forceData = extractForces(context);
double energy = 0;
map<string, double> globalParameters;
for (int i = 0; i < (int) globalParameterNames.size(); i++)
globalParameters[globalParameterNames[i]] = context.getParameter(globalParameterNames[i]);
for (auto& name : globalParameterNames)
globalParameters[name] = context.getParameter(name);
ReferenceCustomAngleIxn customAngle(energyExpression, forceExpression, parameterNames, globalParameters, energyParamDerivExpressions);
if (usePeriodic)
customAngle.setPeriodic(extractBoxVectors(context));
......@@ -870,8 +870,8 @@ double ReferenceCalcCustomTorsionForceKernel::execute(ContextImpl& context, bool
vector<Vec3>& forceData = extractForces(context);
double energy = 0;
map<string, double> globalParameters;
for (int i = 0; i < (int) globalParameterNames.size(); i++)
globalParameters[globalParameterNames[i]] = context.getParameter(globalParameterNames[i]);
for (auto& name : globalParameterNames)
globalParameters[name] = context.getParameter(name);
ReferenceCustomTorsionIxn customTorsion(energyExpression, forceExpression, parameterNames, globalParameters, energyParamDerivExpressions);
if (usePeriodic)
customTorsion.setPeriodic(extractBoxVectors(context));
......@@ -1165,8 +1165,8 @@ void ReferenceCalcCustomNonbondedForceKernel::initialize(const System& system, c
// Delete the custom functions.
for (map<string, Lepton::CustomFunction*>::iterator iter = functions.begin(); iter != functions.end(); iter++)
delete iter->second;
for (auto& function : functions)
delete function.second;
// Record information for the long range correction.
......@@ -1208,11 +1208,11 @@ double ReferenceCalcCustomNonbondedForceKernel::execute(ContextImpl& context, bo
if (interactionGroups.size() > 0)
ixn.setInteractionGroups(interactionGroups);
bool globalParamsChanged = false;
for (int i = 0; i < (int) globalParameterNames.size(); i++) {
double value = context.getParameter(globalParameterNames[i]);
if (globalParamValues[globalParameterNames[i]] != value)
for (auto& name : globalParameterNames) {
double value = context.getParameter(name);
if (globalParamValues[name] != value)
globalParamsChanged = true;
globalParamValues[globalParameterNames[i]] = value;
globalParamValues[name] = value;
}
if (useSwitchingFunction)
ixn.setUseSwitchingFunction(switchingDistance);
......@@ -1459,8 +1459,8 @@ void ReferenceCalcCustomGBForceKernel::initialize(const System& system, const Cu
// Delete the custom functions.
for (map<string, Lepton::CustomFunction*>::iterator iter = functions.begin(); iter != functions.end(); iter++)
delete iter->second;
for (auto& function : functions)
delete function.second;
}
double ReferenceCalcCustomGBForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
......@@ -1477,8 +1477,8 @@ double ReferenceCalcCustomGBForceKernel::execute(ContextImpl& context, bool incl
ixn.setUseCutoff(nonbondedCutoff, *neighborList);
}
map<string, double> globalParameters;
for (int i = 0; i < (int) globalParameterNames.size(); i++)
globalParameters[globalParameterNames[i]] = context.getParameter(globalParameterNames[i]);
for (auto& name : globalParameterNames)
globalParameters[name] = context.getParameter(name);
vector<double> energyParamDerivValues(energyParamDerivNames.size()+1, 0.0);
ixn.calculateIxn(numParticles, posData, particleParamArray, exclusions, globalParameters, forceData, includeEnergy ? &energy : NULL, &energyParamDerivValues[0]);
map<string, double>& energyParamDerivs = extractEnergyParameterDerivatives(context);
......@@ -1593,8 +1593,8 @@ double ReferenceCalcCustomExternalForceKernel::execute(ContextImpl& context, boo
boxVectors = extractBoxVectors(context);
double energy = 0;
map<string, double> globalParameters;
for (int i = 0; i < (int) globalParameterNames.size(); i++)
globalParameters[globalParameterNames[i]] = context.getParameter(globalParameterNames[i]);
for (auto& name : globalParameterNames)
globalParameters[name] = context.getParameter(name);
ReferenceCustomExternalIxn force(energyExpression, forceExpressionX, forceExpressionY, forceExpressionZ, parameterNames, globalParameters);
for (int i = 0; i < numParticles; ++i)
force.calculateForce(particles[i], posData, particleParamArray[i], forceData, includeEnergy ? &energy : NULL);
......@@ -1699,8 +1699,8 @@ void ReferenceCalcCustomHbondForceKernel::initialize(const System& system, const
// Delete the custom functions.
for (map<string, Lepton::CustomFunction*>::iterator iter = functions.begin(); iter != functions.end(); iter++)
delete iter->second;
for (auto& function : functions)
delete function.second;
}
double ReferenceCalcCustomHbondForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
......@@ -1710,8 +1710,8 @@ double ReferenceCalcCustomHbondForceKernel::execute(ContextImpl& context, bool i
ixn->setPeriodic(extractBoxVectors(context));
double energy = 0;
map<string, double> globalParameters;
for (int i = 0; i < (int) globalParameterNames.size(); i++)
globalParameters[globalParameterNames[i]] = context.getParameter(globalParameterNames[i]);
for (auto& name : globalParameterNames)
globalParameters[name] = context.getParameter(name);
ixn->calculatePairIxn(posData, donorParamArray, acceptorParamArray, exclusions, globalParameters, forceData, includeEnergy ? &energy : NULL);
return energy;
}
......@@ -1803,8 +1803,8 @@ void ReferenceCalcCustomCentroidBondForceKernel::initialize(const System& system
// Delete the custom functions.
for (map<string, Lepton::CustomFunction*>::iterator iter = functions.begin(); iter != functions.end(); iter++)
delete iter->second;
for (auto& function : functions)
delete function.second;
}
double ReferenceCalcCustomCentroidBondForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
......@@ -1812,8 +1812,8 @@ double ReferenceCalcCustomCentroidBondForceKernel::execute(ContextImpl& context,
vector<Vec3>& forceData = extractForces(context);
double energy = 0;
map<string, double> globalParameters;
for (int i = 0; i < (int) globalParameterNames.size(); i++)
globalParameters[globalParameterNames[i]] = context.getParameter(globalParameterNames[i]);
for (auto& name : globalParameterNames)
globalParameters[name] = context.getParameter(name);
if (usePeriodic)
ixn->setPeriodic(extractBoxVectors(context));
vector<double> energyParamDerivValues(energyParamDerivNames.size()+1, 0.0);
......@@ -1893,8 +1893,8 @@ void ReferenceCalcCustomCompoundBondForceKernel::initialize(const System& system
// Delete the custom functions.
for (map<string, Lepton::CustomFunction*>::iterator iter = functions.begin(); iter != functions.end(); iter++)
delete iter->second;
for (auto& function : functions)
delete function.second;
}
double ReferenceCalcCustomCompoundBondForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
......@@ -1902,8 +1902,8 @@ double ReferenceCalcCustomCompoundBondForceKernel::execute(ContextImpl& context,
vector<Vec3>& forceData = extractForces(context);
double energy = 0;
map<string, double> globalParameters;
for (int i = 0; i < (int) globalParameterNames.size(); i++)
globalParameters[globalParameterNames[i]] = context.getParameter(globalParameterNames[i]);
for (auto& name : globalParameterNames)
globalParameters[name] = context.getParameter(name);
if (usePeriodic)
ixn->setPeriodic(extractBoxVectors(context));
vector<double> energyParamDerivValues(energyParamDerivNames.size()+1, 0.0);
......@@ -1966,8 +1966,8 @@ double ReferenceCalcCustomManyParticleForceKernel::execute(ContextImpl& context,
vector<Vec3>& forceData = extractForces(context);
double energy = 0;
map<string, double> globalParameters;
for (int i = 0; i < (int) globalParameterNames.size(); i++)
globalParameters[globalParameterNames[i]] = context.getParameter(globalParameterNames[i]);
for (auto& name : globalParameterNames)
globalParameters[name] = context.getParameter(name);
if (nonbondedMethod == CutoffPeriodic) {
Vec3* boxVectors = extractBoxVectors(context);
double minAllowedSize = 2*cutoffDistance;
......@@ -2232,8 +2232,8 @@ void ReferenceIntegrateCustomStepKernel::initialize(const System& system, const
for (int i = 0; i < numParticles; ++i)
masses[i] = system.getParticleMass(i);
perDofValues.resize(integrator.getNumPerDofVariables());
for (int i = 0; i < (int) perDofValues.size(); i++)
perDofValues[i].resize(numParticles);
for (auto& values : perDofValues)
values.resize(numParticles);
// Create the computation objects.
......
platforms/reference/src/SimTKReference/ReferenceAndersenThermostat.cpp
View file @ 6ed5bc4e
......@@ -66,13 +66,12 @@ using namespace OpenMM;
double temperature, double collisionFrequency, double stepSize) const {
const double collisionProbability = 1.0f - exp(-collisionFrequency*stepSize);
for (int i = 0; i < (int) atomGroups.size(); ++i) {
for (auto& group : atomGroups) {
if (SimTKOpenMMUtilities::getUniformlyDistributedRandomNumber() < collisionProbability) {
// A collision occurred, so set the velocities to new values chosen from a Boltzmann distribution.
for (int j = 0; j < (int) atomGroups[i].size(); j++) {
int atom = atomGroups[i][j];
for (int atom : group) {
if (atomMasses[atom] != 0) {
const double velocityScale = static_cast<double>(sqrt(BOLTZ*temperature/atomMasses[atom]));
atomVelocities[atom][0] = velocityScale*SimTKOpenMMUtilities::getNormallyDistributedRandomNumber();
......
platforms/reference/src/SimTKReference/ReferenceCCMAAlgorithm.cpp
View file @ 6ed5bc4e
......@@ -127,8 +127,8 @@ ReferenceCCMAAlgorithm::ReferenceCCMAAlgorithm(int numberOfAtoms,
// We didn't find one, so look for an angle force field term.
const vector<int>& angleCandidates = atomAngles[atomb];
for (vector<int>::const_iterator iter = angleCandidates.begin(); iter != angleCandidates.end(); iter++) {
const AngleInfo& angle = angles[*iter];
for (int candidate : angleCandidates) {
const AngleInfo& angle = angles[candidate];
if ((angle.atom1 == atoma && angle.atom3 == atomc) || (angle.atom3 == atoma && angle.atom1 == atomc)) {
matrix[j].push_back(pair<int, double>(k, scale*cos(angle.angle)));
break;
......@@ -145,8 +145,7 @@ ReferenceCCMAAlgorithm::ReferenceCCMAAlgorithm(int numberOfAtoms,
vector<double> matrixValue;
for (int i = 0; i < numberOfConstraints; i++) {
matrixRowStart.push_back(matrixValue.size());
for (int j = 0; j < (int) matrix[i].size(); j++) {
pair<int, double> element = matrix[i][j];
for (auto& element : matrix[i]) {
matrixColIndex.push_back(element.first);
matrixValue.push_back(element.second);
}
......@@ -292,10 +291,8 @@ void ReferenceCCMAAlgorithm::applyConstraints(vector<Vec3>& atomCoordinates,
if (_matrix.size() > 0) {
for (int i = 0; i < _numberOfConstraints; i++) {
double sum = 0.0;
for (int j = 0; j < (int) _matrix[i].size(); j++) {
pair<int, double> element = _matrix[i][j];
for (auto& element : _matrix[i])
sum += element.second*constraintDelta[element.first];
}
tempDelta[i] = sum;
}
constraintDelta = tempDelta;
......
platforms/reference/src/SimTKReference/ReferenceCustomAngleIxn.cpp
View file @ 6ed5bc4e
......@@ -47,10 +47,10 @@ ReferenceCustomAngleIxn::ReferenceCustomAngleIxn(const Lepton::CompiledExpressio
expressionSet.registerExpression(this->energyParamDerivExpressions[i]);
thetaIndex = expressionSet.getVariableIndex("theta");
numParameters = parameterNames.size();
for (int i = 0; i < (int) numParameters; i++)
angleParamIndex.push_back(expressionSet.getVariableIndex(parameterNames[i]));
for (map<string, double>::const_iterator iter = globalParameters.begin(); iter != globalParameters.end(); ++iter)
expressionSet.setVariable(expressionSet.getVariableIndex(iter->first), iter->second);
for (auto& param : parameterNames)
angleParamIndex.push_back(expressionSet.getVariableIndex(param));
for (auto& param : globalParameters)
expressionSet.setVariable(expressionSet.getVariableIndex(param.first), param.second);
}
/**---------------------------------------------------------------------------------------
......
platforms/reference/src/SimTKReference/ReferenceCustomBondIxn.cpp
View file @ 6ed5bc4e
......@@ -48,10 +48,10 @@ ReferenceCustomBondIxn::ReferenceCustomBondIxn(const Lepton::CompiledExpression&
expressionSet.registerExpression(this->energyParamDerivExpressions[i]);
rIndex = expressionSet.getVariableIndex("r");
numParameters = parameterNames.size();
for (int i = 0; i < (int) numParameters; i++)
bondParamIndex.push_back(expressionSet.getVariableIndex(parameterNames[i]));
for (map<string, double>::const_iterator iter = globalParameters.begin(); iter != globalParameters.end(); ++iter)
expressionSet.setVariable(expressionSet.getVariableIndex(iter->first), iter->second);
for (auto& param : parameterNames)
bondParamIndex.push_back(expressionSet.getVariableIndex(param));
for (auto& param : globalParameters)
expressionSet.setVariable(expressionSet.getVariableIndex(param.first), param.second);
}
/**---------------------------------------------------------------------------------------
......
platforms/reference/src/SimTKReference/ReferenceCustomCentroidBondIxn.cpp
View file @ 6ed5bc4e
......@@ -56,12 +56,12 @@ ReferenceCustomCentroidBondIxn::ReferenceCustomCentroidBondIxn(int numGroupsPerB
positionTerms.push_back(ReferenceCustomCentroidBondIxn::PositionTermInfo(yname.str(), i, 1, energyExpression.differentiate(yname.str()).createCompiledExpression()));
positionTerms.push_back(ReferenceCustomCentroidBondIxn::PositionTermInfo(zname.str(), i, 2, energyExpression.differentiate(zname.str()).createCompiledExpression()));
}
for (map<string, vector<int> >::const_iterator iter = distances.begin(); iter != distances.end(); ++iter)
distanceTerms.push_back(ReferenceCustomCentroidBondIxn::DistanceTermInfo(iter->first, iter->second, energyExpression.differentiate(iter->first).createCompiledExpression()));
for (map<string, vector<int> >::const_iterator iter = angles.begin(); iter != angles.end(); ++iter)
angleTerms.push_back(ReferenceCustomCentroidBondIxn::AngleTermInfo(iter->first, iter->second, energyExpression.differentiate(iter->first).createCompiledExpression()));
for (map<string, vector<int> >::const_iterator iter = dihedrals.begin(); iter != dihedrals.end(); ++iter)
dihedralTerms.push_back(ReferenceCustomCentroidBondIxn::DihedralTermInfo(iter->first, iter->second, energyExpression.differentiate(iter->first).createCompiledExpression()));
for (auto& term : distances)
distanceTerms.push_back(ReferenceCustomCentroidBondIxn::DistanceTermInfo(term.first, term.second, energyExpression.differentiate(term.first).createCompiledExpression()));
for (auto& term : angles)
angleTerms.push_back(ReferenceCustomCentroidBondIxn::AngleTermInfo(term.first, term.second, energyExpression.differentiate(term.first).createCompiledExpression()));
for (auto& term : dihedrals)
dihedralTerms.push_back(ReferenceCustomCentroidBondIxn::DihedralTermInfo(term.first, term.second, energyExpression.differentiate(term.first).createCompiledExpression()));
for (int i = 0; i < positionTerms.size(); i++) {
expressionSet.registerExpression(positionTerms[i].forceExpression);
positionTerms[i].index = expressionSet.getVariableIndex(positionTerms[i].name);
......@@ -108,8 +108,8 @@ void ReferenceCustomCentroidBondIxn::calculatePairIxn(vector<Vec3>& atomCoordina
// Compute the forces on groups.
for (map<string, double>::const_iterator iter = globalParameters.begin(); iter != globalParameters.end(); ++iter)
expressionSet.setVariable(expressionSet.getVariableIndex(iter->first), iter->second);
for (auto& param : globalParameters)
expressionSet.setVariable(expressionSet.getVariableIndex(param.first), param.second);
vector<Vec3> groupForces(numGroups);
int numBonds = bondGroups.size();
for (int bond = 0; bond < numBonds; bond++) {
......@@ -131,23 +131,18 @@ void ReferenceCustomCentroidBondIxn::calculateOneIxn(int bond, vector<Vec3>& gro
// Compute all of the variables the energy can depend on.
const vector<int>& groups = bondGroups[bond];
for (int i = 0; i < (int) positionTerms.size(); i++) {
const PositionTermInfo& term = positionTerms[i];
for (auto& term : positionTerms)
expressionSet.setVariable(term.index, groupCenters[groups[term.group]][term.component]);
}
for (int i = 0; i < (int) distanceTerms.size(); i++) {
const DistanceTermInfo& term = distanceTerms[i];
for (auto& term : distanceTerms) {
computeDelta(groups[term.g1], groups[term.g2], term.delta, groupCenters);
expressionSet.setVariable(term.index, term.delta[ReferenceForce::RIndex]);
}
for (int i = 0; i < (int) angleTerms.size(); i++) {
const AngleTermInfo& term = angleTerms[i];
for (auto& term : angleTerms) {
computeDelta(groups[term.g1], groups[term.g2], term.delta1, groupCenters);
computeDelta(groups[term.g3], groups[term.g2], term.delta2, groupCenters);
expressionSet.setVariable(term.index, computeAngle(term.delta1, term.delta2));
}
for (int i = 0; i < (int) dihedralTerms.size(); i++) {
const DihedralTermInfo& term = dihedralTerms[i];
for (auto& term : dihedralTerms) {
computeDelta(groups[term.g2], groups[term.g1], term.delta1, groupCenters);
computeDelta(groups[term.g2], groups[term.g3], term.delta2, groupCenters);
computeDelta(groups[term.g4], groups[term.g3], term.delta3, groupCenters);
......@@ -158,15 +153,12 @@ void ReferenceCustomCentroidBondIxn::calculateOneIxn(int bond, vector<Vec3>& gro
// Apply forces based on individual particle coordinates.
for (int i = 0; i < (int) positionTerms.size(); i++) {
const PositionTermInfo& term = positionTerms[i];
for (auto& term : positionTerms)
forces[groups[term.group]][term.component] -= term.forceExpression.evaluate();
}
// Apply forces based on distances.
for (int i = 0; i < (int) distanceTerms.size(); i++) {
const DistanceTermInfo& term = distanceTerms[i];
for (auto& term : distanceTerms) {
double dEdR = term.forceExpression.evaluate()/(term.delta[ReferenceForce::RIndex]);
for (int i = 0; i < 3; i++) {
double force = -dEdR*term.delta[i];
......@@ -177,8 +169,7 @@ void ReferenceCustomCentroidBondIxn::calculateOneIxn(int bond, vector<Vec3>& gro
// Apply forces based on angles.
for (int i = 0; i < (int) angleTerms.size(); i++) {
const AngleTermInfo& term = angleTerms[i];
for (auto& term : angleTerms) {
double dEdTheta = term.forceExpression.evaluate();
double thetaCross[ReferenceForce::LastDeltaRIndex];
SimTKOpenMMUtilities::crossProductVector3(term.delta1, term.delta2, thetaCross);
......@@ -204,8 +195,7 @@ void ReferenceCustomCentroidBondIxn::calculateOneIxn(int bond, vector<Vec3>& gro
// Apply forces based on dihedrals.
for (int i = 0; i < (int) dihedralTerms.size(); i++) {
const DihedralTermInfo& term = dihedralTerms[i];
for (auto& term : dihedralTerms) {
double dEdTheta = term.forceExpression.evaluate();
double internalF[4][3];
double forceFactors[4];
......
platforms/reference/src/SimTKReference/ReferenceCustomCompoundBondIxn.cpp
View file @ 6ed5bc4e
......@@ -61,12 +61,12 @@ ReferenceCustomCompoundBondIxn::ReferenceCustomCompoundBondIxn(int numParticlesP
particleTerms.push_back(ReferenceCustomCompoundBondIxn::ParticleTermInfo(yname.str(), i, 1, energyExpression.differentiate(yname.str()).createCompiledExpression()));
particleTerms.push_back(ReferenceCustomCompoundBondIxn::ParticleTermInfo(zname.str(), i, 2, energyExpression.differentiate(zname.str()).createCompiledExpression()));
}
for (map<string, vector<int> >::const_iterator iter = distances.begin(); iter != distances.end(); ++iter)
distanceTerms.push_back(ReferenceCustomCompoundBondIxn::DistanceTermInfo(iter->first, iter->second, energyExpression.differentiate(iter->first).createCompiledExpression()));
for (map<string, vector<int> >::const_iterator iter = angles.begin(); iter != angles.end(); ++iter)
angleTerms.push_back(ReferenceCustomCompoundBondIxn::AngleTermInfo(iter->first, iter->second, energyExpression.differentiate(iter->first).createCompiledExpression()));
for (map<string, vector<int> >::const_iterator iter = dihedrals.begin(); iter != dihedrals.end(); ++iter)
dihedralTerms.push_back(ReferenceCustomCompoundBondIxn::DihedralTermInfo(iter->first, iter->second, energyExpression.differentiate(iter->first).createCompiledExpression()));
for (auto& term : distances)
distanceTerms.push_back(ReferenceCustomCompoundBondIxn::DistanceTermInfo(term.first, term.second, energyExpression.differentiate(term.first).createCompiledExpression()));
for (auto& term : angles)
angleTerms.push_back(ReferenceCustomCompoundBondIxn::AngleTermInfo(term.first, term.second, energyExpression.differentiate(term.first).createCompiledExpression()));
for (auto& term : dihedrals)
dihedralTerms.push_back(ReferenceCustomCompoundBondIxn::DihedralTermInfo(term.first, term.second, energyExpression.differentiate(term.first).createCompiledExpression()));
for (int i = 0; i < particleTerms.size(); i++) {
expressionSet.registerExpression(particleTerms[i].forceExpression);
particleTerms[i].index = expressionSet.getVariableIndex(particleTerms[i].name);
......@@ -119,8 +119,8 @@ void ReferenceCustomCompoundBondIxn::setPeriodic(OpenMM::Vec3* vectors) {
void ReferenceCustomCompoundBondIxn::calculatePairIxn(vector<Vec3>& atomCoordinates, double** bondParameters,
const map<string, double>& globalParameters, vector<Vec3>& forces,
double* totalEnergy, double* energyParamDerivs) {
for (map<string, double>::const_iterator iter = globalParameters.begin(); iter != globalParameters.end(); ++iter)
expressionSet.setVariable(expressionSet.getVariableIndex(iter->first), iter->second);
for (auto& param : globalParameters)
expressionSet.setVariable(expressionSet.getVariableIndex(param.first), param.second);
int numBonds = bondAtoms.size();
for (int bond = 0; bond < numBonds; bond++) {
for (int i = 0; i < numParameters; i++)
......@@ -146,23 +146,18 @@ void ReferenceCustomCompoundBondIxn::calculateOneIxn(int bond, vector<Vec3>& ato
// Compute all of the variables the energy can depend on.
const vector<int>& atoms = bondAtoms[bond];
for (int i = 0; i < (int) particleTerms.size(); i++) {
const ParticleTermInfo& term = particleTerms[i];
for (auto& term : particleTerms)
expressionSet.setVariable(term.index, atomCoordinates[atoms[term.atom]][term.component]);
}
for (int i = 0; i < (int) distanceTerms.size(); i++) {
const DistanceTermInfo& term = distanceTerms[i];
for (auto& term : distanceTerms) {
computeDelta(atoms[term.p1], atoms[term.p2], term.delta, atomCoordinates);
expressionSet.setVariable(term.index, term.delta[ReferenceForce::RIndex]);
}
for (int i = 0; i < (int) angleTerms.size(); i++) {
const AngleTermInfo& term = angleTerms[i];
for (auto& term : angleTerms) {
computeDelta(atoms[term.p1], atoms[term.p2], term.delta1, atomCoordinates);
computeDelta(atoms[term.p3], atoms[term.p2], term.delta2, atomCoordinates);
expressionSet.setVariable(term.index, computeAngle(term.delta1, term.delta2));
}
for (int i = 0; i < (int) dihedralTerms.size(); i++) {
const DihedralTermInfo& term = dihedralTerms[i];
for (auto& term : dihedralTerms) {
computeDelta(atoms[term.p2], atoms[term.p1], term.delta1, atomCoordinates);
computeDelta(atoms[term.p2], atoms[term.p3], term.delta2, atomCoordinates);
computeDelta(atoms[term.p4], atoms[term.p3], term.delta3, atomCoordinates);
......@@ -173,15 +168,12 @@ void ReferenceCustomCompoundBondIxn::calculateOneIxn(int bond, vector<Vec3>& ato
// Apply forces based on individual particle coordinates.
for (int i = 0; i < (int) particleTerms.size(); i++) {
const ParticleTermInfo& term = particleTerms[i];
for (auto& term : particleTerms)
forces[atoms[term.atom]][term.component] -= term.forceExpression.evaluate();
}
// Apply forces based on distances.
for (int i = 0; i < (int) distanceTerms.size(); i++) {
const DistanceTermInfo& term = distanceTerms[i];
for (auto& term : distanceTerms) {
double dEdR = term.forceExpression.evaluate()/(term.delta[ReferenceForce::RIndex]);
for (int i = 0; i < 3; i++) {
double force = -dEdR*term.delta[i];
......@@ -192,8 +184,7 @@ void ReferenceCustomCompoundBondIxn::calculateOneIxn(int bond, vector<Vec3>& ato
// Apply forces based on angles.
for (int i = 0; i < (int) angleTerms.size(); i++) {
const AngleTermInfo& term = angleTerms[i];
for (auto& term : angleTerms) {
double dEdTheta = term.forceExpression.evaluate();
double thetaCross[ReferenceForce::LastDeltaRIndex];
SimTKOpenMMUtilities::crossProductVector3(term.delta1, term.delta2, thetaCross);
......@@ -219,8 +210,7 @@ void ReferenceCustomCompoundBondIxn::calculateOneIxn(int bond, vector<Vec3>& ato
// Apply forces based on dihedrals.
for (int i = 0; i < (int) dihedralTerms.size(); i++) {
const DihedralTermInfo& term = dihedralTerms[i];
for (auto& term : dihedralTerms) {
double dEdTheta = term.forceExpression.evaluate();
double internalF[4][3];
double forceFactors[4];
......
platforms/reference/src/SimTKReference/ReferenceCustomDynamics.cpp
View file @ 6ed5bc4e
......@@ -175,8 +175,8 @@ ExpressionTreeNode ReferenceCustomDynamics::replaceDerivFunctions(const Expressi
}
else {
vector<ExpressionTreeNode> children;
for (int i = 0; i < (int) node.getChildren().size(); i++)
children.push_back(replaceDerivFunctions(node.getChildren()[i], context));
for (auto& child : node.getChildren())
children.push_back(replaceDerivFunctions(child, context));
return ExpressionTreeNode(op.clone(), children);
}
}
......@@ -204,8 +204,8 @@ void ReferenceCustomDynamics::update(ContextImpl& context, int numberOfAtoms, ve
initialize(context, masses, globals);
int numSteps = stepType.size();
globals.insert(context.getParameters().begin(), context.getParameters().end());
for (map<string, double>::const_iterator iter = globals.begin(); iter != globals.end(); ++iter)
expressionSet.setVariable(expressionSet.getVariableIndex(iter->first), iter->second);
for (auto& global : globals)
expressionSet.setVariable(expressionSet.getVariableIndex(global.first), global.second);
oldPos = atomCoordinates;
// Loop over steps and execute them.
......@@ -276,8 +276,8 @@ void ReferenceCustomDynamics::update(ContextImpl& context, int numberOfAtoms, ve
recordChangedParameters(context, globals);
context.updateContextState();
globals.insert(context.getParameters().begin(), context.getParameters().end());
for (map<string, double>::const_iterator iter = globals.begin(); iter != globals.end(); ++iter)
expressionSet.setVariable(expressionSet.getVariableIndex(iter->first), iter->second);
for (auto& global : globals)
expressionSet.setVariable(expressionSet.getVariableIndex(global.first), global.second);
break;
}
case CustomIntegrator::IfBlockStart: {
......@@ -355,10 +355,10 @@ bool ReferenceCustomDynamics::evaluateCondition(int step) {
* Check which context parameters have changed and register them with the context.
*/
void ReferenceCustomDynamics::recordChangedParameters(OpenMM::ContextImpl& context, std::map<std::string, double>& globals) {
for (map<string, double>::const_iterator iter = context.getParameters().begin(); iter != context.getParameters().end(); ++iter) {
string name = iter->first;
for (auto& param : context.getParameters()) {
string name = param.first;
double value = globals[name];
if (value != iter->second)
if (value != param.second)
context.setParameter(name, globals[name]);
}
}
......@@ -385,8 +385,8 @@ double ReferenceCustomDynamics::computeKineticEnergy(OpenMM::ContextImpl& contex
if (invalidatesForces.size() == 0)
initialize(context, masses, globals);
globals.insert(context.getParameters().begin(), context.getParameters().end());
for (map<string, double>::const_iterator iter = globals.begin(); iter != globals.end(); ++iter)
expressionSet.setVariable(expressionSet.getVariableIndex(iter->first), iter->second);
for (auto& global : globals)
expressionSet.setVariable(expressionSet.getVariableIndex(global.first), global.second);
if (kineticEnergyNeedsForce) {
energy = context.calcForcesAndEnergy(true, true, -1);
forcesAreValid = true;
......
platforms/reference/src/SimTKReference/ReferenceCustomExternalIxn.cpp
View file @ 6ed5bc4e
......@@ -57,17 +57,17 @@ ReferenceCustomExternalIxn::ReferenceCustomExternalIxn(const Lepton::CompiledExp
forceZY = ReferenceForce::getVariablePointer(this->forceExpressionZ, "y");
forceZZ = ReferenceForce::getVariablePointer(this->forceExpressionZ, "z");
numParameters = parameterNames.size();
for (int i = 0; i < (int) numParameters; i++) {
energyParams.push_back(ReferenceForce::getVariablePointer(this->energyExpression, parameterNames[i]));
forceXParams.push_back(ReferenceForce::getVariablePointer(this->forceExpressionX, parameterNames[i]));
forceYParams.push_back(ReferenceForce::getVariablePointer(this->forceExpressionY, parameterNames[i]));
forceZParams.push_back(ReferenceForce::getVariablePointer(this->forceExpressionZ, parameterNames[i]));
for (auto& param : parameterNames) {
energyParams.push_back(ReferenceForce::getVariablePointer(this->energyExpression, param));
forceXParams.push_back(ReferenceForce::getVariablePointer(this->forceExpressionX, param));
forceYParams.push_back(ReferenceForce::getVariablePointer(this->forceExpressionY, param));
forceZParams.push_back(ReferenceForce::getVariablePointer(this->forceExpressionZ, param));
}
for (map<string, double>::const_iterator iter = globalParameters.begin(); iter != globalParameters.end(); ++iter) {
ReferenceForce::setVariable(ReferenceForce::getVariablePointer(this->energyExpression, iter->first), iter->second);
ReferenceForce::setVariable(ReferenceForce::getVariablePointer(this->forceExpressionX, iter->first), iter->second);
ReferenceForce::setVariable(ReferenceForce::getVariablePointer(this->forceExpressionY, iter->first), iter->second);
ReferenceForce::setVariable(ReferenceForce::getVariablePointer(this->forceExpressionZ, iter->first), iter->second);
for (auto& param : globalParameters) {
ReferenceForce::setVariable(ReferenceForce::getVariablePointer(this->energyExpression, param.first), param.second);
ReferenceForce::setVariable(ReferenceForce::getVariablePointer(this->forceExpressionX, param.first), param.second);
ReferenceForce::setVariable(ReferenceForce::getVariablePointer(this->forceExpressionY, param.first), param.second);
ReferenceForce::setVariable(ReferenceForce::getVariablePointer(this->forceExpressionZ, param.first), param.second);
}
}
......
platforms/reference/src/SimTKReference/ReferenceCustomGBIxn.cpp
View file @ 6ed5bc4e
......@@ -84,19 +84,19 @@ ReferenceCustomGBIxn::ReferenceCustomGBIxn(const vector<Lepton::CompiledExpressi
xIndex = expressionSet.getVariableIndex("x");
yIndex = expressionSet.getVariableIndex("y");
zIndex = expressionSet.getVariableIndex("z");
for (int i = 0; i < (int) parameterNames.size(); i++) {
paramIndex.push_back(expressionSet.getVariableIndex(parameterNames[i]));
for (auto& param : parameterNames) {
paramIndex.push_back(expressionSet.getVariableIndex(param));
for (int j = 1; j < 3; j++) {
stringstream name;
name << parameterNames[i] << j;
name << param << j;
particleParamIndex.push_back(expressionSet.getVariableIndex(name.str()));
}
}
for (int i = 0; i < (int) valueNames.size(); i++) {
valueIndex.push_back(expressionSet.getVariableIndex(valueNames[i]));
for (auto& value : valueNames) {
valueIndex.push_back(expressionSet.getVariableIndex(value));
for (int j = 1; j < 3; j++) {
stringstream name;
name << valueNames[i] << j;
name << value << j;
particleValueIndex.push_back(expressionSet.getVariableIndex(name.str()));
}
}
......@@ -153,8 +153,8 @@ ReferenceCustomGBIxn::~ReferenceCustomGBIxn() {
void ReferenceCustomGBIxn::calculateIxn(int numberOfAtoms, vector<Vec3>& atomCoordinates, double** atomParameters,
const vector<set<int> >& exclusions, map<string, double>& globalParameters, vector<Vec3>& forces,
double* totalEnergy, double* energyParamDerivs) {
for (map<string, double>::const_iterator iter = globalParameters.begin(); iter != globalParameters.end(); ++iter)
expressionSet.setVariable(expressionSet.getVariableIndex(iter->first), iter->second);
for (auto& param : globalParameters)
expressionSet.setVariable(expressionSet.getVariableIndex(param.first), param.second);
// Initialize arrays for storing values.
......@@ -225,8 +225,7 @@ void ReferenceCustomGBIxn::calculateParticlePairValue(int index, int numAtoms, v
if (cutoff) {
// Loop over all pairs in the neighbor list.
for (int i = 0; i < (int) neighborList->size(); i++) {
OpenMM::AtomPair pair = (*neighborList)[i];
for (auto& pair : *neighborList) {
if (useExclusions && exclusions[pair.first].find(pair.second) != exclusions[pair.first].end())
continue;
calculateOnePairValue(index, pair.first, pair.second, atomCoordinates, atomParameters);
......@@ -306,8 +305,7 @@ void ReferenceCustomGBIxn::calculateParticlePairEnergyTerm(int index, int numAto
if (cutoff) {
// Loop over all pairs in the neighbor list.
for (int i = 0; i < (int) neighborList->size(); i++) {
OpenMM::AtomPair pair = (*neighborList)[i];
for (auto& pair : *neighborList) {
if (useExclusions && exclusions[pair.first].find(pair.second) != exclusions[pair.first].end())
continue;
calculateOnePairEnergyTerm(index, pair.first, pair.second, atomCoordinates, atomParameters, forces, totalEnergy, energyParamDerivs);
......@@ -377,8 +375,7 @@ void ReferenceCustomGBIxn::calculateChainRuleForces(int numAtoms, vector<Vec3>&
if (cutoff) {
// Loop over all pairs in the neighbor list.
for (int i = 0; i < (int) neighborList->size(); i++) {
OpenMM::AtomPair pair = (*neighborList)[i];
for (auto& pair : *neighborList) {
bool isExcluded = (exclusions[pair.first].find(pair.second) != exclusions[pair.first].end());
calculateOnePairChainRule(pair.first, pair.second, atomCoordinates, atomParameters, forces, isExcluded);
calculateOnePairChainRule(pair.second, pair.first, atomCoordinates, atomParameters, forces, isExcluded);
......
platforms/reference/src/SimTKReference/ReferenceCustomHbondIxn.cpp
View file @ 6ed5bc4e
......@@ -49,12 +49,12 @@ ReferenceCustomHbondIxn::ReferenceCustomHbondIxn(const vector<vector<int> >& don
const map<string, vector<int> >& distances, const map<string, vector<int> >& angles, const map<string, vector<int> >& dihedrals) :
cutoff(false), periodic(false), donorAtoms(donorAtoms), acceptorAtoms(acceptorAtoms), energyExpression(energyExpression.createProgram()),
donorParamNames(donorParameterNames), acceptorParamNames(acceptorParameterNames) {
for (map<string, vector<int> >::const_iterator iter = distances.begin(); iter != distances.end(); ++iter)
distanceTerms.push_back(ReferenceCustomHbondIxn::DistanceTermInfo(iter->first, iter->second, energyExpression.differentiate(iter->first).optimize().createProgram()));
for (map<string, vector<int> >::const_iterator iter = angles.begin(); iter != angles.end(); ++iter)
angleTerms.push_back(ReferenceCustomHbondIxn::AngleTermInfo(iter->first, iter->second, energyExpression.differentiate(iter->first).optimize().createProgram()));
for (map<string, vector<int> >::const_iterator iter = dihedrals.begin(); iter != dihedrals.end(); ++iter)
dihedralTerms.push_back(ReferenceCustomHbondIxn::DihedralTermInfo(iter->first, iter->second, energyExpression.differentiate(iter->first).optimize().createProgram()));
for (auto& term : distances)
distanceTerms.push_back(ReferenceCustomHbondIxn::DistanceTermInfo(term.first, term.second, energyExpression.differentiate(term.first).optimize().createProgram()));
for (auto& term : angles)
angleTerms.push_back(ReferenceCustomHbondIxn::AngleTermInfo(term.first, term.second, energyExpression.differentiate(term.first).optimize().createProgram()));
for (auto& term : dihedrals)
dihedralTerms.push_back(ReferenceCustomHbondIxn::DihedralTermInfo(term.first, term.second, energyExpression.differentiate(term.first).optimize().createProgram()));
}
/**---------------------------------------------------------------------------------------
......
platforms/reference/src/SimTKReference/ReferenceCustomManyParticleIxn.cpp
View file @ 6ed5bc4e
......@@ -60,8 +60,8 @@ ReferenceCustomManyParticleIxn::ReferenceCustomManyParticleIxn(const CustomManyP
// Delete the custom functions.
for (map<string, Lepton::CustomFunction*>::iterator iter = functions.begin(); iter != functions.end(); iter++)
delete iter->second;
for (auto& function : functions)
delete function.second;
// Differentiate the energy to get expressions for the force.
......@@ -80,12 +80,12 @@ ReferenceCustomManyParticleIxn::ReferenceCustomManyParticleIxn(const CustomManyP
particleParamNames[i].push_back(paramname.str());
}
}
for (map<string, vector<int> >::const_iterator iter = distances.begin(); iter != distances.end(); ++iter)
distanceTerms.push_back(ReferenceCustomManyParticleIxn::DistanceTermInfo(iter->first, iter->second, energyExpr.differentiate(iter->first).optimize().createProgram()));
for (map<string, vector<int> >::const_iterator iter = angles.begin(); iter != angles.end(); ++iter)
angleTerms.push_back(ReferenceCustomManyParticleIxn::AngleTermInfo(iter->first, iter->second, energyExpr.differentiate(iter->first).optimize().createProgram()));
for (map<string, vector<int> >::const_iterator iter = dihedrals.begin(); iter != dihedrals.end(); ++iter)
dihedralTerms.push_back(ReferenceCustomManyParticleIxn::DihedralTermInfo(iter->first, iter->second, energyExpr.differentiate(iter->first).optimize().createProgram()));
for (auto& term : distances)
distanceTerms.push_back(ReferenceCustomManyParticleIxn::DistanceTermInfo(term.first, term.second, energyExpr.differentiate(term.first).optimize().createProgram()));
for (auto& term : angles)
angleTerms.push_back(ReferenceCustomManyParticleIxn::AngleTermInfo(term.first, term.second, energyExpr.differentiate(term.first).optimize().createProgram()));
for (auto& term : dihedrals)
dihedralTerms.push_back(ReferenceCustomManyParticleIxn::DihedralTermInfo(term.first, term.second, energyExpr.differentiate(term.first).optimize().createProgram()));
// Record exclusions.
......@@ -192,23 +192,18 @@ void ReferenceCustomManyParticleIxn::calculateOneIxn(const vector<int>& particle
// Compute all of the variables the energy can depend on.
for (int i = 0; i < (int) particleTerms.size(); i++) {
const ParticleTermInfo& term = particleTerms[i];
for (auto& term : particleTerms)
variables[term.name] = atomCoordinates[permutedParticles[term.atom]][term.component];
}
for (int i = 0; i < (int) distanceTerms.size(); i++) {
const DistanceTermInfo& term = distanceTerms[i];
for (auto& term : distanceTerms) {
computeDelta(permutedParticles[term.p1], permutedParticles[term.p2], term.delta, atomCoordinates);
variables[term.name] = term.delta[ReferenceForce::RIndex];
}
for (int i = 0; i < (int) angleTerms.size(); i++) {
const AngleTermInfo& term = angleTerms[i];
for (auto& term : angleTerms) {
computeDelta(permutedParticles[term.p1], permutedParticles[term.p2], term.delta1, atomCoordinates);
computeDelta(permutedParticles[term.p3], permutedParticles[term.p2], term.delta2, atomCoordinates);
variables[term.name] = computeAngle(term.delta1, term.delta2);
}
for (int i = 0; i < (int) dihedralTerms.size(); i++) {
const DihedralTermInfo& term = dihedralTerms[i];
for (auto& term : dihedralTerms) {
computeDelta(permutedParticles[term.p2], permutedParticles[term.p1], term.delta1, atomCoordinates);
computeDelta(permutedParticles[term.p2], permutedParticles[term.p3], term.delta2, atomCoordinates);
computeDelta(permutedParticles[term.p4], permutedParticles[term.p3], term.delta3, atomCoordinates);
......@@ -219,15 +214,12 @@ void ReferenceCustomManyParticleIxn::calculateOneIxn(const vector<int>& particle
// Apply forces based on individual particle coordinates.
for (int i = 0; i < (int) particleTerms.size(); i++) {
const ParticleTermInfo& term = particleTerms[i];
for (auto& term : particleTerms)
forces[permutedParticles[term.atom]][term.component] -= term.forceExpression.evaluate(variables);
}
// Apply forces based on distances.
for (int i = 0; i < (int) distanceTerms.size(); i++) {
const DistanceTermInfo& term = distanceTerms[i];
for (auto& term : distanceTerms) {
double dEdR = term.forceExpression.evaluate(variables)/(term.delta[ReferenceForce::RIndex]);
for (int i = 0; i < 3; i++) {
double force = -dEdR*term.delta[i];
......@@ -238,8 +230,7 @@ void ReferenceCustomManyParticleIxn::calculateOneIxn(const vector<int>& particle
// Apply forces based on angles.
for (int i = 0; i < (int) angleTerms.size(); i++) {
const AngleTermInfo& term = angleTerms[i];
for (auto& term : angleTerms) {
double dEdTheta = term.forceExpression.evaluate(variables);
double thetaCross[ReferenceForce::LastDeltaRIndex];
SimTKOpenMMUtilities::crossProductVector3(term.delta1, term.delta2, thetaCross);
......@@ -265,8 +256,7 @@ void ReferenceCustomManyParticleIxn::calculateOneIxn(const vector<int>& particle
// Apply forces based on dihedrals.
for (int i = 0; i < (int) dihedralTerms.size(); i++) {
const DihedralTermInfo& term = dihedralTerms[i];
for (auto& term : dihedralTerms) {
double dEdTheta = term.forceExpression.evaluate(variables);
double internalF[4][3];
double forceFactors[4];
......
platforms/reference/src/SimTKReference/ReferenceCustomNonbondedIxn.cpp
View file @ 6ed5bc4e
......@@ -53,10 +53,10 @@ ReferenceCustomNonbondedIxn::ReferenceCustomNonbondedIxn(const Lepton::CompiledE
for (int i = 0; i < this->energyParamDerivExpressions.size(); i++)
expressionSet.registerExpression(this->energyParamDerivExpressions[i]);
rIndex = expressionSet.getVariableIndex("r");
for (int i = 0; i < (int) paramNames.size(); i++) {
for (auto& param : paramNames) {
for (int j = 1; j < 3; j++) {
stringstream name;
name << paramNames[i] << j;
name << param << j;
particleParamIndex.push_back(expressionSet.getVariableIndex(name.str()));
}
}
......@@ -159,14 +159,14 @@ void ReferenceCustomNonbondedIxn::calculatePairIxn(int numberOfAtoms, vector<Vec
double* fixedParameters, const map<string, double>& globalParameters, vector<Vec3>& forces,
double* energyByAtom, double* totalEnergy, double* energyParamDerivs) {
for (map<string, double>::const_iterator iter = globalParameters.begin(); iter != globalParameters.end(); ++iter)
expressionSet.setVariable(expressionSet.getVariableIndex(iter->first), iter->second);
for (auto& param : globalParameters)
expressionSet.setVariable(expressionSet.getVariableIndex(param.first), param.second);
if (interactionGroups.size() > 0) {
// The user has specified interaction groups, so compute only the requested interactions.
for (int group = 0; group < (int) interactionGroups.size(); group++) {
const set<int>& set1 = interactionGroups[group].first;
const set<int>& set2 = interactionGroups[group].second;
for (auto& group : interactionGroups) {
const set<int>& set1 = group.first;
const set<int>& set2 = group.second;
for (set<int>::const_iterator atom1 = set1.begin(); atom1 != set1.end(); ++atom1) {
for (set<int>::const_iterator atom2 = set2.begin(); atom2 != set2.end(); ++atom2) {
if (*atom1 == *atom2 || exclusions[*atom1].find(*atom2) != exclusions[*atom1].end())
......@@ -185,8 +185,7 @@ void ReferenceCustomNonbondedIxn::calculatePairIxn(int numberOfAtoms, vector<Vec
else if (cutoff) {
// We are using a cutoff, so get the interactions from the neighbor list.
for (int i = 0; i < (int) neighborList->size(); i++) {
OpenMM::AtomPair pair = (*neighborList)[i];
for (auto& pair : *neighborList) {
for (int j = 0; j < (int) paramNames.size(); j++) {
expressionSet.setVariable(particleParamIndex[j*2], atomParameters[pair.first][j]);
expressionSet.setVariable(particleParamIndex[j*2+1], atomParameters[pair.second][j]);
......
platforms/reference/src/SimTKReference/ReferenceCustomTorsionIxn.cpp
View file @ 6ed5bc4e
......@@ -47,10 +47,10 @@ ReferenceCustomTorsionIxn::ReferenceCustomTorsionIxn(const Lepton::CompiledExpre
expressionSet.registerExpression(this->energyParamDerivExpressions[i]);
thetaIndex = expressionSet.getVariableIndex("theta");
numParameters = parameterNames.size();
for (int i = 0; i < (int) numParameters; i++)
torsionParamIndex.push_back(expressionSet.getVariableIndex(parameterNames[i]));
for (map<string, double>::const_iterator iter = globalParameters.begin(); iter != globalParameters.end(); ++iter)
expressionSet.setVariable(expressionSet.getVariableIndex(iter->first), iter->second);
for (auto& param : parameterNames)
torsionParamIndex.push_back(expressionSet.getVariableIndex(param));
for (auto& param : globalParameters)
expressionSet.setVariable(expressionSet.getVariableIndex(param.first), param.second);
}
/**---------------------------------------------------------------------------------------
......
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