--- trunk/OOPSE/libmdtools/ForceFields.cpp 2003/04/07 21:16:35 472 +++ trunk/OOPSE/libmdtools/ForceFields.cpp 2003/10/29 20:41:39 843 @@ -1,5 +1,10 @@ -#include +#include +using namespace std; + + +#include + #ifdef IS_MPI #include #endif // is_mpi @@ -21,16 +26,8 @@ void ForceFields::calcRcut( void ){ //calc rCut and rList - entry_plug->rCut = 2.5 * bigSigma; - if(entry_plug->rCut > (entry_plug->box_x / 2.0)) - entry_plug->rCut = entry_plug->box_x / 2.0; - if(entry_plug->rCut > (entry_plug->box_y / 2.0)) - entry_plug->rCut = entry_plug->box_y / 2.0; - if(entry_plug->rCut > (entry_plug->box_z / 2.0)) - entry_plug->rCut = entry_plug->box_z / 2.0; - - entry_plug->rList = entry_plug->rCut + 1.0; - + entry_plug->setDefaultRcut( 2.5 * bigSigma ); + } void ForceFields::doForces( int calcPot, int calcStress ){ @@ -39,72 +36,52 @@ void ForceFields::doForces( int calcPot, int calcStres double* frc; double* pos; double* trq; - double* tau; double* A; double* u_l;; - DirectionalAtom* dAtom; + SimState* config; - double ut[3]; - - //u_l = new double[entry_plug->n_atoms]; - short int passedCalcPot = (short int)calcPot; short int passedCalcStress = (short int)calcStress; - // forces are zeroed here, before any are acumulated. + // forces are zeroed here, before any are accumulated. // NOTE: do not rezero the forces in Fortran. - for(i=0; in_atoms; i++){ - entry_plug->atoms[i]->zeroForces(); - -// if( entry_plug->atoms[i]->isDirectional() ){ -// dAtom = (DirectionalAtom *)entry_plug->atoms[i]; -// dAtom->getU(ut); - - -// if(dAtom->getIndex()== 1){ -// std::cerr << "atom 2's u_l = " << ut[0] << ", " << ut[1] -// << ", " << ut[2] << "\n"; -// } -// } - + entry_plug->atoms[i]->zeroForces(); } for(i=0; in_mol; i++ ){ entry_plug->molecules[i].calcForces(); } - frc = Atom::getFrcArray(); - pos = Atom::getPosArray(); - trq = Atom::getTrqArray(); - A = Atom::getAmatArray(); - u_l = Atom::getUlArray(); - tau = entry_plug->tau; - + config = entry_plug->getConfiguration(); + frc = config->getFrcArray(); + pos = config->getPosArray(); + trq = config->getTrqArray(); + A = config->getAmatArray(); + u_l = config->getUlArray(); + isError = 0; entry_plug->lrPot = 0.0; - tau[0] = 0.0; - tau[1] = 0.0; - tau[2] = 0.0; - + for (i=0; i<9; i++) { + entry_plug->tau[i] = 0.0; + } - + fortranForceLoop( pos, A, u_l, frc, trq, - tau, + entry_plug->tau, &(entry_plug->lrPot), &passedCalcPot, &passedCalcStress, &isError ); - // delete[] u_l; if( isError ){ sprintf( painCave.errMsg, @@ -118,6 +95,7 @@ void ForceFields::doForces( int calcPot, int calcStres "returned from the force calculation.\n" ); MPIcheckPoint(); #endif // is_mpi + }