broyden.F90 4.08 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13
MODULE m_broyden
  USE m_juDFT
  !################################################################
  !     IMIX = 3 : BROYDEN'S FIRST METHOD
  !     IMIX = 5 : BROYDEN'S SECOND METHOD
  !     IMIX = 7 : GENERALIZED ANDERSEN METHOD
  !     sm   : input charge density of iteration m
  !            afterwards update rho(m+1)
  !     fm   : output minus input charge density of iteration m
  !     sm1  : input charge density of iteration m-1
  !     fm1   : output minus inputcharge density of iteration m-1
  !################################################################
CONTAINS
Daniel Wortmann's avatar
Daniel Wortmann committed
14
  SUBROUTINE broyden(input,fm,sm)
15
    USE m_types
16
    USE m_broyd_io
Daniel Wortmann's avatar
Daniel Wortmann committed
17
    USE m_types_mixvector
18
    IMPLICIT NONE
19

Daniel Wortmann's avatar
Daniel Wortmann committed
20 21 22
    TYPE(t_input),INTENT(IN)        :: input
    TYPE(t_mixvector),INTENT(IN)    :: fm
    TYPE(t_mixvector),INTENT(INOUT) :: sm
23

Daniel Wortmann's avatar
Daniel Wortmann committed
24 25
    
    
26
    ! Local Scalars
27
    INTEGER         :: i,it,k,nit,iread,nmaph, mit
28
    REAL            :: bm,dfivi,fmvm,smnorm,vmnorm,alphan
Daniel Wortmann's avatar
Daniel Wortmann committed
29
    LOGICAL         :: l_exist
30 31 32
    REAL, PARAMETER :: one=1.0, zero=0.0

    ! Local Arrays
Daniel Wortmann's avatar
Daniel Wortmann committed
33 34
    REAL,ALLOCATABLE  :: am(:)
    TYPE(t_mixvector) :: fm1,sm1,ui,um,vi,vm
35

Daniel Wortmann's avatar
Daniel Wortmann committed
36
    
37 38
    dfivi = zero

Daniel Wortmann's avatar
Daniel Wortmann committed
39 40 41 42 43 44 45
    CALL fm1%alloc()
    CALL sm1%alloc()
    CALL ui%alloc()
    CALL um%alloc()
    CALL vi%alloc()
    CALL vm%alloc()
    
46
    ALLOCATE (am(input%maxiter+1))
47
    am  = 0.0
48
    l_exist = initBroydenHistory(input,hybrid,nmap) ! returns true if there already exists a Broyden history
Daniel Wortmann's avatar
Daniel Wortmann committed
49
    IF (l_exist) THEN
50 51
       ! load input charge density (sm1) and difference of 
       ! in and out charge densities (fm1) from previous iteration (m-1)
52

Daniel Wortmann's avatar
Daniel Wortmann committed
53
       CALL readLastIterInAndDiffDen(mit,alphan,sm1,fm1)
54
       IF (ABS(input%alpha-alphan) > 0.0001) THEN
55
          WRITE (6,*) 'mixing parameter has been changed; reset'
Daniel Wortmann's avatar
Daniel Wortmann committed
56
          WRITE (6,*) 'broyden algorithm or set alpha to',alphan
57 58 59 60 61
          CALL juDFT_error("mixing parameter (input) changed", calledby ="broyden")
       END IF

       ! generate F_m   - F_(m-1)  ... sm1
       !      and rho_m - rho_(m-1) .. fm1
Daniel Wortmann's avatar
Daniel Wortmann committed
62 63 64 65
       sm1 = sm - sm1
       fm1 = fm - fm1
    ELSE
       mit=1
66
    END IF
67 68

    ! save F_m and rho_m for next iteration
69 70
    nit = mit +1
    IF (nit > input%maxiter+1) nit = 1
71
    CALL writeLastIterInAndDiffDen(hybrid,nmap,nit,input%alpha,sm,fm)
72

Daniel Wortmann's avatar
Daniel Wortmann committed
73
    IF (.NOT.l_exist) THEN 
74 75
       !     update for rho for mit=1 is straight mixing
       !     sm = sm + alpha*fm
Daniel Wortmann's avatar
Daniel Wortmann committed
76
       sm=sm+input%alpha*fm
77 78 79
    ELSE
       !     |vi> = w|vi> 
       !     loop to generate um : um = sm1 + alpha*fm1 - \sum <fm1|w|vi> ui
Daniel Wortmann's avatar
Daniel Wortmann committed
80 81
       um = input%alpha * fm1 + sm1
       
82 83
       iread = MIN(mit-1,input%maxiter+1)
       DO it = 2, iread
84 85 86
          CALL readUVec(input,hybrid,nmap,it-mit,mit,ui)
          CALL readVVec(input,hybrid,nmap,it-mit,mit,dfivi,vi)

Daniel Wortmann's avatar
Daniel Wortmann committed
87
          am(it) = vi.dot.fm
88
          ! calculate um(:) = -am(it)*ui(:) + um(:)
Daniel Wortmann's avatar
Daniel Wortmann committed
89
          um=um-am(it)*ui
90 91
          WRITE(6,FMT='(5x,"<vi|w|Fm> for it",i2,5x,f10.6)')it,am(it) 
       END DO
92

Daniel Wortmann's avatar
Daniel Wortmann committed
93
       IF (input%imix.EQ.7) THEN
94 95 96 97 98 99
          !****************************************
          !      generalized anderson method
          !****************************************

          ! calculate vm = alpha*wfm1 -\sum <fm1|w|vi> <fi1|w|vi><vi|
          ! convolute fm1 with the metrik and store in vm
Daniel Wortmann's avatar
Daniel Wortmann committed
100
          vm=fm1%apply_metric()
101
          DO it = 2,iread
102
             CALL readVVec(input,hybrid,nmap,it-mit,mit,dfivi,vi)
103
             ! calculate vm(:) = -am(it)*dfivi*vi(:) + vm
Daniel Wortmann's avatar
Daniel Wortmann committed
104
             vm=vm-am(it)*dfivi*vi
105
          END DO
106

Daniel Wortmann's avatar
Daniel Wortmann committed
107
          vmnorm=fm1.dot.vm
108
          ! if (vmnorm.lt.tol_10) stop
109

110
          ! calculate vm(:) = (1.0/vmnorm)*vm(:)
Daniel Wortmann's avatar
Daniel Wortmann committed
111 112
          vm=(1.0/vmnorm)*vm
     
113
          ! save dfivi(mit) for next iteration
114
          dfivi = vmnorm
Daniel Wortmann's avatar
Daniel Wortmann committed
115 116
       ELSE
          CALL judft_error("Only generalized Anderson implemented")
117
       END IF
118 119 120

       ! write um,vm and dfivi on file broyd.?

121 122
       CALL writeUVec(input,hybrid,nmap,mit,um)
       CALL writeVVec(input,hybrid,nmap,mit,dfivi,vm)
123 124 125

       ! update rho(m+1)
       ! calculate <fm|w|vm>
Daniel Wortmann's avatar
Daniel Wortmann committed
126
       fmvm = vm.dot.fm
127
       ! calculate sm(:) = (1.0-fmvm)*um(:) + sm
Daniel Wortmann's avatar
Daniel Wortmann committed
128
       sm=sm+(one-fmvm)*um
129
    END IF
Daniel Wortmann's avatar
Daniel Wortmann committed
130
 
131 132
  END SUBROUTINE broyden
END MODULE m_broyden