ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/mdtools/md_code/Thermo.cpp
(Generate patch)

Comparing trunk/mdtools/md_code/Thermo.cpp (file contents):
Revision 221 by chuckv, Thu Jan 2 20:14:08 2003 UTC vs.
Revision 252 by chuckv, Tue Jan 28 22:16:55 2003 UTC

# Line 1 | Line 1
1   #include <cmath>
2 + #include <iostream>
3 +
4 +
5 + #ifdef IS_MPI
6   #include <mpi++.h>
7 + #endif //is_mpi
8  
9   #include "Thermo.hpp"
10   #include "SRI.hpp"
11   #include "LRI.hpp"
12   #include "Integrator.hpp"
13  
14 + #define BASE_SEED 123456789
15  
16 + Thermo::Thermo( SimInfo* the_entry_plug ) {
17 +  entry_plug = the_entry_plug;
18 +  int baseSeed = BASE_SEED;
19 +  gaussStream = new gaussianSPRNG( baseSeed );
20 + }
21 +
22 + Thermo::~Thermo(){
23 +  delete gaussStream;
24 + }
25 +
26   double Thermo::getKinetic(){
27  
28    const double e_convert = 4.184E-4; // convert kcal/mol -> (amu A^2)/fs^2
# Line 51 | Line 67 | double Thermo::getKinetic(){
67   #ifdef IS_MPI
68    MPI_COMM_WORLD.Allreduce(&kinetic,&kinetic_global,1,MPI_DOUBLE,MPI_SUM);
69    kinetic = kinetic_global;
70 < #endif
70 > #endif //is_mpi
71  
72    kinetic = kinetic * 0.5 / e_convert;
73  
# Line 70 | Line 86 | double Thermo::getPotential(){
86  
87    potential = 0.0;
88    potential_global = 0.0;
89 <  potential += entry_plug->longRange->get_potential();;
89 >  potential += entry_plug->lrPot;
90  
91    // std::cerr << "long range potential: " << potential << "\n";
92    for( el=0; el<nSRI; el++ ){
# Line 82 | Line 98 | double Thermo::getPotential(){
98   #ifdef IS_MPI
99    MPI_COMM_WORLD.Allreduce(&potential,&potential_global,1,MPI_DOUBLE,MPI_SUM);
100    potential = potential_global;
101 < #endif
101 > #endif // is_mpi
102  
103    return potential;
104   }
# Line 148 | Line 164 | void Thermo::velocitize() {
164    ndf = ndfRaw - n_constraints - 3;
165    kebar = kb * temperature * (double)ndf / ( 2.0 * (double)ndfRaw );
166    
167 +  printf("Entered Velocitize\n");
168    for(vr = 0; vr < n_atoms; vr++){
169      
170      // uses equipartition theory to solve for vbar in angstrom/fs
# Line 167 | Line 184 | void Thermo::velocitize() {
184         single processor version of the code, but SPRNG should
185         still be preferred for consistency.
186      */
187 +
188   #ifdef IS_MPI
189   #error "SPRNG random number generator must be used for MPI"
190   #else
191   #warning "Using drand48 for random number generation"
192 < #endif  
192 > #endif  // is_mpi
193 >
194      x = drand48();
195      y = drand48();
196      vx = vbar * sqrt( -2.0 * log(x)) * cos(2 * M_PI * y);
# Line 183 | Line 202 | void Thermo::velocitize() {
202      x = drand48();
203      y = drand48();
204      vz = vbar * sqrt( -2.0 * log(x)) * cos(2 * M_PI * y);
205 < #endif
205 >    printf("Setting new velocities vx: %f\n",vx);
206 > #endif // use_spring
207  
208   #ifdef USE_SPRNG
209 <    vx = vbar * entry_plug->gaussStream->getGaussian();
210 <    vy = vbar * entry_plug->gaussStream->getGaussian();
211 <    vz = vbar * entry_plug->gaussStream->getGaussian();
212 < #endif
209 >    vx = vbar * gaussStream->getGaussian();
210 >    vy = vbar * gaussStream->getGaussian();
211 >    vz = vbar * gaussStream->getGaussian();
212 > #endif // use_spring
213  
214      atoms[vr]->set_vx( vx );
215      atoms[vr]->set_vy( vy );
# Line 238 | Line 258 | void Thermo::velocitize() {
258        if( atoms[i]->isDirectional() ){
259          
260          dAtom = (DirectionalAtom *)atoms[i];
261 +
262 + #ifndef USE_SPRNG
263 +
264   #ifdef IS_MPI
265   #error "SPRNG random number generator must be used for MPI"
266 < #else
266 > #else  // is_mpi
267   #warning "Using drand48 for random number generation"
268 < #endif  
268 > #endif   // is_MPI
269          
270          vbar = sqrt( 2.0 * kebar * dAtom->getIxx() );
271          x = drand48();
# Line 258 | Line 281 | void Thermo::velocitize() {
281          x = drand48();
282          y = drand48();
283          jz = vbar * sqrt( -2.0 * log(x)) * cos(2 * M_PI * y);
284 < #endif
285 < #ifdef USE_SPRNG
284 >
285 > #else //use_sprng
286 >
287          vbar = sqrt( 2.0 * kebar * dAtom->getIxx() );
288 <        jx = vbar * entry_plug->gaussStream->getGaussian();
288 >        jx = vbar * gaussStream->getGaussian();
289  
290          vbar = sqrt( 2.0 * kebar * dAtom->getIyy() );
291 <        jy = vbar * entry_plug->gaussStream->getGaussian();
291 >        jy = vbar * gaussStream->getGaussian();
292  
293          vbar = sqrt( 2.0 * kebar * dAtom->getIzz() );
294 <        jz = vbar * entry_plug->gaussStream->getGaussian();
295 < #endif
294 >        jz = vbar * gaussStream->getGaussian();
295 > #endif //use_sprng
296          
297          dAtom->setJx( jx );
298          dAtom->setJy( jy );

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines