calculating the direct Coulomb potential a.k.a. Hartree potential or classical electron-electron repulsion potential
v_J(x) = INT{ rho(x') / |x-x'| } d^3x' where x and x' indicate a vector in 3D.
Basically, its complexity is O(N^2), where N is # of atoms. When molecular grids are used, it becomes O(N_MG^2), where N_MG is # of MG. Here, we aim for linear-scaling, i.e., ~O(N_MG) or ~O(N), by utilizing (a) multicenter decomposition, (b) spherical harmonics expansion, and (c) multipole expansion techniques, (d) screening w.r.t. distance.
(a) Multicenter decomposition: implemented in potential_J via rho_A rho(x) = SUM_A rho_A(x) where rho_A(x) = w_A(x) rho(x) and SUM_A w_A(x) = 1 for all x. complexity: O(N_MG^2) ==> O(N_NG) O(N) O(N_AG)
(b) Spherical harmonics expansion: implemented in potential_J_atom via rho_lm rho_A(x) = SUM_{l=0}^{inf} SUM_{m=-l}^{l} rho_A,lm(r) Y_l^m(th,phi) where rho_A,lm(r) = INT d{Omega} rho_A(x) Y_l^m*(th,phi) complexity: O(N_MG) O(N) O(N_AG) ==> O(N_MG) O(N) O(N_r) O(lmax^2) Practically, it helps a lot because lmax can be quickly trucated. Truncation of lmax can be controlled by absolute accuracy, "-EPS_abs".
(c) Multipole expansion: implemented in potential_J_atom by finding rmax 1/|x-x'| = SUM_{l=0}^{inf}( r_<^l / r_>^{l+1} ) P_l( cos(gamma) ) where cos(gamma) = cos(th)cos(th') - sin(th)sin(th')cos(phi-phi') complexity: O(N_MG) O(N) O(N_r) O(lmax^2) ------- ------ ==> O(N) O(lmax^2) [ O(N_r^2) + O(N_r) + O(N_MG) ] depending on how much MGs are close to A (r of MG < rmax). The O(N_MG) at the end is trivial, just assigning a value rather than performing integrals, so it would be greatly helpful for a large system. Choosing rmax can be controlled by relative accuracy, "-EPS_rel".
(d) Screening with respect to the distance if r_A < r_J (cut-off radius), then just use a monopole (l=0 and m=0).
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| type(param), | intent(in) | :: | p | |||
| type(grid3d), | intent(in) | :: | g | |||
| type(atom_set), | intent(in) | :: | a(:) | |||
| real(kind=long), | intent(in) | :: | rho(:) | |||
| real(kind=long), | intent(out), | optional | :: | q_a(:) | ||
| real(kind=long), | intent(out), | optional | :: | dipole_a(:,:) |