(git:e7e05ae)
qs_commutators.F
Go to the documentation of this file.
1 !--------------------------------------------------------------------------------------------------!
2 ! CP2K: A general program to perform molecular dynamics simulations !
3 ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 ! !
5 ! SPDX-License-Identifier: GPL-2.0-or-later !
6 !--------------------------------------------------------------------------------------------------!
7 
8 ! **************************************************************************************************
9 !> \brief Calculation of commutator [H,r] matrices
10 !> \par History
11 !> JGH: [7.2016]
12 !> \author Juerg Hutter
13 ! **************************************************************************************************
17  USE cp_control_types, ONLY: dft_control_type
20  USE dbcsr_api, ONLY: dbcsr_create,&
21  dbcsr_p_type,&
22  dbcsr_set
23  USE kinds, ONLY: dp
24  USE qs_environment_types, ONLY: get_qs_env,&
25  qs_environment_type
26  USE qs_kind_types, ONLY: qs_kind_type
27  USE qs_neighbor_list_types, ONLY: neighbor_list_set_p_type
28 
29 !$ USE OMP_LIB, ONLY: omp_get_max_threads, omp_get_thread_num
30 #include "./base/base_uses.f90"
31 
32  IMPLICIT NONE
33 
34  PRIVATE
35 
36 ! *** Global parameters ***
37 
38  CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'qs_commutators'
39 
40 ! *** Public subroutines ***
41 
42  PUBLIC :: build_com_hr_matrix
43 
44 CONTAINS
45 
46 ! **************************************************************************************************
47 !> \brief Calculation of the [H,r] commutators matrices over Cartesian Gaussian functions.
48 !> \param qs_env ...
49 !> \param matrix_hr ...
50 !> \date 26.07.2016
51 !> \par History
52 !> \author JGH
53 !> \version 1.0
54 ! **************************************************************************************************
55  SUBROUTINE build_com_hr_matrix(qs_env, matrix_hr)
56 
57  TYPE(qs_environment_type), POINTER :: qs_env
58  TYPE(dbcsr_p_type), DIMENSION(:), OPTIONAL, &
59  POINTER :: matrix_hr
60 
61  CHARACTER(len=*), PARAMETER :: routinen = 'build_com_hr_matrix'
62 
63  INTEGER :: handle, ir
64  REAL(kind=dp) :: eps_ppnl
65  TYPE(dbcsr_p_type), DIMENSION(:, :), POINTER :: matrix_s
66  TYPE(dft_control_type), POINTER :: dft_control
67  TYPE(neighbor_list_set_p_type), DIMENSION(:), &
68  POINTER :: sab_orb, sap_ppnl
69  TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
70 
71  CALL timeset(routinen, handle)
72 
73  NULLIFY (sab_orb, sap_ppnl)
74  CALL get_qs_env(qs_env=qs_env, sab_orb=sab_orb, sap_ppnl=sap_ppnl)
75  !
76  CALL get_qs_env(qs_env=qs_env, qs_kind_set=qs_kind_set, dft_control=dft_control)
77  eps_ppnl = dft_control%qs_control%eps_ppnl
78  !
79  CALL get_qs_env(qs_env=qs_env, matrix_s_kp=matrix_s)
80  cpassert(.NOT. ASSOCIATED(matrix_hr))
81  CALL dbcsr_allocate_matrix_set(matrix_hr, 3)
82  DO ir = 1, 3
83  ALLOCATE (matrix_hr(ir)%matrix)
84  CALL dbcsr_create(matrix_hr(ir)%matrix, template=matrix_s(1, 1)%matrix, &
85  name="COMMUTATOR")
86  CALL cp_dbcsr_alloc_block_from_nbl(matrix_hr(ir)%matrix, sab_orb)
87  CALL dbcsr_set(matrix_hr(ir)%matrix, 0.0_dp)
88  END DO
89 
90  CALL build_com_tr_matrix(matrix_hr, qs_kind_set, "ORB", sab_orb)
91  CALL build_com_rpnl(matrix_hr, qs_kind_set, sab_orb, sap_ppnl, eps_ppnl)
92 
93  CALL timestop(handle)
94 
95  END SUBROUTINE build_com_hr_matrix
96 
97 END MODULE qs_commutators
98 
Calculation of commutator of kinetic energy and position operator.
subroutine, public build_com_tr_matrix(matrix_tr, qs_kind_set, basis_type, sab_nl)
Calculation of commutator [T,r] over Cartesian Gaussian functions.
Calculation of the non-local pseudopotential contribution to the core Hamiltonian <a|V(non-local)|b> ...
subroutine, public build_com_rpnl(matrix_rv, qs_kind_set, sab_orb, sap_ppnl, eps_ppnl)
...
Defines control structures, which contain the parameters and the settings for the DFT-based calculati...
DBCSR operations in CP2K.
Defines the basic variable types.
Definition: kinds.F:23
integer, parameter, public dp
Definition: kinds.F:34
Calculation of commutator [H,r] matrices.
subroutine, public build_com_hr_matrix(qs_env, matrix_hr)
Calculation of the [H,r] commutators matrices over Cartesian Gaussian functions.
subroutine, public get_qs_env(qs_env, atomic_kind_set, qs_kind_set, cell, super_cell, cell_ref, use_ref_cell, kpoints, dft_control, mos, sab_orb, sab_all, qmmm, qmmm_periodic, sac_ae, sac_ppl, sac_lri, sap_ppnl, sab_vdw, sab_scp, sap_oce, sab_lrc, sab_se, sab_xtbe, sab_tbe, sab_core, sab_xb, sab_xtb_nonbond, sab_almo, sab_kp, sab_kp_nosym, particle_set, energy, force, matrix_h, matrix_h_im, matrix_ks, matrix_ks_im, matrix_vxc, run_rtp, rtp, matrix_h_kp, matrix_h_im_kp, matrix_ks_kp, matrix_ks_im_kp, matrix_vxc_kp, kinetic_kp, matrix_s_kp, matrix_w_kp, matrix_s_RI_aux_kp, matrix_s, matrix_s_RI_aux, matrix_w, matrix_p_mp2, matrix_p_mp2_admm, rho, rho_xc, pw_env, ewald_env, ewald_pw, active_space, mpools, input, para_env, blacs_env, scf_control, rel_control, kinetic, qs_charges, vppl, rho_core, rho_nlcc, rho_nlcc_g, ks_env, ks_qmmm_env, wf_history, scf_env, local_particles, local_molecules, distribution_2d, dbcsr_dist, molecule_kind_set, molecule_set, subsys, cp_subsys, oce, local_rho_set, rho_atom_set, task_list, task_list_soft, rho0_atom_set, rho0_mpole, rhoz_set, ecoul_1c, rho0_s_rs, rho0_s_gs, do_kpoints, has_unit_metric, requires_mo_derivs, mo_derivs, mo_loc_history, nkind, natom, nelectron_total, nelectron_spin, efield, neighbor_list_id, linres_control, xas_env, virial, cp_ddapc_env, cp_ddapc_ewald, outer_scf_history, outer_scf_ihistory, x_data, et_coupling, dftb_potential, results, se_taper, se_store_int_env, se_nddo_mpole, se_nonbond_env, admm_env, lri_env, lri_density, exstate_env, ec_env, dispersion_env, gcp_env, vee, rho_external, external_vxc, mask, mp2_env, bs_env, kg_env, WannierCentres, atprop, ls_scf_env, do_transport, transport_env, v_hartree_rspace, s_mstruct_changed, rho_changed, potential_changed, forces_up_to_date, mscfg_env, almo_scf_env, gradient_history, variable_history, embed_pot, spin_embed_pot, polar_env, mos_last_converged, rhs)
Get the QUICKSTEP environment.
Define the quickstep kind type and their sub types.
Definition: qs_kind_types.F:23
Define the neighbor list data types and the corresponding functionality.