Commit 99e8c98f authored by Daniel Wortmann's avatar Daniel Wortmann

Bugfix in magnetic distance

parent bd43f992
......@@ -98,7 +98,7 @@ contains
maxiter=merge(1,input%maxiter,input%imix==0)
CALL mixing_history(input%imix,maxiter,inden,outden,sm,fsm,it)
CALL distance(mpi%irank,cell%vol,input%jspins,fsm(it),sm(it),inDen%iter,outDen,results,fsm_Mag)
CALL distance(mpi%irank,cell%vol,input%jspins,fsm(it),inDen,outDen,results,fsm_Mag)
! KERKER PRECONDITIONER
IF( input%preconditioning_param /= 0 ) call kerker(field, DIMENSION, mpi, &
......
......@@ -5,16 +5,16 @@
!--------------------------------------------------------------------------------
module m_distance
contains
SUBROUTINE distance(irank,vol,jspins,fsm,sm,iter,outden,results,fsm_mag)
SUBROUTINE distance(irank,vol,jspins,fsm,inden,outden,results,fsm_mag)
use m_types
use m_types_mixvector
use m_xmlOutput
implicit none
integer,intent(in) :: irank,jspins,iter
integer,intent(in) :: irank,jspins
real,intent(in) :: vol
type(t_mixvector),INTENT(IN) :: fsm,sm
TYPE(t_potden),INTENT(IN) :: outden
type(t_mixvector),INTENT(IN) :: fsm
TYPE(t_potden),INTENT(IN) :: inden,outden
TYPE(t_results),INTENT(INOUT) :: results
type(t_mixvector),INTENT(OUT) :: fsm_mag
......@@ -28,7 +28,8 @@ contains
CALL fsm_mag%alloc()
! calculate Magnetisation-difference
CALL fsm_mag%from_density(outden,swapspin=.TRUE.)
fsm_mag=fsm_mag-sm
call fmMet%from_density(inden,swapspin=.true.)
fsm_mag=fsm_mag-fmMet
ENDIF
! Apply metric w to fsm and store in fmMet: w |fsm>
fmMet=fsm%apply_metric()
......@@ -39,8 +40,7 @@ contains
END DO
IF (SIZE(outden%pw,2)>2) dist(6) = fsm%multiply_dot_mask(fmMet,(/.TRUE.,.TRUE.,.TRUE.,.FALSE./),3)
IF (jspins.EQ.2) THEN
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(3) = fsMet%multiply_dot_mask(fsm_mag,(/.true.,.true.,.true.,.false./),1)
dist(4) = dist(1) + dist(2) + 2.0e0*dist(3)
dist(5) = dist(1) + dist(2) - 2.0e0*dist(3)
ENDIF
......@@ -56,10 +56,10 @@ contains
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))
WRITE ( 6,FMT=7900) js,inDen%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))
IF (abs(dist(6))>1E-15) WRITE (6,FMT=7900) 3,inDen%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)| +/_
......@@ -69,8 +69,8 @@ contains
(/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))
WRITE ( 6,FMT=8000) inDen%iter,1000*SQRT(ABS(dist(4)/vol))
WRITE ( 6,FMT=8010) inDen%iter,1000*SQRT(ABS(dist(5)/vol))
!dist/vol should always be >= 0 ,
!but for dist=0 numerically you might obtain dist/vol < 0
......
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