(git:ed6f26b)
Loading...
Searching...
No Matches
negf_atom_map.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2025 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Map atoms between various force environments.
10!> \author Sergey Chulkov
11! **************************************************************************************************
12
13MODULE negf_atom_map
14 USE atomic_kind_types, ONLY: get_atomic_kind
15 USE cell_types, ONLY: cell_type,&
16 real_to_scaled
17 USE kinds, ONLY: default_string_length,&
18 dp
19 USE particle_types, ONLY: particle_type
20 USE qs_kind_types, ONLY: get_qs_kind,&
21 qs_kind_type
22 USE qs_subsys_types, ONLY: qs_subsys_get,&
23 qs_subsys_type
24#include "./base/base_uses.f90"
25
26 IMPLICIT NONE
27 PRIVATE
28
29 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'negf_atom_map'
30 LOGICAL, PARAMETER, PRIVATE :: debug_this_module = .true.
31
32 PUBLIC :: negf_atom_map_type, negf_map_atomic_indices
33
34! **************************************************************************************************
35!> \brief Structure that maps the given atom in the sourse FORCE_EVAL section with another atom
36!> from the target FORCE_EVAL section.
37! **************************************************************************************************
38 TYPE negf_atom_map_type
39 !> atomic index within the target FORCE_EVAL
40 INTEGER :: iatom = -1
41 !> cell replica
42 INTEGER, DIMENSION(3) :: cell = -1
43 END TYPE negf_atom_map_type
44
45 PRIVATE :: qs_kind_group_type, qs_kind_groups_create, qs_kind_groups_release
46
47! **************************************************************************************************
48!> \brief List of equivalent atoms.
49! **************************************************************************************************
50 TYPE qs_kind_group_type
51 !> atomic symbol
52 CHARACTER(len=2) :: element_symbol = ""
53 !> number of atoms of this kind in 'atom_list'
54 INTEGER :: natoms = -1
55 !> number of spherical Gaussian functions per atom
56 INTEGER :: nsgf = -1
57 !> list of atomic indices
58 INTEGER, ALLOCATABLE, DIMENSION(:) :: atom_list
59 !> atomic coordinates [3 x natoms]
60 REAL(kind=dp), ALLOCATABLE, DIMENSION(:, :) :: r
61 END TYPE qs_kind_group_type
62
63CONTAINS
64! **************************************************************************************************
65!> \brief Map atoms in the cell 'subsys_device' listed in 'atom_list' with the corresponding
66!> atoms in the cell 'subsys_contact'.
67!> \param atom_map list of atoms in the cell 'subsys_contact' (initialised on exit)
68!> \param atom_list atomic indices of selected atoms in the cell 'subsys_device'
69!> \param subsys_device QuickStep subsystem of the device force environment
70!> \param subsys_contact QuickStep subsystem of the contact force environment
71!> \param eps_geometry accuracy in mapping atoms based on their Cartesian coordinates
72!> \par History
73!> * 08.2017 created [Sergey Chulkov]
74! **************************************************************************************************
75 SUBROUTINE negf_map_atomic_indices(atom_map, atom_list, subsys_device, subsys_contact, eps_geometry)
76 TYPE(negf_atom_map_type), DIMENSION(:), &
77 INTENT(out) :: atom_map
78 INTEGER, DIMENSION(:), INTENT(in) :: atom_list
79 TYPE(qs_subsys_type), POINTER :: subsys_device, subsys_contact
80 REAL(kind=dp), INTENT(in) :: eps_geometry
81
82 CHARACTER(len=*), PARAMETER :: routinen = 'negf_map_atomic_indices'
83
84 CHARACTER(len=2) :: element_device
85 CHARACTER(len=default_string_length) :: atom_str
86 INTEGER :: atom_index_device, handle, iatom, &
87 iatom_kind, ikind, ikind_contact, &
88 natoms, nkinds_contact, nsgf_device
89 REAL(kind=dp), DIMENSION(3) :: coords, coords_error, coords_scaled
90 TYPE(cell_type), POINTER :: cell_contact
91 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set_contact, particle_set_device
92 TYPE(qs_kind_group_type), ALLOCATABLE, &
93 DIMENSION(:) :: kind_groups_contact
94 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set_contact, qs_kind_set_device
95
96 CALL timeset(routinen, handle)
97
98 natoms = SIZE(atom_map, 1)
99 cpassert(SIZE(atom_list) == natoms)
100
101 CALL qs_subsys_get(subsys_device, particle_set=particle_set_device, qs_kind_set=qs_kind_set_device)
102 CALL qs_subsys_get(subsys_contact, cell=cell_contact, particle_set=particle_set_contact, qs_kind_set=qs_kind_set_contact)
103
104 CALL qs_kind_groups_create(kind_groups_contact, particle_set_contact, qs_kind_set_contact)
105 nkinds_contact = SIZE(kind_groups_contact)
106
107 DO iatom = 1, natoms
108 atom_index_device = atom_list(iatom)
109 CALL get_atomic_kind(particle_set_device(atom_index_device)%atomic_kind, kind_number=ikind)
110 CALL get_qs_kind(qs_kind_set_device(ikind), element_symbol=element_device, nsgf=nsgf_device)
111
112 atom_map(iatom)%iatom = 0
113
114 iterate_kind: DO ikind_contact = 1, nkinds_contact
115 ! looking for an equivalent atomic kind (based on the element symbol and the number of atomic basis functions)
116 IF (kind_groups_contact(ikind_contact)%element_symbol == element_device .AND. &
117 kind_groups_contact(ikind_contact)%nsgf == nsgf_device) THEN
118
119 ! loop over matching atoms
120 DO iatom_kind = 1, kind_groups_contact(ikind_contact)%natoms
121 coords(1:3) = particle_set_device(atom_index_device)%r(1:3) - &
122 kind_groups_contact(ikind_contact)%r(1:3, iatom_kind)
123
124 CALL real_to_scaled(coords_scaled, coords, cell_contact)
125 coords_error = coords_scaled - real(nint(coords_scaled), kind=dp)
126
127 IF (sqrt(dot_product(coords_error, coords_error)) < eps_geometry) THEN
128 atom_map(iatom)%iatom = kind_groups_contact(ikind_contact)%atom_list(iatom_kind)
129 atom_map(iatom)%cell = nint(coords_scaled)
130 EXIT iterate_kind
131 END IF
132 END DO
133 END IF
134 END DO iterate_kind
135
136 IF (atom_map(iatom)%iatom == 0) THEN
137 ! atom has not been found in the corresponding force_env
138 WRITE (atom_str, '(A2,3(1X,F11.6))') element_device, particle_set_device(atom_index_device)%r
139
140 CALL cp_abort(__location__, &
141 "Unable to map the atom ("//trim(atom_str)//") onto the atom from the corresponding FORCE_EVAL section")
142 END IF
143 END DO
144
145 CALL qs_kind_groups_release(kind_groups_contact)
146
147 CALL timestop(handle)
148 END SUBROUTINE negf_map_atomic_indices
149
150! **************************************************************************************************
151!> \brief Group particles from 'particle_set' according to their atomic (QS) kind.
152!> \param kind_groups kind groups that will be created
153!> \param particle_set list of particles
154!> \param qs_kind_set list of QS kinds
155!> \par History
156!> * 08.2017 created [Sergey Chulkov]
157!> \note used within the subroutine negf_map_atomic_indices() in order to map atoms in
158!> a linear scalling fashion
159! **************************************************************************************************
160 SUBROUTINE qs_kind_groups_create(kind_groups, particle_set, qs_kind_set)
161 TYPE(qs_kind_group_type), ALLOCATABLE, &
162 DIMENSION(:), INTENT(inout) :: kind_groups
163 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
164 TYPE(qs_kind_type), DIMENSION(:), POINTER :: qs_kind_set
165
166 CHARACTER(len=*), PARAMETER :: routinen = 'qs_kind_groups_create'
167
168 INTEGER :: handle, iatom, ikind, natoms, nkinds
169
170 CALL timeset(routinen, handle)
171
172 natoms = SIZE(particle_set)
173 nkinds = 0
174
175 DO iatom = 1, natoms
176 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
177 IF (nkinds < ikind) nkinds = ikind
178 END DO
179
180 ALLOCATE (kind_groups(nkinds))
181
182 DO ikind = 1, nkinds
183 kind_groups(ikind)%natoms = 0
184 CALL get_qs_kind(qs_kind_set(ikind), element_symbol=kind_groups(ikind)%element_symbol, nsgf=kind_groups(ikind)%nsgf)
185 END DO
186
187 DO iatom = 1, natoms
188 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
189 kind_groups(ikind)%natoms = kind_groups(ikind)%natoms + 1
190 END DO
191
192 DO ikind = 1, nkinds
193 ALLOCATE (kind_groups(ikind)%atom_list(kind_groups(ikind)%natoms))
194 ALLOCATE (kind_groups(ikind)%r(3, kind_groups(ikind)%natoms))
195
196 kind_groups(ikind)%natoms = 0
197 END DO
198
199 DO iatom = 1, natoms
200 CALL get_atomic_kind(particle_set(iatom)%atomic_kind, kind_number=ikind)
201 kind_groups(ikind)%natoms = kind_groups(ikind)%natoms + 1
202
203 kind_groups(ikind)%atom_list(kind_groups(ikind)%natoms) = iatom
204 kind_groups(ikind)%r(1:3, kind_groups(ikind)%natoms) = particle_set(iatom)%r(1:3)
205 END DO
206
207 CALL timestop(handle)
208 END SUBROUTINE qs_kind_groups_create
209
210! **************************************************************************************************
211!> \brief Release groups of particles.
212!> \param kind_groups kind groups to release
213!> \par History
214!> * 08.2017 created [Sergey Chulkov]
215! **************************************************************************************************
216 SUBROUTINE qs_kind_groups_release(kind_groups)
217 TYPE(qs_kind_group_type), ALLOCATABLE, &
218 DIMENSION(:), INTENT(inout) :: kind_groups
219
220 CHARACTER(len=*), PARAMETER :: routinen = 'qs_kind_groups_release'
221
222 INTEGER :: handle, ikind
223
224 CALL timeset(routinen, handle)
225
226 IF (ALLOCATED(kind_groups)) THEN
227 DO ikind = SIZE(kind_groups), 1, -1
228 IF (ALLOCATED(kind_groups(ikind)%atom_list)) DEALLOCATE (kind_groups(ikind)%atom_list)
229 IF (ALLOCATED(kind_groups(ikind)%r)) DEALLOCATE (kind_groups(ikind)%r)
230 END DO
231
232 DEALLOCATE (kind_groups)
233 END IF
234
235 CALL timestop(handle)
236 END SUBROUTINE qs_kind_groups_release
237END MODULE negf_atom_map
atomic_kind_types
Define the atomic kind types and their sub types.
Definition atomic_kind_types.F:23
atomic_kind_types::get_atomic_kind
subroutine, public get_atomic_kind(atomic_kind, fist_potential, element_symbol, name, mass, kind_number, natom, atom_list, rcov, rvdw, z, qeff, apol, cpol, mm_radius, shell, shell_active, damping)
Get attributes of an atomic kind.
Definition atomic_kind_types.F:138
cell_types
Handles all functions related to the CELL.
Definition cell_types.F:15
cell_types::real_to_scaled
subroutine, public real_to_scaled(s, r, cell)
Transform real to scaled cell coordinates. s=h_inv*r.
Definition cell_types.F:486
kinds
Defines the basic variable types.
Definition kinds.F:23
kinds::dp
integer, parameter, public dp
Definition kinds.F:34
kinds::default_string_length
integer, parameter, public default_string_length
Definition kinds.F:57
negf_atom_map
Map atoms between various force environments.
Definition negf_atom_map.F:13
negf_atom_map::negf_map_atomic_indices
subroutine, public negf_map_atomic_indices(atom_map, atom_list, subsys_device, subsys_contact, eps_geometry)
Map atoms in the cell 'subsys_device' listed in 'atom_list' with the corresponding atoms in the cell ...
Definition negf_atom_map.F:76
particle_types
Define the data structure for the particle information.
Definition particle_types.F:19
qs_kind_types
Define the quickstep kind type and their sub types.
Definition qs_kind_types.F:23
qs_kind_types::get_qs_kind
subroutine, public get_qs_kind(qs_kind, basis_set, basis_type, ncgf, nsgf, all_potential, tnadd_potential, gth_potential, sgp_potential, upf_potential, se_parameter, dftb_parameter, xtb_parameter, dftb3_param, zatom, zeff, elec_conf, mao, lmax_dftb, alpha_core_charge, ccore_charge, core_charge, core_charge_radius, paw_proj_set, paw_atom, hard_radius, hard0_radius, max_rad_local, covalent_radius, vdw_radius, gpw_type_forced, harmonics, max_iso_not0, max_s_harm, grid_atom, ngrid_ang, ngrid_rad, lmax_rho0, dft_plus_u_atom, l_of_dft_plus_u, n_of_dft_plus_u, u_minus_j, u_of_dft_plus_u, j_of_dft_plus_u, alpha_of_dft_plus_u, beta_of_dft_plus_u, j0_of_dft_plus_u, occupation_of_dft_plus_u, dispersion, bs_occupation, magnetization, no_optimize, addel, laddel, naddel, orbitals, max_scf, eps_scf, smear, u_ramping, u_minus_j_target, eps_u_ramping, init_u_ramping_each_scf, reltmat, ghost, floating, name, element_symbol, pao_basis_size, pao_model_file, pao_potentials, pao_descriptors, nelec)
Get attributes of an atomic kind.
Definition qs_kind_types.F:454
qs_subsys_types
types that represent a quickstep subsys
Definition qs_subsys_types.F:12
qs_subsys_types::qs_subsys_get
subroutine, public qs_subsys_get(subsys, atomic_kinds, atomic_kind_set, particles, particle_set, local_particles, molecules, molecule_set, molecule_kinds, molecule_kind_set, local_molecules, para_env, colvar_p, shell_particles, core_particles, gci, multipoles, natom, nparticle, ncore, nshell, nkind, atprop, virial, results, cell, cell_ref, use_ref_cell, energy, force, qs_kind_set, cp_subsys, nelectron_total, nelectron_spin)
...
Definition qs_subsys_types.F:134
cell_types::cell_type
Type defining parameters related to the simulation cell.
Definition cell_types.F:55
negf_atom_map::negf_atom_map_type
Structure that maps the given atom in the sourse FORCE_EVAL section with another atom from the target...
Definition negf_atom_map.F:38
particle_types::particle_type
Definition particle_types.F:35
qs_kind_types::qs_kind_type
Provides all information about a quickstep kind.
Definition qs_kind_types.F:171
qs_subsys_types::qs_subsys_type
Definition qs_subsys_types.F:56