54#include "./base/base_uses.f90"
62 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'qs_kinetic'
68 INTEGER,
DIMENSION(1:56),
SAVE :: ndod = [0, 1, 1, 1, 0, 0, 0, 0, 0, 0, 1, 1, 1, &
69 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, &
70 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, &
71 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1]
100 basis_type, sab_nl, calculate_forces, matrix_p, matrixkp_p, &
101 ext_kpoints, eps_filter, nderivative)
107 POINTER :: matrixkp_t
108 CHARACTER(LEN=*),
INTENT(IN),
OPTIONAL :: matrix_name
109 CHARACTER(LEN=*),
INTENT(IN) :: basis_type
112 LOGICAL,
INTENT(IN),
OPTIONAL :: calculate_forces
113 TYPE(
dbcsr_type),
OPTIONAL,
POINTER :: matrix_p
115 POINTER :: matrixkp_p
116 TYPE(
kpoint_type),
OPTIONAL,
POINTER :: ext_kpoints
117 REAL(kind=
dp),
INTENT(IN),
OPTIONAL :: eps_filter
118 INTEGER,
INTENT(IN),
OPTIONAL :: nderivative
124 CALL build_kinetic_matrix_low(ks_env, matrix_t, matrixkp_t, matrix_name, basis_type, &
125 sab_nl, calculate_forces, matrix_p, matrixkp_p, &
126 ext_kpoints, eps_filter, natom, nderivative)
147 SUBROUTINE build_kinetic_matrix_low(ks_env, matrix_t, matrixkp_t, matrix_name, basis_type, &
148 sab_nl, calculate_forces, matrix_p, matrixkp_p, &
149 ext_kpoints, eps_filter, natom, nderivative)
155 POINTER :: matrixkp_t
156 CHARACTER(LEN=*),
INTENT(IN),
OPTIONAL :: matrix_name
157 CHARACTER(LEN=*),
INTENT(IN) :: basis_type
160 LOGICAL,
INTENT(IN),
OPTIONAL :: calculate_forces
161 TYPE(
dbcsr_type),
OPTIONAL,
POINTER :: matrix_p
163 POINTER :: matrixkp_p
164 TYPE(
kpoint_type),
OPTIONAL,
POINTER :: ext_kpoints
165 REAL(kind=
dp),
INTENT(IN),
OPTIONAL :: eps_filter
166 INTEGER,
INTENT(IN) :: natom
167 INTEGER,
INTENT(IN),
OPTIONAL :: nderivative
169 CHARACTER(len=*),
PARAMETER :: routinen =
'build_kinetic_matrix_low'
171 INTEGER :: atom_a, handle, i, iatom, ic, icol, ikind, irow, iset, jatom, jkind, jset, ldsab, &
172 maxder, ncoa, ncob, nder, nimg, nkind, nseta, nsetb, sgfa, sgfb, slot
173 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: atom_of_kind, kind_of
174 INTEGER,
DIMENSION(3) :: cell
175 INTEGER,
DIMENSION(:),
POINTER :: la_max, la_min, lb_max, lb_min, npgfa, &
177 INTEGER,
DIMENSION(:, :),
POINTER :: first_sgfa, first_sgfb
178 INTEGER,
DIMENSION(:, :, :),
POINTER :: cell_to_index
179 LOGICAL :: do_forces, do_symmetric, dokp, found, &
180 trans, use_cell_mapping, use_virial
181 REAL(kind=
dp) :: f, f0, ff, rab2, tab
182 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :) :: pab, qab
183 REAL(kind=
dp),
ALLOCATABLE,
DIMENSION(:, :, :) :: kab
184 REAL(kind=
dp),
DIMENSION(3) :: force_a, rab
185 REAL(kind=
dp),
DIMENSION(3, 3) :: pv_thread
186 REAL(kind=
dp),
DIMENSION(3, natom) :: force_thread
187 REAL(kind=
dp),
DIMENSION(:),
POINTER :: set_radius_a, set_radius_b
188 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: p_block, rpgfa, rpgfb, scon_a, scon_b, &
190 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: dab
192 TYPE(
block_p_type),
ALLOCATABLE,
DIMENSION(:) :: k_block
199 TYPE(
qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
210 CALL timeset(routinen, handle)
212 NULLIFY (atomic_kind_set, qs_kind_set, p_block, dft_control)
214 dft_control=dft_control, &
215 atomic_kind_set=atomic_kind_set, &
216 qs_kind_set=qs_kind_set)
219 IF (
PRESENT(matrix_t))
THEN
221 use_cell_mapping = .false.
222 nimg = dft_control%nimages
223 ELSEIF (
PRESENT(matrixkp_t))
THEN
225 IF (
PRESENT(ext_kpoints))
THEN
226 IF (
ASSOCIATED(ext_kpoints))
THEN
228 use_cell_mapping = (
SIZE(cell_to_index) > 1)
229 nimg =
SIZE(ext_kpoints%index_to_cell, 2)
231 use_cell_mapping = .false.
235 CALL get_ks_env(ks_env=ks_env, kpoints=kpoints)
237 use_cell_mapping = (
SIZE(cell_to_index) > 1)
238 nimg = dft_control%nimages
244 nkind =
SIZE(atomic_kind_set)
247 IF (
PRESENT(calculate_forces)) do_forces = calculate_forces
250 IF (
PRESENT(nderivative)) nder = nderivative
254 cpassert(
SIZE(sab_nl) > 0)
258 ALLOCATE (basis_set_list(nkind))
263 CALL create_sab_matrix(ks_env, matrixkp_t, matrix_name, basis_set_list, basis_set_list, &
264 sab_nl, do_symmetric)
267 CALL create_sab_matrix(ks_env, matrix_t, matrix_name, basis_set_list, basis_set_list, &
268 sab_nl, do_symmetric)
274 CALL get_ks_env(ks_env=ks_env, force=force, virial=virial)
275 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
278 cpassert(
PRESENT(matrixkp_p))
280 IF (
PRESENT(matrix_p))
THEN
282 ELSEIF (
PRESENT(matrixkp_p))
THEN
283 matp => matrixkp_p(1, 1)%matrix
285 cpabort(
"Missing density matrix")
290 force_thread = 0.0_dp
317 ALLOCATE (kab(ldsab, ldsab, maxder), qab(ldsab, ldsab))
319 ALLOCATE (dab(ldsab, ldsab, 3), pab(ldsab, ldsab))
322 ALLOCATE (k_block(maxder))
324 NULLIFY (k_block(i)%block)
328 DO slot = 1, sab_nl(1)%nl_size
330 ikind = sab_nl(1)%nlist_task(slot)%ikind
331 jkind = sab_nl(1)%nlist_task(slot)%jkind
332 iatom = sab_nl(1)%nlist_task(slot)%iatom
333 jatom = sab_nl(1)%nlist_task(slot)%jatom
334 cell(:) = sab_nl(1)%nlist_task(slot)%cell(:)
335 rab(1:3) = sab_nl(1)%nlist_task(slot)%r(1:3)
337 basis_set_a => basis_set_list(ikind)%gto_basis_set
338 IF (.NOT.
ASSOCIATED(basis_set_a)) cycle
339 basis_set_b => basis_set_list(jkind)%gto_basis_set
340 IF (.NOT.
ASSOCIATED(basis_set_b)) cycle
346 first_sgfa => basis_set_a%first_sgf
347 la_max => basis_set_a%lmax
348 la_min => basis_set_a%lmin
349 npgfa => basis_set_a%npgf
350 nseta = basis_set_a%nset
351 nsgfa => basis_set_a%nsgf_set
352 rpgfa => basis_set_a%pgf_radius
353 set_radius_a => basis_set_a%set_radius
354 scon_a => basis_set_a%scon
355 zeta => basis_set_a%zet
357 first_sgfb => basis_set_b%first_sgf
358 lb_max => basis_set_b%lmax
359 lb_min => basis_set_b%lmin
360 npgfb => basis_set_b%npgf
361 nsetb = basis_set_b%nset
362 nsgfb => basis_set_b%nsgf_set
363 rpgfb => basis_set_b%pgf_radius
364 set_radius_b => basis_set_b%set_radius
365 scon_b => basis_set_b%scon
366 zetb => basis_set_b%zet
368 IF (use_cell_mapping)
THEN
369 ic = cell_to_index(cell(1), cell(2), cell(3))
375 IF (do_symmetric)
THEN
376 IF (iatom <= jatom)
THEN
384 IF (iatom == jatom) f0 = 1.0_dp
394 row=irow, col=icol, block=k_block(1)%block, found=found)
398 NULLIFY (k_block(i)%block)
400 row=irow, col=icol, block=k_block(i)%block, found=found)
409 row=irow, col=icol, block=p_block, found=found)
415 block=p_block, found=found)
420 rab2 = rab(1)*rab(1) + rab(2)*rab(2) + rab(3)*rab(3)
422 trans = do_symmetric .AND. (iatom > jatom)
426 ncoa = npgfa(iset)*(
ncoset(la_max(iset)) -
ncoset(la_min(iset) - 1))
427 sgfa = first_sgfa(1, iset)
431 IF (set_radius_a(iset) + set_radius_b(jset) < tab) cycle
436 ncob = npgfb(jset)*(
ncoset(lb_max(jset)) -
ncoset(lb_min(jset) - 1))
437 sgfb = first_sgfb(1, jset)
439 IF (do_forces .AND.
ASSOCIATED(p_block) .AND. ((iatom /= jatom) .OR. use_virial))
THEN
442 CALL block_add(
"OUT", kab(:, :, 1), nsgfa(iset), nsgfb(jset), p_block, sgfa, sgfb, trans=trans)
443 CALL decontraction(kab(:, :, 1), pab, scon_a(:, sgfa:), ncoa, nsgfa(iset), scon_b(:, sgfb:), ncob, nsgfb(jset), &
446 CALL kinetic(la_max(iset), la_min(iset), npgfa(iset), rpgfa(:, iset), zeta(:, iset), &
447 lb_max(jset), lb_min(jset), npgfb(jset), rpgfb(:, jset), zetb(:, jset), &
448 rab, kab(:, :, 1), dab)
450 force_thread(:, iatom) = force_thread(:, iatom) + ff*force_a(:)
451 force_thread(:, jatom) = force_thread(:, jatom) - ff*force_a(:)
458 CALL kinetic(la_max(iset), la_min(iset), npgfa(iset), rpgfa(:, iset), zeta(:, iset), &
459 lb_max(jset), lb_min(jset), npgfb(jset), rpgfb(:, jset), zetb(:, jset), &
460 rab, kab=kab(:, :, 1))
461 ELSE IF (nder == 1)
THEN
462 CALL kinetic(la_max(iset), la_min(iset), npgfa(iset), rpgfa(:, iset), zeta(:, iset), &
463 lb_max(jset), lb_min(jset), npgfb(jset), rpgfb(:, jset), zetb(:, jset), &
464 rab, kab=kab(:, :, 1), dab=kab(:, :, 2:4))
469 IF (ndod(i) == 1 .AND. trans) f = -1.0_dp
471 CALL contraction(kab(:, :, i), qab, ca=scon_a(:, sgfa:), na=ncoa, ma=nsgfa(iset), &
472 cb=scon_b(:, sgfb:), nb=ncob, mb=nsgfb(jset), fscale=f, &
475 CALL block_add(
"IN", qab, nsgfa(iset), nsgfb(jset), k_block(i)%block, sgfa, sgfb, trans=trans)
482 DEALLOCATE (kab, qab)
483 IF (do_forces)
DEALLOCATE (pab, dab)
502 atom_a = atom_of_kind(iatom)
503 ikind = kind_of(iatom)
504 force(ikind)%kinetic(:, atom_a) = force(ikind)%kinetic(:, atom_a) + force_thread(:, iatom)
508 IF (do_forces .AND. use_virial)
THEN
509 virial%pv_ekinetic = virial%pv_ekinetic + pv_thread
510 virial%pv_virial = virial%pv_virial + pv_thread
516 IF (
PRESENT(eps_filter))
THEN
522 IF (
PRESENT(eps_filter))
THEN
527 DEALLOCATE (basis_set_list)
529 CALL timestop(handle)
531 END SUBROUTINE build_kinetic_matrix_low
Set of routines to: Contract integrals over primitive Gaussians Decontract (density) matrices Trace m...
Calculation of the kinetic energy integrals over Cartesian Gaussian-type functions.
subroutine, public kinetic(la_max, la_min, npgfa, rpgfa, zeta, lb_max, lb_min, npgfb, rpgfb, zetb, rab, kab, dab)
Calculation of the two-center kinetic energy integrals [a|T|b] over Cartesian Gaussian-type functions...
Define the atomic kind types and their sub types.
subroutine, public get_atomic_kind_set(atomic_kind_set, atom_of_kind, kind_of, natom_of_kind, maxatom, natom, nshell, fist_potential_present, shell_present, shell_adiabatic, shell_check_distance, damping_present)
Get attributes of an atomic kind set.
collect pointers to a block of reals
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
subroutine, public dbcsr_get_block_p(matrix, row, col, block, found, row_size, col_size)
...
subroutine, public dbcsr_filter(matrix, eps)
...
subroutine, public dbcsr_finalize(matrix)
...
DBCSR operations in CP2K.
Defines the basic variable types.
integer, parameter, public int_8
integer, parameter, public dp
Types and basic routines needed for a kpoint calculation.
subroutine, public get_kpoint_info(kpoint, kp_scheme, nkp_grid, kp_shift, symmetry, verbose, full_grid, use_real_wfn, eps_geo, parallel_group_size, kp_range, nkp, xkp, wkp, para_env, blacs_env_all, para_env_kp, para_env_inter_kp, blacs_env, kp_env, kp_aux_env, mpools, iogrp, nkp_groups, kp_dist, cell_to_index, index_to_cell, sab_nl, sab_nl_nosym)
Retrieve information from a kpoint environment.
Provides Cartesian and spherical orbital pointers and indices.
integer, dimension(:), allocatable, public ncoset
Some utility functions for the calculation of integrals.
subroutine, public basis_set_list_setup(basis_set_list, basis_type, qs_kind_set)
Set up an easy accessible list of the basis sets for all kinds.
Define the quickstep kind type and their sub types.
Calculation of kinetic energy matrix and forces.
subroutine, public build_kinetic_matrix(ks_env, matrix_t, matrixkp_t, matrix_name, basis_type, sab_nl, calculate_forces, matrix_p, matrixkp_p, ext_kpoints, eps_filter, nderivative)
Calculation of the kinetic energy matrix over Cartesian Gaussian functions.
subroutine, public get_ks_env(ks_env, v_hartree_rspace, s_mstruct_changed, rho_changed, exc_accint, potential_changed, forces_up_to_date, complex_ks, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, kinetic, matrix_s, matrix_s_ri_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_ri_aux_kp, matrix_ks_im_kp, rho, rho_xc, vppl, xcint_weights, rho_core, rho_nlcc, rho_nlcc_g, vee, neighbor_list_id, sab_orb, sab_all, sac_ae, sac_ppl, sac_lri, sap_ppnl, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_pp, sab_xtb_nonbond, sab_vdw, sab_scp, sab_almo, sab_kp, sab_kp_nosym, sab_cneo, task_list, task_list_soft, kpoints, do_kpoints, atomic_kind_set, qs_kind_set, cell, cell_ref, use_ref_cell, particle_set, energy, force, local_particles, local_molecules, molecule_kind_set, molecule_set, subsys, cp_subsys, virial, results, atprop, nkind, natom, dft_control, dbcsr_dist, distribution_2d, pw_env, para_env, blacs_env, nelectron_total, nelectron_spin)
...
Define the neighbor list data types and the corresponding functionality.
subroutine, public get_neighbor_list_set_p(neighbor_list_sets, nlist, symmetric)
Return the components of the first neighbor list set.
Calculation of overlap matrix, its derivatives and forces.
pure subroutine, public virial_pair_force(pv_virial, f0, force, rab)
Computes the contribution to the stress tensor from two-body pair-wise forces.
Provides all information about an atomic kind.
Contains information about kpoints.
Provides all information about a quickstep kind.
calculation environment to calculate the ks matrix, holds all the needed vars. assumes that the core ...