(git:7641cd9)
Loading...
Searching...
No Matches
negf_control_types.F
Go to the documentation of this file.
1!--------------------------------------------------------------------------------------------------!
2! CP2K: A general program to perform molecular dynamics simulations !
3! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4! !
5! SPDX-License-Identifier: GPL-2.0-or-later !
6!--------------------------------------------------------------------------------------------------!
7
8! **************************************************************************************************
9!> \brief Input control types for NEGF based quantum transport calculations
10! **************************************************************************************************
11
12MODULE negf_control_types
13 USE cp_subsys_types, ONLY: cp_subsys_get,&
14 cp_subsys_type
15 USE input_constants, ONLY: negf_run
16 USE input_section_types, ONLY: section_vals_get,&
17 section_vals_get_subs_vals,&
18 section_vals_type,&
19 section_vals_val_get
20 USE kinds, ONLY: default_string_length,&
21 dp
22 USE mathconstants, ONLY: pi
23 USE molecule_kind_types, ONLY: get_molecule_kind,&
24 molecule_kind_type
25 USE molecule_types, ONLY: get_molecule,&
26 molecule_type
27 USE negf_alloc_types, ONLY: negf_allocatable_ivector
28 USE particle_types, ONLY: particle_type
29 USE physcon, ONLY: kelvin
30 USE string_utilities, ONLY: integer_to_string
31 USE util, ONLY: sort
32#include "./base/base_uses.f90"
33
34 IMPLICIT NONE
35 PRIVATE
36
37 CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'negf_control_types'
38 LOGICAL, PARAMETER, PRIVATE :: debug_this_module = .true.
39
40 PUBLIC :: negf_control_type, negf_control_contact_type
41 PUBLIC :: negf_control_create, negf_control_release, read_negf_control
42
43! **************************************************************************************************
44!> \brief Input parameters related to a single contact.
45!> \author Sergey Chulkov
46! **************************************************************************************************
47 TYPE negf_control_contact_type
48 !> atoms belonging to bulk and screening regions
49 INTEGER, ALLOCATABLE, DIMENSION(:) :: atomlist_bulk, atomlist_screening
50 !> atoms belonging to the primary and secondary bulk unit cells
51 TYPE(negf_allocatable_ivector), ALLOCATABLE, &
52 DIMENSION(:) :: atomlist_cell
53 !> index of the sub_force_env which should be used for bulk calculation
54 INTEGER :: force_env_index = -1
55 !> contact Fermi level needs to be computed
56 LOGICAL :: compute_fermi_level = .false.
57 !> to refine contact Fermi level using NEGF
58 LOGICAL :: refine_fermi_level = .false.
59 !> to shift energies to common zero level
60 LOGICAL :: shift_fermi_level = .false.
61 !> to read/write H and S from/to file
62 LOGICAL :: read_write_hs = .false.
63 !> if restart from files is really done
64 LOGICAL :: is_restart = .false.
65 !> Fermi level or starting Fermi level
66 REAL(kind=dp) :: fermi_level = -1.0_dp
67 !> Fermi level shifted to the common zero-energy level
68 REAL(kind=dp) :: fermi_level_shifted = -1.0_dp
69 !> temperature [in a.u.]
70 REAL(kind=dp) :: temperature = -1.0_dp
71 !> applied electric potential
72 REAL(kind=dp) :: v_external = 0.0_dp
73 END TYPE negf_control_contact_type
74
75! **************************************************************************************************
76!> \brief Input parameters related to the NEGF run.
77!> \author Sergey Chulkov
78! **************************************************************************************************
79 TYPE negf_control_type
80 !> input options for every contact
81 TYPE(negf_control_contact_type), ALLOCATABLE, &
82 DIMENSION(:) :: contacts
83 !> atoms belonging to the scattering region
84 INTEGER, ALLOCATABLE, DIMENSION(:) :: atomlist_s
85 !> atoms belonging to the scattering region as well as atoms belonging to
86 !> screening regions of all the contacts
87 INTEGER, ALLOCATABLE, DIMENSION(:) :: atomlist_s_screening
88 !> to read/write H and S from/to file
89 LOGICAL :: read_write_hs = .false.
90 !> to update the atomic Hamiltonian during NEGF self-consistent cycle
91 LOGICAL :: update_hs = .true.
92 !> if NEGF SCF id restart from saved files
93 LOGICAL :: restart_scf = .true.
94 !> if dft of entire system is done
95 LOGICAL :: is_dft_entire = .false.
96 !> if restart from files is really done
97 LOGICAL :: is_restart = .false.
98 !> the common restart file projectname-negf.restart is written if any of is_restart is .TRUE.
99 LOGICAL :: write_common_restart_file = .false.
100 !> do not keep contact self-energy matrices
101 LOGICAL :: disable_cache = .false.
102 !> convergence criteria for adaptive integration methods
103 REAL(kind=dp) :: conv_density = -1.0_dp
104 !> convergence criteria for iterative Lopez-Sancho algorithm
105 REAL(kind=dp) :: conv_green = -1.0_dp
106 !> convergence criteria for self-consistent iterations
107 REAL(kind=dp) :: conv_scf = -1.0_dp
108 !> accuracy in mapping atoms between different force environments
109 REAL(kind=dp) :: eps_geometry = -1.0_dp
110 !> applied bias [in a.u.]
111 REAL(kind=dp) :: v_bias = -1.0_dp
112 !> integration lower bound [in a.u.]
113 REAL(kind=dp) :: energy_lbound = -1.0_dp
114 !> infinitesimal offset along the imaginary axis [in a.u.]
115 REAL(kind=dp) :: eta = -1.0_dp
116 !> initial guess to determine the actual Fermi level of bulk contacts [in a.u.]
117 REAL(kind=dp) :: homo_lumo_gap = -1.0_dp
118 !> number of residuals (poles of the Fermi function)
119 INTEGER :: delta_npoles = -1
120 !> offset along the x-axis away from the poles of the Fermi function [in units of kT]
121 INTEGER :: gamma_kt = -1
122 !> integration method
123 INTEGER :: integr_method = -1
124 !> minimal number of grid points along the closed contour
125 INTEGER :: integr_min_points = -1
126 !> maximal number of grid points along the closed contour
127 INTEGER :: integr_max_points = -1
128 !> maximal number of SCF iterations
129 INTEGER :: max_scf = -1
130 !> minimal number of MPI processes to be used to compute Green's function per energy point
131 INTEGER :: nprocs = -1
132 !> shift in Hartree potential [in a.u.]
133 REAL(kind=dp) :: v_shift = -1.0_dp
134 !> initial offset to determine the correct shift in Hartree potential [in a.u.]
135 REAL(kind=dp) :: v_shift_offset = -1.0_dp
136 !> maximal number of iteration to determine the shift in Hartree potential
137 INTEGER :: v_shift_maxiters = -1
138 END TYPE negf_control_type
139
140 PRIVATE :: read_negf_atomlist
141
142CONTAINS
143
144! **************************************************************************************************
145!> \brief allocate control options for Non-equilibrium Green's Function calculation
146!> \param negf_control an object to create
147!> \par History
148!> * 02.2017 created [Sergey Chulkov]
149! **************************************************************************************************
150 SUBROUTINE negf_control_create(negf_control)
151 TYPE(negf_control_type), POINTER :: negf_control
152
153 CHARACTER(len=*), PARAMETER :: routinen = 'negf_control_create'
154
155 INTEGER :: handle
156
157 cpassert(.NOT. ASSOCIATED(negf_control))
158 CALL timeset(routinen, handle)
159
160 ALLOCATE (negf_control)
161
162 CALL timestop(handle)
163 END SUBROUTINE negf_control_create
164
165! **************************************************************************************************
166!> \brief release memory allocated for NEGF control options
167!> \param negf_control an object to release
168!> \par History
169!> * 02.2017 created [Sergey Chulkov]
170! **************************************************************************************************
171 SUBROUTINE negf_control_release(negf_control)
172 TYPE(negf_control_type), POINTER :: negf_control
173
174 CHARACTER(len=*), PARAMETER :: routinen = 'negf_control_release'
175
176 INTEGER :: handle, i, j
177
178 CALL timeset(routinen, handle)
179
180 IF (ASSOCIATED(negf_control)) THEN
181 IF (ALLOCATED(negf_control%atomlist_S)) DEALLOCATE (negf_control%atomlist_S)
182 IF (ALLOCATED(negf_control%atomlist_S_screening)) DEALLOCATE (negf_control%atomlist_S_screening)
183
184 IF (ALLOCATED(negf_control%contacts)) THEN
185 DO i = SIZE(negf_control%contacts), 1, -1
186 IF (ALLOCATED(negf_control%contacts(i)%atomlist_bulk)) &
187 DEALLOCATE (negf_control%contacts(i)%atomlist_bulk)
188
189 IF (ALLOCATED(negf_control%contacts(i)%atomlist_screening)) &
190 DEALLOCATE (negf_control%contacts(i)%atomlist_screening)
191
192 IF (ALLOCATED(negf_control%contacts(i)%atomlist_cell)) THEN
193 DO j = SIZE(negf_control%contacts(i)%atomlist_cell), 1, -1
194 IF (ALLOCATED(negf_control%contacts(i)%atomlist_cell(j)%vector)) &
195 DEALLOCATE (negf_control%contacts(i)%atomlist_cell(j)%vector)
196 END DO
197 DEALLOCATE (negf_control%contacts(i)%atomlist_cell)
198 END IF
199 END DO
200
201 DEALLOCATE (negf_control%contacts)
202 END IF
203
204 DEALLOCATE (negf_control)
205 END IF
206
207 CALL timestop(handle)
208 END SUBROUTINE negf_control_release
209
210! **************************************************************************************************
211!> \brief Read NEGF input parameters.
212!> \param negf_control NEGF control parameters
213!> \param input root input section
214!> \param subsys subsystem environment
215! **************************************************************************************************
216 SUBROUTINE read_negf_control(negf_control, input, subsys)
217 TYPE(negf_control_type), POINTER :: negf_control
218 TYPE(section_vals_type), POINTER :: input
219 TYPE(cp_subsys_type), POINTER :: subsys
220
221 CHARACTER(len=*), PARAMETER :: routinen = 'read_negf_control'
222
223 CHARACTER(len=default_string_length) :: contact_id_str, eta_current_str, eta_max_str, &
224 npoles_current_str, npoles_min_str, temp_current_str, temp_min_str
225 INTEGER :: delta_npoles_min, handle, i2_rep, i_rep, &
226 n2_rep, n_rep, natoms_current, &
227 natoms_total, run_type
228 INTEGER, ALLOCATABLE, DIMENSION(:) :: inds
229 LOGICAL :: do_negf, is_explicit
230 REAL(kind=dp) :: eta_max, temp_current, temp_min
231 TYPE(section_vals_type), POINTER :: cell_section, contact_section, &
232 negf_section, region_section, &
233 subsection
234
235 CALL timeset(routinen, handle)
236
237 CALL section_vals_val_get(input, "GLOBAL%RUN_TYPE", i_val=run_type)
238 do_negf = run_type == negf_run
239
240 negf_section => section_vals_get_subs_vals(input, "NEGF")
241
242 contact_section => section_vals_get_subs_vals(negf_section, "CONTACT")
243 CALL section_vals_get(contact_section, n_repetition=n_rep, explicit=is_explicit)
244 IF ((.NOT. is_explicit) .AND. do_negf) THEN
245 CALL cp_abort(__location__, &
246 "At least one contact is needed for NEGF calculation.")
247 END IF
248
249 ALLOCATE (negf_control%contacts(n_rep))
250 DO i_rep = 1, n_rep
251 region_section => section_vals_get_subs_vals(contact_section, "SCREENING_REGION", i_rep_section=i_rep)
252 CALL section_vals_get(region_section, explicit=is_explicit)
253
254 IF ((.NOT. is_explicit) .AND. do_negf) THEN
255 WRITE (contact_id_str, '(I11)') i_rep
256 CALL cp_abort(__location__, &
257 "The screening region must be defined for the contact "//trim(adjustl(contact_id_str))//".")
258 END IF
259
260 IF (is_explicit) THEN
261 CALL read_negf_atomlist(negf_control%contacts(i_rep)%atomlist_screening, region_section, 1, subsys)
262 END IF
263
264 region_section => section_vals_get_subs_vals(contact_section, "BULK_REGION", i_rep_section=i_rep)
265
266 CALL section_vals_get(region_section, explicit=is_explicit)
267
268 IF ((.NOT. is_explicit) .AND. do_negf) THEN
269 WRITE (contact_id_str, '(I11)') i_rep
270 CALL cp_abort(__location__, &
271 "The bulk region must be defined for the contact "//trim(adjustl(contact_id_str))//".")
272 END IF
273
274 IF (is_explicit) THEN
275 CALL read_negf_atomlist(negf_control%contacts(i_rep)%atomlist_bulk, region_section, 1, subsys)
276 END IF
277
278 CALL section_vals_val_get(contact_section, "FORCE_EVAL_SECTION", &
279 i_val=negf_control%contacts(i_rep)%force_env_index, &
280 i_rep_section=i_rep)
281
282 cell_section => section_vals_get_subs_vals(region_section, "CELL")
283 CALL section_vals_get(cell_section, n_repetition=n2_rep, explicit=is_explicit)
284
285 IF (((.NOT. is_explicit) .OR. n2_rep /= 2) .AND. negf_control%contacts(i_rep)%force_env_index <= 0 .AND. do_negf) THEN
286 WRITE (contact_id_str, '(I11)') i_rep
287 CALL cp_abort(__location__, &
288 "You must either provide indices of atoms belonging to two adjacent bulk unit cells "// &
289 "(BULK_REGION/CELL) for the contact, or the index of the FORCE_EVAL section (FORCE_EVAL_SECTION) "// &
290 "which will be used to construct Kohn-Sham matrix for the bulk contact "// &
291 trim(adjustl(contact_id_str))//".")
292 END IF
293
294 IF (is_explicit .AND. n2_rep > 0) THEN
295 ALLOCATE (negf_control%contacts(i_rep)%atomlist_cell(n2_rep))
296
297 DO i2_rep = 1, n2_rep
298 CALL read_negf_atomlist(negf_control%contacts(i_rep)%atomlist_cell(i2_rep)%vector, cell_section, i2_rep, subsys)
299 END DO
300 END IF
301
302 CALL section_vals_val_get(contact_section, "REFINE_FERMI_LEVEL", &
303 l_val=negf_control%contacts(i_rep)%refine_fermi_level, &
304 i_rep_section=i_rep)
305
306 CALL section_vals_val_get(contact_section, "FERMI_LEVEL", &
307 r_val=negf_control%contacts(i_rep)%fermi_level, &
308 i_rep_section=i_rep, explicit=is_explicit)
309 IF (.NOT. is_explicit) negf_control%contacts(i_rep)%refine_fermi_level = .false.
310 negf_control%contacts(i_rep)%compute_fermi_level = (.NOT. is_explicit) .OR. &
311 negf_control%contacts(i_rep)%refine_fermi_level
312
313 CALL section_vals_val_get(contact_section, "FERMI_LEVEL_SHIFTED", &
314 r_val=negf_control%contacts(i_rep)%fermi_level_shifted, &
315 i_rep_section=i_rep, explicit=is_explicit)
316 IF (is_explicit) negf_control%contacts(i_rep)%shift_fermi_level = .true.
317
318 CALL section_vals_val_get(contact_section, "TEMPERATURE", &
319 r_val=negf_control%contacts(i_rep)%temperature, &
320 i_rep_section=i_rep)
321 IF (negf_control%contacts(i_rep)%temperature <= 0.0_dp) THEN
322 CALL cp_abort(__location__, "Electronic temperature must be > 0")
323 END IF
324
325 CALL section_vals_val_get(contact_section, "ELECTRIC_POTENTIAL", &
326 r_val=negf_control%contacts(i_rep)%v_external, &
327 i_rep_section=i_rep)
328
329 subsection => section_vals_get_subs_vals(contact_section, "RESTART", i_rep_section=i_rep)
330
331 CALL section_vals_val_get(subsection, "READ_WRITE_HS", &
332 l_val=negf_control%contacts(i_rep)%read_write_HS, &
333 explicit=is_explicit)
334 IF (is_explicit) negf_control%contacts(i_rep)%read_write_HS = .true.
335
336 END DO
337
338 region_section => section_vals_get_subs_vals(negf_section, "SCATTERING_REGION")
339 CALL section_vals_get(region_section, explicit=is_explicit)
340 IF (is_explicit) THEN
341 CALL read_negf_atomlist(negf_control%atomlist_S, region_section, 1, subsys)
342 END IF
343
344 subsection => section_vals_get_subs_vals(negf_section, "SCATTERING_REGION%RESTART")
345 CALL section_vals_val_get(subsection, "READ_WRITE_HS", &
346 l_val=negf_control%read_write_HS, &
347 explicit=is_explicit)
348 IF (is_explicit) negf_control%read_write_HS = .true.
349
350 CALL section_vals_val_get(negf_section, "DISABLE_CACHE", l_val=negf_control%disable_cache)
351
352 CALL section_vals_val_get(negf_section, "EPS_DENSITY", r_val=negf_control%conv_density)
353 CALL section_vals_val_get(negf_section, "EPS_GREEN", r_val=negf_control%conv_green)
354 CALL section_vals_val_get(negf_section, "EPS_SCF", r_val=negf_control%conv_scf)
355
356 CALL section_vals_val_get(negf_section, "EPS_GEO", r_val=negf_control%eps_geometry)
357
358 CALL section_vals_val_get(negf_section, "ENERGY_LBOUND", r_val=negf_control%energy_lbound)
359 CALL section_vals_val_get(negf_section, "ETA", r_val=negf_control%eta)
360 CALL section_vals_val_get(negf_section, "HOMO_LUMO_GAP", r_val=negf_control%homo_lumo_gap)
361 CALL section_vals_val_get(negf_section, "DELTA_NPOLES", i_val=negf_control%delta_npoles)
362 CALL section_vals_val_get(negf_section, "GAMMA_KT", i_val=negf_control%gamma_kT)
363
364 CALL section_vals_val_get(negf_section, "INTEGRATION_METHOD", i_val=negf_control%integr_method)
365 CALL section_vals_val_get(negf_section, "INTEGRATION_MIN_POINTS", i_val=negf_control%integr_min_points)
366 CALL section_vals_val_get(negf_section, "INTEGRATION_MAX_POINTS", i_val=negf_control%integr_max_points)
367
368 IF (negf_control%integr_max_points < negf_control%integr_min_points) &
369 negf_control%integr_max_points = negf_control%integr_min_points
370
371 CALL section_vals_val_get(negf_section, "MAX_SCF", i_val=negf_control%max_scf)
372
373 CALL section_vals_val_get(negf_section, "NPROC_POINT", i_val=negf_control%nprocs)
374
375 CALL section_vals_val_get(negf_section, "V_SHIFT", r_val=negf_control%v_shift)
376 CALL section_vals_val_get(negf_section, "V_SHIFT_OFFSET", r_val=negf_control%v_shift_offset)
377 CALL section_vals_val_get(negf_section, "V_SHIFT_MAX_ITERS", i_val=negf_control%v_shift_maxiters)
378
379 CALL section_vals_val_get(negf_section, "SCF%UPDATE_HS", l_val=negf_control%update_HS)
380 CALL section_vals_val_get(negf_section, "SCF%RESTART_SCF", l_val=negf_control%restart_scf)
381
382 ! check consistency
383 IF (negf_control%eta < 0.0_dp) THEN
384 CALL cp_abort(__location__, "ETA must be >= 0")
385 END IF
386
387 IF (n_rep > 0) THEN
388 delta_npoles_min = nint(0.5_dp*(negf_control%eta/(pi*maxval(negf_control%contacts(:)%temperature)) + 1.0_dp))
389 ELSE
390 delta_npoles_min = 1
391 END IF
392
393 IF (negf_control%delta_npoles < delta_npoles_min) THEN
394 IF (n_rep > 0) THEN
395 eta_max = real(2*negf_control%delta_npoles - 1, kind=dp)*pi*maxval(negf_control%contacts(:)%temperature)
396 temp_current = maxval(negf_control%contacts(:)%temperature)*kelvin
397 temp_min = negf_control%eta/(pi*real(2*negf_control%delta_npoles - 1, kind=dp))*kelvin
398
399 WRITE (eta_current_str, '(ES11.4E2)') negf_control%eta
400 WRITE (eta_max_str, '(ES11.4E2)') eta_max
401 WRITE (npoles_current_str, '(I11)') negf_control%delta_npoles
402 WRITE (npoles_min_str, '(I11)') delta_npoles_min
403 WRITE (temp_current_str, '(F11.3)') temp_current
404 WRITE (temp_min_str, '(F11.3)') temp_min
405
406 CALL cp_abort(__location__, &
407 "Parameter DELTA_NPOLES must be at least "//trim(adjustl(npoles_min_str))// &
408 " (instead of "//trim(adjustl(npoles_current_str))// &
409 ") for given TEMPERATURE ("//trim(adjustl(temp_current_str))// &
410 " K) and ETA ("//trim(adjustl(eta_current_str))// &
411 "). Alternatively you can increase TEMPERATURE above "//trim(adjustl(temp_min_str))// &
412 " K, or decrease ETA below "//trim(adjustl(eta_max_str))// &
413 ". Please keep in mind that very tight ETA may result in dramatical precision loss"// &
414 " due to inversion of ill-conditioned matrices.")
415 ELSE
416 ! no leads have been defined, so calculation will abort anyway
417 negf_control%delta_npoles = delta_npoles_min
418 END IF
419 END IF
420
421 ! expand scattering region by adding atoms from contact screening regions
422 n_rep = SIZE(negf_control%contacts)
423 IF (ALLOCATED(negf_control%atomlist_S)) THEN
424 natoms_total = SIZE(negf_control%atomlist_S)
425 ELSE
426 natoms_total = 0
427 END IF
428
429 DO i_rep = 1, n_rep
430 IF (ALLOCATED(negf_control%contacts(i_rep)%atomlist_screening)) THEN
431 IF (ALLOCATED(negf_control%contacts(i_rep)%atomlist_screening)) &
432 natoms_total = natoms_total + SIZE(negf_control%contacts(i_rep)%atomlist_screening)
433 END IF
434 END DO
435
436 IF (natoms_total > 0) THEN
437 ALLOCATE (negf_control%atomlist_S_screening(natoms_total))
438 IF (ALLOCATED(negf_control%atomlist_S)) THEN
439 natoms_total = SIZE(negf_control%atomlist_S)
440 negf_control%atomlist_S_screening(1:natoms_total) = negf_control%atomlist_S(1:natoms_total)
441 ELSE
442 natoms_total = 0
443 END IF
444
445 DO i_rep = 1, n_rep
446 IF (ALLOCATED(negf_control%contacts(i_rep)%atomlist_screening)) THEN
447 natoms_current = SIZE(negf_control%contacts(i_rep)%atomlist_screening)
448
449 negf_control%atomlist_S_screening(natoms_total + 1:natoms_total + natoms_current) = &
450 negf_control%contacts(i_rep)%atomlist_screening(1:natoms_current)
451
452 natoms_total = natoms_total + natoms_current
453 END IF
454 END DO
455
456 ! sort and remove duplicated atoms
457 ALLOCATE (inds(natoms_total))
458 CALL sort(negf_control%atomlist_S_screening, natoms_total, inds)
459 DEALLOCATE (inds)
460
461 natoms_current = 1
462 DO i_rep = natoms_current + 1, natoms_total
463 IF (negf_control%atomlist_S_screening(i_rep) /= negf_control%atomlist_S_screening(natoms_current)) THEN
464 natoms_current = natoms_current + 1
465 negf_control%atomlist_S_screening(natoms_current) = negf_control%atomlist_S_screening(i_rep)
466 END IF
467 END DO
468
469 IF (natoms_current < natoms_total) THEN
470 CALL move_alloc(negf_control%atomlist_S_screening, inds)
471
472 ALLOCATE (negf_control%atomlist_S_screening(natoms_current))
473 negf_control%atomlist_S_screening(1:natoms_current) = inds(1:natoms_current)
474 DEALLOCATE (inds)
475 END IF
476 END IF
477
478 IF (do_negf .AND. SIZE(negf_control%contacts) > 2) THEN
479 CALL cp_abort(__location__, &
480 "General case (> 2 contacts) has not been implemented yet")
481 END IF
482
483 CALL timestop(handle)
484 END SUBROUTINE read_negf_control
485
486! **************************************************************************************************
487!> \brief Read region-specific list of atoms.
488!> \param atomlist list of atoms
489!> \param input_section input section which contains 'LIST' and 'MOLNAME' keywords
490!> \param i_rep_section repetition index of the input_section
491!> \param subsys subsystem environment
492! **************************************************************************************************
493 SUBROUTINE read_negf_atomlist(atomlist, input_section, i_rep_section, subsys)
494 INTEGER, ALLOCATABLE, DIMENSION(:), INTENT(out) :: atomlist
495 TYPE(section_vals_type), POINTER :: input_section
496 INTEGER, INTENT(in) :: i_rep_section
497 TYPE(cp_subsys_type), POINTER :: subsys
498
499 CHARACTER(len=*), PARAMETER :: routinen = 'read_negf_atomlist'
500
501 CHARACTER(len=default_string_length) :: index_str, natoms_str
502 CHARACTER(len=default_string_length), &
503 DIMENSION(:), POINTER :: cptr
504 INTEGER :: first_atom, handle, iatom, ikind, imol, iname, irep, last_atom, natoms_current, &
505 natoms_max, natoms_total, nkinds, nmols, nnames, nrep_list, nrep_molname
506 INTEGER, ALLOCATABLE, DIMENSION(:) :: inds
507 INTEGER, DIMENSION(:), POINTER :: iptr
508 LOGICAL :: is_list, is_molname
509 TYPE(molecule_kind_type), DIMENSION(:), POINTER :: molecule_kind_set
510 TYPE(molecule_kind_type), POINTER :: molecule_kind
511 TYPE(molecule_type), DIMENSION(:), POINTER :: molecule_set
512 TYPE(molecule_type), POINTER :: molecule
513 TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
514
515 CALL timeset(routinen, handle)
516
517 CALL cp_subsys_get(subsys, particle_set=particle_set, &
518 molecule_set=molecule_set, &
519 molecule_kind_set=molecule_kind_set)
520 natoms_max = SIZE(particle_set)
521 nkinds = SIZE(molecule_kind_set)
522
523 CALL section_vals_val_get(input_section, "LIST", i_rep_section=i_rep_section, &
524 n_rep_val=nrep_list, explicit=is_list)
525 CALL section_vals_val_get(input_section, "MOLNAME", i_rep_section=i_rep_section, &
526 n_rep_val=nrep_molname, explicit=is_molname)
527
528 ! compute the number of atoms in the NEGF region, and check the validity of given atomic indices
529 natoms_total = 0
530 IF (is_list .AND. nrep_list > 0) THEN
531 DO irep = 1, nrep_list
532 CALL section_vals_val_get(input_section, "LIST", i_rep_section=i_rep_section, i_rep_val=irep, i_vals=iptr)
533
534 natoms_current = SIZE(iptr)
535 DO iatom = 1, natoms_current
536 IF (iptr(iatom) > natoms_max) THEN
537 CALL integer_to_string(iptr(iatom), index_str)
538 CALL integer_to_string(natoms_max, natoms_str)
539 CALL cp_abort(__location__, &
540 "NEGF: Atomic index "//trim(index_str)//" given in section "// &
541 trim(input_section%section%name)//" exceeds the maximum number of atoms ("// &
542 trim(natoms_str)//").")
543 END IF
544 END DO
545
546 natoms_total = natoms_total + natoms_current
547 END DO
548 END IF
549
550 IF (is_molname .AND. nrep_molname > 0) THEN
551 DO irep = 1, nrep_molname
552 CALL section_vals_val_get(input_section, "MOLNAME", i_rep_section=i_rep_section, i_rep_val=irep, c_vals=cptr)
553 nnames = SIZE(cptr)
554
555 DO iname = 1, nnames
556 DO ikind = 1, nkinds
557 IF (molecule_kind_set(ikind)%name == cptr(iname)) EXIT
558 END DO
559
560 IF (ikind <= nkinds) THEN
561 molecule_kind => molecule_kind_set(ikind)
562 CALL get_molecule_kind(molecule_kind, nmolecule=nmols, molecule_list=iptr)
563
564 DO imol = 1, nmols
565 molecule => molecule_set(iptr(imol))
566 CALL get_molecule(molecule, first_atom=first_atom, last_atom=last_atom)
567 natoms_current = last_atom - first_atom + 1
568 natoms_total = natoms_total + natoms_current
569 END DO
570 ELSE
571 CALL cp_abort(__location__, &
572 "NEGF: A molecule with the name '"//trim(cptr(iname))//"' mentioned in section "// &
573 trim(input_section%section%name)//" has not been defined. Note that names are case sensitive.")
574 END IF
575 END DO
576 END DO
577 END IF
578
579 ! create a list of atomic indices
580 IF (natoms_total > 0) THEN
581 ALLOCATE (atomlist(natoms_total))
582
583 natoms_total = 0
584
585 IF (is_list .AND. nrep_list > 0) THEN
586 DO irep = 1, nrep_list
587 CALL section_vals_val_get(input_section, "LIST", i_rep_section=i_rep_section, i_rep_val=irep, i_vals=iptr)
588
589 natoms_current = SIZE(iptr)
590 atomlist(natoms_total + 1:natoms_total + natoms_current) = iptr(1:natoms_current)
591 natoms_total = natoms_total + natoms_current
592 END DO
593 END IF
594
595 IF (is_molname .AND. nrep_molname > 0) THEN
596 DO irep = 1, nrep_molname
597 CALL section_vals_val_get(input_section, "MOLNAME", i_rep_section=i_rep_section, i_rep_val=irep, c_vals=cptr)
598 nnames = SIZE(cptr)
599
600 DO iname = 1, nnames
601 DO ikind = 1, nkinds
602 IF (molecule_kind_set(ikind)%name == cptr(iname)) EXIT
603 END DO
604
605 IF (ikind <= nkinds) THEN
606 molecule_kind => molecule_kind_set(ikind)
607 CALL get_molecule_kind(molecule_kind, nmolecule=nmols, molecule_list=iptr)
608
609 DO imol = 1, nmols
610 molecule => molecule_set(iptr(imol))
611 CALL get_molecule(molecule, first_atom=first_atom, last_atom=last_atom)
612
613 DO natoms_current = first_atom, last_atom
614 natoms_total = natoms_total + 1
615 atomlist(natoms_total) = natoms_current
616 END DO
617 END DO
618 END IF
619 END DO
620 END DO
621 END IF
622
623 ! remove duplicated atoms
624 ALLOCATE (inds(natoms_total))
625 CALL sort(atomlist, natoms_total, inds)
626 DEALLOCATE (inds)
627
628 natoms_current = 1
629 DO iatom = natoms_current + 1, natoms_total
630 IF (atomlist(iatom) /= atomlist(natoms_current)) THEN
631 natoms_current = natoms_current + 1
632 atomlist(natoms_current) = atomlist(iatom)
633 END IF
634 END DO
635
636 IF (natoms_current < natoms_total) THEN
637 CALL move_alloc(atomlist, inds)
638
639 ALLOCATE (atomlist(natoms_current))
640 atomlist(1:natoms_current) = inds(1:natoms_current)
641 DEALLOCATE (inds)
642 END IF
643 END IF
644
645 CALL timestop(handle)
646 END SUBROUTINE read_negf_atomlist
647END MODULE negf_control_types
cp_subsys_types
types that represent a subsys, i.e. a part of the system
Definition cp_subsys_types.F:16
cp_subsys_types::cp_subsys_get
subroutine, public cp_subsys_get(subsys, ref_count, 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)
returns information about various attributes of the given subsys
Definition cp_subsys_types.F:364
input_constants
collects all constants needed in input so that they can be used without circular dependencies
Definition input_constants.F:17
input_constants::negf_run
integer, parameter, public negf_run
Definition input_constants.F:124
input_section_types
objects that represent the structure of input sections and the data contained in an input section
Definition input_section_types.F:15
input_section_types::section_vals_get_subs_vals
recursive type(section_vals_type) function, pointer, public section_vals_get_subs_vals(section_vals, subsection_name, i_rep_section, can_return_null)
returns the values of the requested subsection
Definition input_section_types.F:731
input_section_types::section_vals_get
subroutine, public section_vals_get(section_vals, ref_count, n_repetition, n_subs_vals_rep, section, explicit)
returns various attributes about the section_vals
Definition input_section_types.F:704
input_section_types::section_vals_val_get
subroutine, public section_vals_val_get(section_vals, keyword_name, i_rep_section, i_rep_val, n_rep_val, val, l_val, i_val, r_val, c_val, l_vals, i_vals, r_vals, c_vals, explicit)
returns the requested value
Definition input_section_types.F:1047
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
mathconstants
Definition of mathematical constants and functions.
Definition mathconstants.F:16
mathconstants::pi
real(kind=dp), parameter, public pi
Definition mathconstants.F:110
molecule_kind_types
Define the molecule kind structure types and the corresponding functionality.
Definition molecule_kind_types.F:16
molecule_kind_types::get_molecule_kind
subroutine, public get_molecule_kind(molecule_kind, atom_list, bond_list, bend_list, ub_list, impr_list, opbend_list, colv_list, fixd_list, g3x3_list, g4x6_list, vsite_list, torsion_list, shell_list, name, mass, charge, kind_number, natom, nbend, nbond, nub, nimpr, nopbend, nconstraint, nconstraint_fixd, nfixd, ncolv, ng3x3, ng4x6, nvsite, nfixd_restraint, ng3x3_restraint, ng4x6_restraint, nvsite_restraint, nrestraints, nmolecule, nsgf, nshell, ntorsion, molecule_list, nelectron, nelectron_alpha, nelectron_beta, bond_kind_set, bend_kind_set, ub_kind_set, impr_kind_set, opbend_kind_set, torsion_kind_set, molname_generated)
Get informations about a molecule kind.
Definition molecule_kind_types.F:512
molecule_types
Define the data structure for the molecule information.
Definition molecule_types.F:16
molecule_types::get_molecule
subroutine, public get_molecule(molecule, molecule_kind, lmi, lci, lg3x3, lg4x6, lcolv, first_atom, last_atom, first_shell, last_shell)
Get components from a molecule data set.
Definition molecule_types.F:326
negf_alloc_types
Allocatable vectors for NEGF based quantum transport calculations.
Definition negf_alloc_types.F:12
negf_control_types
Input control types for NEGF based quantum transport calculations.
Definition negf_control_types.F:12
negf_control_types::negf_control_create
subroutine, public negf_control_create(negf_control)
allocate control options for Non-equilibrium Green's Function calculation
Definition negf_control_types.F:151
negf_control_types::read_negf_control
subroutine, public read_negf_control(negf_control, input, subsys)
Read NEGF input parameters.
Definition negf_control_types.F:217
negf_control_types::negf_control_release
subroutine, public negf_control_release(negf_control)
release memory allocated for NEGF control options
Definition negf_control_types.F:172
particle_types
Define the data structure for the particle information.
Definition particle_types.F:19
physcon
Definition of physical constants:
Definition physcon.F:68
physcon::kelvin
real(kind=dp), parameter, public kelvin
Definition physcon.F:165
string_utilities
Utilities for string manipulations.
Definition string_utilities.F:16
string_utilities::integer_to_string
subroutine, public integer_to_string(inumber, string)
Converts an integer number to a string. The WRITE statement will return an error message,...
Definition string_utilities.F:3239
util
All kind of helpful little routines.
Definition util.F:14
cp_subsys_types::cp_subsys_type
represents a system: atoms, molecules, their pos,vel,...
Definition cp_subsys_types.F:89
input_section_types::section_vals_type
stores the values of a section
Definition input_section_types.F:127
molecule_kind_types::molecule_kind_type
Definition molecule_kind_types.F:151
molecule_types::molecule_type
Definition molecule_types.F:159
negf_alloc_types::negf_allocatable_ivector
Allocatable 1-D integer vector.
Definition negf_alloc_types.F:27
negf_control_types::negf_control_contact_type
Input parameters related to a single contact.
Definition negf_control_types.F:47
negf_control_types::negf_control_type
Input parameters related to the NEGF run.
Definition negf_control_types.F:79
particle_types::particle_type
Definition particle_types.F:35