distance.F90 3.73 KB
Newer Older
1 2 3 4 5 6 7
!--------------------------------------------------------------------------------
! 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_distance
contains
Daniel Wortmann's avatar
Daniel Wortmann committed
8
  SUBROUTINE distance(irank,vol,jspins,fsm,sm,iter,outden,results,fsm_mag)
9 10 11
    use m_types
    use m_types_mixvector
    use m_xmlOutput
12
 
13 14 15 16 17 18
    implicit none
    integer,intent(in)             :: irank,jspins,iter
    real,intent(in)                :: vol
    type(t_mixvector),INTENT(IN)   :: fsm,sm
    TYPE(t_potden),INTENT(IN)      :: outden
    TYPE(t_results),INTENT(INOUT)  :: results
Daniel Wortmann's avatar
Daniel Wortmann committed
19
    type(t_mixvector),INTENT(OUT)   :: fsm_mag
Daniel Wortmann's avatar
Daniel Wortmann committed
20
    
21 22
    integer         ::js
    REAL            :: dist(6) !1:up,2:down,3:spinoff,4:total,5:magnet,6:noco
Daniel Wortmann's avatar
Daniel Wortmann committed
23
    TYPE(t_mixvector)::fmMet
24
    character(len=100)::attributes(2)
25
    
26
    CALL fmMet%alloc()
Daniel Wortmann's avatar
Daniel Wortmann committed
27 28 29 30 31 32
    IF (jspins==2) THEN
       CALL fsm_mag%alloc()
       ! calculate Magnetisation-difference
       CALL fsm_mag%from_density(outden,swapspin=.TRUE.)
       fsm_mag=fsm_mag-sm
    ENDIF
33
    ! Apply metric w to fsm and store in fmMet:  w |fsm>
34 35 36
    fmMet=fsm%apply_metric()
  
    dist(:) = 0.0
Daniel Wortmann's avatar
Daniel Wortmann committed
37
    DO js = 1,jspins
38 39
       dist(js) = fsm%multiply_dot_mask(fmMet,(/.true.,.true.,.true.,.false./),js)
    END DO
Daniel Wortmann's avatar
Daniel Wortmann committed
40
    IF (SIZE(outden%pw,2)>2) dist(6) = fsm%multiply_dot_mask(fmMet,(/.TRUE.,.TRUE.,.TRUE.,.FALSE./),3)
Daniel Wortmann's avatar
Daniel Wortmann committed
41
    IF (jspins.EQ.2) THEN
42 43 44 45
       dist(3) = fsm_mag%multiply_dot_mask(fmMet,(/.true.,.true.,.true.,.false./),1)+&
            fsm_mag%multiply_dot_mask(fmMet,(/.true.,.true.,.true.,.false./),2)
       dist(4) = dist(1) + dist(2) + 2.0e0*dist(3)
       dist(5) = dist(1) + dist(2) - 2.0e0*dist(3)
Daniel Wortmann's avatar
Daniel Wortmann committed
46 47
    ENDIF
    
48 49
    results%last_distance=maxval(1000*SQRT(ABS(dist/vol)))
    if (irank>1) return
50

51 52 53
    !calculate the distance of charge densities for each spin
    CALL openXMLElement('densityConvergence',(/'units'/),(/'me/bohr^3'/))
    
Daniel Wortmann's avatar
Daniel Wortmann committed
54
    DO js = 1,jspins         
55 56 57 58 59 60 61 62 63 64 65 66
       attributes = ''
       WRITE(attributes(1),'(i0)') js
       WRITE(attributes(2),'(f20.10)') 1000*SQRT(ABS(dist(js)/vol))
       CALL writeXMLElementForm('chargeDensity',(/'spin    ','distance'/),attributes,reshape((/4,8,1,20/),(/2,2/)))
       WRITE ( 6,FMT=7900) js,iter,1000*SQRT(ABS(dist(js)/vol))
    END DO
    
    IF (abs(dist(6))>1E-15) WRITE (6,FMT=7900) 3,iter,1000*SQRT(ABS(dist(6)/vol))
    
    !calculate the distance of total charge and spin density
    !|rho/m(o) - rho/m(i)| = |rh1(o) -rh1(i)|+ |rh2(o) -rh2(i)| +/_
    !                        +/_2<rh2(o) -rh2(i)|rh1(o) -rh1(i)>
Daniel Wortmann's avatar
Daniel Wortmann committed
67
    IF (jspins.EQ.2) THEN
68 69 70 71 72 73 74 75 76 77 78 79
       CALL writeXMLElementFormPoly('overallChargeDensity',(/'distance'/),&
            (/1000*SQRT(ABS(dist(4)/vol))/),reshape((/10,20/),(/1,2/)))
       CALL writeXMLElementFormPoly('spinDensity',(/'distance'/),&
            (/1000*SQRT(ABS(dist(5)/vol))/),reshape((/19,20/),(/1,2/)))
       WRITE ( 6,FMT=8000) iter,1000*SQRT(ABS(dist(4)/vol))
       WRITE ( 6,FMT=8010) iter,1000*SQRT(ABS(dist(5)/vol))
       
       !dist/vol should always be >= 0 ,
       !but for dist=0 numerically you might obtain dist/vol < 0
       !(e.g. when calculating non-magnetic systems with jspins=2).
    END IF
    CALL closeXMLElement('densityConvergence')
80 81


82 83 84 85 86 87 88
7900  FORMAT (/,'---->    distance of charge densities for spin ',i2,'                 it=',i5,':',f13.6,' me/bohr**3')
8000 FORMAT (/,'---->    distance of charge densities for it=',i5,':', f13.6,' me/bohr**3')
8010 FORMAT (/,'---->    distance of spin densities for it=',i5,':', f13.6,' me/bohr**3')
8020 FORMAT (4d25.14)
8030  FORMAT (10i10)
  end SUBROUTINE distance
end module m_distance