| 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; |
| 156 |
< |
Xiorr *= convertConstant; |
| 157 |
< |
|
| 158 |
< |
|
| 156 |
> |
Xiott *= OOPSEConstant::viscoConvert; |
| 157 |
> |
Xiotr *= OOPSEConstant::viscoConvert; |
| 158 |
> |
Xiorr *= OOPSEConstant::viscoConvert; |
| 159 |
|
|
| 160 |
|
Mat3x3d tmp; |
| 161 |
|
Mat3x3d tmpInv; |
| 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 |
|
|
| 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 |
< |
|
| 326 |
< |
const RealType convertConstant = 6.023; //convert poise.angstrom to amu/fs |
| 327 |
< |
Xitt *= convertConstant; |
| 328 |
< |
Xitr *= convertConstant; |
| 329 |
< |
Xirr *= convertConstant; |
| 325 |
> |
// add the volume correction here: |
| 326 |
> |
Xirr += (6.0 * viscosity * volume) * I; |
| 327 |
|
|
| 328 |
< |
RealType kt = OOPSEConstant::kB * temperature; |
| 328 |
> |
Xitt *= OOPSEConstant::viscoConvert; |
| 329 |
> |
Xitr *= OOPSEConstant::viscoConvert; |
| 330 |
> |
Xirr *= OOPSEConstant::viscoConvert; |
| 331 |
|
|
| 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 |
| 336 |
|
Mat3x3d Dotr; //translation-rotation couplingl diffusion tensor at arbitrary origin O |
| 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; |
| 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; |