| 553 |
|
real(kind=dp) :: Ri, Ri3, Ri6, Ri7, Ri12, Ri13, R126, R137, prefactor |
| 554 |
|
real(kind=dp) :: chipoalphap2, chioalpha2, ec, epsnot |
| 555 |
|
real(kind=dp) :: drdotudx, drdotudy, drdotudz |
| 556 |
+ |
real(kind=dp) :: drdotudux, drdotuduy, drdotuduz |
| 557 |
|
real(kind=dp) :: ljeps, ljsigma, sigmaratio, sigmaratioi |
| 558 |
|
integer :: ljt1, ljt2, atid1, atid2, gbt1, gbt2 |
| 559 |
|
logical :: gb_first |
| 628 |
|
ljsigma = getSigma(atid1) |
| 629 |
|
ljeps = getEpsilon(atid1) |
| 630 |
|
endif |
| 631 |
< |
|
| 631 |
> |
|
| 632 |
> |
write(*,*) 'd u', dx, dy, dz, ul(1), ul(2), ul(3) |
| 633 |
> |
|
| 634 |
|
rdotu = (dx*ul(1)+dy*ul(2)+dz*ul(3))*ri |
| 635 |
< |
|
| 635 |
> |
|
| 636 |
|
drdotudx = ul(1)*ri-rdotu*dx*ri*ri |
| 637 |
|
drdotudy = ul(2)*ri-rdotu*dy*ri*ri |
| 638 |
|
drdotudz = ul(3)*ri-rdotu*dz*ri*ri |
| 639 |
+ |
drdotudux = drdx |
| 640 |
+ |
drdotuduy = drdy |
| 641 |
+ |
drdotuduz = drdz |
| 642 |
+ |
|
| 643 |
|
|
| 644 |
|
moom = 1.0d0 / gb_mu |
| 645 |
|
mum1 = gb_mu-1 |
| 669 |
|
|
| 670 |
|
mess = 1-rdotu*rdotu*chioalpha2 |
| 671 |
|
sab = 1.0d0/dsqrt(mess) |
| 672 |
+ |
|
| 673 |
+ |
write(*,*) 's', sc, sab, rdotu, chioalpha2 |
| 674 |
|
dsabdct = sc*sab*sab*sab*rdotu*chioalpha2 |
| 675 |
|
|
| 676 |
|
eab = 1-chipoalphap2*rdotu*rdotu |
| 677 |
|
eabf = enot*eab**gb_mu |
| 678 |
+ |
|
| 679 |
+ |
write(*,*) 'e', enot, chipoalphap2, gb_mu, rdotu, eab, mum1 |
| 680 |
+ |
|
| 681 |
|
depmudct = -2*enot*chipoalphap2*gb_mu*rdotu*eab**mum1 |
| 682 |
|
|
| 683 |
|
BigR = (r - sab*sc + sc)/sc |
| 684 |
|
dBigRdx = (drdx -dsabdct*drdotudx)/sc |
| 685 |
|
dBigRdy = (drdy -dsabdct*drdotudy)/sc |
| 686 |
|
dBigRdz = (drdz -dsabdct*drdotudz)/sc |
| 687 |
< |
dBigRdux = (-dsabdct*drdx)/sc |
| 688 |
< |
dBigRduy = (-dsabdct*drdy)/sc |
| 689 |
< |
dBigRduz = (-dsabdct*drdz)/sc |
| 687 |
> |
dBigRdux = (-dsabdct*drdotudux)/sc |
| 688 |
> |
dBigRduy = (-dsabdct*drdotuduy)/sc |
| 689 |
> |
dBigRduz = (-dsabdct*drdotuduz)/sc |
| 690 |
|
|
| 691 |
+ |
write(*,*) 'ds dep', dsabdct, depmudct |
| 692 |
+ |
write(*,*) 'drdotudu', drdotudux, drdotuduy, drdotuduz |
| 693 |
|
depmudx = depmudct*drdotudx |
| 694 |
|
depmudy = depmudct*drdotudy |
| 695 |
|
depmudz = depmudct*drdotudz |
| 696 |
< |
depmudux = depmudct*drdx |
| 697 |
< |
depmuduy = depmudct*drdy |
| 698 |
< |
depmuduz = depmudct*drdz |
| 696 |
> |
depmudux = depmudct*drdotudux |
| 697 |
> |
depmuduy = depmudct*drdotuduy |
| 698 |
> |
depmuduz = depmudct*drdotuduz |
| 699 |
|
|
| 700 |
|
Ri = 1.0d0/BigR |
| 701 |
|
Ri3 = Ri*Ri*Ri |
| 711 |
|
dUdx = prefactor*(eabf*R137*dBigRdx + R126*depmudx) |
| 712 |
|
dUdy = prefactor*(eabf*R137*dBigRdy + R126*depmudy) |
| 713 |
|
dUdz = prefactor*(eabf*R137*dBigRdz + R126*depmudz) |
| 714 |
< |
write(*,*) 'p', prefactor, eabf, r137,dbigrdux, depmudux, r126 |
| 714 |
> |
write(*,*) 'dRdu', dbigrdux, dbigrduy, dbigrduz |
| 715 |
> |
write(*,*) 'dEdu', depmudux, depmuduy, depmuduz |
| 716 |
|
dUdux = prefactor*(eabf*R137*dBigRdux + R126*depmudux) |
| 717 |
|
dUduy = prefactor*(eabf*R137*dBigRduy + R126*depmuduy) |
| 718 |
|
dUduz = prefactor*(eabf*R137*dBigRduz + R126*depmuduz) |
| 752 |
|
t(1,atom1) = t(1,atom1) + ul(2)*dUduz - ul(3)*dUduy |
| 753 |
|
t(2,atom1) = t(2,atom1) + ul(3)*dUdux - ul(1)*dUduz |
| 754 |
|
t(3,atom1) = t(3,atom1) + ul(1)*dUduy - ul(2)*dUdux |
| 755 |
< |
write(*,*) t(1,atom1), t(2,atom1), t(3,atom1) |
| 755 |
> |
write(*,*) 'T1', t(1,atom1), t(2,atom1), t(3,atom1) |
| 756 |
|
else |
| 757 |
|
t(1,atom2) = t(1,atom2) + ul(2)*dUduz - ul(3)*dUduy |
| 758 |
|
t(2,atom2) = t(2,atom2) + ul(3)*dUdux - ul(1)*dUduz |
| 759 |
|
t(3,atom2) = t(3,atom2) + ul(1)*dUduy - ul(2)*dUdux |
| 760 |
|
|
| 761 |
< |
write(*,*) t(1,atom2), t(2,atom2), t(3,atom2) |
| 761 |
> |
write(*,*) 'T2', t(1,atom2), t(2,atom2), t(3,atom2) |
| 762 |
|
endif |
| 763 |
|
|
| 764 |
|
#endif |
