17 USE dbcsr_api,
ONLY: dbcsr_add,&
19 dbcsr_iterator_blocks_left,&
20 dbcsr_iterator_next_block,&
21 dbcsr_iterator_start,&
39 neighbor_list_iterator_p_type,&
41 neighbor_list_set_p_type
47 #include "./base/base_uses.f90"
53 CHARACTER(len=*),
PARAMETER,
PRIVATE :: moduleN =
'qs_dftb3_methods'
73 sap_int, calculate_forces, just_energy)
75 TYPE(qs_environment_type),
POINTER :: qs_env
76 TYPE(dbcsr_p_type),
DIMENSION(:, :),
POINTER :: ks_matrix
77 TYPE(qs_rho_type),
POINTER :: rho
78 REAL(
dp),
DIMENSION(:),
INTENT(in) :: mcharge
79 TYPE(qs_energy_type),
POINTER :: energy
80 REAL(
dp),
DIMENSION(:),
INTENT(in) :: xgamma, zeff
81 TYPE(sap_int_type),
DIMENSION(:),
POINTER :: sap_int
82 LOGICAL,
INTENT(in) :: calculate_forces, just_energy
84 CHARACTER(len=*),
PARAMETER :: routinen =
'build_dftb3_diagonal'
86 INTEGER :: atom_i, atom_j, blk, handle, i, ia, iac, &
87 iatom, ic, icol, ikind, irow, is, &
88 jatom, jkind, natom, nimg, nkind
89 INTEGER,
ALLOCATABLE,
DIMENSION(:) :: atom_of_kind, kind_of
90 INTEGER,
DIMENSION(3) :: cellind
91 INTEGER,
DIMENSION(:, :, :),
POINTER :: cell_to_index
92 LOGICAL :: found, use_virial
93 REAL(kind=
dp) :: dr, eb3, eloc, fi, gmij, ua, ui, uj
94 REAL(kind=
dp),
DIMENSION(3) :: fij, rij
95 REAL(kind=
dp),
DIMENSION(:, :),
POINTER :: dsblock, ksblock, pblock, sblock
96 REAL(kind=
dp),
DIMENSION(:, :, :),
POINTER :: dsint
97 TYPE(atomic_kind_type),
DIMENSION(:),
POINTER :: atomic_kind_set
98 TYPE(atprop_type),
POINTER :: atprop
99 TYPE(cell_type),
POINTER :: cell
100 TYPE(dbcsr_iterator_type) :: iter
101 TYPE(dbcsr_p_type),
DIMENSION(:, :),
POINTER :: matrix_p, matrix_s
102 TYPE(distribution_1d_type),
POINTER :: local_particles
103 TYPE(kpoint_type),
POINTER :: kpoints
104 TYPE(mp_para_env_type),
POINTER :: para_env
105 TYPE(neighbor_list_iterator_p_type), &
106 DIMENSION(:),
POINTER :: nl_iterator
107 TYPE(neighbor_list_set_p_type),
DIMENSION(:), &
109 TYPE(particle_type),
DIMENSION(:),
POINTER :: particle_set
110 TYPE(qs_force_type),
DIMENSION(:),
POINTER :: force
111 TYPE(qs_kind_type),
DIMENSION(:),
POINTER :: qs_kind_set
112 TYPE(virial_type),
POINTER :: virial
114 CALL timeset(routinen, handle)
118 CALL get_qs_env(qs_env=qs_env, atomic_kind_set=atomic_kind_set, &
119 qs_kind_set=qs_kind_set, atprop=atprop)
121 kind_of=kind_of, atom_of_kind=atom_of_kind)
124 CALL get_qs_env(qs_env=qs_env, local_particles=local_particles)
125 DO ikind = 1,
SIZE(local_particles%n_el)
127 DO ia = 1, local_particles%n_el(ikind)
128 iatom = local_particles%list(ikind)%array(ia)
129 eloc = -1.0_dp/6.0_dp*ua*mcharge(iatom)**3
131 IF (atprop%energy)
THEN
133 eloc = -0.5_dp*eloc - 0.25_dp*ua*zeff(ikind)*mcharge(iatom)**2
134 atprop%atecoul(iatom) = atprop%atecoul(iatom) + eloc
138 CALL get_qs_env(qs_env=qs_env, para_env=para_env)
139 CALL para_env%sum(eb3)
143 IF (calculate_forces)
THEN
144 CALL get_qs_env(qs_env=qs_env, matrix_s_kp=matrix_s, natom=natom, force=force, &
145 cell=cell, virial=virial, particle_set=particle_set)
147 use_virial = virial%pv_availability .AND. (.NOT. virial%pv_numer)
149 IF (
SIZE(matrix_p, 1) == 2)
THEN
150 DO ic = 1,
SIZE(matrix_p, 2)
151 CALL dbcsr_add(matrix_p(1, ic)%matrix, matrix_p(2, ic)%matrix, &
152 alpha_scalar=1.0_dp, beta_scalar=1.0_dp)
156 nimg =
SIZE(matrix_p, 2)
157 NULLIFY (cell_to_index)
160 CALL get_qs_env(qs_env=qs_env, kpoints=kpoints)
165 CALL dbcsr_iterator_start(iter, matrix_s(1, 1)%matrix)
166 DO WHILE (dbcsr_iterator_blocks_left(iter))
167 CALL dbcsr_iterator_next_block(iter, irow, icol, sblock, blk)
168 ikind = kind_of(irow)
169 atom_i = atom_of_kind(irow)
171 jkind = kind_of(icol)
172 atom_j = atom_of_kind(icol)
175 gmij = -0.5_dp*(ui*mcharge(irow)**2 + uj*mcharge(icol)**2)
178 CALL dbcsr_get_block_p(matrix=matrix_p(1, 1)%matrix, &
179 row=irow, col=icol, block=pblock, found=found)
183 CALL dbcsr_get_block_p(matrix=matrix_s(1 + i, 1)%matrix, &
184 row=irow, col=icol, block=dsblock, found=found)
186 fi = -gmij*sum(pblock*dsblock)
187 force(ikind)%rho_elec(i, atom_i) = force(ikind)%rho_elec(i, atom_i) + fi
188 force(jkind)%rho_elec(i, atom_j) = force(jkind)%rho_elec(i, atom_j) - fi
192 CALL dbcsr_iterator_stop(iter)
195 cpassert(
ASSOCIATED(sap_int))
199 iac = ikind + nkind*(jkind - 1)
200 IF (.NOT.
ASSOCIATED(sap_int(iac)%alist)) cycle
203 DO ia = 1, sap_int(iac)%nalist
204 IF (.NOT.
ASSOCIATED(sap_int(iac)%alist(ia)%clist)) cycle
205 iatom = sap_int(iac)%alist(ia)%aatom
206 DO ic = 1, sap_int(iac)%alist(ia)%nclist
207 jatom = sap_int(iac)%alist(ia)%clist(ic)%catom
208 rij = sap_int(iac)%alist(ia)%clist(ic)%rac
209 dr = sqrt(sum(rij(:)**2))
210 IF (dr > 1.e-6_dp)
THEN
211 dsint => sap_int(iac)%alist(ia)%clist(ic)%acint
212 gmij = -0.5_dp*(ui*mcharge(iatom)**2 + uj*mcharge(jatom)**2)
213 icol = max(iatom, jatom)
214 irow = min(iatom, jatom)
216 CALL dbcsr_get_block_p(matrix=matrix_p(1, 1)%matrix, &
217 row=irow, col=icol, block=pblock, found=found)
220 IF (irow == iatom)
THEN
221 fij(i) = -gmij*sum(pblock*dsint(:, :, i))
223 fij(i) = -gmij*sum(transpose(pblock)*dsint(:, :, i))
227 IF (iatom == jatom) fi = 0.5_dp
229 IF (atprop%stress)
THEN
241 CALL get_qs_env(qs_env=qs_env, sab_orb=n_list)
245 iatom=iatom, jatom=jatom, r=rij, cell=cellind)
247 dr = sqrt(sum(rij**2))
248 IF (iatom == jatom .AND. dr < 1.0e-6_dp) cycle
250 icol = max(iatom, jatom)
251 irow = min(iatom, jatom)
253 ic = cell_to_index(cellind(1), cellind(2), cellind(3))
256 ikind = kind_of(iatom)
257 atom_i = atom_of_kind(iatom)
259 jkind = kind_of(jatom)
260 atom_j = atom_of_kind(jatom)
263 gmij = -0.5_dp*(ui*mcharge(iatom)**2 + uj*mcharge(jatom)**2)
266 CALL dbcsr_get_block_p(matrix=matrix_p(1, ic)%matrix, &
267 row=irow, col=icol, block=pblock, found=found)
271 CALL dbcsr_get_block_p(matrix=matrix_s(1 + i, ic)%matrix, &
272 row=irow, col=icol, block=dsblock, found=found)
274 IF (irow == iatom)
THEN
275 fi = -gmij*sum(pblock*dsblock)
277 fi = gmij*sum(pblock*dsblock)
279 force(ikind)%rho_elec(i, atom_i) = force(ikind)%rho_elec(i, atom_i) + fi
280 force(jkind)%rho_elec(i, atom_j) = force(jkind)%rho_elec(i, atom_j) - fi
285 IF (iatom == jatom) fi = 0.5_dp
287 IF (atprop%stress)
THEN
298 IF (
SIZE(matrix_p, 1) == 2)
THEN
299 DO ic = 1,
SIZE(matrix_p, 2)
300 CALL dbcsr_add(matrix_p(1, ic)%matrix, matrix_p(2, ic)%matrix, &
301 alpha_scalar=1.0_dp, beta_scalar=-1.0_dp)
307 IF (.NOT. just_energy)
THEN
308 CALL get_qs_env(qs_env=qs_env, matrix_s_kp=matrix_s, natom=natom)
310 nimg =
SIZE(ks_matrix, 2)
311 NULLIFY (cell_to_index)
314 CALL get_qs_env(qs_env=qs_env, kpoints=kpoints)
320 CALL dbcsr_iterator_start(iter, matrix_s(1, 1)%matrix)
321 DO WHILE (dbcsr_iterator_blocks_left(iter))
322 CALL dbcsr_iterator_next_block(iter, irow, icol, sblock, blk)
323 ikind = kind_of(irow)
325 jkind = kind_of(icol)
327 gmij = -0.5_dp*(ui*mcharge(irow)**2 + uj*mcharge(icol)**2)
328 DO is = 1,
SIZE(ks_matrix, 1)
330 CALL dbcsr_get_block_p(matrix=ks_matrix(is, 1)%matrix, &
331 row=irow, col=icol, block=ksblock, found=found)
333 ksblock = ksblock - 0.5_dp*gmij*sblock
336 CALL dbcsr_iterator_stop(iter)
339 CALL get_qs_env(qs_env=qs_env, sab_orb=n_list)
343 iatom=iatom, jatom=jatom, r=rij, cell=cellind)
345 icol = max(iatom, jatom)
346 irow = min(iatom, jatom)
348 ic = cell_to_index(cellind(1), cellind(2), cellind(3))
351 ikind = kind_of(iatom)
353 jkind = kind_of(jatom)
355 gmij = -0.5_dp*(ui*mcharge(iatom)**2 + uj*mcharge(jatom)**2)
358 CALL dbcsr_get_block_p(matrix=matrix_s(1, ic)%matrix, &
359 row=irow, col=icol, block=sblock, found=found)
361 DO is = 1,
SIZE(ks_matrix, 1)
363 CALL dbcsr_get_block_p(matrix=ks_matrix(is, ic)%matrix, &
364 row=irow, col=icol, block=ksblock, found=found)
366 ksblock = ksblock - 0.5_dp*gmij*sblock
375 CALL timestop(handle)
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.
Holds information on atomic properties.
Handles all functions related to the CELL.
stores a lists of integer that are local to a processor. The idea is that these integers represent ob...
Defines the basic variable types.
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.
Interface to the message passing library MPI.
Define the data structure for the particle information.
Calculation of DFTB3 Terms.
subroutine, public build_dftb3_diagonal(qs_env, ks_matrix, rho, mcharge, energy, xgamma, zeff, sap_int, calculate_forces, just_energy)
...
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.
Define the neighbor list data types and the corresponding functionality.
subroutine, public neighbor_list_iterator_create(iterator_set, nl, search, nthread)
Neighbor list iterator functions.
subroutine, public neighbor_list_iterator_release(iterator_set)
...
integer function, public neighbor_list_iterate(iterator_set, mepos)
...
subroutine, public get_iterator_info(iterator_set, mepos, ikind, jkind, nkind, ilist, nlist, inode, nnode, iatom, jatom, r, cell)
...
superstucture that hold various representations of the density and keeps track of which ones are vali...
subroutine, public qs_rho_get(rho_struct, rho_ao, rho_ao_im, rho_ao_kp, rho_ao_im_kp, rho_r, drho_r, rho_g, drho_g, tau_r, tau_g, rho_r_valid, drho_r_valid, rho_g_valid, drho_g_valid, tau_r_valid, tau_g_valid, tot_rho_r, tot_rho_g, rho_r_sccs, soft_valid, complex_rho_ao)
returns info about the density described by this object. If some representation is not available an e...
General overlap type integrals containers.
pure subroutine, public virial_pair_force(pv_virial, f0, force, rab)
Computes the contribution to the stress tensor from two-body pair-wise forces.