Commit 2b5c9ee1 authored by Daniel Wortmann's avatar Daniel Wortmann

Bugfixes to make new mixer work, also added IO again.

parent 755fc78d
......@@ -59,7 +59,7 @@ kpoints/tetcon.f kpoints/kvecon.f init/boxdim.f global/radsra.f math/intgr.F glo
)
set(inpgen_F90 ${inpgen_F90} global/constants.f90 io/xsf_io.f90
eigen/orthoglo.F90 juDFT/usage_data.F90 math/ylm4.F90
eigen/orthoglo.F90 juDFT/usage_data.F90 math/ylm4.F90 mpi/mpi_bc_tool.F90
global/sort.f90 global/chkmt.f90 inpgen/inpgen.f90 inpgen/set_inp.f90 inpgen/inpgen_help.f90 io/rw_inp.f90 juDFT/juDFT.F90 global/find_enpara.f90
inpgen/closure.f90 inpgen/inpgen_arguments.F90
juDFT/info.F90 juDFT/stop.F90 juDFT/args.F90 juDFT/time.F90 juDFT/init.F90 juDFT/sysinfo.F90 juDFT/string.f90 io/w_inpXML.f90 kpoints/julia.f90 global/utility.F90
......
......@@ -20,7 +20,6 @@ io/cdnpot_io_common.F90
io/cdn_io.F90
io/pot_io.F90
io/banddos_io.F90
io/broyd_io.F90
io/rw_inp.f90
io/rw_noco.f90
io/r_inpXML.F90
......
......@@ -41,7 +41,7 @@ CONTAINS
USE m_fleur_init
USE m_optional
USE m_cdn_io
USE m_broyd_io
USE m_mixing_history
USE m_qfix
USE m_vgen
USE m_writexcstuff
......@@ -212,7 +212,7 @@ CONTAINS
cell,oneD,enpara,results,sym,xcpot,vTot,iter,iterHF)
END SELECT
IF(hybrid%l_calhf) THEN
CALL system("rm broyd*")
call mixing_history_reset(mpi)
iter = 0
END IF
ENDIF
......@@ -412,8 +412,9 @@ CONTAINS
! mix input and output densities
CALL timestart("mixing")
CALL mix(field2,dimension,mpi,stars,atoms,sphhar,vacuum,input,&
sym,cell,noco,oneD,archiveType,inDen,outDen,results)
CALL mix_charge(field2,DIMENSION,mpi,(iter==input%itmax.OR.judft_was_argument("-mix_io")),&
stars,atoms,sphhar,vacuum,input,&
sym,cell,noco,oneD,archiveType,inDen,outDen,results)
CALL timestop("mixing")
IF(mpi%irank == 0) THEN
......@@ -485,7 +486,7 @@ CONTAINS
CLOSE(2)
PRINT *,"qfix set to F"
ENDIF
CALL resetBroydenHistory()
call mixing_history_reset(mpi)
ENDIF
CALL juDFT_end(" GEO new inp.xml created ! ",mpi%irank)
END SUBROUTINE priv_geo_end
......
......@@ -13,7 +13,7 @@ MODULE m_fleur_arguments
CHARACTER(len=200) :: values
END TYPE t_fleur_param
INTEGER,PARAMETER:: no_params=24
INTEGER,PARAMETER:: no_params=25
TYPE(t_fleur_param) :: fleur_param(no_params)=(/&
!Input options
t_fleur_param(0,"-toXML","Convert an old 'inp' file into the new XML format",""),&
......@@ -58,6 +58,7 @@ MODULE m_fleur_arguments
t_fleur_param(0,"-trace","Try to generate a stacktrace in case of an error",""),&
t_fleur_param(0,"-debugtime","Write the start/stop of all timers to the console",""),&
!Output
t_fleur_param(0,"-mix_io","Do not store mixing history in memory but do IO in each iteration",""),&
t_fleur_param(0,"-no_out","Do not open the 'out' file but write to stdout",""),&
t_fleur_param(0,"-genEnpara","Generate an 'enpara' file for the energy parameters",""),&
t_fleur_param(0,"-kpts_gw","add alternative k point set for GW in all outputs for the XML input file",""),&
......
......@@ -28,7 +28,7 @@
USE m_setupMPI
USE m_cdn_io
USE m_fleur_info
USE m_broyd_io
USE m_mixing_history
USE m_checks
USE m_prpqfftmap
USE m_writeOutHeader
......@@ -531,7 +531,7 @@
CALL results%init(dimension,input,atoms,kpts,noco)
IF (mpi%irank.EQ.0) THEN
IF(input%gw.NE.0) CALL resetBroydenHistory()
IF(input%gw.NE.0) CALL mixing_history_reset(mpi)
CALL setStartingDensity(noco%l_noco)
END IF
......
......@@ -15,9 +15,9 @@ MODULE m_mix
contains
SUBROUTINE mix( field, DIMENSION, mpi, &
stars, atoms, sphhar, vacuum, input, sym, cell, noco, &
oneD, archiveType, inDen, outDen, results )
SUBROUTINE mix_charge( field, DIMENSION, mpi, l_writehistory,&
stars, atoms, sphhar, vacuum, input, sym, cell, noco, &
oneD, archiveType, inDen, outDen, results )
use m_juDFT
use m_constants
......@@ -38,17 +38,18 @@ contains
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_stars),TARGET,INTENT(in) :: stars
TYPE(t_cell),TARGET,INTENT(in) :: cell
TYPE(t_sphhar),TARGET,INTENT(in) :: sphhar
type(t_field), intent(inout) :: field
type(t_dimension), intent(in) :: dimension
type(t_mpi), intent(in) :: mpi
type(t_atoms), intent(inout) :: atoms !n_u is modified temporarily
TYPE(t_atoms),TARGET,INTENT(in) :: atoms
type(t_potden), intent(inout) :: outDen
type(t_results), intent(inout) :: results
type(t_potden), intent(inout) :: inDen
integer, intent(in) :: archiveType
LOGICAL, INTENT(IN) :: l_writehistory
real :: fix
type(t_potden) :: resDen, vYukawa
......@@ -57,68 +58,84 @@ contains
LOGICAL :: l_densitymatrix
INTEGER :: it,maxiter
MPI0_a: IF( mpi%irank == 0 ) THEN
!determine type of mixing:
!imix=0:straight, imix=o broyden first, imix=5:broyden second
!imix=:generalozed anderson mixing
select case( input%imix )
case( 0 )
WRITE( 6, fmt='(a,2f10.5)' ) 'STRAIGHT MIXING',input%alpha
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
case( 3 )
!determine type of mixing:
!imix=0:straight, imix=o broyden first, imix=5:broyden second
!imix=:generalozed anderson mixing
select case( input%imix )
case( 0 )
WRITE( 6, fmt='(a,2f10.5)' ) 'STRAIGHT MIXING',input%alpha
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
case( 3 )
write( 6, fmt='(a,f10.5)' ) 'BROYDEN FIRST MIXING',input%alpha
case( 5 )
case( 5 )
write( 6, fmt='(a,f10.5)' ) 'BROYDEN SECOND MIXING',input%alpha
case( 7 )
case( 7 )
write( 6, fmt='(a,f10.5)' ) 'ANDERSON GENERALIZED',input%alpha
case default
case default
call juDFT_error( "mix: input%imix =/= 0,3,5,7 ", calledby ="mix" )
end select
if ( input%jspins == 2 .and. input%imix /= 0 ) then
write( 6, '(''WARNING : for QUASI-NEWTON METHODS SPINF=1'')' )
end if
ENDIF MPI0_a
l_densitymatrix=.FALSE.
IF (atoms%n_u>0) THEN
l_densitymatrix=.NOT.input%ldaulinmix
CALL u_mix(input,atoms,inDen%mmpMat,outDen%mmpMat)
IF (ALL(inDen%mmpMat==0.0)) THEN
l_densitymatrix=.FALSE.
inDen%mmpMat=outDen%mmpMat
if (mpi%irank.ne.0) inden%mmpmat=0.0
ENDIF
ENDIF
CALL mixvector_init(mpi%mpi_comm,l_densitymatrix,oneD,input,vacuum,noco,sym,stars,cell,sphhar,atoms)
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),inDen,outDen,results,fsm_Mag)
end select
if ( input%jspins == 2 .and. input%imix /= 0 ) then
write( 6, '(''WARNING : for QUASI-NEWTON METHODS SPINF=1'')' )
end if
ENDIF MPI0_a
l_densitymatrix=.FALSE.
IF (atoms%n_u>0) THEN
l_densitymatrix=.NOT.input%ldaulinmix
IF (mpi%irank==0) CALL u_mix(input,atoms,inDen%mmpMat,outDen%mmpMat)
IF (ALL(inDen%mmpMat==0.0)) THEN
l_densitymatrix=.FALSE.
inDen%mmpMat=outDen%mmpMat
if (mpi%irank.ne.0) inden%mmpmat=0.0
ENDIF
ENDIF
CALL mixvector_init(mpi%mpi_comm,l_densitymatrix,oneD,input,vacuum,noco,sym,stars,cell,sphhar,atoms)
CALL mixing_history_open(mpi,input%maxiter)
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),inDen,outDen,results,fsm_Mag)
! KERKER PRECONDITIONER
IF( input%preconditioning_param /= 0 ) call kerker(field, DIMENSION, mpi, &
stars, atoms, sphhar, vacuum, input, sym, cell, noco, &
oneD, inDen, outDen, fsm(it) )
IF( input%preconditioning_param /= 0 ) THEN
CALL kerker(field, DIMENSION, mpi, &
stars, atoms, sphhar, vacuum, input, sym, cell, noco, &
oneD, inDen, outDen, fsm(it) )
!Store modified density in history
CALL mixing_history_store(fsm(it))
END IF
!mixing of the densities
if(input%imix==0.or.it==1) CALL stmix(atoms,input,noco,fsm(it),fsm_mag,sm(it))
!if(it>1.and.input%imix==9) CALL pulay(input%alpha,fsm,sm)
if(it>1.and.(input%imix==3.or.input%imix==5.or.input%imix==7)) Call broyden(input%alpha,fsm,sm)
!initiatlize mixed density and extract it
!extracte mixed density
inDen%pw=0.0;inDen%mt=0.0
IF (ALLOCATED(inDen%vacz)) inden%vacz=0.0
IF (ALLOCATED(inDen%vacxy)) inden%vacxy=0.0
IF (ALLOCATED(inDen%mmpMat).AND.l_densitymatrix) inden%mmpMat=0.0
CALL sm(it)%to_density(inDen)
IF (atoms%n_u>0.AND..NOT.l_densitymatrix.AND..NOT.input%ldaulinmix) THEN
!No density matrix was present
!but is now created...
inden%mmpMAT=outden%mmpMat
CALL mixing_history_reset(mpi)
CALL mixvector_reset()
ENDIF
!fix charge of the new density
CALL qfix(mpi,stars,atoms,sym,vacuum, sphhar,input,cell,oneD,inDen,noco%l_noco,.FALSE.,.FALSE., fix)
IF (mpi%irank==0) CALL qfix(mpi,stars,atoms,sym,vacuum, sphhar,input,cell,oneD,inDen,noco%l_noco,.FALSE.,.FALSE., fix)
IF(vacuum%nvac.EQ.1) THEN
inDen%vacz(:,2,:) = inDen%vacz(:,1,:)
......@@ -128,30 +145,25 @@ contains
inDen%vacxy(:,:,2,:) = inDen%vacxy(:,:,1,:)
END IF
END IF
IF (atoms%n_u>0.AND..NOT.l_densitymatrix.AND..NOT.input%ldaulinmix) THEN
!No density matrix was present
!but is now created...
CALL mixing_history_reset()
CALL mixvector_reset()
ENDIF
!write out mixed density
CALL writeDensity(stars,vacuum,atoms,cell,sphhar,input,sym,oneD,archiveType,CDN_INPUT_DEN_const,&
IF (mpi%irank==0) CALL writeDensity(stars,vacuum,atoms,cell,sphhar,input,sym,oneD,archiveType,CDN_INPUT_DEN_const,&
1,results%last_distance,results%ef,.TRUE.,inDen)
#ifdef CPP_HDF
IF (judft_was_argument("-last_extra")) THEN
IF (mpi%irank==0.and.judft_was_argument("-last_extra")) THEN
CALL system("rm cdn_last.hdf")
CALL writeDensity(stars,vacuum,atoms,cell,sphhar,input,sym,oneD,archiveType,CDN_INPUT_DEN_const,&
1,results%last_distance,results%ef,.TRUE.,inDen,'cdn_last')
END IF
#endif
inDen%iter = inDen%iter + 1
IF (l_writehistory.AND.input%imix.NE.0) CALL mixing_history_close(mpi)
END SUBROUTINE mix
END SUBROUTINE mix_charge
END MODULE m_mix
......@@ -14,9 +14,9 @@ contains
integer,intent(in) :: irank,jspins
real,intent(in) :: vol
type(t_mixvector),INTENT(IN) :: fsm
TYPE(t_potden),INTENT(IN) :: inden,outden
TYPE(t_potden),INTENT(INOUT) :: inden,outden
TYPE(t_results),INTENT(INOUT) :: results
type(t_mixvector),INTENT(OUT) :: fsm_mag
type(t_mixvector),INTENT(OUT) :: fsm_mag
integer ::js
REAL :: dist(6) !1:up,2:down,3:spinoff,4:total,5:magnet,6:noco
......@@ -46,8 +46,7 @@ contains
ENDIF
results%last_distance=maxval(1000*SQRT(ABS(dist/vol)))
if (irank>1) return
if (irank>0) return
!calculate the distance of charge densities for each spin
CALL openXMLElement('densityConvergence',(/'units'/),(/'me/bohr^3'/))
......
......@@ -4,79 +4,78 @@
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------
MODULE m_kerker
contains
CONTAINS
SUBROUTINE kerker(field, DIMENSION, mpi, &
stars, atoms, sphhar, vacuum, input, sym, cell, noco, &
oneD, inDen, outDen, precon_v )
!Implementation of the Kerker preconditioner by M.Hinzen
USE m_vgen_coulomb
USE m_VYukawaFilm
USE m_juDFT
USE m_qfix
USE m_types
USE m_types_mixvector
USE m_constants
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_field), INTENT(inout) :: field
TYPE(t_dimension), INTENT(in) :: DIMENSION
TYPE(t_mpi), INTENT(in) :: mpi
TYPE(t_atoms), INTENT(in) :: atoms
TYPE(t_potden), INTENT(inout) :: outDen
TYPE(t_potden), INTENT(in) :: inDen
TYPE(t_mixvector), INTENT(INOUT) :: precon_v
!Locals
type(t_potden) :: resDen, vYukawa
real :: fix
integer :: lh,n
CALL resDen%init( stars, atoms, sphhar, vacuum, noco, input%jspins, POTDEN_TYPE_DEN )
CALL vYukawa%init( stars, atoms, sphhar, vacuum, noco, input%jspins, 4 )
MPI0_b: IF( mpi%irank == 0 ) THEN
CALL resDen%subPotDen( outDen, inDen )
IF( input%jspins == 2 ) CALL resDen%SpinsToChargeAndMagnetisation()
END IF MPI0_b
SUBROUTINE kerker(field, DIMENSION, mpi, &
stars, atoms, sphhar, vacuum, input, sym, cell, noco, &
oneD, inDen, outDen, precon_v )
!Implementation of the Kerker preconditioner by M.Hinzen
USE m_vgen_coulomb
USE m_VYukawaFilm
USE m_juDFT
USE m_qfix
USE m_types
USE m_types_mixvector
USE m_constants
USE m_mpi_bc_potden
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_field), INTENT(inout) :: field
TYPE(t_dimension), INTENT(in) :: DIMENSION
TYPE(t_mpi), INTENT(in) :: mpi
TYPE(t_atoms), INTENT(in) :: atoms
TYPE(t_potden), INTENT(inout) :: outDen
TYPE(t_potden), INTENT(in) :: inDen
TYPE(t_mixvector), INTENT(INOUT) :: precon_v
!Locals
type(t_potden) :: resDen, vYukawa
real :: fix
integer :: lh,n
CALL resDen%init( stars, atoms, sphhar, vacuum, noco, input%jspins, POTDEN_TYPE_DEN )
CALL vYukawa%init( stars, atoms, sphhar, vacuum, noco, input%jspins, 4 )
MPI0_b: IF( mpi%irank == 0 ) THEN
CALL resDen%subPotDen( outDen, inDen )
IF( input%jspins == 2 ) CALL resDen%SpinsToChargeAndMagnetisation()
END IF MPI0_b
#ifdef CPP_MPI
CALL mpi_bc_potden( mpi, stars, sphhar, atoms, input, vacuum, oneD, noco, resDen )
CALL mpi_bc_potden( mpi, stars, sphhar, atoms, input, vacuum, oneD, noco, resDen )
#endif
IF ( .NOT. input%film ) THEN
CALL vgen_coulomb( 1, mpi, DIMENSION, oneD, input, field, vacuum, sym, stars, cell, &
sphhar, atoms, resDen, vYukawa )
ELSE
vYukawa%iter = resDen%iter
CALL VYukawaFilm( stars, vacuum, cell, sym, input, mpi, atoms, sphhar, DIMENSION, oneD, resDen, &
vYukawa )
END IF
MPI0_c: IF( mpi%irank == 0 ) THEN
IF( input%preconditioning_param /= 0 ) THEN
resDen%pw(1:stars%ng3,1) = resDen%pw(1:stars%ng3,1) - input%preconditioning_param ** 2 / fpi_const * vYukawa%pw(1:stars%ng3,1)
DO n = 1, atoms%ntype
DO lh = 0, sphhar%nlhd
resDen%mt(1:atoms%jri(n),lh,n,1) = resDen%mt(1:atoms%jri(n),lh,n,1) &
- input%preconditioning_param ** 2 / fpi_const &
* vYukawa%mt(1:atoms%jri(n),lh,n,1) * atoms%rmsh(1:atoms%jri(n),n) ** 2
END DO
END DO
resDen%vacz = resDen%vacz - input%preconditioning_param ** 2 / fpi_const * vYukawa%vacz
resDen%vacxy = resDen%vacxy - input%preconditioning_param ** 2 / fpi_const * vYukawa%vacxy
IF( input%jspins == 2 ) CALL resDen%ChargeAndMagnetisationToSpins()
! fix the preconditioned density
CALL outDen%addPotDen( resDen, inDen )
CALL qfix(mpi,stars, atoms, sym, vacuum, sphhar, input, cell, oneD, outDen, noco%l_noco, .FALSE., .TRUE., fix )
CALL resDen%subPotDen( outDen, inDen )
call precon_v%from_density(resden)
END IF
! end of preconditioner
END IF MPI0_c
END SUBROUTINE kerker
end MODULE m_kerker
IF ( .NOT. input%film ) THEN
CALL vgen_coulomb( 1, mpi, DIMENSION, oneD, input, field, vacuum, sym, stars, cell, &
sphhar, atoms, resDen, vYukawa )
ELSE
vYukawa%iter = resDen%iter
CALL VYukawaFilm( stars, vacuum, cell, sym, input, mpi, atoms, sphhar, DIMENSION, oneD, resDen, &
vYukawa )
END IF
MPI0_c: IF( mpi%irank == 0 ) THEN
resDen%pw(1:stars%ng3,1) = resDen%pw(1:stars%ng3,1) - input%preconditioning_param ** 2 / fpi_const * vYukawa%pw(1:stars%ng3,1)
DO n = 1, atoms%ntype
DO lh = 0, sphhar%nlhd
resDen%mt(1:atoms%jri(n),lh,n,1) = resDen%mt(1:atoms%jri(n),lh,n,1) &
- input%preconditioning_param ** 2 / fpi_const &
* vYukawa%mt(1:atoms%jri(n),lh,n,1) * atoms%rmsh(1:atoms%jri(n),n) ** 2
END DO
END DO
resDen%vacz = resDen%vacz - input%preconditioning_param ** 2 / fpi_const * vYukawa%vacz
resDen%vacxy = resDen%vacxy - input%preconditioning_param ** 2 / fpi_const * vYukawa%vacxy
IF( input%jspins == 2 ) CALL resDen%ChargeAndMagnetisationToSpins()
! fix the preconditioned density
CALL outDen%addPotDen( resDen, inDen )
CALL qfix(mpi,stars, atoms, sym, vacuum, sphhar, input, cell, oneD, outDen, noco%l_noco, .FALSE., .TRUE., fix )
! CALL resDen%subPotDen( outDen, inDen )
! CALL precon_v%from_density(resden)
END IF MPI0_c
! end of preconditioner
END SUBROUTINE kerker
END MODULE m_kerker
......@@ -4,17 +4,74 @@
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------
MODULE m_mixing_history
use m_types_mixvector
implicit none
integer:: iter_stored=0
type(t_mixvector),allocatable::sm_store(:),fsm_store(:)
contains
USE m_types_mixvector
IMPLICIT NONE
PRIVATE
INTEGER:: iter_stored=0
TYPE(t_mixvector),ALLOCATABLE::sm_store(:),fsm_store(:)
PUBLIC :: mixing_history,mixing_history_reset,mixing_history_store
PUBLIC :: mixing_history_open,mixing_history_close
CONTAINS
SUBROUTINE mixing_history_open(mpi,maxiter)
USE m_types,ONLY:t_mpi
INTEGER,INTENT(IN) :: maxiter
TYPE(t_mpi),INTENT(in):: mpi
CHARACTER(len=20):: filename
LOGICAL :: l_fileexist
INTEGER :: n
IF (iter_stored>0) RETURN ! History in memory found, no need to do IO
IF (mpi%isize>1) THEN
WRITE(filename,'(a,i0)') "mixing_history.",mpi%irank
ELSE
filename="mixing_history"
ENDIF
INQUIRE(file=filename,exist=l_fileexist)
IF (.NOT.l_fileexist) RETURN !No previous data
OPEN(888,file=filename,status='old',form='unformatted')
READ(888) iter_stored
IF (.NOT.ALLOCATED(sm_store)) ALLOCATE(sm_store(maxiter),fsm_store(maxiter))
DO n=1,MIN(iter_stored,maxiter)
READ(888) sm_store(n)
READ(888) fsm_store(n)
ENDDO
CLOSE(888)
END SUBROUTINE mixing_history_open
SUBROUTINE mixing_history_close(mpi)
USE m_types,ONLY:t_mpi
TYPE(t_mpi),INTENT(in):: mpi
CHARACTER(len=20):: filename
INTEGER :: n
IF (iter_stored==0) RETURN ! Nothing found to be stored
IF (mpi%isize>1) THEN
WRITE(filename,'(a,i0)') "mixing_history.",mpi%irank
ELSE
filename="mixing_history"
ENDIF
OPEN(888,file=filename,form='unformatted',status='replace')
WRITE(888) iter_stored
DO n=1,iter_stored
WRITE(888) sm_store(n)
WRITE(888) fsm_store(n)
ENDDO
CLOSE(888)
DEALLOCATE(sm_store,fsm_store)
iter_stored=0
END SUBROUTINE mixing_history_close
SUBROUTINE mixing_history(imix,maxiter,inden,outden,sm,fsm,it)
use m_types
USE m_types
implicit none
INTEGER,INTENT(in)::imix,maxiter
type(t_potden),intent(in)::inden,outden
type(t_potden),intent(inout)::inden,outden
type(t_mixvector),ALLOCATABLE::sm(:),fsm(:)
INTEGER,INTENT(out)::it
......@@ -47,9 +104,17 @@ contains
endif
end subroutine mixing_history
SUBROUTINE mixing_history_reset()
SUBROUTINE mixing_history_reset(mpi)
USE m_types,ONLY:t_mpi
IMPLICIT NONE
TYPE(t_mpi),INTENT(in)::mpi
iter_stored=0
IF (mpi%irank==0) CALL system('rm mixing_history*')
END SUBROUTINE mixing_history_reset
SUBROUTINE mixing_history_store(fsm)
IMPLICIT NONE
TYPE(t_mixvector),INTENT(IN)::fsm
IF (iter_stored>0) fsm_store(iter_stored)=fsm
END SUBROUTINE mixing_history_store
end MODULE m_mixing_history
......@@ -7,10 +7,12 @@ MODULE m_types_mixvector
!TODO!!!
! LDA+U
! Noco (third spin)
use m_types
implicit none
#ifdef CPP_MPI
include 'mpif.h'
#endif
PRIVATE
!Here we store the pointers used for metric
TYPE(t_stars),POINTER :: stars
......@@ -49,6 +51,10 @@ MODULE m_types_mixvector
PROCEDURE :: to_density=>mixvector_to_density
PROCEDURE :: apply_metric=>mixvector_metric
PROCEDURE :: multiply_dot_mask
PROCEDURE :: read_unformatted
PROCEDURE :: write_unformatted
GENERIC :: READ(UNFORMATTED) =>read_unformatted
GENERIC :: WRITE(UNFORMATTED) =>write_unformatted
END TYPE t_mixvector
INTERFACE OPERATOR (*)
......@@ -70,6 +76,35 @@ MODULE m_types_mixvector
CONTAINS
SUBROUTINE READ_unformatted(this,unit,iostat,iomsg)
IMPLICIT NONE
CLASS(t_mixvector),INTENT(INOUT)::this
INTEGER,INTENT(IN)::unit
INTEGER,INTENT(OUT)::iostat
CHARACTER(len=*),INTENT(INOUT)::iomsg
CALL this%alloc()
IF (pw_here) READ(unit) this%vec_pw
IF (mt_here) READ(unit) this%vec_mt
IF (vac_here) READ(unit) this%vec_vac
IF (misc_here) READ(unit) this%vec_misc
END SUBROUTINE READ_unformatted
SUBROUTINE write_unformatted(this,unit,iostat,iomsg)
IMPLICIT NONE
CLASS(t_mixvector),INTENT(IN)::this
INTEGER,INTENT(IN)::unit
INTEGER,INTENT(OUT)::iostat
CHARACTER(len=*),INTENT(INOUT)::iomsg
IF (pw_here) WRITE(unit) this%vec_pw
IF (mt_here) WRITE(unit) this%vec_mt
IF (vac_here) WRITE(unit) this%vec_vac
IF (misc_here) WRITE(unit) this%vec_misc
END SUBROUTINE write_unformatted
SUBROUTINE mixvector_reset()
IMPLICIT NONE
atoms=>NULL()
......@@ -369,19 +404,18 @@ CONTAINS
SUBROUTINE init_storage_mpi(mpi_comm)
IMPLICIT NONE
INTEGER,INTENT(in):: mpi_comm
INTEGER :: irank,isize,err,js,new_comm
mix_mpi_comm=mpi_comm
#ifdef CPP_MPI
INCLUDE 'mpif.h'
INTEGER :: irank,isize,err,js,new_comm
CALL mpi_comm_rank(mpi_comm,irank,err)
CALL mpi_comm_size(mpi_comm,isize,err)
IF (isize==1) RETURN !No parallelization
js=MERGE(jspins,3,.NOT.l_noco)!distribute spins
js=MAX(js,isize)
CALL MPI_COMM_SPLIT(mpi_comm,MOD(irank,js),irank,new_comm)
spin_here=(/MOD(irank,js)==0,MOD(irank,js)==1,(isize=2.AND.irank==0).or.MOD(irank,js)==2/)
js=MIN(js,isize)
CALL MPI_COMM_SPLIT(mpi_comm,MOD(irank,js),irank,new_comm,err)
spin_here=(/MOD(irank,js)==0,MOD(irank,js)==1,(isize==2.AND.irank==0).or.MOD(irank,js)==2/)
CALL mpi_comm_rank(new_comm,irank,err)
CALL mpi_comm_size(new_comm,isize,err)
......@@ -391,20 +425,21 @@ CONTAINS
IF(isize==1) return !No further parallelism
!Split off the pw-part
pw_here=(irank==0)
mt_here=(irank>1)
vac_here=vac_here.AND.(irank>1)
misc_here=misc_here.AND.(irank>1)
mt_here=(irank>0)
vac_here=vac_here.AND.(irank>0)
misc_here=misc_here.AND.(irank>0)
isize=isize-1
irank=irank-1
mt_rank=irank
mt_size=isize
IF(isize==1.OR.irank<0) RETURN !No further parallelism
IF (vac_here.OR.misc_here) THEN !split off-vacuum&misc part
vac_here=vac_here.AND.(irank==0)
misc_here=misc_here.AND.(irank==0)
mt_here=(irank>1)
mt_here=(irank>0)
isize=isize-1
irank=irank-1
ENDIF
IF(isize==1.OR.irank<0) RETURN !No further parallelism