Hi All,
I have a question on the esp printing for ff02 types of force fields.
Here is the piece of code from subroutine esp() in "sander.f" in amber10.
+++++++++++++++++++++++++++
do jn = 1,nesp
e_x = zero
e_y = zero
e_z = zero
e_q = zero
read(dat_unit,'(1x,4e16.0)')esp_qm,xb_esp,yb_esp,zb_esp
x_esp = xb_esp * BOHRS_TO_A
y_esp = yb_esp * BOHRS_TO_A
z_esp = zb_esp * BOHRS_TO_A
do kn = 1,natom
dist = (sqrt((x(1,kn)-x_esp)**2 + &
(x(2,kn)-y_esp)**2 + &
(x(3,kn)-z_esp)**2))
dist3 = dist**3
e_x = e_x - mom_ind(1,kn )*(x(1,kn)-x_esp)/dist3
e_y = e_y - mom_ind(2,kn )*(x(2,kn)-y_esp)/dist3
e_z = e_z - mom_ind(3,kn )*(x(3,kn)-z_esp)/dist3
end do
e_x = e_x * BOHRS_TO_A * INV_AMBER_ELECTROSTATIC
e_y = e_y * BOHRS_TO_A * INV_AMBER_ELECTROSTATIC
e_z = e_z * BOHRS_TO_A * INV_AMBER_ELECTROSTATIC
e_q = e_q * BOHRS_TO_A * INV_AMBER_ELECTROSTATIC
esp_new = e_x + e_y + e_z
write(new_unit, '(1x,4e16.7)')esp_new &
,xb_esp,yb_esp,zb_esp
write(minus_new_unit,'(1x,4e16.7)')esp_qm-esp_new &
,xb_esp,yb_esp,zb_esp
end do
++++++++++++++++++++++++++++
My question is related to the three lines after finishing the potential at
field point (xb_esp,yb_esp,zb_esp) from all induced dipoles:
e_x = e_x * BOHRS_TO_A * INV_AMBER_ELECTROSTATIC
e_y = e_y * BOHRS_TO_A * INV_AMBER_ELECTROSTATIC
e_z = e_z * BOHRS_TO_A * INV_AMBER_ELECTROSTATIC
The code first computes the potential in the electron/Angstrom unit (after
multiplied by INV_AMBER_ELECTROSTATIC), and is the process to dump the
potential out after converting it to the QM unit. So why not multiple the
potential by A_TO_BOHRS but by BOHRS_TO_A?
All the best,
Ray
==========================================
Ray Luo, Ph.D.
Associate Professor
Dept Molecular Biology & Biochemistry
University of California, Irvine, CA 92697
USPS: PO Box 3900 Email: rluo.uci.edu
Phones: (949) 824-9528, 9562
Web: http://rayl0.bio.uci.edu/
==========================================
Received on Wed Oct 15 2008 - 05:12:13 PDT