Commit fd2b2b92 authored by Daniel Wortmann's avatar Daniel Wortmann

More coding...

parent 4f9b2987
This diff is collapsed.
......@@ -11,79 +11,46 @@ MODULE m_broyden
! fm1 : output minus inputcharge density of iteration m-1
!################################################################
CONTAINS
SUBROUTINE broyden(cell,stars,atoms,vacuum,sphhar,input,noco,oneD,sym,&
hybrid,mmap,nmaph,mapmt,mapvac2,nmap,fm,sm,lpot)
#include"cpp_double.h"
USE m_metric
SUBROUTINE broyden(input,fm,sm)
USE m_types
USE m_broyd_io
USE m_types_mixvector
IMPLICIT NONE
TYPE(t_oneD),INTENT(IN) :: oneD
TYPE(t_input),INTENT(IN) :: input
TYPE(t_vacuum),INTENT(IN) :: vacuum
TYPE(t_noco),INTENT(IN) :: noco
TYPE(t_sym),INTENT(IN) :: sym
TYPE(t_stars),INTENT(IN) :: stars
TYPE(t_cell),INTENT(IN) :: cell
TYPE(t_sphhar),INTENT(IN) :: sphhar
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_hybrid),INTENT(IN) :: hybrid
! Scalar Arguments
INTEGER, INTENT (IN) :: mmap,nmap
INTEGER, INTENT (IN) :: mapmt,mapvac2
LOGICAL,OPTIONAL,INTENT(IN) :: lpot
! Array Arguments
REAL, INTENT (IN) :: fm(nmap)
REAL, INTENT (INOUT) :: sm(nmap)
TYPE(t_input),INTENT(IN) :: input
TYPE(t_mixvector),INTENT(IN) :: fm
TYPE(t_mixvector),INTENT(INOUT) :: sm
! Local Scalars
INTEGER :: i,it,k,nit,iread,nmaph, mit
REAL :: bm,dfivi,fmvm,smnorm,vmnorm,alphan
LOGICAL :: l_pot, l_exist
LOGICAL :: l_exist
REAL, PARAMETER :: one=1.0, zero=0.0
! Local Arrays
REAL, ALLOCATABLE :: am(:)
REAL, ALLOCATABLE :: fm1(:),sm1(:),ui(:),um(:),vi(:),vm(:)
! External Functions
REAL CPP_BLAS_sdot
EXTERNAL CPP_BLAS_sdot
! External Subroutines
EXTERNAL CPP_BLAS_saxpy,CPP_BLAS_sscal
REAL,ALLOCATABLE :: am(:)
TYPE(t_mixvector) :: fm1,sm1,ui,um,vi,vm
dfivi = zero
l_pot = .FALSE.
IF (PRESENT(lpot)) l_pot = lpot
ALLOCATE (fm1(mmap),sm1(mmap),ui(mmap),um(mmap),vi(mmap),vm(mmap))
CALL fm1%alloc()
CALL sm1%alloc()
CALL ui%alloc()
CALL um%alloc()
CALL vi%alloc()
CALL vm%alloc()
ALLOCATE (am(input%maxiter+1))
fm1 = 0.0
sm1 = 0.0
ui = 0.0
um = 0.0
vi = 0.0
vm = 0.0
am = 0.0
mit = 0
l_exist = initBroydenHistory(input,hybrid,nmap) ! returns true if there already exists a Broyden history
IF(.NOT.l_exist) mit = 1
IF (mit.NE.1) THEN
IF (l_exist) THEN
! load input charge density (sm1) and difference of
! in and out charge densities (fm1) from previous iteration (m-1)
CALL readLastIterInAndDiffDen(hybrid,nmap,mit,alphan,sm1(:nmap),fm1(:nmap))
CALL readLastIterInAndDiffDen(mit,alphan,sm1,fm1)
IF (ABS(input%alpha-alphan) > 0.0001) THEN
WRITE (6,*) 'mixing parameter has been changed; reset'
WRITE (6,*) 'broyden algorithm or set alpha to',alphan
......@@ -92,8 +59,10 @@ CONTAINS
! generate F_m - F_(m-1) ... sm1
! and rho_m - rho_(m-1) .. fm1
sm1(1:nmap) = sm(1:nmap) - sm1(1:nmap)
fm1(1:nmap) = fm(1:nmap) - fm1(1:nmap)
sm1 = sm - sm1
fm1 = fm - fm1
ELSE
mit=1
END IF
! save F_m and rho_m for next iteration
......@@ -101,95 +70,50 @@ CONTAINS
IF (nit > input%maxiter+1) nit = 1
CALL writeLastIterInAndDiffDen(hybrid,nmap,nit,input%alpha,sm,fm)
IF (mit.EQ.1) THEN
IF (.NOT.l_exist) THEN
! update for rho for mit=1 is straight mixing
! sm = sm + alpha*fm
CALL CPP_BLAS_saxpy(nmap,input%alpha,fm,1,sm,1)
sm=sm+input%alpha*fm
ELSE
! |vi> = w|vi>
! loop to generate um : um = sm1 + alpha*fm1 - \sum <fm1|w|vi> ui
um(:nmap) = input%alpha * fm1(:nmap) + sm1(:nmap)
um = input%alpha * fm1 + sm1
iread = MIN(mit-1,input%maxiter+1)
DO it = 2, iread
CALL readUVec(input,hybrid,nmap,it-mit,mit,ui)
CALL readVVec(input,hybrid,nmap,it-mit,mit,dfivi,vi)
am(it) = CPP_BLAS_sdot(nmap,vi,1,fm1,1)
am(it) = vi.dot.fm
! calculate um(:) = -am(it)*ui(:) + um(:)
CALL CPP_BLAS_saxpy(nmap,-am(it),ui,1,um,1)
um=um-am(it)*ui
WRITE(6,FMT='(5x,"<vi|w|Fm> for it",i2,5x,f10.6)')it,am(it)
END DO
IF (input%imix.EQ.3) THEN
!****************************************
! broyden's first method
!****************************************
! convolute drho(m) with the metric: |fm1> = w|sm1>
CALL metric(cell,atoms,vacuum,sphhar,input,noco,stars,sym,oneD,&
mmap,nmaph,mapmt,mapvac2,sm1,fm1,l_pot)
! calculate the norm of sm1 : <sm1|w|sm1>
smnorm = CPP_BLAS_sdot(nmap,sm1,1,fm1,1)
! generate vm = alpha*sm1 - \sum <ui|w|sm1> vi
vm(:) = input%alpha * fm1(:)
DO it = 2,iread
CALL readUVec(input,hybrid,nmap,it-mit,mit,ui)
CALL readVVec(input,hybrid,nmap,it-mit,mit,dfivi,vi)
bm = CPP_BLAS_sdot(nmap,ui,1,fm1,1)
! calculate vm(:) = -bm*vi(:) + vm
CALL CPP_BLAS_saxpy(nmap,-bm,vi,1,vm,1)
!write(6,FMT='(5x,"<ui|w|Fm> for it",i2,5x,f10.6)') it, bm
END DO
! complete evaluation of vm
! vmnorm = <um|w|sm1>-<sm1|w|sm1>
vmnorm = CPP_BLAS_sdot(nmap,fm1,1,um,1) - smnorm
! if (vmnorm.lt.tol_10) stop
CALL CPP_BLAS_sscal(nmap,one/vmnorm,vm,1)
ELSE IF (input%imix.EQ.5) THEN
!****************************************
! broyden's second method
!****************************************
! multiply fm1 with metric matrix and store in vm: w |fm1>
CALL metric(cell,atoms,vacuum,sphhar,input,noco,stars,sym,oneD,&
mmap,nmaph,mapmt,mapvac2,fm1,vm,l_pot)
! calculate the norm of fm1 and normalize vm it: vm = wfm1 / <fm1|w|fm1>
vmnorm = one / CPP_BLAS_sdot(nmap,fm1,1,vm,1)
CALL CPP_BLAS_sscal(nmap,vmnorm,vm,1)
ELSE IF (input%imix.EQ.7) THEN
IF (input%imix.EQ.7) THEN
!****************************************
! generalized anderson method
!****************************************
! calculate vm = alpha*wfm1 -\sum <fm1|w|vi> <fi1|w|vi><vi|
! convolute fm1 with the metrik and store in vm
CALL metric(cell,atoms,vacuum,sphhar,input,noco,stars,sym,oneD,&
mmap,nmaph,mapmt,mapvac2,fm1,vm,l_pot)
vm=fm1%apply_metric()
DO it = 2,iread
CALL readVVec(input,hybrid,nmap,it-mit,mit,dfivi,vi)
! calculate vm(:) = -am(it)*dfivi*vi(:) + vm
CALL CPP_BLAS_saxpy(nmap,-am(it)*dfivi,vi,1,vm,1)
vm=vm-am(it)*dfivi*vi
END DO
vmnorm = CPP_BLAS_sdot(nmap,fm1,1,vm,1)
vmnorm=fm1.dot.vm
! if (vmnorm.lt.tol_10) stop
! calculate vm(:) = (1.0/vmnorm)*vm(:)
CALL CPP_BLAS_sscal(nmap,one/vmnorm,vm,1)
vm=(1.0/vmnorm)*vm
! save dfivi(mit) for next iteration
dfivi = vmnorm
ELSE
CALL judft_error("Only generalized Anderson implemented")
END IF
! write um,vm and dfivi on file broyd.?
......@@ -199,12 +123,10 @@ CONTAINS
! update rho(m+1)
! calculate <fm|w|vm>
fmvm = CPP_BLAS_sdot(nmap,vm,1,fm,1)
fmvm = vm.dot.fm
! calculate sm(:) = (1.0-fmvm)*um(:) + sm
CALL CPP_BLAS_saxpy(nmap,one-fmvm,um,1,sm,1)
sm=sm+(one-fmvm)*um
END IF
DEALLOCATE (fm1,sm1,ui,um,vi,vm,am)
END SUBROUTINE broyden
END MODULE m_broyden
......@@ -5,11 +5,17 @@
!--------------------------------------------------------------------------------
module m_distance
contains
subroutine distance(fsm)
SUBROUTINE distance(irank,fsm)
real :: dist(6)
REAL :: dist(6) !1:up,2:down,3:total,4:
CALL fsm_up%alloc()
IF (input%jspin
fsm_up=
dist(:) = 0.0
......
......@@ -15,68 +15,34 @@ MODULE m_stmix
CONTAINS
SUBROUTINE stmix(&
& atoms,input,noco,&
& nmap,nmaph,fsm,sm)
& fsm,fsm_mag,sm)
USE m_types_mixvector
USE m_types
IMPLICIT NONE
TYPE(t_input),INTENT(IN) :: input
TYPE(t_noco),INTENT(IN) :: noco
TYPE(t_atoms),INTENT(IN) :: atoms
! ..
! .. Scalar Arguments ..
INTEGER, INTENT (IN) :: nmaph,nmap
! ..
! .. Array Arguments ..
REAL fsm(:),sm(:)
TYPE(t_input),INTENT(IN) :: input
TYPE(t_noco),INTENT(IN) :: noco
TYPE(t_atoms),INTENT(IN) :: atoms
TYPE(t_mixvector),INTENT(IN) :: fsm,fsm_mag
TYPE(t_mixvector),INTENT(INOUT) :: sm
! ..
! .. Local Scalars ..
INTEGER imap
REAL,PARAMETER:: tol_6=1.0e-6
! ..
!
WRITE (6,FMT='(a)') 'STRAIGHT MIXING'
IF (input%jspins.EQ.1) WRITE (6,FMT='(a,2f10.5)')&
& 'charge density mixing parameter:',input%alpha
IF (input%jspins.EQ.2) WRITE (6,FMT='(a,2f10.5)')&
& 'spin density mixing parameter:',input%alpha*input%spinf
IF ( ABS(input%spinf-1.0e0).LE.tol_6 .OR. input%jspins.EQ.1 ) THEN
! --> perform simple mixing
!
! sm1 = sm + alpha * F(sm)
sm = sm + input%alpha*fsm
sm(:nmap) = sm(:nmap) + input%alpha*fsm(:nmap)
IF ( ABS(input%spinf-1.0e0).LE.tol_6 .OR. input%jspins.EQ.1 ) THEN
! Done with sm1 = sm + alpha * F(sm)
!No spin
RETURN
ELSE
sm = sm + input%alpha/2.0*(input%spinf-1.0)*fsm_mag
! -->perform simple mixing with the mixing parameters
! for charge and spin
!
! sm1+/_ = (sm+/_) + alpha* F(sm)
! +/-0.5alpha(spinf-1)( F(sm+) + F(sm-) )
DO imap = 1,nmaph
sm(imap) = sm(imap) + input%alpha*fsm(imap) &
& + input%alpha/2.0*(input%spinf-1.0)*(fsm(imap) - fsm(imap+nmaph))
ENDDO
DO imap = nmaph+1,2*nmaph
sm(imap) = sm(imap) + input%alpha*fsm(imap) &
& + input%alpha/2.0*(input%spinf-1.0)*(fsm(imap) - fsm(imap-nmaph))
ENDDO
IF (noco%l_noco) THEN
DO imap = 2*nmaph+1, nmap - 98*input%jspins*atoms%n_u
sm(imap) = sm(imap) + input%alpha*input%spinf*fsm(imap)
ENDDO
ENDIF
IF ( atoms%n_u > 0 ) THEN
DO imap = nmap - 98*input%jspins*atoms%n_u + 1, nmap - 98*atoms%n_u
sm(imap) = sm(imap) + input%alpha*fsm(imap) &
& + input%alpha/2.0*(input%spinf-1.0)*(fsm(imap) - fsm(imap+98*atoms%n_u))
ENDDO
DO imap = nmap - 98*atoms%n_u + 1, nmap
sm(imap) = sm(imap) + input%alpha*fsm(imap) &
& + input%alpha/2.0*(input%spinf-1.0)*(fsm(imap) - fsm(imap-98*atoms%n_u))
ENDDO
ENDIF
END IF
END SUBROUTINE stmix
......
This diff is collapsed.
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment