kpoints.f90 3.71 KB
Newer Older
Daniel Wortmann's avatar
Daniel Wortmann committed
1 2 3 4 5 6 7 8
!--------------------------------------------------------------------------------
! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
! This file is part of FLEUR and available as free software under the conditions
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------

module m_kpoints
contains
9
  subroutine kpoints(oneD,sym,cell,input,noco,banddos,kpts,l_kpts)
10 11 12 13 14
    USE m_juDFT
    USE m_types
    USE m_julia
    USE m_kptgen_hybrid
    USE m_od_kptsgen
15
    USE m_unfold_band_kpts
Daniel Wortmann's avatar
Daniel Wortmann committed
16
  
Daniel Wortmann's avatar
Daniel Wortmann committed
17 18 19 20 21 22 23 24 25
    implicit none
    TYPE(t_input),INTENT(IN)   :: input
    TYPE(t_sym),INTENT(IN)     :: sym
    TYPE(t_oneD),INTENT(IN)     :: oneD
    TYPE(t_cell),INTENT(IN)     :: cell
    TYPE(t_banddos),INTENT(IN)  :: banddos
    TYPE(t_noco),INTENT(IN)     :: noco
    TYPE(t_kpts),INTENT(INOUT)  :: kpts
    LOGICAL,INTENT(IN)          :: l_kpts
26 27
    TYPE(t_kpts)                :: p_kpts
    TYPE(t_cell)                :: p_cell
Daniel Wortmann's avatar
Daniel Wortmann committed
28

29 30
    TYPE(t_sym) :: sym_hlp

31 32 33 34 35 36
    IF (input%l_wann) THEN
       IF (kpts%specificationType.NE.2) THEN
          CALL juDFT_error('l_wann only with kPointMesh', calledby = 'kpoints')
       END IF
    END IF

37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73

    IF (banddos%unfoldband) THEN
     CALL unfold_band_kpts(banddos,p_cell,cell,p_kpts,kpts)
     IF (.NOT.l_kpts) THEN
       IF (.NOT.oneD%odd%d1) THEN
          IF (input%l_wann) THEN
             sym_hlp=sym
             sym_hlp%nop=1
             sym_hlp%nop2=1
             CALL kptgen_hybrid(p_kpts,sym_hlp%invs,noco%l_soc,sym_hlp%nop,sym_hlp%mrot,sym_hlp%tau)
          ELSE IF (.FALSE.) THEN !this was used to generate q-points in jij case
             sym_hlp=sym
             sym_hlp%nop=1
             sym_hlp%nop2=1
             CALL julia(sym_hlp,p_cell,input,noco,banddos,p_kpts,.FALSE.,.TRUE.)
          ELSE IF (kpts%l_gamma.and.(banddos%ndir.eq.0)) THEN
             CALL kptgen_hybrid(p_kpts,sym%invs,noco%l_soc,sym%nop,sym%mrot,sym%tau)
          ELSE
             CALL julia(sym,p_cell,input,noco,banddos,p_kpts,.FALSE.,.TRUE.)
          END IF
       ELSE
          CALL juDFT_error('Error: No kpoint set generation for 1D systems yet!', calledby = 'kpoints')
          CALL od_kptsgen (kpts%nkpt)
       END IF
      END IF

     !Rescale weights and kpoints

     p_kpts%wtkpt(:) = p_kpts%wtkpt(:) / sum(p_kpts%wtkpt)
     p_kpts%bk(:,:) = p_kpts%bk(:,:) / p_kpts%posScale
     p_kpts%posScale = 1.0
     IF (p_kpts%nkpt3(3).EQ.0) p_kpts%nkpt3(3) = 1

    CALL find_supercell_kpts(banddos,p_cell,cell,p_kpts,kpts)

    ELSE
     IF (.NOT.l_kpts) THEN
74
       IF (.NOT.oneD%odd%d1) THEN
75 76 77 78 79
          IF (input%l_wann) THEN
             sym_hlp=sym
             sym_hlp%nop=1
             sym_hlp%nop2=1
             CALL kptgen_hybrid(kpts,sym_hlp%invs,noco%l_soc,sym_hlp%nop,sym_hlp%mrot,sym_hlp%tau)
80
          ELSE IF (.FALSE.) THEN !this was used to generate q-points in jij case
81 82 83 84
             sym_hlp=sym
             sym_hlp%nop=1
             sym_hlp%nop2=1
             CALL julia(sym_hlp,cell,input,noco,banddos,kpts,.FALSE.,.TRUE.)
85
          ELSE IF (kpts%l_gamma.and.(banddos%ndir.eq.0)) THEN
86 87 88 89 90
             CALL kptgen_hybrid(kpts,sym%invs,noco%l_soc,sym%nop,sym%mrot,sym%tau)
          ELSE
             CALL julia(sym,cell,input,noco,banddos,kpts,.FALSE.,.TRUE.)
          END IF
       ELSE
91
          CALL juDFT_error('Error: No kpoint set generation for 1D systems yet!', calledby = 'kpoints')
92 93
          CALL od_kptsgen (kpts%nkpt)
       END IF
94
     END IF
95

96 97 98 99 100 101
    !Rescale weights and kpoints

    kpts%wtkpt(:) = kpts%wtkpt(:) / sum(kpts%wtkpt)
    kpts%bk(:,:) = kpts%bk(:,:) / kpts%posScale
    kpts%posScale = 1.0
    IF (kpts%nkpt3(3).EQ.0) kpts%nkpt3(3) = 1
102
  END IF
Daniel Wortmann's avatar
Daniel Wortmann committed
103 104 105

end subroutine kpoints
end module m_kpoints