392 |
|
constERFC = derfc(dampingAlpha*defaultCutoff) |
393 |
|
invRootPi = 0.56418958354775628695d0 |
394 |
|
alphaPi = 2*dampingAlpha*invRootPi |
395 |
< |
|
395 |
> |
|
396 |
|
haveDWAconstants = .true. |
397 |
|
endif |
398 |
|
endif |
448 |
|
real (kind=dp) :: dudx, dudy, dudz |
449 |
|
real (kind=dp) :: scale, sc2, bigR |
450 |
|
real (kind=dp) :: varERFC, varEXP |
451 |
+ |
real (kind=dp) :: limScale |
452 |
|
|
453 |
|
if (.not.allocated(ElectrostaticMap)) then |
454 |
|
call handleError("electrostatic", "no ElectrostaticMap was present before first call of do_electrostatic_pair!") |
472 |
|
!! some variables we'll need independent of electrostatic type: |
473 |
|
|
474 |
|
riji = 1.0d0 / rij |
475 |
< |
|
475 |
> |
|
476 |
|
xhat = d(1) * riji |
477 |
|
yhat = d(2) * riji |
478 |
|
zhat = d(3) * riji |
618 |
|
vpair = vpair + vterm |
619 |
|
epot = epot + sw*vterm |
620 |
|
|
621 |
< |
dudr = -sw*pre11*q_i*q_j * (riji*riji*riji - rcuti2*rcuti) |
621 |
> |
dudr = -sw*pre11*q_i*q_j * (riji*riji-rcuti2)*riji |
622 |
|
|
623 |
|
dudx = dudx + dudr * d(1) |
624 |
|
dudy = dudy + dudr * d(2) |
632 |
|
vpair = vpair + vterm |
633 |
|
epot = epot + sw*vterm |
634 |
|
|
635 |
< |
dudr = -sw*pre11*q_i*q_j * ( riji*(varERFC*riji*riji & |
636 |
< |
+ alphaPi*varEXP) & |
637 |
< |
- rcuti*(constERFC*rcuti2 & |
638 |
< |
+ alphaPi*constEXP) ) |
635 |
> |
dudr = -sw*pre11*q_i*q_j * ( riji*((varERFC*riji*riji & |
636 |
> |
+ alphaPi*varEXP) & |
637 |
> |
- (constERFC*rcuti2 & |
638 |
> |
+ alphaPi*constEXP)) ) |
639 |
|
|
640 |
|
dudx = dudx + dudr * d(1) |
641 |
|
dudy = dudy + dudr * d(2) |