TABLE OF CONTENTS


ABINIT/gspot_transgrid_and_pack [ Functions ]

[ Top ] [ Functions ]

NAME

 gspot_transgrid_and_pack

FUNCTION

  Set up local potential vlocal on the coarse FFT mesh with proper dimensioning from vtrial given on the fine mesh.
  Also take into account the spin.

INPUTS

  isppol=Spin polarization.
  usepaw=1 if PAW
  paral_kgb: 1 if paral_kgb
  nfft=(effective) number of FFT grid points on the coarse mesh (for this processor)
  ngfft(18)contain all needed information about 3D FFT, for the coarse FFT mesh. see ~abinit/doc/variables/vargs.htm#ngfft
  nfftf=(effective) number of FFT grid points on the fine mesh (for this processor)
  nvloc==1 if nspden <=2, nvloc==4 for nspden==4,
  nspden=Number of spin density components.
  ncomp=Number of extra components in vtrial and vlocal (e.g. 1 if LDA/GGA pot, 4 for Meta-GGA, etc).
  pawfgr <type(pawfgr_type)>=fine grid parameters and related data
  vtrial(nfftf,nspden)=INPUT potential Vtrial(r).
  mpi_enreg=information about MPI parallelization

OUTPUT

  vlocal(n4,n5,n6,nvloc,ncomp): Potential on the coarse grid.

SIDE EFFECTS

SOURCE

1977 subroutine gspot_transgrid_and_pack(isppol, usepaw, paral_kgb,  nfft, ngfft, nfftf, &
1978                                     nspden, nvloc, ncomp, pawfgr, mpi_enreg, vtrial, vlocal)
1979 
1980 !Arguments -------------------------------
1981  integer,intent(in) :: isppol, nspden, ncomp, usepaw, paral_kgb, nfft, nfftf, nvloc
1982  type(pawfgr_type), intent(in) :: pawfgr
1983  type(MPI_type), intent(in) :: mpi_enreg
1984 !arrays
1985  integer,intent(in) :: ngfft(18)
1986  real(dp),intent(inout) :: vtrial(nfftf, nspden, ncomp)
1987  real(dp),intent(out) :: vlocal(ngfft(4), ngfft(5), ngfft(6), nvloc, ncomp)
1988 
1989 
1990 !Local variables-------------------------------
1991 !scalars
1992  integer :: n1,n2,n3,n4,n5,n6,ispden,ic
1993  real(dp) :: rhodum(1)
1994  real(dp),allocatable :: cgrvtrial(:,:), vlocal_tmp(:,:,:)
1995 
1996 ! *************************************************************************
1997 
1998  ! Coarse mesh.
1999  n1=ngfft(1); n2=ngfft(2); n3=ngfft(3)
2000  n4=ngfft(4); n5=ngfft(5); n6=ngfft(6)
2001 
2002  ! Set up local potential vlocal with proper dimensioning, from vtrial
2003  ! Also take into account the spin.
2004  if (nspden /= 4) then
2005    if (usepaw == 0 .or. pawfgr%usefinegrid == 0) then
2006      ! Fine mesh == Coarse mesh.
2007      do ic=1,ncomp
2008        call fftpac(isppol,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,vtrial(:,:,ic),vlocal(:,:,:,:,ic),2)
2009      end do
2010    else
2011      ! Transfer from fine mesh to coarse and then pack data
2012      ABI_MALLOC(cgrvtrial,(nfft,nspden))
2013      do ic=1,ncomp
2014        call transgrid(1,mpi_enreg,nspden,-1,0,0,paral_kgb,pawfgr,&
2015                       rhodum,rhodum,cgrvtrial,vtrial(:,:,ic))
2016        call fftpac(isppol,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,&
2017                    cgrvtrial,vlocal(:,:,:,:,ic),2)
2018      end do
2019      ABI_FREE(cgrvtrial)
2020    end if
2021  else
2022    ! nspden == 4. replace isppol by loop over ispden.
2023    ABI_MALLOC(vlocal_tmp, (n4,n5,n6))
2024    if (usepaw == 0 .or. pawfgr%usefinegrid == 0) then
2025      ! Fine mesh == Coarse mesh.
2026      do ic=1,ncomp
2027        do ispden=1,nspden
2028          call fftpac(ispden,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,vtrial(:,:,ic),vlocal_tmp,2)
2029          vlocal(:,:,:,ispden,ic) = vlocal_tmp(:,:,:)
2030        end do
2031      end do
2032    else
2033      ! Transfer from fine mesh to coarse and then pack data
2034      ABI_MALLOC(cgrvtrial,(nfft,nspden))
2035      do ic=1,ncomp
2036        call transgrid(1,mpi_enreg,nspden,-1,0,0,paral_kgb,pawfgr,rhodum,rhodum,cgrvtrial,vtrial(:,:,ic))
2037        do ispden=1,nspden
2038          call fftpac(ispden,mpi_enreg,nspden,n1,n2,n3,n4,n5,n6,ngfft,cgrvtrial,vlocal_tmp,2)
2039          vlocal(:,:,:,ispden,ic) = vlocal_tmp(:,:,:)
2040        end do
2041      end do
2042      ABI_FREE(cgrvtrial)
2043    end if
2044    ABI_FREE(vlocal_tmp)
2045  end if ! nspden
2046 
2047 end subroutine gspot_transgrid_and_pack

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.

TODO

  All array pointers in H datatypes should be declared as contiguous for efficient reasons
  (well, here performance is critical)
  Client code should make sure they always point contiguous targets.

COPYRIGHT

 Copyright (C) 2009-2022 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 .

SOURCE

24 #if defined HAVE_CONFIG_H
25 #include "config.h"
26 #endif
27 
28 #include "abi_common.h"
29 
30 module m_hamiltonian
31 
32  use defs_basis
33  use m_abicore
34  use m_errors
35  use m_xmpi
36 
37  use defs_datatypes,      only : pseudopotential_type
38  use defs_abitypes,       only : MPI_type
39  use m_copy,              only : addr_copy
40  use m_geometry,          only : metric
41  use m_pawtab,            only : pawtab_type
42  use m_pawfgr,            only : pawfgr_type
43  use m_fftcore,           only : sphereboundary
44  use m_fft,               only : fftpac
45  use m_fourier_interpol,  only : transgrid
46  use m_pawcprj,           only : pawcprj_getdim
47  use m_paw_ij,            only : paw_ij_type
48  use m_paral_atom,        only : get_my_atmtab, free_my_atmtab
49  use m_electronpositron,  only : electronpositron_type, electronpositron_calctype
50  use m_kg,                only : ph1d3d, getph
51  use m_fock,              only : fock_common_type, fock_BZ_type, fock_ACE_type, fock_type
52 
53 #if defined HAVE_GPU_CUDA
54  use m_manage_cuda
55 #endif
56 
57 #if defined HAVE_FC_ISO_C_BINDING
58  use, intrinsic :: iso_c_binding, only : c_ptr,c_loc,c_f_pointer,c_int32_t
59 #endif
60 
61  implicit none
62 
63  private
64 
65  public :: pawdij2ekb
66  public :: pawdij2e1kb
67  public :: gspot_transgrid_and_pack  ! Set up local potential vlocal on the coarse FFT mesh from vtrial on the fine mesh.
68 
69  ! These constants select how H is applied in reciprocal space
70  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.

SOURCE

1210 subroutine copy_hamiltonian(gs_hamk_out,gs_hamk_in)
1211 
1212 !Arguments ------------------------------------
1213  type(gs_hamiltonian_type),intent(in),target :: gs_hamk_in
1214  type(gs_hamiltonian_type),intent(out),target :: gs_hamk_out
1215 
1216 !Local variables-------------------------------
1217  integer :: tmp2i(5)
1218 #if defined HAVE_FC_ISO_C_BINDING
1219  type(C_PTR) :: ham_ptr
1220 #endif
1221 
1222 ! *************************************************************************
1223 
1224  DBG_ENTER("COLL")
1225 
1226 !@gs_hamiltonian_type
1227 
1228  gs_hamk_out%dimekb1 = gs_hamk_in%dimekb1
1229  gs_hamk_out%dimekb2 = gs_hamk_in%dimekb2
1230  gs_hamk_out%dimekbq = gs_hamk_in%dimekbq
1231  gs_hamk_out%istwf_k = gs_hamk_in%istwf_k
1232  gs_hamk_out%istwf_kp = gs_hamk_in%istwf_kp
1233  gs_hamk_out%lmnmax = gs_hamk_in%lmnmax
1234  gs_hamk_out%matblk = gs_hamk_in%matblk
1235  gs_hamk_out%mgfft = gs_hamk_in%mgfft
1236  gs_hamk_out%mpsang = gs_hamk_in%mpsang
1237  gs_hamk_out%mpssoang = gs_hamk_in%mpssoang
1238  gs_hamk_out%natom = gs_hamk_in%natom
1239  gs_hamk_out%nfft = gs_hamk_in%nfft
1240  gs_hamk_out%npw_k = gs_hamk_in%npw_k
1241  gs_hamk_out%npw_kp = gs_hamk_in%npw_kp
1242  gs_hamk_out%npw_fft_k = gs_hamk_in%npw_fft_k
1243  gs_hamk_out%npw_fft_kp = gs_hamk_in%npw_fft_kp
1244  gs_hamk_out%nspinor = gs_hamk_in%nspinor
1245  gs_hamk_out%nsppol = gs_hamk_in%nsppol
1246  gs_hamk_out%ntypat = gs_hamk_in%ntypat
1247  gs_hamk_out%nvloc = gs_hamk_in%nvloc
1248  gs_hamk_out%n4 = gs_hamk_in%n4
1249  gs_hamk_out%n5 = gs_hamk_in%n5
1250  gs_hamk_out%n6 = gs_hamk_in%n6
1251  gs_hamk_out%use_gpu_cuda = gs_hamk_in%use_gpu_cuda
1252  gs_hamk_out%usecprj = gs_hamk_in%usecprj
1253  gs_hamk_out%usepaw = gs_hamk_in%usepaw
1254  gs_hamk_out%useylm = gs_hamk_in%useylm
1255  gs_hamk_out%ngfft = gs_hamk_in%ngfft
1256  gs_hamk_out%nloalg = gs_hamk_in%nloalg
1257  gs_hamk_out%ucvol = gs_hamk_in%ucvol
1258  gs_hamk_out%gmet = gs_hamk_in%gmet
1259  gs_hamk_out%gprimd = gs_hamk_in%gprimd
1260  gs_hamk_out%kpt_k = gs_hamk_in%kpt_k
1261  gs_hamk_out%kpt_kp = gs_hamk_in%kpt_kp
1262 
1263  ABI_MALLOC(gs_hamk_out%atindx,(gs_hamk_out%natom))
1264  gs_hamk_out%atindx = gs_hamk_in%atindx
1265  ABI_MALLOC(gs_hamk_out%atindx1,(gs_hamk_out%natom))
1266  gs_hamk_out%atindx1 = gs_hamk_in%atindx1
1267  ABI_MALLOC(gs_hamk_out%dimcprj,(gs_hamk_out%natom*gs_hamk_out%usepaw))
1268  if (gs_hamk_out%usepaw==1) gs_hamk_out%dimcprj = gs_hamk_in%dimcprj
1269  ABI_MALLOC(gs_hamk_out%typat,(gs_hamk_out%natom))
1270  gs_hamk_out%typat = gs_hamk_in%typat
1271  ABI_MALLOC(gs_hamk_out%gbound_k,(2*gs_hamk_out%mgfft+8,2))
1272  gs_hamk_out%gbound_k = gs_hamk_in%gbound_k
1273  ABI_MALLOC(gs_hamk_out%indlmn,(6,gs_hamk_out%lmnmax,gs_hamk_out%ntypat))
1274  gs_hamk_out%indlmn = gs_hamk_in%indlmn
1275  ABI_MALLOC(gs_hamk_out%nattyp,(gs_hamk_out%ntypat))
1276  gs_hamk_out%nattyp = gs_hamk_in%nattyp
1277  ABI_MALLOC(gs_hamk_out%nucdipmom,(3,gs_hamk_out%natom))
1278  gs_hamk_out%nucdipmom = gs_hamk_in%nucdipmom
1279  ABI_MALLOC(gs_hamk_out%phkxred,(2,gs_hamk_out%natom))
1280  gs_hamk_out%phkxred = gs_hamk_in%phkxred
1281  ABI_MALLOC(gs_hamk_out%ph1d,(2,3*(2*gs_hamk_out%mgfft+1)*gs_hamk_out%natom))
1282  gs_hamk_out%ph1d = gs_hamk_in%ph1d
1283  ABI_MALLOC(gs_hamk_out%pspso,(gs_hamk_out%ntypat))
1284  gs_hamk_out%pspso = gs_hamk_in%pspso
1285  tmp2i(1:5)=shape(gs_hamk_in%ekb_spin)
1286  ABI_MALLOC(gs_hamk_out%ekb_spin,(tmp2i(1),tmp2i(2),tmp2i(3),tmp2i(4),tmp2i(5)))
1287  gs_hamk_out%ekb_spin = gs_hamk_in%ekb_spin
1288  gs_hamk_out%ekb => gs_hamk_out%ekb_spin(:,:,:,:,1)
1289  tmp2i(1:2)=shape(gs_hamk_in%sij)
1290  ABI_MALLOC(gs_hamk_out%sij,(tmp2i(1),tmp2i(2)))
1291  gs_hamk_out%sij = gs_hamk_in%sij
1292 
1293  if (associated(gs_hamk_in%gbound_kp,gs_hamk_in%gbound_k)) then
1294    gs_hamk_out%gbound_kp => gs_hamk_out%gbound_k
1295  else
1296    ABI_MALLOC(gs_hamk_out%gbound_kp,(2,gs_hamk_out%natom))
1297    gs_hamk_out%gbound_kp = gs_hamk_in%gbound_kp
1298  end if
1299  if (associated(gs_hamk_in%phkpxred,gs_hamk_in%phkxred)) then
1300    gs_hamk_out%phkpxred => gs_hamk_out%phkxred
1301  else
1302    ABI_MALLOC(gs_hamk_out%phkpxred,(2,gs_hamk_out%natom))
1303    gs_hamk_out%phkpxred = gs_hamk_in%phkpxred
1304  end if
1305 
1306  call addr_copy(gs_hamk_in%xred,gs_hamk_out%xred)
1307  call addr_copy(gs_hamk_in%vectornd,gs_hamk_out%vectornd)
1308  call addr_copy(gs_hamk_in%vlocal,gs_hamk_out%vlocal)
1309  call addr_copy(gs_hamk_in%vxctaulocal,gs_hamk_out%vxctaulocal)
1310  call addr_copy(gs_hamk_in%kinpw_k,gs_hamk_out%kinpw_k)
1311  call addr_copy(gs_hamk_in%kinpw_kp,gs_hamk_out%kinpw_kp)
1312  call addr_copy(gs_hamk_in%kg_k,gs_hamk_out%kg_k)
1313  call addr_copy(gs_hamk_in%kg_kp,gs_hamk_out%kg_kp)
1314  call addr_copy(gs_hamk_in%kpg_k,gs_hamk_out%kpg_k)
1315  call addr_copy(gs_hamk_in%kpg_kp,gs_hamk_out%kpg_kp)
1316  call addr_copy(gs_hamk_in%ffnl_k,gs_hamk_out%ffnl_k)
1317  call addr_copy(gs_hamk_in%ffnl_kp,gs_hamk_out%ffnl_kp)
1318  call addr_copy(gs_hamk_in%ph3d_k,gs_hamk_out%ph3d_k)
1319  call addr_copy(gs_hamk_in%ph3d_kp,gs_hamk_out%ph3d_kp)
1320 
1321 !For pointers to structured datatypes, have to copy the address
1322 !manually because there is no generic addr_copy function for that
1323  if (associated(gs_hamk_in%fockcommon)) then
1324 #if defined HAVE_FC_ISO_C_BINDING
1325    ham_ptr=c_loc(gs_hamk_in%fockcommon)
1326    call c_f_pointer(ham_ptr,gs_hamk_out%fockcommon)
1327 #else
1328    gs_hamk_out%fockcommon=transfer(gs_hamk_in%fockcommon,gs_hamk_out%fockcommon)
1329 #endif
1330  else
1331    nullify(gs_hamk_out%fockcommon)
1332  end if
1333  if (associated(gs_hamk_in%fockbz)) then
1334 #if defined HAVE_FC_ISO_C_BINDING
1335    ham_ptr=c_loc(gs_hamk_in%fockbz)
1336    call c_f_pointer(ham_ptr,gs_hamk_out%fockbz)
1337 #else
1338    gs_hamk_out%fockbz=transfer(gs_hamk_in%fockbz,gs_hamk_out%fockbz)
1339 #endif
1340  else
1341    nullify(gs_hamk_out%fockbz)
1342  end if
1343  if (associated(gs_hamk_in%fockACE_k)) then
1344 #if defined HAVE_FC_ISO_C_BINDING
1345    ham_ptr=c_loc(gs_hamk_in%fockACE_k)
1346    call c_f_pointer(ham_ptr,gs_hamk_out%fockACE_k)
1347 #else
1348    gs_hamk_out%fockACE_k=transfer(gs_hamk_in%fockACE_k,gs_hamk_out%fockACE_k)
1349 #endif
1350  else
1351    nullify(gs_hamk_out%fockACE_k)
1352  end if
1353 
1354  DBG_EXIT("COLL")
1355 
1356 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.

SOURCE

543 subroutine destroy_hamiltonian(Ham)
544 
545 !Arguments ------------------------------------
546 !scalars
547  class(gs_hamiltonian_type),intent(inout),target :: Ham
548 
549 ! *************************************************************************
550 
551  DBG_ENTER("COLL")
552 
553 !@gs_hamiltonian_type
554 
555 ! Integer Pointers
556  if (associated(Ham%gbound_kp,Ham%gbound_k)) then
557    nullify(Ham%gbound_kp)
558  else if (associated(Ham%gbound_kp)) then
559    ABI_FREE(Ham%gbound_kp)
560  end if
561 
562 ! Integer arrays
563  ABI_SFREE(Ham%atindx)
564  ABI_SFREE(Ham%atindx1)
565  ABI_SFREE(Ham%dimcprj)
566  ABI_SFREE(Ham%gbound_k)
567  ABI_SFREE(Ham%indlmn)
568  ABI_SFREE(Ham%nattyp)
569  ABI_SFREE(Ham%pspso)
570  ABI_SFREE(Ham%typat)
571 
572 ! Real Pointers
573  if (associated(Ham%phkpxred,Ham%phkxred)) then
574    nullify(Ham%phkpxred)
575  else if (associated(Ham%phkpxred)) then
576    ABI_FREE(Ham%phkpxred)
577  end if
578  ABI_SFREE(Ham%phkxred)
579  if (associated(Ham%ekb)) nullify(Ham%ekb)
580  if (associated(Ham%vectornd)) nullify(Ham%vectornd)
581  if (associated(Ham%vlocal)) nullify(Ham%vlocal)
582  if (associated(Ham%vxctaulocal)) nullify(Ham%vxctaulocal)
583  if (associated(Ham%xred)) nullify(Ham%xred)
584  if (associated(Ham%kinpw_k)) nullify(Ham%kinpw_k)
585  if (associated(Ham%kinpw_kp)) nullify(Ham%kinpw_kp)
586  if (associated(Ham%kg_k)) nullify(Ham%kg_k)
587  if (associated(Ham%kg_kp)) nullify(Ham%kg_kp)
588  if (associated(Ham%kpg_k)) nullify(Ham%kpg_k)
589  if (associated(Ham%kpg_kp)) nullify(Ham%kpg_kp)
590  if (associated(Ham%ffnl_k)) nullify(Ham%ffnl_k)
591  if (associated(Ham%ffnl_kp)) nullify(Ham%ffnl_kp)
592  if (associated(Ham%ph3d_k)) nullify(Ham%ph3d_k)
593  if (associated(Ham%ph3d_kp)) nullify(Ham%ph3d_kp)
594 
595 ! Real arrays
596  ABI_SFREE(Ham%ekb_spin)
597  ABI_SFREE(Ham%sij)
598  ABI_SFREE(Ham%nucdipmom)
599  ABI_SFREE(Ham%ph1d)
600 
601 ! Structured datatype pointers
602  if (associated(Ham%fockcommon)) nullify(Ham%fockcommon)
603  if (associated(Ham%fockACE_k)) nullify(Ham%fockACE_k)
604  if (associated(Ham%fockbz)) nullify(Ham%fockbz)
605 #if defined HAVE_GPU_CUDA
606  if(Ham%use_gpu_cuda==1) then
607    call gpu_finalize_ham_data()
608  end if
609 #endif
610 
611  DBG_EXIT("COLL")
612 
613 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.

SOURCE

1449 subroutine destroy_rf_hamiltonian(rf_Ham)
1450 
1451 !Arguments ------------------------------------
1452 !scalars
1453  class(rf_hamiltonian_type),intent(inout) :: rf_Ham
1454 
1455 ! *************************************************************************
1456 
1457  DBG_ENTER("COLL")
1458 
1459 !@rf_hamiltonian_type
1460 
1461 ! Real arrays
1462  ABI_SFREE(rf_Ham%e1kbfr_spin)
1463  ABI_SFREE(rf_Ham%e1kbsc_spin)
1464 
1465 ! Real pointers
1466  if (associated(rf_Ham%dkinpw_k)) nullify(rf_Ham%dkinpw_k)
1467  if (associated(rf_Ham%dkinpw_kp)) nullify(rf_Ham%dkinpw_kp)
1468  if (associated(rf_Ham%ddkinpw_k)) nullify(rf_Ham%ddkinpw_k)
1469  if (associated(rf_Ham%ddkinpw_kp)) nullify(rf_Ham%ddkinpw_kp)
1470  if (associated(rf_Ham%vectornd)) nullify(rf_Ham%vectornd)
1471  if (associated(rf_Ham%vlocal1)) nullify(rf_Ham%vlocal1)
1472  if (associated(rf_Ham%e1kbfr)) nullify(rf_Ham%e1kbfr)
1473  if (associated(rf_Ham%e1kbsc)) nullify(rf_Ham%e1kbsc)
1474 
1475  DBG_EXIT("COLL")
1476 
1477 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

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

SOURCE

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

SOURCE

1512 subroutine init_rf_hamiltonian(cplex,gs_Ham,ipert,rf_Ham,&
1513 &          comm_atom,mpi_atmtab,mpi_spintab,paw_ij1,has_e1kbsc) ! optional arguments
1514 
1515 !Arguments ------------------------------------
1516 !scalars
1517  integer,intent(in) :: cplex,ipert
1518  integer,intent(in),optional :: comm_atom
1519  logical,intent(in),optional :: has_e1kbsc
1520  type(gs_hamiltonian_type),intent(in) :: gs_Ham
1521  type(rf_hamiltonian_type),intent(inout),target :: rf_Ham
1522 !arrays
1523  integer,optional,intent(in)  :: mpi_atmtab(:),mpi_spintab(2)
1524  type(paw_ij_type),optional,intent(in) :: paw_ij1(:)
1525 
1526 !Local variables-------------------------------
1527 !scalars
1528  integer :: cplex_dij1,isp,jsp,my_comm_atom,my_nsppol
1529  logical :: has_e1kbsc_
1530 !arrays
1531  integer :: my_spintab(2)
1532  real(dp),allocatable,target :: e1kb_tmp(:,:,:,:)
1533 
1534 ! *************************************************************************
1535 
1536  DBG_ENTER("COLL")
1537 
1538 !@rf_hamiltonian_type
1539 
1540 !Manage optional parameters
1541  has_e1kbsc_=.false.;if (present(has_e1kbsc)) has_e1kbsc_=has_e1kbsc
1542  my_comm_atom=xmpi_comm_self;if (present(comm_atom)) my_comm_atom=comm_atom
1543  my_spintab=0;my_spintab(1:gs_Ham%nsppol)=1;if(present(mpi_spintab)) my_spintab=mpi_spintab
1544  my_nsppol=count(my_spintab==1)
1545 
1546  rf_Ham%cplex    =cplex
1547  rf_Ham%n4       =gs_Ham%n4
1548  rf_Ham%n5       =gs_Ham%n5
1549  rf_Ham%n6       =gs_Ham%n6
1550  rf_Ham%nvloc    =gs_Ham%nvloc
1551  rf_Ham%nsppol   =gs_Ham%nsppol
1552  rf_Ham%nspinor  =gs_Ham%nspinor
1553 
1554  rf_Ham%dime1kb1=0
1555  rf_Ham%dime1kb2=gs_Ham%dimekb2
1556  if (gs_Ham%usepaw==1.and.ipert/=gs_Ham%natom+1.and.ipert/=gs_Ham%natom+10) then
1557    cplex_dij1=1;if ((gs_Ham%nspinor==2).or.any(abs(gs_Ham%nucdipmom)>tol8)) cplex_dij1=2
1558    rf_Ham%dime1kb1=cplex_dij1*(gs_Ham%lmnmax*(gs_Ham%lmnmax+1))/2
1559  end if
1560 
1561   ! Allocate the arrays of the 1st-order Hamiltonian
1562   ! We preload here 1st-order non-local factors in order to avoid
1563   ! a communication over atoms inside the loop over spins.
1564  if (gs_Ham%usepaw==1.and.rf_Ham%dime1kb1>0) then
1565    if ((ipert>=1.and.ipert<=gs_Ham%natom).or.ipert==gs_Ham%natom+2.or.&
1566         ipert==gs_Ham%natom+3.or.ipert==gs_Ham%natom+4.or.ipert==gs_Ham%natom+11) then
1567 
1568      ABI_MALLOC(rf_Ham%e1kbfr_spin,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex,my_nsppol))
1569      rf_Ham%e1kbfr_spin=zero
1570      if (has_e1kbsc_) then
1571        ABI_MALLOC(rf_Ham%e1kbsc_spin,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex,my_nsppol))
1572        rf_Ham%e1kbsc_spin=zero
1573      end if
1574 
1575      if (present(paw_ij1)) then
1576 
1577        if (my_nsppol<rf_Ham%nsppol) then
1578          ABI_MALLOC(e1kb_tmp,(rf_Ham%dime1kb1,rf_Ham%dime1kb2,rf_Ham%nspinor**2,cplex))
1579        end if
1580 
1581 !      === Frozen term
1582        jsp=0
1583        do isp=1,rf_Ham%nsppol
1584          if (my_spintab(isp)==1) then
1585            jsp=jsp+1 ; rf_Ham%e1kbfr => rf_Ham%e1kbfr_spin(:,:,:,:,jsp)
1586          else
1587            rf_Ham%e1kbfr => e1kb_tmp
1588          end if
1589          if (present(mpi_atmtab)) then
1590            call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbfr=rf_Ham%e1kbfr,mpi_atmtab=mpi_atmtab)
1591          else
1592            call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbfr=rf_Ham%e1kbfr)
1593          end if
1594        end do
1595 
1596 !      === Self-consistent term
1597        if (has_e1kbsc_) then
1598          jsp=0
1599          do isp=1,rf_Ham%nsppol
1600            if (my_spintab(isp)==1) then
1601              jsp=jsp+1 ; rf_Ham%e1kbsc => rf_Ham%e1kbsc_spin(:,:,:,:,jsp)
1602            else
1603              rf_Ham%e1kbsc => e1kb_tmp
1604            end if
1605            if (present(mpi_atmtab)) then
1606              call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbsc=rf_Ham%e1kbsc,mpi_atmtab=mpi_atmtab)
1607            else
1608              call pawdij2e1kb(paw_ij1,isp,my_comm_atom,e1kbsc=rf_Ham%e1kbsc)
1609            end if
1610          end do
1611        end if
1612 
1613        if (my_nsppol<rf_Ham%nsppol) then
1614          ABI_FREE(e1kb_tmp)
1615        end if
1616 
1617      end if
1618    end if
1619  end if
1620 
1621  if (.not.allocated(rf_Ham%e1kbfr_spin)) then
1622    ABI_MALLOC(rf_Ham%e1kbfr_spin,(0,0,0,0,0))
1623  end if
1624  if (.not.allocated(rf_Ham%e1kbsc_spin)) then
1625    ABI_MALLOC(rf_Ham%e1kbsc_spin,(0,0,0,0,0))
1626  end if
1627  nullify(rf_Ham%e1kbfr)
1628  nullify(rf_Ham%e1kbsc)
1629 
1630  DBG_EXIT("COLL")
1631 
1632 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

SOURCE

 916 subroutine load_k_hamiltonian(ham,ffnl_k,fockACE_k,gbound_k,istwf_k,kinpw_k,&
 917                               kg_k,kpg_k,kpt_k,npw_k,npw_fft_k,ph3d_k,&
 918                               compute_gbound,compute_ph3d)
 919 
 920 !Arguments ------------------------------------
 921 !scalars
 922  integer,intent(in),optional :: npw_k,npw_fft_k,istwf_k
 923  logical,intent(in),optional :: compute_gbound,compute_ph3d
 924  class(gs_hamiltonian_type),intent(inout),target :: ham
 925 !arrays
 926  integer,intent(in),optional,target :: gbound_k(:,:),kg_k(:,:)
 927  real(dp),intent(in),optional :: kpt_k(3)
 928  real(dp),intent(in),optional,target :: ffnl_k(:,:,:,:),kinpw_k(:),kpg_k(:,:),ph3d_k(:,:,:)
 929  type(fock_ACE_type),intent(in),optional,target :: fockACE_k
 930 
 931 !Local variables-------------------------------
 932 !scalars
 933  integer :: iat,iatom
 934  logical :: compute_gbound_
 935  real(dp) :: arg
 936  !character(len=500) :: msg
 937 
 938 ! *************************************************************************
 939 
 940  DBG_ENTER("COLL")
 941 
 942 !@gs_hamiltonian_type
 943 
 944 !k-dependent scalars
 945  if (present(kpt_k)) then
 946    ham%kpt_k(:)  = kpt_k(:)
 947    ham%kpt_kp(:) = kpt_k(:)
 948  end if
 949  if (present(istwf_k)) then
 950    ham%istwf_k  = istwf_k
 951    ham%istwf_kp = istwf_k
 952  end if
 953  if (present(npw_k)) then
 954    ham%npw_k  = npw_k
 955    ham%npw_kp = npw_k
 956  end if
 957  if (present(npw_fft_k)) then
 958    ham%npw_fft_k  = npw_fft_k
 959    ham%npw_fft_kp = npw_fft_k
 960  else if (present(npw_k)) then
 961    ham%npw_fft_k  = npw_k
 962    ham%npw_fft_kp = npw_k
 963  end if
 964 
 965 !Pointers to k-dependent quantitites
 966  if (present(kinpw_k)) then
 967    ham%kinpw_k  => kinpw_k
 968    ham%kinpw_kp => kinpw_k
 969  end if
 970  if (present(kg_k)) then
 971    ham%kg_k  => kg_k
 972    ham%kg_kp => kg_k
 973  end if
 974  if (present(kpg_k)) then
 975    ham%kpg_k  => kpg_k
 976    ham%kpg_kp => kpg_k
 977  end if
 978  if (present(ffnl_k)) then
 979    ham%ffnl_k  => ffnl_k
 980    ham%ffnl_kp => ffnl_k
 981  end if
 982  if (present(ph3d_k)) then
 983    ham%ph3d_k  => ph3d_k
 984    ham%ph3d_kp => ph3d_k
 985  end if
 986  if (present(fockACE_k)) then
 987    ham%fockACE_k  => fockACE_k
 988  end if
 989 !Compute exp(i.k.R) for each atom
 990  if (present(kpt_k)) then
 991    if (associated(Ham%phkpxred).and.(.not.associated(Ham%phkpxred,Ham%phkxred))) then
 992      ABI_FREE(Ham%phkpxred)
 993    end if
 994    ABI_SFREE(ham%phkxred)
 995    ABI_MALLOC(ham%phkxred,(2,ham%natom))
 996    do iat=1,ham%natom
 997      iatom=ham%atindx(iat)
 998      arg=two_pi*DOT_PRODUCT(kpt_k,ham%xred(:,iat))
 999      ham%phkxred(1,iatom)=DCOS(arg)
1000      ham%phkxred(2,iatom)=DSIN(arg)
1001    end do
1002    ham%phkpxred => ham%phkxred
1003  end if
1004 
1005 !Compute or copy G sphere boundary at k+g
1006  compute_gbound_=.false.;if (present(compute_gbound)) compute_gbound_=compute_gbound
1007  if (present(gbound_k)) compute_gbound_=.true.
1008  if (compute_gbound_) then
1009    if (associated(Ham%gbound_kp,Ham%gbound_k)) then
1010      nullify(Ham%gbound_kp)
1011    else if (associated(Ham%gbound_kp)) then
1012      ABI_FREE(Ham%gbound_kp)
1013    end if
1014    ABI_SFREE(ham%gbound_k)
1015  end if
1016  if (.not.allocated(ham%gbound_k)) then
1017    ABI_MALLOC(ham%gbound_k,(2*ham%mgfft+8,2))
1018    ham%gbound_k(:,:)=0
1019    ham%gbound_kp => ham%gbound_k
1020  end if
1021  if (compute_gbound_) then
1022    if (present(gbound_k)) then
1023      ham%gbound_k(:,:)=gbound_k(:,:)
1024    else
1025      if (.not.associated(ham%kg_k)) then
1026        ABI_BUG('Something is missing for gbound_k computation!')
1027      end if
1028      !write(std_out,*)"About to call sphereboundary"
1029      !write(std_out,*)"size(kg_k), npw_k, mgfft",size(ham%kg_k, dim=2), ham%npw_k, ham%mgfft
1030      call sphereboundary(ham%gbound_k,ham%istwf_k,ham%kg_k,ham%mgfft,ham%npw_k)
1031    end if
1032    ham%gbound_kp => ham%gbound_k
1033  end if
1034 
1035 !Compute 3D structure factors for each atom at k+g
1036  if (present(compute_ph3d).and.present(ph3d_k)) then
1037    if (compute_ph3d.and.ham%nloalg(2)>0) then
1038      if ((.not.allocated(ham%phkxred)).or.(.not.associated(ham%kg_k)).or.&
1039          (.not.associated(ham%ph3d_k))) then
1040        ABI_BUG('Something is missing for ph3d_k computation!')
1041      end if
1042      call ph1d3d(1,ham%natom,ham%kg_k,ham%matblk,ham%natom,ham%npw_k,ham%ngfft(1),&
1043                  ham%ngfft(2),ham%ngfft(3),ham%phkxred,ham%ph1d,ham%ph3d_k)
1044    end if
1045  end if
1046 
1047  DBG_EXIT("COLL")
1048 
1049 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.

SOURCE

1727 subroutine load_k_rf_hamiltonian(rf_Ham,dkinpw_k,ddkinpw_k,npw_k)
1728 
1729 !Arguments ------------------------------------
1730 !scalars
1731  integer,intent(in),optional :: npw_k
1732  class(rf_hamiltonian_type),intent(inout),target :: rf_Ham
1733 !arrays
1734  real(dp),intent(in),optional,target :: dkinpw_k(:),ddkinpw_k(:)
1735 
1736 !Local variables-------------------------------
1737 
1738 ! *************************************************************************
1739 
1740  DBG_ENTER("COLL")
1741 
1742 !@gs_hamiltonian_type
1743 
1744 !k-dependent scalars
1745  if (present(npw_k)) then
1746    rf_Ham%npw_k  = npw_k
1747    rf_Ham%npw_kp = npw_k
1748  end if
1749 
1750 !Pointers to k-dependent quantitites
1751  if (present(dkinpw_k)) then
1752    rf_Ham%dkinpw_k  => dkinpw_k
1753    rf_Ham%dkinpw_kp => dkinpw_k
1754  end if
1755  if (present(ddkinpw_k)) then
1756    rf_Ham%ddkinpw_k  => ddkinpw_k
1757    rf_Ham%ddkinpw_kp => ddkinpw_k
1758  end if
1759 
1760  DBG_EXIT("COLL")
1761 
1762 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

SOURCE

1084 subroutine load_kprime_hamiltonian(ham,ffnl_kp,gbound_kp,istwf_kp,kinpw_kp,&
1085                                    kg_kp,kpg_kp,kpt_kp,npw_kp,npw_fft_kp,&
1086                                    ph3d_kp,compute_gbound,compute_ph3d)
1087 
1088 !Arguments ------------------------------------
1089 !scalars
1090  integer,intent(in),optional :: npw_kp,npw_fft_kp,istwf_kp
1091  logical,intent(in),optional :: compute_gbound,compute_ph3d
1092  class(gs_hamiltonian_type),intent(inout),target :: ham
1093 !arrays
1094  integer,intent(in),optional,target :: gbound_kp(:,:),kg_kp(:,:)
1095  real(dp),intent(in),optional :: kpt_kp(3)
1096  real(dp),intent(in),optional,target :: ffnl_kp(:,:,:,:),kinpw_kp(:),kpg_kp(:,:),ph3d_kp(:,:,:)
1097 
1098 !Local variables-------------------------------
1099 !scalars
1100  integer :: iat,iatom
1101  logical :: compute_gbound_
1102  real(dp) :: arg
1103  !character(len=500) :: msg
1104 
1105 ! *************************************************************************
1106 
1107  DBG_ENTER("COLL")
1108 
1109 !@gs_hamiltonian_type
1110 
1111 !k-dependent scalars
1112  if (present(kpt_kp))   ham%kpt_kp(:)= kpt_kp(:)
1113  if (present(istwf_kp)) ham%istwf_kp = istwf_kp
1114  if (present(npw_kp))   ham%npw_kp   = npw_kp
1115  if (present(npw_fft_kp)) then
1116     ham%npw_fft_kp = npw_fft_kp
1117  else if (present(npw_kp)) then
1118     ham%npw_fft_kp = npw_kp
1119  end if
1120 
1121 !Pointers to k-dependent quantitites
1122  if (present(kinpw_kp)) ham%kinpw_kp => kinpw_kp
1123  if (present(kg_kp))    ham%kg_kp    => kg_kp
1124  if (present(kpg_kp))   ham%kpg_kp   => kpg_kp
1125  if (present(ffnl_kp))  ham%ffnl_kp  => ffnl_kp
1126  if (present(ph3d_kp))  ham%ph3d_kp  => ph3d_kp
1127 
1128 !Compute exp(i.k^prime.R) for each atom
1129  if (present(kpt_kp)) then
1130    if (associated(ham%phkpxred,ham%phkxred)) then
1131      nullify(ham%phkpxred)
1132    else if (associated(ham%phkpxred)) then
1133      ABI_FREE(ham%phkpxred)
1134    end if
1135    ABI_MALLOC(ham%phkpxred,(2,ham%natom))
1136    do iat=1,ham%natom
1137      iatom=ham%atindx(iat)
1138      arg=two_pi*DOT_PRODUCT(kpt_kp,ham%xred(:,iat))
1139      ham%phkpxred(1,iatom)=DCOS(arg)
1140      ham%phkpxred(2,iatom)=DSIN(arg)
1141    end do
1142  end if
1143 
1144 !Compute or copy G sphere boundary at k^prime+g
1145  compute_gbound_=.false.
1146  if (present(kpt_kp).and.present(compute_gbound)) compute_gbound_=compute_gbound
1147  if (present(gbound_kp)) compute_gbound_=.true.
1148  if (compute_gbound_) then
1149    if (associated(ham%gbound_kp,ham%gbound_k)) then
1150      nullify(ham%gbound_kp)
1151    else if (associated(ham%gbound_kp)) then
1152      ABI_FREE(ham%gbound_kp)
1153    end if
1154    if (present(gbound_kp)) then
1155      ham%gbound_kp(:,:)=gbound_kp(:,:)
1156    else
1157      if (.not.associated(ham%kg_kp)) then
1158        ABI_BUG('Something is missing for gbound_kp computation!')
1159      end if
1160      ABI_MALLOC(ham%gbound_kp,(2*ham%mgfft+8,2))
1161      call sphereboundary(ham%gbound_kp,ham%istwf_kp,ham%kg_kp,ham%mgfft,ham%npw_kp)
1162    end if
1163  end if
1164 
1165 !Compute 3D structure factors for each atom at k^prime+g
1166  if (present(compute_ph3d).and.present(ph3d_kp)) then
1167    if (compute_ph3d.and.ham%nloalg(2)>0) then
1168      if ((.not.associated(ham%phkpxred)).or.(.not.associated(ham%kg_kp)).or.&
1169 &        (.not.associated(ham%ph3d_kp))) then
1170        ABI_BUG('Something is missing for ph3d_kp computation!')
1171      end if
1172      call ph1d3d(1,ham%natom,ham%kg_kp,ham%matblk,ham%natom,ham%npw_kp,ham%ngfft(1),&
1173 &                ham%ngfft(2),ham%ngfft(3),ham%phkpxred,ham%ph1d,ham%ph3d_kp)
1174    end if
1175  end if
1176 
1177  DBG_EXIT("COLL")
1178 
1179 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
  [vectornd(n4,n5,n6,nvloc,3)]=optional, vector potential of nuclear magnetic dipoles in real space
  [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.

SOURCE

1382 subroutine load_spin_hamiltonian(Ham,isppol,vectornd,vlocal,vxctaulocal,with_nonlocal)
1383 
1384 !Arguments ------------------------------------
1385 !scalars
1386  integer,intent(in) :: isppol
1387  logical,optional,intent(in) :: with_nonlocal
1388  class(gs_hamiltonian_type),intent(inout),target :: Ham
1389 !arrays
1390  real(dp),optional,intent(in),target :: vectornd(:,:,:,:,:)
1391  real(dp),optional,intent(in),target :: vlocal(:,:,:,:),vxctaulocal(:,:,:,:,:)
1392 
1393 !Local variables-------------------------------
1394 !scalars
1395  integer :: jsppol
1396 
1397 ! *************************************************************************
1398 
1399  DBG_ENTER("COLL")
1400 
1401  !@gs_hamiltonian_type
1402  if (present(vlocal)) then
1403    ABI_CHECK(size(vlocal)==Ham%n4*Ham%n5*Ham%n6*Ham%nvloc, "Wrong vlocal")
1404    Ham%vlocal => vlocal
1405  end if
1406  if (present(vxctaulocal)) then
1407    ABI_CHECK(size(vxctaulocal)==Ham%n4*Ham%n5*Ham%n6*Ham%nvloc*4, "Wrong vxctaulocal")
1408    Ham%vxctaulocal => vxctaulocal
1409  end if
1410  if (present(vectornd)) then
1411    ABI_CHECK(size(vectornd)==Ham%n4*Ham%n5*Ham%n6*Ham%nvloc*3, "Wrong vectornd")
1412    Ham%vectornd => vectornd
1413  end if
1414 
1415  ! Retrieve non-local factors for this spin component
1416  if (present(with_nonlocal)) then
1417    if (with_nonlocal) then
1418      jsppol=min(isppol,size(Ham%ekb_spin,5))
1419      if (jsppol>0) Ham%ekb => Ham%ekb_spin(:,:,:,:,jsppol)
1420    end if
1421  end if
1422 
1423  ! Update enl and sij on GPU
1424 #if defined HAVE_GPU_CUDA
1425  if (Ham%use_gpu_cuda==1) then
1426    call gpu_update_ham_data(Ham%ekb(:,:,:,1),size(Ham%ekb),Ham%sij,size(Ham%sij),Ham%gprimd,size(Ham%gprimd))
1427  end if
1428 #endif
1429 
1430  DBG_EXIT("COLL")
1431 
1432 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
  [vectornd(n4,n5,n6,nvloc)]=optional, vector potential of nuclear magnetic dipoles in real space in
   ddk direction idir
  [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.

SOURCE

1657 subroutine load_spin_rf_hamiltonian(rf_Ham,isppol,vectornd,vlocal1,with_nonlocal)
1658 
1659 !Arguments ------------------------------------
1660 !scalars
1661  integer,intent(in) :: isppol
1662  logical,optional,intent(in) :: with_nonlocal
1663  class(rf_hamiltonian_type),intent(inout),target :: rf_Ham
1664 !arrays
1665  real(dp),optional,target,intent(in) :: vlocal1(:,:,:,:)
1666  real(dp),optional,target,intent(in) :: vectornd(:,:,:,:)
1667 
1668 !Local variables-------------------------------
1669 !scalars
1670  integer :: jsppol
1671 
1672 ! *************************************************************************
1673 
1674  DBG_ENTER("COLL")
1675 
1676 !@rf_hamiltonian_type
1677 
1678  if (present(vlocal1)) then
1679    ABI_CHECK(size(vlocal1)==rf_Ham%cplex*rf_Ham%n4*rf_Ham%n5*rf_Ham%n6*rf_Ham%nvloc,"Wrong vlocal1")
1680    rf_Ham%vlocal1 => vlocal1
1681  end if
1682 
1683  if (present(vectornd)) then
1684    ABI_CHECK(size(vectornd)==rf_Ham%cplex*rf_Ham%n4*rf_Ham%n5*rf_Ham%n6*rf_Ham%nvloc,"Wrong vectornd")
1685    rf_Ham%vectornd => vectornd
1686  end if
1687 
1688  ! Retrieve non-local factors for this spin component
1689  if (present(with_nonlocal)) then
1690    if (with_nonlocal) then
1691      if (size(rf_Ham%e1kbfr_spin)>0) then
1692        jsppol=min(isppol,size(rf_Ham%e1kbfr_spin,5))
1693        if (jsppol>0) rf_Ham%e1kbfr => rf_Ham%e1kbfr_spin(:,:,:,:,jsppol)
1694      end if
1695      if (size(rf_Ham%e1kbsc_spin)>0) then
1696        jsppol=min(isppol,size(rf_Ham%e1kbsc_spin,5))
1697        if (jsppol>0) rf_Ham%e1kbsc => rf_Ham%e1kbsc_spin(:,:,:,:,jsppol)
1698      end if
1699    end if
1700  end if
1701 
1702  DBG_EXIT("COLL")
1703 
1704 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

SOURCE

1859 subroutine pawdij2e1kb(paw_ij1,isppol,comm_atom,mpi_atmtab,e1kbfr,e1kbsc)
1860 
1861 !Arguments ------------------------------------
1862 !scalars
1863  integer,intent(in) :: isppol,comm_atom
1864 !arrays
1865  integer,intent(in),optional,target :: mpi_atmtab(:)
1866  real(dp),optional,intent(out) :: e1kbfr(:,:,:,:),e1kbsc(:,:,:,:)
1867  type(paw_ij_type),intent(in) :: paw_ij1(:)
1868 
1869 !Local variables-------------------------------
1870 !scalars
1871  integer :: dimdij1,dime1kb1,dime1kb3,dime1kb4,iatom,iatom_tot,ierr,isp,ispden
1872  integer :: my_natom,natom,qphase
1873  logical :: my_atmtab_allocated,paral_atom
1874 !arrays
1875  integer,pointer :: my_atmtab(:)
1876 
1877 ! *************************************************************************
1878 
1879  DBG_ENTER("COLL")
1880 
1881  if ((.not.present(e1kbfr)).and.(.not.present(e1kbsc))) return
1882  if (present(e1kbfr)) then
1883    e1kbfr=zero ; natom=size(e1kbfr,2)
1884  end if
1885  if (present(e1kbsc)) then
1886    e1kbsc=zero ; natom=size(e1kbsc,2)
1887  end if
1888 
1889 !Set up parallelism over atoms
1890  my_natom=size(paw_ij1) ; paral_atom=(xmpi_comm_size(comm_atom)>1)
1891  nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
1892  call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
1893 
1894 !Retrieve 1st-order PAW Dij coefficients for this spin component (frozen)
1895  if (my_natom>0.and.present(e1kbfr)) then
1896    if (allocated(paw_ij1(1)%dijfr)) then
1897      dime1kb1=size(e1kbfr,1) ; dime1kb3=size(e1kbfr,3) ; dime1kb4=size(e1kbfr,4)
1898      ABI_CHECK(paw_ij1(1)%qphase==dime1kb4,'BUG in pawdij2e1kb (1)!')
1899      do ispden=1,dime1kb3
1900        isp=isppol;if (dime1kb3==4) isp=ispden
1901        do iatom=1,my_natom
1902          iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom)
1903          qphase=paw_ij1(iatom)%qphase
1904          dimdij1=paw_ij1(iatom)%cplex_dij*paw_ij1(iatom)%lmn2_size
1905          ABI_CHECK(dimdij1<=dime1kb1,'BUG: size of paw_ij1%dij>dime1kb1!')
1906          e1kbfr(1:dimdij1,iatom_tot,ispden,1)=paw_ij1(iatom)%dijfr(1:dimdij1,isp)
1907          if (qphase==2) e1kbfr(1:dimdij1,iatom_tot,ispden,2)=paw_ij1(iatom)%dijfr(dimdij1+1:2*dimdij1,isp)
1908        end do
1909      end do
1910    end if
1911  end if
1912 
1913 !Retrieve 1st-order PAW Dij coefficients for this spin component (self-consistent)
1914  if (my_natom>0.and.present(e1kbsc)) then
1915    if (allocated(paw_ij1(1)%dijfr).and.allocated(paw_ij1(1)%dij)) then
1916      dime1kb1=size(e1kbsc,1) ; dime1kb3=size(e1kbsc,3) ; dime1kb4=size(e1kbsc,4)
1917      ABI_CHECK(paw_ij1(1)%qphase==dime1kb4,'BUG in pawdij2e1kb (1)!')
1918      do ispden=1,dime1kb3
1919        isp=isppol;if (dime1kb3==4) isp=ispden
1920        do iatom=1,my_natom
1921          iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom)
1922          qphase=paw_ij1(iatom)%qphase
1923          dimdij1=paw_ij1(iatom)%cplex_dij*paw_ij1(iatom)%lmn2_size
1924          ABI_CHECK(dimdij1<=dime1kb1,'BUG: size of paw_ij1%dij>dime1kb1!')
1925          e1kbsc(1:dimdij1,iatom_tot,ispden,1)=paw_ij1(iatom)%dij  (1:dimdij1,isp) &
1926 &                                            -paw_ij1(iatom)%dijfr(1:dimdij1,isp)
1927          if (qphase==2) e1kbsc(1:dimdij1,iatom_tot,ispden,2)=paw_ij1(iatom)%dij  (dimdij1+1:2*dimdij1,isp) &
1928 &                                                           -paw_ij1(iatom)%dijfr(dimdij1+1:2*dimdij1,isp)
1929        end do
1930      end do
1931    end if
1932  end if
1933 
1934  ! Communication in case of distribution over atomic sites
1935  if (paral_atom) then
1936    if (present(e1kbfr)) call xmpi_sum(e1kbfr,comm_atom,ierr)
1937    if (present(e1kbsc)) call xmpi_sum(e1kbsc,comm_atom,ierr)
1938  end if
1939 
1940 !Destroy atom table used for parallelism
1941  call free_my_atmtab(my_atmtab,my_atmtab_allocated)
1942 
1943  DBG_EXIT("COLL")
1944 
1945 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

SOURCE

1781 subroutine pawdij2ekb(ekb,paw_ij,isppol,comm_atom,mpi_atmtab)
1782 
1783 !Arguments ------------------------------------
1784 !scalars
1785  integer,intent(in) :: isppol,comm_atom
1786 !arrays
1787  integer,intent(in),optional,target :: mpi_atmtab(:)
1788  real(dp),intent(out) :: ekb(:,:,:,:)
1789  type(paw_ij_type),intent(in) :: paw_ij(:)
1790 
1791 !Local variables-------------------------------
1792 !scalars
1793  integer :: dimdij,dimekb1,dimekb3,dimekb4,iatom,iatom_tot,ierr,ii,isp,ispden,my_natom,natom,qphase
1794  logical :: my_atmtab_allocated,paral_atom
1795 !arrays
1796  integer,pointer :: my_atmtab(:)
1797 
1798 ! *************************************************************************
1799 
1800  DBG_ENTER("COLL")
1801 
1802  ekb=zero
1803 
1804 !Set up parallelism over atoms
1805  natom=size(ekb,2); my_natom=size(paw_ij)
1806  paral_atom=(xmpi_comm_size(comm_atom)>1)
1807  nullify(my_atmtab);if (present(mpi_atmtab)) my_atmtab => mpi_atmtab
1808  call get_my_atmtab(comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
1809 
1810  !Retrieve PAW Dij coefficients for this spin component
1811  if (my_natom>0) then
1812    if (allocated(paw_ij(1)%dij)) then
1813      dimekb1=size(ekb,1) ; dimekb3=size(ekb,3) ; dimekb4=size(ekb,4)
1814      qphase=paw_ij(1)%qphase
1815      ABI_CHECK(qphase<=dimekb4,'paw_ij%qphase>dimekb4!')
1816      do ii=1,qphase
1817        do ispden=1,dimekb3
1818          isp=isppol; if (dimekb3==4) isp=ispden
1819          do iatom=1,my_natom
1820            iatom_tot=iatom;if (paral_atom) iatom_tot=my_atmtab(iatom)
1821            dimdij=paw_ij(iatom)%cplex_dij*paw_ij(iatom)%lmn2_size
1822            ABI_CHECK(dimdij<=dimekb1,'Size of paw_ij%dij>dimekb1!')
1823            ekb(1:dimdij,iatom_tot,ispden,ii)=paw_ij(iatom)%dij(1+(ii-1)*dimdij:ii*dimdij,isp)
1824          end do
1825        end do
1826      end do
1827    end if
1828  end if
1829 
1830 !Communication in case of distribution over atomic sites
1831  if (paral_atom) then
1832    call xmpi_sum(ekb,comm_atom,ierr)
1833  end if
1834 
1835 !Destroy atom table used for parallelism
1836  call free_my_atmtab(my_atmtab,my_atmtab_allocated)
1837 
1838  DBG_EXIT("COLL")
1839 
1840 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

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