TABLE OF CONTENTS


ABINIT/m_hamiltonian [ Modules ]

[ Top ] [ Modules ]

NAME

 m_hamiltonian

FUNCTION

  This module provides the definition of the gs_hamiltonian_type and rf_hamiltonian_type
  datastructures used in the "getghc" and "getgh1c" routines to apply the Hamiltonian (or
  its derivative) on a wavefunction.
  Methods to initialize or destroy the objects are defined here.
  It also defines the ddiago_ctl_type structures datatype used to control the algorithm
  used in ks_ddiago for performing the direct diagonalization of the KS Hamiltonian.

COPYRIGHT

 Copyright (C) 2009-2018 ABINIT group (MG, MT)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

CHILDREN

SOURCE

26 #if defined HAVE_CONFIG_H
27 #include "config.h"
28 #endif
29 
30 #include "abi_common.h"
31 
32 module m_hamiltonian
33 
34  use defs_basis
35  use defs_datatypes
36  use defs_abitypes
37  use m_abicore
38  use m_errors
39  use m_xmpi
40 
41  use m_copy,              only : addr_copy
42  use m_geometry,          only : metric
43  use m_pawtab,            only : pawtab_type
44  use m_fftcore,           only : sphereboundary
45  use m_pawcprj,           only : pawcprj_getdim
46  use m_paw_ij,            only : paw_ij_type
47  use m_paral_atom,        only : get_my_atmtab, free_my_atmtab
48  use m_electronpositron,  only : electronpositron_type, electronpositron_calctype
49  use m_kg,                only : ph1d3d, getph
50  use m_fock,              only : fock_common_type, fock_BZ_type, fock_ACE_type, fock_type
51 
52 #if defined HAVE_FC_ISO_C_BINDING
53  use iso_c_binding, only : c_ptr,c_loc,c_f_pointer
54 #endif
55 
56  implicit none
57 
58  private
59 
60  public ::  pawdij2ekb
61  public ::  pawdij2e1kb
62 
63 !These constants select how H is applied in reciprocal space
64  integer,parameter,public :: KPRIME_H_K=1, K_H_KPRIME=2, K_H_K=3, KPRIME_H_KPRIME=4

m_hamiltonian/copy_hamiltonian [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  copy_hamiltonian

INPUTS

  gs_hamk_in<gs_hamiltonian_type>=Structured datatype completely initialized,
                                  to be copied.

FUNCTION

  Copy a gs_hamiltonian_type variable (gs_hamk_in) in another (gs_hamk_out).
  In contrast to an assignment statement (gs_hamk_out=gs_hamk_in), this
  subroutine allocate memory space for the pointers contained in the data
  structure (gs_hamk_out) and copy the content of the corresponding memory
  space of gs_hamk_in in it. In contrast, the assignment statement would
  only associate the pointers of gs_hamk_out to the same memory space than
  the corresponding ones in gs_hamk_in. This can cause trouble if one data
  structure is destroyed before a reading/writing statement for the other
  structure, causing access to unallocated memory space (silently, without
  segmentation fault being generated).

OUTPUT

  gs_hamk_out<gs_hamiltonian_type>=Structured datatype containing separate
                                   copies of all data of gs_hamk_in upon exit.

PARENTS

      gwls_hamiltonian

CHILDREN

      destroy_mpi_enreg,initmpi_seq,kpgsph,wrtout

SOURCE

1297 subroutine copy_hamiltonian(gs_hamk_out,gs_hamk_in)
1298 
1299 
1300 !This section has been created automatically by the script Abilint (TD).
1301 !Do not modify the following lines by hand.
1302 #undef ABI_FUNC
1303 #define ABI_FUNC 'copy_hamiltonian'
1304 !End of the abilint section
1305 
1306 implicit none
1307 
1308 !Arguments ------------------------------------
1309  type(gs_hamiltonian_type),intent(in),target :: gs_hamk_in
1310  type(gs_hamiltonian_type),intent(out),target :: gs_hamk_out
1311 
1312 !Local variables-------------------------------
1313  integer :: tmp2i(5)
1314 #if defined HAVE_FC_ISO_C_BINDING
1315  type(C_PTR) :: ham_ptr
1316 #endif
1317 
1318 ! *************************************************************************
1319 
1320  DBG_ENTER("COLL")
1321 
1322 !@gs_hamiltonian_type
1323 
1324  gs_hamk_out%dimekb1 = gs_hamk_in%dimekb1
1325  gs_hamk_out%dimekb2 = gs_hamk_in%dimekb2
1326  gs_hamk_out%dimekbq = gs_hamk_in%dimekbq
1327  gs_hamk_out%istwf_k = gs_hamk_in%istwf_k
1328  gs_hamk_out%istwf_kp = gs_hamk_in%istwf_kp
1329  gs_hamk_out%lmnmax = gs_hamk_in%lmnmax
1330  gs_hamk_out%matblk = gs_hamk_in%matblk
1331  gs_hamk_out%mgfft = gs_hamk_in%mgfft
1332  gs_hamk_out%mpsang = gs_hamk_in%mpsang
1333  gs_hamk_out%mpssoang = gs_hamk_in%mpssoang
1334  gs_hamk_out%natom = gs_hamk_in%natom
1335  gs_hamk_out%nfft = gs_hamk_in%nfft
1336  gs_hamk_out%npw_k = gs_hamk_in%npw_k
1337  gs_hamk_out%npw_kp = gs_hamk_in%npw_kp
1338  gs_hamk_out%npw_fft_k = gs_hamk_in%npw_fft_k
1339  gs_hamk_out%npw_fft_kp = gs_hamk_in%npw_fft_kp
1340  gs_hamk_out%nspinor = gs_hamk_in%nspinor
1341  gs_hamk_out%nsppol = gs_hamk_in%nsppol
1342  gs_hamk_out%ntypat = gs_hamk_in%ntypat
1343  gs_hamk_out%nvloc = gs_hamk_in%nvloc
1344  gs_hamk_out%n4 = gs_hamk_in%n4
1345  gs_hamk_out%n5 = gs_hamk_in%n5
1346  gs_hamk_out%n6 = gs_hamk_in%n6
1347  gs_hamk_out%use_gpu_cuda = gs_hamk_in%use_gpu_cuda
1348  gs_hamk_out%usecprj = gs_hamk_in%usecprj
1349  gs_hamk_out%usepaw = gs_hamk_in%usepaw
1350  gs_hamk_out%useylm = gs_hamk_in%useylm
1351  gs_hamk_out%ngfft = gs_hamk_in%ngfft
1352  gs_hamk_out%nloalg = gs_hamk_in%nloalg
1353  gs_hamk_out%ucvol = gs_hamk_in%ucvol
1354  gs_hamk_out%gmet = gs_hamk_in%gmet
1355  gs_hamk_out%gprimd = gs_hamk_in%gprimd
1356  gs_hamk_out%kpt_k = gs_hamk_in%kpt_k
1357  gs_hamk_out%kpt_kp = gs_hamk_in%kpt_kp
1358 
1359  ABI_ALLOCATE(gs_hamk_out%atindx,(gs_hamk_out%natom))
1360  gs_hamk_out%atindx = gs_hamk_in%atindx
1361  ABI_ALLOCATE(gs_hamk_out%atindx1,(gs_hamk_out%natom))
1362  gs_hamk_out%atindx1 = gs_hamk_in%atindx1
1363  ABI_ALLOCATE(gs_hamk_out%dimcprj,(gs_hamk_out%natom*gs_hamk_out%usepaw))
1364  if (gs_hamk_out%usepaw==1) gs_hamk_out%dimcprj = gs_hamk_in%dimcprj
1365  ABI_ALLOCATE(gs_hamk_out%typat,(gs_hamk_out%natom))
1366  gs_hamk_out%typat = gs_hamk_in%typat
1367  ABI_ALLOCATE(gs_hamk_out%gbound_k,(2*gs_hamk_out%mgfft+8,2))
1368  gs_hamk_out%gbound_k = gs_hamk_in%gbound_k
1369  ABI_ALLOCATE(gs_hamk_out%indlmn,(6,gs_hamk_out%lmnmax,gs_hamk_out%ntypat))
1370  gs_hamk_out%indlmn = gs_hamk_in%indlmn
1371  ABI_ALLOCATE(gs_hamk_out%nattyp,(gs_hamk_out%ntypat))
1372  gs_hamk_out%nattyp = gs_hamk_in%nattyp
1373  ABI_ALLOCATE(gs_hamk_out%nucdipmom,(3,gs_hamk_out%natom))
1374  gs_hamk_out%nucdipmom = gs_hamk_in%nucdipmom
1375  ABI_ALLOCATE(gs_hamk_out%phkxred,(2,gs_hamk_out%natom))
1376  gs_hamk_out%phkxred = gs_hamk_in%phkxred
1377  ABI_ALLOCATE(gs_hamk_out%ph1d,(2,3*(2*gs_hamk_out%mgfft+1)*gs_hamk_out%natom))
1378  gs_hamk_out%ph1d = gs_hamk_in%ph1d
1379  ABI_ALLOCATE(gs_hamk_out%pspso,(gs_hamk_out%ntypat))
1380  gs_hamk_out%pspso = gs_hamk_in%pspso
1381  tmp2i(1:5)=shape(gs_hamk_in%ekb_spin)
1382  ABI_ALLOCATE(gs_hamk_out%ekb_spin,(tmp2i(1),tmp2i(2),tmp2i(3),tmp2i(4),tmp2i(5)))
1383  gs_hamk_out%ekb_spin = gs_hamk_in%ekb_spin
1384  gs_hamk_out%ekb => gs_hamk_out%ekb_spin(:,:,:,:,1)
1385  tmp2i(1:2)=shape(gs_hamk_in%sij)
1386  ABI_ALLOCATE(gs_hamk_out%sij,(tmp2i(1),tmp2i(2)))
1387  gs_hamk_out%sij = gs_hamk_in%sij
1388 
1389  if (associated(gs_hamk_in%gbound_kp,gs_hamk_in%gbound_k)) then
1390    gs_hamk_out%gbound_kp => gs_hamk_out%gbound_k
1391  else
1392    ABI_ALLOCATE(gs_hamk_out%gbound_kp,(2,gs_hamk_out%natom))
1393    gs_hamk_out%gbound_kp = gs_hamk_in%gbound_kp
1394  end if
1395  if (associated(gs_hamk_in%phkpxred,gs_hamk_in%phkxred)) then
1396    gs_hamk_out%phkpxred => gs_hamk_out%phkxred
1397  else
1398    ABI_ALLOCATE(gs_hamk_out%phkpxred,(2,gs_hamk_out%natom))
1399    gs_hamk_out%phkpxred = gs_hamk_in%phkpxred
1400  end if
1401 
1402  call addr_copy(gs_hamk_in%xred,gs_hamk_out%xred)
1403  call addr_copy(gs_hamk_in%vlocal,gs_hamk_out%vlocal)
1404  call addr_copy(gs_hamk_in%vxctaulocal,gs_hamk_out%vxctaulocal)
1405  call addr_copy(gs_hamk_in%kinpw_k,gs_hamk_out%kinpw_k)
1406  call addr_copy(gs_hamk_in%kinpw_kp,gs_hamk_out%kinpw_kp)
1407  call addr_copy(gs_hamk_in%kg_k,gs_hamk_out%kg_k)
1408  call addr_copy(gs_hamk_in%kg_kp,gs_hamk_out%kg_kp)
1409  call addr_copy(gs_hamk_in%kpg_k,gs_hamk_out%kpg_k)
1410  call addr_copy(gs_hamk_in%kpg_kp,gs_hamk_out%kpg_kp)
1411  call addr_copy(gs_hamk_in%ffnl_k,gs_hamk_out%ffnl_k)
1412  call addr_copy(gs_hamk_in%ffnl_kp,gs_hamk_out%ffnl_kp)
1413  call addr_copy(gs_hamk_in%ph3d_k,gs_hamk_out%ph3d_k)
1414  call addr_copy(gs_hamk_in%ph3d_kp,gs_hamk_out%ph3d_kp)
1415 
1416 !For pointers to structured datatypes, have to copy the address
1417 !manually because there is no generic addr_copy function for that
1418  if (associated(gs_hamk_in%fockcommon)) then
1419 #if defined HAVE_FC_ISO_C_BINDING
1420    ham_ptr=c_loc(gs_hamk_in%fockcommon)
1421    call c_f_pointer(ham_ptr,gs_hamk_out%fockcommon)
1422 #else
1423    gs_hamk_out%fockcommon=transfer(gs_hamk_in%fockcommon,gs_hamk_out%fockcommon)
1424 #endif
1425  else
1426    nullify(gs_hamk_out%fockcommon)
1427  end if
1428  if (associated(gs_hamk_in%fockbz)) then
1429 #if defined HAVE_FC_ISO_C_BINDING
1430    ham_ptr=c_loc(gs_hamk_in%fockbz)
1431    call c_f_pointer(ham_ptr,gs_hamk_out%fockbz)
1432 #else
1433    gs_hamk_out%fockbz=transfer(gs_hamk_in%fockbz,gs_hamk_out%fockbz)
1434 #endif
1435  else
1436    nullify(gs_hamk_out%fockbz)
1437  end if
1438  if (associated(gs_hamk_in%fockACE_k)) then
1439 #if defined HAVE_FC_ISO_C_BINDING
1440    ham_ptr=c_loc(gs_hamk_in%fockACE_k)
1441    call c_f_pointer(ham_ptr,gs_hamk_out%fockACE_k)
1442 #else
1443    gs_hamk_out%fockACE_k=transfer(gs_hamk_in%fockACE_k,gs_hamk_out%fockACE_k)
1444 #endif
1445  else
1446    nullify(gs_hamk_out%fockACE_k)
1447  end if
1448 
1449  DBG_EXIT("COLL")
1450 
1451 end subroutine copy_hamiltonian

m_hamiltonian/destroy_hamiltonian [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  destroy_hamiltonian

FUNCTION

  Clean and destroy gs_hamiltonian_type datastructure

SIDE EFFECTS

  Ham<gs_hamiltonian_type>=All dynamic memory defined in the structure is deallocated.

PARENTS

      d2frnl,dfpt_nselt,dfpt_nstdy,dfpt_nstpaw,dfpt_rhofermi,dfpt_vtorho
      dfptnl_resp,energy,fock2ACE,forstrnps,gwls_hamiltonian,ks_ddiago,m_gkk
      m_io_kss,m_phgamma,m_phpi,m_shirley,m_sigmaph,nonlop_test,vtorho

CHILDREN

      destroy_mpi_enreg,initmpi_seq,kpgsph,wrtout

SOURCE

526 subroutine destroy_hamiltonian(Ham)
527 
528 !Arguments ------------------------------------
529 !scalars
530 
531 !This section has been created automatically by the script Abilint (TD).
532 !Do not modify the following lines by hand.
533 #undef ABI_FUNC
534 #define ABI_FUNC 'destroy_hamiltonian'
535 !End of the abilint section
536 
537  type(gs_hamiltonian_type),intent(inout),target :: Ham
538 
539 ! *************************************************************************
540 
541  DBG_ENTER("COLL")
542 
543 !@gs_hamiltonian_type
544 
545 ! Integer Pointers
546  if (associated(Ham%gbound_kp,Ham%gbound_k))   then
547    nullify(Ham%gbound_kp)
548  else if (associated(Ham%gbound_kp)) then
549    ABI_DEALLOCATE(Ham%gbound_kp)
550  end if
551 
552 ! Integer arrays
553  if (allocated(Ham%atindx))  then
554    ABI_DEALLOCATE(Ham%atindx)
555  end if
556  if (allocated(Ham%atindx1))  then
557    ABI_DEALLOCATE(Ham%atindx1)
558  end if
559  if (allocated(Ham%dimcprj))  then
560    ABI_DEALLOCATE(Ham%dimcprj)
561  end if
562  if (allocated(Ham%gbound_k))  then
563    ABI_DEALLOCATE(Ham%gbound_k)
564  end if
565  if (allocated(Ham%indlmn ))  then
566    ABI_DEALLOCATE(Ham%indlmn)
567  end if
568  if (allocated(Ham%nattyp))  then
569    ABI_DEALLOCATE(Ham%nattyp)
570  end if
571  if (allocated(Ham%pspso))  then
572    ABI_DEALLOCATE(Ham%pspso)
573  end if
574  if (allocated(Ham%typat))  then
575    ABI_DEALLOCATE(Ham%typat)
576  end if
577 
578 ! Real Pointers
579  if (associated(Ham%phkpxred,Ham%phkxred))   then
580    nullify(Ham%phkpxred)
581  else if (associated(Ham%phkpxred)) then
582    ABI_DEALLOCATE(Ham%phkpxred)
583  end if
584  if (allocated(Ham%phkxred))   then
585    ABI_DEALLOCATE(Ham%phkxred)
586  end if
587  if (associated(Ham%ekb)) nullify(Ham%ekb)
588  if (associated(Ham%vlocal)) nullify(Ham%vlocal)
589  if (associated(Ham%vxctaulocal)) nullify(Ham%vxctaulocal)
590  if (associated(Ham%xred)) nullify(Ham%xred)
591  if (associated(Ham%kinpw_k)) nullify(Ham%kinpw_k)
592  if (associated(Ham%kinpw_kp)) nullify(Ham%kinpw_kp)
593  if (associated(Ham%kg_k)) nullify(Ham%kg_k)
594  if (associated(Ham%kg_kp)) nullify(Ham%kg_kp)
595  if (associated(Ham%kpg_k)) nullify(Ham%kpg_k)
596  if (associated(Ham%kpg_kp)) nullify(Ham%kpg_kp)
597  if (associated(Ham%ffnl_k)) nullify(Ham%ffnl_k)
598  if (associated(Ham%ffnl_kp)) nullify(Ham%ffnl_kp)
599  if (associated(Ham%ph3d_k)) nullify(Ham%ph3d_k)
600  if (associated(Ham%ph3d_kp)) nullify(Ham%ph3d_kp)
601 
602 ! Real arrays
603  if (allocated(Ham%ekb_spin))   then
604    ABI_DEALLOCATE(Ham%ekb_spin)
605  end if
606  if (allocated(Ham%sij))   then
607    ABI_DEALLOCATE(Ham%sij)
608  end if
609  if (allocated(Ham%nucdipmom)) then
610    ABI_DEALLOCATE(Ham%nucdipmom)
611  end if
612  if (allocated(Ham%ph1d))   then
613    ABI_DEALLOCATE(Ham%ph1d)
614  end if
615 
616 ! Complex arrays
617  if(allocated(Ham%nucdipmom_k)) then
618    ABI_DEALLOCATE(Ham%nucdipmom_k)
619  end if
620 
621 ! Structured datatype pointers
622  if (associated(Ham%fockcommon)) nullify(Ham%fockcommon)
623  if (associated(Ham%fockACE_k)) nullify(Ham%fockACE_k)
624  if (associated(Ham%fockbz)) nullify(Ham%fockbz)
625 #if defined HAVE_GPU_CUDA
626  if(Ham%use_gpu_cuda==1) then
627    call gpu_finalize_ham_data()
628  end if
629 #endif
630 
631  DBG_EXIT("COLL")
632 
633 end subroutine destroy_hamiltonian

m_hamiltonian/destroy_rf_hamiltonian [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  destroy_rf_hamiltonian

FUNCTION

  Clean and destroy rf_hamiltonian_type datastructure

SIDE EFFECTS

  rf_Ham<rf_hamiltonian_type>=All dynamic memory defined in the structure is deallocated.

PARENTS

      dfpt_nstpaw,dfpt_nstwf,dfpt_rhofermi,dfpt_vtorho,m_gkk,m_phgamma,m_phpi
      m_sigmaph

CHILDREN

      destroy_mpi_enreg,initmpi_seq,kpgsph,wrtout

SOURCE

1563 subroutine destroy_rf_hamiltonian(rf_Ham)
1564 
1565 !Arguments ------------------------------------
1566 !scalars
1567 
1568 !This section has been created automatically by the script Abilint (TD).
1569 !Do not modify the following lines by hand.
1570 #undef ABI_FUNC
1571 #define ABI_FUNC 'destroy_rf_hamiltonian'
1572 !End of the abilint section
1573 
1574  type(rf_hamiltonian_type),intent(inout) :: rf_Ham
1575 
1576 ! *************************************************************************
1577 
1578  DBG_ENTER("COLL")
1579 
1580 !@rf_hamiltonian_type
1581 
1582 ! Real arrays
1583  if (allocated(rf_Ham%e1kbfr_spin))   then
1584    ABI_DEALLOCATE(rf_Ham%e1kbfr_spin)
1585  end if
1586  if (allocated(rf_Ham%e1kbsc_spin))   then
1587    ABI_DEALLOCATE(rf_Ham%e1kbsc_spin)
1588  end if
1589 
1590 ! Real pointers
1591  if (associated(rf_Ham%dkinpw_k)) nullify(rf_Ham%dkinpw_k)
1592  if (associated(rf_Ham%dkinpw_kp)) nullify(rf_Ham%dkinpw_kp)
1593  if (associated(rf_Ham%ddkinpw_k)) nullify(rf_Ham%ddkinpw_k)
1594  if (associated(rf_Ham%ddkinpw_kp)) nullify(rf_Ham%ddkinpw_kp)
1595  if (associated(rf_Ham%vlocal1)) nullify(rf_Ham%vlocal1)
1596  if (associated(rf_Ham%e1kbfr)) nullify(rf_Ham%e1kbfr)
1597  if (associated(rf_Ham%e1kbsc)) nullify(rf_Ham%e1kbsc)
1598 
1599  DBG_EXIT("COLL")
1600 
1601 end subroutine destroy_rf_hamiltonian

m_hamiltonian/gs_hamiltonian_type [ Types ]

[ Top ] [ m_hamiltonian ] [ Types ]

NAME

 gs_hamiltonian_type

FUNCTION

 This datastructure contains the information about one Hamiltonian,
 needed in the "getghc" routine, that apply the Hamiltonian
 on a wavefunction.
 The Hamiltonian is expressed in reciprocal space:
       H_k^prime,k = exp(-i.k^prime.r^prime) H exp(i.k.r)
 In most cases k=k^prime and the k^prime objects are simply pointers
 to k objects.

SOURCE

 84  type,public :: gs_hamiltonian_type
 85 
 86 ! ===== Integer scalars
 87 
 88   integer :: dimekb1
 89    ! First dimension of Ekb
 90    ! Same as psps%dimekb
 91    ! ->Norm conserving : Max. number of Kleinman-Bylander energies
 92    !                     for each atom type
 93    !                     dimekb1=lnmax
 94    ! ->PAW : Max. number of Dij coefficients connecting projectors
 95    !                     for each atom
 96    !                     dimekb1=cplex_dij*lmnmax*(lmnmax+1)/2
 97 
 98   integer :: dimekb2
 99    ! Second dimension of Ekb
100    ! ->Norm conserving psps: dimekb2=ntypat
101    ! ->PAW                 : dimekb2=natom
102 
103   integer :: dimekbq
104    ! Fourth dimension of Ekb
105    ! 2 if Ekb factors contain a exp(-iqR) phase, 1 otherwise
106 
107   integer :: istwf_k
108    ! option parameter that describes the storage of wfs at k
109 
110   integer :: istwf_kp
111    ! option parameter that describes the storage of wfs at k^prime
112 
113   integer :: lmnmax
114    ! Maximum number of different l,m,n components over all types of psps.
115    ! same as dtset%lmnmax
116 
117   integer :: matblk
118    ! dimension of the array ph3d
119 
120   integer :: mgfft
121    ! maximum size for 1D FFTs (same as dtset%mgfft)
122 
123   integer :: mpsang
124    ! Highest angular momentum of non-local projectors over all type of psps.
125    ! shifted by 1 : for all local psps, mpsang=0; for largest s, mpsang=1,
126    ! for largest p, mpsang=2; for largest d, mpsang=3; for largest f, mpsang=4
127    ! This gives also the number of non-local "channels"
128    ! same as psps%mpsang
129 
130   integer :: mpssoang
131    ! Maximum number of channels, including those for treating the spin-orbit coupling
132    ! For NC pseudopotentials only:
133    !   when mpspso=1, mpssoang=mpsang
134    !   when mpspso=2, mpssoang=2*mpsang-1
135    ! For PAW: same as mpsang
136    ! same as psps%mpssoang
137 
138   integer :: natom
139    ! The number of atoms for this dataset; same as dtset%natom
140 
141   integer :: nfft
142    ! number of FFT grid points same as dtset%nfft
143 
144   integer :: npw_k
145    ! number of plane waves at k
146    ! In case of band-FFT parallelism, npw_k is the number of plane waves
147    ! processed by current proc
148 
149   integer :: npw_fft_k
150    ! number of plane waves at k used to apply Hamiltonian when band-FFT
151    ! parallelism is activated (i.e. data are distributed in the "FFT" configuration)
152 
153   integer :: npw_kp
154    ! number of plane waves at k^prime
155    ! In case of band-FFT parallelism, npw_kp is the number of plane waves
156    ! processed by current proc
157 
158   integer :: npw_fft_kp
159    ! number of plane waves at k^prime used to apply Hamiltonian when band-FFT
160    ! parallelism is activated (i.e. data are distributed in the "FFT" configuration)
161 
162   integer :: nspinor
163    ! Number of spinorial components
164 
165   integer :: nsppol
166    ! Total number of spin components (1=non-polarized, 2=polarized)
167 
168   integer :: ntypat
169    ! Number of types of pseudopotentials same as dtset%ntypat
170 
171   integer :: nvloc
172    ! Number of components of vloc
173    ! usually, nvloc=1, except in the non-collinear magnetism case, where nvloc=4
174 
175   integer :: n4,n5,n6
176    ! same as ngfft(4:6)
177 
178   integer :: use_gpu_cuda
179   ! governs wheter we do the hamiltonian calculation on gpu or not
180 
181   integer :: usecprj
182    ! usecprj= 1 if cprj projected WF are stored in memory
183    !        = 0 if they are to be computed on the fly
184 
185   integer :: usepaw
186    ! if usepaw=0 , use norm-conserving psps part of the code
187    ! is usepaw=1 , use paw part of the code
188 
189   integer :: useylm
190    ! governs the way the nonlocal operator is to be applied:
191    !   1=using Ylm, 0=using Legendre polynomials
192 
193 ! ===== Integer arrays
194 
195   integer, allocatable :: atindx(:)
196    ! atindx(natom)
197    ! index table for atoms (see gstate.f)
198 
199   integer, allocatable :: atindx1(:)
200    ! atindx1(natom)
201    ! index table for atoms, inverse of atindx (see gstate.f)
202 
203   integer, allocatable :: dimcprj(:)
204    ! dimcprj(natom*usepaw)=dimensions of array cprj
205    ! dimcprj(ia)=cprj(ia,:)%nlmn
206    ! atoms are ordered by atom-type
207 
208   integer, allocatable :: gbound_k(:,:)
209    ! gbound_k(2*mgfft+8,2)
210    ! G sphere boundary, for each plane wave at k
211 
212   integer, allocatable :: indlmn(:,:,:)
213    ! indlmn(6,lmnmax,ntypat)
214    ! For each type of psp,
215    ! array giving l,m,n,lm,ln,spin for i=ln  (if useylm=0)
216    !                                or i=lmn (if useylm=1)
217 
218   integer, allocatable :: nattyp(:)
219    ! nattyp(ntypat)
220    ! # of atoms of each type
221 
222   integer :: ngfft(18)
223    ! ngfft(1:3)=integer fft box dimensions
224    ! ngfft(4:6)=integer fft box dimensions, might be augmented for CPU speed
225    ! ngfft(7)=fftalg
226    ! ngfft(8)=fftalg
227 
228   integer :: nloalg(3)
229    ! governs the choice of the algorithm for non-local operator same as dtset%nloalg
230 
231   integer, allocatable :: pspso(:)
232    ! pspso(ntypat)
233    ! For each type of psp, 1 if no spin-orbit component is taken
234    ! into account, 2 if a spin-orbit component is used
235    ! Revelant for NC-psps and PAW
236 
237   integer, allocatable :: typat(:)
238    ! typat(natom)
239    ! type of each atom
240 
241 ! integer, allocatable :: indpw_k(:,:)
242    ! indpw_k(4,npw_fft_k)
243    ! array which gives fft box index for given basis sphere
244 
245 ! Integer pointers
246 
247   integer, ABI_CONTIGUOUS pointer :: gbound_kp(:,:) => null()
248    ! gbound_kp(2*mgfft+8,2)
249    ! G sphere boundary, for each plane wave at k^prime
250 
251   integer, pointer :: kg_k(:,:) => null()
252    ! kg_k(3,npw_fft_k)
253    ! G vector coordinates with respect to reciprocal lattice translations
254    ! at k
255 
256   integer, pointer :: kg_kp(:,:) => null()
257    ! kg_kp(3,npw_fft_kp)
258    ! G vector coordinates with respect to reciprocal lattice translations
259    ! at k^prime
260 
261 ! ===== Real scalars
262 
263   real(dp) :: ucvol
264    ! unit cell volume (Bohr**3)
265 
266 ! ===== Real arrays
267 
268   real(dp), allocatable :: ekb_spin(:,:,:,:,:)
269    ! ekb_spin(dimekb1,dimekb2,nspinor**2,dimekbq,my_nsppol)
270    ! Contains the values of ekb array for all spins treated by current process
271    ! See ekb description ; ekb is pointer to ekb_spin(:,:,:,:,my_isppol)
272 
273   real(dp), allocatable :: sij(:,:)
274    ! sij(dimekb1,ntypat*usepaw) = overlap matrix for paw calculation
275 
276   real(dp) :: gmet(3,3)
277    ! reciprocal space metric tensor in Bohr**-2
278 
279   real(dp) :: gprimd(3,3)
280    ! dimensional reciprocal space primitive translations (Bohr^-1)
281 
282   real(dp) :: kpt_k(3)
283    ! dimensionless k point coordinates wrt reciprocal lattice vectors
284 
285   real(dp) :: kpt_kp(3)
286    ! dimensionless k^prime point coordinates wrt reciprocal lattice vectors
287 
288   real(dp), allocatable :: nucdipmom(:,:)
289    ! nucdipmom(3,natom)
290    ! nuclear dipole moments at each atomic position
291 
292   real(dp), allocatable :: ph1d(:,:)
293    ! ph1d(2,3*(2*mgfft+1)*natom)
294    ! 1-dim phase arrays for structure factor (see getph.f).
295 
296   real(dp), allocatable :: phkxred(:,:)
297    ! phkxred(2,natom)
298    ! phase factors exp(2 pi k.xred) at k
299 
300 ! ===== Complex arrays
301 
302   complex(dpc), allocatable :: nucdipmom_k(:)
303    ! nucdipmom_k(npw_k*(npw_k+1)/2)
304    ! nuclear dipole moment Hamiltonian in reciprocal space, stored as
305    ! lower triangular part of Hermitian matrix
306 
307 ! ===== Real pointers
308 
309   real(dp), ABI_CONTIGUOUS pointer :: ekb(:,:,:,:) => null()
310    ! ekb(dimekb1,dimekb2,nspinor**2,dimekbq)
311    !  ->Norm conserving : (Real) Kleinman-Bylander energies (hartree)
312    !          for number of basis functions (l,n) (lnmax)
313    !          and number of atom types (ntypat)
314    !          dimekb1=lnmax ; dimekb2=ntypat ; dimekbq=1
315    !  ->PAW : (Real, symmetric) Frozen part of Dij coefficients
316    !                            to connect projectors
317    !          for number of basis functions (l,m,n) (lmnmax)
318    !          and number of atom (natom)
319    !          dimekb1=lmnmax*(lmnmax+1)/2 ; dimekb2=natom ; dimekbq=1
320    ! ekb is spin dependent in the case of PAW calculations.
321    ! For each spin component, ekb points to ekb_spin(:,:,:,:,my_isppol)
322    ! dimekbq=2 if Ekb factors contain a exp(-iqR) phase, dimekbq=1 otherwise
323    ! About the non-local factors symmetry:
324    !   - The lower triangular part of the Dij matrix can be deduced from the upper one
325    !     with the following relation: D^s2s1_ji = (D^s1s2_ij)^*
326    !     where s1,s2 are spinor components
327 
328   real(dp), pointer :: ffnl_k(:,:,:,:) => null()
329    ! ffnl_k(npw_fft_k,2,dimffnl_k,ntypat)
330    ! nonlocal form factors at k
331 
332   real(dp), pointer :: ffnl_kp(:,:,:,:) => null()
333    ! ffnl_kp(npw_fft_kp,2,dimffnl_kp,ntypat)
334    ! nonlocal form factors at k_prime
335 
336   real(dp), pointer :: kinpw_k(:) => null()
337    ! kinpw_k(npw_fft_k)
338    ! (modified) kinetic energy for each plane wave at k
339 
340   real(dp), pointer :: kinpw_kp(:) => null()
341    ! kinpw_kp(npw_fft_kp)
342    ! (modified) kinetic energy for each plane wave at k^prime
343 
344   real(dp), pointer :: kpg_k(:,:) => null()
345    ! kpg_k(3,npw_fft_k)
346    ! k+G vector coordinates at k
347 
348   real(dp), pointer :: kpg_kp(:,:) => null()
349    ! kpg_kp(3,npw_fft_kp)
350    ! k^prime+G vector coordinates at k^prime
351 
352   real(dp), ABI_CONTIGUOUS pointer :: phkpxred(:,:) => null()
353    ! phkpxred(2,natom)
354    ! phase factors exp(2 pi k^prime.xred) at k^prime
355 
356   real(dp), pointer :: ph3d_k(:,:,:) => null()
357    ! ph3d_k(2,npw_fft_k,matblk)
358    ! 3-dim structure factors, for each atom and plane wave at k
359 
360   real(dp), pointer :: ph3d_kp(:,:,:) => null()
361    ! ph3d_kp(2,npw_fft_kp,matblk)
362    ! 3-dim structure factors, for each atom and plane wave at k^prime
363 
364   real(dp), pointer :: vlocal(:,:,:,:) => null()
365    ! vlocal(n4,n5,n6,nvloc)
366    ! local potential in real space, on the augmented fft grid
367 
368   real(dp), pointer :: vxctaulocal(:,:,:,:,:) => null()
369    ! vxctaulocal(n4,n5,n6,nvloc,4)
370    ! derivative of XC energy density with respect to kinetic energy density,
371    ! in real space, on the augmented fft grid
372 
373   real(dp), ABI_CONTIGUOUS pointer :: xred(:,:) => null()
374    ! xred(3,natom)
375    ! reduced coordinates of atoms (dimensionless)
376 
377 ! ===== Structured datatype pointers
378 
379   type(fock_common_type), pointer :: fockcommon => null()
380    ! fock
381    ! common quantities needed to calculate Fock exact exchange
382 
383   type(fock_BZ_type), pointer :: fockbz => null()
384    ! fock
385    ! total brillouin zone quantities needed to calculate Fock exact exchange
386 
387   type(fock_ACE_type), pointer :: fockACE_k => null()
388    ! fock
389    ! ACE quantities needed to calculate Fock exact exchange in the ACE context
390 
391  end type gs_hamiltonian_type
392 
393 
394  public :: init_hamiltonian         ! Initialize the GS Hamiltonian
395  public :: destroy_hamiltonian      ! Free the memory in the GS Hamiltonian
396  public :: load_spin_hamiltonian    ! Setup of the spin-dependent part of the GS Hamiltonian
397  public :: load_k_hamiltonian       ! Setup of the k-dependent part of the GS Hamiltonian
398  public :: load_kprime_hamiltonian  ! Setup of the k^prime-dependent part of the GS Hamiltonian
399  public :: copy_hamiltonian         ! Copy the object

m_hamiltonian/init_hamiltonian [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  init_hamiltonian

FUNCTION

  Creation method for the gs_hamiltonian_type structure.
  It allocates memory and initializes all quantities that do not depend on the k-point or spin.

INPUTS

  [comm_atom]=optional, MPI communicator over atoms
  [fockcommon <type(fock_common_type)>]= common quantities to calculate Fock exact exchange
  [fockbz <type(fock_BZ_type)>]= quantities to calculate Fock exact exchange in the total BZ
  natom=Number of atoms in the unit cell.
  nfft=Number of FFT grid points (for this processors).
  nspinor=Number of spinorial components
  nsppol=1 for unpolarized, 2 for spin-polarized
  nspden=Number of spin density components.
  mgfft=Maximum size for 1D FFTs i.e., MAXVAL(ngfft(1:3))
  [mpi_atmtab(:)]=optional, indexes of the atoms treated by current proc
  [mpi_spintab(2)]=optional, flags defining the spin(s) treated be current process:
                   mpi_spintab(1)=1 if non-polarized or spin-up treated
                   mpi_spintab(2)=1 if polarized and spin-dn treated
  psps<pseudopotential_type>=structure datatype gathering data on the pseudopotentials.
  [electronpositron<electronpositron_type>]=Structured datatype storing data for the
    electron-positron two-component DFT (optional).
  ngfft(18)=integer array with FFT box dimensions and other information on FFTs, for the FINE rectangular grid.
  nloalg(3)=governs the choice of the algorithm for non-local operator
  [nucdipmom(3,natom)]= (optional) array of nuclear dipole moments at atomic sites
  [ph1d(2,3*(2*mgfft+1)*natom)]=1-dimensions phase arrays for structure factor (see getph.f).
              Optional, recalculated inside the routine if not present in input.
  rprimd(3,3)=Direct lattice vectors in Bohr.
  typat(natom)=Type of each atom.
  [usecprj]=flag use only for PAW; 1 if cprj datastructure is allocated
  xred(3,natom)=Reduced coordinates of the atoms.
  pawtab(ntypat*psps%usepaw)<pawtab_type>=PAW TABulated data initialized at start.
  [paw_ij(:) <type(paw_ij_type)>]=optional, paw arrays given on (i,j) channels

SIDE EFFECTS

  Ham<gs_hamiltonian_type>=Structured datatype almost completely initialized:
   * Basic variables and dimensions are transfered to the structure.
   * All pointers are allocated with correct dimensions.
   * Quantities that do not depend on the k-point or spin are initialized.

PARENTS

      d2frnl,dfpt_nselt,dfpt_nstdy,dfpt_nstpaw,dfpt_rhofermi,dfpt_vtorho
      dfptnl_resp,energy,fock2ACE,forstrnps,ks_ddiago,m_gkk,m_io_kss
      m_phgamma,m_phpi,m_shirley,m_sigmaph,nonlop_test,vtorho

CHILDREN

      destroy_mpi_enreg,initmpi_seq,kpgsph,wrtout

SOURCE

691 subroutine init_hamiltonian(ham,Psps,pawtab,nspinor,nsppol,nspden,natom,typat,&
692 &                           xred,nfft,mgfft,ngfft,rprimd,nloalg,&
693 &                           ph1d,usecprj,comm_atom,mpi_atmtab,mpi_spintab,paw_ij,&  ! optional
694 &                           electronpositron,fock,nucdipmom,use_gpu_cuda)           ! optional
695 
696 
697 !This section has been created automatically by the script Abilint (TD).
698 !Do not modify the following lines by hand.
699 #undef ABI_FUNC
700 #define ABI_FUNC 'init_hamiltonian'
701 !End of the abilint section
702 
703  implicit none
704 
705 !Arguments ------------------------------------
706 !scalars
707  integer,intent(in) :: nfft,natom,nspinor,nsppol,nspden,mgfft
708  integer,optional,intent(in) :: comm_atom,usecprj,use_gpu_cuda
709  type(gs_hamiltonian_type),intent(inout),target :: ham
710  type(electronpositron_type),optional,pointer :: electronpositron
711  type(fock_type),optional,pointer :: fock
712  type(pseudopotential_type),intent(in) :: psps
713 !arrays
714  integer,intent(in) :: ngfft(18),nloalg(3),typat(natom)
715  integer,optional,intent(in)  :: mpi_atmtab(:),mpi_spintab(2)
716  real(dp),intent(in) :: rprimd(3,3)
717  real(dp),intent(in),target :: xred(3,natom)
718  real(dp),optional,intent(in) :: nucdipmom(3,natom),ph1d(2,3*(2*mgfft+1)*natom)
719  type(pawtab_type),intent(in) :: pawtab(psps%ntypat*psps%usepaw)
720  type(paw_ij_type),optional,intent(in) :: paw_ij(:)
721 
722 !Local variables-------------------------------
723 !scalars
724  integer :: my_comm_atom,my_nsppol,itypat,iat,ilmn,indx,isp,cplex,cplex_dij,jsp
725  real(dp) :: ucvol
726 !arrays
727  integer :: my_spintab(2)
728  real(dp) :: gmet(3,3),gprimd(3,3),rmet(3,3)
729  real(dp),allocatable,target :: ekb_tmp(:,:,:,:)
730 
731 ! *************************************************************************
732 
733  DBG_ENTER("COLL")
734 
735  !@gs_hamiltonian_type
736 
737 !Manage optional parameters
738  my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
739  my_spintab=0;my_spintab(1:nsppol)=1;if (present(mpi_spintab)) my_spintab(1:2)=mpi_spintab(1:2)
740  my_nsppol=count(my_spintab==1)
741 
742  call metric(gmet,gprimd,-1,rmet,rprimd,ucvol)
743 
744  ABI_CHECK(mgfft==MAXVAL(ngfft(1:3)),"Wrong mgfft")
745 
746 !Allocate the arrays of the Hamiltonian whose dimensions do not depend on k
747  ABI_ALLOCATE(ham%atindx,(natom))
748  ABI_ALLOCATE(ham%atindx1,(natom))
749  ABI_ALLOCATE(ham%typat,(natom))
750  ABI_ALLOCATE(ham%indlmn,(6,psps%lmnmax,psps%ntypat))
751  ABI_ALLOCATE(ham%nattyp,(psps%ntypat))
752  ABI_ALLOCATE(ham%nucdipmom,(3,natom))
753  ABI_ALLOCATE(ham%ph1d,(2,3*(2*mgfft+1)*natom))
754  ABI_ALLOCATE(ham%pspso,(psps%ntypat))
755 
756 !Initialize most of the Hamiltonian
757  indx=1
758  do itypat=1,psps%ntypat
759    ham%nattyp(itypat)=0
760    do iat=1,natom
761      if (typat(iat)==itypat) then
762        ham%atindx (iat )=indx
763        ham%atindx1(indx)=iat
764        indx=indx+1
765        ham%nattyp(itypat)=ham%nattyp(itypat)+1
766      end if
767    end do
768  end do
769 
770  ham%gmet(:,:)  =gmet(:,:)
771  ham%gprimd(:,:)=gprimd(:,:)
772  ham%indlmn(:,:,:)=psps%indlmn(:,:,:)
773  ham%lmnmax     =psps%lmnmax
774  ham%mgfft      =mgfft
775  ham%mpsang     =psps%mpsang
776  ham%mpssoang   =psps%mpssoang
777  ham%natom      =natom
778  ham%nfft       =nfft
779  ham%ngfft(:)   =ngfft(:)
780  ham%nloalg(:)  =nloalg(:)
781  ham%matblk=min(NLO_MINCAT,maxval(ham%nattyp)); if (nloalg(2)>0) ham%matblk=natom
782  ham%nsppol     =nsppol
783  ham%nspinor    =nspinor
784  ham%ntypat     =psps%ntypat
785  ham%typat      =typat(1:natom)
786  ham%nvloc=1; if(nspden==4)ham%nvloc=4
787  ham%n4         =ngfft(4)
788  ham%n5         =ngfft(5)
789  ham%n6         =ngfft(6)
790  ham%usepaw     =psps%usepaw
791  ham%ucvol      =ucvol
792  ham%useylm     =psps%useylm
793  ham%use_gpu_cuda=0 ; if(PRESENT(use_gpu_cuda)) ham%use_gpu_cuda=use_gpu_cuda
794 
795  ham%pspso(:)   =psps%pspso(1:psps%ntypat)
796  if (psps%usepaw==1) then
797    do itypat=1,psps%ntypat
798      ham%pspso(itypat)=1+pawtab(itypat)%usespnorb
799    end do
800  end if
801 
802  if (present(nucdipmom)) then
803    ham%nucdipmom(:,:) = nucdipmom(:,:)
804  else
805    ham%nucdipmom(:,:) = zero
806  end if
807 
808  ham%xred => xred
809 
810  if (present(fock)) then
811    if (associated(fock)) then
812      ham%fockcommon => fock%fock_common
813      if (fock%fock_common%use_ACE==0) ham%fockbz=>fock%fock_BZ
814    end if
815  end if
816 
817  if (present(ph1d)) then
818    ham%ph1d(:,:) = ph1d(:,:)
819  else ! Recalculate structure factor phases
820    call getph(ham%atindx,natom,ngfft(1),ngfft(2),ngfft(3),ham%ph1d,xred)
821  end if
822 
823  if (ham%usepaw==1) then
824    ham%usecprj=0;if (present(usecprj)) ham%usecprj=usecprj
825    ABI_ALLOCATE(ham%dimcprj,(natom))
826    !Be carefull cprj are ordered by atom type (used in non-local operator)
827    call pawcprj_getdim(ham%dimcprj,natom,ham%nattyp,ham%ntypat,ham%typat,pawtab,'O')
828  else
829    ham%usecprj=0
830    ABI_ALLOCATE(ham%dimcprj,(0))
831  end if
832 
833 ! ===========================
834 ! ==== Non-local factors ====
835 ! ===========================
836 
837 
838  if (ham%usepaw==0) then ! Norm-conserving: use constant Kleimann-Bylander energies.
839    ham%dimekb1=psps%dimekb
840    ham%dimekb2=psps%ntypat
841    ham%dimekbq=1
842    ABI_ALLOCATE(ham%ekb_spin,(psps%dimekb,psps%ntypat,nspinor**2,1,1))
843    ham%ekb => ham%ekb_spin(:,:,:,:,1)
844    ABI_ALLOCATE(ham%sij,(0,0))
845    ham%ekb(:,:,1,1)=psps%ekb(:,:)
846    if (nspinor==2) then
847      ham%ekb(:,:,2,1)=psps%ekb(:,:)
848      ham%ekb(:,:,3:4,1)=zero
849    end if
850    if (PRESENT(electronpositron)) then
851      if (electronpositron_calctype(electronpositron)==1) ham%ekb(:,:,:,:)=-ham%ekb(:,:,:,:)
852    end if
853 
854 !  Update enl on GPU (will do it later for PAW)
855 #if defined HAVE_GPU_CUDA
856    if (ham%use_gpu_cuda==1) then
857      call gpu_update_ham_data(ham%ekb(:,:,:,1),size(ham%ekb),ham%sij,size(ham%sij),ham%gprimd,size(ham%gprimd))
858    end if
859 #endif
860 
861  else ! PAW: store overlap coefficients (spin non dependent) and Dij coefficients (spin dependent)
862    cplex_dij=1
863    if (present(paw_ij)) then
864      if (size(paw_ij)>0) cplex_dij=paw_ij(1)%cplex_dij
865    end if
866    if ((nspinor==2).or.any(abs(ham%nucdipmom)>tol8)) cplex_dij=2
867    ham%dimekb1=psps%dimekb*cplex_dij
868    ham%dimekb2=natom
869    ham%dimekbq=1
870    if (present(paw_ij)) then
871      if (size(paw_ij)>0) ham%dimekbq=paw_ij(1)%cplex_rf
872    end if
873    ABI_ALLOCATE(ham%sij,(ham%dimekb1,psps%ntypat))
874    do itypat=1,psps%ntypat
875      if (cplex_dij==1) then
876        ham%sij(1:pawtab(itypat)%lmn2_size,itypat)=pawtab(itypat)%sij(:)
877      else
878        do ilmn=1,pawtab(itypat)%lmn2_size
879          ham%sij(2*ilmn-1,itypat)=pawtab(itypat)%sij(ilmn)
880          ham%sij(2*ilmn  ,itypat)=zero
881        end do
882      end if
883      if (cplex_dij*pawtab(itypat)%lmn2_size<ham%dimekb1) then
884        ham%sij(cplex_dij*pawtab(itypat)%lmn2_size+1:ham%dimekb1,itypat)=zero
885      end if
886    end do
887    !We preload here PAW non-local factors in order to avoid a communication over atoms
888    ! inside the loop over spins.
889    ABI_ALLOCATE(ham%ekb_spin,(ham%dimekb1,ham%dimekb2,nspinor**2,ham%dimekbq,my_nsppol))
890    ham%ekb_spin=zero
891    if (present(paw_ij)) then
892      if (my_nsppol<ham%nsppol) then
893        ABI_ALLOCATE(ekb_tmp,(ham%dimekb1,ham%dimekb2,nspinor**2,ham%dimekbq))
894      end if
895      jsp=0
896      do isp=1,ham%nsppol
897        if (my_spintab(isp)==1) then
898          jsp=jsp+1 ; ham%ekb => ham%ekb_spin(:,:,:,:,jsp)
899        else
900          ham%ekb => ekb_tmp
901        end if
902        if (present(mpi_atmtab)) then
903          call pawdij2ekb(ham%ekb,paw_ij,isp,my_comm_atom,mpi_atmtab=mpi_atmtab)
904        else
905          call pawdij2ekb(ham%ekb,paw_ij,isp,my_comm_atom)
906        end if
907      end do
908      if (my_nsppol<ham%nsppol) then
909        ABI_DEALLOCATE(ekb_tmp)
910      end if
911    end if
912    nullify(ham%ekb)
913  end if
914 
915  DBG_EXIT("COLL")
916 
917 end subroutine init_hamiltonian

m_hamiltonian/init_rf_hamiltonian [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  init_rf_hamiltonian

FUNCTION

  Creation method for the rf_hamiltonian_type structure.
  It allocates memory and initializes all quantities that do not depend on the k-point or spin.

INPUTS

  [comm_atom]=optional, MPI communicator over atoms
  cplex_paw=1 if all on-site PAW quantities are real (GS), 2 if they are complex (RF)
  gs_Ham<gs_hamiltonian_type>=Structured datatype containing data for ground-state Hamiltonian at (k+q)
  [has_e1kbsc]=optional, true if rf_Ham%e1kbsc has to be initialized.
               e1kbsc contains the self-consistent 1st-order PAW Dij coefficients (depending on VHxc^(1))
  ipert=index of perturbation
  [mpi_atmtab(:)]=optional, indexes of the atoms treated by current proc
  [mpi_spintab(2)]=optional, flags defining the spin(s) treated be current process:
                    mpi_spintab(1)=1 if non-polarized or spin-up treated
                    mpi_spintab(2)=1 if polarized and spin-dn treated
  [paw_ij1(:)<paw_ij_type>]=Various 1st-order arrays given on (i,j) (partial waves)
                            channels (paw_ij1%dij and paw_ij1%difr only used here).

SIDE EFFECTS

  rf_Ham<rf_hamiltonian_type>=Structured datatype almost completely initialized:
   * Basic variables and dimensions are transfered to the structure.
   * All pointers are allocated with correct dimensions.
   * Quantities that do not depend on the k-point or spin are initialized.

PARENTS

      dfpt_nstpaw,dfpt_nstwf,dfpt_rhofermi,dfpt_vtorho,m_gkk,m_phgamma,m_phpi
      m_sigmaph

CHILDREN

      destroy_mpi_enreg,initmpi_seq,kpgsph,wrtout

SOURCE

1643 subroutine init_rf_hamiltonian(cplex,gs_Ham,ipert,rf_Ham,&
1644 &          comm_atom,mpi_atmtab,mpi_spintab,paw_ij1,has_e1kbsc) ! optional arguments
1645 
1646 
1647 !This section has been created automatically by the script Abilint (TD).
1648 !Do not modify the following lines by hand.
1649 #undef ABI_FUNC
1650 #define ABI_FUNC 'init_rf_hamiltonian'
1651 !End of the abilint section
1652 
1653  implicit none
1654 
1655 !Arguments ------------------------------------
1656 !scalars
1657  integer,intent(in) :: cplex,ipert
1658  integer,intent(in),optional :: comm_atom
1659  logical,intent(in),optional :: has_e1kbsc
1660  type(gs_hamiltonian_type),intent(in) :: gs_Ham
1661  type(rf_hamiltonian_type),intent(inout),target :: rf_Ham
1662 !arrays
1663  integer,optional,intent(in)  :: mpi_atmtab(:),mpi_spintab(2)
1664  type(paw_ij_type),optional,intent(in) :: paw_ij1(:)
1665 
1666 !Local variables-------------------------------
1667 !scalars
1668  integer :: cplex_dij1,isp,jsp,my_comm_atom,my_nsppol
1669  logical :: has_e1kbsc_
1670 !arrays
1671  integer :: my_spintab(2)
1672  real(dp),allocatable,target :: e1kb_tmp(:,:,:,:)
1673 
1674 ! *************************************************************************
1675 
1676  DBG_ENTER("COLL")
1677 
1678 !@rf_hamiltonian_type
1679 
1680 !Manage optional parameters
1681  has_e1kbsc_=.false.;if (present(has_e1kbsc)) has_e1kbsc_=has_e1kbsc
1682  my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
1683  my_spintab=0;my_spintab(1:gs_Ham%nsppol)=1;if(present(mpi_spintab)) my_spintab=mpi_spintab
1684  my_nsppol=count(my_spintab==1)
1685 
1686  rf_Ham%cplex    =cplex
1687 
1688  rf_Ham%n4       =gs_Ham%n4
1689  rf_Ham%n5       =gs_Ham%n5
1690  rf_Ham%n6       =gs_Ham%n6
1691  rf_Ham%nvloc    =gs_Ham%nvloc
1692  rf_Ham%nsppol   =gs_Ham%nsppol
1693  rf_Ham%nspinor  =gs_Ham%nspinor
1694 
1695  rf_Ham%dime1kb1=0
1696  rf_Ham%dime1kb2=gs_Ham%dimekb2
1697  if (gs_Ham%usepaw==1.and.ipert/=gs_Ham%natom+1.and.ipert/=gs_Ham%natom+10) then
1698    cplex_dij1=1;if ((gs_Ham%nspinor==2).or.any(abs(gs_Ham%nucdipmom)>tol8)) cplex_dij1=2
1699    rf_Ham%dime1kb1=cplex_dij1*(gs_Ham%lmnmax*(gs_Ham%lmnmax+1))/2
1700  end if
1701 
1702 !Allocate the arrays of the 1st-order Hamiltonian
1703 !  We preload here 1st-order non-local factors in order to avoid
1704 !  a communication over atoms inside the loop over spins.
1705  if (gs_Ham%usepaw==1.and.rf_Ham%dime1kb1>0) then
1706    if ((ipert>=1.and.ipert<=gs_Ham%natom).or.ipert==gs_Ham%natom+2.or.&
1707 &    ipert==gs_Ham%natom+3.or.ipert==gs_Ham%natom+4.or.ipert==gs_Ham%natom+11) then
1708 
1709      ABI_ALLOCATE(rf_Ham%e1kbfr_spin,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex,my_nsppol))
1710      rf_Ham%e1kbfr_spin=zero
1711      if (has_e1kbsc_) then
1712        ABI_ALLOCATE(rf_Ham%e1kbsc_spin,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex,my_nsppol))
1713        rf_Ham%e1kbsc_spin=zero
1714      end if
1715 
1716      if (present(paw_ij1)) then
1717 
1718        if (my_nsppol<rf_Ham%nsppol) then
1719          ABI_ALLOCATE(e1kb_tmp,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex))
1720        end if
1721 
1722 !      === Frozen term
1723        jsp=0
1724        do isp=1,rf_Ham%nsppol
1725          if (my_spintab(isp)==1) then
1726            jsp=jsp+1 ; rf_Ham%e1kbfr => rf_Ham%e1kbfr_spin(:,:,:,:,jsp)
1727          else
1728            rf_Ham%e1kbfr => e1kb_tmp
1729          end if
1730          if (present(mpi_atmtab)) then
1731            call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbfr=rf_Ham%e1kbfr,mpi_atmtab=mpi_atmtab)
1732          else
1733            call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbfr=rf_Ham%e1kbfr)
1734          end if
1735        end do
1736 
1737 !      === Self-consistent term
1738        if (has_e1kbsc_) then
1739          jsp=0
1740          do isp=1,rf_Ham%nsppol
1741            if (my_spintab(isp)==1) then
1742              jsp=jsp+1 ; rf_Ham%e1kbsc => rf_Ham%e1kbsc_spin(:,:,:,:,jsp)
1743            else
1744              rf_Ham%e1kbsc => e1kb_tmp
1745            end if
1746            if (present(mpi_atmtab)) then
1747              call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbsc=rf_Ham%e1kbsc,mpi_atmtab=mpi_atmtab)
1748            else
1749              call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbsc=rf_Ham%e1kbsc)
1750            end if
1751          end do
1752        end if
1753 
1754        if (my_nsppol<rf_Ham%nsppol) then
1755          ABI_DEALLOCATE(e1kb_tmp)
1756        end if
1757 
1758      end if
1759    end if
1760  end if
1761 
1762  if (.not.allocated(rf_Ham%e1kbfr_spin)) then
1763    ABI_ALLOCATE(rf_Ham%e1kbfr_spin,(0,0,0,0,0))
1764  end if
1765  if (.not.allocated(rf_Ham%e1kbsc_spin)) then
1766    ABI_ALLOCATE(rf_Ham%e1kbsc_spin,(0,0,0,0,0))
1767  end if
1768  nullify(rf_Ham%e1kbfr)
1769  nullify(rf_Ham%e1kbsc)
1770 
1771  DBG_EXIT("COLL")
1772 
1773 end subroutine init_rf_hamiltonian

m_hamiltonian/load_k_hamiltonian [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  load_k_hamiltonian

FUNCTION

  Setup of the k-dependent part of the Hamiltonian H_k_k^prime

INPUTS

  [compute_gbound]=flag. if true, G sphere boundary is computed here
  [compute_ph3d]=flag. if true, 3D structure factors are computed here (only if nloalg(1)>0)
  [gbound_k]=G sphere boundary (not compatible with compute_gbound=TRUE)
  [ffnl_k]=nonlocal form factors on basis sphere
  [istwf_k]=parameter that describes the storage of wfs
  [kinpw_k]=(modified) kinetic energy for each plane wave
  [kg_k]=planewave reduced coordinates in basis sphere (g vectors)
  [kpg_k]=(k+g) vectors in reciprocal space
  [kpt_k]=k point coordinates
  [npw_k]=number of plane waves (processed by current proc when band-FFT parallelism is on)
  [npw_fft_k]=number of plane waves used to apply Hamiltonian (in the "FFT" configuration)
  [ph3d_k]=3-dim structure factors, for each atom and plane wave

SIDE EFFECTS

  ham<gs_hamiltonian_type>=structured datatype completed with k-dependent quantitites.
          Quantities at k^prime are set equal to quantities at k.
    k-dependent scalars and pointers associated
    phkxred=exp(.k.xred) for each atom
    [ham%gbound_k]=G sphere boundary, for each plane wave
    [ham%ph3d_k]=3-dim structure factors, for each atom and plane wave

PARENTS

      d2frnl,dfpt_nsteltwf,dfpt_nstpaw,dfpt_nstwf,dfpt_rhofermi,dfptnl_resp
      energy,fock2ACE,forstrnps,getgh1c,gwls_hamiltonian,ks_ddiago,m_io_kss
      m_shirley,nonlop_test,vtorho,wfd_vnlpsi

CHILDREN

      destroy_mpi_enreg,initmpi_seq,kpgsph,wrtout

SOURCE

 961 subroutine load_k_hamiltonian(ham,ffnl_k,fockACE_k,gbound_k,istwf_k,kinpw_k,&
 962 &                             kg_k,kpg_k,kpt_k,nucdipmom_k,npw_k,npw_fft_k,ph3d_k,&
 963 &                             compute_gbound,compute_ph3d)
 964 
 965 
 966 !This section has been created automatically by the script Abilint (TD).
 967 !Do not modify the following lines by hand.
 968 #undef ABI_FUNC
 969 #define ABI_FUNC 'load_k_hamiltonian'
 970 !End of the abilint section
 971 
 972  implicit none
 973 
 974 !Arguments ------------------------------------
 975 !scalars
 976  integer,intent(in),optional :: npw_k,npw_fft_k,istwf_k
 977  logical,intent(in),optional :: compute_gbound,compute_ph3d
 978  type(gs_hamiltonian_type),intent(inout),target :: ham
 979 !arrays
 980  integer,intent(in),optional,target :: gbound_k(:,:),kg_k(:,:)
 981  real(dp),intent(in),optional :: kpt_k(3)
 982  real(dp),intent(in),optional,target :: ffnl_k(:,:,:,:),kinpw_k(:),kpg_k(:,:),ph3d_k(:,:,:)
 983  complex(dpc),intent(in),optional :: nucdipmom_k(:)
 984  type(fock_ACE_type),intent(in),optional,target :: fockACE_k
 985 
 986 !Local variables-------------------------------
 987 !scalars
 988  integer :: iat,iatom
 989  logical :: compute_gbound_
 990  real(dp) :: arg
 991  !character(len=500) :: msg
 992 
 993 ! *************************************************************************
 994 
 995  DBG_ENTER("COLL")
 996 
 997 !@gs_hamiltonian_type
 998 
 999 !k-dependent scalars
1000  if (present(kpt_k)) then
1001    ham%kpt_k(:)  = kpt_k(:)
1002    ham%kpt_kp(:) = kpt_k(:)
1003  end if
1004  if (present(istwf_k)) then
1005    ham%istwf_k  = istwf_k
1006    ham%istwf_kp = istwf_k
1007  end if
1008  if (present(npw_k)) then
1009    ham%npw_k  = npw_k
1010    ham%npw_kp = npw_k
1011  end if
1012  if (present(npw_fft_k)) then
1013    ham%npw_fft_k  = npw_fft_k
1014    ham%npw_fft_kp = npw_fft_k
1015  else if (present(npw_k)) then
1016    ham%npw_fft_k  = npw_k
1017    ham%npw_fft_kp = npw_k
1018  end if
1019 
1020  ! k-dependend complex quantities
1021   if (present(nucdipmom_k)) then
1022    ham%nucdipmom_k(:) = nucdipmom_k(:)
1023   end if
1024 
1025 !Pointers to k-dependent quantitites
1026  if (present(kinpw_k)) then
1027    ham%kinpw_k  => kinpw_k
1028    ham%kinpw_kp => kinpw_k
1029  end if
1030  if (present(kg_k)) then
1031    ham%kg_k  => kg_k
1032    ham%kg_kp => kg_k
1033  end if
1034  if (present(kpg_k)) then
1035    ham%kpg_k  => kpg_k
1036    ham%kpg_kp => kpg_k
1037  end if
1038  if (present(ffnl_k)) then
1039    ham%ffnl_k  => ffnl_k
1040    ham%ffnl_kp => ffnl_k
1041  end if
1042  if (present(ph3d_k)) then
1043    ham%ph3d_k  => ph3d_k
1044    ham%ph3d_kp => ph3d_k
1045  end if
1046  if (present(fockACE_k)) then
1047    ham%fockACE_k  => fockACE_k
1048  end if
1049 !Compute exp(i.k.R) for each atom
1050  if (present(kpt_k)) then
1051    if (associated(Ham%phkpxred).and.(.not.associated(Ham%phkpxred,Ham%phkxred))) then
1052      ABI_DEALLOCATE(Ham%phkpxred)
1053    end if
1054    if (allocated(ham%phkxred)) then
1055      ABI_DEALLOCATE(ham%phkxred)
1056    end if
1057    ABI_ALLOCATE(ham%phkxred,(2,ham%natom))
1058    do iat=1,ham%natom
1059      iatom=ham%atindx(iat)
1060      arg=two_pi*DOT_PRODUCT(kpt_k,ham%xred(:,iat))
1061      ham%phkxred(1,iatom)=DCOS(arg)
1062      ham%phkxred(2,iatom)=DSIN(arg)
1063    end do
1064    ham%phkpxred => ham%phkxred
1065  end if
1066 
1067 !Compute or copy G sphere boundary at k+g
1068  compute_gbound_=.false.;if (present(compute_gbound)) compute_gbound_=compute_gbound
1069  if (present(gbound_k)) compute_gbound_=.true.
1070  if (compute_gbound_) then
1071    if (associated(Ham%gbound_kp,Ham%gbound_k)) then
1072      nullify(Ham%gbound_kp)
1073    else if (associated(Ham%gbound_kp)) then
1074      ABI_DEALLOCATE(Ham%gbound_kp)
1075    end if
1076    if (allocated(ham%gbound_k)) then
1077      ABI_DEALLOCATE(ham%gbound_k)
1078    end if
1079  end if
1080  if (.not.allocated(ham%gbound_k)) then
1081    ABI_ALLOCATE(ham%gbound_k,(2*ham%mgfft+8,2))
1082    ham%gbound_k(:,:)=0
1083    ham%gbound_kp => ham%gbound_k
1084  end if
1085  if (compute_gbound_) then
1086    if (present(gbound_k)) then
1087      ham%gbound_k(:,:)=gbound_k(:,:)
1088    else
1089      if (.not.associated(ham%kg_k)) then
1090        MSG_BUG('Something is missing for gbound_k computation!')
1091      end if
1092      !write(std_out,*)"About to call sphereboundary"
1093      !write(std_out,*)"size(kg_k), npw_k, mgfft",size(ham%kg_k, dim=2), ham%npw_k, ham%mgfft
1094      call sphereboundary(ham%gbound_k,ham%istwf_k,ham%kg_k,ham%mgfft,ham%npw_k)
1095    end if
1096    ham%gbound_kp => ham%gbound_k
1097  end if
1098 
1099 !Compute 3D structure factors for each atom at k+g
1100  if (present(compute_ph3d).and.present(ph3d_k)) then
1101    if (compute_ph3d.and.ham%nloalg(2)>0) then
1102      if ((.not.allocated(ham%phkxred)).or.(.not.associated(ham%kg_k)).or.&
1103 &        (.not.associated(ham%ph3d_k))) then
1104        MSG_BUG('Something is missing for ph3d_k computation!')
1105      end if
1106      call ph1d3d(1,ham%natom,ham%kg_k,ham%matblk,ham%natom,ham%npw_k,ham%ngfft(1),&
1107 &                ham%ngfft(2),ham%ngfft(3),ham%phkxred,ham%ph1d,ham%ph3d_k)
1108    end if
1109  end if
1110 
1111  DBG_EXIT("COLL")
1112 
1113 end subroutine load_k_hamiltonian

m_hamiltonian/load_k_rf_hamiltonian [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  load_k_rf_hamiltonian

FUNCTION

  Setup of the k-dependent part of the 1st- and 2nd- order Hamiltonian

INPUTS

  [dkinpw_k]=1st derivative of the (modified) kinetic energy for each plane wave
  [ddkinpw_k]=2nd derivative of the (modified) kinetic energy for each plane wave
  [npw_k]=number of plane waves

SIDE EFFECTS

  rf_Ham<rf_hamiltonian_type>=structured datatype completed with k-dependent quantitites.
          Quantities at k^prime are set equal to quantities at k.

PARENTS

      dfpt_nstpaw,dfpt_nstwf,dfpt_rhofermi,getgh1c

CHILDREN

      destroy_mpi_enreg,initmpi_seq,kpgsph,wrtout

SOURCE

1881 subroutine load_k_rf_hamiltonian(rf_Ham,dkinpw_k,ddkinpw_k,npw_k)
1882 
1883 
1884 !This section has been created automatically by the script Abilint (TD).
1885 !Do not modify the following lines by hand.
1886 #undef ABI_FUNC
1887 #define ABI_FUNC 'load_k_rf_hamiltonian'
1888 !End of the abilint section
1889 
1890  implicit none
1891 
1892 !Arguments ------------------------------------
1893 !scalars
1894  integer,intent(in),optional :: npw_k
1895  type(rf_hamiltonian_type),intent(inout),target :: rf_Ham
1896 !arrays
1897  real(dp),intent(in),optional,target :: dkinpw_k(:),ddkinpw_k(:)
1898 
1899 !Local variables-------------------------------
1900 
1901 ! *************************************************************************
1902 
1903  DBG_ENTER("COLL")
1904 
1905 !@gs_hamiltonian_type
1906 
1907 !k-dependent scalars
1908  if (present(npw_k)) then
1909    rf_Ham%npw_k  = npw_k
1910    rf_Ham%npw_kp = npw_k
1911  end if
1912 
1913 !Pointers to k-dependent quantitites
1914  if (present(dkinpw_k)) then
1915    rf_Ham%dkinpw_k  => dkinpw_k
1916    rf_Ham%dkinpw_kp => dkinpw_k
1917  end if
1918  if (present(ddkinpw_k)) then
1919    rf_Ham%ddkinpw_k  => ddkinpw_k
1920    rf_Ham%ddkinpw_kp => ddkinpw_k
1921  end if
1922 
1923  DBG_EXIT("COLL")
1924 
1925 end subroutine load_k_rf_hamiltonian

m_hamiltonian/load_kprime_hamiltonian [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  load_kprime_hamiltonian

FUNCTION

  Setup of the k^prime-dependent part of the Hamiltonian H_k_k^prime

INPUTS

  [compute_gbound]=flag. if true, G sphere boundary is computed here
  [compute_ph3d]=flag. if true, 3D structure factors are computed here (only if nloalg(2)>0)
  [gbound_kp]=G sphere boundary (not compatible with compute_gbound=TRUE)
  [ffnl_kp]=nonlocal form factors on basis sphere
  [istwf_kp]=parameter that describes the storage of wfs
  [kinpw_kp]=(modified) kinetic energy for each plane wave
  [kg_kp]=planewave reduced coordinates in basis sphere (g vectors)
  [kpg_kp]=(k+g) vectors in reciprocal space
  [kpt_kp]=k point coordinates
  [npw_kp]=number of plane waves (processed by current proc when band-FFT parallelism is on)
  [npw_fft_kp]=number of plane waves used to apply Hamiltonian (in the "FFT" configuration)
  [ph3d_kp]=3-dim structure factors, for each atom and plane wave

SIDE EFFECTS

  ham<gs_hamiltonian_type>=structured datatype completed with k^prime-dependent quantitites.
    k^prime-dependent scalars and pointers associated
    phkpxred=exp(.k^prime.xred) for each atom
    [ham%gbound_kp]=G sphere boundary, for each plane wave
    [ham%ph3d_kp]=3-dim structure factors at k^prime at k, for each atom and plane wave

PARENTS

      dfpt_nstpaw,dfpt_nstwf,dfpt_rhofermi,fock_getghc,getgh1c

CHILDREN

      destroy_mpi_enreg,initmpi_seq,kpgsph,wrtout

SOURCE

1154 subroutine load_kprime_hamiltonian(ham,ffnl_kp,gbound_kp,istwf_kp,kinpw_kp,&
1155 &                                  kg_kp,kpg_kp,kpt_kp,npw_kp,npw_fft_kp,&
1156 &                                  ph3d_kp,compute_gbound,compute_ph3d)
1157 
1158 
1159 !This section has been created automatically by the script Abilint (TD).
1160 !Do not modify the following lines by hand.
1161 #undef ABI_FUNC
1162 #define ABI_FUNC 'load_kprime_hamiltonian'
1163 !End of the abilint section
1164 
1165  implicit none
1166 
1167 !Arguments ------------------------------------
1168 !scalars
1169  integer,intent(in),optional :: npw_kp,npw_fft_kp,istwf_kp
1170  logical,intent(in),optional :: compute_gbound,compute_ph3d
1171  type(gs_hamiltonian_type),intent(inout),target :: ham
1172 !arrays
1173  integer,intent(in),optional,target :: gbound_kp(:,:),kg_kp(:,:)
1174  real(dp),intent(in),optional :: kpt_kp(3)
1175  real(dp),intent(in),optional,target :: ffnl_kp(:,:,:,:),kinpw_kp(:),kpg_kp(:,:),ph3d_kp(:,:,:)
1176 
1177 !Local variables-------------------------------
1178 !scalars
1179  integer :: iat,iatom
1180  logical :: compute_gbound_
1181  real(dp) :: arg
1182  character(len=100) :: msg
1183 
1184 ! *************************************************************************
1185 
1186  DBG_ENTER("COLL")
1187 
1188 !@gs_hamiltonian_type
1189 
1190 !k-dependent scalars
1191  if (present(kpt_kp))   ham%kpt_kp(:)= kpt_kp(:)
1192  if (present(istwf_kp)) ham%istwf_kp = istwf_kp
1193  if (present(npw_kp))   ham%npw_kp   = npw_kp
1194  if (present(npw_fft_kp)) then
1195     ham%npw_fft_kp = npw_fft_kp
1196  else if (present(npw_kp)) then
1197     ham%npw_fft_kp = npw_kp
1198  end if
1199 
1200 !Pointers to k-dependent quantitites
1201  if (present(kinpw_kp)) ham%kinpw_kp => kinpw_kp
1202  if (present(kg_kp))    ham%kg_kp    => kg_kp
1203  if (present(kpg_kp))   ham%kpg_kp   => kpg_kp
1204  if (present(ffnl_kp))  ham%ffnl_kp  => ffnl_kp
1205  if (present(ph3d_kp))  ham%ph3d_kp  => ph3d_kp
1206 
1207 !Compute exp(i.k^prime.R) for each atom
1208  if (present(kpt_kp)) then
1209    if (associated(ham%phkpxred,ham%phkxred)) then
1210      nullify(ham%phkpxred)
1211    else if (associated(ham%phkpxred)) then
1212      ABI_DEALLOCATE(ham%phkpxred)
1213    end if
1214    ABI_ALLOCATE(ham%phkpxred,(2,ham%natom))
1215    do iat=1,ham%natom
1216      iatom=ham%atindx(iat)
1217      arg=two_pi*DOT_PRODUCT(kpt_kp,ham%xred(:,iat))
1218      ham%phkpxred(1,iatom)=DCOS(arg)
1219      ham%phkpxred(2,iatom)=DSIN(arg)
1220    end do
1221  end if
1222 
1223 !Compute or copy G sphere boundary at k^prime+g
1224  compute_gbound_=.false.
1225  if (present(kpt_kp).and.present(compute_gbound)) compute_gbound_=compute_gbound
1226  if (present(gbound_kp)) compute_gbound_=.true.
1227  if (compute_gbound_) then
1228    if (associated(ham%gbound_kp,ham%gbound_k)) then
1229      nullify(ham%gbound_kp)
1230    else if (associated(ham%gbound_kp)) then
1231      ABI_DEALLOCATE(ham%gbound_kp)
1232    end if
1233    if (present(gbound_kp)) then
1234      ham%gbound_kp(:,:)=gbound_kp(:,:)
1235    else
1236      if (.not.associated(ham%kg_kp)) then
1237        msg='Something is missing for gbound_kp computation!'
1238        MSG_BUG(msg)
1239      end if
1240      ABI_ALLOCATE(ham%gbound_kp,(2*ham%mgfft+8,2))
1241      call sphereboundary(ham%gbound_kp,ham%istwf_kp,ham%kg_kp,ham%mgfft,ham%npw_kp)
1242    end if
1243  end if
1244 
1245 !Compute 3D structure factors for each atom at k^prime+g
1246  if (present(compute_ph3d).and.present(ph3d_kp)) then
1247    if (compute_ph3d.and.ham%nloalg(2)>0) then
1248      if ((.not.associated(ham%phkpxred)).or.(.not.associated(ham%kg_kp)).or.&
1249 &        (.not.associated(ham%ph3d_kp))) then
1250        msg='Something is missing for ph3d_kp computation!'
1251        MSG_BUG(msg)
1252      end if
1253      call ph1d3d(1,ham%natom,ham%kg_kp,ham%matblk,ham%natom,ham%npw_kp,ham%ngfft(1),&
1254 &                ham%ngfft(2),ham%ngfft(3),ham%phkpxred,ham%ph1d,ham%ph3d_kp)
1255    end if
1256  end if
1257 
1258  DBG_EXIT("COLL")
1259 
1260 end subroutine load_kprime_hamiltonian

m_hamiltonian/load_spin_hamiltonian [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  load_spin_hamiltonian

INPUTS

  isppol=index of current spin
  [vlocal(n4,n5,n6,nvloc)]=optional, local potential in real space
  [vxctaulocal(n4,n5,n6,nvloc,4)]=optional, derivative of XC energy density with respect
                                  to kinetic energy density in real space
  [with_nonlocal]=optional, true if non-local factors have to be loaded

FUNCTION

  Setup of the spin-dependent part of the GS Hamiltonian.

SIDE EFFECTS

  Ham<gs_hamiltonian_type>=Structured datatype initialization phase:
   * Quantities that depend spin are initialized.

PARENTS

      d2frnl,dfpt_nstdy,dfpt_nstpaw,dfpt_rhofermi,dfpt_vtorho,energy,fock2ACE
      forstrnps,ks_ddiago,m_gkk,m_io_kss,m_phgamma,m_phpi,m_shirley,m_sigmaph
      nonlop_test,vtorho

CHILDREN

      destroy_mpi_enreg,initmpi_seq,kpgsph,wrtout

SOURCE

1484 subroutine load_spin_hamiltonian(Ham,isppol,vlocal,vxctaulocal,with_nonlocal)
1485 
1486 
1487 !This section has been created automatically by the script Abilint (TD).
1488 !Do not modify the following lines by hand.
1489 #undef ABI_FUNC
1490 #define ABI_FUNC 'load_spin_hamiltonian'
1491 !End of the abilint section
1492 
1493  implicit none
1494 
1495 !Arguments ------------------------------------
1496 !scalars
1497  integer,intent(in) :: isppol
1498  logical,optional,intent(in) :: with_nonlocal
1499  type(gs_hamiltonian_type),intent(inout),target :: Ham
1500 !arrays
1501  real(dp),optional,intent(in),target :: vlocal(:,:,:,:),vxctaulocal(:,:,:,:,:)
1502 
1503 !Local variables-------------------------------
1504 !scalars
1505  integer :: jsppol
1506 
1507 ! *************************************************************************
1508 
1509  DBG_ENTER("COLL")
1510 
1511 !@gs_hamiltonian_type
1512 
1513  if (present(vlocal)) then
1514    ABI_CHECK(size(vlocal)==Ham%n4*Ham%n5*Ham%n6*Ham%nvloc,"Wrong vlocal")
1515    Ham%vlocal => vlocal
1516  end if
1517  if (present(vxctaulocal)) then
1518    ABI_CHECK(size(vxctaulocal)==Ham%n4*Ham%n5*Ham%n6*Ham%nvloc*4,"Wrong vxctaulocal")
1519    Ham%vxctaulocal => vxctaulocal
1520  end if
1521 
1522 !Retrieve non-local factors for this spin component
1523  if (present(with_nonlocal)) then
1524    if (with_nonlocal) then
1525      jsppol=min(isppol,size(Ham%ekb_spin,5))
1526      if (jsppol>0) Ham%ekb => Ham%ekb_spin(:,:,:,:,jsppol)
1527    end if
1528  end if
1529 
1530 !Update enl and sij on GPU
1531 #if defined HAVE_GPU_CUDA
1532  if (Ham%use_gpu_cuda==1) then
1533    call gpu_update_ham_data(Ham%ekb(:,:,:,1),size(Ham%ekb),Ham%sij,size(Ham%sij),Ham%gprimd,size(Ham%gprimd))
1534  end if
1535 #endif
1536 
1537  DBG_EXIT("COLL")
1538 
1539 end subroutine load_spin_hamiltonian

m_hamiltonian/load_spin_rf_hamiltonian [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  load_spin_rf_hamiltonian

FUNCTION

  Setup of the spin-dependent part of the 1st- and 2nd- order Hamiltonian.

INPUTS

  isppol=index of current spin
  [vlocal1(cplex*n4,n5,n6,nvloc)]=optional, 1st-order local potential in real space
  [with_nonlocal]=optional, true if non-local factors have to be loaded

SIDE EFFECTS

  rf_Ham<rf_hamiltonian_type>=Structured datatype initialization phase:
   * Quantities that depend on spin are initialized.

PARENTS

      dfpt_rhofermi,dfpt_vtorho,m_gkk,m_phgamma,m_phpi,m_sigmaph

CHILDREN

      destroy_mpi_enreg,initmpi_seq,kpgsph,wrtout

SOURCE

1802 subroutine load_spin_rf_hamiltonian(rf_Ham,isppol,vlocal1,with_nonlocal)
1803 
1804 
1805 !This section has been created automatically by the script Abilint (TD).
1806 !Do not modify the following lines by hand.
1807 #undef ABI_FUNC
1808 #define ABI_FUNC 'load_spin_rf_hamiltonian'
1809 !End of the abilint section
1810 
1811  implicit none
1812 
1813 !Arguments ------------------------------------
1814 !scalars
1815  integer,intent(in) :: isppol
1816  logical,optional,intent(in) :: with_nonlocal
1817  type(rf_hamiltonian_type),intent(inout),target :: rf_Ham
1818 !arrays
1819  real(dp),optional,target,intent(in) :: vlocal1(:,:,:,:)
1820 
1821 !Local variables-------------------------------
1822 !scalars
1823  integer :: jsppol
1824 
1825 ! *************************************************************************
1826 
1827  DBG_ENTER("COLL")
1828 
1829 !@rf_hamiltonian_type
1830 
1831  if (present(vlocal1)) then
1832    ABI_CHECK(size(vlocal1)==rf_Ham%cplex*rf_Ham%n4*rf_Ham%n5*rf_Ham%n6*rf_Ham%nvloc,"Wrong vlocal1")
1833    rf_Ham%vlocal1 => vlocal1
1834  end if
1835 
1836 !Retrieve non-local factors for this spin component
1837  if (present(with_nonlocal)) then
1838    if (with_nonlocal) then
1839      if (size(rf_Ham%e1kbfr_spin)>0) then
1840        jsppol=min(isppol,size(rf_Ham%e1kbfr_spin,5))
1841        if (jsppol>0) rf_Ham%e1kbfr => rf_Ham%e1kbfr_spin(:,:,:,:,jsppol)
1842      end if
1843      if (size(rf_Ham%e1kbsc_spin)>0) then
1844        jsppol=min(isppol,size(rf_Ham%e1kbsc_spin,5))
1845        if (jsppol>0) rf_Ham%e1kbsc => rf_Ham%e1kbsc_spin(:,:,:,:,jsppol)
1846      end if
1847    end if
1848  end if
1849 
1850  DBG_EXIT("COLL")
1851 
1852 end subroutine load_spin_rf_hamiltonian

m_hamiltonian/pawdij2e1kb [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  pawdij2e1kb

FUNCTION

  Transfer PAW Dij (on-site RF Hamiltonian) values
  from paw_ij datastructure to e1kb array

INPUTS

OUTPUT

PARENTS

      d2frnl,dfpt_nstpaw,m_hamiltonian

CHILDREN

      destroy_mpi_enreg,initmpi_seq,kpgsph,wrtout

SOURCE

2043 subroutine pawdij2e1kb(paw_ij1,isppol,comm_atom,mpi_atmtab,e1kbfr,e1kbsc)
2044 
2045 
2046 !This section has been created automatically by the script Abilint (TD).
2047 !Do not modify the following lines by hand.
2048 #undef ABI_FUNC
2049 #define ABI_FUNC 'pawdij2e1kb'
2050 !End of the abilint section
2051 
2052  implicit none
2053 
2054 !Arguments ------------------------------------
2055 !scalars
2056  integer,intent(in) :: isppol,comm_atom
2057 !arrays
2058  integer,intent(in),optional,target :: mpi_atmtab(:)
2059  real(dp),optional,intent(out) :: e1kbfr(:,:,:,:),e1kbsc(:,:,:,:)
2060  type(paw_ij_type),intent(in) :: paw_ij1(:)
2061 
2062 !Local variables-------------------------------
2063 !scalars
2064  integer :: cplex_rf,dimdij1,dime1kb1,dime1kb3,dime1kb4,iatom,iatom_tot,ierr,isp,ispden
2065  integer :: my_natom,natom
2066  logical :: my_atmtab_allocated,paral_atom
2067 !arrays
2068  integer,pointer :: my_atmtab(:)
2069 
2070 ! *************************************************************************
2071 
2072  DBG_ENTER("COLL")
2073 
2074  if ((.not.present(e1kbfr)).and.(.not.present(e1kbsc))) return
2075  if (present(e1kbfr)) then
2076    e1kbfr=zero ; natom=size(e1kbfr,2)
2077  end if
2078  if (present(e1kbsc)) then
2079    e1kbsc=zero ; natom=size(e1kbsc,2)
2080  end if
2081 
2082 !Set up parallelism over atoms
2083  my_natom=size(paw_ij1) ; paral_atom=(xmpi_comm_size(comm_atom)>1)
2084  nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
2085  call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
2086 
2087 !Retrieve 1st-order PAW Dij coefficients for this spin component (frozen)
2088  if (my_natom>0.and.present(e1kbfr)) then
2089    if (allocated(paw_ij1(1)%dijfr)) then
2090      dime1kb1=size(e1kbfr,1) ; dime1kb3=size(e1kbfr,3) ; dime1kb4=size(e1kbfr,4)
2091      ABI_CHECK(paw_ij1(1)%cplex_rf==dime1kb4,'BUG in pawdij2e1kb (1)!')
2092      do ispden=1,dime1kb3
2093        isp=isppol;if (dime1kb3==4) isp=ispden
2094        do iatom=1,my_natom
2095          iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom)
2096          cplex_rf=paw_ij1(iatom)%cplex_rf
2097          dimdij1=paw_ij1(iatom)%cplex_dij*paw_ij1(iatom)%lmn2_size
2098          ABI_CHECK(dimdij1<=dime1kb1,'BUG: size of paw_ij1%dij>dime1kb1!')
2099          e1kbfr(1:dimdij1,iatom_tot,ispden,1)=paw_ij1(iatom)%dijfr(1:dimdij1,isp)
2100          if (cplex_rf==2) e1kbfr(1:dimdij1,iatom_tot,ispden,2)=paw_ij1(iatom)%dijfr(dimdij1+1:2*dimdij1,isp)
2101        end do
2102      end do
2103    end if
2104  end if
2105 
2106 !Retrieve 1st-order PAW Dij coefficients for this spin component (self-consistent)
2107  if (my_natom>0.and.present(e1kbsc)) then
2108    if (allocated(paw_ij1(1)%dijfr).and.allocated(paw_ij1(1)%dij)) then
2109      dime1kb1=size(e1kbsc,1) ; dime1kb3=size(e1kbsc,3) ; dime1kb4=size(e1kbsc,4)
2110      ABI_CHECK(paw_ij1(1)%cplex_rf==dime1kb4,'BUG in pawdij2e1kb (1)!')
2111      do ispden=1,dime1kb3
2112        isp=isppol;if (dime1kb3==4) isp=ispden
2113        do iatom=1,my_natom
2114          iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom)
2115          cplex_rf=paw_ij1(iatom)%cplex_rf
2116          dimdij1=paw_ij1(iatom)%cplex_dij*paw_ij1(iatom)%lmn2_size
2117          ABI_CHECK(dimdij1<=dime1kb1,'BUG: size of paw_ij1%dij>dime1kb1!')
2118          e1kbsc(1:dimdij1,iatom_tot,ispden,1)=paw_ij1(iatom)%dij  (1:dimdij1,isp) &
2119 &                                            -paw_ij1(iatom)%dijfr(1:dimdij1,isp)
2120          if (cplex_rf==2) e1kbsc(1:dimdij1,iatom_tot,ispden,2)=paw_ij1(iatom)%dij  (dimdij1+1:2*dimdij1,isp) &
2121 &                                                             -paw_ij1(iatom)%dijfr(dimdij1+1:2*dimdij1,isp)
2122        end do
2123      end do
2124    end if
2125  end if
2126 
2127 !Communication in case of distribution over atomic sites
2128  if (paral_atom) then
2129    if (present(e1kbfr)) then
2130      call xmpi_sum(e1kbfr,comm_atom,ierr)
2131    end if
2132    if (present(e1kbsc)) then
2133      call xmpi_sum(e1kbsc,comm_atom,ierr)
2134    end if
2135  end if
2136 
2137 !Destroy atom table used for parallelism
2138  call free_my_atmtab(my_atmtab,my_atmtab_allocated)
2139 
2140  DBG_EXIT("COLL")
2141 
2142 end subroutine pawdij2e1kb

m_hamiltonian/pawdij2ekb [ Functions ]

[ Top ] [ m_hamiltonian ] [ Functions ]

NAME

  pawdij2ekb

FUNCTION

  Transfer PAW Dij (on-site GS Hamiltonian) values
  from paw_ij datastructure to ekb array

INPUTS

OUTPUT

PARENTS

      m_hamiltonian

CHILDREN

      destroy_mpi_enreg,initmpi_seq,kpgsph,wrtout

SOURCE

1950 subroutine pawdij2ekb(ekb,paw_ij,isppol,comm_atom,mpi_atmtab)
1951 
1952 
1953 !This section has been created automatically by the script Abilint (TD).
1954 !Do not modify the following lines by hand.
1955 #undef ABI_FUNC
1956 #define ABI_FUNC 'pawdij2ekb'
1957 !End of the abilint section
1958 
1959  implicit none
1960 
1961 !Arguments ------------------------------------
1962 !scalars
1963  integer,intent(in) :: isppol,comm_atom
1964 !arrays
1965  integer,intent(in),optional,target :: mpi_atmtab(:)
1966  real(dp),intent(out) :: ekb(:,:,:,:)
1967  type(paw_ij_type),intent(in) :: paw_ij(:)
1968 
1969 !Local variables-------------------------------
1970 !scalars
1971  integer :: cplex,dimdij,dimekb1,dimekb3,dimekb4,iatom,iatom_tot,ierr,ii,isp,ispden,my_natom,natom
1972  logical :: my_atmtab_allocated,paral_atom
1973 !arrays
1974  integer,pointer :: my_atmtab(:)
1975 
1976 ! *************************************************************************
1977 
1978  DBG_ENTER("COLL")
1979 
1980  ekb=zero
1981 
1982 !Set up parallelism over atoms
1983  natom=size(ekb,2); my_natom=size(paw_ij)
1984  paral_atom=(xmpi_comm_size(comm_atom)>1)
1985  nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
1986  call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
1987 
1988  !Retrieve PAW Dij coefficients for this spin component
1989  if (my_natom>0) then
1990    if (allocated(paw_ij(1)%dij)) then
1991      dimekb1=size(ekb,1) ; dimekb3=size(ekb,3) ; dimekb4=size(ekb,4)
1992      cplex=paw_ij(1)%cplex_rf
1993      ABI_CHECK(cplex<=dimekb4,'paw_ij%cplex_rf>dimekb4!')
1994      do ii=1,cplex
1995        do ispden=1,dimekb3
1996          isp=isppol; if (dimekb3==4) isp=ispden
1997          do iatom=1,my_natom
1998            iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom)
1999            dimdij=paw_ij(iatom)%cplex_dij*paw_ij(iatom)%lmn2_size
2000            ABI_CHECK(dimdij<=dimekb1,'Size of paw_ij%dij>dimekb1!')
2001            ekb(1:dimdij,iatom_tot,ispden,ii)=paw_ij(iatom)%dij(1+(ii-1)*dimdij:ii*dimdij,isp)
2002          end do
2003        end do
2004      end do
2005    end if
2006  end if
2007 
2008 !Communication in case of distribution over atomic sites
2009  if (paral_atom) then
2010    call xmpi_sum(ekb,comm_atom,ierr)
2011  end if
2012 
2013 !Destroy atom table used for parallelism
2014  call free_my_atmtab(my_atmtab,my_atmtab_allocated)
2015 
2016  DBG_EXIT("COLL")
2017 
2018 end subroutine pawdij2ekb

m_hamiltonian/rf_hamiltonian_type [ Types ]

[ Top ] [ m_hamiltonian ] [ Types ]

NAME

 rf_hamiltonian_type

FUNCTION

 This datastructure contains few data about one 1st-order Hamiltonian,
 needed in the "getgh1c" routine, that apply the 1st-order Hamiltonian
 on a wavefunction.

SOURCE

415  type,public :: rf_hamiltonian_type
416 
417 ! ===== Integer scalars
418 
419   integer :: cplex
420    ! if 1, real space 1-order functions on FFT grid are REAL; if 2, COMPLEX
421 
422   integer :: dime1kb1
423    ! First dimension of E1kb, derivative of Ekb with respect to a perturbation
424 
425   integer :: dime1kb2
426    ! Second dimension of E1kb, derivative of Ekb with respect to a perturbation
427    ! NCPP: dime1kb2=ntypat, PAW: dime1kb2=natom
428 
429   integer :: npw_k
430    ! number of plane waves at k
431 
432   integer :: npw_kp
433    ! number of plane waves at k^prime
434 
435   integer:: nspinor
436    ! Number of spinorial components
437 
438   integer :: nsppol
439    ! Total number of spin components (1=non-polarized, 2=polarized)
440 
441   integer :: nvloc
442    ! Number of components of vloc
443    ! usually, nvloc=1, except in the non-collinear magnetism case, where nvloc=4
444 
445   integer :: n4,n5,n6
446    ! same as ngfft(4:6)
447 
448 ! ===== Real arrays
449 
450   real(dp), allocatable :: e1kbfr_spin(:,:,:,:,:)
451    ! e1kbfr_spin(dimekb1,dimekb2,nspinor**2,cplex,my_nsppol)
452    ! Contains the values of e1kbfr array for all spins treated by current process
453    ! See e1kbfr description ; e1kbfr is pointer to e1kbfr_spin(:,:,:,:,isppol)
454 
455   real(dp), allocatable :: e1kbsc_spin(:,:,:,:,:)
456    ! e1kbsc_spin(dimekb1,dimekb2,nspinor**2,cplex,my_nsppol)
457    ! Contains the values of e1kbsc array for all spins treated by current process
458    ! See e1kbsc description ; e1kbsc is pointer to e1kbsc_spin(:,:,:,:,isppol)
459 
460 ! ===== Real pointers
461 
462   real(dp), pointer :: dkinpw_k(:) => null()
463    ! dkinpw_k(npw_k)
464    ! 1st derivative of the (modified) kinetic energy for each plane wave at k
465 
466   real(dp), pointer :: dkinpw_kp(:) => null()
467    ! dkinpw_kp(npw_kp)
468    ! 1st derivative of the (modified) kinetic energy for each plane wave at k^prime
469 
470   real(dp), pointer :: ddkinpw_k(:) => null()
471    ! ddkinpw_k(npw_k)
472    ! 2nd derivative of the (modified) kinetic energy for each plane wave at k
473 
474   real(dp), pointer :: ddkinpw_kp(:) => null()
475    ! ddkinpw_kp(npw_kp)
476    ! 2nd derivative of the (modified) kinetic energy for each plane wave at k^prime
477 
478   real(dp), pointer :: e1kbfr(:,:,:,:) => null()
479    ! Frozen part of 1st derivative of ekb for the considered perturbation
480    ! (part not depending on VHxc^(1))
481    ! e1kbfr(dime1kb1,dime1kb2,nspinor**2,cplex)
482    ! For each spin component, e1kbfr points to e1kbfr_spin(:,:,:,:,my_isppol)
483 
484   real(dp), ABI_CONTIGUOUS pointer :: e1kbsc(:,:,:,:) => null()
485    ! Self-consistent 1st derivative of ekb for the considered perturbation
486    ! (part depending only on self-consistent VHxc^(1))
487    ! e1kbsc(dime1kb1,dime1kb2,nspinor**2,cplex)
488    ! For each spin component, e1kbfr points to e1kbfr_spin(:,:,:,:,my_isppol)
489 
490   real(dp), pointer :: vlocal1(:,:,:,:) => null()
491    ! vlocal1(cplex*n4,n5,n6,nvloc)
492    ! 1st-order local potential in real space, on the augmented fft grid
493 
494  end type rf_hamiltonian_type
495 
496  public :: init_rf_hamiltonian      ! Initialize the RF Hamiltonian
497  public :: destroy_rf_hamiltonian   ! Free the memory in the RF Hamiltonian
498  public :: load_spin_rf_hamiltonian ! Setup of the spin-dependent part of the RF Hamiltonian.
499  public :: load_k_rf_hamiltonian    ! Setup of the k-dependent part of the RF Hamiltonian