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

Continuing refactoring of the reference platform

parent 0843c5f3
......@@ -121,9 +121,9 @@ static vector<RealVec>& extractForces(ContextImpl& context) {
return *((vector<RealVec>*) data->forces);
}
static RealOpenMM* extractBoxSize(ContextImpl& context) {
static RealVec& extractBoxSize(ContextImpl& context) {
ReferencePlatform::PlatformData* data = reinterpret_cast<ReferencePlatform::PlatformData*>(context.getPlatformData());
return (RealOpenMM*) data->periodicBoxSize;
return *(RealVec*) data->periodicBoxSize;
}
static void findAnglesForCCMA(const System& system, vector<ReferenceCCMAAlgorithm::AngleInfo>& angles) {
......@@ -215,14 +215,14 @@ void ReferenceUpdateStateDataKernel::getForces(ContextImpl& context, std::vector
}
void ReferenceUpdateStateDataKernel::getPeriodicBoxVectors(ContextImpl& context, Vec3& a, Vec3& b, Vec3& c) const {
RealOpenMM* box = extractBoxSize(context);
RealVec& box = extractBoxSize(context);
a = Vec3(box[0], 0, 0);
b = Vec3(0, box[1], 0);
c = Vec3(0, 0, box[2]);
}
void ReferenceUpdateStateDataKernel::setPeriodicBoxVectors(ContextImpl& context, const Vec3& a, const Vec3& b, const Vec3& c) const {
RealOpenMM* box = extractBoxSize(context);
RealVec& box = extractBoxSize(context);
box[0] = (RealOpenMM) a[0];
box[1] = (RealOpenMM) b[1];
box[2] = (RealOpenMM) c[2];
......@@ -668,7 +668,7 @@ double ReferenceCalcNonbondedForceKernel::execute(ContextImpl& context, bool inc
bool ewald = (nonbondedMethod == Ewald);
bool pme = (nonbondedMethod == PME);
if (nonbondedMethod != NoCutoff) {
computeNeighborListVoxelHash(*neighborList, numParticles, posData, exclusions, (periodic || ewald || pme) ? extractBoxSize(context) : NULL, nonbondedCutoff, 0.0);
computeNeighborListVoxelHash(*neighborList, numParticles, posData, exclusions, extractBoxSize(context), periodic || ewald || pme, nonbondedCutoff, 0.0);
clj.setUseCutoff(nonbondedCutoff, *neighborList, rfDielectric);
}
if (periodic || ewald || pme)
......@@ -682,7 +682,7 @@ double ReferenceCalcNonbondedForceKernel::execute(ContextImpl& context, bool inc
ReferenceLJCoulomb14 nonbonded14;
refBondForce.calculateForce(num14, bonded14IndexArray, posData, bonded14ParamArray, forceData, includeEnergy ? &energy : NULL, nonbonded14);
if (periodic || ewald || pme) {
RealOpenMM* boxSize = extractBoxSize(context);
RealVec& boxSize = extractBoxSize(context);
energy += dispersionCoefficient/(boxSize[0]*boxSize[1]*boxSize[2]);
}
return energy;
......@@ -799,7 +799,7 @@ double ReferenceCalcCustomNonbondedForceKernel::execute(ContextImpl& context, bo
ReferenceCustomNonbondedIxn ixn(energyExpression, forceExpression, parameterNames);
bool periodic = (nonbondedMethod == CutoffPeriodic);
if (nonbondedMethod != NoCutoff) {
computeNeighborListVoxelHash(*neighborList, numParticles, posData, exclusions, periodic ? extractBoxSize(context) : NULL, nonbondedCutoff, 0.0);
computeNeighborListVoxelHash(*neighborList, numParticles, posData, exclusions, extractBoxSize(context), periodic, nonbondedCutoff, 0.0);
ixn.setUseCutoff(nonbondedCutoff, *neighborList);
}
if (periodic)
......@@ -885,10 +885,10 @@ void ReferenceCalcGBVIForceKernel::initialize(const System& system, const GBVIFo
double ReferenceCalcGBVIForceKernel::execute(ContextImpl& context, bool includeForces, bool includeEnergy) {
vector<RealVec>& posData = extractPositions(context);
RealOpenMM* bornRadii = new RealOpenMM[context.getSystem().getNumParticles()];
vector<RealOpenMM> bornRadii(context.getSystem().getNumParticles());
if (isPeriodic)
gbvi->getGBVIParameters()->setPeriodic(extractBoxSize(context));
gbvi->computeBornRadii(posData, bornRadii, NULL );
gbvi->computeBornRadii(posData, bornRadii );
if (includeForces) {
vector<RealVec>& forceData = extractForces(context);
gbvi->computeBornForces(bornRadii, posData, &charges[0], forceData);
......@@ -896,7 +896,6 @@ double ReferenceCalcGBVIForceKernel::execute(ContextImpl& context, bool includeF
RealOpenMM energy = 0.0;
if (includeEnergy)
energy = gbvi->computeBornEnergy(bornRadii ,posData, &charges[0]);
delete[] bornRadii;
return static_cast<double>(energy);
}
......@@ -1029,7 +1028,7 @@ double ReferenceCalcCustomGBForceKernel::execute(ContextImpl& context, bool incl
if (periodic)
ixn.setPeriodic(extractBoxSize(context));
if (nonbondedMethod != NoCutoff) {
computeNeighborListVoxelHash(*neighborList, numParticles, posData, exclusions, periodic ? extractBoxSize(context) : NULL, nonbondedCutoff, 0.0);
computeNeighborListVoxelHash(*neighborList, numParticles, posData, exclusions, extractBoxSize(context), periodic, nonbondedCutoff, 0.0);
ixn.setUseCutoff(nonbondedCutoff, *neighborList);
}
map<string, double> globalParameters;
......@@ -1525,7 +1524,7 @@ void ReferenceApplyMonteCarloBarostatKernel::scaleCoordinates(ContextImpl& conte
if (barostat == NULL)
barostat = new ReferenceMonteCarloBarostat(context.getSystem().getNumParticles(), context.getMolecules());
vector<RealVec>& posData = extractPositions(context);
RealOpenMM* boxSize = extractBoxSize(context);
RealVec& boxSize = extractBoxSize(context);
barostat->applyBarostat(posData, boxSize, scale);
}
......
......@@ -88,13 +88,12 @@ ReferencePlatform::PlatformData::PlatformData(int numParticles) : time(0.0), ste
positions = new vector<RealVec>(numParticles);
velocities = new vector<RealVec>(numParticles);
forces = new vector<RealVec>(numParticles);
periodicBoxSize = new RealOpenMM[3];
periodicBoxSize = new RealVec();
}
ReferencePlatform::PlatformData::~PlatformData() {
delete (vector<RealVec>*) positions;
delete (vector<RealVec>*) velocities;
delete (vector<RealVec>*) forces;
RealOpenMM* periodicBoxSize = (RealOpenMM*) this->periodicBoxSize;
delete[] periodicBoxSize;
delete (RealVec*) periodicBoxSize;
}
......@@ -122,7 +122,7 @@ ReferenceCustomGBIxn::~ReferenceCustomGBIxn( ){
--------------------------------------------------------------------------------------- */
void ReferenceCustomGBIxn::setPeriodic( RealOpenMM* boxSize ) {
void ReferenceCustomGBIxn::setPeriodic( RealVec& boxSize ) {
if (cutoff) {
assert(boxSize[0] >= 2.0*cutoffDistance);
......
......@@ -267,7 +267,7 @@ class ReferenceCustomGBIxn {
--------------------------------------------------------------------------------------- */
void setPeriodic( RealOpenMM* boxSize );
void setPeriodic( OpenMM::RealVec& boxSize );
/**---------------------------------------------------------------------------------------
......
......@@ -90,7 +90,7 @@ void ReferenceCustomHbondIxn::setUseCutoff(RealOpenMM distance) {
--------------------------------------------------------------------------------------- */
void ReferenceCustomHbondIxn::setPeriodic(RealOpenMM* boxSize) {
void ReferenceCustomHbondIxn::setPeriodic(RealVec& boxSize) {
assert(cutoff);
assert(boxSize[0] >= 2.0*cutoffDistance);
assert(boxSize[1] >= 2.0*cutoffDistance);
......
......@@ -114,7 +114,7 @@ class ReferenceCustomHbondIxn : public ReferenceBondIxn {
--------------------------------------------------------------------------------------- */
void setPeriodic(RealOpenMM* boxSize);
void setPeriodic(OpenMM::RealVec& boxSize);
/**---------------------------------------------------------------------------------------
......
......@@ -110,7 +110,7 @@ ReferenceCustomNonbondedIxn::~ReferenceCustomNonbondedIxn( ){
--------------------------------------------------------------------------------------- */
int ReferenceCustomNonbondedIxn::setPeriodic( RealOpenMM* boxSize ) {
int ReferenceCustomNonbondedIxn::setPeriodic( RealVec& boxSize ) {
assert(cutoff);
assert(boxSize[0] >= 2.0*cutoffDistance);
......
......@@ -110,7 +110,7 @@ class ReferenceCustomNonbondedIxn {
--------------------------------------------------------------------------------------- */
int setPeriodic( RealOpenMM* boxSize );
int setPeriodic( OpenMM::RealVec& boxSize );
/**---------------------------------------------------------------------------------------
......
......@@ -107,7 +107,7 @@ ReferenceLJCoulombIxn::~ReferenceLJCoulombIxn( ){
--------------------------------------------------------------------------------------- */
int ReferenceLJCoulombIxn::setPeriodic( RealOpenMM* boxSize ) {
int ReferenceLJCoulombIxn::setPeriodic( RealVec& boxSize ) {
assert(cutoff);
assert(boxSize[0] >= 2.0*cutoffDistance);
......
......@@ -117,7 +117,7 @@ class ReferenceLJCoulombIxn : public ReferencePairIxn {
--------------------------------------------------------------------------------------- */
int setPeriodic( RealOpenMM* boxSize );
int setPeriodic( OpenMM::RealVec& boxSize );
/**---------------------------------------------------------------------------------------
......
......@@ -62,7 +62,7 @@ ReferenceMonteCarloBarostat::~ReferenceMonteCarloBarostat( ) {
--------------------------------------------------------------------------------------- */
void ReferenceMonteCarloBarostat::applyBarostat(vector<RealVec>& atomPositions, RealOpenMM* boxSize, RealOpenMM scale) {
void ReferenceMonteCarloBarostat::applyBarostat(vector<RealVec>& atomPositions, const RealVec& boxSize, RealOpenMM scale) {
int numAtoms = savedAtomPositions[0].size();
for (int i = 0; i < numAtoms; i++)
for (int j = 0; j < 3; j++)
......
......@@ -66,7 +66,7 @@ class ReferenceMonteCarloBarostat {
--------------------------------------------------------------------------------------- */
void applyBarostat(std::vector<OpenMM::RealVec>& atomPositions, RealOpenMM* boxSize, RealOpenMM scale);
void applyBarostat(std::vector<OpenMM::RealVec>& atomPositions, const OpenMM::RealVec& boxSize, RealOpenMM scale);
/**---------------------------------------------------------------------------------------
......
......@@ -18,9 +18,9 @@ static double periodicDifference(double val1, double val2, double period) {
}
// squared distance between two points
static double compPairDistanceSquared(const RealVec& pos1, const RealVec& pos2, const RealOpenMM* periodicBoxSize) {
static double compPairDistanceSquared(const RealVec& pos1, const RealVec& pos2, const RealVec& periodicBoxSize, bool usePeriodic) {
double dx, dy, dz;
if (periodicBoxSize == NULL) {
if (!usePeriodic) {
dx = pos2[0] - pos1[0];
dy = pos2[1] - pos1[1];
dz = pos2[2] - pos1[2];
......@@ -40,7 +40,8 @@ void OPENMM_EXPORT computeNeighborListNaive(
int nAtoms,
const AtomLocationList& atomLocations,
const vector<set<int> >& exclusions,
const RealOpenMM* periodicBoxSize,
const RealVec& periodicBoxSize,
bool usePeriodic,
double maxDistance,
double minDistance,
bool reportSymmetricPairs
......@@ -55,7 +56,7 @@ void OPENMM_EXPORT computeNeighborListNaive(
{
for (AtomIndex atomJ = atomI + 1; atomJ < (AtomIndex) nAtoms; ++atomJ)
{
double pairDistanceSquared = compPairDistanceSquared(atomLocations[atomI], atomLocations[atomJ], periodicBoxSize);
double pairDistanceSquared = compPairDistanceSquared(atomLocations[atomI], atomLocations[atomJ], periodicBoxSize, usePeriodic);
if ( (pairDistanceSquared <= maxDistanceSquared) && (pairDistanceSquared >= minDistanceSquared))
if (exclusions[atomI].find(atomJ) == exclusions[atomI].end())
{
......@@ -95,9 +96,9 @@ typedef std::vector< VoxelItem > Voxel;
class VoxelHash
{
public:
VoxelHash(double vsx, double vsy, double vsz, const RealOpenMM* periodicBoxSize) :
voxelSizeX(vsx), voxelSizeY(vsy), voxelSizeZ(vsz), periodicBoxSize(periodicBoxSize) {
if (periodicBoxSize != NULL) {
VoxelHash(double vsx, double vsy, double vsz, const RealVec& periodicBoxSize, bool usePeriodic) :
voxelSizeX(vsx), voxelSizeY(vsy), voxelSizeZ(vsz), periodicBoxSize(periodicBoxSize), usePeriodic(usePeriodic) {
if (usePeriodic) {
nx = (int) floor(periodicBoxSize[0]/voxelSizeX+0.5);
ny = (int) floor(periodicBoxSize[1]/voxelSizeY+0.5);
nz = (int) floor(periodicBoxSize[2]/voxelSizeZ+0.5);
......@@ -115,7 +116,7 @@ public:
VoxelIndex getVoxelIndex(const RealVec& location) const {
double xperiodic, yperiodic, zperiodic;
if (periodicBoxSize == NULL) {
if (!usePeriodic) {
xperiodic = location[0];
yperiodic = location[1];
zperiodic = location[2];
......@@ -162,7 +163,7 @@ public:
int lastx = centerVoxelIndex.x+dIndexX;
int lasty = centerVoxelIndex.y+dIndexY;
int lastz = centerVoxelIndex.z+dIndexZ;
if (periodicBoxSize != NULL) {
if (usePeriodic) {
lastx = min(lastx, centerVoxelIndex.x-dIndexX+nx-1);
lasty = min(lasty, centerVoxelIndex.y-dIndexY+ny-1);
lastz = min(lastz, centerVoxelIndex.z-dIndexZ+nz-1);
......@@ -174,7 +175,7 @@ public:
for (int z = centerVoxelIndex.z - dIndexZ; z <= lastz; ++z)
{
VoxelIndex voxelIndex(x, y, z);
if (periodicBoxSize != NULL) {
if (usePeriodic) {
voxelIndex.x = (x+nx)%nx;
voxelIndex.y = (y+ny)%ny;
voxelIndex.z = (z+nz)%nz;
......@@ -192,7 +193,7 @@ public:
// Ignore exclusions.
if (exclusions[atomI].find(atomJ) != exclusions[atomI].end()) continue;
double dSquared = compPairDistanceSquared(locationI, locationJ, periodicBoxSize);
double dSquared = compPairDistanceSquared(locationI, locationJ, periodicBoxSize, usePeriodic);
if (dSquared > maxDistanceSquared) continue;
if (dSquared < minDistanceSquared) continue;
......@@ -208,7 +209,8 @@ public:
private:
double voxelSizeX, voxelSizeY, voxelSizeZ;
int nx, ny, nz;
const RealOpenMM* periodicBoxSize;
const RealVec& periodicBoxSize;
const bool usePeriodic;
std::map<VoxelIndex, Voxel> voxelMap;
};
......@@ -219,7 +221,8 @@ void OPENMM_EXPORT computeNeighborListVoxelHash(
int nAtoms,
const AtomLocationList& atomLocations,
const vector<set<int> >& exclusions,
const RealOpenMM* periodicBoxSize,
const RealVec& periodicBoxSize,
bool usePeriodic,
double maxDistance,
double minDistance,
bool reportSymmetricPairs
......@@ -228,14 +231,14 @@ void OPENMM_EXPORT computeNeighborListVoxelHash(
neighborList.clear();
double edgeSizeX, edgeSizeY, edgeSizeZ;
if (periodicBoxSize == NULL)
if (!usePeriodic)
edgeSizeX = edgeSizeY = edgeSizeZ = maxDistance; // TODO - adjust this as needed
else {
edgeSizeX = periodicBoxSize[0]/floor(periodicBoxSize[0]/maxDistance);
edgeSizeY = periodicBoxSize[1]/floor(periodicBoxSize[1]/maxDistance);
edgeSizeZ = periodicBoxSize[2]/floor(periodicBoxSize[2]/maxDistance);
}
VoxelHash voxelHash(edgeSizeX, edgeSizeY, edgeSizeZ, periodicBoxSize);
VoxelHash voxelHash(edgeSizeX, edgeSizeY, edgeSizeZ, periodicBoxSize, usePeriodic);
for (AtomIndex atomJ = 0; atomJ < (AtomIndex) nAtoms; ++atomJ) // use "j", because j > i for pairs
{
// 1) Find other atoms that are close to this one
......
......@@ -22,7 +22,8 @@ void OPENMM_EXPORT computeNeighborListNaive(
int nAtoms,
const AtomLocationList& atomLocations,
const std::vector<std::set<int> >& exclusions,
const RealOpenMM* periodicBoxSize,
const RealVec& periodicBoxSize,
bool usePeriodic,
double maxDistance,
double minDistance = 0.0,
bool reportSymmetricPairs = false
......@@ -36,7 +37,8 @@ void OPENMM_EXPORT computeNeighborListVoxelHash(
int nAtoms,
const AtomLocationList& atomLocations,
const std::vector<std::set<int> >& exclusions,
const RealOpenMM* periodicBoxSize,
const RealVec& periodicBoxSize,
bool usePeriodic,
double maxDistance,
double minDistance = 0.0,
bool reportSymmetricPairs = false
......
......@@ -149,7 +149,7 @@ int CpuGBVI::setGBVIParameters( GBVIParameters* gbviParameters ){
#define GBVIDebug 0
int CpuGBVI::computeBornRadii( vector<RealVec>& atomCoordinates, RealOpenMM* bornRadii, RealOpenMM* chain ){
int CpuGBVI::computeBornRadii( vector<RealVec>& atomCoordinates, vector<RealOpenMM>& bornRadii ){
// ---------------------------------------------------------------------------------------
......@@ -438,7 +438,7 @@ RealOpenMM CpuGBVI::Sgb( RealOpenMM t ){
--------------------------------------------------------------------------------------- */
RealOpenMM CpuGBVI::computeBornEnergy( const RealOpenMM* bornRadii, vector<RealVec>& atomCoordinates,
RealOpenMM CpuGBVI::computeBornEnergy( const vector<RealOpenMM>& bornRadii, vector<RealVec>& atomCoordinates,
const RealOpenMM* partialCharges ){
// ---------------------------------------------------------------------------------------
......@@ -462,10 +462,6 @@ RealOpenMM CpuGBVI::computeBornEnergy( const RealOpenMM* bornRadii, vector<RealV
const RealOpenMM* atomicRadii = gbviParameters->getAtomicRadii();
const RealOpenMM* gammaParameters = gbviParameters->getGammaParameters();
if( bornRadii == NULL ){
bornRadii = getBornRadii();
}
#if( GBVIDebug == 1 )
FILE* logFile = stderr;
(void) fprintf( logFile, "\n%s\n", methodName );
......@@ -557,7 +553,7 @@ partialCharges[atomJ]*Sgb( t )/deltaR[ReferenceForce::RIndex];
--------------------------------------------------------------------------------------- */
int CpuGBVI::computeBornForces( const RealOpenMM* bornRadii, vector<RealVec>& atomCoordinates,
int CpuGBVI::computeBornForces( const std::vector<RealOpenMM>& bornRadii, vector<RealVec>& atomCoordinates,
const RealOpenMM* partialCharges, std::vector<OpenMM::RealVec>& inputForces){
// ---------------------------------------------------------------------------------------
......@@ -587,10 +583,6 @@ int CpuGBVI::computeBornForces( const RealOpenMM* bornRadii, vector<RealVec>& at
const RealOpenMM* atomicRadii = gbviParameters->getAtomicRadii();
const RealOpenMM* gammaParameters = gbviParameters->getGammaParameters();
if( bornRadii == NULL ){
bornRadii = getBornRadii();
}
// ---------------------------------------------------------------------------------------
// constants
......@@ -601,8 +593,6 @@ int CpuGBVI::computeBornForces( const RealOpenMM* bornRadii, vector<RealVec>& at
// set energy/forces to zero
const unsigned int arraySzInBytes = sizeof( RealOpenMM )*numberOfAtoms;
RealOpenMM** forces = new RealOpenMM*[numberOfAtoms];
RealOpenMM* block = new RealOpenMM[numberOfAtoms*3];
memset( block, 0, sizeof( RealOpenMM )*numberOfAtoms*3 );
......@@ -612,8 +602,8 @@ int CpuGBVI::computeBornForces( const RealOpenMM* bornRadii, vector<RealVec>& at
blockPtr += 3;
}
RealOpenMM* bornForces = getBornForce();
memset( bornForces, 0, arraySzInBytes );
vector<RealOpenMM>& bornForces = getBornForce();
bornForces.assign(numberOfAtoms, 0.0);
// ---------------------------------------------------------------------------------------
......
......@@ -101,8 +101,7 @@ class CpuGBVI : public CpuImplicitSolvent {
--------------------------------------------------------------------------------------- */
int computeBornRadii( std::vector<OpenMM::RealVec>& atomCoordinates, RealOpenMM* bornRadii,
RealOpenMM* gbviChain = NULL );
int computeBornRadii( std::vector<OpenMM::RealVec>& atomCoordinates, std::vector<RealOpenMM>& bornRadii );
/**---------------------------------------------------------------------------------------
......@@ -232,7 +231,7 @@ class CpuGBVI : public CpuImplicitSolvent {
--------------------------------------------------------------------------------------- */
RealOpenMM computeBornEnergy( const RealOpenMM* bornRadii, std::vector<OpenMM::RealVec>& atomCoordinates,
RealOpenMM computeBornEnergy( const std::vector<RealOpenMM>& bornRadii, std::vector<OpenMM::RealVec>& atomCoordinates,
const RealOpenMM* partialCharges );
/**---------------------------------------------------------------------------------------
......@@ -248,7 +247,7 @@ class CpuGBVI : public CpuImplicitSolvent {
--------------------------------------------------------------------------------------- */
int computeBornForces( const RealOpenMM* bornRadii, std::vector<OpenMM::RealVec>& atomCoordinates,
int computeBornForces( const std::vector<RealOpenMM>& bornRadii, std::vector<OpenMM::RealVec>& atomCoordinates,
const RealOpenMM* partialCharges, std::vector<OpenMM::RealVec>& inputForces );
/**---------------------------------------------------------------------------------------
......
......@@ -90,18 +90,7 @@ CpuImplicitSolvent::~CpuImplicitSolvent( ){
// ---------------------------------------------------------------------------------------
//#ifdef UseGromacsMalloc
// save_free( "_bornForce", __FILE__, __LINE__, _bornForce );
//#else
// delete[] _bornForce;
//#endif
delete _implicitSolventParameters;
delete[] _bornForce;
delete[] _bornRadii;
delete[] _tempBornRadii;
}
/**---------------------------------------------------------------------------------------
......@@ -118,10 +107,6 @@ void CpuImplicitSolvent::_initializeDataMembers( void ){
// ---------------------------------------------------------------------------------------
_bornRadii = NULL;
_tempBornRadii = NULL;
_bornForce = NULL;
_includeAceApproximation = 0;
_forceConversionFactor = (RealOpenMM) 1.0;
......@@ -490,7 +475,7 @@ int CpuImplicitSolvent::incrementForceCallIndex( void ){
--------------------------------------------------------------------------------------- */
RealOpenMM* CpuImplicitSolvent::getBornForce( void ){
vector<RealOpenMM>& CpuImplicitSolvent::getBornForce( void ){
// ---------------------------------------------------------------------------------------
......@@ -498,8 +483,8 @@ RealOpenMM* CpuImplicitSolvent::getBornForce( void ){
// ---------------------------------------------------------------------------------------
if( _bornForce == NULL ){
_bornForce = new RealOpenMM[_implicitSolventParameters->getNumberOfAtoms()];
if( _bornForce.size() == 0 ){
_bornForce.resize(_implicitSolventParameters->getNumberOfAtoms());
}
return _bornForce;
......@@ -514,7 +499,7 @@ RealOpenMM* CpuImplicitSolvent::getBornForce( void ){
--------------------------------------------------------------------------------------- */
RealOpenMM* CpuImplicitSolvent::getBornRadii( void ){
vector<RealOpenMM>& CpuImplicitSolvent::getBornRadii( void ){
// ---------------------------------------------------------------------------------------
......@@ -522,9 +507,8 @@ RealOpenMM* CpuImplicitSolvent::getBornRadii( void ){
// ---------------------------------------------------------------------------------------
if( _bornRadii == NULL ){
_bornRadii = new RealOpenMM[_implicitSolventParameters->getNumberOfAtoms()];
_bornRadii[0] = 0.0;
if( _bornRadii.size() == 0 ){
_bornRadii.resize(_implicitSolventParameters->getNumberOfAtoms(), 0.0);
}
return _bornRadii;
}
......@@ -537,7 +521,7 @@ RealOpenMM* CpuImplicitSolvent::getBornRadii( void ){
--------------------------------------------------------------------------------------- */
const RealOpenMM* CpuImplicitSolvent::getBornRadiiConst( void ) const {
const vector<RealOpenMM>& CpuImplicitSolvent::getBornRadiiConst( void ) const {
// ---------------------------------------------------------------------------------------
......@@ -557,7 +541,7 @@ const RealOpenMM* CpuImplicitSolvent::getBornRadiiConst( void ) const {
--------------------------------------------------------------------------------------- */
RealOpenMM* CpuImplicitSolvent::getBornRadiiTemp( void ){
vector<RealOpenMM>& CpuImplicitSolvent::getBornRadiiTemp( void ){
// ---------------------------------------------------------------------------------------
......@@ -565,8 +549,8 @@ RealOpenMM* CpuImplicitSolvent::getBornRadiiTemp( void ){
// ---------------------------------------------------------------------------------------
if( _tempBornRadii == NULL ){
_tempBornRadii = new RealOpenMM[_implicitSolventParameters->getNumberOfAtoms()];
if( _tempBornRadii.size() == 0 ){
_tempBornRadii.resize(_implicitSolventParameters->getNumberOfAtoms(), 0.0);
}
return _tempBornRadii;
}
......@@ -584,8 +568,7 @@ RealOpenMM* CpuImplicitSolvent::getBornRadiiTemp( void ){
--------------------------------------------------------------------------------------- */
int CpuImplicitSolvent::computeBornRadii( vector<RealVec>& atomCoordinates, RealOpenMM* bornRadii,
RealOpenMM* obcChain ){
int CpuImplicitSolvent::computeBornRadii( vector<RealVec>& atomCoordinates, vector<RealOpenMM>& bornRadii ){
// ---------------------------------------------------------------------------------------
......@@ -647,7 +630,7 @@ int CpuImplicitSolvent::computeImplicitSolventForces( vector<RealVec>& atomCoord
// after first iteration Born radii are updated in force calculation (computeBornEnergyForces())
// unless updateBornRadii is set
RealOpenMM* bornRadii = getBornRadii();
vector<RealOpenMM>& bornRadii = getBornRadii();
if( updateBornRadii || bornRadii[0] < (RealOpenMM) 0.0001 || callId == 1 ){
computeBornRadii( atomCoordinates, bornRadii );
}
......@@ -674,7 +657,7 @@ int CpuImplicitSolvent::computeImplicitSolventForces( vector<RealVec>& atomCoord
--------------------------------------------------------------------------------------- */
int CpuImplicitSolvent::computeBornEnergyForces( RealOpenMM* bornRadii,
int CpuImplicitSolvent::computeBornEnergyForces( vector<RealOpenMM>& bornRadii,
vector<RealVec>& atomCoordinates,
const RealOpenMM* partialCharges,
vector<RealVec>& forces ){
......@@ -708,8 +691,8 @@ int CpuImplicitSolvent::computeBornEnergyForces( RealOpenMM* bornRadii,
--------------------------------------------------------------------------------------- */
int CpuImplicitSolvent::computeAceNonPolarForce( const ImplicitSolventParameters* implicitSolventParameters,
const RealOpenMM* bornRadii, RealOpenMM* energy,
RealOpenMM* forces ) const {
const vector<RealOpenMM>& bornRadii, RealOpenMM* energy,
vector<RealOpenMM>& forces ) const {
// ---------------------------------------------------------------------------------------
......
......@@ -26,6 +26,7 @@
#define __CpuImplicitSolvent_H__
#include "ImplicitSolventParameters.h"
#include <vector>
// ---------------------------------------------------------------------------------------
......@@ -71,12 +72,12 @@ class OPENMM_EXPORT CpuImplicitSolvent {
// work arrays
RealOpenMM* _bornForce;
std::vector<RealOpenMM> _bornForce;
// Born radii and force
RealOpenMM* _bornRadii;
RealOpenMM* _tempBornRadii;
std::vector<RealOpenMM> _bornRadii;
std::vector<RealOpenMM> _tempBornRadii;
// convert units for energy/force
......@@ -107,7 +108,7 @@ class OPENMM_EXPORT CpuImplicitSolvent {
--------------------------------------------------------------------------------------- */
RealOpenMM* getBornForce( void );
std::vector<RealOpenMM>& getBornForce( void );
/**---------------------------------------------------------------------------------------
......@@ -118,7 +119,7 @@ class OPENMM_EXPORT CpuImplicitSolvent {
--------------------------------------------------------------------------------------- */
RealOpenMM* getBornRadiiTemp( void );
std::vector<RealOpenMM>& getBornRadiiTemp( void );
/**---------------------------------------------------------------------------------------
......@@ -330,7 +331,7 @@ class OPENMM_EXPORT CpuImplicitSolvent {
--------------------------------------------------------------------------------------- */
RealOpenMM* getBornRadii( void );
std::vector<RealOpenMM>& getBornRadii( void );
/**---------------------------------------------------------------------------------------
......@@ -340,7 +341,7 @@ class OPENMM_EXPORT CpuImplicitSolvent {
--------------------------------------------------------------------------------------- */
const RealOpenMM* getBornRadiiConst( void ) const;
const std::vector<RealOpenMM>& getBornRadiiConst( void ) const;
/**---------------------------------------------------------------------------------------
......@@ -372,8 +373,7 @@ class OPENMM_EXPORT CpuImplicitSolvent {
--------------------------------------------------------------------------------------- */
virtual int computeBornRadii( std::vector<OpenMM::RealVec>& atomCoordinates, RealOpenMM* bornRadii,
RealOpenMM* obcChain = NULL );
virtual int computeBornRadii( std::vector<OpenMM::RealVec>& atomCoordinates, std::vector<RealOpenMM>& bornRadii );
/**---------------------------------------------------------------------------------------
......@@ -388,7 +388,7 @@ class OPENMM_EXPORT CpuImplicitSolvent {
--------------------------------------------------------------------------------------- */
virtual int computeBornEnergyForces( RealOpenMM* bornRadii, std::vector<OpenMM::RealVec>& atomCoordinates,
virtual int computeBornEnergyForces( std::vector<RealOpenMM>& bornRadii, std::vector<OpenMM::RealVec>& atomCoordinates,
const RealOpenMM* partialCharges, std::vector<OpenMM::RealVec>& forces );
/**---------------------------------------------------------------------------------------
......@@ -405,8 +405,8 @@ class OPENMM_EXPORT CpuImplicitSolvent {
--------------------------------------------------------------------------------------- */
int computeAceNonPolarForce( const ImplicitSolventParameters* implicitSolventParameters,
const RealOpenMM* bornRadii, RealOpenMM* energy,
RealOpenMM* forces ) const;
const std::vector<RealOpenMM>& bornRadii, RealOpenMM* energy,
std::vector<RealOpenMM>& forces ) const;
/**---------------------------------------------------------------------------------------
......
......@@ -71,13 +71,6 @@ CpuObc::~CpuObc( ){
// static const char* methodName = "\nCpuObc::~CpuObc";
// ---------------------------------------------------------------------------------------
//if( _obcParameters != NULL ){
// delete _obcParameters;
//}
delete[] _obcChain;
delete[] _obcChainTemp;
}
/**---------------------------------------------------------------------------------------
......@@ -95,8 +88,6 @@ void CpuObc::_initializeObcDataMembers( void ){
// ---------------------------------------------------------------------------------------
_obcParameters = NULL;
_obcChain = NULL;
_obcChainTemp = NULL;
}
/**---------------------------------------------------------------------------------------
......@@ -149,7 +140,7 @@ int CpuObc::setObcParameters( ObcParameters* obcParameters ){
--------------------------------------------------------------------------------------- */
RealOpenMM* CpuObc::getObcChain( void ){
vector<RealOpenMM>& CpuObc::getObcChain( void ){
// ---------------------------------------------------------------------------------------
......@@ -157,8 +148,8 @@ RealOpenMM* CpuObc::getObcChain( void ){
// ---------------------------------------------------------------------------------------
if( _obcChain == NULL ){
_obcChain = new RealOpenMM[_obcParameters->getNumberOfAtoms()];
if( _obcChain.size() == 0 ){
_obcChain.resize(_obcParameters->getNumberOfAtoms());
}
return _obcChain;
}
......@@ -171,7 +162,7 @@ RealOpenMM* CpuObc::getObcChain( void ){
--------------------------------------------------------------------------------------- */
RealOpenMM* CpuObc::getObcChainConst( void ) const {
const vector<RealOpenMM>& CpuObc::getObcChainConst( void ) const {
// ---------------------------------------------------------------------------------------
......@@ -191,7 +182,7 @@ RealOpenMM* CpuObc::getObcChainConst( void ) const {
--------------------------------------------------------------------------------------- */
RealOpenMM* CpuObc::getObcChainTemp( void ){
vector<RealOpenMM>& CpuObc::getObcChainTemp( void ){
// ---------------------------------------------------------------------------------------
......@@ -199,9 +190,6 @@ RealOpenMM* CpuObc::getObcChainTemp( void ){
// ---------------------------------------------------------------------------------------
if( _obcChainTemp == NULL ){
_obcChainTemp = new RealOpenMM[_obcParameters->getNumberOfAtoms()];
}
return _obcChainTemp;
}
......@@ -219,7 +207,7 @@ RealOpenMM* CpuObc::getObcChainTemp( void ){
--------------------------------------------------------------------------------------- */
int CpuObc::computeBornRadii( vector<RealVec>& atomCoordinates, RealOpenMM* bornRadii, RealOpenMM* obcChain ){
int CpuObc::computeBornRadii( vector<RealVec>& atomCoordinates, vector<RealOpenMM>& bornRadii ){
// ---------------------------------------------------------------------------------------
......@@ -239,9 +227,7 @@ int CpuObc::computeBornRadii( vector<RealVec>& atomCoordinates, RealOpenMM* born
int numberOfAtoms = obcParameters->getNumberOfAtoms();
RealOpenMM* atomicRadii = obcParameters->getAtomicRadii();
const RealOpenMM* scaledRadiusFactor = obcParameters->getScaledRadiusFactors();
if( !obcChain ){
obcChain = getObcChain();
}
vector<RealOpenMM>& obcChain = getObcChain();
RealOpenMM dielectricOffset = obcParameters->getDielectricOffset();
RealOpenMM alphaObc = obcParameters->getAlphaObc();
......@@ -361,7 +347,7 @@ if( logFile ){
--------------------------------------------------------------------------------------- */
int CpuObc::computeBornEnergyForces( RealOpenMM* bornRadii, vector<RealVec>& atomCoordinates,
int CpuObc::computeBornEnergyForces( vector<RealOpenMM>& bornRadii, vector<RealVec>& atomCoordinates,
const RealOpenMM* partialCharges, vector<RealVec>& inputForces ){
// ---------------------------------------------------------------------------------------
......@@ -382,10 +368,6 @@ int CpuObc::computeBornEnergyForces( RealOpenMM* bornRadii, vector<RealVec>& ato
const ObcParameters* obcParameters = getObcParameters();
const int numberOfAtoms = obcParameters->getNumberOfAtoms();
if( bornRadii == NULL ){
bornRadii = getBornRadii();
}
// ---------------------------------------------------------------------------------------
// constants
......@@ -414,7 +396,6 @@ int CpuObc::computeBornEnergyForces( RealOpenMM* bornRadii, vector<RealVec>& ato
// set energy/forces to zero
RealOpenMM obcEnergy = zero;
const unsigned int arraySzInBytes = sizeof( RealOpenMM )*numberOfAtoms;
RealOpenMM** forces = (RealOpenMM**) malloc( sizeof( RealOpenMM* )*numberOfAtoms );
RealOpenMM* block = (RealOpenMM*) malloc( sizeof( RealOpenMM )*numberOfAtoms*3 );
......@@ -425,8 +406,8 @@ int CpuObc::computeBornEnergyForces( RealOpenMM* bornRadii, vector<RealVec>& ato
blockPtr += 3;
}
RealOpenMM* bornForces = getBornForce();
memset( bornForces, 0, arraySzInBytes );
vector<RealOpenMM>& bornForces = getBornForce();
bornForces.assign(numberOfAtoms, 0.0);
// ---------------------------------------------------------------------------------------
......@@ -518,13 +499,13 @@ int CpuObc::computeBornEnergyForces( RealOpenMM* bornRadii, vector<RealVec>& ato
// initialize Born radii & ObcChain temp arrays -- contain values
// used in next iteration
RealOpenMM* bornRadiiTemp = getBornRadiiTemp();
memset( bornRadiiTemp, 0, arraySzInBytes );
vector<RealOpenMM>& bornRadiiTemp = getBornRadiiTemp();
bornRadiiTemp.assign(numberOfAtoms, 0.0);
RealOpenMM* obcChainTemp = getObcChainTemp();
memset( obcChainTemp, 0, arraySzInBytes );
vector<RealOpenMM>& obcChainTemp = getObcChainTemp();
obcChainTemp.assign(numberOfAtoms, 0.0);
RealOpenMM* obcChain = getObcChain();
vector<RealOpenMM>& obcChain = getObcChain();
const RealOpenMM* atomicRadii = obcParameters->getAtomicRadii();
const RealOpenMM alphaObc = obcParameters->getAlphaObc();
......@@ -644,8 +625,8 @@ int CpuObc::computeBornEnergyForces( RealOpenMM* bornRadii, vector<RealVec>& ato
// copy new Born radii and obcChain values into permanent array
memcpy( bornRadii, bornRadiiTemp, arraySzInBytes );
memcpy( obcChain, obcChainTemp, arraySzInBytes );
bornRadii = bornRadiiTemp;
obcChain = obcChainTemp;
free( (char*) block );
free( (char*) forces );
......
......@@ -40,8 +40,8 @@ class CpuObc : public CpuImplicitSolvent {
// arrays containing OBC chain derivative
RealOpenMM* _obcChain;
RealOpenMM* _obcChainTemp;
std::vector<RealOpenMM> _obcChain;
std::vector<RealOpenMM> _obcChainTemp;
// initialize data members (more than
// one constructor, so centralize intialization here)
......@@ -101,8 +101,8 @@ class CpuObc : public CpuImplicitSolvent {
--------------------------------------------------------------------------------------- */
RealOpenMM* getObcChain( void );
RealOpenMM* getObcChainConst( void ) const;
std::vector<RealOpenMM>& getObcChain( void );
const std::vector<RealOpenMM>& getObcChainConst( void ) const;
/**---------------------------------------------------------------------------------------
......@@ -113,7 +113,7 @@ class CpuObc : public CpuImplicitSolvent {
--------------------------------------------------------------------------------------- */
RealOpenMM* getObcChainTemp( void );
std::vector<RealOpenMM>& getObcChainTemp( void );
/**---------------------------------------------------------------------------------------
......@@ -128,8 +128,7 @@ class CpuObc : public CpuImplicitSolvent {
--------------------------------------------------------------------------------------- */
int computeBornRadii( std::vector<OpenMM::RealVec>& atomCoordinates, RealOpenMM* bornRadii,
RealOpenMM* obcChain = NULL );
int computeBornRadii( std::vector<OpenMM::RealVec>& atomCoordinates, std::vector<RealOpenMM>& bornRadii );
/**---------------------------------------------------------------------------------------
......@@ -144,7 +143,7 @@ class CpuObc : public CpuImplicitSolvent {
--------------------------------------------------------------------------------------- */
int computeBornEnergyForces( RealOpenMM* bornRadii, std::vector<OpenMM::RealVec>& atomCoordinates,
int computeBornEnergyForces( std::vector<RealOpenMM>& bornRadii, std::vector<OpenMM::RealVec>& atomCoordinates,
const RealOpenMM* partialCharges, std::vector<OpenMM::RealVec>& forces );
/**---------------------------------------------------------------------------------------
......
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