40#include "./base/base_uses.f90"
44 LOGICAL,
PRIVATE,
PARAMETER :: debug_this_module = .true.
45 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'cp_ddapc_types'
54 REAL(kind=
dp) :: c0 = 0.0_dp
55 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: ami => null()
56 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: md => null()
57 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: mr => null()
58 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: mt => null()
59 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: ms => null()
60 REAL(kind=
dp),
POINTER,
DIMENSION(:, :) :: gfunc => null()
61 REAL(kind=
dp),
POINTER,
DIMENSION(:) :: w => null()
66 LOGICAL :: do_decoupling = .false.
67 LOGICAL :: do_qmmm_periodic_decpl = .false.
68 LOGICAL :: do_solvation = .false.
69 LOGICAL :: do_property = .false.
70 LOGICAL :: do_restraint = .false.
72 TYPE(
pw_pool_type),
POINTER :: pw_pool_qm => null(), pw_pool_mm => null()
73 TYPE(
pw_grid_type),
POINTER :: pw_grid_qm => null(), pw_grid_mm => null()
97 particle_set, radii, cell, super_cell, rho_tot_g, gcut, iw2, Vol, &
103 REAL(kind=
dp),
DIMENSION(:),
POINTER :: radii
104 TYPE(
cell_type),
POINTER :: cell, super_cell
106 REAL(kind=
dp),
INTENT(IN) :: gcut
107 INTEGER,
INTENT(IN) :: iw2
108 REAL(kind=
dp),
INTENT(IN) :: vol
111 CHARACTER(len=*),
PARAMETER :: routinen =
'cp_ddapc_create'
116 CALL timeset(routinen, handle)
117 NULLIFY (cp_ddapc_env%AmI, &
122 cp_ddapc_env%gfunc, &
125 CALL ddapc_eval_gfunc(cp_ddapc_env%gfunc, cp_ddapc_env%w, gcut, rho_tot_g, radii)
128 cp_ddapc_env%gfunc, &
136 IF (cp_ddapc_ewald%do_qmmm_periodic_decpl .OR. &
137 cp_ddapc_ewald%do_decoupling)
THEN
141 param_section => cp_ddapc_ewald%ewald_section
150 IF (cp_ddapc_ewald%do_qmmm_periodic_decpl .OR. cp_ddapc_ewald%do_decoupling)
THEN
151 ALLOCATE (cp_ddapc_env%Mt(
SIZE(cp_ddapc_env%Md, 1),
SIZE(cp_ddapc_env%Md, 2)))
152 IF (cp_ddapc_ewald%do_decoupling)
THEN
154 cp_ddapc_env%Mt = cp_ddapc_env%Md
158 CALL ewald_ddapc_pot(cp_para_env, cp_ddapc_ewald%coeff_mm, -1.0_dp, super_cell, param_section, &
159 particle_set, cp_ddapc_env%Mr, radii)
160 cp_ddapc_env%Mt = cp_ddapc_env%Md + cp_ddapc_env%Mr
164 IF (cp_ddapc_ewald%do_solvation)
THEN
168 particle_set, cp_ddapc_env%Ms, radii)
170 CALL timestop(handle)
183 IF (
ASSOCIATED(cp_ddapc_env%AmI))
THEN
184 DEALLOCATE (cp_ddapc_env%AmI)
186 IF (
ASSOCIATED(cp_ddapc_env%Mt))
THEN
187 DEALLOCATE (cp_ddapc_env%Mt)
189 IF (
ASSOCIATED(cp_ddapc_env%Md))
THEN
190 DEALLOCATE (cp_ddapc_env%Md)
192 IF (
ASSOCIATED(cp_ddapc_env%Mr))
THEN
193 DEALLOCATE (cp_ddapc_env%Mr)
195 IF (
ASSOCIATED(cp_ddapc_env%Ms))
THEN
196 DEALLOCATE (cp_ddapc_env%Ms)
198 IF (
ASSOCIATED(cp_ddapc_env%gfunc))
THEN
199 DEALLOCATE (cp_ddapc_env%gfunc)
201 IF (
ASSOCIATED(cp_ddapc_env%w))
THEN
202 DEALLOCATE (cp_ddapc_env%w)
220 force_env_section, subsys_section, para_env)
222 LOGICAL,
INTENT(IN) :: qmmm_decoupl
227 INTEGER :: my_val, npts(3), output_unit
228 INTEGER,
DIMENSION(:),
POINTER :: ngrids
229 LOGICAL :: analyt, decoupling, &
230 do_qmmm_periodic_decpl, do_restraint, &
231 do_restraintb, do_solvation
232 REAL(kind=
dp) :: hmat(3, 3)
233 REAL(kind=
dp),
DIMENSION(:),
POINTER :: gx, gy, gz, lg
234 TYPE(
cell_type),
POINTER :: dummy_cell, mm_cell
236 TYPE(
section_vals_type),
POINTER :: cell_section, grid_print_section, multipole_section, &
237 poisson_section, printc_section, qmmm_per_section, restraint_section, restraint_sectionb, &
242 cpassert(.NOT.
ASSOCIATED(cp_ddapc_ewald))
243 ALLOCATE (cp_ddapc_ewald)
244 NULLIFY (cp_ddapc_ewald%pw_grid_mm, &
245 cp_ddapc_ewald%pw_grid_qm, &
246 cp_ddapc_ewald%ewald_section, &
247 cp_ddapc_ewald%pw_pool_mm, &
248 cp_ddapc_ewald%pw_pool_qm, &
249 cp_ddapc_ewald%coeff_mm, &
250 cp_ddapc_ewald%coeff_qm)
251 NULLIFY (multipole_section)
259 "PROPERTIES%ET_COUPLING%DDAPC_RESTRAINT_A")
264 do_qmmm_periodic_decpl = qmmm_decoupl
265 cp_ddapc_ewald%do_solvation = do_solvation
266 cp_ddapc_ewald%do_qmmm_periodic_decpl = do_qmmm_periodic_decpl
268 cp_ddapc_ewald%do_restraint = do_restraint .OR. do_restraintb
270 IF (do_qmmm_periodic_decpl .AND. decoupling)
THEN
272 IF (output_unit > 0) &
273 WRITE (output_unit,
'(T2,"WARNING",A)') &
274 "A calculation with the QMMM periodic model has been requested.", &
275 "The explicit POISSON section in DFT section will be IGNORED.", &
276 "QM Electrostatic controlled only by the PERIODIC section in QMMM section"
289 cp_ddapc_ewald%do_decoupling = decoupling
290 IF (cp_ddapc_ewald%do_qmmm_periodic_decpl)
THEN
294 cp_ddapc_ewald%ewald_section => multipole_section
295 IF (cp_ddapc_ewald%do_decoupling .OR. cp_ddapc_ewald%do_qmmm_periodic_decpl)
THEN
298 IF (.NOT. analyt)
THEN
302 NULLIFY (lg, gx, gy, gz)
304 CALL eval_lg(multipole_section, hmat, qm_cell%deth, lg, gx, gy, gz)
306 CALL setup_ewald_spline(pw_grid=cp_ddapc_ewald%pw_grid_qm, pw_pool=cp_ddapc_ewald%pw_pool_qm, &
307 coeff=cp_ddapc_ewald%coeff_qm, lg=lg, gx=gx, gy=gy, gz=gz, hmat=hmat, npts=npts, &
308 param_section=multipole_section, tag=
"ddapc", &
309 print_section=grid_print_section)
314 IF (cp_ddapc_ewald%do_qmmm_periodic_decpl)
THEN
315 NULLIFY (mm_cell, dummy_cell)
317 CALL read_cell(mm_cell, dummy_cell, cell_section=cell_section, para_env=para_env)
319 CALL eval_lg(multipole_section, hmat, mm_cell%deth, lg, gx, gy, gz)
321 CALL setup_ewald_spline(pw_grid=cp_ddapc_ewald%pw_grid_mm, pw_pool=cp_ddapc_ewald%pw_pool_mm, &
322 coeff=cp_ddapc_ewald%coeff_mm, lg=lg, gx=gx, gy=gy, gz=gz, hmat=hmat, npts=npts, &
323 param_section=multipole_section, tag=
"ddapc", &
324 print_section=grid_print_section)
349 SUBROUTINE eval_lg(multipole_section, hmat, deth, LG, gx, gy, gz)
351 REAL(kind=
dp),
INTENT(IN) :: hmat(3, 3), deth
352 REAL(kind=
dp),
DIMENSION(:),
POINTER :: lg, gx, gy, gz
354 INTEGER :: i, k1, k2, k3, n_rep, ndim, nmax1, &
356 REAL(kind=
dp) :: alpha, eps,
fac, fs, fvec(3), galpha, &
357 gsq, gsqi, rcut, tol, tol1
359 rcut = min(hmat(1, 1), hmat(2, 2), hmat(3, 3))/2.0_dp
363 eps = min(abs(eps), 0.5_dp)
364 tol = sqrt(abs(log(eps*rcut)))
365 alpha = sqrt(abs(log(eps*rcut*tol)))/rcut
366 galpha = 1.0_dp/(4.0_dp*alpha*alpha)
367 tol1 = sqrt(-log(eps*rcut*(2.0_dp*tol*alpha)**2))
368 nmax1 = nint(0.25_dp + hmat(1, 1)*alpha*tol1/
pi)
369 nmax2 = nint(0.25_dp + hmat(2, 2)*alpha*tol1/
pi)
370 nmax3 = nint(0.25_dp + hmat(3, 3)*alpha*tol1/
pi)
372 fvec = 2.0_dp*
pi/(/hmat(1, 1), hmat(2, 2), hmat(3, 3)/)
373 ndim = (nmax1 + 1)*(2*nmax2 + 1)*(2*nmax3 + 1) - 1
381 DO k2 = -nmax2, nmax2
382 DO k3 = -nmax3, nmax3
383 IF (k1 == 0 .AND. k2 == 0 .AND. k3 == 0) cycle
385 fs = 2.0_dp;
IF (k1 == 0) fs = 1.0_dp
386 gx(i) = fvec(1)*real(k1, kind=
dp)
387 gy(i) = fvec(2)*real(k2, kind=
dp)
388 gz(i) = fvec(3)*real(k3, kind=
dp)
389 gsq = gx(i)*gx(i) + gy(i)*gy(i) + gz(i)*gz(i)
391 lg(i) =
fac*gsqi*exp(-galpha*gsq)
396 END SUBROUTINE eval_lg
408 IF (
ASSOCIATED(cp_ddapc_ewald))
THEN
409 IF (
ASSOCIATED(cp_ddapc_ewald%coeff_qm))
THEN
410 CALL cp_ddapc_ewald%pw_pool_qm%give_back_pw(cp_ddapc_ewald%coeff_qm)
411 DEALLOCATE (cp_ddapc_ewald%coeff_qm)
413 IF (
ASSOCIATED(cp_ddapc_ewald%coeff_mm))
THEN
414 CALL cp_ddapc_ewald%pw_pool_mm%give_back_pw(cp_ddapc_ewald%coeff_mm)
415 DEALLOCATE (cp_ddapc_ewald%coeff_mm)
417 IF (
ASSOCIATED(cp_ddapc_ewald%pw_pool_qm))
THEN
419 cpassert(.NOT.
ASSOCIATED(cp_ddapc_ewald%pw_pool_qm))
421 IF (
ASSOCIATED(cp_ddapc_ewald%pw_pool_mm))
THEN
423 cpassert(.NOT.
ASSOCIATED(cp_ddapc_ewald%pw_pool_mm))
425 IF (
ASSOCIATED(cp_ddapc_ewald%pw_grid_qm))
THEN
427 cpassert(.NOT.
ASSOCIATED(cp_ddapc_ewald%pw_grid_qm))
429 IF (
ASSOCIATED(cp_ddapc_ewald%pw_grid_mm))
THEN
431 cpassert(.NOT.
ASSOCIATED(cp_ddapc_ewald%pw_grid_mm))
433 DEALLOCATE (cp_ddapc_ewald)
Handles all functions related to the CELL.
recursive subroutine, public read_cell(cell, cell_ref, use_ref_cell, cell_section, check_for_ref, para_env)
...
Handles all functions related to the CELL.
subroutine, public cell_release(cell)
releases the given cell (see doc/ReferenceCounting.html)
contains information regarding the decoupling/recoupling method of Bloechl
subroutine, public solvation_ddapc_pot(solvation_section, particle_set, m, radii)
Evaluates the electrostatic potential due to a simple solvation model Spherical cavity in a dieletric...
recursive subroutine, public ewald_ddapc_pot(cp_para_env, coeff, factor, cell, multipole_section, particle_set, m, radii)
Evaluates the Ewald term E2 and E3 energy term for the decoupling/coupling of periodic images.
subroutine, public ddapc_eval_ami(gami, c0, gfunc, w, particle_set, gcut, rho_tot_g, radii, iw, vol)
Computes the inverse AmI of the Am matrix.
subroutine, public ddapc_eval_gfunc(gfunc, w, gcut, rho_tot_g, radii)
...
contains information regarding the decoupling/recoupling method of Bloechl
subroutine, public cp_ddapc_release(cp_ddapc_env)
...
subroutine, public cp_ddapc_ewald_create(cp_ddapc_ewald, qmmm_decoupl, qm_cell, force_env_section, subsys_section, para_env)
...
subroutine, public cp_ddapc_ewald_release(cp_ddapc_ewald)
...
subroutine, public cp_ddapc_create(cp_para_env, cp_ddapc_env, cp_ddapc_ewald, particle_set, radii, cell, super_cell, rho_tot_g, gcut, iw2, vol, force_env_section)
...
various routines to log and control the output. The idea is that decisions about where to log should ...
integer function, public cp_logger_get_default_io_unit(logger)
returns the unit nr for the ionode (-1 on all other processors) skips as well checks if the procs cal...
type(cp_logger_type) function, pointer, public cp_get_default_logger()
returns the default logger
routines to handle the output, The idea is to remove the decision of wheter to output and what to out...
logical function, public cp_printkey_is_on(iteration_info, print_key)
returns true if the printlevel activates this printkey does not look if this iteration it should be p...
Setting up the Spline coefficients used to Interpolate the G-Term in Ewald sums.
subroutine, public setup_ewald_spline(pw_grid, pw_pool, coeff, lg, gx, gy, gz, hmat, npts, param_section, tag, print_section)
Setup of the G-space Ewald Term Spline Coefficients.
Defines the basic variable types.
integer, parameter, public dp
Definition of mathematical constants and functions.
real(kind=dp), parameter, public pi
real(kind=dp), dimension(0:maxfac), parameter, public fac
Interface to the message passing library MPI.
Define the data structure for the particle information.
This module defines the grid data type and some basic operations on it.
subroutine, public pw_grid_release(pw_grid)
releases the given pw grid
functions related to the poisson solver on regular grids
integer, parameter, public pw_poisson_multipole
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...
subroutine, public pw_pool_release(pool)
releases the given pool (see cp2k/doc/ReferenceCounting.html)
Type defining parameters related to the simulation cell.
type of a logger, at the moment it contains just a print level starting at which level it should be l...
stores all the informations relevant to an mpi environment
Manages a pool of grids (to be used for example as tmp objects), but can also be used to instantiate ...