97 |
|
|
98 |
|
end subroutine createMomentMap |
99 |
|
|
100 |
< |
subroutine accumulate_rf(atom1, atom2, rij, u_l, taper) |
100 |
> |
subroutine accumulate_rf(atom1, atom2, rij, eFrame, taper) |
101 |
|
|
102 |
|
integer, intent(in) :: atom1, atom2 |
103 |
|
real (kind = dp), intent(in) :: rij |
104 |
< |
real (kind = dp), dimension(3,nLocal) :: u_l |
104 |
> |
real (kind = dp), dimension(9,nLocal) :: eFrame |
105 |
|
|
106 |
|
integer :: me1, me2 |
107 |
|
real (kind = dp), intent(in) :: taper |
128 |
|
|
129 |
|
#ifdef IS_MPI |
130 |
|
me1 = atid_Row(atom1) |
131 |
< |
ul1(1) = u_l_Row(1,atom1) |
132 |
< |
ul1(2) = u_l_Row(2,atom1) |
133 |
< |
ul1(3) = u_l_Row(3,atom1) |
131 |
> |
ul1(1) = eFrame_Row(7,atom1) |
132 |
> |
ul1(2) = eFrame_Row(8,atom1) |
133 |
> |
ul1(3) = eFrame_Row(9,atom1) |
134 |
|
|
135 |
|
me2 = atid_Col(atom2) |
136 |
< |
ul2(1) = u_l_Col(1,atom2) |
137 |
< |
ul2(2) = u_l_Col(2,atom2) |
138 |
< |
ul2(3) = u_l_Col(3,atom2) |
136 |
> |
ul2(1) = eFrame_Col(7,atom2) |
137 |
> |
ul2(2) = eFrame_Col(8,atom2) |
138 |
> |
ul2(3) = eFrame_Col(9,atom2) |
139 |
|
#else |
140 |
|
me1 = atid(atom1) |
141 |
< |
ul1(1) = u_l(1,atom1) |
142 |
< |
ul1(2) = u_l(2,atom1) |
143 |
< |
ul1(3) = u_l(3,atom1) |
141 |
> |
ul1(1) = eFrame(7,atom1) |
142 |
> |
ul1(2) = eFrame(8,atom1) |
143 |
> |
ul1(3) = eFrame(9,atom1) |
144 |
|
|
145 |
|
me2 = atid(atom2) |
146 |
< |
ul2(1) = u_l(1,atom2) |
147 |
< |
ul2(2) = u_l(2,atom2) |
148 |
< |
ul2(3) = u_l(3,atom2) |
146 |
> |
ul2(1) = eFrame(7,atom2) |
147 |
> |
ul2(2) = eFrame(8,atom2) |
148 |
> |
ul2(3) = eFrame(9,atom2) |
149 |
|
#endif |
150 |
|
|
151 |
|
mu1 = MomentMap(me1)%dipole_moment |
173 |
|
return |
174 |
|
end subroutine accumulate_rf |
175 |
|
|
176 |
< |
subroutine accumulate_self_rf(atom1, mu1, u_l) |
176 |
> |
subroutine accumulate_self_rf(atom1, mu1, eFrame) |
177 |
|
|
178 |
|
integer, intent(in) :: atom1 |
179 |
|
real(kind=dp), intent(in) :: mu1 |
180 |
< |
real(kind=dp), dimension(3,nLocal) :: u_l |
180 |
> |
real(kind=dp), dimension(9,nLocal) :: eFrame |
181 |
|
|
182 |
|
!! should work for both MPI and non-MPI version since this is not pairwise. |
183 |
< |
rf(1,atom1) = rf(1,atom1) + u_l(1,atom1)*mu1 |
184 |
< |
rf(2,atom1) = rf(2,atom1) + u_l(2,atom1)*mu1 |
185 |
< |
rf(3,atom1) = rf(3,atom1) + u_l(3,atom1)*mu1 |
183 |
> |
rf(1,atom1) = rf(1,atom1) + eFrame(7,atom1)*mu1 |
184 |
> |
rf(2,atom1) = rf(2,atom1) + eFrame(8,atom1)*mu1 |
185 |
> |
rf(3,atom1) = rf(3,atom1) + eFrame(9,atom1)*mu1 |
186 |
|
|
187 |
|
return |
188 |
|
end subroutine accumulate_self_rf |
189 |
|
|
190 |
< |
subroutine reaction_field_final(a1, mu1, u_l, rfpot, t, do_pot) |
190 |
> |
subroutine reaction_field_final(a1, mu1, eFrame, rfpot, t, do_pot) |
191 |
|
|
192 |
|
integer, intent(in) :: a1 |
193 |
|
real (kind=dp), intent(in) :: mu1 |
194 |
|
real (kind=dp), intent(inout) :: rfpot |
195 |
|
logical, intent(in) :: do_pot |
196 |
< |
real (kind = dp), dimension(3,nLocal) :: u_l |
196 |
> |
real (kind = dp), dimension(9,nLocal) :: eFrame |
197 |
|
real (kind = dp), dimension(3,nLocal) :: t |
198 |
|
|
199 |
|
integer :: localError |
217 |
|
|
218 |
|
! The torque contribution is dipole cross reaction_field |
219 |
|
|
220 |
< |
t(1,a1) = t(1,a1) + pre*mu1*(u_l(2,a1)*rf(3,a1) - u_l(3,a1)*rf(2,a1)) |
221 |
< |
t(2,a1) = t(2,a1) + pre*mu1*(u_l(3,a1)*rf(1,a1) - u_l(1,a1)*rf(3,a1)) |
222 |
< |
t(3,a1) = t(3,a1) + pre*mu1*(u_l(1,a1)*rf(2,a1) - u_l(2,a1)*rf(1,a1)) |
220 |
> |
t(1,a1) = t(1,a1) + pre*mu1*(eFrame(8,a1)*rf(3,a1) - eFrame(9,a1)*rf(2,a1)) |
221 |
> |
t(2,a1) = t(2,a1) + pre*mu1*(eFrame(9,a1)*rf(1,a1) - eFrame(7,a1)*rf(3,a1)) |
222 |
> |
t(3,a1) = t(3,a1) + pre*mu1*(eFrame(7,a1)*rf(2,a1) - eFrame(8,a1)*rf(1,a1)) |
223 |
|
|
224 |
|
! the potential contribution is -1/2 dipole dot reaction_field |
225 |
|
|
226 |
|
if (do_pot) then |
227 |
|
rfpot = rfpot - 0.5d0 * pre * mu1 * & |
228 |
< |
(rf(1,a1)*u_l(1,a1) + rf(2,a1)*u_l(2,a1) + rf(3,a1)*u_l(3,a1)) |
228 |
> |
(rf(1,a1)*eFrame(7,a1) + rf(2,a1)*eFrame(8,a1) + rf(3,a1)*eFrame(9,a1)) |
229 |
|
endif |
230 |
|
|
231 |
|
return |
232 |
|
end subroutine reaction_field_final |
233 |
|
|
234 |
< |
subroutine rf_correct_forces(atom1, atom2, d, rij, u_l, taper, f, fpair) |
234 |
> |
subroutine rf_correct_forces(atom1, atom2, d, rij, eFrame, taper, f, fpair) |
235 |
|
|
236 |
|
integer, intent(in) :: atom1, atom2 |
237 |
|
real(kind=dp), dimension(3), intent(in) :: d |
238 |
|
real(kind=dp), intent(in) :: rij, taper |
239 |
< |
real( kind = dp ), dimension(3,nLocal) :: u_l |
239 |
> |
real( kind = dp ), dimension(9,nLocal) :: eFrame |
240 |
|
real( kind = dp ), dimension(3,nLocal) :: f |
241 |
|
real( kind = dp ), dimension(3), intent(inout) :: fpair |
242 |
|
|
274 |
|
|
275 |
|
#ifdef IS_MPI |
276 |
|
me1 = atid_Row(atom1) |
277 |
< |
ul1(1) = u_l_Row(1,atom1) |
278 |
< |
ul1(2) = u_l_Row(2,atom1) |
279 |
< |
ul1(3) = u_l_Row(3,atom1) |
277 |
> |
ul1(1) = eFrame_Row(7,atom1) |
278 |
> |
ul1(2) = eFrame_Row(8,atom1) |
279 |
> |
ul1(3) = eFrame_Row(9,atom1) |
280 |
|
|
281 |
|
me2 = atid_Col(atom2) |
282 |
< |
ul2(1) = u_l_Col(1,atom2) |
283 |
< |
ul2(2) = u_l_Col(2,atom2) |
284 |
< |
ul2(3) = u_l_Col(3,atom2) |
282 |
> |
ul2(1) = eFrame_Col(7,atom2) |
283 |
> |
ul2(2) = eFrame_Col(8,atom2) |
284 |
> |
ul2(3) = eFrame_Col(9,atom2) |
285 |
|
#else |
286 |
|
me1 = atid(atom1) |
287 |
< |
ul1(1) = u_l(1,atom1) |
288 |
< |
ul1(2) = u_l(2,atom1) |
289 |
< |
ul1(3) = u_l(3,atom1) |
287 |
> |
ul1(1) = eFrame(7,atom1) |
288 |
> |
ul1(2) = eFrame(8,atom1) |
289 |
> |
ul1(3) = eFrame(9,atom1) |
290 |
|
|
291 |
|
me2 = atid(atom2) |
292 |
< |
ul2(1) = u_l(1,atom2) |
293 |
< |
ul2(2) = u_l(2,atom2) |
294 |
< |
ul2(3) = u_l(3,atom2) |
292 |
> |
ul2(1) = eFrame(7,atom2) |
293 |
> |
ul2(2) = eFrame(8,atom2) |
294 |
> |
ul2(3) = eFrame(9,atom2) |
295 |
|
#endif |
296 |
|
|
297 |
|
mu1 = MomentMap(me1)%dipole_moment |