IffGit has a new shared runner for building Docker images in GitLab CI. Visit https://iffgit.fz-juelich.de/examples/ci-docker-in-docker for more details.

kpoints.f90 2.99 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
28
    TYPE(t_kpts)                :: tmp_kpts
Daniel Wortmann's avatar
Daniel Wortmann committed
29

30
31
    TYPE(t_sym) :: sym_hlp

32
33
34
35
36
37
    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

38
    IF (.NOT.l_kpts) THEN
39
       IF (.NOT.oneD%odd%d1) THEN
40
41
42
43
          IF (input%l_wann) THEN
             sym_hlp=sym
             sym_hlp%nop=1
             sym_hlp%nop2=1
44
             CALL kptgen_hybrid(input,cell,sym_hlp,kpts,noco%l_soc)
45
          ELSE IF (.FALSE.) THEN !this was used to generate q-points in jij case
46
47
48
49
             sym_hlp=sym
             sym_hlp%nop=1
             sym_hlp%nop2=1
             CALL julia(sym_hlp,cell,input,noco,banddos,kpts,.FALSE.,.TRUE.)
50
          ELSE IF (kpts%l_gamma.and.(banddos%ndir.eq.0)) THEN
51
             CALL kptgen_hybrid(input,cell,sym,kpts,noco%l_soc)
52
          ELSE
53
54
             IF (banddos%unfoldband) THEN
               CALL unfold_band_kpts(banddos,p_cell,cell,p_kpts,kpts)
55
 	       CALL julia(sym,cell,input,noco,banddos,kpts,.FALSE.,.TRUE.)
56
57
58
59
60
               CALL julia(sym,p_cell,input,noco,banddos,p_kpts,.FALSE.,.TRUE.)
               CALL find_supercell_kpts(banddos,p_cell,cell,p_kpts,kpts)
             ELSE
               CALL julia(sym,cell,input,noco,banddos,kpts,.FALSE.,.TRUE.)
             END IF
61
62
          END IF
       ELSE
63
          CALL juDFT_error('Error: No kpoint set generation for 1D systems yet!', calledby = 'kpoints')
64
65
          CALL od_kptsgen (kpts%nkpt)
       END IF
66

67
68
69
70
71
72
73
74
75
76
77
78
       !Rescale weights and kpoints
       IF (.not.banddos%unfoldband) THEN
          kpts%wtkpt(:) = kpts%wtkpt(:) / sum(kpts%wtkpt)
       END IF
       kpts%bk(:,:) = kpts%bk(:,:) / kpts%posScale
       kpts%posScale = 1.0
       IF (kpts%nkpt3(3).EQ.0) kpts%nkpt3(3) = 1
    ELSE
       IF (banddos%unfoldband) THEN
          CALL unfold_band_kpts(banddos,p_cell,cell,p_kpts,kpts)
          CALL find_supercell_kpts(banddos,p_cell,cell,p_kpts,kpts)
       END IF
79
    END IF
Daniel Wortmann's avatar
Daniel Wortmann committed
80
81
82

end subroutine kpoints
end module m_kpoints