| 35 |
|
* |
| 36 |
|
* [1] Meineke, et al., J. Comp. Chem. 26, 252-271 (2005). |
| 37 |
|
* [2] Fennell & Gezelter, J. Chem. Phys. 124, 234104 (2006). |
| 38 |
< |
* [3] Sun, Lin & Gezelter, J. Chem. Phys. 128, 24107 (2008). |
| 38 |
> |
* [3] Sun, Lin & Gezelter, J. Chem. Phys. 128, 234107 (2008). |
| 39 |
|
* [4] Kuang & Gezelter, J. Chem. Phys. 133, 164101 (2010). |
| 40 |
|
* [5] Vardeman, Stocker & Gezelter, J. Chem. Theory Comput. 7, 834 (2011). |
| 41 |
|
*/ |
| 263 |
|
// sin(theta) is near 0.0: |
| 264 |
|
|
| 265 |
|
RealType sti2 = 1.0 - cti*cti; |
| 266 |
+ |
RealType proji; |
| 267 |
+ |
Vector3d dcpidr, dcpidu, dspidr, dspidu; |
| 268 |
|
if (fabs(sti2) < 1.0e-12) { |
| 269 |
< |
RealType proji = sqrt(r * 1.0e-12); |
| 270 |
< |
Vector3d dcpidx(1.0 / proji, |
| 271 |
< |
0.0, |
| 272 |
< |
|
| 273 |
< |
dcpidx = 1.0_dp / proji |
| 274 |
< |
dcpidy = 0.0_dp |
| 275 |
< |
dcpidux = xi / proji |
| 276 |
< |
dcpiduy = 0.0_dp |
| 277 |
< |
dspidx = 0.0_dp |
| 278 |
< |
dspidy = 1.0_dp / proji |
| 279 |
< |
dspidux = 0.0_dp |
| 280 |
< |
dspiduy = yi / proji |
| 281 |
< |
else |
| 282 |
< |
proji = sqrt(xi2 + yi2) |
| 283 |
< |
proji3 = proji*proji*proji |
| 284 |
< |
dcpidx = 1.0_dp / proji - xi2 / proji3 |
| 285 |
< |
dcpidy = - xi * yi / proji3 |
| 286 |
< |
dcpidux = xi / proji - (xi2 * xi) / proji3 |
| 287 |
< |
dcpiduy = - (xi * yi2) / proji3 |
| 288 |
< |
dspidx = - xi * yi / proji3 |
| 289 |
< |
dspidy = 1.0_dp / proji - yi2 / proji3 |
| 288 |
< |
dspidux = - (yi * xi2) / proji3 |
| 289 |
< |
dspiduy = yi / proji - (yi2 * yi) / proji3 |
| 290 |
< |
endif |
| 291 |
< |
|
| 292 |
< |
cpi = xi / proji |
| 293 |
< |
dcpidz = 0.0_dp |
| 294 |
< |
dcpiduz = 0.0_dp |
| 269 |
> |
proji = sqrt(r * 1.0e-12); |
| 270 |
> |
dcpidr = Vector3d( 1.0 / proji, 0.0, 0.0); |
| 271 |
> |
dcpidu = Vector3d( xi / proji, 0.0, 0.0); |
| 272 |
> |
dspidr = Vector3d( 0.0, 1.0 / proji, 0.0); |
| 273 |
> |
dspidu = Vector3d (0.0, yi / proji, 0.0); |
| 274 |
> |
} else { |
| 275 |
> |
proji = sqrt(xi2 + yi2); |
| 276 |
> |
RealType proji3 = proji*proji*proji; |
| 277 |
> |
dcpidr = Vector3d(1.0_dp / proji - xi2 / proji3, |
| 278 |
> |
- xi * yi / proji3, |
| 279 |
> |
0.0); |
| 280 |
> |
dcpidu = Vector3d( xi / proji - (xi2 * xi) / proji3, |
| 281 |
> |
- (xi * yi2) / proji3, |
| 282 |
> |
0.0); |
| 283 |
> |
dspidr = Vector3d(- xi * yi / proji3, |
| 284 |
> |
1.0_dp / proji - yi2 / proji3, |
| 285 |
> |
0.0); |
| 286 |
> |
dspidu = Vector3d(- (yi * xi2) / proji3, |
| 287 |
> |
yi / proji - (yi2 * yi) / proji3, |
| 288 |
> |
0.0); |
| 289 |
> |
} |
| 290 |
|
|
| 291 |
< |
spi = yi / proji |
| 292 |
< |
dspidz = 0.0_dp |
| 293 |
< |
dspiduz = 0.0_dp |
| 291 |
> |
cpi = xi / proji; |
| 292 |
> |
dcpidr.z() = 0.0; |
| 293 |
> |
dcpidu.z() = 0.0; |
| 294 |
> |
|
| 295 |
> |
spi = yi / proji; |
| 296 |
> |
dspidr.z() = 0.0; |
| 297 |
> |
dspidu.z() = 0.0; |
| 298 |
|
|
| 299 |
|
|
| 301 |
– |
|
| 302 |
– |
|
| 300 |
|
RealType sigma0 = mixer.sigma0; |
| 301 |
|
RealType dw = mixer.dw; |
| 302 |
|
RealType eps0 = mixer.eps0; |
