u_setup.f90 4.66 KB
Newer Older
1 2 3 4 5 6
!--------------------------------------------------------------------------------
! 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.
!--------------------------------------------------------------------------------

7 8 9 10 11 12 13 14 15
MODULE m_usetup
  USE m_juDFT
  !-------------------------------------------------------------------+
  ! Sets up the quantities needed for the LDA+U subroutines:          |
  !     radial integrals: us,dus,uds,duds                             |
  !     overlap of dot u: ddn                                         |
  !     potential matrix: vs_mmp                                      |
  !     total energy contribution: e_ldau                             |
  !                                                  G.B. Oct. 2000   |
Daniel Wortmann's avatar
Daniel Wortmann committed
16 17
  !                                                                   |
  !     Extension to multiple U per atom type  G.M. 2017              |
18 19
  !-------------------------------------------------------------------+
CONTAINS
20
  SUBROUTINE u_setup(sym,atoms,sphhar, input,el,inDen,pot,mpi,results)
21 22 23 24 25
    USE m_umtx
    USE m_uj2f
    USE m_nmat_rot
    USE m_vmmp
    USE m_types
26
    USE m_constants
27
    USE m_cdn_io
28
    IMPLICIT NONE
29
    TYPE(t_sym),INTENT(IN)          :: sym
30 31 32 33
    TYPE(t_results),INTENT(INOUT)   :: results
    TYPE(t_mpi),INTENT(IN)          :: mpi
    TYPE(t_input),INTENT(IN)        :: input
    TYPE(t_sphhar),INTENT(IN)       :: sphhar
34
    TYPE(t_atoms),INTENT(IN)        :: atoms
35
    TYPE(t_potden),INTENT(IN)       :: inDen
36
    TYPE(t_potden),INTENT(INOUT)    :: pot
37

38
    REAL,    INTENT(IN)           :: el(0:,:,:) !(0:atoms%lmaxd,ntype,jspd)
39
    ! ... Local Variables ...
Daniel Wortmann's avatar
Daniel Wortmann committed
40
    INTEGER itype,ispin,j,k,l,jspin,urec,i_u
41
    INTEGER noded,nodeu,ios
42 43
    REAL wronk
    CHARACTER*8 l_type*2,l_form*9
44
    REAL f(atoms%jmtd,2),g(atoms%jmtd,2),zero(atoms%n_u)
45 46
    REAL f0(atoms%n_u,input%jspins),f2(atoms%n_u,input%jspins),f4(atoms%n_u,input%jspins),f6(atoms%n_u,input%jspins)
    REAL, ALLOCATABLE :: u(:,:,:,:,:,:)
47
    COMPLEX, ALLOCATABLE :: n_mmp(:,:,:,:)
48 49 50
    !
    ! look, whether density matrix exists already:
    !
51
    IF (ANY(inDen%mmpMat(:,:,:,:).NE.0.0).AND.atoms%n_u>0) THEN
52

53
       ! calculate slater integrals from u and j
54 55
       CALL uj2f(input%jspins,atoms,f0,f2,f4,f6)

56
       ! set up e-e- interaction matrix
57
       ALLOCATE (u(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,&
58 59
                   -lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,MAX(1,atoms%n_u),input%jspins))
       ALLOCATE (n_mmp(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,MAX(1,atoms%n_u),input%jspins))
60
       n_mmp(:,:,:,:) = inDen%mmpMat(:,:,:,:)
61 62 63 64 65
       DO ispin = 1, 1 ! input%jspins
          f0(:,1) = (f0(:,1) + f0(:,input%jspins) ) / 2
          f2(:,1) = (f2(:,1) + f2(:,input%jspins) ) / 2
          f4(:,1) = (f4(:,1) + f4(:,input%jspins) ) / 2
          f6(:,1) = (f6(:,1) + f6(:,input%jspins) ) / 2
66 67
          CALL umtx(atoms,f0(1,ispin),f2(1,ispin),f4(1,ispin),f6(1,ispin),&
                    u(-lmaxU_const,-lmaxU_const,-lmaxU_const,-lmaxU_const,1,ispin))
Daniel Wortmann's avatar
Daniel Wortmann committed
68
       END DO
69

70
       ! check for possible rotation of n_mmp
71 72 73
       zero = 0.0
       CALL nmat_rot(atoms%lda_u(:)%phi,Atoms%lda_u(:)%theta,zero,3,atoms%n_u,input%jspins,atoms%lda_u%l,n_mmp)
       
74
       ! calculate potential matrix and total energy correction
75 76
       CALL v_mmp(sym,atoms,input%jspins,n_mmp,u,f0,f2,pot%mmpMat,results)

77 78 79
       IF (mpi%irank.EQ.0) THEN
          DO jspin = 1,input%jspins
             WRITE (6,'(a7,i3)') 'spin #',jspin
Daniel Wortmann's avatar
Daniel Wortmann committed
80 81 82 83 84 85
             DO i_u = 1, atoms%n_u
                itype = atoms%lda_u(i_u)%atomType
                l = atoms%lda_u(i_u)%l
                WRITE (l_type,'(i2)') 2*(2*l+1)
                l_form = '('//l_type//'f12.7)'
                WRITE (6,'(a20,i3)') 'n-matrix for atom # ',itype
86
                WRITE (6,l_form) ((n_mmp(k,j,i_u,jspin),k=-l,l),j=-l,l)
Daniel Wortmann's avatar
Daniel Wortmann committed
87 88 89 90 91
                WRITE (6,'(a20,i3)') 'V-matrix for atom # ',itype
                IF (atoms%lda_u(i_u)%l_amf) THEN
                   WRITE (6,*) 'using the around-mean-field limit '
                ELSE
                   WRITE (6,*) 'using the atomic limit of LDA+U '
92
                ENDIF
93
                WRITE (6,l_form) ((pot%mmpMat(k,j,i_u,jspin),k=-l,l),j=-l,l)
Daniel Wortmann's avatar
Daniel Wortmann committed
94 95
             END DO
          END DO
96 97
          WRITE (6,*) results%e_ldau
       ENDIF
98
       DEALLOCATE (u,n_mmp)
99 100 101 102
    ELSE
       IF (mpi%irank.EQ.0) THEN
          WRITE (*,*) 'no density matrix found ... skipping LDA+U'
       ENDIF
103
       pot%mmpMat(:,:,:,:) = CMPLX(0.0,0.0)
104 105 106 107 108 109
       results%e_ldau = 0.0
    ENDIF

    RETURN
  END SUBROUTINE u_setup
END MODULE m_usetup