Skip to content

GitLab

Commit f72062af authored by peastman's avatar peastman
Browse files

Deleted obsolete debugging code

parent 3d69185a
Hide whitespace changes
Inline Side-by-side
Showing with 0 additions and 267 deletions
platforms/reference/include/ReferenceGBVI.h
View file @ f72062af
......@@ -289,26 +289,6 @@ class ReferenceGBVI {
GBVIParameters* gbviParameters,
RealOpenMM* bornRadius, RealOpenMM* switchDeriviative);
/**---------------------------------------------------------------------------------------
Print GB/VI parameters, radii, forces, ...
@param atomCoordinates atomic coordinates
@param partialCharges partial charges
@param bornRadii Born radii (may be empty)
@param bornForces Born forces (may be empty)
@param forces forces (may be empty)
@param idString id string (who is calling)
@param log log file
--------------------------------------------------------------------------------------- */
void printGbvi(const std::vector<OpenMM::RealVec>& atomCoordinates, const std::vector<RealOpenMM>& partialCharges,
const std::vector<RealOpenMM>& bornRadii,
const std::vector<RealOpenMM>& bornForces,
const std::vector<OpenMM::RealVec>& forces,
const std::string& idString, FILE* log);
};
} // namespace OpenMM
......
platforms/reference/include/ReferenceObc.h
View file @ f72062af
......@@ -157,27 +157,6 @@ class ReferenceObc {
RealOpenMM computeBornEnergyForces(const std::vector<OpenMM::RealVec>& atomCoordinates,
const std::vector<RealOpenMM>& partialCharges, std::vector<OpenMM::RealVec>& forces);
/**---------------------------------------------------------------------------------------
Print Obc parameters, radii, forces, ...
@param atomCoordinates atomic coordinates
@param partialCharges partial charges
@param bornRadii Born radii (may be empty)
@param bornForces Born forces (may be empty)
@param forces forces (may be empty)
@param idString id string (who is calling)
@param log log file
--------------------------------------------------------------------------------------- */
void printObc(const std::vector<OpenMM::RealVec>& atomCoordinates,
const std::vector<RealOpenMM>& partialCharges,
const std::vector<RealOpenMM>& bornRadii,
const std::vector<RealOpenMM>& bornForces,
const std::vector<OpenMM::RealVec>& forces,
const std::string& idString, FILE* log);
};
......
platforms/reference/src/SimTKReference/ReferenceGBVI.cpp
View file @ f72062af
......@@ -720,8 +720,6 @@ void ReferenceGBVI::computeBornForces(std::vector<RealVec>& atomCoordinates, con
}
}
//printGbvi(atomCoordinates, partialCharges, bornRadii, bornForces, forces, "GBVI: Post loop2", stderr);
// convert from cal to Joule & apply prefactor tau = (1/diel_solute - 1/diel_solvent)
RealOpenMM conversion = static_cast<RealOpenMM>(gbviParameters->getTau());
......@@ -733,93 +731,6 @@ void ReferenceGBVI::computeBornForces(std::vector<RealVec>& atomCoordinates, con
}
/**---------------------------------------------------------------------------------------
Print GB/VI parameters, radii, forces, ...
@param atomCoordinates atomic coordinates
@param partialCharges partial charges
@param bornRadii Born radii (may be empty)
@param bornForces Born forces (may be empty)
@param forces forces (may be empty)
@param idString id string (who is calling)
@param log log file
--------------------------------------------------------------------------------------- */
void ReferenceGBVI::printGbvi(const std::vector<OpenMM::RealVec>& atomCoordinates, const vector<RealOpenMM>& partialCharges,
const vector<RealOpenMM>& bornRadii,
const vector<RealOpenMM>& bornForces,
const std::vector<OpenMM::RealVec>& forces,
const std::string& idString, FILE* log) {
// ---------------------------------------------------------------------------------------
const GBVIParameters* gbviParameters = getGBVIParameters();
const int numberOfAtoms = gbviParameters->getNumberOfAtoms();
const vector<RealOpenMM>& atomicRadii = gbviParameters->getAtomicRadii();
const vector<RealOpenMM>& gammaParameters = gbviParameters->getGammaParameters();
// ---------------------------------------------------------------------------------------
// constants
const RealOpenMM preFactor = 2.0*gbviParameters->getElectricConstant();
// ---------------------------------------------------------------------------------------
const vector<RealOpenMM>& scaledRadii = gbviParameters->getScaledRadii();
const vector<RealOpenMM>& switchDeriviative = getSwitchDeriviative();
RealOpenMM tau = static_cast<RealOpenMM>(gbviParameters->getTau());
int useComparisonFormat = 1;
(void) fprintf(log, "Reference Gbvi %s atoms=%d\n", idString.c_str(), numberOfAtoms);
(void) fprintf(log, " tau %15.7e\n", tau);
(void) fprintf(log, " scaleMethod %d (QuinticEnum=%d)\n",
_gbviParameters->getBornRadiusScalingMethod(), GBVIParameters::QuinticSpline);
(void) fprintf(log, " preFactor %15.7e)\n", preFactor);
if (useComparisonFormat) {
(void) fprintf(log, " br bF swd r scR tau*gamma q)\n");
for (unsigned int atomI = 0; atomI < static_cast<unsigned int>(numberOfAtoms); atomI++) {
(void) fprintf(log, "%6d ", atomI);
if (bornRadii.size() > atomI) {
(void) fprintf(log, "%15.7e ", bornRadii[atomI]);
}
if (bornForces.size() > atomI) {
(void) fprintf(log, "%15.7e ", tau*bornForces[atomI]);
}
(void) fprintf(log, " %15.7e %15.7e %15.7e %15.7e %15.7e",
switchDeriviative[atomI],
atomicRadii[atomI],
scaledRadii[atomI],
tau*gammaParameters[atomI],
partialCharges[atomI]);
(void) fprintf(log, "\n");
}
} else {
for (unsigned int atomI = 0; atomI < static_cast<unsigned int>(numberOfAtoms); atomI++) {
(void) fprintf(log, "%6d r=%15.7e rSc=%15.7e swd=%15.7e tau*gam=%15.7e q=%15.7e", atomI,
atomicRadii[atomI],
scaledRadii[atomI],
switchDeriviative[atomI],
tau*gammaParameters[atomI],
partialCharges[atomI]);
if (bornRadii.size() > atomI) {
(void) fprintf(log, " bR=%15.7e", bornRadii[atomI]);
}
if (bornForces.size() > atomI) {
(void) fprintf(log, " tau*bF=%15.7e", tau*bornForces[atomI]);
}
(void) fprintf(log, "\n");
}
}
return;
}
/**---------------------------------------------------------------------------------------
Use double precision
......
platforms/reference/src/SimTKReference/ReferenceObc.cpp
View file @ f72062af
......@@ -494,90 +494,5 @@ RealOpenMM ReferenceObc::computeBornEnergyForces(const vector<RealVec>& atomCoor
}
//printObc(atomCoordinates, partialCharges, bornRadii, bornForces, inputForces, "Obc Post loop2", stderr);
return obcEnergy;
}
/**---------------------------------------------------------------------------------------
Print Obc parameters, radii, forces, ...
@param atomCoordinates atomic coordinates
@param partialCharges partial charges
@param bornRadii Born radii (may be empty)
@param bornForces Born forces (may be empty)
@param forces forces (may be empty)
@param idString id string (who is calling)
@param log log file
--------------------------------------------------------------------------------------- */
void ReferenceObc::printObc(const std::vector<OpenMM::RealVec>& atomCoordinates,
const vector<RealOpenMM>& partialCharges,
const vector<RealOpenMM>& bornRadii,
const vector<RealOpenMM>& bornForces,
const std::vector<OpenMM::RealVec>& forces,
const std::string& idString, FILE* log) {
// ---------------------------------------------------------------------------------------
const ObcParameters* obcParameters = getObcParameters();
const int numberOfAtoms = obcParameters->getNumberOfAtoms();
const vector<RealOpenMM>& atomicRadii = obcParameters->getAtomicRadii();
const RealOpenMM preFactor = 2.0*obcParameters->getElectricConstant();
const vector<RealOpenMM>& obcChain = getObcChain();
const vector<RealOpenMM>& scaledRadiusFactor = obcParameters->getScaledRadiusFactors();
const RealOpenMM alphaObc = obcParameters->getAlphaObc();
const RealOpenMM betaObc = obcParameters->getBetaObc();
const RealOpenMM gammaObc = obcParameters->getGammaObc();
const int comparisonFormat = 1;
// ---------------------------------------------------------------------------------------
(void) fprintf(log, "Reference Obc %s atoms=%d\n", idString.c_str(), numberOfAtoms);
if (comparisonFormat) {
(void) fprintf(log, "Reference Obc %s atoms=%d Chain/Radii/Force\n", idString.c_str(), numberOfAtoms);
for (unsigned int atomI = 0; atomI < static_cast<unsigned int>(numberOfAtoms); atomI++) {
(void) fprintf(log, "%6d ", atomI);
if (obcChain.size() > atomI) {
(void) fprintf(log, " %15.7e", obcChain[atomI]);
}
if (bornRadii.size() > atomI) {
(void) fprintf(log, " %15.7e", bornRadii[atomI]);
}
if (bornForces.size() > atomI) {
(void) fprintf(log, " %15.7e", bornForces[atomI]);
}
(void) fprintf(log, " %15.7e %6.3f", atomicRadii[atomI], partialCharges[atomI]);
(void) fprintf(log, "\n");
}
} else {
(void) fprintf(log, "Reference Obc %s atoms=%d\n", idString.c_str(), numberOfAtoms);
(void) fprintf(log, " preFactor %15.7e\n", preFactor);
(void) fprintf(log, " alpha %15.7e\n", alphaObc);
(void) fprintf(log, " beta %15.7e\n", betaObc);
(void) fprintf(log, " gamma %15.7e\n", gammaObc);
for (unsigned int atomI = 0; atomI < static_cast<unsigned int>(numberOfAtoms); atomI++) {
(void) fprintf(log, "%6d r=%15.7e q=%6.3f", atomI,
atomicRadii[atomI], partialCharges[atomI]);
if (obcChain.size() > atomI) {
(void) fprintf(log, " bChn=%15.7e", obcChain[atomI]);
}
if (bornRadii.size() > atomI) {
(void) fprintf(log, " bR=%15.7e", bornRadii[atomI]);
}
if (bornForces.size() > atomI) {
(void) fprintf(log, " bF=%15.7e", bornForces[atomI]);
}
(void) fprintf(log, "\n");
}
}
return;
}
platforms/reference/tests/TestReferenceCustomGBForce.cpp
View file @ f72062af
......@@ -674,10 +674,6 @@ static void findScaledRadii(GBVIForce& gbviForce, std::vector<double> & scaledRa
// load 1-2 atom pairs along w/ bond distance using HarmonicBondForce & constraints
// numberOfBonds < 1, indicating they were not set by the user
if (numberOfBonds < 1 && numberOfParticles > 1) {
(void) fprintf(stderr, "Warning: no covalent bonds set for GB/VI force!\n");
}
std::vector< std::vector<int> > bondIndices;
bondIndices.resize(numberOfBonds);
......@@ -742,9 +738,6 @@ static void findScaledRadii(GBVIForce& gbviForce, std::vector<double> & scaledRa
gbviForce.getParticleParameters(j, charge, radiusJ, gamma);
if ( bonded12[j].size() == 0) {
if (numberOfParticles > 1) {
(void) fprintf(stderr, "Warning GBVIForceImpl::findScaledRadii atom %d has no covalent bonds; using atomic radius=%.3f.\n", j, radiusJ);
}
scaledRadiusJ = radiusJ;
// errors++;
}
......@@ -783,7 +776,6 @@ static void findScaledRadii(GBVIForce& gbviForce, std::vector<double> & scaledRa
scaledRadiusJ = 0.0;
}
}
//(void) fprintf(stderr, "scaledRadii %d %12.4f\n", j, scaledRadiusJ);
scaledRadii[j] = scaledRadiusJ;
}
......@@ -793,26 +785,6 @@ static void findScaledRadii(GBVIForce& gbviForce, std::vector<double> & scaledRa
if (errors) {
throw OpenMMException("GBVIForceImpl::findScaledRadii errors -- aborting");
}
#if GBVIDebug
(void) fprintf(stderr, " R q gamma scaled radii no. bnds\n");
double totalQ = 0.0;
for (int i = 0; i < (int) scaledRadii.size(); i++) {
double charge;
double gamma;
double radiusI;
gbviForce.getParticleParameters(i, charge, radiusI, gamma);
totalQ += charge;
(void) fprintf(stderr, "%4d %14.5e %14.5e %14.5e %14.5e %d\n", i, radiusI, charge, gamma, scaledRadii[i], (int) bonded12[i].size());
}
(void) fprintf(stderr, "Total charge=%e\n", totalQ);
(void) fflush(stderr);
#endif
#undef GBVIDebug
}
// load parameters from gbviForce to customGbviForce
......
platforms/reference/tests/TestReferenceGBVIForce.cpp
View file @ f72062af
......@@ -53,7 +53,6 @@ using namespace std;
const double TOL = 1e-5;
void testSingleParticle() {
const int log = 0;
ReferencePlatform platform;
System system;
system.addParticle(2.0);
......@@ -85,10 +84,6 @@ void testSingleParticle() {
double expectedE = (bornEnergy+nonpolarEnergy);
double obtainedE = state.getPotentialEnergy();
double diff = fabs((obtainedE - expectedE)/expectedE);
if (log) {
(void) fprintf(stderr, "testSingleParticle expected=%14.6e obtained=%14.6e diff=%14.6e breakdown:[%14.6e %14.6e]\n",
expectedE, obtainedE, diff, bornEnergy, nonpolarEnergy);
}
ASSERT_EQUAL_TOL((bornEnergy+nonpolarEnergy), state.getPotentialEnergy(), 0.01);
}
......@@ -97,7 +92,6 @@ void testEnergyEthane(int applyBornRadiiScaling) {
ReferencePlatform platform;
const int numParticles = 8;
const int log = 0;
System system;
LangevinIntegrator integrator(0, 0.1, 0.01);
......@@ -139,9 +133,6 @@ void testEnergyEthane(int applyBornRadiiScaling) {
NonbondedForce* nonbonded = new NonbondedForce();
nonbonded->setNonbondedMethod(NonbondedForce::NoCutoff);
if (log) {
(void) fprintf(stderr, "Applying GB/VI\n");
}
GBVIForce* forceField = new GBVIForce();
if (applyBornRadiiScaling) {
forceField->setBornRadiusScalingMethod(GBVIForce::QuinticSpline);
......@@ -212,9 +203,6 @@ void testEnergyEthane(int applyBornRadiiScaling) {
context.setPositions(positions);
State state = context.getState(State::Forces | State::Energy);
if (log) {
(void) fprintf(stderr, "Energy %.4e\n", state.getPotentialEnergy());
}
// Take a small step in the direction of the energy gradient.
......@@ -222,18 +210,12 @@ void testEnergyEthane(int applyBornRadiiScaling) {
double forceSum[3] = { 0.0, 0.0, 0.0 };
for (int i = 0; i < numParticles; ++i) {
Vec3 f = state.getForces()[i];
if (log) {
(void) fprintf(stderr, "F %d [%14.6e %14.6e %14.6e]\n", i, f[0], f[1], f[2]);
}
norm += f[0]*f[0] + f[1]*f[1] + f[2]*f[2];
forceSum[0] += f[0];
forceSum[1] += f[1];
forceSum[2] += f[2];
}
norm = std::sqrt(norm);
if (log) {
(void) fprintf(stderr, "Fsum [%14.6e %14.6e %14.6e] norm=%14.6e\n", forceSum[0], forceSum[1], forceSum[2], norm);
}
const double delta = 1e-4;
double step = delta/norm;
......@@ -246,12 +228,6 @@ void testEnergyEthane(int applyBornRadiiScaling) {
State state2 = context.getState(State::Energy);
if (log) {
double deltaE = fabs(state.getPotentialEnergy() - state2.getPotentialEnergy())/delta;
double diff = (deltaE - norm)/norm;
(void) fprintf(stderr, "Energies %.8e %.8e deltaE=%14.7e %14.7e diff=%14.7e\n", state.getPotentialEnergy(), state2.getPotentialEnergy(), deltaE, norm, diff);
}
// See whether the potential energy changed by the expected amount.
ASSERT_EQUAL_TOL(norm, (state2.getPotentialEnergy()-state.getPotentialEnergy())/delta, 0.01)
......
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