ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/applications/hydrodynamics/ApproximationModel.cpp
(Generate patch)

Comparing trunk/src/applications/hydrodynamics/ApproximationModel.cpp (file contents):
Revision 981 by gezelter, Mon Jun 5 18:24:45 2006 UTC vs.
Revision 1322 by gezelter, Mon Nov 24 22:01:56 2008 UTC

# Line 72 | Line 72 | namespace oopse {
72    bool ApproximationModel::calcHydroProps(Shape* shape, RealType viscosity, RealType temperature) {
73      
74      bool ret = true;
75 <    HydroProp* cr;
76 <    HydroProp* cd;
75 >    HydroProp* cr = new HydroProp();
76 >    HydroProp* cd = new HydroProp();
77      calcHydroPropsAtCR(beads_, viscosity, temperature, cr);
78 <    //calcHydroPropsAtCD(beads_, viscosity, temperature, cd);
78 >    calcHydroPropsAtCD(beads_, viscosity, temperature, cd);
79      setCR(cr);
80      setCD(cd);
81 –    
81      return true;    
82    }
83    
# Line 146 | Line 145 | namespace oopse {
145          
146          Xiott += Cij;
147          Xiotr += U[i] * Cij;
148 <        //Xiorr += -U[i] * Cij * U[j] + (6 * viscosity * volume) * I;    
148 >        // uncorrected here.  Volume correction is added after we assemble Xiorr
149          Xiorr += -U[i] * Cij * U[j];
150        }
151      }
152 +
153 +    // add the volume correction
154 +    Xiorr += (6.0 * viscosity * volume) * I;    
155      
156 <    const RealType convertConstant = 6.023; //convert poise.angstrom to amu/fs
157 <    Xiott *= convertConstant;
158 <    Xiotr *= convertConstant;
157 <    Xiorr *= convertConstant;
158 <    
159 <    
156 >    Xiott *= OOPSEConstant::viscoConvert;
157 >    Xiotr *= OOPSEConstant::viscoConvert;
158 >    Xiorr *= OOPSEConstant::viscoConvert;
159      
160      Mat3x3d tmp;
161      Mat3x3d tmpInv;
# Line 204 | Line 203 | namespace oopse {
203      Dr6x6.getSubMatrix(0, 3, Drrt);
204      Dr6x6.getSubMatrix(3, 0, Drtr);
205      Dr6x6.getSubMatrix(3, 3, Drrr);
206 <    RealType kt = OOPSEConstant::kB * temperature ;
206 >    RealType kt = OOPSEConstant::kb * temperature ; // in kcal mol^-1
207      Drtt *= kt;
208      Drrt *= kt;
209      Drtr *= kt;
210      Drrr *= kt;
211 <    Xirtt *= OOPSEConstant::kb * temperature;
212 <    Xirtr *= OOPSEConstant::kb * temperature;
213 <    Xirrr *= OOPSEConstant::kb * temperature;
211 >    //Xirtt *= OOPSEConstant::kb * temperature;
212 >    //Xirtr *= OOPSEConstant::kb * temperature;
213 >    //Xirrr *= OOPSEConstant::kb * temperature;
214      
215      Mat6x6d Xi, D;
216  
# Line 255 | Line 254 | namespace oopse {
254      return true;
255   }
256    
257 <  bool ApproximationModel::calcHydroPropsAtCD(std::vector<BeadParam>& beads, RealType viscosity, RealType temperature, HydroProp* cr) {
257 >  bool ApproximationModel::calcHydroPropsAtCD(std::vector<BeadParam>& beads, RealType viscosity, RealType temperature, HydroProp* cd) {
258      
259      int nbeads = beads.size();
260      DynamicRectMatrix<RealType> B(3*nbeads, 3*nbeads);
# Line 319 | Line 318 | namespace oopse {
318              
319          Xitt += Cij;
320          Xitr += U[i] * Cij;
321 <            //Xirr += -U[i] * Cij * U[j] + (6 * viscosity * volume) * I;            
321 >        // uncorrected here.  Volume correction is added after we assemble Xiorr
322          Xirr += -U[i] * Cij * U[j];
323        }
324      }
325 +    // add the volume correction here:
326 +    Xirr += (6.0 * viscosity * volume) * I;    
327      
328 <    const RealType convertConstant = 6.023; //convert poise.angstrom to amu/fs
329 <    Xitt *= convertConstant;
330 <    Xitr *= convertConstant;
330 <    Xirr *= convertConstant;
328 >    Xitt *= OOPSEConstant::viscoConvert;
329 >    Xitr *= OOPSEConstant::viscoConvert;
330 >    Xirr *= OOPSEConstant::viscoConvert;
331      
332 <    RealType kt = OOPSEConstant::kB * temperature;
332 >    RealType kt = OOPSEConstant::kb * temperature; // in kcal mol^-1
333      
334      Mat3x3d Dott; //translational diffusion tensor at arbitrary origin O
335      Mat3x3d Dorr; //rotational diffusion tensor at arbitrary origin O
# Line 407 | Line 407 | namespace oopse {
407  
408      Mat6x6d Xi, D;
409  
410 <    cr->setCOR(rod);
410 >    cd->setCOR(rod);
411  
412 <    cr->setXi(Xid);
412 >    cd->setXi(Xid);
413  
414      D.setSubMatrix(0, 0, Ddtt);
415      D.setSubMatrix(0, 3, Ddtr);
416      D.setSubMatrix(3, 0, Ddtr);
417      D.setSubMatrix(3, 3, Ddrr);
418  
419 <    cr->setD(D);
419 >    cd->setD(D);
420  
421      std::cout << "viscosity = " << viscosity << std::endl;
422      std::cout << "temperature = " << temperature << std::endl;
423      std::cout << "center of diffusion :" << std::endl;
424      std::cout << rod << std::endl;
425      std::cout << "diffusion tensor at center of diffusion " << std::endl;
426 <    std::cout << "translation(A^2/fs) :" << std::endl;
426 >    std::cout << "translation(A^2 / fs) :" << std::endl;
427      std::cout << Ddtt << std::endl;
428 <    std::cout << "translation-rotation(A^3/fs):" << std::endl;
428 >    std::cout << "translation-rotation(A / fs):" << std::endl;
429      std::cout << Ddtr << std::endl;
430 <    std::cout << "rotation(A^4/fs):" << std::endl;
430 >    std::cout << "rotation(fs^-1):" << std::endl;
431      std::cout << Ddrr << std::endl;
432  
433      std::cout << "resistance tensor at center of diffusion " << std::endl;
# Line 443 | Line 443 | namespace oopse {
443      Xid.getSubMatrix(3, 3, Xidrr);
444  
445      std::cout << Xidtt << std::endl;
446 <    std::cout << "rotation-translation (kcal*fs*mol^-1*Ang^-3):" << std::endl;
446 >    std::cout << "rotation-translation (kcal*fs*mol^-1*Ang^-1):" << std::endl;
447      std::cout << Xidrt << std::endl;
448 <    std::cout << "translation-rotation(kcal*fs*mol^-1*Ang^-3):" << std::endl;
448 >    std::cout << "translation-rotation(kcal*fs*mol^-1*Ang^-1):" << std::endl;
449      std::cout << Xidtr << std::endl;
450 <    std::cout << "rotation(kcal*fs*mol^-1*Ang^-4):" << std::endl;
450 >    std::cout << "rotation(kcal*fs*mol^-1):" << std::endl;
451      std::cout << Xidrr << std::endl;
452  
453      return true;

Diff Legend

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