Branch data Line data Source code
1 : : !!****m* ABINIT/m_hamiltonian
2 : : !! NAME
3 : : !! m_hamiltonian
4 : : !!
5 : : !! FUNCTION
6 : : !! This module provides the definition of the gs_hamiltonian_type and of the rf_hamiltonian_type
7 : : !! datastructures used in the "getghc" and "getgh1c" routines to apply the Hamiltonian (or
8 : : !! its derivative) on a wavefunction. Methods to initialize or destroy the objects are defined here.
9 : : !!
10 : : !! TODO
11 : : !! All array pointers in H datatypes should be declared as contiguous for efficiency reasons.
12 : : !! (well, here performance is critical). Client code should make sure they always point contiguous targets.
13 : : !!
14 : : !! COPYRIGHT
15 : : !! Copyright (C) 2009-2026 ABINIT group (MG, MT)
16 : : !! This file is distributed under the terms of the
17 : : !! GNU General Public License, see ~abinit/COPYING
18 : : !! or http://www.gnu.org/copyleft/gpl.txt .
19 : : !!
20 : : !! SOURCE
21 : :
22 : : #if defined HAVE_CONFIG_H
23 : : #include "config.h"
24 : : #endif
25 : :
26 : : #include "abi_common.h"
27 : :
28 : : module m_hamiltonian
29 : :
30 : : use iso_fortran_env, only : int32, int64, real32, real64
31 : :
32 : : use defs_basis
33 : : use m_abicore
34 : : use m_errors
35 : : use m_xmpi
36 : :
37 : : use m_fstrings, only : sjoin, itoa, ktoa, yesno
38 : : use defs_datatypes, only : pseudopotential_type
39 : : use defs_abitypes, only : MPI_type
40 : : use m_dtset, only : dataset_type
41 : : use m_copy, only : addr_copy
42 : : use m_crystal, only : crystal_t
43 : : use m_geometry, only : metric
44 : : use m_pawtab, only : pawtab_type
45 : : use m_pawfgr, only : pawfgr_type
46 : : use m_fftcore, only : sphereboundary
47 : : use m_fft, only : fftpac
48 : : use m_fourier_interpol, only : transgrid
49 : : use m_pawcprj, only : pawcprj_getdim
50 : : use m_paw_ij, only : paw_ij_type
51 : : use m_paral_atom, only : get_my_atmtab, free_my_atmtab
52 : : use m_electronpositron, only : electronpositron_type, electronpositron_calctype
53 : : use m_kg, only : ph1d3d, getph, mkkin, mkkpg
54 : : use m_fock, only : fock_common_type, fock_BZ_type, fock_ACE_type, fock_type
55 : : use m_mkffnl, only : mkffnl_objs
56 : : use m_initylmg, only : initylmg_k
57 : : use m_ompgpu_fourwf, only : free_ompgpu_fourwf
58 : :
59 : : #if defined HAVE_GPU_CUDA
60 : : use m_manage_cuda
61 : : #endif
62 : :
63 : : #if defined HAVE_FC_ISO_C_BINDING
64 : : use, intrinsic :: iso_c_binding, only : c_ptr,c_loc,c_f_pointer,c_int32_t,c_int64_t,c_size_t
65 : : #endif
66 : :
67 : : #if defined HAVE_GPU && defined HAVE_YAKL
68 : : use gator_mod
69 : : #endif
70 : :
71 : : implicit none
72 : :
73 : : private
74 : :
75 : : public :: pawdij2ekb
76 : : public :: pawdij2e1kb
77 : : public :: gspot_transgrid_and_pack ! Set up local potential vlocal on the coarse FFT mesh from vtrial on the fine mesh.
78 : :
79 : : ! These constants select how H_{k',k} is applied in reciprocal space
80 : : integer,parameter,public :: KPRIME_H_K=1, K_H_KPRIME=2, K_H_K=3, KPRIME_H_KPRIME=4
81 : : !!***
82 : :
83 : : !----------------------------------------------------------------------
84 : :
85 : : !!****t* m_hamiltonian/gs_hamiltonian_type
86 : : !! NAME
87 : : !! gs_hamiltonian_type
88 : : !!
89 : : !! FUNCTION
90 : : !! This datastructure contains the information about one Hamiltonian,
91 : : !! needed in the "getghc" routine, that applies the Hamiltonian on a wavefunction.
92 : : !! The Hamiltonian is expressed in reciprocal space:
93 : : !!
94 : : !! H_k^prime,k = exp(-i.k^prime.r^prime) H exp(i.k.r)
95 : : !!
96 : : !! In most cases k = k^prime and the k^prime objects are simply pointers to k objects.
97 : : !!
98 : : !! SOURCE
99 : :
100 : : type,public :: gs_hamiltonian_type
101 : :
102 : : ! ===== Integer scalars
103 : :
104 : : integer :: dimekb1 = -1
105 : : ! First dimension of Ekb
106 : : ! Same as psps%dimekb
107 : : ! ->Norm conserving : Max. number of Kleinman-Bylander energies
108 : : ! for each atom type
109 : : ! dimekb1=lnmax
110 : : ! ->PAW : Max. number of Dij coefficients connecting projectors
111 : : ! for each atom
112 : : ! dimekb1=cplex_dij*lmnmax*(lmnmax+1)/2
113 : :
114 : : integer :: dimekb2 = -1
115 : : ! Second dimension of Ekb
116 : : ! ->Norm conserving psps: dimekb2=ntypat
117 : : ! ->PAW : dimekb2=natom
118 : :
119 : : integer :: dimekbq = -1
120 : : ! Fourth dimension of Ekb
121 : : ! 2 if Ekb factors contain a exp(-iqR) phase, 1 otherwise
122 : :
123 : : integer :: ispin_gbt = -1
124 : : ! 1, GBT spin up channel (uses k - q/2)
125 : : ! 2, GBT spin down channel (uses k + q/2)
126 : :
127 : : integer :: istwf_k = -1
128 : : ! option parameter that describes the storage of wfs at k
129 : :
130 : : integer :: istwf_kp = -1
131 : : ! option parameter that describes the storage of wfs at k^prime
132 : :
133 : : integer :: lmnmax = -1
134 : : ! Maximum number of different l,m,n components over all types of psps.
135 : : ! same as dtset%lmnmax
136 : :
137 : : integer :: matblk = -1
138 : : ! dimension of the array ph3d
139 : :
140 : : integer :: mgfft = -1
141 : : ! maximum size for 1D FFTs (same as dtset%mgfft)
142 : :
143 : : integer :: mpsang = -1
144 : : ! Highest angular momentum of non-local projectors over all type of psps.
145 : : ! shifted by 1 : for all local psps, mpsang=0; for largest s, mpsang=1,
146 : : ! for largest p, mpsang=2; for largest d, mpsang=3; for largest f, mpsang=4
147 : : ! This gives also the number of non-local "channels"
148 : : ! same as psps%mpsang
149 : :
150 : : integer :: mpssoang = -1
151 : : ! Maximum number of channels, including those for treating the spin-orbit coupling
152 : : ! For NC pseudopotentials only:
153 : : ! when mpspso=1, mpssoang=mpsang
154 : : ! when mpspso=2, mpssoang=2*mpsang-1
155 : : ! For PAW: same as mpsang
156 : : ! same as psps%mpssoang
157 : :
158 : : integer :: natom = -1
159 : : ! The number of atoms for this dataset; same as dtset%natom
160 : :
161 : : integer :: nfft = -1
162 : : ! number of FFT grid points same as dtset%nfft
163 : :
164 : : integer :: nfft_blocks = 1
165 : : ! number of blocks fourwf computation is divided into
166 : :
167 : : integer :: npw_k = -1
168 : : ! number of plane waves at k
169 : : ! In case of band-FFT parallelism, npw_k is the number of plane waves
170 : : ! processed by current proc
171 : :
172 : : integer :: npw_fft_k = -1
173 : : ! number of plane waves at k used to apply Hamiltonian when band-FFT
174 : : ! parallelism is activated (i.e. data are distributed in the "FFT" configuration)
175 : :
176 : : integer :: npw_kp = -1
177 : : ! number of plane waves at k^prime
178 : : ! In case of band-FFT parallelism, npw_kp is the number of plane waves
179 : : ! processed by current proc
180 : :
181 : : integer :: npw_fft_kp = -1
182 : : ! number of plane waves at k^prime used to apply Hamiltonian when band-FFT
183 : : ! parallelism is activated (i.e. data are distributed in the "FFT" configuration)
184 : :
185 : : integer :: nspinor = -1
186 : : ! Number of spinorial components
187 : :
188 : : integer :: nsppol = -1
189 : : ! Total number of spin components (1=non-polarized, 2=polarized)
190 : :
191 : : integer :: ntypat = -1
192 : : ! Number of types of pseudopotentials same as dtset%ntypat
193 : :
194 : : integer :: nvloc = -1
195 : : ! Number of components of vloc
196 : : ! usually, nvloc=1, except in the non-collinear magnetism case, where nvloc=4
197 : :
198 : : integer :: n4 = -1, n5 = -1, n6 = -1
199 : : ! same as ngfft(4:6)
200 : :
201 : : integer :: gpu_option = -1
202 : : ! Governs the choice of the GPU implementation:
203 : : ! = 0 ==> do not use GPU
204 : : ! > 0 ==> see defs_basis.F90 to have the list of possible GPU implementations
205 : :
206 : : integer :: usecprj = -1
207 : : ! usecprj= 1 if cprj projected WF are stored in memory
208 : : ! = 0 if they are to be computed on the fly
209 : :
210 : : integer :: usepaw = -1
211 : : ! if usepaw=0 , use norm-conserving psps part of the code
212 : : ! is usepaw=1 , use paw part of the code
213 : :
214 : : integer :: useylm = -1
215 : : ! governs the way the nonlocal operator is to be applied:
216 : : ! 1=using Ylm, 0=using Legendre polynomials
217 : :
218 : : integer :: use_gbt = 0
219 : : ! 0, use normal non-collinear calculation
220 : : ! 1, use spin spiral calculation
221 : : ! 2, use spin spiral with z-component SOC
222 : :
223 : : integer :: zora = 0
224 : : ! zora=0: no zora terms. zora=1: use available zora terms
225 : : ! currently this is limited to nuclear dipole moment terms,
226 : : ! although pawspnorb should eventually be included here as well
227 : :
228 : : ! ===== Integer arrays
229 : :
230 : : #if defined HAVE_GPU && defined HAVE_YAKL
231 : : integer(c_int32_t), contiguous, pointer :: atindx(:) => null()
232 : : #else
233 : : integer, allocatable :: atindx(:)
234 : : #endif
235 : : ! atindx(natom)
236 : : ! index table for atoms (see gstate.f)
237 : :
238 : : #if defined HAVE_GPU && defined HAVE_YAKL
239 : : integer(c_int32_t), contiguous, pointer :: atindx1(:) => null()
240 : : #else
241 : : integer, allocatable :: atindx1(:)
242 : : #endif
243 : : ! atindx1(natom)
244 : : ! index table for atoms, inverse of atindx (see gstate.f)
245 : :
246 : : integer, allocatable :: dimcprj(:)
247 : : ! dimcprj(natom*usepaw)=dimensions of array cprj
248 : : ! dimcprj(ia)=cprj(ia,:)%nlmn
249 : : ! atoms are ordered by atom-type
250 : :
251 : : integer, allocatable :: gbound_k(:,:)
252 : : ! gbound_k(2*mgfft+8,2)
253 : : ! G sphere boundary, for each plane wave at k
254 : :
255 : : #if defined HAVE_GPU && defined HAVE_YAKL
256 : : integer(c_int32_t), contiguous, pointer :: indlmn(:,:,:) => null()
257 : : #else
258 : : integer(c_int32_t), allocatable :: indlmn(:,:,:)
259 : : #endif
260 : : ! indlmn(6,lmnmax,ntypat)
261 : : ! For each type of psp,
262 : : ! array giving l,m,n,lm,ln,spin for i=ln (if useylm=0)
263 : : ! or i=lmn (if useylm=1)
264 : :
265 : : #if defined HAVE_GPU && defined HAVE_YAKL
266 : : integer(c_int32_t), contiguous, pointer :: nattyp(:) => null()
267 : : #else
268 : : integer, allocatable :: nattyp(:)
269 : : #endif
270 : : ! nattyp(ntypat)
271 : : ! # of atoms of each type
272 : :
273 : : integer :: ngfft(18)
274 : : ! ngfft(1:3)=integer fft box dimensions
275 : : ! ngfft(4:6)=integer fft box dimensions, might be augmented for CPU speed
276 : : ! ngfft(7)=fftalg
277 : : ! ngfft(8)=fftalg
278 : :
279 : : integer :: nloalg(3)
280 : : ! governs the choice of the algorithm for non-local operator same as dtset%nloalg
281 : :
282 : : integer, allocatable :: pspso(:)
283 : : ! pspso(ntypat)
284 : : ! For each type of psp, 1 if no spin-orbit component is taken
285 : : ! into account, 2 if a spin-orbit component is used
286 : : ! Relevant for NC-psps and PAW.
287 : :
288 : : #if defined HAVE_GPU && defined HAVE_YAKL
289 : : integer(c_int32_t), contiguous, pointer :: typat(:) => null()
290 : : #else
291 : : integer, allocatable :: typat(:)
292 : : #endif
293 : : ! typat(natom)
294 : : ! type of each atom
295 : :
296 : : ! integer, allocatable :: indpw_k(:,:)
297 : : ! indpw_k(4,npw_fft_k)
298 : : ! array which gives fft box index for given basis sphere
299 : :
300 : : ! Integer pointers
301 : :
302 : : integer, contiguous, pointer :: gbound_kp(:,:) => null()
303 : : ! gbound_kp(2*mgfft+8,2)
304 : : ! G sphere boundary, for each plane wave at k^prime
305 : :
306 : : #if defined HAVE_GPU && defined HAVE_YAKL
307 : : integer(int32), contiguous, pointer :: kg_k(:,:) => null()
308 : : #else
309 : : integer, contiguous, pointer :: kg_k(:,:) => null()
310 : : #endif
311 : : ! kg_k(3,npw_fft_k)
312 : : ! G vector coordinates with respect to reciprocal lattice translations
313 : : ! at k
314 : :
315 : : integer, contiguous, pointer :: kg_kp(:,:) => null()
316 : : ! kg_kp(3,npw_fft_kp)
317 : : ! G vector coordinates with respect to reciprocal lattice translations
318 : : ! at k^prime
319 : :
320 : : ! ===== Real scalars
321 : :
322 : : real(dp) :: ucvol = -one
323 : : ! unit cell volume (Bohr**3)
324 : :
325 : : ! ===== Real arrays
326 : :
327 : : real(dp), allocatable :: ekb_spin(:,:,:,:,:)
328 : : ! ekb_spin(dimekb1,dimekb2,nspinor**2,dimekbq,my_nsppol)
329 : : ! Contains the values of ekb array for all spins treated by current process
330 : : ! See ekb description ; ekb is pointer to ekb_spin(:,:,:,:,my_isppol)
331 : :
332 : : real(dp), allocatable :: sij(:,:)
333 : : ! sij(dimekb1,ntypat*usepaw) = overlap matrix for paw calculation
334 : :
335 : : real(dp) :: gmet(3,3)
336 : : ! reciprocal space metric tensor in Bohr**-2
337 : :
338 : : real(dp) :: gprimd(3,3)
339 : : ! dimensional reciprocal space primitive translations (Bohr^-1)
340 : :
341 : : real(dp) :: kpt_k(3)
342 : : ! dimensionless k point coordinates wrt reciprocal lattice vectors
343 : :
344 : : real(dp) :: kpt_kp(3)
345 : : ! dimensionless k^prime point coordinates wrt reciprocal lattice vectors
346 : :
347 : : real(dp) :: spinaxis(3)
348 : : ! spin quantization axis
349 : :
350 : : real(dp), allocatable :: nucdipmom(:,:)
351 : : ! nucdipmom(3,natom)
352 : : ! nuclear dipole moments at each atomic position
353 : :
354 : : #if defined HAVE_GPU && defined HAVE_YAKL
355 : : real(c_double), contiguous, pointer :: ph1d(:,:) => null()
356 : : #else
357 : : real(dp), allocatable :: ph1d(:,:)
358 : : #endif
359 : : ! ph1d(2,3*(2*mgfft+1)*natom)
360 : : ! 1-dim phase arrays for structure factor (see getph.f).
361 : :
362 : : real(dp), allocatable :: phkxred(:,:)
363 : : ! phkxred(2,natom)
364 : : ! phase factors exp(2 pi k.xred) at k
365 : :
366 : : ! ===== Real pointers
367 : :
368 : : real(dp), contiguous, pointer :: ekb(:,:,:,:) => null()
369 : : ! ekb(dimekb1,dimekb2,nspinor**2,dimekbq)
370 : : ! ->Norm conserving : (Real) Kleinman-Bylander energies (hartree)
371 : : ! for number of basis functions (l,n) (lnmax)
372 : : ! and number of atom types (ntypat)
373 : : ! dimekb1=lnmax ; dimekb2=ntypat ; dimekbq=1
374 : : ! ->PAW : (Real, symmetric) Frozen part of Dij coefficients
375 : : ! to connect projectors
376 : : ! for number of basis functions (l,m,n) (lmnmax)
377 : : ! and number of atom (natom)
378 : : ! dimekb1=lmnmax*(lmnmax+1)/2 ; dimekb2=natom ; dimekbq=1
379 : : ! ekb is spin dependent in the case of PAW calculations.
380 : : ! For each spin component, ekb points to ekb_spin(:,:,:,:,my_isppol)
381 : : ! dimekbq=2 if Ekb factors contain a exp(-iqR) phase, dimekbq=1 otherwise
382 : : ! About the non-local factors symmetry:
383 : : ! - The lower triangular part of the Dij matrix can be deduced from the upper one
384 : : ! with the following relation: D^s2s1_ji = (D^s1s2_ij)^*
385 : : ! where s1,s2 are spinor components
386 : :
387 : : real(dp), contiguous, pointer :: ffnl_k(:,:,:,:) => null()
388 : : ! ffnl_k(npw_fft_k,2,dimffnl_k,ntypat)
389 : : ! nonlocal form factors at k
390 : :
391 : : real(dp), contiguous, pointer :: ffnl_kp(:,:,:,:) => null()
392 : : ! ffnl_kp(npw_fft_kp,2,dimffnl_kp,ntypat)
393 : : ! nonlocal form factors at k_prime
394 : :
395 : : real(dp), contiguous, pointer :: kinpw_k(:) => null()
396 : : ! kinpw_k(npw_fft_k)
397 : : ! (modified) kinetic energy for each plane wave at k
398 : : ! CAVEAT: In band mode, this array is NOT EQUIVALENT to kinpw(npw_k)
399 : :
400 : : real(dp), contiguous, pointer :: kinpw_kp(:) => null()
401 : : ! kinpw_kp(npw_fft_kp)
402 : : ! (modified) kinetic energy for each plane wave at k^prime
403 : :
404 : : real(dp), contiguous, pointer :: kpg_k(:,:) => null()
405 : : ! kpg_k(3,npw_fft_k)
406 : : ! k+G vector coordinates at k
407 : :
408 : : real(dp), contiguous, pointer :: kpg_kp(:,:) => null()
409 : : ! kpg_kp(3,npw_fft_kp)
410 : : ! k^prime+G vector coordinates at k^prime
411 : :
412 : : real(dp), contiguous, pointer :: phkpxred(:,:) => null()
413 : : ! phkpxred(2,natom)
414 : : ! phase factors exp(2 pi k^prime.xred) at k^prime
415 : :
416 : : real(dp), contiguous, pointer :: ph3d_k(:,:,:) => null()
417 : : ! ph3d_k(2,npw_fft_k,matblk)
418 : : ! 3-dim structure factors, for each atom and plane wave at k
419 : :
420 : : real(dp), contiguous, pointer :: ph3d_kp(:,:,:) => null()
421 : : ! ph3d_kp(2,npw_fft_kp,matblk)
422 : : ! 3-dim structure factors, for each atom and plane wave at k^prime
423 : :
424 : : real(dp), contiguous, pointer :: vectornd(:,:,:,:,:) => null()
425 : : ! vectornd(n4,n5,n6,nvloc,3)
426 : : ! vector potential of nuclear magnetic dipoles
427 : : ! in real space, on the augmented fft grid
428 : :
429 : : real(dp), contiguous, pointer :: vlocal(:,:,:,:) => null()
430 : : ! vlocal(n4,n5,n6,nvloc)
431 : : ! local potential in real space, on the augmented fft grid
432 : :
433 : : real(dp), contiguous, pointer :: vxctaulocal(:,:,:,:,:) => null()
434 : : ! vxctaulocal(n4,n5,n6,nvloc,4)
435 : : ! derivative of XC energy density with respect to kinetic energy density,
436 : : ! in real space, on the augmented fft grid
437 : :
438 : : real(dp), contiguous, pointer :: xred(:,:) => null()
439 : : ! xred(3,natom)
440 : : ! reduced coordinates of atoms (dimensionless)
441 : :
442 : : real(dp), allocatable :: fofr_work(:,:,:,:)
443 : : ! (2,n4,n5,n6,ndat)
444 : : ! Buffer used in getgh1c when calling fourwf to compute <r|vlocal1|u_nk> for n =1, ndat.
445 : : ! It is automatically allocated/reallocated by alloc_fofr according on the input ndat.
446 : : ! In the case of gpu_option = 2, this is the buffer that is mapped to the GPU but only when ndat changes.
447 : :
448 : : ! ===== Structured datatype pointers
449 : :
450 : : type(fock_common_type), pointer :: fockcommon => null()
451 : : ! common quantities needed to calculate Fock exact exchange
452 : :
453 : : type(fock_BZ_type), pointer :: fockbz => null()
454 : : ! total brillouin zone quantities needed to calculate Fock exact exchange
455 : :
456 : : type(fock_ACE_type), pointer :: fockACE_k => null()
457 : : ! ACE quantities needed to calculate Fock exact exchange in the ACE context
458 : :
459 : : contains
460 : :
461 : : procedure :: init => gsham_init
462 : : ! Initialize the GS Hamiltonian
463 : :
464 : : procedure :: free => gsham_free
465 : : ! Free the memory in the GS Hamiltonian
466 : :
467 : : procedure :: load_spin => gsham_load_spin
468 : : ! Setup of the spin-dependent part of the GS Hamiltonian
469 : :
470 : : procedure :: load_k => gsham_load_k
471 : : ! Setup of the k-dependent part of the GS Hamiltonian
472 : :
473 : : procedure :: load_kprime => gsham_load_kprime
474 : : ! Setup of the k^prime-dependent part of the GS Hamiltonian
475 : :
476 : : procedure :: eph_setup_k => gsham_eph_setup_k
477 : : ! Simplified interface to load either k or kprime in the case of e-ph calculations.
478 : :
479 : : procedure :: copy => gsham_copy
480 : : ! Copy the object
481 : :
482 : : procedure :: print => gsham_print
483 : : ! Print the object
484 : :
485 : : procedure :: alloc_fofr_work => gsham_alloc_fofr_work
486 : : ! Allocate work space array before calling fourwf for ndat bands and map it to GPU
487 : :
488 : : end type gs_hamiltonian_type
489 : : !!***
490 : :
491 : : !----------------------------------------------------------------------
492 : :
493 : : !!****t* m_hamiltonian/rf_hamiltonian_type
494 : : !! NAME
495 : : !! rf_hamiltonian_type
496 : : !!
497 : : !! FUNCTION
498 : : !! This datastructure contains few data about one 1st-order Hamiltonian,
499 : : !! needed in the "getgh1c" routine, that apply the 1st-order Hamiltonian
500 : : !! on a wavefunction.
501 : : !!
502 : : !! SOURCE
503 : :
504 : : type,public :: rf_hamiltonian_type
505 : :
506 : : ! ===== Integer scalars
507 : :
508 : : integer :: cplex = -1
509 : : ! if 1, real space 1-order functions on FFT grid are REAL; if 2, COMPLEX
510 : :
511 : : integer :: dime1kb1 = -1
512 : : ! First dimension of E1kb, derivative of Ekb with respect to a perturbation
513 : :
514 : : integer :: dime1kb2 = -1
515 : : ! Second dimension of E1kb, derivative of Ekb with respect to a perturbation
516 : : ! NCPP: dime1kb2=ntypat, PAW: dime1kb2=natom
517 : :
518 : : integer :: npw_k = -1
519 : : ! number of plane waves at k
520 : :
521 : : integer :: npw_kp = -1
522 : : ! number of plane waves at k^prime
523 : :
524 : : integer:: nspinor = -1
525 : : ! Number of spinorial components
526 : :
527 : : integer :: nsppol = -1
528 : : ! Total number of spin components (1=non-polarized, 2=polarized)
529 : :
530 : : integer :: nvloc = -1
531 : : ! Number of components of vloc
532 : : ! usually, nvloc=1, except in the non-collinear magnetism case, where nvloc=4
533 : :
534 : : integer :: n4 = -1, n5 = -1, n6 = -1
535 : : ! same as ngfft(4:6)
536 : :
537 : : integer :: zora = 0
538 : : ! determines zora term use. 0 means nonrelativistic, 1 means use zora terms
539 : : ! where available (currently only nuclear dipole moment terms)
540 : :
541 : : ! ===== Real arrays
542 : :
543 : : real(dp), allocatable :: e1kbfr_spin(:,:,:,:,:)
544 : : ! e1kbfr_spin(dimekb1,dimekb2,nspinor**2,cplex,my_nsppol)
545 : : ! Contains the values of e1kbfr array for all spins treated by current process
546 : : ! See e1kbfr description; e1kbfr is pointer to e1kbfr_spin(:,:,:,:,isppol)
547 : :
548 : : real(dp), allocatable :: e1kbsc_spin(:,:,:,:,:)
549 : : ! e1kbsc_spin(dimekb1,dimekb2,nspinor**2,cplex,my_nsppol)
550 : : ! Contains the values of e1kbsc array for all spins treated by current process
551 : : ! See e1kbsc description; e1kbsc is pointer to e1kbsc_spin(:,:,:,:,isppol)
552 : :
553 : : ! ===== Real pointers
554 : :
555 : : real(dp), contiguous, pointer :: dkinpw_k(:) => null()
556 : : ! dkinpw_k(npw_k)
557 : : ! 1st derivative of the (modified) kinetic energy for each plane wave at k
558 : :
559 : : real(dp), contiguous, pointer :: dkinpw_kp(:) => null()
560 : : ! dkinpw_kp(npw_kp)
561 : : ! 1st derivative of the (modified) kinetic energy for each plane wave at k^prime
562 : :
563 : : real(dp), contiguous, pointer :: ddkinpw_k(:) => null()
564 : : ! ddkinpw_k(npw_k)
565 : : ! 2nd derivative of the (modified) kinetic energy for each plane wave at k
566 : :
567 : : real(dp), contiguous, pointer :: ddkinpw_kp(:) => null()
568 : : ! ddkinpw_kp(npw_kp)
569 : : ! 2nd derivative of the (modified) kinetic energy for each plane wave at k^prime
570 : :
571 : : real(dp), contiguous, pointer :: e1kbfr(:,:,:,:) => null()
572 : : ! Frozen part of 1st derivative of ekb for the considered perturbation
573 : : ! (part not depending on VHxc^(1))
574 : : ! e1kbfr(dime1kb1,dime1kb2,nspinor**2,cplex)
575 : : ! For each spin component, e1kbfr points to e1kbfr_spin(:,:,:,:,my_isppol)
576 : :
577 : : real(dp), contiguous, pointer :: e1kbsc(:,:,:,:) => null()
578 : : ! Self-consistent 1st derivative of ekb for the considered perturbation
579 : : ! (part depending only on self-consistent VHxc^(1))
580 : : ! e1kbsc(dime1kb1,dime1kb2,nspinor**2,cplex)
581 : : ! For each spin component, e1kbfr points to e1kbfr_spin(:,:,:,:,my_isppol)
582 : :
583 : : real(dp), contiguous, pointer :: vectornd(:,:,:,:) => null()
584 : : ! vectornd(n4,n5,n6,nvloc)
585 : : ! vector potential of nuclear magnetic dipoles
586 : : ! in real space, on the augmented fft grid, in direction idir
587 : : ! (the ddk pert direction)
588 : :
589 : : real(dp), contiguous, pointer :: vlocal1(:,:,:,:) => null()
590 : : ! vlocal1(cplex*n4,n5,n6,nvloc)
591 : : ! 1st-order local potential in real space, on the augmented fft grid
592 : :
593 : : real(dp), contiguous, pointer :: vxctaulocal(:,:,:,:,:) => null()
594 : : ! vxctaulocal(n4,n5,n6,nvloc,4)
595 : : ! derivative of XC energy density with respect to kinetic energy density,
596 : : ! in real space, on the augmented fft grid
597 : :
598 : : contains
599 : : procedure :: init => rfham_init
600 : : ! Initialize the RF Hamiltonian
601 : :
602 : : procedure :: free => rfham_free
603 : : ! Free the memory in the RF Hamiltonian
604 : :
605 : : procedure :: load_spin => rfham_load_spin
606 : : ! Setup of the spin-dependent part of the RF Hamiltonian.
607 : :
608 : : procedure :: load_k => rfham_load_k
609 : : ! Setup of the k-dependent part of the RF Hamiltonian
610 : :
611 : : end type rf_hamiltonian_type
612 : : !!***
613 : :
614 : : contains !===========================================================
615 : :
616 : : !----------------------------------------------------------------------
617 : :
618 : : !!****f* m_hamiltonian/gsham_free
619 : : !! NAME
620 : : !! gsham_free
621 : : !!
622 : : !! FUNCTION
623 : : !! Clean and destroy gs_hamiltonian_type datastructure
624 : : !!
625 : : !! SOURCE
626 : :
627 : 79605 : subroutine gsham_free(Ham)
628 : :
629 : : !Arguments ------------------------------------
630 : : class(gs_hamiltonian_type),target,intent(inout) :: Ham
631 : :
632 : : !Local variables-------------------------------
633 : 79605 : real(dp), contiguous, pointer :: fofr_work_ptr(:,:,:,:)
634 : : ! *************************************************************************
635 : :
636 : : DBG_ENTER("COLL")
637 : :
638 : : !@gs_hamiltonian_type
639 : :
640 : : ! Integer Pointers
641 : 79605 : if (associated(Ham%gbound_kp,Ham%gbound_k)) then
642 : 59239 : nullify(Ham%gbound_kp)
643 [ + + ]: 20366 : else if (associated(Ham%gbound_kp)) then
644 : 19969 : ABI_FREE(Ham%gbound_kp)
645 : : end if
646 : :
647 : : ! Integer arrays
648 [ + - ]: 79605 : if(Ham%gpu_option == ABI_GPU_KOKKOS) then
649 : : #if defined HAVE_GPU && defined HAVE_YAKL
650 : : ABI_SFREE_MANAGED(Ham%atindx)
651 : : ABI_SFREE_MANAGED(Ham%atindx1)
652 : : ABI_SFREE_MANAGED(Ham%typat)
653 : : ABI_SFREE_MANAGED(Ham%indlmn)
654 : : ABI_SFREE_MANAGED(Ham%nattyp)
655 : : #endif
656 : : else
657 [ - + ]: 79605 : ABI_FREE(Ham%atindx)
658 [ - + ]: 79605 : ABI_FREE(Ham%atindx1)
659 [ - + ]: 79605 : ABI_FREE(Ham%typat)
660 [ - + ]: 79605 : ABI_FREE(Ham%indlmn)
661 [ - + ]: 79605 : ABI_FREE(Ham%nattyp)
662 : : end if
663 [ + + ]: 79605 : ABI_SFREE(Ham%gbound_k)
664 [ + - ]: 79605 : ABI_SFREE(Ham%pspso)
665 [ + - ]: 79605 : ABI_SFREE(Ham%dimcprj)
666 : :
667 : : ! Real Pointers
668 [ + + ]: 79605 : if (associated(Ham%phkpxred,Ham%phkxred)) then
669 : 56038 : nullify(Ham%phkpxred)
670 [ + + ]: 23567 : else if (associated(Ham%phkpxred)) then
671 : 23170 : ABI_FREE(Ham%phkpxred)
672 : : end if
673 [ + + ]: 79605 : ABI_SFREE(Ham%phkxred)
674 [ + + ]: 79605 : if (associated(Ham%ekb)) nullify(Ham%ekb)
675 [ + + ]: 79605 : if (associated(Ham%vectornd)) nullify(Ham%vectornd)
676 [ + + ]: 79605 : if (associated(Ham%vlocal)) nullify(Ham%vlocal)
677 [ + + ]: 79605 : if (associated(Ham%vxctaulocal)) nullify(Ham%vxctaulocal)
678 [ + - ]: 79605 : if (associated(Ham%xred)) nullify(Ham%xred)
679 [ + + ]: 79605 : if (associated(Ham%kinpw_k)) nullify(Ham%kinpw_k)
680 [ + + ]: 79605 : if (associated(Ham%kinpw_kp)) nullify(Ham%kinpw_kp)
681 [ + + ]: 79605 : if (associated(Ham%kg_k)) nullify(Ham%kg_k)
682 [ + + ]: 79605 : if (associated(Ham%kg_kp)) nullify(Ham%kg_kp)
683 [ + + ]: 79605 : if (associated(Ham%kpg_k)) nullify(Ham%kpg_k)
684 [ + + ]: 79605 : if (associated(Ham%kpg_kp)) nullify(Ham%kpg_kp)
685 [ + + ]: 79605 : if (associated(Ham%ffnl_k)) nullify(Ham%ffnl_k)
686 [ + + ]: 79605 : if (associated(Ham%ffnl_kp)) nullify(Ham%ffnl_kp)
687 [ + + ]: 79605 : if (associated(Ham%ph3d_k)) nullify(Ham%ph3d_k)
688 [ + + ]: 79605 : if (associated(Ham%ph3d_kp)) nullify(Ham%ph3d_kp)
689 : :
690 : :
691 : : ! Real arrays
692 [ + - ]: 79605 : ABI_SFREE(Ham%ekb_spin)
693 [ + - ]: 79605 : ABI_SFREE(Ham%sij)
694 [ + - ]: 79605 : ABI_SFREE(Ham%nucdipmom)
695 [ + - ]: 79605 : if(Ham%gpu_option == ABI_GPU_KOKKOS) then
696 : : #if defined HAVE_GPU && defined HAVE_YAKL
697 : : ABI_SFREE_MANAGED(Ham%ph1d)
698 : : #endif
699 : : else
700 [ - + ]: 79605 : ABI_FREE(Ham%ph1d)
701 : : end if
702 : :
703 [ - + ]: 79605 : if(Ham%gpu_option==ABI_GPU_OPENMP) then
704 : 0 : call free_ompgpu_fourwf()
705 : : end if
706 : :
707 : : ! Structured datatype pointers
708 [ + + ]: 79605 : if (associated(Ham%fockcommon)) nullify(Ham%fockcommon)
709 [ + + ]: 79605 : if (associated(Ham%fockACE_k)) nullify(Ham%fockACE_k)
710 [ + + ]: 79605 : if (associated(Ham%fockbz)) nullify(Ham%fockbz)
711 : : #if defined HAVE_GPU_CUDA
712 : : if(Ham%gpu_option==ABI_GPU_LEGACY .or. Ham%gpu_option==ABI_GPU_KOKKOS) then
713 : : call gpu_finalize_ham_data()
714 : : end if
715 : : #endif
716 : :
717 : : if (Ham%gpu_option==ABI_GPU_OPENMP) then
718 : : fofr_work_ptr => Ham%fofr_work
719 : : #ifdef HAVE_OPENMP_OFFLOAD
720 : : !$OMP TARGET EXIT DATA MAP(delete:fofr_work_ptr) IF (Ham%gpu_option==ABI_GPU_OPENMP)
721 : : #endif
722 : : end if
723 [ - + ]: 79605 : ABI_SFREE(Ham%fofr_work)
724 : :
725 : : DBG_EXIT("COLL")
726 : :
727 [ + + ]: 79605 : end subroutine gsham_free
728 : : !!***
729 : :
730 : : !----------------------------------------------------------------------
731 : :
732 : : !!****f* m_hamiltonian/gsham_init
733 : : !! NAME
734 : : !! gsham_init
735 : : !!
736 : : !! FUNCTION
737 : : !! Creation method for the gs_hamiltonian_type structure.
738 : : !! It allocates memory and initializes all quantities that do not depend on the k-point or spin.
739 : : !!
740 : : !! INPUTS
741 : : !! [comm_atom]=optional, MPI communicator over atoms
742 : : !! [fock <type(fock_type)>]= common quantities to calculate Fock exact exchange
743 : : !! natom=Number of atoms in the unit cell.
744 : : !! nfft=Number of FFT grid points (for this processors).
745 : : !! nspinor=Number of spinorial components
746 : : !! nsppol=1 for unpolarized, 2 for spin-polarized
747 : : !! nspden=Number of spin density components.
748 : : !! mgfft=Maximum size for 1D FFTs i.e., MAXVAL(ngfft(1:3))
749 : : !! [mpi_atmtab(:)]=optional, indexes of the atoms treated by current proc
750 : : !! [mpi_spintab(2)]=optional, flags defining the spin(s) treated be current process:
751 : : !! mpi_spintab(1)=1 if non-polarized or spin-up treated
752 : : !! mpi_spintab(2)=1 if polarized and spin-dn treated
753 : : !! psps<pseudopotential_type>=structure datatype gathering data on the pseudopotentials.
754 : : !! [electronpositron<electronpositron_type>]=Structured datatype storing data for the
755 : : !! electron-positron two-component DFT (optional).
756 : : !! ngfft(18)=integer array with FFT box dimensions and other information on FFTs, for the FINE rectangular grid.
757 : : !! nloalg(3)=governs the choice of the algorithm for non-local operator
758 : : !! [nucdipmom(3,natom)]= (optional) array of nuclear dipole moments at atomic sites
759 : : !! [ph1d(2,3*(2*mgfft+1)*natom)]=1-dimensions phase arrays for structure factor (see getph.f).
760 : : !! Optional, recalculated inside the routine if not present in input.
761 : : !! rprimd(3,3)=Direct lattice vectors in Bohr.
762 : : !! typat(natom)=Type of each atom.
763 : : !! [usecprj]=flag use only for PAW; 1 if cprj datastructure is allocated
764 : : !! [gpu_option] = GPU implementation to use, i.e. cuda, openMP, ... (0=not using GPU)
765 : : !! xred(3,natom)=Reduced coordinates of the atoms.
766 : : !! pawtab(ntypat*psps%usepaw)<pawtab_type>=PAW TABulated data initialized at start.
767 : : !! [paw_ij(:) <type(paw_ij_type)>]=optional, paw arrays given on (i,j) channels
768 : : !!
769 : : !! SIDE EFFECTS
770 : : !! Ham<gs_hamiltonian_type>=Structured datatype almost completely initialized:
771 : : !! * Basic variables and dimensions are transferred to the structure.
772 : : !! * All pointers are allocated with correct dimensions.
773 : : !! * Quantities that do not depend on the k-point or spin are initialized.
774 : : !!
775 : : !! SOURCE
776 : :
777 : 79597 : subroutine gsham_init(ham,Psps,pawtab,nspinor,nsppol,nspden,natom,typat,&
778 : 79597 : xred,nfft,mgfft,ngfft,rprimd,nloalg,&
779 [ + + + + : 79597 : ph1d,usecprj,comm_atom,mpi_atmtab,mpi_spintab,paw_ij,& ! optional
- + + + +
+ - + ]
780 : 75576 : electronpositron,fock,nucdipmom,gpu_option,spinaxis,& ! optional
781 : : use_gbt,zora,nfft_blocks) ! optional
782 : :
783 : : !Arguments ------------------------------------
784 : : !scalars
785 : : class(gs_hamiltonian_type),intent(inout),target :: ham
786 : : integer,intent(in) :: nfft,natom,nspinor,nsppol,nspden,mgfft
787 : : integer,optional,intent(in) :: comm_atom,usecprj,gpu_option,use_gbt,zora,nfft_blocks
788 : : type(electronpositron_type),optional,pointer :: electronpositron
789 : : type(fock_type),optional,pointer :: fock
790 : : type(pseudopotential_type),intent(in) :: psps
791 : : !arrays
792 : : integer,intent(in) :: ngfft(18),nloalg(3),typat(natom)
793 : : integer,optional,intent(in) :: mpi_atmtab(:),mpi_spintab(2)
794 : : real(dp),intent(in) :: rprimd(3,3)
795 : : real(dp),intent(in),target :: xred(3,natom)
796 : : real(dp),optional,intent(in) :: nucdipmom(3,natom),ph1d(2,3*(2*mgfft+1)*natom),spinaxis(3)
797 : : type(pawtab_type),intent(in) :: pawtab(psps%ntypat*psps%usepaw)
798 : : type(paw_ij_type),optional,intent(in) :: paw_ij(:)
799 : :
800 : : !Local variables-------------------------------
801 : : !scalars
802 : : integer :: my_comm_atom,my_nsppol,my_zora,itypat,iat,ilmn,indx,isp
803 : : integer :: cplex_dij,jsp,l_gpu_option
804 : : real(dp) :: ucvol
805 : : !arrays
806 : : integer :: my_spintab(2)
807 : : real(dp) :: gmet(3,3),gprimd(3,3),rmet(3,3)
808 : 79597 : real(dp),allocatable,target :: ekb_tmp(:,:,:,:)
809 : : ! *************************************************************************
810 : :
811 : : DBG_ENTER("COLL")
812 : :
813 : : !@gs_hamiltonian_type
814 : :
815 : : !Manage optional parameters
816 [ + + ]: 79597 : my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
817 [ + + + + ]: 166580 : my_spintab=0;my_spintab(1:nsppol)=1;if (present(mpi_spintab)) my_spintab(1:2)=mpi_spintab(1:2)
818 [ + + + + ]: 238791 : my_nsppol=count(my_spintab==1)
819 [ + + ]: 79597 : l_gpu_option=ABI_GPU_DISABLED; if(present(gpu_option)) l_gpu_option=gpu_option
820 [ - + ]: 79597 : my_zora=0; if (present(zora)) my_zora=zora
821 : :
822 : 79597 : ham%nfft_blocks=1;
823 [ - + ]: 79597 : if (present(nfft_blocks) .and. l_gpu_option==ABI_GPU_OPENMP) then
824 : 0 : ham%nfft_blocks=nfft_blocks
825 : : end if
826 [ + + ]: 79597 : ham%use_gbt = 0; if (present(use_gbt)) ham%use_gbt = use_gbt
827 [ + + + + : 460870 : ham%spinaxis = zero; if (present(spinaxis)) ham%spinaxis = spinaxis
+ + ]
828 : :
829 : 79597 : call metric(gmet,gprimd,-1,rmet,rprimd,ucvol)
830 : :
831 [ + + - + : 318388 : ABI_CHECK_IEQ(mgfft, MAXVAL(ngfft(1:3)), "Wrong mgfft")
- - ]
832 : :
833 : : !Allocate the arrays of the Hamiltonian whose dimensions do not depend on k
834 [ + - ]: 79597 : if(l_gpu_option == ABI_GPU_KOKKOS) then
835 : : #if defined HAVE_GPU && defined HAVE_YAKL
836 : : ABI_MALLOC_MANAGED(ham%atindx,(/natom/))
837 : : ABI_MALLOC_MANAGED(ham%atindx1,(/natom/))
838 : : ABI_MALLOC_MANAGED(ham%typat,(/natom/))
839 : : ABI_MALLOC_MANAGED(ham%indlmn,(/6,psps%lmnmax,psps%ntypat/))
840 : : ABI_MALLOC_MANAGED(ham%nattyp,(/psps%ntypat/))
841 : : ABI_MALLOC_MANAGED(ham%ph1d,(/2,3*(2*mgfft+1)*natom/))
842 : : #endif
843 : : else
844 [ - + + - : 238791 : ABI_MALLOC(ham%atindx,(natom))
+ - + - -
+ - + -
+ ]
845 [ - + + - : 159194 : ABI_MALLOC(ham%atindx1,(natom))
- + - + -
+ ]
846 [ - + + - : 159194 : ABI_MALLOC(ham%typat,(natom))
- + - + -
+ ]
847 [ - + + - : 318388 : ABI_MALLOC(ham%indlmn,(6,psps%lmnmax,psps%ntypat))
- + + - +
- + - - +
- + - + ]
848 [ - + + - : 238791 : ABI_MALLOC(ham%nattyp,(psps%ntypat))
+ - + - -
+ - + -
+ ]
849 [ + - - + : 238791 : ABI_MALLOC(ham%ph1d,(2,3*(2*mgfft+1)*natom))
+ - + - +
- - + - +
- + ]
850 : : end if
851 : :
852 [ - + + - : 238791 : ABI_MALLOC(ham%pspso,(psps%ntypat))
+ - + - -
+ - + -
+ ]
853 [ - + + - : 238791 : ABI_MALLOC(ham%nucdipmom,(3,natom))
+ - + - -
+ - + -
+ ]
854 : :
855 : : !Initialize most of the Hamiltonian
856 : 79597 : indx=1
857 [ + + ]: 194862 : do itypat=1,psps%ntypat
858 : 115265 : ham%nattyp(itypat)=0
859 [ + + ]: 501157 : do iat=1,natom
860 [ + + ]: 421560 : if (typat(iat)==itypat) then
861 : 204838 : ham%atindx (iat )=indx
862 : 204838 : ham%atindx1(indx)=iat
863 : 204838 : indx=indx+1
864 : 204838 : ham%nattyp(itypat)=ham%nattyp(itypat)+1
865 : : end if
866 : : end do
867 : : end do
868 : :
869 [ + + + + ]: 1034761 : ham%gmet(:,:) =gmet(:,:)
870 [ + + + + ]: 1034761 : ham%gprimd(:,:)=gprimd(:,:)
871 [ + + + + : 4698872 : ham%indlmn(:,:,:)=psps%indlmn(:,:,:)
+ + ]
872 : 79597 : ham%ispin_gbt =-1
873 : 79597 : ham%lmnmax =psps%lmnmax
874 : 79597 : ham%mgfft =mgfft
875 : 79597 : ham%mpsang =psps%mpsang
876 : 79597 : ham%mpssoang =psps%mpssoang
877 : 79597 : ham%natom =natom
878 : 79597 : ham%nfft =nfft
879 [ + + ]: 1512343 : ham%ngfft(:) =ngfft(:)
880 [ + + ]: 318388 : ham%nloalg(:) =nloalg(:)
881 [ + + + - ]: 194862 : ham%matblk=min(NLO_MINCAT,maxval(ham%nattyp)); if (nloalg(2)>0) ham%matblk=natom
882 : 79597 : ham%nsppol =nsppol
883 : 79597 : ham%nspinor =nspinor
884 : 79597 : ham%ntypat =psps%ntypat
885 [ + - - + : 364032 : ham%typat =typat(1:natom)
- - + + ]
886 [ + + ]: 79597 : ham%nvloc=1; if(nspden==4)ham%nvloc=4
887 : 79597 : ham%n4 =ngfft(4)
888 : 79597 : ham%n5 =ngfft(5)
889 : 79597 : ham%n6 =ngfft(6)
890 : 79597 : ham%usepaw =psps%usepaw
891 : 79597 : ham%ucvol =ucvol
892 : 79597 : ham%useylm =psps%useylm
893 : 79597 : ham%zora =my_zora
894 [ + + ]: 79597 : ham%gpu_option=ABI_GPU_DISABLED ; if(PRESENT(gpu_option)) ham%gpu_option=gpu_option
895 : :
896 [ + + ]: 194862 : ham%pspso(:) =psps%pspso(1:psps%ntypat)
897 [ + + ]: 79597 : if (psps%usepaw==1) then
898 [ + + ]: 57817 : do itypat=1,psps%ntypat
899 : 57817 : ham%pspso(itypat)=1+pawtab(itypat)%usespnorb
900 : : end do
901 : : end if
902 : :
903 [ + + ]: 79597 : if (present(nucdipmom)) then
904 [ + + + + ]: 860436 : ham%nucdipmom(:,:) = nucdipmom(:,:)
905 : : else
906 [ + + + + ]: 38513 : ham%nucdipmom(:,:) = zero
907 : : end if
908 : :
909 : 79597 : ham%xred => xred
910 : :
911 [ + + ]: 79597 : if (present(fock)) then
912 [ + + ]: 55047 : if (associated(fock)) then
913 : 2735 : ham%fockcommon => fock%fock_common
914 [ + + ]: 2735 : if (fock%fock_common%use_ACE==0) ham%fockbz=>fock%fock_BZ
915 : : end if
916 : : end if
917 : :
918 [ + + ]: 79597 : if (present(ph1d)) then
919 [ + + + + ]: 91195810 : ham%ph1d(:,:) = ph1d(:,:)
920 : : else ! Recalculate structure factor phases
921 : 264 : call getph(ham%atindx,natom,ngfft(1),ngfft(2),ngfft(3),ham%ph1d,xred)
922 : : end if
923 : :
924 [ + + ]: 79597 : if (ham%usepaw==1) then
925 [ + + ]: 22882 : ham%usecprj=0;if (present(usecprj)) ham%usecprj=usecprj
926 [ - + + - : 68646 : ABI_MALLOC(ham%dimcprj,(natom))
+ - + - -
+ - + -
+ ]
927 : : !Be careful cprj are ordered by atom type (used in non-local operator)
928 : 22882 : call pawcprj_getdim(ham%dimcprj,natom,ham%nattyp,ham%ntypat,ham%typat,pawtab,'O')
929 : : else
930 : 56715 : ham%usecprj=0
931 [ - + - + ]: 56715 : ABI_MALLOC(ham%dimcprj,(0))
932 : : end if
933 : :
934 : : ! ===========================
935 : : ! ==== Non-local factors ====
936 : : ! ===========================
937 : :
938 [ + + ]: 79597 : if (ham%usepaw==0) then
939 : : ! Norm-conserving: use constant Kleimann-Bylander energies.
940 : : ! nspinor ** 2 is a fake dimension here in the sense that
941 : : ! the KB energies for the scalar part and the SOC part are packed in the firs dimension (dimekb).
942 : : ! In nonlop_pl, ekb are accessed using ekb(iln,itypat,ispinor) where iln runs over all projects (scalar + SOC)
943 : : ! The ispinor index is irrelevant as ekb(:,:,1) = ekb(:,:,2). See nonlop_pl
944 : 56715 : ham%dimekb1=psps%dimekb
945 : 56715 : ham%dimekb2=psps%ntypat
946 : 56715 : ham%dimekbq=1
947 [ - + + - : 283575 : ABI_MALLOC(ham%ekb_spin,(psps%dimekb,psps%ntypat,nspinor**2,1,1))
- + + - -
+ + - + -
+ - - + -
+ - + ]
948 : 56715 : ham%ekb => ham%ekb_spin(:,:,:,:,1)
949 [ - + - + ]: 56715 : ABI_MALLOC(ham%sij,(0,0))
950 [ + + + + ]: 382737 : ham%ekb(:,:,1,1)=psps%ekb(:,:)
951 [ + + ]: 56715 : if (nspinor==2) then
952 [ + + + + ]: 37216 : ham%ekb(:,:,2,1)=psps%ekb(:,:)
953 [ + + + + : 77173 : ham%ekb(:,:,3:4,1)=zero
+ + ]
954 : : end if
955 [ + + ]: 56715 : if (PRESENT(electronpositron)) then
956 [ + + + + : 36903 : if (electronpositron_calctype(electronpositron)==1) ham%ekb(:,:,:,:)=-ham%ekb(:,:,:,:)
+ + + + +
+ ]
957 : : end if
958 : :
959 : : ! Update enl on GPU (will do it later for PAW)
960 : : #if defined HAVE_GPU_CUDA
961 : : if (ham%gpu_option==ABI_GPU_LEGACY .or. ham%gpu_option==ABI_GPU_KOKKOS) then
962 : : call gpu_update_ham_data(&
963 : : & ham%ekb(:,:,:,1), INT(size(ham%ekb), c_int64_t), &
964 : : & ham%sij, INT(size(ham%sij), c_int64_t), &
965 : : & ham%gprimd, INT(size(ham%gprimd),c_int64_t))
966 : : end if
967 : : #endif
968 : :
969 : : else
970 : : ! PAW: store overlap coefficients (spin non dependent) and Dij coefficients (spin dependent)
971 : 22882 : cplex_dij=1
972 [ + - ]: 22882 : if (present(paw_ij)) then
973 [ + + ]: 22882 : if (size(paw_ij)>0) cplex_dij=paw_ij(1)%cplex_dij
974 : : end if
975 [ + + + + : 272028 : if ((nspinor==2).or.any(abs(ham%nucdipmom)>tol8)) cplex_dij=2
+ + + + +
+ ]
976 : 22882 : ham%dimekb1=psps%dimekb*cplex_dij
977 : 22882 : ham%dimekb2=natom
978 : 22882 : ham%dimekbq=1
979 [ + - ]: 22882 : if (present(paw_ij)) then
980 [ + + ]: 22882 : if (size(paw_ij)>0) ham%dimekbq=paw_ij(1)%qphase
981 : : end if
982 [ - + + - : 91528 : ABI_MALLOC(ham%sij,(ham%dimekb1,psps%ntypat))
- + + - +
- + - - +
- + - + ]
983 [ + + ]: 57817 : do itypat=1,psps%ntypat
984 [ + + ]: 34935 : if (cplex_dij==1) then
985 [ + + ]: 1700145 : ham%sij(1:pawtab(itypat)%lmn2_size,itypat)=pawtab(itypat)%sij(:)
986 : : else
987 [ + + ]: 411505 : do ilmn=1,pawtab(itypat)%lmn2_size
988 : 406390 : ham%sij(2*ilmn-1,itypat)=pawtab(itypat)%sij(ilmn)
989 : 411505 : ham%sij(2*ilmn ,itypat)=zero
990 : : end do
991 : : end if
992 [ + + ]: 57817 : if (cplex_dij*pawtab(itypat)%lmn2_size<ham%dimekb1) then
993 [ + + ]: 228967 : ham%sij(cplex_dij*pawtab(itypat)%lmn2_size+1:ham%dimekb1,itypat)=zero
994 : : end if
995 : : end do
996 : : ! We preload here PAW non-local factors in order to avoid a communication over atoms
997 : : ! inside the loop over spins.
998 [ - + + - : 160174 : ABI_MALLOC(ham%ekb_spin,(ham%dimekb1,ham%dimekb2,nspinor**2,ham%dimekbq,my_nsppol))
- + + - -
+ + - - +
+ - + + +
- + - + +
- + - + -
+ ]
999 [ + + + + : 8221893 : ham%ekb_spin=zero
+ + + + +
+ ]
1000 [ - + ]: 22882 : if (present(paw_ij)) then
1001 [ + + ]: 22882 : if (my_nsppol<ham%nsppol) then
1002 [ - + + - : 3168 : ABI_MALLOC(ekb_tmp,(ham%dimekb1,ham%dimekb2,nspinor**2,ham%dimekbq))
- + + - -
+ + - - +
+ - + - +
- - + -
+ ]
1003 : : end if
1004 : 22882 : jsp=0
1005 [ + + ]: 48550 : do isp=1,ham%nsppol
1006 [ + + ]: 25668 : if (my_spintab(isp)==1) then
1007 : 25112 : jsp=jsp+1 ; ham%ekb => ham%ekb_spin(:,:,:,:,jsp)
1008 : : else
1009 : 556 : ham%ekb => ekb_tmp
1010 : : end if
1011 [ + + ]: 48550 : if (present(mpi_atmtab)) then
1012 : 7556 : call pawdij2ekb(ham%ekb,paw_ij,isp,my_comm_atom,mpi_atmtab=mpi_atmtab)
1013 : : else
1014 : 18112 : call pawdij2ekb(ham%ekb,paw_ij,isp,my_comm_atom)
1015 : : end if
1016 : : end do
1017 [ + + ]: 22882 : if (my_nsppol<ham%nsppol) then
1018 [ - + ]: 528 : ABI_FREE(ekb_tmp)
1019 : : end if
1020 : : end if
1021 : 22882 : nullify(ham%ekb)
1022 : : end if
1023 : :
1024 : : DBG_EXIT("COLL")
1025 : :
1026 [ + + + + : 234770 : end subroutine gsham_init
- + ]
1027 : : !!***
1028 : :
1029 : : !!****f* m_hamiltonian/gsham_load_k
1030 : : !! NAME
1031 : : !! gsham_load_k
1032 : : !!
1033 : : !! FUNCTION
1034 : : !! Setup of the k-dependent part of the Hamiltonian H_k_k^prime
1035 : : !!
1036 : : !! INPUTS
1037 : : !! [compute_gbound]=flag. if true, G sphere boundary is computed here
1038 : : !! [compute_ph3d]=flag. if true, 3D structure factors are computed here (only if nloalg(1)>0)
1039 : : !! [gbound_k]=G sphere boundary (not compatible with compute_gbound=TRUE)
1040 : : !! [ffnl_k]=nonlocal form factors on basis sphere
1041 : : !! [istwf_k]=parameter that describes the storage of wfs
1042 : : !! [kinpw_k]=(modified) kinetic energy for each plane wave
1043 : : !! [kg_k]=planewave reduced coordinates in basis sphere (g vectors)
1044 : : !! [kpg_k]=(k+g) vectors in reciprocal space
1045 : : !! [kpt_k]=k point coordinates
1046 : : !! [npw_k]=number of plane waves (processed by current proc when band-FFT parallelism is on)
1047 : : !! [npw_fft_k]=number of plane waves used to apply Hamiltonian (in the "FFT" configuration)
1048 : : !! [ph3d_k]=3-dim structure factors, for each atom and plane wave
1049 : : !!
1050 : : !! SIDE EFFECTS
1051 : : !! ham<gs_hamiltonian_type>=structured datatype completed with k-dependent quantities.
1052 : : !! Quantities at k^prime are set equal to quantities at k.
1053 : : !! k-dependent scalars and pointers associated
1054 : : !! phkxred=exp(.k.xred) for each atom
1055 : : !! [ham%gbound_k]=G sphere boundary, for each plane wave
1056 : : !! [ham%ph3d_k]=3-dim structure factors, for each atom and plane wave
1057 : : !!
1058 : : !! SOURCE
1059 : :
1060 [ + + + - : 9551097 : subroutine gsham_load_k(ham,ffnl_k,fockACE_k,gbound_k,istwf_k,kinpw_k,&
- + + + +
- - + + +
+ - - + ]
1061 [ + + + - : 6367398 : kg_k,kpg_k,kpt_k,npw_k,npw_fft_k,ph3d_k,&
- + + + +
- - + + +
+ - - + ]
1062 : : compute_gbound,compute_ph3d)
1063 : :
1064 : : !Arguments ------------------------------------
1065 : : !scalars
1066 : : class(gs_hamiltonian_type),intent(inout),target :: ham
1067 : : integer,intent(in),optional :: npw_k,npw_fft_k,istwf_k
1068 : : logical,intent(in),optional :: compute_gbound,compute_ph3d
1069 : : !arrays
1070 : : integer,intent(in),optional,target :: gbound_k(:,:),kg_k(:,:)
1071 : : real(dp),intent(in),optional :: kpt_k(3)
1072 : : real(dp),intent(in),optional,target :: ffnl_k(:,:,:,:),kinpw_k(:),kpg_k(:,:),ph3d_k(:,:,:)
1073 : : type(fock_ACE_type),intent(in),optional,target :: fockACE_k
1074 : :
1075 : : !Local variables-------------------------------
1076 : : !scalars
1077 : : integer :: iat,iatom
1078 : : logical :: compute_gbound_
1079 : : real(dp) :: arg
1080 : : !character(len=500) :: msg
1081 : : ! *************************************************************************
1082 : :
1083 : : DBG_ENTER("COLL")
1084 : :
1085 : : !@gs_hamiltonian_type
1086 : :
1087 : : !k-dependent scalars
1088 [ + + ]: 3183699 : if (present(kpt_k)) then
1089 [ + + ]: 7486960 : ham%kpt_k(:) = kpt_k(:)
1090 [ + + ]: 7486960 : ham%kpt_kp(:) = kpt_k(:)
1091 : : end if
1092 [ + + ]: 3183699 : if (present(istwf_k)) then
1093 : 1871740 : ham%istwf_k = istwf_k
1094 : 1871740 : ham%istwf_kp = istwf_k
1095 : : end if
1096 [ + + ]: 3183699 : if (present(npw_k)) then
1097 : 1871740 : ham%npw_k = npw_k
1098 : 1871740 : ham%npw_kp = npw_k
1099 : : end if
1100 [ + + ]: 3183699 : if (present(npw_fft_k)) then
1101 : 46250 : ham%npw_fft_k = npw_fft_k
1102 : 46250 : ham%npw_fft_kp = npw_fft_k
1103 [ + + ]: 3137449 : else if (present(npw_k)) then
1104 : 1871740 : ham%npw_fft_k = npw_k
1105 : 1871740 : ham%npw_fft_kp = npw_k
1106 : : end if
1107 : :
1108 : : !Pointers to k-dependent quantities
1109 [ + + ]: 3183699 : if (present(kinpw_k)) then
1110 : 507674 : ham%kinpw_k => kinpw_k
1111 : 507674 : ham%kinpw_kp => kinpw_k
1112 : : end if
1113 [ + + ]: 3183699 : if (present(kg_k)) then
1114 : 1917998 : ham%kg_k => kg_k
1115 : 1917998 : ham%kg_kp => kg_k
1116 : : end if
1117 [ + + ]: 3183699 : if (present(kpg_k)) then
1118 : 1907620 : ham%kpg_k => kpg_k
1119 : 1907620 : ham%kpg_kp => kpg_k
1120 : : end if
1121 [ + + ]: 3183699 : if (present(ffnl_k)) then
1122 : 1976866 : ham%ffnl_k => ffnl_k
1123 : 1976866 : ham%ffnl_kp => ffnl_k
1124 : : end if
1125 [ + + ]: 3183699 : if (present(ph3d_k)) then
1126 : 1798057 : ham%ph3d_k => ph3d_k
1127 : 1798057 : ham%ph3d_kp => ph3d_k
1128 : : end if
1129 [ + + ]: 3183699 : if (present(fockACE_k)) then
1130 : 6063 : ham%fockACE_k => fockACE_k
1131 : : end if
1132 : : !Compute exp(i.k.R) for each atom
1133 [ + + ]: 3183699 : if (present(kpt_k)) then
1134 [ + + + + ]: 1871740 : if (associated(Ham%phkpxred).and.(.not.associated(Ham%phkpxred,Ham%phkxred))) then
1135 : 1357224 : ABI_FREE(Ham%phkpxred)
1136 : : end if
1137 [ + + ]: 1871740 : ABI_SFREE(ham%phkxred)
1138 [ - + + - : 5615220 : ABI_MALLOC(ham%phkxred,(2,ham%natom))
+ - + - -
+ - + -
+ ]
1139 [ + + ]: 5486657 : do iat=1,ham%natom
1140 : 3614917 : iatom=ham%atindx(iat)
1141 [ + + ]: 14459668 : arg=two_pi*DOT_PRODUCT(kpt_k,ham%xred(:,iat))
1142 : 3614917 : ham%phkxred(1,iatom)=DCOS(arg)
1143 : 5486657 : ham%phkxred(2,iatom)=DSIN(arg)
1144 : : end do
1145 : 1871740 : ham%phkpxred => ham%phkxred
1146 : : end if
1147 : :
1148 : : ! Compute or copy G sphere boundary at k+g
1149 [ + + ]: 3183699 : compute_gbound_=.false.;if (present(compute_gbound)) compute_gbound_=compute_gbound
1150 [ + + ]: 3183699 : if (present(gbound_k)) compute_gbound_=.true.
1151 [ + + ]: 3137449 : if (compute_gbound_) then
1152 [ + + ]: 1696167 : if (associated(Ham%gbound_kp,Ham%gbound_k)) then
1153 : 378404 : nullify(Ham%gbound_kp)
1154 [ + + ]: 1317763 : else if (associated(Ham%gbound_kp)) then
1155 : 1254876 : ABI_FREE(Ham%gbound_kp)
1156 : : end if
1157 [ + + ]: 1696167 : ABI_SFREE(ham%gbound_k)
1158 : : end if
1159 [ + + ]: 3183699 : if (.not.allocated(ham%gbound_k)) then
1160 [ + - - + : 6849920 : ABI_MALLOC(ham%gbound_k,(2*ham%mgfft+8,2))
+ - + - +
- + - - +
- + - + ]
1161 [ + + + + ]: 139520412 : ham%gbound_k(:,:)=0
1162 : 1712480 : ham%gbound_kp => ham%gbound_k
1163 : : end if
1164 [ + + ]: 3183699 : if (compute_gbound_) then
1165 [ + + ]: 1696167 : if (present(gbound_k)) then
1166 [ + + + + ]: 4388306 : ham%gbound_k(:,:)=gbound_k(:,:)
1167 : : else
1168 [ - + ]: 1649917 : if (.not.associated(ham%kg_k)) then
1169 : 0 : ABI_BUG('Something is missing for gbound_k computation!')
1170 : : end if
1171 : : !write(std_out,*)"About to call sphereboundary"
1172 : : !write(std_out,*)"size(kg_k), npw_k, mgfft",size(ham%kg_k, dim=2), ham%npw_k, ham%mgfft
1173 [ + - - - : 1649917 : call sphereboundary(ham%gbound_k,ham%istwf_k,ham%kg_k,ham%mgfft,ham%npw_k)
- - - - -
+ ]
1174 : : end if
1175 : 1696167 : ham%gbound_kp => ham%gbound_k
1176 : : end if
1177 : :
1178 : : ! Compute 3D structure factors for each atom at k+g
1179 [ + + ]: 3183699 : if (present(compute_ph3d).and.present(ph3d_k)) then
1180 [ + + + - ]: 1751120 : if (compute_ph3d.and.ham%nloalg(2)>0) then
1181 [ + - + - : 1718160 : if ((.not.allocated(ham%phkxred)).or.(.not.associated(ham%kg_k)).or.&
- + ]
1182 : : (.not.associated(ham%ph3d_k))) then
1183 : 0 : ABI_BUG('Something is missing for ph3d_k computation!')
1184 : : end if
1185 : : call ph1d3d(1,ham%natom,ham%kg_k,ham%matblk,ham%natom,ham%npw_k,ham%ngfft(1),&
1186 [ + - - - : 1718160 : ham%ngfft(2),ham%ngfft(3),ham%phkxred,ham%ph1d,ham%ph3d_k)
- - - - +
- - - - +
- + - - -
- - - ]
1187 : : end if
1188 : : end if
1189 : :
1190 : : DBG_EXIT("COLL")
1191 : :
1192 : 3183699 : end subroutine gsham_load_k
1193 : : !!***
1194 : :
1195 : : !----------------------------------------------------------------------
1196 : :
1197 : : !!****f* m_hamiltonian/gsham_load_kprime
1198 : : !! NAME
1199 : : !! gsham_load_kprime
1200 : : !!
1201 : : !! FUNCTION
1202 : : !! Setup of the k^prime-dependent part of the Hamiltonian H_k_k^prime
1203 : : !!
1204 : : !! INPUTS
1205 : : !! [compute_gbound]=flag. if true, G sphere boundary is computed here
1206 : : !! [compute_ph3d]=flag. if true, 3D structure factors are computed here (only if nloalg(2)>0)
1207 : : !! [gbound_kp]=G sphere boundary (not compatible with compute_gbound=TRUE)
1208 : : !! [ffnl_kp]=nonlocal form factors on basis sphere
1209 : : !! [istwf_kp]=parameter that describes the storage of wfs
1210 : : !! [kinpw_kp]=(modified) kinetic energy for each plane wave
1211 : : !! [kg_kp]=planewave reduced coordinates in basis sphere (g vectors)
1212 : : !! [kpg_kp]=(k+g) vectors in reciprocal space
1213 : : !! [kpt_kp]=k point coordinates
1214 : : !! [npw_kp]=number of plane waves (processed by current proc when band-FFT parallelism is on)
1215 : : !! [npw_fft_kp]=number of plane waves used to apply Hamiltonian (in the "FFT" configuration)
1216 : : !! [ph3d_kp]=3-dim structure factors, for each atom and plane wave
1217 : : !!
1218 : : !! SIDE EFFECTS
1219 : : !! ham<gs_hamiltonian_type>=structured datatype completed with k^prime-dependent quantities.
1220 : : !! k^prime-dependent scalars and pointers associated
1221 : : !! phkpxred=exp(.k^prime.xred) for each atom
1222 : : !! [ham%gbound_kp]=G sphere boundary, for each plane wave
1223 : : !! [ham%ph3d_kp]=3-dim structure factors at k^prime at k, for each atom and plane wave
1224 : : !!
1225 : : !! SOURCE
1226 : :
1227 [ + + + - : 6121458 : subroutine gsham_load_kprime(ham,ffnl_kp,gbound_kp,istwf_kp,kinpw_kp,&
- + + + +
- - + - +
- - - - ]
1228 [ + + + - : 2040486 : kg_kp,kpg_kp,kpt_kp,npw_kp,npw_fft_kp,&
- + + + +
- - + ]
1229 [ + + + - : 2040486 : ph3d_kp,compute_gbound,compute_ph3d)
- + ]
1230 : :
1231 : : !Arguments ------------------------------------
1232 : : !scalars
1233 : : class(gs_hamiltonian_type),intent(inout),target :: ham
1234 : : integer,intent(in),optional :: npw_kp,npw_fft_kp,istwf_kp
1235 : : logical,intent(in),optional :: compute_gbound,compute_ph3d
1236 : : !arrays
1237 : : integer,intent(in),optional,target :: gbound_kp(:,:),kg_kp(:,:)
1238 : : real(dp),intent(in),optional :: kpt_kp(3)
1239 : : real(dp),intent(in),optional,target :: ffnl_kp(:,:,:,:),kinpw_kp(:),kpg_kp(:,:),ph3d_kp(:,:,:)
1240 : :
1241 : : !Local variables-------------------------------
1242 : : !scalars
1243 : : integer :: iat,iatom
1244 : : logical :: compute_gbound_
1245 : : real(dp) :: arg
1246 : : !character(len=500) :: msg
1247 : : ! *************************************************************************
1248 : :
1249 : : DBG_ENTER("COLL")
1250 : :
1251 : : !@gs_hamiltonian_type
1252 : :
1253 : : !k-dependent scalars
1254 [ + + + + ]: 6869217 : if (present(kpt_kp)) ham%kpt_kp(:)= kpt_kp(:)
1255 [ + + ]: 2040486 : if (present(istwf_kp)) ham%istwf_kp = istwf_kp
1256 [ + + ]: 2040486 : if (present(npw_kp)) ham%npw_kp = npw_kp
1257 [ - + ]: 2040486 : if (present(npw_fft_kp)) then
1258 : 0 : ham%npw_fft_kp = npw_fft_kp
1259 [ + + ]: 2040486 : else if (present(npw_kp)) then
1260 : 1608177 : ham%npw_fft_kp = npw_kp
1261 : : end if
1262 : :
1263 : : !Pointers to k-dependent quantities
1264 [ + + ]: 2040486 : if (present(kinpw_kp)) ham%kinpw_kp => kinpw_kp
1265 [ + + ]: 2040486 : if (present(kg_kp)) ham%kg_kp => kg_kp
1266 [ + + ]: 2040486 : if (present(kpg_kp)) ham%kpg_kp => kpg_kp
1267 [ + + ]: 2040486 : if (present(ffnl_kp)) ham%ffnl_kp => ffnl_kp
1268 [ + + ]: 2040486 : if (present(ph3d_kp)) ham%ph3d_kp => ph3d_kp
1269 : :
1270 : : !Compute exp(i.k^prime.R) for each atom
1271 [ + + ]: 2040486 : if (present(kpt_kp)) then
1272 [ + + ]: 1609577 : if (associated(ham%phkpxred,ham%phkxred)) then
1273 : 1380394 : nullify(ham%phkpxred)
1274 [ + - ]: 229183 : else if (associated(ham%phkpxred)) then
1275 : 229183 : ABI_FREE(ham%phkpxred)
1276 : : end if
1277 [ - + + - : 4828731 : ABI_MALLOC(ham%phkpxred,(2,ham%natom))
+ - + - -
+ - + ]
1278 [ + + ]: 4745605 : do iat=1,ham%natom
1279 : 3136028 : iatom=ham%atindx(iat)
1280 [ + + ]: 12544112 : arg=two_pi*DOT_PRODUCT(kpt_kp,ham%xred(:,iat))
1281 : 3136028 : ham%phkpxred(1,iatom)=DCOS(arg)
1282 : 4745605 : ham%phkpxred(2,iatom)=DSIN(arg)
1283 : : end do
1284 : : end if
1285 : :
1286 : : !Compute or copy G sphere boundary at k^prime+g
1287 : 2040486 : compute_gbound_=.false.
1288 [ + + ]: 2040486 : if (present(kpt_kp).and.present(compute_gbound)) compute_gbound_=compute_gbound
1289 [ + - ]: 2040486 : if (present(gbound_kp)) compute_gbound_=.true.
1290 [ + + ]: 2040486 : if (compute_gbound_) then
1291 [ + + ]: 1278613 : if (associated(ham%gbound_kp,ham%gbound_k)) then
1292 : 1274845 : nullify(ham%gbound_kp)
1293 [ + - ]: 3768 : else if (associated(ham%gbound_kp)) then
1294 : 3768 : ABI_FREE(ham%gbound_kp)
1295 : : end if
1296 [ - + ]: 1278613 : if (present(gbound_kp)) then
1297 [ # # # # : 0 : ABI_MALLOC(ham%gbound_kp,(size(gbound_kp,1),size(gbound_kp,2)))
# # # # #
# # # # #
# # ]
1298 [ # # # # ]: 0 : ham%gbound_kp(:,:)=gbound_kp(:,:)
1299 : : else
1300 [ - + ]: 1278613 : if (.not.associated(ham%kg_kp)) then
1301 : 0 : ABI_BUG('Something is missing for gbound_kp computation!')
1302 : : end if
1303 [ + - - + : 5114452 : ABI_MALLOC(ham%gbound_kp,(2*ham%mgfft+8,2))
+ - + - +
- + - - +
- + ]
1304 [ - + - - : 1278613 : call sphereboundary(ham%gbound_kp,ham%istwf_kp,ham%kg_kp,ham%mgfft,ham%npw_kp)
+ - - - -
- - - - +
- - - - -
+ ]
1305 : : end if
1306 : : end if
1307 : :
1308 : : !Compute 3D structure factors for each atom at k^prime+g
1309 [ + + ]: 2040486 : if (present(compute_ph3d).and.present(ph3d_kp)) then
1310 [ + + + - ]: 603822 : if (compute_ph3d.and.ham%nloalg(2)>0) then
1311 [ + - + - : 586452 : if ((.not.associated(ham%phkpxred)).or.(.not.associated(ham%kg_kp)).or.&
- + ]
1312 : : & (.not.associated(ham%ph3d_kp))) then
1313 : 0 : ABI_BUG('Something is missing for ph3d_kp computation!')
1314 : : end if
1315 : : call ph1d3d(1,ham%natom,ham%kg_kp,ham%matblk,ham%natom,ham%npw_kp,ham%ngfft(1),&
1316 [ + - - - : 586452 : & ham%ngfft(2),ham%ngfft(3),ham%phkpxred,ham%ph1d,ham%ph3d_kp)
- - - - +
- - - - -
- - + - -
- - + - +
- + - - -
- - - ]
1317 : : end if
1318 : : end if
1319 : :
1320 : : DBG_EXIT("COLL")
1321 : :
1322 : 2040486 : end subroutine gsham_load_kprime
1323 : : !!***
1324 : :
1325 : : !----------------------------------------------------------------------
1326 : :
1327 : : !!****f* m_hamiltonian/gsham_eph_setup_k
1328 : : !! NAME
1329 : : !! gsham_eph_setup_k
1330 : : !!
1331 : : !! FUNCTION
1332 : : !! Simplified interface to load either k or kprime in the case of e-ph calculations.
1333 : : !!
1334 : : !! INPUTS
1335 : : !! which_k= "k" to load k, "kq" to load kprime
1336 : : !! See load_k or load_kprime for the meaning of arguments.
1337 : : !!
1338 : : !! SOURCE
1339 : :
1340 : 157025 : subroutine gsham_eph_setup_k(gs_ham, which_k, kk, istwf_k, npw_k, kg_k, dtset, cryst, psps, & ! in
1341 : : nkpg_k, kpg_k, ffnl_k, kinpw_k, ph3d_k, comm) ! out
1342 : :
1343 : : !Arguments ------------------------------------
1344 : : !scalars
1345 : : class(gs_hamiltonian_type),intent(inout) :: gs_ham
1346 : : character(len=*),intent(in) :: which_k
1347 : : type(dataset_type),intent(in) :: dtset
1348 : : type(crystal_t),intent(in) :: cryst
1349 : : type(pseudopotential_type),intent(in) :: psps
1350 : : integer,intent(in) :: istwf_k, npw_k, comm
1351 : : integer,intent(out) :: nkpg_k
1352 : : !arrays
1353 : : real(dp),intent(in) :: kk(3)
1354 : : integer,intent(in) :: kg_k(3,npw_k)
1355 : : real(dp),allocatable,intent(out) :: kpg_k(:,:), ffnl_k(:,:,:,:), kinpw_k(:), ph3d_k(:,:,:)
1356 : :
1357 : : !Local variables-------------------------------
1358 : : !scalars
1359 : : integer,parameter :: ider0 = 0, idir0 = 0, optder0 = 0
1360 : : !arrays
1361 : : real(dp) :: ylmgr_k_dum(1,1,1)
1362 : : real(dp),allocatable :: ylm_k(:,:)
1363 : : ! *************************************************************************
1364 : :
1365 : : ! Compute k+G vectors
1366 : 157025 : nkpg_k = 3 * dtset%nloalg(3)
1367 [ - + + - : 628100 : ABI_MALLOC(kpg_k, (npw_k, nkpg_k))
+ - - - +
- - + - +
- + - + ]
1368 [ - + ]: 157025 : if (nkpg_k > 0) call mkkpg(kg_k, kpg_k, kk, nkpg_k, npw_k)
1369 : :
1370 : : ! Spherical Harmonics at k for useylm == 1.
1371 [ - + + - : 628100 : ABI_MALLOC(ylm_k, (npw_k, psps%mpsang**2 * psps%useylm))
+ + + - +
- + + - +
- + ]
1372 [ + + ]: 157025 : if (psps%useylm == 1) call initylmg_k(npw_k, psps%mpsang, optder0, cryst%rprimd, cryst%gprimd, kk, kg_k, ylm_k, ylmgr_k_dum)
1373 : :
1374 : : ! Compute nonlocal form factors ffnl_k at (k+G)
1375 [ - + + - : 785125 : ABI_MALLOC(ffnl_k, (npw_k, 1, psps%lmnmax, psps%ntypat))
- + + - -
+ + - + -
+ - - + -
+ - + ]
1376 : 157025 : call mkffnl_objs(cryst, psps, 1, ffnl_k, ider0, idir0, kg_k, kpg_k, kk, nkpg_k, npw_k, ylm_k, ylmgr_k_dum, comm=comm)
1377 [ - + ]: 157025 : ABI_FREE(ylm_k)
1378 : :
1379 : : ! Compute (1/2) (2 Pi)**2 (kG)**2:
1380 [ - + + - : 134657935 : ABI_CALLOC(kinpw_k, (npw_k))
+ - + - -
+ - + - +
+ + ]
1381 : 157025 : call mkkin(dtset%ecut, dtset%ecutsm, dtset%effmass_free, cryst%gmet, kg_k, kinpw_k, kk, npw_k, 0, 0)
1382 : :
1383 [ - + + - : 628100 : ABI_MALLOC(ph3d_k, (2, npw_k, gs_ham%matblk))
- + + - +
- + - - +
- + - + ]
1384 : :
1385 : : ! Load the k dependent parts of the Hamiltonian
1386 : 76649 : select case (which_k)
1387 : : case ("k")
1388 : : call gs_ham%load_k(kpt_k=kk, npw_k=npw_k, istwf_k=istwf_k, kg_k=kg_k, kpg_k=kpg_k, kinpw_k=kinpw_k, &
1389 : 76649 : ph3d_k=ph3d_k, ffnl_k=ffnl_k, compute_ph3d=.true., compute_gbound=.true.)
1390 : : case ("kq")
1391 : : call gs_ham%load_kprime(kpt_kp=kk, npw_kp=npw_k, istwf_kp=istwf_k, kg_kp=kg_k, kpg_kp=kpg_k, kinpw_kp=kinpw_k, &
1392 : 80376 : ph3d_kp=ph3d_k, ffnl_kp=ffnl_k, compute_ph3d=.true., compute_gbound=.true.)
1393 : : case default
1394 [ + + - ]: 157025 : ABI_ERROR(sjoin("Invalid value for which_k:", which_k))
1395 : : end select
1396 : :
1397 : 157025 : end subroutine gsham_eph_setup_k
1398 : : !!***
1399 : :
1400 : : !----------------------------------------------------------------------
1401 : :
1402 : : !!****f* m_hamiltonian/gsham_copy
1403 : : !! NAME
1404 : : !! gsham_copy
1405 : : !!
1406 : : !! INPUTS
1407 : : !! gs_hamk_in<gs_hamiltonian_type>=Structured datatype completely initialized, to be copied.
1408 : : !!
1409 : : !! FUNCTION
1410 : : !! Copy a gs_hamiltonian_type variable (gs_hamk_in) in another (gs_hamk_out).
1411 : : !! In contrast to an assignment statement (gs_hamk_out=gs_hamk_in), this
1412 : : !! subroutine allocate memory space for the pointers contained in the data
1413 : : !! structure (gs_hamk_out) and copy the content of the corresponding memory
1414 : : !! space of gs_hamk_in in it. In contrast, the assignment statement would
1415 : : !! only associate the pointers of gs_hamk_out to the same memory space than
1416 : : !! the corresponding ones in gs_hamk_in. This can cause trouble if one data
1417 : : !! structure is destroyed before a reading/writing statement for the other
1418 : : !! structure, causing access to unallocated memory space (silently, without
1419 : : !! segmentation fault being generated).
1420 : : !!
1421 : : !! OUTPUT
1422 : : !! gs_hamk_out<gs_hamiltonian_type>=Structured datatype containing separate
1423 : : !! copies of all data of gs_hamk_in upon exit.
1424 : : !!
1425 : : !! SOURCE
1426 : :
1427 [ + - ]: 8 : subroutine gsham_copy(gs_hamk_in, gs_hamk_out)
1428 : :
1429 : : !Arguments ------------------------------------
1430 : : class(gs_hamiltonian_type),intent(in),target :: gs_hamk_in
1431 : : class(gs_hamiltonian_type),intent(out),target :: gs_hamk_out
1432 : :
1433 : : !Local variables-------------------------------
1434 : : integer :: tmp2i(5)
1435 : : type(C_PTR) :: ham_ptr
1436 : : ! *************************************************************************
1437 : :
1438 : : DBG_ENTER("COLL")
1439 : :
1440 : : !@gs_hamiltonian_type
1441 : :
1442 : 8 : gs_hamk_out%dimekb1 = gs_hamk_in%dimekb1
1443 : 8 : gs_hamk_out%dimekb2 = gs_hamk_in%dimekb2
1444 : 8 : gs_hamk_out%dimekbq = gs_hamk_in%dimekbq
1445 : 8 : gs_hamk_out%istwf_k = gs_hamk_in%istwf_k
1446 : 8 : gs_hamk_out%istwf_kp = gs_hamk_in%istwf_kp
1447 : 8 : gs_hamk_out%lmnmax = gs_hamk_in%lmnmax
1448 : 8 : gs_hamk_out%matblk = gs_hamk_in%matblk
1449 : 8 : gs_hamk_out%mgfft = gs_hamk_in%mgfft
1450 : 8 : gs_hamk_out%mpsang = gs_hamk_in%mpsang
1451 : 8 : gs_hamk_out%mpssoang = gs_hamk_in%mpssoang
1452 : 8 : gs_hamk_out%natom = gs_hamk_in%natom
1453 : 8 : gs_hamk_out%nfft = gs_hamk_in%nfft
1454 : 8 : gs_hamk_out%nfft_blocks = gs_hamk_in%nfft_blocks
1455 : 8 : gs_hamk_out%npw_k = gs_hamk_in%npw_k
1456 : 8 : gs_hamk_out%npw_kp = gs_hamk_in%npw_kp
1457 : 8 : gs_hamk_out%npw_fft_k = gs_hamk_in%npw_fft_k
1458 : 8 : gs_hamk_out%npw_fft_kp = gs_hamk_in%npw_fft_kp
1459 : 8 : gs_hamk_out%nspinor = gs_hamk_in%nspinor
1460 : 8 : gs_hamk_out%nsppol = gs_hamk_in%nsppol
1461 : 8 : gs_hamk_out%ntypat = gs_hamk_in%ntypat
1462 : 8 : gs_hamk_out%nvloc = gs_hamk_in%nvloc
1463 : 8 : gs_hamk_out%n4 = gs_hamk_in%n4
1464 : 8 : gs_hamk_out%n5 = gs_hamk_in%n5
1465 : 8 : gs_hamk_out%n6 = gs_hamk_in%n6
1466 : 8 : gs_hamk_out%gpu_option = gs_hamk_in%gpu_option
1467 [ + + + + ]: 56 : gs_hamk_out%spinaxis = gs_hamk_in%spinaxis
1468 : 8 : gs_hamk_out%usecprj = gs_hamk_in%usecprj
1469 : 8 : gs_hamk_out%usepaw = gs_hamk_in%usepaw
1470 : 8 : gs_hamk_out%useylm = gs_hamk_in%useylm
1471 : 8 : gs_hamk_out%use_gbt = gs_hamk_in%use_gbt
1472 : 8 : gs_hamk_out%zora = gs_hamk_in%zora
1473 [ + + + + ]: 296 : gs_hamk_out%ngfft = gs_hamk_in%ngfft
1474 [ + + + + ]: 56 : gs_hamk_out%nloalg = gs_hamk_in%nloalg
1475 : 8 : gs_hamk_out%ucvol = gs_hamk_in%ucvol
1476 [ + + + + : 200 : gs_hamk_out%gmet = gs_hamk_in%gmet
+ + + + ]
1477 [ + + + + : 200 : gs_hamk_out%gprimd = gs_hamk_in%gprimd
+ + + + ]
1478 : 8 : gs_hamk_out%ispin_gbt = gs_hamk_in%ispin_gbt
1479 [ + + + + ]: 56 : gs_hamk_out%kpt_k = gs_hamk_in%kpt_k
1480 [ + + + + ]: 56 : gs_hamk_out%kpt_kp = gs_hamk_in%kpt_kp
1481 : :
1482 [ - + + - : 24 : ABI_MALLOC(gs_hamk_out%atindx,(gs_hamk_out%natom))
+ - + - -
+ - + -
+ ]
1483 [ + - + + : 96 : gs_hamk_out%atindx = gs_hamk_in%atindx
+ - - + -
- - - - -
+ + ]
1484 [ - + + - : 24 : ABI_MALLOC(gs_hamk_out%atindx1,(gs_hamk_out%natom))
+ - + - -
+ - + -
+ ]
1485 [ + - + + : 96 : gs_hamk_out%atindx1 = gs_hamk_in%atindx1
+ - - + -
- - - - -
+ + ]
1486 [ + - - - : 24 : ABI_MALLOC(gs_hamk_out%dimcprj,(gs_hamk_out%natom*gs_hamk_out%usepaw))
+ - - + -
+ - + -
+ ]
1487 [ - + - - : 8 : if (gs_hamk_out%usepaw==1) gs_hamk_out%dimcprj = gs_hamk_in%dimcprj
- - - - -
- - - - -
- - - - ]
1488 [ - + + - : 24 : ABI_MALLOC(gs_hamk_out%typat,(gs_hamk_out%natom))
+ - + - -
+ - + -
+ ]
1489 [ + - + + : 96 : gs_hamk_out%typat = gs_hamk_in%typat
+ - - + -
- - - - -
+ + ]
1490 [ + - - + : 32 : ABI_MALLOC(gs_hamk_out%gbound_k,(2*gs_hamk_out%mgfft+8,2))
+ - + - +
- + - - +
- + - + ]
1491 [ + - + + : 1456 : gs_hamk_out%gbound_k = gs_hamk_in%gbound_k
+ + + - +
- - + - -
- - - - -
- - - + +
+ + ]
1492 [ - + + - : 32 : ABI_MALLOC(gs_hamk_out%indlmn,(6,gs_hamk_out%lmnmax,gs_hamk_out%ntypat))
- + + - +
- + - - +
- + - + ]
1493 [ + - + + : 720 : gs_hamk_out%indlmn = gs_hamk_in%indlmn
+ + + + +
- + - + -
- + - - -
- - - - -
- - - - -
- + + + +
+ + ]
1494 [ - + + - : 24 : ABI_MALLOC(gs_hamk_out%nattyp,(gs_hamk_out%ntypat))
+ - + - -
+ - + -
+ ]
1495 [ + - + + : 48 : gs_hamk_out%nattyp = gs_hamk_in%nattyp
+ - - + -
- - - - -
+ + ]
1496 [ - + + - : 24 : ABI_MALLOC(gs_hamk_out%nucdipmom,(3,gs_hamk_out%natom))
+ - + - -
+ - + -
+ ]
1497 [ + - + + : 336 : gs_hamk_out%nucdipmom = gs_hamk_in%nucdipmom
+ + + - +
- - + - -
- - - - -
- - - + +
+ + ]
1498 [ - + + - : 24 : ABI_MALLOC(gs_hamk_out%phkxred,(2,gs_hamk_out%natom))
+ - + - -
+ - + -
+ ]
1499 [ + - + + : 256 : gs_hamk_out%phkxred = gs_hamk_in%phkxred
+ + + - +
- - + - -
- - - - -
- - - + +
+ + ]
1500 [ + - - + : 24 : ABI_MALLOC(gs_hamk_out%ph1d,(2,3*(2*gs_hamk_out%mgfft+1)*gs_hamk_out%natom))
+ - + - +
- - + - +
- + ]
1501 [ + - + + : 26656 : gs_hamk_out%ph1d = gs_hamk_in%ph1d
+ + + - +
- - + - -
- - - - -
- - - + +
+ + ]
1502 [ - + + - : 24 : ABI_MALLOC(gs_hamk_out%pspso,(gs_hamk_out%ntypat))
+ - + - -
+ - + -
+ ]
1503 [ + - + + : 48 : gs_hamk_out%pspso = gs_hamk_in%pspso
+ - - + -
- - - - -
+ + ]
1504 [ + + ]: 48 : tmp2i(1:5)=shape(gs_hamk_in%ekb_spin)
1505 [ - + + - : 56 : ABI_MALLOC(gs_hamk_out%ekb_spin,(tmp2i(1),tmp2i(2),tmp2i(3),tmp2i(4),tmp2i(5)))
- + + - -
+ + - - +
+ - - + +
- + - + -
- + - + -
+ ]
1506 [ + - + + : 192 : gs_hamk_out%ekb_spin = gs_hamk_in%ekb_spin
+ + + + +
+ + + + -
+ - + - +
- + - - +
- - - - -
- - - - -
- - - - -
- - - - -
- - + + +
+ + + + +
+ + ]
1507 : 8 : gs_hamk_out%ekb => gs_hamk_out%ekb_spin(:,:,:,:,1)
1508 [ + + ]: 24 : tmp2i(1:2)=shape(gs_hamk_in%sij)
1509 [ + - - - : 32 : ABI_MALLOC(gs_hamk_out%sij,(tmp2i(1),tmp2i(2)))
+ - - - +
- - + - +
- + - + ]
1510 [ - + - + : 16 : gs_hamk_out%sij = gs_hamk_in%sij
- - + - +
- - + - -
- - - - -
- - - - +
- - ]
1511 : :
1512 [ + - ]: 8 : if (associated(gs_hamk_in%gbound_kp,gs_hamk_in%gbound_k)) then
1513 : 8 : gs_hamk_out%gbound_kp => gs_hamk_out%gbound_k
1514 : : else
1515 [ # # # # : 0 : ABI_MALLOC(gs_hamk_out%gbound_kp,(2,gs_hamk_out%natom))
# # # # #
# # # ]
1516 [ # # # # : 0 : gs_hamk_out%gbound_kp = gs_hamk_in%gbound_kp
# # # # #
# ]
1517 : : end if
1518 [ + - ]: 8 : if (associated(gs_hamk_in%phkpxred,gs_hamk_in%phkxred)) then
1519 : 8 : gs_hamk_out%phkpxred => gs_hamk_out%phkxred
1520 : : else
1521 [ # # # # : 0 : ABI_MALLOC(gs_hamk_out%phkpxred,(2,gs_hamk_out%natom))
# # # # #
# # # ]
1522 [ # # # # : 0 : gs_hamk_out%phkpxred = gs_hamk_in%phkpxred
# # # # #
# ]
1523 : : end if
1524 : :
1525 : 8 : call addr_copy(gs_hamk_in%xred,gs_hamk_out%xred)
1526 : 8 : call addr_copy(gs_hamk_in%vectornd,gs_hamk_out%vectornd)
1527 : 8 : call addr_copy(gs_hamk_in%vlocal,gs_hamk_out%vlocal)
1528 : 8 : call addr_copy(gs_hamk_in%vxctaulocal,gs_hamk_out%vxctaulocal)
1529 : 8 : call addr_copy(gs_hamk_in%kinpw_k,gs_hamk_out%kinpw_k)
1530 : 8 : call addr_copy(gs_hamk_in%kinpw_kp,gs_hamk_out%kinpw_kp)
1531 : 8 : call addr_copy(gs_hamk_in%kg_k,gs_hamk_out%kg_k)
1532 : 8 : call addr_copy(gs_hamk_in%kg_kp,gs_hamk_out%kg_kp)
1533 : 8 : call addr_copy(gs_hamk_in%kpg_k,gs_hamk_out%kpg_k)
1534 : 8 : call addr_copy(gs_hamk_in%kpg_kp,gs_hamk_out%kpg_kp)
1535 : 8 : call addr_copy(gs_hamk_in%ffnl_k,gs_hamk_out%ffnl_k)
1536 : 8 : call addr_copy(gs_hamk_in%ffnl_kp,gs_hamk_out%ffnl_kp)
1537 : 8 : call addr_copy(gs_hamk_in%ph3d_k,gs_hamk_out%ph3d_k)
1538 : 8 : call addr_copy(gs_hamk_in%ph3d_kp,gs_hamk_out%ph3d_kp)
1539 : :
1540 : : !For pointers to structured datatypes, have to copy the address
1541 : : !manually because there is no generic addr_copy function for that
1542 [ - + ]: 8 : if (associated(gs_hamk_in%fockcommon)) then
1543 : : #if defined HAVE_FC_ISO_C_BINDING
1544 : 0 : ham_ptr=c_loc(gs_hamk_in%fockcommon)
1545 : 0 : call c_f_pointer(ham_ptr,gs_hamk_out%fockcommon)
1546 : : #else
1547 : : gs_hamk_out%fockcommon=transfer(gs_hamk_in%fockcommon,gs_hamk_out%fockcommon)
1548 : : #endif
1549 : : else
1550 : 8 : nullify(gs_hamk_out%fockcommon)
1551 : : end if
1552 [ - + ]: 8 : if (associated(gs_hamk_in%fockbz)) then
1553 : : #if defined HAVE_FC_ISO_C_BINDING
1554 : 0 : ham_ptr=c_loc(gs_hamk_in%fockbz)
1555 : 0 : call c_f_pointer(ham_ptr,gs_hamk_out%fockbz)
1556 : : #else
1557 : : gs_hamk_out%fockbz=transfer(gs_hamk_in%fockbz,gs_hamk_out%fockbz)
1558 : : #endif
1559 : : else
1560 : 8 : nullify(gs_hamk_out%fockbz)
1561 : : end if
1562 [ - + ]: 8 : if (associated(gs_hamk_in%fockACE_k)) then
1563 : : #if defined HAVE_FC_ISO_C_BINDING
1564 : 0 : ham_ptr=c_loc(gs_hamk_in%fockACE_k)
1565 : 0 : call c_f_pointer(ham_ptr,gs_hamk_out%fockACE_k)
1566 : : #else
1567 : : gs_hamk_out%fockACE_k=transfer(gs_hamk_in%fockACE_k,gs_hamk_out%fockACE_k)
1568 : : #endif
1569 : : else
1570 : 8 : nullify(gs_hamk_out%fockACE_k)
1571 : : end if
1572 : :
1573 : : DBG_EXIT("COLL")
1574 : :
1575 : 8 : end subroutine gsham_copy
1576 : : !!***
1577 : :
1578 : : !----------------------------------------------------------------------
1579 : :
1580 : : !!****f* m_hamiltonian/gsham_load_spin
1581 : : !! NAME
1582 : : !! gsham_load_spin
1583 : : !!
1584 : : !! INPUTS
1585 : : !! isppol=index of current spin
1586 : : !! [vectornd(n4,n5,n6,nvloc,3)]=optional, vector potential of nuclear magnetic dipoles in real space
1587 : : !! [vlocal(n4,n5,n6,nvloc)]=optional, local potential in real space
1588 : : !! [vxctaulocal(n4,n5,n6,nvloc,4)]=optional, derivative of XC energy density with respect
1589 : : !! to kinetic energy density in real space
1590 : : !! [with_nonlocal]=optional, true if non-local factors have to be loaded
1591 : : !!
1592 : : !! FUNCTION
1593 : : !! Setup of the spin-dependent part of the GS Hamiltonian.
1594 : : !!
1595 : : !! SIDE EFFECTS
1596 : : !! Ham<gs_hamiltonian_type>=Structured datatype initialization phase:
1597 : : !! * Quantities that depend spin are initialized.
1598 : : !!
1599 : : !! SOURCE
1600 : :
1601 [ + + + - : 200445 : subroutine gsham_load_spin(Ham,isppol,vectornd,vlocal,vxctaulocal,with_nonlocal)
- + + + +
- - + + +
+ - - + ]
1602 : :
1603 : : !Arguments ------------------------------------
1604 : : !scalars
1605 : : class(gs_hamiltonian_type),intent(inout),target :: Ham
1606 : : integer,intent(in) :: isppol
1607 : : logical,optional,intent(in) :: with_nonlocal
1608 : : !arrays
1609 : : real(dp),optional,intent(in),target :: vectornd(:,:,:,:,:)
1610 : : real(dp),optional,intent(in),target :: vlocal(:,:,:,:),vxctaulocal(:,:,:,:,:)
1611 : :
1612 : : !Local variables-------------------------------
1613 : : integer :: jsppol
1614 : : ! *************************************************************************
1615 : :
1616 : : DBG_ENTER("COLL")
1617 : :
1618 : : !@gs_hamiltonian_type
1619 [ + + ]: 200445 : if (present(vlocal)) then
1620 [ + + - + : 925640 : ABI_CHECK_IEQ(size(vlocal), Ham%n4*Ham%n5*Ham%n6*Ham%nvloc, "Wrong vlocal")
- - ]
1621 : 185128 : Ham%vlocal => vlocal
1622 : : end if
1623 [ + + ]: 200445 : if (present(vxctaulocal)) then
1624 [ + + - + : 4560 : ABI_CHECK_IEQ(size(vxctaulocal), Ham%n4*Ham%n5*Ham%n6*Ham%nvloc*4, "Wrong vxctaulocal")
- - ]
1625 : 760 : Ham%vxctaulocal => vxctaulocal
1626 : : end if
1627 [ + + ]: 200445 : if (present(vectornd)) then
1628 [ + + - + : 6402 : ABI_CHECK_IEQ(size(vectornd), Ham%n4*Ham%n5*Ham%n6*Ham%nvloc*3, "Wrong vectornd")
- - ]
1629 : 1067 : Ham%vectornd => vectornd
1630 : : end if
1631 : :
1632 : : ! Retrieve non-local factors for this spin component
1633 [ + + ]: 200445 : if (present(with_nonlocal)) then
1634 [ + - ]: 198618 : if (with_nonlocal) then
1635 : 198618 : jsppol=min(isppol,size(Ham%ekb_spin,5))
1636 [ + + ]: 198618 : if (jsppol>0) Ham%ekb => Ham%ekb_spin(:,:,:,:,jsppol)
1637 : : end if
1638 : : end if
1639 : :
1640 : : ! Update enl and sij on GPU
1641 : : #if defined HAVE_GPU_CUDA
1642 : : if (Ham%gpu_option==ABI_GPU_LEGACY .or. Ham%gpu_option==ABI_GPU_KOKKOS) then
1643 : : call gpu_update_ham_data(&
1644 : : & Ham%ekb(:,:,:,1), INT(size(Ham%ekb), c_int64_t), &
1645 : : & Ham%sij, INT(size(Ham%sij), c_int64_t), &
1646 : : & Ham%gprimd, INT(size(Ham%gprimd),c_int64_t))
1647 : : end if
1648 : : #endif
1649 : :
1650 : : DBG_EXIT("COLL")
1651 : :
1652 : 200445 : end subroutine gsham_load_spin
1653 : : !!***
1654 : :
1655 : : !!****f* m_hamiltonian/gsham_alloc_fofr_work
1656 : : !! NAME
1657 : : !! gsham_alloc_fofr_work
1658 : : !!
1659 : : !! FUNCTION
1660 : : !!
1661 : : !! INPUTS
1662 : : !!
1663 : : !! SOURCE
1664 : :
1665 : 0 : subroutine gsham_alloc_fofr_work(gs_ham, ndat)
1666 : :
1667 : : !Arguments ------------------------------------
1668 : : class(gs_hamiltonian_type),target,intent(inout) :: gs_ham
1669 : : integer,intent(in) :: ndat
1670 : :
1671 : : !Local variables-------------------------------
1672 : 0 : real(dp), contiguous, pointer :: fofr_work_ptr(:,:,:,:)
1673 : : ! *************************************************************************
1674 : :
1675 : 0 : if (.not. allocated(gs_ham%fofr_work)) then
1676 : : !print *, "first allocation"
1677 : : ! First allocation on CPU and GPU.
1678 [ # # # # : 0 : ABI_MALLOC(gs_ham%fofr_work, (2, gs_ham%n4, gs_ham%n5, gs_ham%n6*ndat))
# # # # #
# # # # #
# # # # #
# # # ]
1679 : 0 : fofr_work_ptr => gs_ham%fofr_work
1680 : : #ifdef HAVE_OPENMP_OFFLOAD
1681 : : !$OMP TARGET ENTER DATA MAP(alloc:fofr_work_ptr) IF (gs_ham%gpu_option==ABI_GPU_OPENMP)
1682 : : #endif
1683 : : end if
1684 : :
1685 : : ! Realloc and remap if the buffer is not large enough.
1686 : : !if (gs_ham%n6*ndat > size(gs_ham%fofr_work, dim=4)) then
1687 : : ! Realloc and remap if buffer size changed.
1688 [ # # ]: 0 : if (gs_ham%n6*ndat /= size(gs_ham%fofr_work, dim=4)) then
1689 : : !print *, "reallocating:", gs_ham%n6*ndat, size(gs_ham%fofr_work, dim=4)
1690 : 0 : fofr_work_ptr => gs_ham%fofr_work
1691 : : #ifdef HAVE_OPENMP_OFFLOAD
1692 : : !$OMP TARGET EXIT DATA MAP(delete:fofr_work_ptr) IF (gs_ham%gpu_option==ABI_GPU_OPENMP)
1693 : : #endif
1694 [ # # # # : 0 : ABI_REMALLOC(gs_ham%fofr_work, (2, gs_ham%n4, gs_ham%n5, gs_ham%n6*ndat))
# # # # #
# # # # #
# # # # #
# # # #
# ]
1695 : 0 : fofr_work_ptr => gs_ham%fofr_work
1696 : : #ifdef HAVE_OPENMP_OFFLOAD
1697 : : !$OMP TARGET ENTER DATA MAP(alloc:fofr_work_ptr) IF (gs_ham%gpu_option==ABI_GPU_OPENMP)
1698 : : #endif
1699 : : end if ! realloc condition.
1700 : :
1701 [ # # ]: 0 : end subroutine gsham_alloc_fofr_work
1702 : : !!***
1703 : :
1704 : : !!****f* m_hamiltonian/gsham_print
1705 : : !! NAME
1706 : : !! gsham_print
1707 : : !!
1708 : : !! FUNCTION
1709 : : !! Print info on the object.
1710 : : !!
1711 : : !! INPUTS
1712 : : !! units=Unit numbers for output.
1713 : : !! header=Header string
1714 : : !! prtvol=verbosity level
1715 : : !!
1716 : : !! SOURCE
1717 : :
1718 : 0 : subroutine gsham_print(gs_ham, units, header, prtvol)
1719 : :
1720 : : !Arguments ------------------------------------
1721 : : class(gs_hamiltonian_type),intent(in) :: gs_ham
1722 : : integer,intent(in) :: units(:), prtvol
1723 : : character(len=*),intent(in) :: header
1724 : :
1725 : : !Local variables-------------------------------
1726 : : character(len=5000) :: msg
1727 : : ! *************************************************************************
1728 : :
1729 [ # # ]: 0 : if (len_trim(header) == 0) then
1730 : 0 : msg = ' ==== Info on the gs_hamiltonian_type object ==== '
1731 : : else
1732 [ # # # # ]: 0 : msg = ' ==== '//trim(adjustl(header))//' ==== '
1733 : : end if
1734 : 0 : call wrtout(units, msg)
1735 : :
1736 : 0 : call wrtout(units, sjoin(" kpt_k:", ktoa(gs_ham%kpt_k)))
1737 : 0 : call wrtout(units, sjoin(" kpt_kp:", ktoa(gs_ham%kpt_kp)))
1738 : 0 : call wrtout(units, sjoin(" npw_k:", itoa(gs_ham%npw_k)))
1739 : 0 : call wrtout(units, sjoin(" npw_kp:", itoa(gs_ham%npw_kp)))
1740 : 0 : call wrtout(units, sjoin(" associated(vlocal):", yesno(associated(gs_ham%vlocal))))
1741 : :
1742 : : if (prtvol > 0) then
1743 : : end if
1744 : :
1745 [ # # ]: 0 : end subroutine gsham_print
1746 : : !!***
1747 : :
1748 : : !----------------------------------------------------------------------
1749 : :
1750 : : !!****f* m_hamiltonian/rfham_free
1751 : : !! NAME
1752 : : !! rfham_free
1753 : : !!
1754 : : !! FUNCTION
1755 : : !! Clean and destroy rf_hamiltonian_type datastructure
1756 : : !!
1757 : : !! SOURCE
1758 : :
1759 : 2156061 : subroutine rfham_free(rf_Ham)
1760 : :
1761 : : !Arguments ------------------------------------
1762 : : class(rf_hamiltonian_type),intent(inout) :: rf_Ham
1763 : : ! *************************************************************************
1764 : :
1765 : : DBG_ENTER("COLL")
1766 : :
1767 : : !@rf_hamiltonian_type
1768 : :
1769 : : ! Real arrays
1770 [ + - ]: 2156061 : ABI_SFREE(rf_Ham%e1kbfr_spin)
1771 [ + - ]: 2156061 : ABI_SFREE(rf_Ham%e1kbsc_spin)
1772 : :
1773 : : ! Real pointers
1774 [ + + ]: 2156061 : if (associated(rf_Ham%dkinpw_k)) nullify(rf_Ham%dkinpw_k)
1775 [ + + ]: 2156061 : if (associated(rf_Ham%dkinpw_kp)) nullify(rf_Ham%dkinpw_kp)
1776 [ + + ]: 2156061 : if (associated(rf_Ham%ddkinpw_k)) nullify(rf_Ham%ddkinpw_k)
1777 [ + + ]: 2156061 : if (associated(rf_Ham%ddkinpw_kp)) nullify(rf_Ham%ddkinpw_kp)
1778 [ + + ]: 2156061 : if (associated(rf_Ham%vectornd)) nullify(rf_Ham%vectornd)
1779 [ + + ]: 2156061 : if (associated(rf_Ham%vlocal1)) nullify(rf_Ham%vlocal1)
1780 [ + + ]: 2156061 : if (associated(rf_Ham%vxctaulocal)) nullify(rf_Ham%vxctaulocal)
1781 [ + + ]: 2156061 : if (associated(rf_Ham%e1kbfr)) nullify(rf_Ham%e1kbfr)
1782 [ + + ]: 2156061 : if (associated(rf_Ham%e1kbsc)) nullify(rf_Ham%e1kbsc)
1783 : :
1784 : : DBG_EXIT("COLL")
1785 : :
1786 : 2156061 : end subroutine rfham_free
1787 : : !!***
1788 : :
1789 : : !----------------------------------------------------------------------
1790 : :
1791 : : !!****f* m_hamiltonian/rfham_init
1792 : : !! NAME
1793 : : !! rfham_init
1794 : : !!
1795 : : !! FUNCTION
1796 : : !! Creation method for the rf_hamiltonian_type structure.
1797 : : !! It allocates memory and initializes all quantities that do not depend on the k-point or spin.
1798 : : !!
1799 : : !! INPUTS
1800 : : !! [comm_atom]=optional, MPI communicator over atoms
1801 : : !! cplex_paw=1 if all on-site PAW quantities are real (GS), 2 if they are complex (RF)
1802 : : !! gs_Ham<gs_hamiltonian_type>=Structured datatype containing data for ground-state Hamiltonian at (k+q)
1803 : : !! [has_e1kbsc]=optional, true if rf_Ham%e1kbsc has to be initialized.
1804 : : !! e1kbsc contains the self-consistent 1st-order PAW Dij coefficients (depending on VHxc^(1))
1805 : : !! ipert=index of perturbation
1806 : : !! [mpi_atmtab(:)]=optional, indexes of the atoms treated by current proc
1807 : : !! [mpi_spintab(2)]=optional, flags defining the spin(s) treated be current process:
1808 : : !! mpi_spintab(1)=1 if non-polarized or spin-up treated
1809 : : !! mpi_spintab(2)=1 if polarized and spin-dn treated
1810 : : !! [paw_ij1(:)<paw_ij_type>]=Various 1st-order arrays given on (i,j) (partial waves)
1811 : : !! channels (paw_ij1%dij and paw_ij1%difr only used here).
1812 : : !!
1813 : : !! SIDE EFFECTS
1814 : : !! rf_Ham<rf_hamiltonian_type>=Structured datatype almost completely initialized:
1815 : : !! * Basic variables and dimensions are transferred to the structure.
1816 : : !! * All pointers are allocated with correct dimensions.
1817 : : !! * Quantities that do not depend on the k-point or spin are initialized.
1818 : : !!
1819 : : !! SOURCE
1820 : :
1821 : 2229761 : subroutine rfham_init(rf_ham, cplex, gs_Ham, ipert,&
1822 [ + + + + : 2229761 : comm_atom, mpi_atmtab, mpi_spintab, paw_ij1, has_e1kbsc) ! optional arguments
- + + + +
+ - + ]
1823 : :
1824 : : !Arguments ------------------------------------
1825 : : !scalars
1826 : : class(rf_hamiltonian_type),intent(inout),target :: rf_Ham
1827 : : integer,intent(in) :: cplex,ipert
1828 : : integer,intent(in),optional :: comm_atom
1829 : : logical,intent(in),optional :: has_e1kbsc
1830 : : type(gs_hamiltonian_type),intent(in) :: gs_Ham
1831 : : !arrays
1832 : : integer,optional,intent(in) :: mpi_atmtab(:),mpi_spintab(2)
1833 : : type(paw_ij_type),optional,intent(in) :: paw_ij1(:)
1834 : :
1835 : : !Local variables-------------------------------
1836 : : !scalars
1837 : : integer :: cplex_dij1,isp,jsp,my_comm_atom,my_nsppol
1838 : : logical :: has_e1kbsc_
1839 : : !arrays
1840 : : integer :: my_spintab(2)
1841 : 2229761 : real(dp),allocatable,target :: e1kb_tmp(:,:,:,:)
1842 : : ! *************************************************************************
1843 : :
1844 : : DBG_ENTER("COLL")
1845 : :
1846 : : !@rf_hamiltonian_type
1847 : :
1848 : : !Manage optional parameters
1849 [ + + ]: 2229761 : has_e1kbsc_=.false.;if (present(has_e1kbsc)) has_e1kbsc_=has_e1kbsc
1850 [ + + ]: 2229761 : my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
1851 [ + + + + ]: 4475094 : my_spintab=0;my_spintab(1:gs_Ham%nsppol)=1;if(present(mpi_spintab)) my_spintab=mpi_spintab
1852 [ + + + + ]: 6689283 : my_nsppol=count(my_spintab==1)
1853 : :
1854 : 2229761 : rf_Ham%cplex =cplex
1855 : 2229761 : rf_Ham%n4 =gs_Ham%n4
1856 : 2229761 : rf_Ham%n5 =gs_Ham%n5
1857 : 2229761 : rf_Ham%n6 =gs_Ham%n6
1858 : 2229761 : rf_Ham%nvloc =gs_Ham%nvloc
1859 : 2229761 : rf_Ham%nsppol =gs_Ham%nsppol
1860 : 2229761 : rf_Ham%nspinor =gs_Ham%nspinor
1861 : 2229761 : rf_Ham%zora =gs_Ham%zora
1862 : :
1863 : 2229761 : rf_Ham%dime1kb1=0
1864 : 2229761 : rf_Ham%dime1kb2=gs_Ham%dimekb2
1865 [ + + + + : 2229761 : if (gs_Ham%usepaw==1.and.ipert/=gs_Ham%natom+1.and.ipert/=gs_Ham%natom+10) then
+ + ]
1866 [ + + + + : 62500 : cplex_dij1=1;if ((gs_Ham%nspinor==2).or.any(abs(gs_Ham%nucdipmom)>tol8)) cplex_dij1=2
+ - + - -
+ ]
1867 : 4192 : rf_Ham%dime1kb1=cplex_dij1*(gs_Ham%lmnmax*(gs_Ham%lmnmax+1))/2
1868 : : end if
1869 : :
1870 : : ! Allocate the arrays of the 1st-order Hamiltonian
1871 : : ! We preload here 1st-order non-local factors in order to avoid
1872 : : ! a communication over atoms inside the loop over spins.
1873 [ + + + + ]: 2229761 : if (gs_Ham%usepaw==1.and.rf_Ham%dime1kb1>0) then
1874 : : if ((ipert>=1.and.ipert<=gs_Ham%natom).or.ipert==gs_Ham%natom+2.or.&
1875 [ + - + + : 4192 : ipert==gs_Ham%natom+3.or.ipert==gs_Ham%natom+4.or.ipert==gs_Ham%natom+11) then
+ + + + +
+ - + ]
1876 : :
1877 [ - + + - : 29344 : ABI_MALLOC(rf_Ham%e1kbfr_spin,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex,my_nsppol))
- + + - -
+ + - - +
+ - + + +
- + - + +
- + - + -
+ ]
1878 [ + + + + : 513811 : rf_Ham%e1kbfr_spin=zero
+ + + + +
+ ]
1879 [ + + ]: 4192 : if (has_e1kbsc_) then
1880 [ - + + - : 20951 : ABI_MALLOC(rf_Ham%e1kbsc_spin,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex,my_nsppol))
- + + - -
+ + - - +
+ - - + +
- + - + -
- + - + -
+ ]
1881 [ + + + + : 496875 : rf_Ham%e1kbsc_spin=zero
+ + + + +
+ ]
1882 : : end if
1883 : :
1884 [ + + ]: 4192 : if (present(paw_ij1)) then
1885 : :
1886 [ - + ]: 3108 : if (my_nsppol<rf_Ham%nsppol) then
1887 [ # # # # : 0 : ABI_MALLOC(e1kb_tmp,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex))
# # # # #
# # # # #
# # # # #
# # # #
# ]
1888 : : end if
1889 : :
1890 : : ! === Frozen term
1891 : 3108 : jsp=0
1892 [ + + ]: 6936 : do isp=1,rf_Ham%nsppol
1893 [ + - ]: 3828 : if (my_spintab(isp)==1) then
1894 : 3828 : jsp=jsp+1 ; rf_Ham%e1kbfr => rf_Ham%e1kbfr_spin(:,:,:,:,jsp)
1895 : : else
1896 : 0 : rf_Ham%e1kbfr => e1kb_tmp
1897 : : end if
1898 [ + + ]: 6936 : if (present(mpi_atmtab)) then
1899 : 1314 : call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbfr=rf_Ham%e1kbfr,mpi_atmtab=mpi_atmtab)
1900 : : else
1901 : 2514 : call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbfr=rf_Ham%e1kbfr)
1902 : : end if
1903 : : end do
1904 : :
1905 : : ! === Self-consistent term
1906 [ + + ]: 3108 : if (has_e1kbsc_) then
1907 : 2993 : jsp=0
1908 [ + + ]: 6706 : do isp=1,rf_Ham%nsppol
1909 [ + - ]: 3713 : if (my_spintab(isp)==1) then
1910 : 3713 : jsp=jsp+1 ; rf_Ham%e1kbsc => rf_Ham%e1kbsc_spin(:,:,:,:,jsp)
1911 : : else
1912 : 0 : rf_Ham%e1kbsc => e1kb_tmp
1913 : : end if
1914 [ + + ]: 6706 : if (present(mpi_atmtab)) then
1915 : 1254 : call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbsc=rf_Ham%e1kbsc,mpi_atmtab=mpi_atmtab)
1916 : : else
1917 : 2459 : call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbsc=rf_Ham%e1kbsc)
1918 : : end if
1919 : : end do
1920 : : end if
1921 : :
1922 [ - + ]: 3108 : if (my_nsppol<rf_Ham%nsppol) then
1923 [ # # ]: 0 : ABI_FREE(e1kb_tmp)
1924 : : end if
1925 : :
1926 : : end if
1927 : : end if
1928 : : end if
1929 : :
1930 [ + + ]: 2229761 : if (.not.allocated(rf_Ham%e1kbfr_spin)) then
1931 [ - + - + ]: 2151869 : ABI_MALLOC(rf_Ham%e1kbfr_spin,(0,0,0,0,0))
1932 : : end if
1933 [ + + ]: 2229761 : if (.not.allocated(rf_Ham%e1kbsc_spin)) then
1934 [ - + - + ]: 2153068 : ABI_MALLOC(rf_Ham%e1kbsc_spin,(0,0,0,0,0))
1935 : : end if
1936 : 2229761 : nullify(rf_Ham%e1kbfr)
1937 : 2229761 : nullify(rf_Ham%e1kbsc)
1938 : :
1939 : : DBG_EXIT("COLL")
1940 : :
1941 [ - + ]: 2229761 : end subroutine rfham_init
1942 : : !!***
1943 : :
1944 : : !----------------------------------------------------------------------
1945 : :
1946 : : !!****f* m_hamiltonian/rfham_load_spin
1947 : : !! NAME
1948 : : !! rfham_load_spin
1949 : : !!
1950 : : !! FUNCTION
1951 : : !! Setup of the spin-dependent part of the 1st- and 2nd- order Hamiltonian.
1952 : : !!
1953 : : !! INPUTS
1954 : : !! isppol=index of current spin
1955 : : !! [vectornd(n4,n5,n6,nvloc)]=optional, vector potential of nuclear magnetic dipoles in real space in
1956 : : !! ddk direction idir
1957 : : !! [vlocal1(cplex*n4,n5,n6,nvloc)]=optional, 1st-order local potential in real space
1958 : : !! [vxctaulocal(n4,n5,n6,nvloc,4)]=optional, deriv of e_XC wrt kin energy, for mGGA
1959 : : !! [with_nonlocal]=optional, true if non-local factors have to be loaded
1960 : : !!
1961 : : !! SIDE EFFECTS
1962 : : !! rf_Ham<rf_hamiltonian_type>=Structured datatype initialization phase:
1963 : : !! * Quantities that depend on spin are initialized.
1964 : : !!
1965 : : !! SOURCE
1966 : :
1967 : 904128 : subroutine rfham_load_spin(rf_Ham, isppol, &
1968 [ + + + - : 904128 : vectornd, vlocal1, vxctaulocal, with_nonlocal) ! optional
- + + + +
- - + + +
+ - - + ]
1969 : :
1970 : : !Arguments ------------------------------------
1971 : : !scalars
1972 : : class(rf_hamiltonian_type),intent(inout),target :: rf_Ham
1973 : : integer,intent(in) :: isppol
1974 : : logical,optional,intent(in) :: with_nonlocal
1975 : : !arrays
1976 : : real(dp),optional,target,intent(in) :: vlocal1(:,:,:,:)
1977 : : real(dp),optional,target,intent(in) :: vectornd(:,:,:,:)
1978 : : real(dp),optional,target,intent(in) :: vxctaulocal(:,:,:,:,:)
1979 : :
1980 : : !Local variables-------------------------------
1981 : : integer :: jsppol
1982 : : ! *************************************************************************
1983 : :
1984 : : DBG_ENTER("COLL")
1985 : :
1986 : : !@rf_hamiltonian_type
1987 : :
1988 [ + + ]: 904128 : if (present(vlocal1)) then
1989 [ + + - + : 4509405 : ABI_CHECK_IEQ(size(vlocal1), rf_Ham%cplex*rf_Ham%n4*rf_Ham%n5*rf_Ham%n6*rf_Ham%nvloc, "Wrong vlocal1")
- - ]
1990 : 901881 : rf_Ham%vlocal1 => vlocal1
1991 : : end if
1992 : :
1993 [ + + ]: 904128 : if (present(vectornd)) then
1994 [ + + - + : 2610 : ABI_CHECK_IEQ(size(vectornd), rf_Ham%cplex*rf_Ham%n4*rf_Ham%n5*rf_Ham%n6*rf_Ham%nvloc, "Wrong vectornd")
- - ]
1995 : 522 : rf_Ham%vectornd => vectornd
1996 : : end if
1997 : :
1998 [ + + ]: 904128 : if (present(vxctaulocal)) then
1999 [ + + - + : 288 : ABI_CHECK_IEQ(size(vxctaulocal), rf_Ham%n4*rf_Ham%n5*rf_Ham%n6*rf_Ham%nvloc*4, "Wrong vxctaulocal")
- - ]
2000 : 48 : rf_Ham%vxctaulocal => vxctaulocal
2001 : : end if
2002 : :
2003 : : ! Retrieve non-local factors for this spin component
2004 [ + + ]: 904128 : if (present(with_nonlocal)) then
2005 [ + - ]: 903311 : if (with_nonlocal) then
2006 [ + + + + ]: 5419866 : if (size(rf_Ham%e1kbfr_spin)>0) then
2007 : 3828 : jsppol=min(isppol,size(rf_Ham%e1kbfr_spin,5))
2008 [ + - ]: 3828 : if (jsppol>0) rf_Ham%e1kbfr => rf_Ham%e1kbfr_spin(:,:,:,:,jsppol)
2009 : : end if
2010 [ + + + + ]: 5419866 : if (size(rf_Ham%e1kbsc_spin)>0) then
2011 : 3713 : jsppol=min(isppol,size(rf_Ham%e1kbsc_spin,5))
2012 [ + - ]: 3713 : if (jsppol>0) rf_Ham%e1kbsc => rf_Ham%e1kbsc_spin(:,:,:,:,jsppol)
2013 : : end if
2014 : : end if
2015 : : end if
2016 : :
2017 : : DBG_EXIT("COLL")
2018 : :
2019 : 904128 : end subroutine rfham_load_spin
2020 : : !!***
2021 : :
2022 : : !----------------------------------------------------------------------
2023 : :
2024 : : !!****f* m_hamiltonian/rfham_load_k
2025 : : !! NAME
2026 : : !! rfham_load_k
2027 : : !!
2028 : : !! FUNCTION
2029 : : !! Setup of the k-dependent part of the 1st- and 2nd- order Hamiltonian
2030 : : !!
2031 : : !! INPUTS
2032 : : !! [dkinpw_k]=1st derivative of the (modified) kinetic energy for each plane wave
2033 : : !! [ddkinpw_k]=2nd derivative of the (modified) kinetic energy for each plane wave
2034 : : !! [npw_k]=number of plane waves
2035 : : !!
2036 : : !! SIDE EFFECTS
2037 : : !! rf_Ham<rf_hamiltonian_type>=structured datatype completed with k-dependent quantities.
2038 : : !! Quantities at k^prime are set equal to quantities at k.
2039 : : !!
2040 : : !! SOURCE
2041 : :
2042 : 1324372 : subroutine rfham_load_k(rf_Ham, &
2043 [ + + + - : 1324372 : dkinpw_k, ddkinpw_k, npw_k) ! optional.
- + + + +
- - + ]
2044 : :
2045 : : !Arguments ------------------------------------
2046 : : !scalars
2047 : : class(rf_hamiltonian_type),intent(inout),target :: rf_Ham
2048 : : integer,intent(in),optional :: npw_k
2049 : : !arrays
2050 : : real(dp),intent(in),optional,target :: dkinpw_k(:),ddkinpw_k(:)
2051 : : ! *************************************************************************
2052 : :
2053 : : DBG_ENTER("COLL")
2054 : :
2055 : : !@gs_hamiltonian_type
2056 : :
2057 : : !k-dependent scalars
2058 [ + + ]: 1324372 : if (present(npw_k)) then
2059 : 1306156 : rf_Ham%npw_k = npw_k
2060 : 1306156 : rf_Ham%npw_kp = npw_k
2061 : : end if
2062 : :
2063 : : !Pointers to k-dependent quantities
2064 [ + + ]: 1324372 : if (present(dkinpw_k)) then
2065 : 1313212 : rf_Ham%dkinpw_k => dkinpw_k
2066 : 1313212 : rf_Ham%dkinpw_kp => dkinpw_k
2067 : : end if
2068 [ + + ]: 1324372 : if (present(ddkinpw_k)) then
2069 : 18216 : rf_Ham%ddkinpw_k => ddkinpw_k
2070 : 18216 : rf_Ham%ddkinpw_kp => ddkinpw_k
2071 : : end if
2072 : :
2073 : : DBG_EXIT("COLL")
2074 : :
2075 : 1324372 : end subroutine rfham_load_k
2076 : : !!***
2077 : :
2078 : : !----------------------------------------------------------------------
2079 : :
2080 : : !!****f* m_hamiltonian/pawdij2ekb
2081 : : !! NAME
2082 : : !! pawdij2ekb
2083 : : !!
2084 : : !! FUNCTION
2085 : : !! Transfer PAW Dij (on-site GS Hamiltonian) values
2086 : : !! from paw_ij datastructure to ekb array
2087 : : !!
2088 : : !! INPUTS
2089 : : !!
2090 : : !! OUTPUT
2091 : : !!
2092 : : !! SOURCE
2093 : :
2094 [ - + - + : 25668 : subroutine pawdij2ekb(ekb,paw_ij,isppol,comm_atom,mpi_atmtab)
+ + + - -
+ ]
2095 : :
2096 : : !Arguments ------------------------------------
2097 : : !scalars
2098 : : integer,intent(in) :: isppol,comm_atom
2099 : : !arrays
2100 : : integer,intent(in),optional,target :: mpi_atmtab(:)
2101 : : real(dp),intent(out) :: ekb(:,:,:,:)
2102 : : type(paw_ij_type),intent(in) :: paw_ij(:)
2103 : :
2104 : : !Local variables-------------------------------
2105 : : !scalars
2106 : : integer :: dimdij,dimekb1,dimekb3,dimekb4,iatom,iatom_tot,ierr,ii,isp,ispden,my_natom,natom,qphase
2107 : : logical :: my_atmtab_allocated,paral_atom
2108 : : !arrays
2109 : 25668 : integer,pointer :: my_atmtab(:)
2110 : : ! *************************************************************************
2111 : :
2112 : : DBG_ENTER("COLL")
2113 : :
2114 [ + + + + : 8422503 : ekb=zero
+ + + + ]
2115 : :
2116 : : !Set up parallelism over atoms
2117 : 25668 : natom=size(ekb,2); my_natom=size(paw_ij)
2118 : 25668 : paral_atom=(xmpi_comm_size(comm_atom)>1)
2119 [ + + ]: 25668 : nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
2120 : 25668 : call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
2121 : :
2122 : : !Retrieve PAW Dij coefficients for this spin component
2123 [ + + ]: 25668 : if (my_natom>0) then
2124 [ + - ]: 23682 : if (allocated(paw_ij(1)%dij)) then
2125 : 23682 : dimekb1=size(ekb,1) ; dimekb3=size(ekb,3) ; dimekb4=size(ekb,4)
2126 : 23682 : qphase=paw_ij(1)%qphase
2127 [ - + ]: 23682 : ABI_CHECK(qphase <= dimekb4, 'paw_ij%qphase>dimekb4!')
2128 [ + + ]: 47364 : do ii=1,qphase
2129 [ + + ]: 78024 : do ispden=1,dimekb3
2130 [ + + ]: 30660 : isp=isppol; if (dimekb3==4) isp=ispden
2131 [ + + ]: 120318 : do iatom=1,my_natom
2132 [ + + ]: 65976 : iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom)
2133 : 65976 : dimdij=paw_ij(iatom)%cplex_dij*paw_ij(iatom)%lmn2_size
2134 [ - + ]: 65976 : ABI_CHECK(dimdij <= dimekb1, 'Size of paw_ij%dij>dimekb1!')
2135 [ + + ]: 6137511 : ekb(1:dimdij,iatom_tot,ispden,ii)=paw_ij(iatom)%dij(1+(ii-1)*dimdij:ii*dimdij,isp)
2136 : : end do
2137 : : end do
2138 : : end do
2139 : : end if
2140 : : end if
2141 : :
2142 : : ! Communication in case of distribution over atomic sites
2143 [ + + ]: 25668 : if (paral_atom) call xmpi_sum(ekb,comm_atom,ierr)
2144 : :
2145 : : ! Destroy atom table used for parallelism
2146 : 25668 : call free_my_atmtab(my_atmtab,my_atmtab_allocated)
2147 : :
2148 : : DBG_EXIT("COLL")
2149 : :
2150 : 25668 : end subroutine pawdij2ekb
2151 : : !!***
2152 : :
2153 : : !----------------------------------------------------------------------
2154 : :
2155 : : !!****f* m_hamiltonian/pawdij2e1kb
2156 : : !! NAME
2157 : : !! pawdij2e1kb
2158 : : !!
2159 : : !! FUNCTION
2160 : : !! Transfer PAW Dij (on-site RF Hamiltonian) values
2161 : : !! from paw_ij datastructure to e1kb array
2162 : : !!
2163 : : !! INPUTS
2164 : : !!
2165 : : !! OUTPUT
2166 : : !!
2167 : : !! SOURCE
2168 : :
2169 [ + + + - : 11192 : subroutine pawdij2e1kb(paw_ij1,isppol,comm_atom,mpi_atmtab,e1kbfr,e1kbsc)
- + + + +
- - + - +
+ + + + -
+ ]
2170 : :
2171 : : !Arguments ------------------------------------
2172 : : !scalars
2173 : : integer,intent(in) :: isppol,comm_atom
2174 : : !arrays
2175 : : integer,intent(in),optional,target :: mpi_atmtab(:)
2176 : : real(dp),optional,intent(out) :: e1kbfr(:,:,:,:),e1kbsc(:,:,:,:)
2177 : : type(paw_ij_type),intent(in) :: paw_ij1(:)
2178 : :
2179 : : !Local variables-------------------------------
2180 : : !scalars
2181 : : integer :: dimdij1,dime1kb1,dime1kb3,dime1kb4,iatom,iatom_tot,ierr,isp,ispden
2182 : : integer :: my_natom,natom,qphase
2183 : : logical :: my_atmtab_allocated,paral_atom
2184 : : !arrays
2185 : 11192 : integer,pointer :: my_atmtab(:)
2186 : : ! *************************************************************************
2187 : :
2188 : : DBG_ENTER("COLL")
2189 : :
2190 : 0 : if ((.not.present(e1kbfr)).and.(.not.present(e1kbsc))) return
2191 [ + + ]: 11192 : if (present(e1kbfr)) then
2192 [ + + + + : 1088748 : e1kbfr=zero ; natom=size(e1kbfr,2)
+ + + + ]
2193 : : end if
2194 [ + + ]: 11192 : if (present(e1kbsc)) then
2195 [ + + + + : 493882 : e1kbsc=zero ; natom=size(e1kbsc,2)
+ + + + ]
2196 : : end if
2197 : :
2198 : : ! Set up parallelism over atoms
2199 : 11192 : my_natom=size(paw_ij1) ; paral_atom=(xmpi_comm_size(comm_atom)>1)
2200 [ + + ]: 11192 : nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
2201 : 11192 : call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
2202 : :
2203 : : ! Retrieve 1st-order PAW Dij coefficients for this spin component (frozen)
2204 [ + + + + ]: 11192 : if (my_natom>0.and.present(e1kbfr)) then
2205 [ + - ]: 7417 : if (allocated(paw_ij1(1)%dijfr)) then
2206 : 7417 : dime1kb1=size(e1kbfr,1) ; dime1kb3=size(e1kbfr,3) ; dime1kb4=size(e1kbfr,4)
2207 [ - + ]: 7417 : ABI_CHECK_IEQ(paw_ij1(1)%qphase, dime1kb4,'BUG in pawdij2e1kb (1)!')
2208 [ + + ]: 14834 : do ispden=1,dime1kb3
2209 [ - + ]: 7417 : isp=isppol;if (dime1kb3==4) isp=ispden
2210 [ + + ]: 30468 : do iatom=1,my_natom
2211 [ + + ]: 15634 : iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom)
2212 : 15634 : qphase=paw_ij1(iatom)%qphase
2213 : 15634 : dimdij1=paw_ij1(iatom)%cplex_dij*paw_ij1(iatom)%lmn2_size
2214 [ - + ]: 15634 : ABI_CHECK(dimdij1<=dime1kb1,'BUG: size of paw_ij1%dij>dime1kb1!')
2215 [ + + ]: 620200 : e1kbfr(1:dimdij1,iatom_tot,ispden,1)=paw_ij1(iatom)%dijfr(1:dimdij1,isp)
2216 [ + + + + ]: 61499 : if (qphase==2) e1kbfr(1:dimdij1,iatom_tot,ispden,2)=paw_ij1(iatom)%dijfr(dimdij1+1:2*dimdij1,isp)
2217 : : end do
2218 : : end do
2219 : : end if
2220 : : end if
2221 : :
2222 : : ! Retrieve 1st-order PAW Dij coefficients for this spin component (self-consistent)
2223 [ + + + + ]: 11192 : if (my_natom>0.and.present(e1kbsc)) then
2224 [ + - + - ]: 3693 : if (allocated(paw_ij1(1)%dijfr).and.allocated(paw_ij1(1)%dij)) then
2225 : 3693 : dime1kb1=size(e1kbsc,1) ; dime1kb3=size(e1kbsc,3) ; dime1kb4=size(e1kbsc,4)
2226 [ - + ]: 3693 : ABI_CHECK_IEQ(paw_ij1(1)%qphase, dime1kb4, 'BUG in pawdij2e1kb (1)!')
2227 [ + + ]: 7386 : do ispden=1,dime1kb3
2228 [ - + ]: 3693 : isp=isppol;if (dime1kb3==4) isp=ispden
2229 [ + + ]: 15057 : do iatom=1,my_natom
2230 [ + + ]: 7671 : iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom)
2231 : 7671 : qphase=paw_ij1(iatom)%qphase
2232 : 7671 : dimdij1=paw_ij1(iatom)%cplex_dij*paw_ij1(iatom)%lmn2_size
2233 [ - + ]: 7671 : ABI_CHECK(dimdij1<=dime1kb1, 'BUG: size of paw_ij1%dij>dime1kb1!')
2234 : : e1kbsc(1:dimdij1,iatom_tot,ispden,1)=paw_ij1(iatom)%dij (1:dimdij1,isp) &
2235 [ + + ]: 302475 : & -paw_ij1(iatom)%dijfr(1:dimdij1,isp)
2236 [ + + ]: 7671 : if (qphase==2) e1kbsc(1:dimdij1,iatom_tot,ispden,2)=paw_ij1(iatom)%dij (dimdij1+1:2*dimdij1,isp) &
2237 [ + + ]: 28245 : & -paw_ij1(iatom)%dijfr(dimdij1+1:2*dimdij1,isp)
2238 : : end do
2239 : : end do
2240 : : end if
2241 : : end if
2242 : :
2243 : : ! Communication in case of distribution over atomic sites
2244 [ + + ]: 11192 : if (paral_atom) then
2245 [ + + ]: 4494 : if (present(e1kbfr)) call xmpi_sum(e1kbfr,comm_atom,ierr)
2246 [ + + ]: 4494 : if (present(e1kbsc)) call xmpi_sum(e1kbsc,comm_atom,ierr)
2247 : : end if
2248 : :
2249 : : ! Destroy atom table used for parallelism
2250 : 11192 : call free_my_atmtab(my_atmtab,my_atmtab_allocated)
2251 : :
2252 : : DBG_EXIT("COLL")
2253 : :
2254 [ - + ]: 11192 : end subroutine pawdij2e1kb
2255 : : !!***
2256 : :
2257 : : !!****f* ABINIT/gspot_transgrid_and_pack
2258 : : !! NAME
2259 : : !! gspot_transgrid_and_pack
2260 : : !!
2261 : : !! FUNCTION
2262 : : !! Set up local potential vlocal on the coarse FFT mesh with proper dimensioning from vtrial given on the fine mesh.
2263 : : !! Also take into account the spin.
2264 : : !!
2265 : : !! INPUTS
2266 : : !! isppol=Spin polarization.
2267 : : !! usepaw=1 if PAW
2268 : : !! paral_kgb: 1 if paral_kgb
2269 : : !! nfft=(effective) number of FFT grid points on the coarse mesh (for this processor)
2270 : : !! ngfft(18)contain all needed information about 3D FFT, for the coarse FFT mesh. see ~abinit/doc/variables/vargs.htm#ngfft
2271 : : !! nfftf=(effective) number of FFT grid points on the fine mesh (for this processor)
2272 : : !! nvloc==1 if nspden <=2, nvloc==4 for nspden==4,
2273 : : !! nspden=Number of spin density components.
2274 : : !! ncomp=Number of extra components in vtrial and vlocal (e.g. 1 if LDA/GGA pot, 4 for Meta-GGA, etc).
2275 : : !! pawfgr <type(pawfgr_type)>=fine grid parameters and related data
2276 : : !! vtrial(nfftf,nspden)=INPUT potential Vtrial(r).
2277 : : !! mpi_enreg=information about MPI parallelization
2278 : : !!
2279 : : !! OUTPUT
2280 : : !! vlocal(n4,n5,n6,nvloc,ncomp): Potential on the coarse grid.
2281 : : !!
2282 : : !! SOURCE
2283 : :
2284 : 56053 : subroutine gspot_transgrid_and_pack(isppol, usepaw, paral_kgb, nfft, ngfft, nfftf, &
2285 : 56053 : nspden, nvloc, ncomp, pawfgr, mpi_enreg, vtrial, vlocal)
2286 : :
2287 : : !Arguments -------------------------------
2288 : : integer,intent(in) :: isppol, nspden, ncomp, usepaw, paral_kgb, nfft, nfftf, nvloc
2289 : : type(pawfgr_type), intent(in) :: pawfgr
2290 : : type(MPI_type), intent(in) :: mpi_enreg
2291 : : !arrays
2292 : : integer,intent(in) :: ngfft(18)
2293 : : real(dp),intent(inout) :: vtrial(nfftf, nspden, ncomp)
2294 : : real(dp),intent(out) :: vlocal(ngfft(4), ngfft(5), ngfft(6), nvloc, ncomp)
2295 : :
2296 : : !Local variables-------------------------------
2297 : : !scalars
2298 : : integer :: n1,n2,n3,n4,n5,n6,ispden,ic
2299 : : real(dp) :: rhodum(1)
2300 : 56053 : real(dp),allocatable :: cgrvtrial(:,:), vlocal_tmp(:,:,:)
2301 : : ! *************************************************************************
2302 : :
2303 : : ! Coarse mesh.
2304 : 56053 : n1=ngfft(1); n2=ngfft(2); n3=ngfft(3)
2305 : 56053 : n4=ngfft(4); n5=ngfft(5); n6=ngfft(6)
2306 : :
2307 : : ! Set up local potential vlocal with proper dimensioning, from vtrial
2308 : : ! Also take into account the spin.
2309 [ + + ]: 56053 : if (nspden /= 4) then
2310 [ + + + + ]: 53523 : if (usepaw == 0 .or. pawfgr%usefinegrid == 0) then
2311 : : ! Fine mesh == Coarse mesh.
2312 [ + + ]: 70239 : do ic=1,ncomp
2313 : 70239 : call fftpac(isppol,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,vtrial(:,:,ic),vlocal(:,:,:,:,ic),2)
2314 : : end do
2315 : : else
2316 : : ! Transfer from fine mesh to coarse and then pack data
2317 [ - + + - : 74996 : ABI_MALLOC(cgrvtrial,(nfft,nspden))
- + + - +
- + - - +
+ - ]
2318 [ + + ]: 41161 : do ic=1,ncomp
2319 : 22412 : call transgrid(1,mpi_enreg,nspden,-1,0,0,paral_kgb,pawfgr, rhodum,rhodum,cgrvtrial,vtrial(:,:,ic))
2320 : 41161 : call fftpac(isppol,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft, cgrvtrial,vlocal(:,:,:,:,ic),2)
2321 : : end do
2322 [ - + ]: 18749 : ABI_FREE(cgrvtrial)
2323 : : end if
2324 : : else
2325 : : ! nspden == 4. replace isppol by loop over ispden.
2326 [ - + + - : 12650 : ABI_MALLOC(vlocal_tmp, (n4,n5,n6))
- + + - -
+ + - + -
+ - - + -
+ ]
2327 [ + + - + ]: 2530 : if (usepaw == 0 .or. pawfgr%usefinegrid == 0) then
2328 : : ! Fine mesh == Coarse mesh.
2329 [ + + ]: 2896 : do ic=1,ncomp
2330 [ + + ]: 8688 : do ispden=1,nspden
2331 : 5792 : call fftpac(ispden,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,vtrial(:,:,ic),vlocal_tmp,2)
2332 [ + + + + : 76604432 : vlocal(:,:,:,ispden,ic) = vlocal_tmp(:,:,:)
+ + ]
2333 : : end do
2334 : : end do
2335 : : else
2336 : : ! Transfer from fine mesh to coarse and then pack data
2337 [ - + + - : 4328 : ABI_MALLOC(cgrvtrial,(nfft,nspden))
- + + - +
- + - - +
+ - ]
2338 [ + + ]: 2164 : do ic=1,ncomp
2339 : 1082 : call transgrid(1,mpi_enreg,nspden,-1,0,0,paral_kgb,pawfgr,rhodum,rhodum,cgrvtrial,vtrial(:,:,ic))
2340 [ + + ]: 6492 : do ispden=1,nspden
2341 : 4328 : call fftpac(ispden,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,cgrvtrial,vlocal_tmp,2)
2342 [ + + + + : 21768666 : vlocal(:,:,:,ispden,ic) = vlocal_tmp(:,:,:)
+ + ]
2343 : : end do
2344 : : end do
2345 [ - + ]: 1082 : ABI_FREE(cgrvtrial)
2346 : : end if
2347 [ - + ]: 2530 : ABI_FREE(vlocal_tmp)
2348 : : end if ! nspden
2349 : :
2350 : 56053 : end subroutine gspot_transgrid_and_pack
2351 : : !!***
2352 : :
2353 [ - + - - : 24 : end module m_hamiltonian
+ + - + +
- - + - +
+ - - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ ]
2354 : : !!***
|
