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 : 79388 : subroutine gsham_free(Ham)
628 : :
629 : : !Arguments ------------------------------------
630 : : class(gs_hamiltonian_type),target,intent(inout) :: Ham
631 : :
632 : : !Local variables-------------------------------
633 : 79388 : real(dp), contiguous, pointer :: fofr_work_ptr(:,:,:,:)
634 : : ! *************************************************************************
635 : :
636 : : DBG_ENTER("COLL")
637 : :
638 : : !@gs_hamiltonian_type
639 : :
640 : : ! Integer Pointers
641 : 79388 : if (associated(Ham%gbound_kp,Ham%gbound_k)) then
642 : 59136 : nullify(Ham%gbound_kp)
643 [ + + ]: 20252 : else if (associated(Ham%gbound_kp)) then
644 : 19855 : ABI_FREE(Ham%gbound_kp)
645 : : end if
646 : :
647 : : ! Integer arrays
648 [ + - ]: 79388 : 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 [ - + ]: 79388 : ABI_FREE(Ham%atindx)
658 [ - + ]: 79388 : ABI_FREE(Ham%atindx1)
659 [ - + ]: 79388 : ABI_FREE(Ham%typat)
660 [ - + ]: 79388 : ABI_FREE(Ham%indlmn)
661 [ - + ]: 79388 : ABI_FREE(Ham%nattyp)
662 : : end if
663 [ + + ]: 79388 : ABI_SFREE(Ham%gbound_k)
664 [ + - ]: 79388 : ABI_SFREE(Ham%pspso)
665 [ + - ]: 79388 : ABI_SFREE(Ham%dimcprj)
666 : :
667 : : ! Real Pointers
668 [ + + ]: 79388 : if (associated(Ham%phkpxred,Ham%phkxred)) then
669 : 55941 : nullify(Ham%phkpxred)
670 [ + + ]: 23447 : else if (associated(Ham%phkpxred)) then
671 : 23050 : ABI_FREE(Ham%phkpxred)
672 : : end if
673 [ + + ]: 79388 : ABI_SFREE(Ham%phkxred)
674 [ + + ]: 79388 : if (associated(Ham%ekb)) nullify(Ham%ekb)
675 [ + + ]: 79388 : if (associated(Ham%vectornd)) nullify(Ham%vectornd)
676 [ + + ]: 79388 : if (associated(Ham%vlocal)) nullify(Ham%vlocal)
677 [ + + ]: 79388 : if (associated(Ham%vxctaulocal)) nullify(Ham%vxctaulocal)
678 [ + - ]: 79388 : if (associated(Ham%xred)) nullify(Ham%xred)
679 [ + + ]: 79388 : if (associated(Ham%kinpw_k)) nullify(Ham%kinpw_k)
680 [ + + ]: 79388 : if (associated(Ham%kinpw_kp)) nullify(Ham%kinpw_kp)
681 [ + + ]: 79388 : if (associated(Ham%kg_k)) nullify(Ham%kg_k)
682 [ + + ]: 79388 : if (associated(Ham%kg_kp)) nullify(Ham%kg_kp)
683 [ + + ]: 79388 : if (associated(Ham%kpg_k)) nullify(Ham%kpg_k)
684 [ + + ]: 79388 : if (associated(Ham%kpg_kp)) nullify(Ham%kpg_kp)
685 [ + + ]: 79388 : if (associated(Ham%ffnl_k)) nullify(Ham%ffnl_k)
686 [ + + ]: 79388 : if (associated(Ham%ffnl_kp)) nullify(Ham%ffnl_kp)
687 [ + + ]: 79388 : if (associated(Ham%ph3d_k)) nullify(Ham%ph3d_k)
688 [ + + ]: 79388 : if (associated(Ham%ph3d_kp)) nullify(Ham%ph3d_kp)
689 : :
690 : :
691 : : ! Real arrays
692 [ + - ]: 79388 : ABI_SFREE(Ham%ekb_spin)
693 [ + - ]: 79388 : ABI_SFREE(Ham%sij)
694 [ + - ]: 79388 : ABI_SFREE(Ham%nucdipmom)
695 [ + - ]: 79388 : 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 [ - + ]: 79388 : ABI_FREE(Ham%ph1d)
701 : : end if
702 : :
703 [ - + ]: 79388 : if(Ham%gpu_option==ABI_GPU_OPENMP) then
704 : 0 : call free_ompgpu_fourwf()
705 : : end if
706 : :
707 : : ! Structured datatype pointers
708 [ + + ]: 79388 : if (associated(Ham%fockcommon)) nullify(Ham%fockcommon)
709 [ + + ]: 79388 : if (associated(Ham%fockACE_k)) nullify(Ham%fockACE_k)
710 [ + + ]: 79388 : 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 [ - + ]: 79388 : ABI_SFREE(Ham%fofr_work)
724 : :
725 : : DBG_EXIT("COLL")
726 : :
727 [ + + ]: 79388 : 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 : 79380 : subroutine gsham_init(ham,Psps,pawtab,nspinor,nsppol,nspden,natom,typat,&
778 : 79380 : xred,nfft,mgfft,ngfft,rprimd,nloalg,&
779 [ + + + + : 79380 : ph1d,usecprj,comm_atom,mpi_atmtab,mpi_spintab,paw_ij,& ! optional
- + + + +
+ - + ]
780 : 75367 : 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),ispin_gbt
807 : : real(dp) :: gmet(3,3),gprimd(3,3),rmet(3,3)
808 : 79380 : real(dp),allocatable,target :: ekb_tmp(:,:,:,:)
809 : : ! *************************************************************************
810 : :
811 : : DBG_ENTER("COLL")
812 : :
813 : : !@gs_hamiltonian_type
814 : :
815 : : !Manage optional parameters
816 [ + + ]: 79380 : my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
817 [ + + + + ]: 166062 : my_spintab=0;my_spintab(1:nsppol)=1;if (present(mpi_spintab)) my_spintab(1:2)=mpi_spintab(1:2)
818 [ + + + + ]: 238140 : my_nsppol=count(my_spintab==1)
819 [ + + ]: 79380 : l_gpu_option=ABI_GPU_DISABLED; if(present(gpu_option)) l_gpu_option=gpu_option
820 [ - + ]: 79380 : my_zora=0; if (present(zora)) my_zora=zora
821 : :
822 : 79380 : ham%nfft_blocks=1;
823 [ - + ]: 79380 : if (present(nfft_blocks) .and. l_gpu_option==ABI_GPU_OPENMP) then
824 : 0 : ham%nfft_blocks=nfft_blocks
825 : : end if
826 [ + + ]: 79380 : ham%use_gbt = 0; if (present(use_gbt)) ham%use_gbt = use_gbt
827 [ + + + + : 459534 : ham%spinaxis = zero; if (present(spinaxis)) ham%spinaxis = spinaxis
+ + ]
828 : :
829 : 79380 : call metric(gmet,gprimd,-1,rmet,rprimd,ucvol)
830 : :
831 [ + + - + : 317520 : 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 [ + - ]: 79380 : 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 [ - + + - : 238140 : ABI_MALLOC(ham%atindx,(natom))
+ - + - -
+ - + -
+ ]
845 [ - + + - : 158760 : ABI_MALLOC(ham%atindx1,(natom))
- + - + -
+ ]
846 [ - + + - : 158760 : ABI_MALLOC(ham%typat,(natom))
- + - + -
+ ]
847 [ - + + - : 317520 : ABI_MALLOC(ham%indlmn,(6,psps%lmnmax,psps%ntypat))
- + + - +
- + - - +
- + - + ]
848 [ - + + - : 238140 : ABI_MALLOC(ham%nattyp,(psps%ntypat))
+ - + - -
+ - + -
+ ]
849 [ + - - + : 238140 : ABI_MALLOC(ham%ph1d,(2,3*(2*mgfft+1)*natom))
+ - + - +
- - + - +
- + ]
850 : : end if
851 : :
852 [ - + + - : 238140 : ABI_MALLOC(ham%pspso,(psps%ntypat))
+ - + - -
+ - + -
+ ]
853 [ - + + - : 238140 : ABI_MALLOC(ham%nucdipmom,(3,natom))
+ - + - -
+ - + -
+ ]
854 : :
855 : : !Initialize most of the Hamiltonian
856 : 79380 : indx=1
857 [ + + ]: 194315 : do itypat=1,psps%ntypat
858 : 114935 : ham%nattyp(itypat)=0
859 [ + + ]: 499686 : do iat=1,natom
860 [ + + ]: 420306 : if (typat(iat)==itypat) then
861 : 204308 : ham%atindx (iat )=indx
862 : 204308 : ham%atindx1(indx)=iat
863 : 204308 : indx=indx+1
864 : 204308 : ham%nattyp(itypat)=ham%nattyp(itypat)+1
865 : : end if
866 : : end do
867 : : end do
868 : :
869 [ + + + + ]: 1031940 : ham%gmet(:,:) =gmet(:,:)
870 [ + + + + ]: 1031940 : ham%gprimd(:,:)=gprimd(:,:)
871 [ + + + + : 4669793 : ham%indlmn(:,:,:)=psps%indlmn(:,:,:)
+ + ]
872 : 79380 : ham%ispin_gbt =ispin_gbt
873 : 79380 : ham%lmnmax =psps%lmnmax
874 : 79380 : ham%mgfft =mgfft
875 : 79380 : ham%mpsang =psps%mpsang
876 : 79380 : ham%mpssoang =psps%mpssoang
877 : 79380 : ham%natom =natom
878 : 79380 : ham%nfft =nfft
879 [ + + ]: 1508220 : ham%ngfft(:) =ngfft(:)
880 [ + + ]: 317520 : ham%nloalg(:) =nloalg(:)
881 [ + + + - ]: 194315 : ham%matblk=min(NLO_MINCAT,maxval(ham%nattyp)); if (nloalg(2)>0) ham%matblk=natom
882 : 79380 : ham%nsppol =nsppol
883 : 79380 : ham%nspinor =nspinor
884 : 79380 : ham%ntypat =psps%ntypat
885 [ + - - + : 363068 : ham%typat =typat(1:natom)
- - + + ]
886 [ + + ]: 79380 : ham%nvloc=1; if(nspden==4)ham%nvloc=4
887 : 79380 : ham%n4 =ngfft(4)
888 : 79380 : ham%n5 =ngfft(5)
889 : 79380 : ham%n6 =ngfft(6)
890 : 79380 : ham%usepaw =psps%usepaw
891 : 79380 : ham%ucvol =ucvol
892 : 79380 : ham%useylm =psps%useylm
893 : 79380 : ham%zora =my_zora
894 [ + + ]: 79380 : ham%gpu_option=ABI_GPU_DISABLED ; if(PRESENT(gpu_option)) ham%gpu_option=gpu_option
895 : :
896 [ + + ]: 194315 : ham%pspso(:) =psps%pspso(1:psps%ntypat)
897 [ + + ]: 79380 : if (psps%usepaw==1) then
898 [ + + ]: 57511 : do itypat=1,psps%ntypat
899 : 57511 : ham%pspso(itypat)=1+pawtab(itypat)%usespnorb
900 : : end do
901 : : end if
902 : :
903 [ + + ]: 79380 : if (present(nucdipmom)) then
904 [ + + + + ]: 858163 : ham%nucdipmom(:,:) = nucdipmom(:,:)
905 : : else
906 [ + + + + ]: 38449 : ham%nucdipmom(:,:) = zero
907 : : end if
908 : :
909 : 79380 : ham%xred => xred
910 : :
911 [ + + ]: 79380 : if (present(fock)) then
912 [ + + ]: 54883 : 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 [ + + ]: 79380 : if (present(ph1d)) then
919 [ + + + + ]: 90925781 : ham%ph1d(:,:) = ph1d(:,:)
920 : : else ! Recalculate structure factor phases
921 : 256 : call getph(ham%atindx,natom,ngfft(1),ngfft(2),ngfft(3),ham%ph1d,xred)
922 : : end if
923 : :
924 [ + + ]: 79380 : if (ham%usepaw==1) then
925 [ + + ]: 22780 : ham%usecprj=0;if (present(usecprj)) ham%usecprj=usecprj
926 [ - + + - : 68340 : ABI_MALLOC(ham%dimcprj,(natom))
+ - + - -
+ - + -
+ ]
927 : : !Be careful cprj are ordered by atom type (used in non-local operator)
928 : 22780 : call pawcprj_getdim(ham%dimcprj,natom,ham%nattyp,ham%ntypat,ham%typat,pawtab,'O')
929 : : else
930 : 56600 : ham%usecprj=0
931 [ - + - + ]: 56600 : ABI_MALLOC(ham%dimcprj,(0))
932 : : end if
933 : :
934 : : ! ===========================
935 : : ! ==== Non-local factors ====
936 : : ! ===========================
937 : :
938 [ + + ]: 79380 : 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 : 56600 : ham%dimekb1=psps%dimekb
945 : 56600 : ham%dimekb2=psps%ntypat
946 : 56600 : ham%dimekbq=1
947 [ - + + - : 283000 : ABI_MALLOC(ham%ekb_spin,(psps%dimekb,psps%ntypat,nspinor**2,1,1))
- + + - -
+ + - + -
+ - - + -
+ - + ]
948 : 56600 : ham%ekb => ham%ekb_spin(:,:,:,:,1)
949 [ - + - + ]: 56600 : ABI_MALLOC(ham%sij,(0,0))
950 [ + + + + ]: 381396 : ham%ekb(:,:,1,1)=psps%ekb(:,:)
951 [ + + ]: 56600 : if (nspinor==2) then
952 [ + + + + ]: 36208 : ham%ekb(:,:,2,1)=psps%ekb(:,:)
953 [ + + + + : 75085 : ham%ekb(:,:,3:4,1)=zero
+ + ]
954 : : end if
955 [ + + ]: 56600 : if (PRESENT(electronpositron)) then
956 [ + + + + : 36796 : 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 : 22780 : cplex_dij=1
972 [ + - ]: 22780 : if (present(paw_ij)) then
973 [ + + ]: 22780 : if (size(paw_ij)>0) cplex_dij=paw_ij(1)%cplex_dij
974 : : end if
975 [ + + + + : 270438 : if ((nspinor==2).or.any(abs(ham%nucdipmom)>tol8)) cplex_dij=2
+ + + + +
+ ]
976 : 22780 : ham%dimekb1=psps%dimekb*cplex_dij
977 : 22780 : ham%dimekb2=natom
978 : 22780 : ham%dimekbq=1
979 [ + - ]: 22780 : if (present(paw_ij)) then
980 [ + + ]: 22780 : if (size(paw_ij)>0) ham%dimekbq=paw_ij(1)%qphase
981 : : end if
982 [ - + + - : 91120 : ABI_MALLOC(ham%sij,(ham%dimekb1,psps%ntypat))
- + + - +
- + - - +
- + - + ]
983 [ + + ]: 57511 : do itypat=1,psps%ntypat
984 [ + + ]: 34731 : if (cplex_dij==1) then
985 [ + + ]: 1685547 : 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 [ + + ]: 57511 : if (cplex_dij*pawtab(itypat)%lmn2_size<ham%dimekb1) then
993 [ + + ]: 221815 : 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 [ - + + - : 159460 : ABI_MALLOC(ham%ekb_spin,(ham%dimekb1,ham%dimekb2,nspinor**2,ham%dimekbq,my_nsppol))
- + + - -
+ + - - +
+ - + + +
- + - + +
- + - + -
+ ]
999 [ + + + + : 8153217 : ham%ekb_spin=zero
+ + + + +
+ ]
1000 [ - + ]: 22780 : if (present(paw_ij)) then
1001 [ + + ]: 22780 : if (my_nsppol<ham%nsppol) then
1002 [ - + + - : 3168 : ABI_MALLOC(ekb_tmp,(ham%dimekb1,ham%dimekb2,nspinor**2,ham%dimekbq))
- + + - -
+ + - - +
+ - + - +
- - + -
+ ]
1003 : : end if
1004 : 22780 : jsp=0
1005 [ + + ]: 48262 : do isp=1,ham%nsppol
1006 [ + + ]: 25482 : if (my_spintab(isp)==1) then
1007 : 24926 : jsp=jsp+1 ; ham%ekb => ham%ekb_spin(:,:,:,:,jsp)
1008 : : else
1009 : 556 : ham%ekb => ekb_tmp
1010 : : end if
1011 [ + + ]: 48262 : if (present(mpi_atmtab)) then
1012 : 7556 : call pawdij2ekb(ham%ekb,paw_ij,isp,my_comm_atom,mpi_atmtab=mpi_atmtab)
1013 : : else
1014 : 17926 : call pawdij2ekb(ham%ekb,paw_ij,isp,my_comm_atom)
1015 : : end if
1016 : : end do
1017 [ + + ]: 22780 : if (my_nsppol<ham%nsppol) then
1018 [ - + ]: 528 : ABI_FREE(ekb_tmp)
1019 : : end if
1020 : : end if
1021 : 22780 : nullify(ham%ekb)
1022 : : end if
1023 : :
1024 : : DBG_EXIT("COLL")
1025 : :
1026 [ + + + + : 234127 : 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 [ + + + - : 9542514 : subroutine gsham_load_k(ham,ffnl_k,fockACE_k,gbound_k,istwf_k,kinpw_k,&
- + + + +
- - + + +
+ - - + ]
1061 [ + + + - : 6361676 : 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 [ + + ]: 3180838 : if (present(kpt_k)) then
1089 [ + + ]: 7478684 : ham%kpt_k(:) = kpt_k(:)
1090 [ + + ]: 7478684 : ham%kpt_kp(:) = kpt_k(:)
1091 : : end if
1092 [ + + ]: 3180838 : if (present(istwf_k)) then
1093 : 1869671 : ham%istwf_k = istwf_k
1094 : 1869671 : ham%istwf_kp = istwf_k
1095 : : end if
1096 [ + + ]: 3180838 : if (present(npw_k)) then
1097 : 1869671 : ham%npw_k = npw_k
1098 : 1869671 : ham%npw_kp = npw_k
1099 : : end if
1100 [ + + ]: 3180838 : 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 [ + + ]: 3134588 : else if (present(npw_k)) then
1104 : 1869671 : ham%npw_fft_k = npw_k
1105 : 1869671 : ham%npw_fft_kp = npw_k
1106 : : end if
1107 : :
1108 : : !Pointers to k-dependent quantities
1109 [ + + ]: 3180838 : if (present(kinpw_k)) then
1110 : 506514 : ham%kinpw_k => kinpw_k
1111 : 506514 : ham%kinpw_kp => kinpw_k
1112 : : end if
1113 [ + + ]: 3180838 : if (present(kg_k)) then
1114 : 1915929 : ham%kg_k => kg_k
1115 : 1915929 : ham%kg_kp => kg_k
1116 : : end if
1117 [ + + ]: 3180838 : if (present(kpg_k)) then
1118 : 1905599 : ham%kpg_k => kpg_k
1119 : 1905599 : ham%kpg_kp => kpg_k
1120 : : end if
1121 [ + + ]: 3180838 : if (present(ffnl_k)) then
1122 : 1974797 : ham%ffnl_k => ffnl_k
1123 : 1974797 : ham%ffnl_kp => ffnl_k
1124 : : end if
1125 [ + + ]: 3180838 : if (present(ph3d_k)) then
1126 : 1795988 : ham%ph3d_k => ph3d_k
1127 : 1795988 : ham%ph3d_kp => ph3d_k
1128 : : end if
1129 [ + + ]: 3180838 : 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 [ + + ]: 3180838 : if (present(kpt_k)) then
1134 [ + + + + ]: 1869671 : if (associated(Ham%phkpxred).and.(.not.associated(Ham%phkpxred,Ham%phkxred))) then
1135 : 1355976 : ABI_FREE(Ham%phkpxred)
1136 : : end if
1137 [ + + ]: 1869671 : ABI_SFREE(ham%phkxred)
1138 [ - + + - : 5609013 : ABI_MALLOC(ham%phkxred,(2,ham%natom))
+ - + - -
+ - + -
+ ]
1139 [ + + ]: 5479202 : do iat=1,ham%natom
1140 : 3609531 : iatom=ham%atindx(iat)
1141 [ + + ]: 14438124 : arg=two_pi*DOT_PRODUCT(kpt_k,ham%xred(:,iat))
1142 : 3609531 : ham%phkxred(1,iatom)=DCOS(arg)
1143 : 5479202 : ham%phkxred(2,iatom)=DSIN(arg)
1144 : : end do
1145 : 1869671 : ham%phkpxred => ham%phkxred
1146 : : end if
1147 : :
1148 : : ! Compute or copy G sphere boundary at k+g
1149 [ + + ]: 3180838 : compute_gbound_=.false.;if (present(compute_gbound)) compute_gbound_=compute_gbound
1150 [ + + ]: 3180838 : if (present(gbound_k)) compute_gbound_=.true.
1151 [ + + ]: 3134588 : if (compute_gbound_) then
1152 [ + + ]: 1694183 : if (associated(Ham%gbound_kp,Ham%gbound_k)) then
1153 : 377832 : nullify(Ham%gbound_kp)
1154 [ + + ]: 1316351 : else if (associated(Ham%gbound_kp)) then
1155 : 1253670 : ABI_FREE(Ham%gbound_kp)
1156 : : end if
1157 [ + + ]: 1694183 : ABI_SFREE(ham%gbound_k)
1158 : : end if
1159 [ + + ]: 3180838 : if (.not.allocated(ham%gbound_k)) then
1160 [ + - - + : 6841940 : ABI_MALLOC(ham%gbound_k,(2*ham%mgfft+8,2))
+ - + - +
- + - - +
- + - + ]
1161 [ + + + + ]: 139264843 : ham%gbound_k(:,:)=0
1162 : 1710485 : ham%gbound_kp => ham%gbound_k
1163 : : end if
1164 [ + + ]: 3180838 : if (compute_gbound_) then
1165 [ + + ]: 1694183 : if (present(gbound_k)) then
1166 [ + + + + ]: 4388306 : ham%gbound_k(:,:)=gbound_k(:,:)
1167 : : else
1168 [ - + ]: 1647933 : 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 [ + - - - : 1647933 : call sphereboundary(ham%gbound_k,ham%istwf_k,ham%kg_k,ham%mgfft,ham%npw_k)
- - - - -
+ ]
1174 : : end if
1175 : 1694183 : ham%gbound_kp => ham%gbound_k
1176 : : end if
1177 : :
1178 : : ! Compute 3D structure factors for each atom at k+g
1179 [ + + ]: 3180838 : if (present(compute_ph3d).and.present(ph3d_k)) then
1180 [ + + + - ]: 1749051 : if (compute_ph3d.and.ham%nloalg(2)>0) then
1181 [ + - + - : 1716091 : 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 [ + - - - : 1716091 : 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 : 3180838 : 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 [ + + + - : 6117354 : subroutine gsham_load_kprime(ham,ffnl_kp,gbound_kp,istwf_kp,kinpw_kp,&
- + + + +
- - + - +
- - - - ]
1228 [ + + + - : 2039118 : kg_kp,kpg_kp,kpt_kp,npw_kp,npw_fft_kp,&
- + + + +
+ - + ]
1229 [ + + + - : 2039118 : 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 [ + + + + ]: 6863745 : if (present(kpt_kp)) ham%kpt_kp(:)= kpt_kp(:)
1255 [ + + ]: 2039118 : if (present(istwf_kp)) ham%istwf_kp = istwf_kp
1256 [ + + ]: 2039118 : if (present(npw_kp)) ham%npw_kp = npw_kp
1257 [ - + ]: 2039118 : if (present(npw_fft_kp)) then
1258 : 0 : ham%npw_fft_kp = npw_fft_kp
1259 [ + + ]: 2039118 : else if (present(npw_kp)) then
1260 : 1607369 : ham%npw_fft_kp = npw_kp
1261 : : end if
1262 : :
1263 : : !Pointers to k-dependent quantities
1264 [ + + ]: 2039118 : if (present(kinpw_kp)) ham%kinpw_kp => kinpw_kp
1265 [ + + ]: 2039118 : if (present(kg_kp)) ham%kg_kp => kg_kp
1266 [ + + ]: 2039118 : if (present(kpg_kp)) ham%kpg_kp => kpg_kp
1267 [ + + ]: 2039118 : if (present(ffnl_kp)) ham%ffnl_kp => ffnl_kp
1268 [ + + ]: 2039118 : if (present(ph3d_kp)) ham%ph3d_kp => ph3d_kp
1269 : :
1270 : : !Compute exp(i.k^prime.R) for each atom
1271 [ + + ]: 2039118 : if (present(kpt_kp)) then
1272 [ + + ]: 1608209 : if (associated(ham%phkpxred,ham%phkxred)) then
1273 : 1379026 : nullify(ham%phkpxred)
1274 [ + - ]: 229183 : else if (associated(ham%phkpxred)) then
1275 : 229183 : ABI_FREE(ham%phkpxred)
1276 : : end if
1277 [ - + + - : 4824627 : ABI_MALLOC(ham%phkpxred,(2,ham%natom))
+ - + - -
+ - + ]
1278 [ + + ]: 4740925 : do iat=1,ham%natom
1279 : 3132716 : iatom=ham%atindx(iat)
1280 [ + + ]: 12530864 : arg=two_pi*DOT_PRODUCT(kpt_kp,ham%xred(:,iat))
1281 : 3132716 : ham%phkpxred(1,iatom)=DCOS(arg)
1282 : 4740925 : 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 : 2039118 : compute_gbound_=.false.
1288 [ + + ]: 2039118 : if (present(kpt_kp).and.present(compute_gbound)) compute_gbound_=compute_gbound
1289 [ + - ]: 2039118 : if (present(gbound_kp)) compute_gbound_=.true.
1290 [ + + ]: 2039118 : if (compute_gbound_) then
1291 [ + + ]: 1277293 : if (associated(ham%gbound_kp,ham%gbound_k)) then
1292 : 1273525 : nullify(ham%gbound_kp)
1293 [ + - ]: 3768 : else if (associated(ham%gbound_kp)) then
1294 : 3768 : ABI_FREE(ham%gbound_kp)
1295 : : end if
1296 [ - + ]: 1277293 : if (present(gbound_kp)) then
1297 [ # # # # ]: 0 : ham%gbound_kp(:,:)=gbound_kp(:,:)
1298 : : else
1299 [ - + ]: 1277293 : if (.not.associated(ham%kg_kp)) then
1300 : 0 : ABI_BUG('Something is missing for gbound_kp computation!')
1301 : : end if
1302 [ + - - + : 5109172 : ABI_MALLOC(ham%gbound_kp,(2*ham%mgfft+8,2))
+ - + - +
- + - - +
- + ]
1303 [ - + - - : 1277293 : call sphereboundary(ham%gbound_kp,ham%istwf_kp,ham%kg_kp,ham%mgfft,ham%npw_kp)
+ - - - -
- - - - +
- - - - -
+ ]
1304 : : end if
1305 : : end if
1306 : :
1307 : : !Compute 3D structure factors for each atom at k^prime+g
1308 [ + + ]: 2039118 : if (present(compute_ph3d).and.present(ph3d_kp)) then
1309 [ + + + - ]: 603246 : if (compute_ph3d.and.ham%nloalg(2)>0) then
1310 [ + - + - : 585876 : if ((.not.associated(ham%phkpxred)).or.(.not.associated(ham%kg_kp)).or.&
- + ]
1311 : : & (.not.associated(ham%ph3d_kp))) then
1312 : 0 : ABI_BUG('Something is missing for ph3d_kp computation!')
1313 : : end if
1314 : : call ph1d3d(1,ham%natom,ham%kg_kp,ham%matblk,ham%natom,ham%npw_kp,ham%ngfft(1),&
1315 [ + - - - : 585876 : & ham%ngfft(2),ham%ngfft(3),ham%phkpxred,ham%ph1d,ham%ph3d_kp)
- - - - +
- - - - -
- - + - -
- - + - +
- + - - -
- - - ]
1316 : : end if
1317 : : end if
1318 : :
1319 : : DBG_EXIT("COLL")
1320 : :
1321 : 2039118 : end subroutine gsham_load_kprime
1322 : : !!***
1323 : :
1324 : : !----------------------------------------------------------------------
1325 : :
1326 : : !!****f* m_hamiltonian/gsham_eph_setup_k
1327 : : !! NAME
1328 : : !! gsham_eph_setup_k
1329 : : !!
1330 : : !! FUNCTION
1331 : : !! Simplified interface to load either k or kprime in the case of e-ph calculations.
1332 : : !!
1333 : : !! INPUTS
1334 : : !! which_k= "k" to load k, "kq" to load kprime
1335 : : !! See load_k or load_kprime for the meaning of arguments.
1336 : : !!
1337 : : !! SOURCE
1338 : :
1339 : 157025 : subroutine gsham_eph_setup_k(gs_ham, which_k, kk, istwf_k, npw_k, kg_k, dtset, cryst, psps, & ! in
1340 : : nkpg_k, kpg_k, ffnl_k, kinpw_k, ph3d_k, comm) ! out
1341 : :
1342 : : !Arguments ------------------------------------
1343 : : !scalars
1344 : : class(gs_hamiltonian_type),intent(inout) :: gs_ham
1345 : : character(len=*),intent(in) :: which_k
1346 : : type(dataset_type),intent(in) :: dtset
1347 : : type(crystal_t),intent(in) :: cryst
1348 : : type(pseudopotential_type),intent(in) :: psps
1349 : : integer,intent(in) :: istwf_k, npw_k, comm
1350 : : integer,intent(out) :: nkpg_k
1351 : : !arrays
1352 : : real(dp),intent(in) :: kk(3)
1353 : : integer,intent(in) :: kg_k(3,npw_k)
1354 : : real(dp),allocatable,intent(out) :: kpg_k(:,:), ffnl_k(:,:,:,:), kinpw_k(:), ph3d_k(:,:,:)
1355 : :
1356 : : !Local variables-------------------------------
1357 : : !scalars
1358 : : integer,parameter :: ider0 = 0, idir0 = 0, optder0 = 0
1359 : : !arrays
1360 : : real(dp) :: ylmgr_k_dum(1,1,1)
1361 : : real(dp),allocatable :: ylm_k(:,:)
1362 : : ! *************************************************************************
1363 : :
1364 : : ! Compute k+G vectors
1365 : 157025 : nkpg_k = 3 * dtset%nloalg(3)
1366 [ - + + - : 628100 : ABI_MALLOC(kpg_k, (npw_k, nkpg_k))
+ - - - +
- - + - +
- + - + ]
1367 [ - + ]: 157025 : if (nkpg_k > 0) call mkkpg(kg_k, kpg_k, kk, nkpg_k, npw_k)
1368 : :
1369 : : ! Spherical Harmonics at k for useylm == 1.
1370 [ - + + - : 628100 : ABI_MALLOC(ylm_k, (npw_k, psps%mpsang**2 * psps%useylm))
+ + + - +
- + + - +
- + ]
1371 [ + + ]: 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)
1372 : :
1373 : : ! Compute nonlocal form factors ffnl_k at (k+G)
1374 [ - + + - : 785125 : ABI_MALLOC(ffnl_k, (npw_k, 1, psps%lmnmax, psps%ntypat))
- + + - -
+ + - + -
+ - - + -
+ - + ]
1375 : 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)
1376 [ - + ]: 157025 : ABI_FREE(ylm_k)
1377 : :
1378 : : ! Compute (1/2) (2 Pi)**2 (kG)**2:
1379 [ - + + - : 134657935 : ABI_CALLOC(kinpw_k, (npw_k))
+ - + - -
+ - + - +
+ + ]
1380 : 157025 : call mkkin(dtset%ecut, dtset%ecutsm, dtset%effmass_free, cryst%gmet, kg_k, kinpw_k, kk, npw_k, 0, 0)
1381 : :
1382 [ - + + - : 628100 : ABI_MALLOC(ph3d_k, (2, npw_k, gs_ham%matblk))
- + + - +
- + - - +
- + - + ]
1383 : :
1384 : : ! Load the k dependent parts of the Hamiltonian
1385 : 76649 : select case (which_k)
1386 : : case ("k")
1387 : : 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, &
1388 : 76649 : ph3d_k=ph3d_k, ffnl_k=ffnl_k, compute_ph3d=.true., compute_gbound=.true.)
1389 : : case ("kq")
1390 : : 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, &
1391 : 80376 : ph3d_kp=ph3d_k, ffnl_kp=ffnl_k, compute_ph3d=.true., compute_gbound=.true.)
1392 : : case default
1393 [ + + - ]: 157025 : ABI_ERROR(sjoin("Invalid value for which_k:", which_k))
1394 : : end select
1395 : :
1396 : 157025 : end subroutine gsham_eph_setup_k
1397 : : !!***
1398 : :
1399 : : !----------------------------------------------------------------------
1400 : :
1401 : : !!****f* m_hamiltonian/gsham_copy
1402 : : !! NAME
1403 : : !! gsham_copy
1404 : : !!
1405 : : !! INPUTS
1406 : : !! gs_hamk_in<gs_hamiltonian_type>=Structured datatype completely initialized, to be copied.
1407 : : !!
1408 : : !! FUNCTION
1409 : : !! Copy a gs_hamiltonian_type variable (gs_hamk_in) in another (gs_hamk_out).
1410 : : !! In contrast to an assignment statement (gs_hamk_out=gs_hamk_in), this
1411 : : !! subroutine allocate memory space for the pointers contained in the data
1412 : : !! structure (gs_hamk_out) and copy the content of the corresponding memory
1413 : : !! space of gs_hamk_in in it. In contrast, the assignment statement would
1414 : : !! only associate the pointers of gs_hamk_out to the same memory space than
1415 : : !! the corresponding ones in gs_hamk_in. This can cause trouble if one data
1416 : : !! structure is destroyed before a reading/writing statement for the other
1417 : : !! structure, causing access to unallocated memory space (silently, without
1418 : : !! segmentation fault being generated).
1419 : : !!
1420 : : !! OUTPUT
1421 : : !! gs_hamk_out<gs_hamiltonian_type>=Structured datatype containing separate
1422 : : !! copies of all data of gs_hamk_in upon exit.
1423 : : !!
1424 : : !! SOURCE
1425 : :
1426 [ + - ]: 8 : subroutine gsham_copy(gs_hamk_in, gs_hamk_out)
1427 : :
1428 : : !Arguments ------------------------------------
1429 : : class(gs_hamiltonian_type),intent(in),target :: gs_hamk_in
1430 : : class(gs_hamiltonian_type),intent(out),target :: gs_hamk_out
1431 : :
1432 : : !Local variables-------------------------------
1433 : : integer :: tmp2i(5)
1434 : : type(C_PTR) :: ham_ptr
1435 : : ! *************************************************************************
1436 : :
1437 : : DBG_ENTER("COLL")
1438 : :
1439 : : !@gs_hamiltonian_type
1440 : :
1441 : 8 : gs_hamk_out%dimekb1 = gs_hamk_in%dimekb1
1442 : 8 : gs_hamk_out%dimekb2 = gs_hamk_in%dimekb2
1443 : 8 : gs_hamk_out%dimekbq = gs_hamk_in%dimekbq
1444 : 8 : gs_hamk_out%istwf_k = gs_hamk_in%istwf_k
1445 : 8 : gs_hamk_out%istwf_kp = gs_hamk_in%istwf_kp
1446 : 8 : gs_hamk_out%lmnmax = gs_hamk_in%lmnmax
1447 : 8 : gs_hamk_out%matblk = gs_hamk_in%matblk
1448 : 8 : gs_hamk_out%mgfft = gs_hamk_in%mgfft
1449 : 8 : gs_hamk_out%mpsang = gs_hamk_in%mpsang
1450 : 8 : gs_hamk_out%mpssoang = gs_hamk_in%mpssoang
1451 : 8 : gs_hamk_out%natom = gs_hamk_in%natom
1452 : 8 : gs_hamk_out%nfft = gs_hamk_in%nfft
1453 : 8 : gs_hamk_out%nfft_blocks = gs_hamk_in%nfft_blocks
1454 : 8 : gs_hamk_out%npw_k = gs_hamk_in%npw_k
1455 : 8 : gs_hamk_out%npw_kp = gs_hamk_in%npw_kp
1456 : 8 : gs_hamk_out%npw_fft_k = gs_hamk_in%npw_fft_k
1457 : 8 : gs_hamk_out%npw_fft_kp = gs_hamk_in%npw_fft_kp
1458 : 8 : gs_hamk_out%nspinor = gs_hamk_in%nspinor
1459 : 8 : gs_hamk_out%nsppol = gs_hamk_in%nsppol
1460 : 8 : gs_hamk_out%ntypat = gs_hamk_in%ntypat
1461 : 8 : gs_hamk_out%nvloc = gs_hamk_in%nvloc
1462 : 8 : gs_hamk_out%n4 = gs_hamk_in%n4
1463 : 8 : gs_hamk_out%n5 = gs_hamk_in%n5
1464 : 8 : gs_hamk_out%n6 = gs_hamk_in%n6
1465 : 8 : gs_hamk_out%gpu_option = gs_hamk_in%gpu_option
1466 [ + + + + ]: 56 : gs_hamk_out%spinaxis = gs_hamk_in%spinaxis
1467 : 8 : gs_hamk_out%usecprj = gs_hamk_in%usecprj
1468 : 8 : gs_hamk_out%usepaw = gs_hamk_in%usepaw
1469 : 8 : gs_hamk_out%useylm = gs_hamk_in%useylm
1470 : 8 : gs_hamk_out%use_gbt = gs_hamk_in%use_gbt
1471 : 8 : gs_hamk_out%zora = gs_hamk_in%zora
1472 [ + + + + ]: 296 : gs_hamk_out%ngfft = gs_hamk_in%ngfft
1473 [ + + + + ]: 56 : gs_hamk_out%nloalg = gs_hamk_in%nloalg
1474 : 8 : gs_hamk_out%ucvol = gs_hamk_in%ucvol
1475 [ + + + + : 200 : gs_hamk_out%gmet = gs_hamk_in%gmet
+ + + + ]
1476 [ + + + + : 200 : gs_hamk_out%gprimd = gs_hamk_in%gprimd
+ + + + ]
1477 : 8 : gs_hamk_out%ispin_gbt = gs_hamk_in%ispin_gbt
1478 [ + + + + ]: 56 : gs_hamk_out%kpt_k = gs_hamk_in%kpt_k
1479 [ + + + + ]: 56 : gs_hamk_out%kpt_kp = gs_hamk_in%kpt_kp
1480 : :
1481 [ - + + - : 24 : ABI_MALLOC(gs_hamk_out%atindx,(gs_hamk_out%natom))
+ - + - -
+ - + -
+ ]
1482 [ + - + + : 96 : gs_hamk_out%atindx = gs_hamk_in%atindx
+ - - + -
- - - - -
+ + ]
1483 [ - + + - : 24 : ABI_MALLOC(gs_hamk_out%atindx1,(gs_hamk_out%natom))
+ - + - -
+ - + -
+ ]
1484 [ + - + + : 96 : gs_hamk_out%atindx1 = gs_hamk_in%atindx1
+ - - + -
- - - - -
+ + ]
1485 [ + - - - : 24 : ABI_MALLOC(gs_hamk_out%dimcprj,(gs_hamk_out%natom*gs_hamk_out%usepaw))
+ - - + -
+ - + -
+ ]
1486 [ - + - - : 8 : if (gs_hamk_out%usepaw==1) gs_hamk_out%dimcprj = gs_hamk_in%dimcprj
- - - - -
- - - - -
- - - - ]
1487 [ - + + - : 24 : ABI_MALLOC(gs_hamk_out%typat,(gs_hamk_out%natom))
+ - + - -
+ - + -
+ ]
1488 [ + - + + : 96 : gs_hamk_out%typat = gs_hamk_in%typat
+ - - + -
- - - - -
+ + ]
1489 [ + - - + : 32 : ABI_MALLOC(gs_hamk_out%gbound_k,(2*gs_hamk_out%mgfft+8,2))
+ - + - +
- + - - +
- + - + ]
1490 [ + - + + : 1456 : gs_hamk_out%gbound_k = gs_hamk_in%gbound_k
+ + + - +
- - + - -
- - - - -
- - - + +
+ + ]
1491 [ - + + - : 32 : ABI_MALLOC(gs_hamk_out%indlmn,(6,gs_hamk_out%lmnmax,gs_hamk_out%ntypat))
- + + - +
- + - - +
- + - + ]
1492 [ + - + + : 720 : gs_hamk_out%indlmn = gs_hamk_in%indlmn
+ + + + +
- + - + -
- + - - -
- - - - -
- - - - -
- + + + +
+ + ]
1493 [ - + + - : 24 : ABI_MALLOC(gs_hamk_out%nattyp,(gs_hamk_out%ntypat))
+ - + - -
+ - + -
+ ]
1494 [ + - + + : 48 : gs_hamk_out%nattyp = gs_hamk_in%nattyp
+ - - + -
- - - - -
+ + ]
1495 [ - + + - : 24 : ABI_MALLOC(gs_hamk_out%nucdipmom,(3,gs_hamk_out%natom))
+ - + - -
+ - + -
+ ]
1496 [ + - + + : 336 : gs_hamk_out%nucdipmom = gs_hamk_in%nucdipmom
+ + + - +
- - + - -
- - - - -
- - - + +
+ + ]
1497 [ - + + - : 24 : ABI_MALLOC(gs_hamk_out%phkxred,(2,gs_hamk_out%natom))
+ - + - -
+ - + -
+ ]
1498 [ + - + + : 256 : gs_hamk_out%phkxred = gs_hamk_in%phkxred
+ + + - +
- - + - -
- - - - -
- - - + +
+ + ]
1499 [ + - - + : 24 : ABI_MALLOC(gs_hamk_out%ph1d,(2,3*(2*gs_hamk_out%mgfft+1)*gs_hamk_out%natom))
+ - + - +
- - + - +
- + ]
1500 [ + - + + : 26656 : gs_hamk_out%ph1d = gs_hamk_in%ph1d
+ + + - +
- - + - -
- - - - -
- - - + +
+ + ]
1501 [ - + + - : 24 : ABI_MALLOC(gs_hamk_out%pspso,(gs_hamk_out%ntypat))
+ - + - -
+ - + -
+ ]
1502 [ + - + + : 48 : gs_hamk_out%pspso = gs_hamk_in%pspso
+ - - + -
- - - - -
+ + ]
1503 [ + + ]: 48 : tmp2i(1:5)=shape(gs_hamk_in%ekb_spin)
1504 [ - + + - : 56 : ABI_MALLOC(gs_hamk_out%ekb_spin,(tmp2i(1),tmp2i(2),tmp2i(3),tmp2i(4),tmp2i(5)))
- + + - -
+ + - - +
+ - - + +
- + - + -
- + - + -
+ ]
1505 [ + - + + : 192 : gs_hamk_out%ekb_spin = gs_hamk_in%ekb_spin
+ + + + +
+ + + + -
+ - + - +
- + - - +
- - - - -
- - - - -
- - - - -
- - - - -
- - + + +
+ + + + +
+ + ]
1506 : 8 : gs_hamk_out%ekb => gs_hamk_out%ekb_spin(:,:,:,:,1)
1507 [ + + ]: 24 : tmp2i(1:2)=shape(gs_hamk_in%sij)
1508 [ + - - - : 32 : ABI_MALLOC(gs_hamk_out%sij,(tmp2i(1),tmp2i(2)))
+ - - - +
- - + - +
- + - + ]
1509 [ - + - + : 16 : gs_hamk_out%sij = gs_hamk_in%sij
- - + - +
- - + - -
- - - - -
- - - - +
- - ]
1510 : :
1511 [ + - ]: 8 : if (associated(gs_hamk_in%gbound_kp,gs_hamk_in%gbound_k)) then
1512 : 8 : gs_hamk_out%gbound_kp => gs_hamk_out%gbound_k
1513 : : else
1514 [ # # # # : 0 : ABI_MALLOC(gs_hamk_out%gbound_kp,(2,gs_hamk_out%natom))
# # # # #
# # # ]
1515 [ # # # # : 0 : gs_hamk_out%gbound_kp = gs_hamk_in%gbound_kp
# # # # #
# ]
1516 : : end if
1517 [ + - ]: 8 : if (associated(gs_hamk_in%phkpxred,gs_hamk_in%phkxred)) then
1518 : 8 : gs_hamk_out%phkpxred => gs_hamk_out%phkxred
1519 : : else
1520 [ # # # # : 0 : ABI_MALLOC(gs_hamk_out%phkpxred,(2,gs_hamk_out%natom))
# # # # #
# # # ]
1521 [ # # # # : 0 : gs_hamk_out%phkpxred = gs_hamk_in%phkpxred
# # # # #
# ]
1522 : : end if
1523 : :
1524 : 8 : call addr_copy(gs_hamk_in%xred,gs_hamk_out%xred)
1525 : 8 : call addr_copy(gs_hamk_in%vectornd,gs_hamk_out%vectornd)
1526 : 8 : call addr_copy(gs_hamk_in%vlocal,gs_hamk_out%vlocal)
1527 : 8 : call addr_copy(gs_hamk_in%vxctaulocal,gs_hamk_out%vxctaulocal)
1528 : 8 : call addr_copy(gs_hamk_in%kinpw_k,gs_hamk_out%kinpw_k)
1529 : 8 : call addr_copy(gs_hamk_in%kinpw_kp,gs_hamk_out%kinpw_kp)
1530 : 8 : call addr_copy(gs_hamk_in%kg_k,gs_hamk_out%kg_k)
1531 : 8 : call addr_copy(gs_hamk_in%kg_kp,gs_hamk_out%kg_kp)
1532 : 8 : call addr_copy(gs_hamk_in%kpg_k,gs_hamk_out%kpg_k)
1533 : 8 : call addr_copy(gs_hamk_in%kpg_kp,gs_hamk_out%kpg_kp)
1534 : 8 : call addr_copy(gs_hamk_in%ffnl_k,gs_hamk_out%ffnl_k)
1535 : 8 : call addr_copy(gs_hamk_in%ffnl_kp,gs_hamk_out%ffnl_kp)
1536 : 8 : call addr_copy(gs_hamk_in%ph3d_k,gs_hamk_out%ph3d_k)
1537 : 8 : call addr_copy(gs_hamk_in%ph3d_kp,gs_hamk_out%ph3d_kp)
1538 : :
1539 : : !For pointers to structured datatypes, have to copy the address
1540 : : !manually because there is no generic addr_copy function for that
1541 [ - + ]: 8 : if (associated(gs_hamk_in%fockcommon)) then
1542 : : #if defined HAVE_FC_ISO_C_BINDING
1543 : 0 : ham_ptr=c_loc(gs_hamk_in%fockcommon)
1544 : 0 : call c_f_pointer(ham_ptr,gs_hamk_out%fockcommon)
1545 : : #else
1546 : : gs_hamk_out%fockcommon=transfer(gs_hamk_in%fockcommon,gs_hamk_out%fockcommon)
1547 : : #endif
1548 : : else
1549 : 8 : nullify(gs_hamk_out%fockcommon)
1550 : : end if
1551 [ - + ]: 8 : if (associated(gs_hamk_in%fockbz)) then
1552 : : #if defined HAVE_FC_ISO_C_BINDING
1553 : 0 : ham_ptr=c_loc(gs_hamk_in%fockbz)
1554 : 0 : call c_f_pointer(ham_ptr,gs_hamk_out%fockbz)
1555 : : #else
1556 : : gs_hamk_out%fockbz=transfer(gs_hamk_in%fockbz,gs_hamk_out%fockbz)
1557 : : #endif
1558 : : else
1559 : 8 : nullify(gs_hamk_out%fockbz)
1560 : : end if
1561 [ - + ]: 8 : if (associated(gs_hamk_in%fockACE_k)) then
1562 : : #if defined HAVE_FC_ISO_C_BINDING
1563 : 0 : ham_ptr=c_loc(gs_hamk_in%fockACE_k)
1564 : 0 : call c_f_pointer(ham_ptr,gs_hamk_out%fockACE_k)
1565 : : #else
1566 : : gs_hamk_out%fockACE_k=transfer(gs_hamk_in%fockACE_k,gs_hamk_out%fockACE_k)
1567 : : #endif
1568 : : else
1569 : 8 : nullify(gs_hamk_out%fockACE_k)
1570 : : end if
1571 : :
1572 : : DBG_EXIT("COLL")
1573 : :
1574 : 8 : end subroutine gsham_copy
1575 : : !!***
1576 : :
1577 : : !----------------------------------------------------------------------
1578 : :
1579 : : !!****f* m_hamiltonian/gsham_load_spin
1580 : : !! NAME
1581 : : !! gsham_load_spin
1582 : : !!
1583 : : !! INPUTS
1584 : : !! isppol=index of current spin
1585 : : !! [vectornd(n4,n5,n6,nvloc,3)]=optional, vector potential of nuclear magnetic dipoles in real space
1586 : : !! [vlocal(n4,n5,n6,nvloc)]=optional, local potential in real space
1587 : : !! [vxctaulocal(n4,n5,n6,nvloc,4)]=optional, derivative of XC energy density with respect
1588 : : !! to kinetic energy density in real space
1589 : : !! [with_nonlocal]=optional, true if non-local factors have to be loaded
1590 : : !!
1591 : : !! FUNCTION
1592 : : !! Setup of the spin-dependent part of the GS Hamiltonian.
1593 : : !!
1594 : : !! SIDE EFFECTS
1595 : : !! Ham<gs_hamiltonian_type>=Structured datatype initialization phase:
1596 : : !! * Quantities that depend spin are initialized.
1597 : : !!
1598 : : !! SOURCE
1599 : :
1600 [ + + + - : 199924 : subroutine gsham_load_spin(Ham,isppol,vectornd,vlocal,vxctaulocal,with_nonlocal)
- + + + +
- - + + +
+ - - + ]
1601 : :
1602 : : !Arguments ------------------------------------
1603 : : !scalars
1604 : : class(gs_hamiltonian_type),intent(inout),target :: Ham
1605 : : integer,intent(in) :: isppol
1606 : : logical,optional,intent(in) :: with_nonlocal
1607 : : !arrays
1608 : : real(dp),optional,intent(in),target :: vectornd(:,:,:,:,:)
1609 : : real(dp),optional,intent(in),target :: vlocal(:,:,:,:),vxctaulocal(:,:,:,:,:)
1610 : :
1611 : : !Local variables-------------------------------
1612 : : integer :: jsppol
1613 : : ! *************************************************************************
1614 : :
1615 : : DBG_ENTER("COLL")
1616 : :
1617 : : !@gs_hamiltonian_type
1618 [ + + ]: 199924 : if (present(vlocal)) then
1619 [ + + - + : 924260 : ABI_CHECK_IEQ(size(vlocal), Ham%n4*Ham%n5*Ham%n6*Ham%nvloc, "Wrong vlocal")
- - ]
1620 : 184852 : Ham%vlocal => vlocal
1621 : : end if
1622 [ + + ]: 199924 : if (present(vxctaulocal)) then
1623 [ + + - + : 4560 : ABI_CHECK_IEQ(size(vxctaulocal), Ham%n4*Ham%n5*Ham%n6*Ham%nvloc*4, "Wrong vxctaulocal")
- - ]
1624 : 760 : Ham%vxctaulocal => vxctaulocal
1625 : : end if
1626 [ + + ]: 199924 : if (present(vectornd)) then
1627 [ + + - + : 6402 : ABI_CHECK_IEQ(size(vectornd), Ham%n4*Ham%n5*Ham%n6*Ham%nvloc*3, "Wrong vectornd")
- - ]
1628 : 1067 : Ham%vectornd => vectornd
1629 : : end if
1630 : :
1631 : : ! Retrieve non-local factors for this spin component
1632 [ + + ]: 199924 : if (present(with_nonlocal)) then
1633 [ + - ]: 198097 : if (with_nonlocal) then
1634 : 198097 : jsppol=min(isppol,size(Ham%ekb_spin,5))
1635 [ + + ]: 198097 : if (jsppol>0) Ham%ekb => Ham%ekb_spin(:,:,:,:,jsppol)
1636 : : end if
1637 : : end if
1638 : :
1639 : : ! Update enl and sij on GPU
1640 : : #if defined HAVE_GPU_CUDA
1641 : : if (Ham%gpu_option==ABI_GPU_LEGACY .or. Ham%gpu_option==ABI_GPU_KOKKOS) then
1642 : : call gpu_update_ham_data(&
1643 : : & Ham%ekb(:,:,:,1), INT(size(Ham%ekb), c_int64_t), &
1644 : : & Ham%sij, INT(size(Ham%sij), c_int64_t), &
1645 : : & Ham%gprimd, INT(size(Ham%gprimd),c_int64_t))
1646 : : end if
1647 : : #endif
1648 : :
1649 : : DBG_EXIT("COLL")
1650 : :
1651 : 199924 : end subroutine gsham_load_spin
1652 : : !!***
1653 : :
1654 : : !!****f* m_hamiltonian/gsham_alloc_fofr_work
1655 : : !! NAME
1656 : : !! gsham_alloc_fofr_work
1657 : : !!
1658 : : !! FUNCTION
1659 : : !!
1660 : : !! INPUTS
1661 : : !!
1662 : : !! SOURCE
1663 : :
1664 : 0 : subroutine gsham_alloc_fofr_work(gs_ham, ndat)
1665 : :
1666 : : !Arguments ------------------------------------
1667 : : class(gs_hamiltonian_type),target,intent(inout) :: gs_ham
1668 : : integer,intent(in) :: ndat
1669 : :
1670 : : !Local variables-------------------------------
1671 : 0 : real(dp), contiguous, pointer :: fofr_work_ptr(:,:,:,:)
1672 : : ! *************************************************************************
1673 : :
1674 : 0 : if (.not. allocated(gs_ham%fofr_work)) then
1675 : : !print *, "first allocation"
1676 : : ! First allocation on CPU and GPU.
1677 [ # # # # : 0 : ABI_MALLOC(gs_ham%fofr_work, (2, gs_ham%n4, gs_ham%n5, gs_ham%n6*ndat))
# # # # #
# # # # #
# # # # #
# # # ]
1678 : 0 : fofr_work_ptr => gs_ham%fofr_work
1679 : : #ifdef HAVE_OPENMP_OFFLOAD
1680 : : !$OMP TARGET ENTER DATA MAP(alloc:fofr_work_ptr) IF (gs_ham%gpu_option==ABI_GPU_OPENMP)
1681 : : #endif
1682 : : end if
1683 : :
1684 : : ! Realloc and remap if the buffer is not large enough.
1685 : : !if (gs_ham%n6*ndat > size(gs_ham%fofr_work, dim=4)) then
1686 : : ! Realloc and remap if buffer size changed.
1687 [ # # ]: 0 : if (gs_ham%n6*ndat /= size(gs_ham%fofr_work, dim=4)) then
1688 : : !print *, "reallocating:", gs_ham%n6*ndat, size(gs_ham%fofr_work, dim=4)
1689 : 0 : fofr_work_ptr => gs_ham%fofr_work
1690 : : #ifdef HAVE_OPENMP_OFFLOAD
1691 : : !$OMP TARGET EXIT DATA MAP(delete:fofr_work_ptr) IF (gs_ham%gpu_option==ABI_GPU_OPENMP)
1692 : : #endif
1693 [ # # # # : 0 : ABI_REMALLOC(gs_ham%fofr_work, (2, gs_ham%n4, gs_ham%n5, gs_ham%n6*ndat))
# # # # #
# # # # #
# # # # #
# # # #
# ]
1694 : 0 : fofr_work_ptr => gs_ham%fofr_work
1695 : : #ifdef HAVE_OPENMP_OFFLOAD
1696 : : !$OMP TARGET ENTER DATA MAP(alloc:fofr_work_ptr) IF (gs_ham%gpu_option==ABI_GPU_OPENMP)
1697 : : #endif
1698 : : end if ! realloc condition.
1699 : :
1700 [ # # ]: 0 : end subroutine gsham_alloc_fofr_work
1701 : : !!***
1702 : :
1703 : : !!****f* m_hamiltonian/gsham_print
1704 : : !! NAME
1705 : : !! gsham_print
1706 : : !!
1707 : : !! FUNCTION
1708 : : !! Print info on the object.
1709 : : !!
1710 : : !! INPUTS
1711 : : !! units=Unit numbers for output.
1712 : : !! header=Header string
1713 : : !! prtvol=verbosity level
1714 : : !!
1715 : : !! SOURCE
1716 : :
1717 : 0 : subroutine gsham_print(gs_ham, units, header, prtvol)
1718 : :
1719 : : !Arguments ------------------------------------
1720 : : class(gs_hamiltonian_type),intent(in) :: gs_ham
1721 : : integer,intent(in) :: units(:), prtvol
1722 : : character(len=*),intent(in) :: header
1723 : :
1724 : : !Local variables-------------------------------
1725 : : character(len=5000) :: msg
1726 : : ! *************************************************************************
1727 : :
1728 [ # # ]: 0 : if (len_trim(header) == 0) then
1729 : 0 : msg = ' ==== Info on the gs_hamiltonian_type object ==== '
1730 : : else
1731 [ # # # # ]: 0 : msg = ' ==== '//trim(adjustl(header))//' ==== '
1732 : : end if
1733 : 0 : call wrtout(units, msg)
1734 : :
1735 : 0 : call wrtout(units, sjoin(" kpt_k:", ktoa(gs_ham%kpt_k)))
1736 : 0 : call wrtout(units, sjoin(" kpt_kp:", ktoa(gs_ham%kpt_kp)))
1737 : 0 : call wrtout(units, sjoin(" npw_k:", itoa(gs_ham%npw_k)))
1738 : 0 : call wrtout(units, sjoin(" npw_kp:", itoa(gs_ham%npw_kp)))
1739 : 0 : call wrtout(units, sjoin(" associated(vlocal):", yesno(associated(gs_ham%vlocal))))
1740 : :
1741 : : if (prtvol > 0) then
1742 : : end if
1743 : :
1744 [ # # ]: 0 : end subroutine gsham_print
1745 : : !!***
1746 : :
1747 : : !----------------------------------------------------------------------
1748 : :
1749 : : !!****f* m_hamiltonian/rfham_free
1750 : : !! NAME
1751 : : !! rfham_free
1752 : : !!
1753 : : !! FUNCTION
1754 : : !! Clean and destroy rf_hamiltonian_type datastructure
1755 : : !!
1756 : : !! SOURCE
1757 : :
1758 : 2156007 : subroutine rfham_free(rf_Ham)
1759 : :
1760 : : !Arguments ------------------------------------
1761 : : class(rf_hamiltonian_type),intent(inout) :: rf_Ham
1762 : : ! *************************************************************************
1763 : :
1764 : : DBG_ENTER("COLL")
1765 : :
1766 : : !@rf_hamiltonian_type
1767 : :
1768 : : ! Real arrays
1769 [ + - ]: 2156007 : ABI_SFREE(rf_Ham%e1kbfr_spin)
1770 [ + - ]: 2156007 : ABI_SFREE(rf_Ham%e1kbsc_spin)
1771 : :
1772 : : ! Real pointers
1773 [ + + ]: 2156007 : if (associated(rf_Ham%dkinpw_k)) nullify(rf_Ham%dkinpw_k)
1774 [ + + ]: 2156007 : if (associated(rf_Ham%dkinpw_kp)) nullify(rf_Ham%dkinpw_kp)
1775 [ + + ]: 2156007 : if (associated(rf_Ham%ddkinpw_k)) nullify(rf_Ham%ddkinpw_k)
1776 [ + + ]: 2156007 : if (associated(rf_Ham%ddkinpw_kp)) nullify(rf_Ham%ddkinpw_kp)
1777 [ + + ]: 2156007 : if (associated(rf_Ham%vectornd)) nullify(rf_Ham%vectornd)
1778 [ + + ]: 2156007 : if (associated(rf_Ham%vlocal1)) nullify(rf_Ham%vlocal1)
1779 [ + + ]: 2156007 : if (associated(rf_Ham%vxctaulocal)) nullify(rf_Ham%vxctaulocal)
1780 [ + + ]: 2156007 : if (associated(rf_Ham%e1kbfr)) nullify(rf_Ham%e1kbfr)
1781 [ + + ]: 2156007 : if (associated(rf_Ham%e1kbsc)) nullify(rf_Ham%e1kbsc)
1782 : :
1783 : : DBG_EXIT("COLL")
1784 : :
1785 : 2156007 : end subroutine rfham_free
1786 : : !!***
1787 : :
1788 : : !----------------------------------------------------------------------
1789 : :
1790 : : !!****f* m_hamiltonian/rfham_init
1791 : : !! NAME
1792 : : !! rfham_init
1793 : : !!
1794 : : !! FUNCTION
1795 : : !! Creation method for the rf_hamiltonian_type structure.
1796 : : !! It allocates memory and initializes all quantities that do not depend on the k-point or spin.
1797 : : !!
1798 : : !! INPUTS
1799 : : !! [comm_atom]=optional, MPI communicator over atoms
1800 : : !! cplex_paw=1 if all on-site PAW quantities are real (GS), 2 if they are complex (RF)
1801 : : !! gs_Ham<gs_hamiltonian_type>=Structured datatype containing data for ground-state Hamiltonian at (k+q)
1802 : : !! [has_e1kbsc]=optional, true if rf_Ham%e1kbsc has to be initialized.
1803 : : !! e1kbsc contains the self-consistent 1st-order PAW Dij coefficients (depending on VHxc^(1))
1804 : : !! ipert=index of perturbation
1805 : : !! [mpi_atmtab(:)]=optional, indexes of the atoms treated by current proc
1806 : : !! [mpi_spintab(2)]=optional, flags defining the spin(s) treated be current process:
1807 : : !! mpi_spintab(1)=1 if non-polarized or spin-up treated
1808 : : !! mpi_spintab(2)=1 if polarized and spin-dn treated
1809 : : !! [paw_ij1(:)<paw_ij_type>]=Various 1st-order arrays given on (i,j) (partial waves)
1810 : : !! channels (paw_ij1%dij and paw_ij1%difr only used here).
1811 : : !!
1812 : : !! SIDE EFFECTS
1813 : : !! rf_Ham<rf_hamiltonian_type>=Structured datatype almost completely initialized:
1814 : : !! * Basic variables and dimensions are transferred to the structure.
1815 : : !! * All pointers are allocated with correct dimensions.
1816 : : !! * Quantities that do not depend on the k-point or spin are initialized.
1817 : : !!
1818 : : !! SOURCE
1819 : :
1820 : 2229707 : subroutine rfham_init(rf_ham, cplex, gs_Ham, ipert,&
1821 [ + + + + : 2229707 : comm_atom, mpi_atmtab, mpi_spintab, paw_ij1, has_e1kbsc) ! optional arguments
- + + + +
+ - + ]
1822 : :
1823 : : !Arguments ------------------------------------
1824 : : !scalars
1825 : : class(rf_hamiltonian_type),intent(inout),target :: rf_Ham
1826 : : integer,intent(in) :: cplex,ipert
1827 : : integer,intent(in),optional :: comm_atom
1828 : : logical,intent(in),optional :: has_e1kbsc
1829 : : type(gs_hamiltonian_type),intent(in) :: gs_Ham
1830 : : !arrays
1831 : : integer,optional,intent(in) :: mpi_atmtab(:),mpi_spintab(2)
1832 : : type(paw_ij_type),optional,intent(in) :: paw_ij1(:)
1833 : :
1834 : : !Local variables-------------------------------
1835 : : !scalars
1836 : : integer :: cplex_dij1,isp,jsp,my_comm_atom,my_nsppol
1837 : : logical :: has_e1kbsc_
1838 : : !arrays
1839 : : integer :: my_spintab(2)
1840 : 2229707 : real(dp),allocatable,target :: e1kb_tmp(:,:,:,:)
1841 : : ! *************************************************************************
1842 : :
1843 : : DBG_ENTER("COLL")
1844 : :
1845 : : !@rf_hamiltonian_type
1846 : :
1847 : : !Manage optional parameters
1848 [ + + ]: 2229707 : has_e1kbsc_=.false.;if (present(has_e1kbsc)) has_e1kbsc_=has_e1kbsc
1849 [ + + ]: 2229707 : my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
1850 [ + + + + ]: 4474938 : my_spintab=0;my_spintab(1:gs_Ham%nsppol)=1;if(present(mpi_spintab)) my_spintab=mpi_spintab
1851 [ + + + + ]: 6689121 : my_nsppol=count(my_spintab==1)
1852 : :
1853 : 2229707 : rf_Ham%cplex =cplex
1854 : 2229707 : rf_Ham%n4 =gs_Ham%n4
1855 : 2229707 : rf_Ham%n5 =gs_Ham%n5
1856 : 2229707 : rf_Ham%n6 =gs_Ham%n6
1857 : 2229707 : rf_Ham%nvloc =gs_Ham%nvloc
1858 : 2229707 : rf_Ham%nsppol =gs_Ham%nsppol
1859 : 2229707 : rf_Ham%nspinor =gs_Ham%nspinor
1860 : 2229707 : rf_Ham%zora =gs_Ham%zora
1861 : :
1862 : 2229707 : rf_Ham%dime1kb1=0
1863 : 2229707 : rf_Ham%dime1kb2=gs_Ham%dimekb2
1864 [ + + + + : 2229707 : if (gs_Ham%usepaw==1.and.ipert/=gs_Ham%natom+1.and.ipert/=gs_Ham%natom+10) then
+ + ]
1865 [ + + + + : 61684 : cplex_dij1=1;if ((gs_Ham%nspinor==2).or.any(abs(gs_Ham%nucdipmom)>tol8)) cplex_dij1=2
+ - + - -
+ ]
1866 : 4144 : rf_Ham%dime1kb1=cplex_dij1*(gs_Ham%lmnmax*(gs_Ham%lmnmax+1))/2
1867 : : end if
1868 : :
1869 : : ! Allocate the arrays of the 1st-order Hamiltonian
1870 : : ! We preload here 1st-order non-local factors in order to avoid
1871 : : ! a communication over atoms inside the loop over spins.
1872 [ + + + + ]: 2229707 : if (gs_Ham%usepaw==1.and.rf_Ham%dime1kb1>0) then
1873 : : if ((ipert>=1.and.ipert<=gs_Ham%natom).or.ipert==gs_Ham%natom+2.or.&
1874 [ + - + + : 4144 : ipert==gs_Ham%natom+3.or.ipert==gs_Ham%natom+4.or.ipert==gs_Ham%natom+11) then
+ + + + +
+ - + ]
1875 : :
1876 [ - + + - : 29008 : ABI_MALLOC(rf_Ham%e1kbfr_spin,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex,my_nsppol))
- + + - -
+ + - - +
+ - + + +
- + - + +
- + - + -
+ ]
1877 [ + + + + : 484083 : rf_Ham%e1kbfr_spin=zero
+ + + + +
+ ]
1878 [ + + ]: 4144 : if (has_e1kbsc_) then
1879 [ - + + - : 20671 : ABI_MALLOC(rf_Ham%e1kbsc_spin,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex,my_nsppol))
- + + - -
+ + - - +
+ - - + +
- + - + -
- + - + -
+ ]
1880 [ + + + + : 467155 : rf_Ham%e1kbsc_spin=zero
+ + + + +
+ ]
1881 : : end if
1882 : :
1883 [ + + ]: 4144 : if (present(paw_ij1)) then
1884 : :
1885 [ - + ]: 3068 : if (my_nsppol<rf_Ham%nsppol) then
1886 [ # # # # : 0 : ABI_MALLOC(e1kb_tmp,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex))
# # # # #
# # # # #
# # # # #
# # # #
# ]
1887 : : end if
1888 : :
1889 : : ! === Frozen term
1890 : 3068 : jsp=0
1891 [ + + ]: 6816 : do isp=1,rf_Ham%nsppol
1892 [ + - ]: 3748 : if (my_spintab(isp)==1) then
1893 : 3748 : jsp=jsp+1 ; rf_Ham%e1kbfr => rf_Ham%e1kbfr_spin(:,:,:,:,jsp)
1894 : : else
1895 : 0 : rf_Ham%e1kbfr => e1kb_tmp
1896 : : end if
1897 [ + + ]: 6816 : if (present(mpi_atmtab)) then
1898 : 1314 : call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbfr=rf_Ham%e1kbfr,mpi_atmtab=mpi_atmtab)
1899 : : else
1900 : 2434 : call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbfr=rf_Ham%e1kbfr)
1901 : : end if
1902 : : end do
1903 : :
1904 : : ! === Self-consistent term
1905 [ + + ]: 3068 : if (has_e1kbsc_) then
1906 : 2953 : jsp=0
1907 [ + + ]: 6586 : do isp=1,rf_Ham%nsppol
1908 [ + - ]: 3633 : if (my_spintab(isp)==1) then
1909 : 3633 : jsp=jsp+1 ; rf_Ham%e1kbsc => rf_Ham%e1kbsc_spin(:,:,:,:,jsp)
1910 : : else
1911 : 0 : rf_Ham%e1kbsc => e1kb_tmp
1912 : : end if
1913 [ + + ]: 6586 : if (present(mpi_atmtab)) then
1914 : 1254 : call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbsc=rf_Ham%e1kbsc,mpi_atmtab=mpi_atmtab)
1915 : : else
1916 : 2379 : call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbsc=rf_Ham%e1kbsc)
1917 : : end if
1918 : : end do
1919 : : end if
1920 : :
1921 [ - + ]: 3068 : if (my_nsppol<rf_Ham%nsppol) then
1922 [ # # ]: 0 : ABI_FREE(e1kb_tmp)
1923 : : end if
1924 : :
1925 : : end if
1926 : : end if
1927 : : end if
1928 : :
1929 [ + + ]: 2229707 : if (.not.allocated(rf_Ham%e1kbfr_spin)) then
1930 [ - + - + ]: 2151863 : ABI_MALLOC(rf_Ham%e1kbfr_spin,(0,0,0,0,0))
1931 : : end if
1932 [ + + ]: 2229707 : if (.not.allocated(rf_Ham%e1kbsc_spin)) then
1933 [ - + - + ]: 2153054 : ABI_MALLOC(rf_Ham%e1kbsc_spin,(0,0,0,0,0))
1934 : : end if
1935 : 2229707 : nullify(rf_Ham%e1kbfr)
1936 : 2229707 : nullify(rf_Ham%e1kbsc)
1937 : :
1938 : : DBG_EXIT("COLL")
1939 : :
1940 [ - + ]: 2229707 : end subroutine rfham_init
1941 : : !!***
1942 : :
1943 : : !----------------------------------------------------------------------
1944 : :
1945 : : !!****f* m_hamiltonian/rfham_load_spin
1946 : : !! NAME
1947 : : !! rfham_load_spin
1948 : : !!
1949 : : !! FUNCTION
1950 : : !! Setup of the spin-dependent part of the 1st- and 2nd- order Hamiltonian.
1951 : : !!
1952 : : !! INPUTS
1953 : : !! isppol=index of current spin
1954 : : !! [vectornd(n4,n5,n6,nvloc)]=optional, vector potential of nuclear magnetic dipoles in real space in
1955 : : !! ddk direction idir
1956 : : !! [vlocal1(cplex*n4,n5,n6,nvloc)]=optional, 1st-order local potential in real space
1957 : : !! [vxctaulocal(n4,n5,n6,nvloc,4)]=optional, deriv of e_XC wrt kin energy, for mGGA
1958 : : !! [with_nonlocal]=optional, true if non-local factors have to be loaded
1959 : : !!
1960 : : !! SIDE EFFECTS
1961 : : !! rf_Ham<rf_hamiltonian_type>=Structured datatype initialization phase:
1962 : : !! * Quantities that depend on spin are initialized.
1963 : : !!
1964 : : !! SOURCE
1965 : :
1966 : 903832 : subroutine rfham_load_spin(rf_Ham, isppol, &
1967 [ + + + - : 903832 : vectornd, vlocal1, vxctaulocal, with_nonlocal) ! optional
- + + + +
- - + + +
+ - - + ]
1968 : :
1969 : : !Arguments ------------------------------------
1970 : : !scalars
1971 : : class(rf_hamiltonian_type),intent(inout),target :: rf_Ham
1972 : : integer,intent(in) :: isppol
1973 : : logical,optional,intent(in) :: with_nonlocal
1974 : : !arrays
1975 : : real(dp),optional,target,intent(in) :: vlocal1(:,:,:,:)
1976 : : real(dp),optional,target,intent(in) :: vectornd(:,:,:,:)
1977 : : real(dp),optional,target,intent(in) :: vxctaulocal(:,:,:,:,:)
1978 : :
1979 : : !Local variables-------------------------------
1980 : : integer :: jsppol
1981 : : ! *************************************************************************
1982 : :
1983 : : DBG_ENTER("COLL")
1984 : :
1985 : : !@rf_hamiltonian_type
1986 : :
1987 [ + + ]: 903832 : if (present(vlocal1)) then
1988 [ + + - + : 4509005 : ABI_CHECK_IEQ(size(vlocal1), rf_Ham%cplex*rf_Ham%n4*rf_Ham%n5*rf_Ham%n6*rf_Ham%nvloc, "Wrong vlocal1")
- - ]
1989 : 901801 : rf_Ham%vlocal1 => vlocal1
1990 : : end if
1991 : :
1992 [ + + ]: 903832 : if (present(vectornd)) then
1993 [ + + - + : 2610 : ABI_CHECK_IEQ(size(vectornd), rf_Ham%cplex*rf_Ham%n4*rf_Ham%n5*rf_Ham%n6*rf_Ham%nvloc, "Wrong vectornd")
- - ]
1994 : 522 : rf_Ham%vectornd => vectornd
1995 : : end if
1996 : :
1997 [ + + ]: 903832 : if (present(vxctaulocal)) then
1998 [ + + - + : 288 : ABI_CHECK_IEQ(size(vxctaulocal), rf_Ham%n4*rf_Ham%n5*rf_Ham%n6*rf_Ham%nvloc*4, "Wrong vxctaulocal")
- - ]
1999 : 48 : rf_Ham%vxctaulocal => vxctaulocal
2000 : : end if
2001 : :
2002 : : ! Retrieve non-local factors for this spin component
2003 [ + + ]: 903832 : if (present(with_nonlocal)) then
2004 [ + - ]: 903015 : if (with_nonlocal) then
2005 [ + + + + ]: 5418090 : if (size(rf_Ham%e1kbfr_spin)>0) then
2006 : 3748 : jsppol=min(isppol,size(rf_Ham%e1kbfr_spin,5))
2007 [ + - ]: 3748 : if (jsppol>0) rf_Ham%e1kbfr => rf_Ham%e1kbfr_spin(:,:,:,:,jsppol)
2008 : : end if
2009 [ + + + + ]: 5418090 : if (size(rf_Ham%e1kbsc_spin)>0) then
2010 : 3633 : jsppol=min(isppol,size(rf_Ham%e1kbsc_spin,5))
2011 [ + - ]: 3633 : if (jsppol>0) rf_Ham%e1kbsc => rf_Ham%e1kbsc_spin(:,:,:,:,jsppol)
2012 : : end if
2013 : : end if
2014 : : end if
2015 : :
2016 : : DBG_EXIT("COLL")
2017 : :
2018 : 903832 : end subroutine rfham_load_spin
2019 : : !!***
2020 : :
2021 : : !----------------------------------------------------------------------
2022 : :
2023 : : !!****f* m_hamiltonian/rfham_load_k
2024 : : !! NAME
2025 : : !! rfham_load_k
2026 : : !!
2027 : : !! FUNCTION
2028 : : !! Setup of the k-dependent part of the 1st- and 2nd- order Hamiltonian
2029 : : !!
2030 : : !! INPUTS
2031 : : !! [dkinpw_k]=1st derivative of the (modified) kinetic energy for each plane wave
2032 : : !! [ddkinpw_k]=2nd derivative of the (modified) kinetic energy for each plane wave
2033 : : !! [npw_k]=number of plane waves
2034 : : !!
2035 : : !! SIDE EFFECTS
2036 : : !! rf_Ham<rf_hamiltonian_type>=structured datatype completed with k-dependent quantities.
2037 : : !! Quantities at k^prime are set equal to quantities at k.
2038 : : !!
2039 : : !! SOURCE
2040 : :
2041 : 1323580 : subroutine rfham_load_k(rf_Ham, &
2042 [ + + + - : 1323580 : dkinpw_k, ddkinpw_k, npw_k) ! optional.
- + + + +
- - + ]
2043 : :
2044 : : !Arguments ------------------------------------
2045 : : !scalars
2046 : : class(rf_hamiltonian_type),intent(inout),target :: rf_Ham
2047 : : integer,intent(in),optional :: npw_k
2048 : : !arrays
2049 : : real(dp),intent(in),optional,target :: dkinpw_k(:),ddkinpw_k(:)
2050 : : ! *************************************************************************
2051 : :
2052 : : DBG_ENTER("COLL")
2053 : :
2054 : : !@gs_hamiltonian_type
2055 : :
2056 : : !k-dependent scalars
2057 [ + + ]: 1323580 : if (present(npw_k)) then
2058 : 1305364 : rf_Ham%npw_k = npw_k
2059 : 1305364 : rf_Ham%npw_kp = npw_k
2060 : : end if
2061 : :
2062 : : !Pointers to k-dependent quantities
2063 [ + + ]: 1323580 : if (present(dkinpw_k)) then
2064 : 1312420 : rf_Ham%dkinpw_k => dkinpw_k
2065 : 1312420 : rf_Ham%dkinpw_kp => dkinpw_k
2066 : : end if
2067 [ + + ]: 1323580 : if (present(ddkinpw_k)) then
2068 : 18216 : rf_Ham%ddkinpw_k => ddkinpw_k
2069 : 18216 : rf_Ham%ddkinpw_kp => ddkinpw_k
2070 : : end if
2071 : :
2072 : : DBG_EXIT("COLL")
2073 : :
2074 : 1323580 : end subroutine rfham_load_k
2075 : : !!***
2076 : :
2077 : : !----------------------------------------------------------------------
2078 : :
2079 : : !!****f* m_hamiltonian/pawdij2ekb
2080 : : !! NAME
2081 : : !! pawdij2ekb
2082 : : !!
2083 : : !! FUNCTION
2084 : : !! Transfer PAW Dij (on-site GS Hamiltonian) values
2085 : : !! from paw_ij datastructure to ekb array
2086 : : !!
2087 : : !! INPUTS
2088 : : !!
2089 : : !! OUTPUT
2090 : : !!
2091 : : !! SOURCE
2092 : :
2093 [ - + - + : 25482 : subroutine pawdij2ekb(ekb,paw_ij,isppol,comm_atom,mpi_atmtab)
+ + + - -
+ ]
2094 : :
2095 : : !Arguments ------------------------------------
2096 : : !scalars
2097 : : integer,intent(in) :: isppol,comm_atom
2098 : : !arrays
2099 : : integer,intent(in),optional,target :: mpi_atmtab(:)
2100 : : real(dp),intent(out) :: ekb(:,:,:,:)
2101 : : type(paw_ij_type),intent(in) :: paw_ij(:)
2102 : :
2103 : : !Local variables-------------------------------
2104 : : !scalars
2105 : : integer :: dimdij,dimekb1,dimekb3,dimekb4,iatom,iatom_tot,ierr,ii,isp,ispden,my_natom,natom,qphase
2106 : : logical :: my_atmtab_allocated,paral_atom
2107 : : !arrays
2108 : 25482 : integer,pointer :: my_atmtab(:)
2109 : : ! *************************************************************************
2110 : :
2111 : : DBG_ENTER("COLL")
2112 : :
2113 [ + + + + : 8353929 : ekb=zero
+ + + + ]
2114 : :
2115 : : !Set up parallelism over atoms
2116 : 25482 : natom=size(ekb,2); my_natom=size(paw_ij)
2117 : 25482 : paral_atom=(xmpi_comm_size(comm_atom)>1)
2118 [ + + ]: 25482 : nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
2119 : 25482 : call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
2120 : :
2121 : : !Retrieve PAW Dij coefficients for this spin component
2122 [ + + ]: 25482 : if (my_natom>0) then
2123 [ + - ]: 23496 : if (allocated(paw_ij(1)%dij)) then
2124 : 23496 : dimekb1=size(ekb,1) ; dimekb3=size(ekb,3) ; dimekb4=size(ekb,4)
2125 : 23496 : qphase=paw_ij(1)%qphase
2126 [ - + ]: 23496 : ABI_CHECK(qphase <= dimekb4, 'paw_ij%qphase>dimekb4!')
2127 [ + + ]: 46992 : do ii=1,qphase
2128 [ + + ]: 77466 : do ispden=1,dimekb3
2129 [ + + ]: 30474 : isp=isppol; if (dimekb3==4) isp=ispden
2130 [ + + ]: 119238 : do iatom=1,my_natom
2131 [ + + ]: 65268 : iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom)
2132 : 65268 : dimdij=paw_ij(iatom)%cplex_dij*paw_ij(iatom)%lmn2_size
2133 [ - + ]: 65268 : ABI_CHECK(dimdij <= dimekb1, 'Size of paw_ij%dij>dimekb1!')
2134 [ + + ]: 6090219 : ekb(1:dimdij,iatom_tot,ispden,ii)=paw_ij(iatom)%dij(1+(ii-1)*dimdij:ii*dimdij,isp)
2135 : : end do
2136 : : end do
2137 : : end do
2138 : : end if
2139 : : end if
2140 : :
2141 : : ! Communication in case of distribution over atomic sites
2142 [ + + ]: 25482 : if (paral_atom) call xmpi_sum(ekb,comm_atom,ierr)
2143 : :
2144 : : ! Destroy atom table used for parallelism
2145 : 25482 : call free_my_atmtab(my_atmtab,my_atmtab_allocated)
2146 : :
2147 : : DBG_EXIT("COLL")
2148 : :
2149 : 25482 : end subroutine pawdij2ekb
2150 : : !!***
2151 : :
2152 : : !----------------------------------------------------------------------
2153 : :
2154 : : !!****f* m_hamiltonian/pawdij2e1kb
2155 : : !! NAME
2156 : : !! pawdij2e1kb
2157 : : !!
2158 : : !! FUNCTION
2159 : : !! Transfer PAW Dij (on-site RF Hamiltonian) values
2160 : : !! from paw_ij datastructure to e1kb array
2161 : : !!
2162 : : !! INPUTS
2163 : : !!
2164 : : !! OUTPUT
2165 : : !!
2166 : : !! SOURCE
2167 : :
2168 [ + + + - : 10984 : subroutine pawdij2e1kb(paw_ij1,isppol,comm_atom,mpi_atmtab,e1kbfr,e1kbsc)
- + + + +
- - + - +
+ + + + -
+ ]
2169 : :
2170 : : !Arguments ------------------------------------
2171 : : !scalars
2172 : : integer,intent(in) :: isppol,comm_atom
2173 : : !arrays
2174 : : integer,intent(in),optional,target :: mpi_atmtab(:)
2175 : : real(dp),optional,intent(out) :: e1kbfr(:,:,:,:),e1kbsc(:,:,:,:)
2176 : : type(paw_ij_type),intent(in) :: paw_ij1(:)
2177 : :
2178 : : !Local variables-------------------------------
2179 : : !scalars
2180 : : integer :: dimdij1,dime1kb1,dime1kb3,dime1kb4,iatom,iatom_tot,ierr,isp,ispden
2181 : : integer :: my_natom,natom,qphase
2182 : : logical :: my_atmtab_allocated,paral_atom
2183 : : !arrays
2184 : 10984 : integer,pointer :: my_atmtab(:)
2185 : : ! *************************************************************************
2186 : :
2187 : : DBG_ENTER("COLL")
2188 : :
2189 : 0 : if ((.not.present(e1kbfr)).and.(.not.present(e1kbsc))) return
2190 [ + + ]: 10984 : if (present(e1kbfr)) then
2191 [ + + + + : 1041260 : e1kbfr=zero ; natom=size(e1kbfr,2)
+ + + + ]
2192 : : end if
2193 [ + + ]: 10984 : if (present(e1kbsc)) then
2194 [ + + + + : 464202 : e1kbsc=zero ; natom=size(e1kbsc,2)
+ + + + ]
2195 : : end if
2196 : :
2197 : : ! Set up parallelism over atoms
2198 : 10984 : my_natom=size(paw_ij1) ; paral_atom=(xmpi_comm_size(comm_atom)>1)
2199 [ + + ]: 10984 : nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
2200 : 10984 : call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
2201 : :
2202 : : ! Retrieve 1st-order PAW Dij coefficients for this spin component (frozen)
2203 [ + + + + ]: 10984 : if (my_natom>0.and.present(e1kbfr)) then
2204 [ + - ]: 7289 : if (allocated(paw_ij1(1)%dijfr)) then
2205 : 7289 : dime1kb1=size(e1kbfr,1) ; dime1kb3=size(e1kbfr,3) ; dime1kb4=size(e1kbfr,4)
2206 [ - + ]: 7289 : ABI_CHECK_IEQ(paw_ij1(1)%qphase, dime1kb4,'BUG in pawdij2e1kb (1)!')
2207 [ + + ]: 14578 : do ispden=1,dime1kb3
2208 [ - + ]: 7289 : isp=isppol;if (dime1kb3==4) isp=ispden
2209 [ + + ]: 29700 : do iatom=1,my_natom
2210 [ + + ]: 15122 : iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom)
2211 : 15122 : qphase=paw_ij1(iatom)%qphase
2212 : 15122 : dimdij1=paw_ij1(iatom)%cplex_dij*paw_ij1(iatom)%lmn2_size
2213 [ - + ]: 15122 : ABI_CHECK(dimdij1<=dime1kb1,'BUG: size of paw_ij1%dij>dime1kb1!')
2214 [ + + ]: 587176 : e1kbfr(1:dimdij1,iatom_tot,ispden,1)=paw_ij1(iatom)%dijfr(1:dimdij1,isp)
2215 [ + + + + ]: 60859 : if (qphase==2) e1kbfr(1:dimdij1,iatom_tot,ispden,2)=paw_ij1(iatom)%dijfr(dimdij1+1:2*dimdij1,isp)
2216 : : end do
2217 : : end do
2218 : : end if
2219 : : end if
2220 : :
2221 : : ! Retrieve 1st-order PAW Dij coefficients for this spin component (self-consistent)
2222 [ + + + + ]: 10984 : if (my_natom>0.and.present(e1kbsc)) then
2223 [ + - + - ]: 3613 : if (allocated(paw_ij1(1)%dijfr).and.allocated(paw_ij1(1)%dij)) then
2224 : 3613 : dime1kb1=size(e1kbsc,1) ; dime1kb3=size(e1kbsc,3) ; dime1kb4=size(e1kbsc,4)
2225 [ - + ]: 3613 : ABI_CHECK_IEQ(paw_ij1(1)%qphase, dime1kb4, 'BUG in pawdij2e1kb (1)!')
2226 [ + + ]: 7226 : do ispden=1,dime1kb3
2227 [ - + ]: 3613 : isp=isppol;if (dime1kb3==4) isp=ispden
2228 [ + + ]: 14577 : do iatom=1,my_natom
2229 [ + + ]: 7351 : iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom)
2230 : 7351 : qphase=paw_ij1(iatom)%qphase
2231 : 7351 : dimdij1=paw_ij1(iatom)%cplex_dij*paw_ij1(iatom)%lmn2_size
2232 [ - + ]: 7351 : ABI_CHECK(dimdij1<=dime1kb1, 'BUG: size of paw_ij1%dij>dime1kb1!')
2233 : : e1kbsc(1:dimdij1,iatom_tot,ispden,1)=paw_ij1(iatom)%dij (1:dimdij1,isp) &
2234 [ + + ]: 281835 : & -paw_ij1(iatom)%dijfr(1:dimdij1,isp)
2235 [ + + ]: 7351 : if (qphase==2) e1kbsc(1:dimdij1,iatom_tot,ispden,2)=paw_ij1(iatom)%dij (dimdij1+1:2*dimdij1,isp) &
2236 [ + + ]: 28165 : & -paw_ij1(iatom)%dijfr(dimdij1+1:2*dimdij1,isp)
2237 : : end do
2238 : : end do
2239 : : end if
2240 : : end if
2241 : :
2242 : : ! Communication in case of distribution over atomic sites
2243 [ + + ]: 10984 : if (paral_atom) then
2244 [ + + ]: 4494 : if (present(e1kbfr)) call xmpi_sum(e1kbfr,comm_atom,ierr)
2245 [ + + ]: 4494 : if (present(e1kbsc)) call xmpi_sum(e1kbsc,comm_atom,ierr)
2246 : : end if
2247 : :
2248 : : ! Destroy atom table used for parallelism
2249 : 10984 : call free_my_atmtab(my_atmtab,my_atmtab_allocated)
2250 : :
2251 : : DBG_EXIT("COLL")
2252 : :
2253 [ - + ]: 10984 : end subroutine pawdij2e1kb
2254 : : !!***
2255 : :
2256 : : !!****f* ABINIT/gspot_transgrid_and_pack
2257 : : !! NAME
2258 : : !! gspot_transgrid_and_pack
2259 : : !!
2260 : : !! FUNCTION
2261 : : !! Set up local potential vlocal on the coarse FFT mesh with proper dimensioning from vtrial given on the fine mesh.
2262 : : !! Also take into account the spin.
2263 : : !!
2264 : : !! INPUTS
2265 : : !! isppol=Spin polarization.
2266 : : !! usepaw=1 if PAW
2267 : : !! paral_kgb: 1 if paral_kgb
2268 : : !! nfft=(effective) number of FFT grid points on the coarse mesh (for this processor)
2269 : : !! ngfft(18)contain all needed information about 3D FFT, for the coarse FFT mesh. see ~abinit/doc/variables/vargs.htm#ngfft
2270 : : !! nfftf=(effective) number of FFT grid points on the fine mesh (for this processor)
2271 : : !! nvloc==1 if nspden <=2, nvloc==4 for nspden==4,
2272 : : !! nspden=Number of spin density components.
2273 : : !! ncomp=Number of extra components in vtrial and vlocal (e.g. 1 if LDA/GGA pot, 4 for Meta-GGA, etc).
2274 : : !! pawfgr <type(pawfgr_type)>=fine grid parameters and related data
2275 : : !! vtrial(nfftf,nspden)=INPUT potential Vtrial(r).
2276 : : !! mpi_enreg=information about MPI parallelization
2277 : : !!
2278 : : !! OUTPUT
2279 : : !! vlocal(n4,n5,n6,nvloc,ncomp): Potential on the coarse grid.
2280 : : !!
2281 : : !! SOURCE
2282 : :
2283 : 55857 : subroutine gspot_transgrid_and_pack(isppol, usepaw, paral_kgb, nfft, ngfft, nfftf, &
2284 : 55857 : nspden, nvloc, ncomp, pawfgr, mpi_enreg, vtrial, vlocal)
2285 : :
2286 : : !Arguments -------------------------------
2287 : : integer,intent(in) :: isppol, nspden, ncomp, usepaw, paral_kgb, nfft, nfftf, nvloc
2288 : : type(pawfgr_type), intent(in) :: pawfgr
2289 : : type(MPI_type), intent(in) :: mpi_enreg
2290 : : !arrays
2291 : : integer,intent(in) :: ngfft(18)
2292 : : real(dp),intent(inout) :: vtrial(nfftf, nspden, ncomp)
2293 : : real(dp),intent(out) :: vlocal(ngfft(4), ngfft(5), ngfft(6), nvloc, ncomp)
2294 : :
2295 : : !Local variables-------------------------------
2296 : : !scalars
2297 : : integer :: n1,n2,n3,n4,n5,n6,ispden,ic
2298 : : real(dp) :: rhodum(1)
2299 : 55857 : real(dp),allocatable :: cgrvtrial(:,:), vlocal_tmp(:,:,:)
2300 : : ! *************************************************************************
2301 : :
2302 : : ! Coarse mesh.
2303 : 55857 : n1=ngfft(1); n2=ngfft(2); n3=ngfft(3)
2304 : 55857 : n4=ngfft(4); n5=ngfft(5); n6=ngfft(6)
2305 : :
2306 : : ! Set up local potential vlocal with proper dimensioning, from vtrial
2307 : : ! Also take into account the spin.
2308 [ + + ]: 55857 : if (nspden /= 4) then
2309 [ + + + + ]: 53391 : if (usepaw == 0 .or. pawfgr%usefinegrid == 0) then
2310 : : ! Fine mesh == Coarse mesh.
2311 [ + + ]: 70169 : do ic=1,ncomp
2312 : 70169 : call fftpac(isppol,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,vtrial(:,:,ic),vlocal(:,:,:,:,ic),2)
2313 : : end do
2314 : : else
2315 : : ! Transfer from fine mesh to coarse and then pack data
2316 [ - + + - : 74608 : ABI_MALLOC(cgrvtrial,(nfft,nspden))
- + + - +
- + - - +
+ - ]
2317 [ + + ]: 40967 : do ic=1,ncomp
2318 : 22315 : call transgrid(1,mpi_enreg,nspden,-1,0,0,paral_kgb,pawfgr, rhodum,rhodum,cgrvtrial,vtrial(:,:,ic))
2319 : 40967 : call fftpac(isppol,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft, cgrvtrial,vlocal(:,:,:,:,ic),2)
2320 : : end do
2321 [ - + ]: 18652 : ABI_FREE(cgrvtrial)
2322 : : end if
2323 : : else
2324 : : ! nspden == 4. replace isppol by loop over ispden.
2325 [ - + + - : 12330 : ABI_MALLOC(vlocal_tmp, (n4,n5,n6))
- + + - -
+ + - + -
+ - - + -
+ ]
2326 [ + + - + ]: 2466 : if (usepaw == 0 .or. pawfgr%usefinegrid == 0) then
2327 : : ! Fine mesh == Coarse mesh.
2328 [ + + ]: 2768 : do ic=1,ncomp
2329 [ + + ]: 8304 : do ispden=1,nspden
2330 : 5536 : call fftpac(ispden,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,vtrial(:,:,ic),vlocal_tmp,2)
2331 [ + + + + : 71371728 : vlocal(:,:,:,ispden,ic) = vlocal_tmp(:,:,:)
+ + ]
2332 : : end do
2333 : : end do
2334 : : else
2335 : : ! Transfer from fine mesh to coarse and then pack data
2336 [ - + + - : 4328 : ABI_MALLOC(cgrvtrial,(nfft,nspden))
- + + - +
- + - - +
+ - ]
2337 [ + + ]: 2164 : do ic=1,ncomp
2338 : 1082 : call transgrid(1,mpi_enreg,nspden,-1,0,0,paral_kgb,pawfgr,rhodum,rhodum,cgrvtrial,vtrial(:,:,ic))
2339 [ + + ]: 6492 : do ispden=1,nspden
2340 : 4328 : call fftpac(ispden,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,cgrvtrial,vlocal_tmp,2)
2341 [ + + + + : 21768666 : vlocal(:,:,:,ispden,ic) = vlocal_tmp(:,:,:)
+ + ]
2342 : : end do
2343 : : end do
2344 [ - + ]: 1082 : ABI_FREE(cgrvtrial)
2345 : : end if
2346 [ - + ]: 2466 : ABI_FREE(vlocal_tmp)
2347 : : end if ! nspden
2348 : :
2349 : 55857 : end subroutine gspot_transgrid_and_pack
2350 : : !!***
2351 : :
2352 [ - + - - : 24 : end module m_hamiltonian
+ + - + +
- - + - +
+ - - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ ]
2353 : : !!***
|
