kpoints.f90 2.78 KB
Newer Older
Daniel Wortmann's avatar
Daniel Wortmann committed
1 2 3 4 5 6 7 8 9
!--------------------------------------------------------------------------------
! 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
  subroutine kpoints(oneD,jij,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
Daniel Wortmann's avatar
Daniel Wortmann committed
15
  
Daniel Wortmann's avatar
Daniel Wortmann committed
16 17 18 19 20 21 22 23 24 25 26
    implicit none
    TYPE(t_input),INTENT(IN)   :: input
    TYPE(t_sym),INTENT(IN)     :: sym
    TYPE(t_oneD),INTENT(IN)     :: oneD
    TYPE(t_Jij),INTENT(IN)      :: Jij
    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

27 28 29 30 31 32 33 34 35 36 37 38 39 40 41
    TYPE(t_sym) :: sym_hlp
    INTEGER     :: i
    REAL        :: recVecLength

    IF (kpts%specificationType.EQ.4) THEN
       DO i = 1, 3
          IF (kpts%kPointDensity(i).LE.0.0) THEN
             CALL juDFT_error('Error: Nonpositive kpointDensity provided', calledby = 'kpoints')
          END IF
          recVecLength = SQRT(cell%bmat(i,1)**2 + cell%bmat(i,2)**2 + cell%bmat(i,3)**2)
          kpts%nkpt3(i) = CEILING(kpts%kPointDensity(i) * recVecLength)
       END DO
       kpts%nkpt = kpts%nkpt3(1) * kpts%nkpt3(2) * kpts%nkpt3(3)
    END IF

42 43 44 45 46 47
    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

48 49
    IF (.NOT.l_kpts) THEN
       IF (.NOT.oneD%odd%d1) THEN
50 51 52 53 54 55
          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)
          ELSE IF ((jij%l_J)) THEN
56 57 58 59
             sym_hlp=sym
             sym_hlp%nop=1
             sym_hlp%nop2=1
             CALL julia(sym_hlp,cell,input,noco,banddos,kpts,.FALSE.,.TRUE.)
60
          ELSE IF (kpts%l_gamma.and.(banddos%ndir.eq.0)) THEN
61 62 63 64 65
             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
66
          CALL juDFT_error('Error: No kpoint set generation for 1D systems yet!', calledby = 'kpoints')
67 68 69
          CALL od_kptsgen (kpts%nkpt)
       END IF
    END IF
70

71 72 73 74 75 76
    !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
Daniel Wortmann's avatar
Daniel Wortmann committed
77 78 79

end subroutine kpoints
end module m_kpoints