- added a subroutine to write the charge density in XML plane by plane

- minor fixes to PW/pw_restart
- Added the possibility to restart using CP from a PW run (at gamma),
  working but still sperimental.


git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@2131 c92efa57-630b-4861-b058-cf58834340f0
This commit is contained in:
cavazzon 2005-08-26 09:28:33 +00:00
parent 3ecea7dd3c
commit 08a4556d85
10 changed files with 660 additions and 201 deletions

View File

@ -1689,4 +1689,3 @@ SUBROUTINE gmeshinfo( )
RETURN
END SUBROUTINE gmeshinfo

View File

@ -28,8 +28,8 @@ MODULE cp_restart
SUBROUTINE cp_writefile( ndw, scradir, ascii, nfi, simtime, acc, nk, xk, wk, &
ht, htm, htvel, gvel, xnhh0, xnhhm, vnhh, taui, cdmi, stau0, &
svel0, staum, svelm, force, vnhp, xnhp0, xnhpm, nhpcl, occ0, &
occm, lambda0, lambdam, xnhe0, xnhem, vnhe, ekincm, mat_z, et, &
c04, cm4, c02, cm2 )
occm, lambda0, lambdam, xnhe0, xnhem, vnhe, ekincm, et, rho, &
c04, cm4, c02, cm2, mat_z )
!
USE iotk_module
USE kinds, ONLY: dbl
@ -37,7 +37,7 @@ MODULE cp_restart
USE control_flags, ONLY: gamma_only, force_pairing
USE io_files, ONLY: iunpun, xmlpun, psfile, pseudo_dir, prefix
USE printout_base, ONLY: title
USE grid_dimensions, ONLY: nr1, nr2, nr3
USE grid_dimensions, ONLY: nr1, nr2, nr3, nr1x, nr2x, nr3l
USE smooth_grid_dimensions, ONLY: nr1s, nr2s, nr3s
USE smallbox_grid_dimensions, ONLY: nr1b, nr2b, nr3b
USE gvecp, ONLY: ngm, ngmt, ecutp, gcutp
@ -51,6 +51,7 @@ MODULE cp_restart
USE funct, ONLY: dft
USE mp, ONLY: mp_sum, mp_bcast
USE parameters, ONLY: nhclm
USE fft_base, ONLY: dfftp
IMPLICIT NONE
INTEGER, INTENT(IN) :: ndw !
@ -88,14 +89,15 @@ MODULE cp_restart
REAL(dbl), INTENT(IN) :: xnhem !
REAL(dbl), INTENT(IN) :: vnhe !
REAL(dbl), INTENT(IN) :: ekincm !
REAL(dbl), INTENT(IN) :: mat_z(:,:,:) !
REAL(dbl), INTENT(IN) :: et(:,:,:) !
REAL(dbl), INTENT(IN) :: rho(:,:) !
COMPLEX(dbl), OPTIONAL, INTENT(IN) :: c04(:,:,:,:) !
COMPLEX(dbl), OPTIONAL, INTENT(IN) :: cm4(:,:,:,:) !
COMPLEX(dbl), OPTIONAL, INTENT(IN) :: c02(:,:) !
COMPLEX(dbl), OPTIONAL, INTENT(IN) :: cm2(:,:) !
REAL(dbl), OPTIONAL, INTENT(IN) :: mat_z(:,:,:) !
!
CHARACTER(LEN=256) :: dirname, filename
CHARACTER(LEN=256) :: dirname, filename, rho_file
CHARACTER(iotk_attlenx) :: attr
CHARACTER(LEN=4) :: cspin
INTEGER :: kunit, ib
@ -107,9 +109,12 @@ MODULE cp_restart
INTEGER, ALLOCATABLE :: ftmp(:,:)
INTEGER, ALLOCATABLE :: ityp(:)
REAL(dbl), ALLOCATABLE :: tau(:,:)
REAL(dbl), ALLOCATABLE :: rhosum(:)
REAL(dbl) :: omega, htm1(3,3), h(3,3)
REAL(dbl) :: a1(3), a2(3), a3(3)
REAL(dbl) :: b1(3), b2(3), b3(3)
REAL(dbl) :: nelec
REAL(dbl) :: scalef
LOGICAL :: lsda
!
@ -144,6 +149,8 @@ MODULE cp_restart
a3 = ht( 3, : )
CALL recips( a1, a2, a3, b1, b2, b3 )
!
scalef = 1.0d0 / SQRT( omega )
!
! Compute array ityp, and tau
!
ALLOCATE( ityp( nat ) )
@ -211,7 +218,7 @@ MODULE cp_restart
! ... IONS
!
CALL write_ions( nsp, nat, atm, ityp, &
psfile, pseudo_dir, amass, tau, iforce, dirname )
psfile, pseudo_dir, amass, tau, iforce, dirname, "Bohr" )
!
! ... PLANE_WAVES
!
@ -251,12 +258,52 @@ MODULE cp_restart
!
CALL iotk_write_begin( iunpun, "CHARGE-DENSITY" )
!
CALL iotk_link( iunpun, "RHO_FILE", TRIM( prefix ) // ".rho", &
CREATE = .FALSE., BINARY = .TRUE., RAW = .TRUE. )
!
! CALL iotk_write_begin( iunpun, "CHARGE-DENSITY", attr = attr )
! CALL iotk_write_dat( iunpun, "RHO", rho )
! CALL iotk_write_end( iunpun, "CHARGE-DENSITY" )
END IF
!
rho_file = TRIM( prefix ) // ".rho"
!
IF( ionode ) THEN
CALL iotk_link( iunpun, "RHO_FILE", rho_file, CREATE = .FALSE., BINARY = .FALSE. )
END IF
!
rho_file = TRIM( dirname ) // '/' // TRIM( rho_file )
!
IF( nspin == 1 ) THEN
!
CALL write_rho_xml( rho_file, rho(:,1), nr1, nr2, nr3, nr1x, nr2x, dfftp%ipp, dfftp%npp )
!
ELSE IF( nspin == 2 ) THEN
!
ALLOCATE( rhosum( SIZE( rho, 1 ) ) )
rhosum = rho(:,1) + rho(:,2)
!
CALL write_rho_xml( rho_file, rhosum, nr1, nr2, nr3, nr1x, nr2x, dfftp%ipp, dfftp%npp )
!
DEALLOCATE( rhosum )
!
rho_file = TRIM( prefix ) // ".rhoup"
!
IF( ionode ) THEN
CALL iotk_link( iunpun, "RHO_FILE", rho_file, CREATE = .FALSE., BINARY = .FALSE. )
END IF
!
rho_file = TRIM( dirname ) // '/' // TRIM( rho_file )
!
CALL write_rho_xml( rho_file, rho(:,1), nr1, nr2, nr3, nr1x, nr2x, dfftp%ipp, dfftp%npp )
!
rho_file = TRIM( prefix ) // ".rhodw"
!
IF( ionode ) THEN
CALL iotk_link( iunpun, "RHO_FILE", rho_file, CREATE = .FALSE., BINARY = .FALSE. )
END IF
!
rho_file = TRIM( dirname ) // '/' // TRIM( rho_file )
!
CALL write_rho_xml( rho_file, rho(:,2), nr1, nr2, nr3, nr1x, nr2x, dfftp%ipp, dfftp%npp )
!
END IF
!
IF( ionode ) THEN
!
CALL iotk_write_end( iunpun, "CHARGE-DENSITY" )
@ -342,12 +389,14 @@ MODULE cp_restart
!
CALL iotk_write_dat( iunpun, "NUMBER_OF_SPIN_COMPONENTS", nspin )
!
nelec = nelt
!
IF( nspin == 2 ) THEN
call iotk_write_attr (attr,"up",nel(1),first=.true.)
call iotk_write_attr (attr,"dw",nel(2))
CALL iotk_write_dat( iunpun, "NUMBER_OF_ELECTRONS", nelt, ATTR = attr )
CALL iotk_write_dat( iunpun, "NUMBER_OF_ELECTRONS", nelec, ATTR = attr )
ELSE
CALL iotk_write_dat( iunpun, "NUMBER_OF_ELECTRONS", nelt )
CALL iotk_write_dat( iunpun, "NUMBER_OF_ELECTRONS", nelec )
END IF
!
IF( nspin == 2 ) THEN
@ -401,46 +450,62 @@ MODULE cp_restart
!
IF ( ionode ) THEN
!
filename = TRIM( wfc_filename( ".", 'evc', ik, ispin ) )
IF( nspin == 1 ) THEN
filename = TRIM( wfc_filename( ".", 'evc', ik ) )
ELSE
filename = TRIM( wfc_filename( ".", 'evc', ik, ispin ) )
END IF
!
CALL iotk_link( iunpun, "wfc", filename, &
CREATE = .FALSE., BINARY = .TRUE., RAW = .TRUE. )
!
filename = TRIM( wfc_filename( dirname, 'evc', ik, ispin ) )
IF( nspin == 1 ) THEN
filename = TRIM( wfc_filename( dirname, 'evc', ik ) )
ELSE
filename = TRIM( wfc_filename( dirname, 'evc', ik, ispin ) )
END IF
END IF
IF( PRESENT( C04 ) ) THEN
CALL write_wfc( iunout, ik, nk, kunit, ispin, nspin, &
c04(:,:,ik,ispin), ngwt, nbnd, ig_l2g, &
ngw, filename )
ngw, filename, scalef )
ELSE IF( PRESENT( C02 ) ) THEN
ib = iupdwn(ispin)
CALL write_wfc( iunout, ik, nk, kunit, ispin, nspin, &
c02(:,ib:), ngwt, nbnd, ig_l2g, &
ngw, filename )
ngw, filename, scalef )
END IF
IF ( ionode ) THEN
!
filename = TRIM( wfc_filename( ".", 'evcm', ik, ispin ) )
IF( nspin == 1 ) THEN
filename = TRIM( wfc_filename( ".", 'evcm', ik ) )
ELSE
filename = TRIM( wfc_filename( ".", 'evcm', ik, ispin ) )
END IF
!
CALL iotk_link( iunpun, "wfcm", filename, &
CREATE = .FALSE., BINARY = .TRUE., RAW = .TRUE. )
!
filename = TRIM( wfc_filename( dirname, 'evcm', ik, ispin ) )
IF( nspin == 1 ) THEN
filename = TRIM( wfc_filename( dirname, 'evcm', ik ) )
ELSE
filename = TRIM( wfc_filename( dirname, 'evcm', ik, ispin ) )
END IF
END IF
IF( PRESENT( cm4 ) ) THEN
CALL write_wfc( iunout, ik, nk, kunit, ispin, nspin, &
cm4(:,:,ik,ispin), ngwt, nbnd, ig_l2g, &
ngw, filename )
ngw, filename, scalef )
ELSE IF( PRESENT( c02 ) ) THEN
ib = iupdwn(ispin)
CALL write_wfc( iunout, ik, nk, kunit, ispin, nspin, &
cm2(:,ib:), ngwt, nbnd, ig_l2g, &
ngw, filename )
ngw, filename, scalef )
END IF
!
@ -458,7 +523,7 @@ MODULE cp_restart
!
cspin = iotk_index(ispin)
!
IF ( ionode ) THEN
IF ( ionode .AND. PRESENT( mat_z ) ) THEN
!
filename = 'mat_z' // cspin
call iotk_link(iunpun,"mat_z" // TRIM( cspin ), filename,create=.true.,binary=.true.,raw=.true.)
@ -508,7 +573,7 @@ MODULE cp_restart
ht, htm, htvel, gvel, xnhh0, xnhhm, vnhh, taui, cdmi, stau0, &
svel0, staum, svelm, force, vnhp, xnhp0, xnhpm, nhpcl, &
occ0, occm, lambda0, lambdam, b1, b2, b3, xnhe0, xnhem, &
vnhe, ekincm, mat_z, tens, c04, cm4, c02, cm2 )
vnhe, ekincm, c04, cm4, c02, cm2, mat_z )
!
USE iotk_module
USE kinds, ONLY: dbl
@ -524,11 +589,12 @@ MODULE cp_restart
USE gvecs, ONLY: ngs, ngst
USE gvecw, ONLY: ngw, ngwt, ecutw
USE electrons_base, ONLY: nspin, nbnd, nelt, nel, nupdwn, iupdwn
USE cell_base, ONLY: ibrav, alat, celldm, symm_type
USE ions_base, ONLY: nsp, nat, na, atm, zv, pmass
USE cell_base, ONLY: ibrav, alat, celldm, symm_type, s_to_r, r_to_s
USE ions_base, ONLY: nsp, nat, na, atm, zv, pmass, sort_tau, atm, ityp, ions_cofmass
USE reciprocal_vectors, ONLY: ngwt, ngw, ig_l2g, mill_l
USE mp, ONLY: mp_sum, mp_bcast
USE parameters, ONLY: nhclm
USE parameters, ONLY: nhclm, ntypx
USE constants, ONLY: eps8, ANGSTROM_AU
IMPLICIT NONE
INTEGER, INTENT(IN) :: ndr ! I/O unit number
@ -569,41 +635,54 @@ MODULE cp_restart
REAL(dbl), INTENT(INOUT) :: xnhem !
REAL(dbl), INTENT(INOUT) :: vnhe !
REAL(dbl), INTENT(INOUT) :: ekincm !
REAL(dbl), INTENT(INOUT) :: mat_z(:,:,:) !
LOGICAL, INTENT(IN) :: tens
COMPLEX(dbl), OPTIONAL, INTENT(INOUT) :: c04(:,:,:,:) !
COMPLEX(dbl), OPTIONAL, INTENT(INOUT) :: cm4(:,:,:,:) !
COMPLEX(dbl), OPTIONAL, INTENT(INOUT) :: c02(:,:) !
COMPLEX(dbl), OPTIONAL, INTENT(INOUT) :: cm2(:,:) !
REAL(dbl), OPTIONAL, INTENT(INOUT) :: mat_z(:,:,:) !
!
CHARACTER(LEN=256) :: dirname, kdirname, filename
CHARACTER(LEN=5) :: kindex
CHARACTER(LEN=4) :: cspin
INTEGER :: strlen
INTEGER :: strlen
INTEGER :: kunit
INTEGER :: k1, k2, k3
INTEGER :: nk1, nk2, nk3
INTEGER :: i, j, ispin, ig, nspin_wfc, ierr, ik
REAL(dbl) :: omega, htm1( 3, 3 ), hinv( 3, 3 ), scalef
LOGICAL :: found
LOGICAL :: tread_cm
INTEGER, ALLOCATABLE :: mill(:,:)
CHARACTER(iotk_attlenx) :: attr
INTEGER :: kunit
INTEGER :: k1, k2, k3
INTEGER :: nk1, nk2, nk3
INTEGER :: i, j, ispin, ig, nspin_wfc, ierr, ik
INTEGER, ALLOCATABLE :: mill(:,:)
REAL(dbl) :: omega, htm1(3,3)
!
! Variables read for testing pourposes
!
INTEGER :: ibrav_
CHARACTER(LEN=256) :: symm_type_
INTEGER :: nat_ , nsp_, na_
INTEGER :: nk_ , ik_ , nt_
LOGICAL :: gamma_only_ , found
REAL(dbl) :: alat_ , a1_ (3), a2_ (3), a3_ (3)
REAL(dbl) :: pmass_ , zv_
REAL(dbl) :: celldm_ ( 6 )
INTEGER :: ispin_ , nspin_ , ngwt_ , nbnd_ , nelt_
INTEGER :: nhpcl_
INTEGER :: ib
INTEGER :: ibrav_
CHARACTER(LEN=9) :: symm_type_
CHARACTER(LEN=3) :: atm_( ntypx )
INTEGER :: nat_ , nsp_, na_
INTEGER :: nk_ , ik_ , nt_
LOGICAL :: gamma_only_
REAL(dbl) :: alat_ , a1_ (3), a2_ (3), a3_ (3)
REAL(dbl) :: pmass_ , zv_
REAL(dbl) :: celldm_ ( 6 )
INTEGER :: ispin_ , nspin_ , ngwt_ , nbnd_
REAL(dbl) :: nelec_
REAL(dbl) :: scalef_
REAL(dbl) :: wk_
INTEGER :: nhpcl_
INTEGER :: ib
REAL(dbl) :: amass_ ( ntypx )
INTEGER, ALLOCATABLE :: ityp_ ( : )
INTEGER, ALLOCATABLE :: isrt_ ( : )
REAL(dbl), ALLOCATABLE :: tau_ ( :, : )
INTEGER, ALLOCATABLE :: if_pos_ ( :, : )
CHARACTER(LEN=256) :: psfile_ ( ntypx )
CHARACTER(LEN=80) :: pos_unit
kunit = 1
found = .FALSE.
dirname = restart_dir( scradir, ndr )
filename = TRIM( dirname ) // '/' // 'restart.xml'
@ -617,37 +696,63 @@ MODULE cp_restart
call errore(" cp_readfile ", " cannot open restart file for reading ", ierr )
ierr = 0
IF( ionode ) THEN
!
call iotk_scan_begin(iunpun,"STATUS" )
call iotk_scan_empty(iunpun,"STEP",attr)
call iotk_scan_attr (attr,"nfi",nfi)
call iotk_scan_dat (iunpun, "TIME", simtime )
call iotk_scan_dat (iunpun, "TITLE", title )
call iotk_scan_end(iunpun,"STATUS")
call iotk_scan_begin(iunpun,"STATUS", found = found )
!
IF( found ) THEN
!
call iotk_scan_empty(iunpun,"STEP",attr)
call iotk_scan_attr (attr,"nfi",nfi)
call iotk_scan_dat (iunpun, "TIME", simtime )
call iotk_scan_dat (iunpun, "TITLE", title )
call iotk_scan_end(iunpun,"STATUS")
!
END IF
!
END IF
!
! Read cell and positions
!
ALLOCATE( tau_ ( 3, nat ) )
ALLOCATE( if_pos_ ( 3, nat ) )
ALLOCATE( ityp_ ( nat ) )
!
IF( ionode ) THEN
!
call read_cell( ibrav_ , symm_type_ , celldm_ , alat_ , &
a1_ , a2_ , a3_ , b1, b2, b3 )
END IF
!
IF( ionode ) THEN
!
call iotk_scan_begin(iunpun,"IONS")
call iotk_scan_dat ( iunpun, "NUMBER_OF_ATOMS", nat_ )
call iotk_scan_dat ( iunpun, "NUMBER_OF_SPECIES", nsp_ )
IF( nsp_ /= nsp .OR. nat_ /= nat ) THEN
ierr = 10
GOTO 100
END IF
call iotk_scan_end(iunpun,"IONS")
call read_ions( nsp_ , nat_ , atm_ , ityp_ , psfile_ , amass_ , tau_ , if_pos_ , pos_unit, ierr )
!
IF( ierr == 0 ) THEN
IF( nsp_ /= nsp .OR. nat_ /= nat ) ierr = 2
DO i = 1, nat
IF( ityp_( i ) /= ityp( i ) ) ierr = 3
END DO
END IF
!
END IF
!
CALL mp_bcast( ierr, ionode_id )
IF( ierr /= 0 ) &
call errore(" cp_readfile ", " cannot read positions from restart file ", ierr )
!
! Read MD timesteps variables
!
IF( ionode ) THEN
!
call iotk_scan_begin( iunpun, "TIMESTEPS", attr, found = found )
!
END IF
!
IF( ionode ) THEN
IF( ionode .AND. found ) THEN
!
call iotk_scan_begin(iunpun,"TIMESTEPS", attr)
call iotk_scan_attr (attr, "nt", nt_ )
!
END IF
!
IF( ionode ) THEN
call iotk_scan_attr ( attr, "nt", nt_ )
!
IF( nt_ > 0 ) THEN
!
@ -688,17 +793,12 @@ MODULE cp_restart
call iotk_scan_end(iunpun,"CELL_NOSE")
!
call iotk_scan_end(iunpun,"STEP0")
!
ELSE
ierr = 40
GOTO 100
END IF
!
END IF
!
IF( ionode ) THEN
!
IF( nt_ > 1 ) THEN
!
call iotk_scan_begin(iunpun,"STEPM")
@ -732,14 +832,76 @@ MODULE cp_restart
!
END IF
!
END IF ! ionode
!
IF( ionode ) THEN
call iotk_scan_end(iunpun,"TIMESTEPS")
!
ELSE IF( ionode ) THEN
!
! MD time steps not found, try to recover from CELL and POSITIONS
!
acc = 0.0d0
!
ALLOCATE( isrt_ ( nat ) )
!
SELECT CASE ( TRIM( pos_unit ) )
CASE ( "alat" )
tau_ = tau_ * alat_
CASE ( "Angstrom" )
tau_ = tau_ * ANGSTROM_AU
CASE DEFAULT
END SELECT
!
CALL sort_tau( taui, isrt_ , tau_ , ityp_ , nat_ , nsp_ )
!
ht( 1, : ) = a1_
ht( 2, : ) = a2_
ht( 3, : ) = a3_
!
CALL invmat( 3, ht, htm1, omega )
!
hinv = TRANSPOSE( htm1 )
!
CALL r_to_s( taui, stau0, na, nsp, hinv )
!
CALL ions_cofmass( taui, amass_ , na, nsp, cdmi )
!
staum = stau0
svel0 = 0.0d0
svelm = 0.0d0
force = 0.0d0
!
htm = ht
htvel = 0.0d0
gvel = 0.0d0
xnhh0 = 0.0d0
vnhh = 0.0d0
xnhhm = 0.0d0
!
xnhe0 = 0.0d0
xnhem = 0.0d0
vnhe = 0.0d0
!
ekincm = 0.0d0
!
xnhp0 = 0.0d0
xnhpm = 0.0d0
vnhp = 0.0d0
!
DEALLOCATE( isrt_ )
END IF ! ionode
DEALLOCATE( tau_ )
DEALLOCATE( if_pos_ )
DEALLOCATE( ityp_ )
!
! Compute the scale factor
!
IF( ionode ) CALL invmat( 3, ht, htm1, omega )
CALL mp_bcast( omega, ionode_id )
scalef = 1.0d0 / SQRT( ABS( omega ) )
!
! Band Structure
!
@ -751,14 +913,14 @@ MODULE cp_restart
call iotk_scan_dat( iunpun, "NUMBER_OF_SPIN_COMPONENTS", nspin_ )
!
IF( nspin_ == 2 ) THEN
call iotk_scan_dat( iunpun, "NUMBER_OF_ELECTRONS", nelt_ , ATTR = attr)
call iotk_scan_dat( iunpun, "NUMBER_OF_ELECTRONS", nelec_ , ATTR = attr)
call iotk_scan_dat( iunpun, "NUMBER_OF_BANDS", nbnd_ , ATTR = attr)
ELSE
call iotk_scan_dat( iunpun, "NUMBER_OF_ELECTRONS", nelt_ )
call iotk_scan_dat( iunpun, "NUMBER_OF_ELECTRONS", nelec_ )
call iotk_scan_dat( iunpun, "NUMBER_OF_BANDS", nbnd_ )
END IF
IF( ( nspin_ /= nspin ) .OR. ( nbnd_ /= nbnd ) .OR. ( nelt_ /= nelt ) ) THEN
IF( ( nspin_ /= nspin ) .OR. ( nbnd_ /= nbnd ) .OR. ( NINT( nelec_ ) /= nelt ) ) THEN
ierr = 30
GOTO 100
END IF
@ -772,53 +934,85 @@ MODULE cp_restart
k_points_loop: DO ik = 1, nk
!
IF ( ionode ) THEN
!
CALL iotk_scan_begin( iunpun, "K-POINT" // TRIM( iotk_index(ik) ) )
!
CALL iotk_scan_dat( iunpun, "WEIGHT", wk_ )
!
END IF
!
DO ispin = 1, nspin
!
cspin = iotk_index( ispin )
!
tread_cm = .TRUE.
!
IF( ionode ) THEN
CALL iotk_scan_dat( iunpun, "OCC" // TRIM( cspin ), occ0(:, ik, ispin ) )
CALL iotk_scan_dat( iunpun, "OCCM" // TRIM( cspin ), occm(:, ik, ispin ), FOUND=found, IERR=ierr )
IF ( .NOT. found ) occm(:, ik, ispin ) = occ0(:, ik, ispin )
occ0(:, ik, ispin ) = occ0(:, ik, ispin ) * wk_
IF ( .NOT. found ) THEN
occm(:, ik, ispin ) = occ0(:, ik, ispin )
tread_cm = .FALSE.
END IF
END IF
!
CALL mp_bcast( tread_cm, ionode_id )
IF ( ionode ) THEN
!
filename = TRIM( wfc_filename( dirname, 'evc', ik, ispin ) )
IF( nspin == 1 ) THEN
filename = TRIM( wfc_filename( dirname, 'evc', ik ) )
ELSE
filename = TRIM( wfc_filename( dirname, 'evc', ik, ispin ) )
END IF
!
END IF
!
IF( PRESENT( c04 ) ) THEN
CALL read_wfc( iunout, ik, nk, kunit, ispin_ , nspin_ , &
c04(:,:,ik,ispin), ngwt_ , nbnd_ , ig_l2g, &
ngw, filename )
ngw, filename, scalef_ )
ELSE IF( PRESENT( c02 ) ) THEN
ib = iupdwn(ispin)
CALL read_wfc( iunout, ik, nk, kunit, ispin_ , nspin_ , &
c02( :, ib: ), ngwt_ , nbnd_ , ig_l2g, &
ngw, filename )
END IF
IF ( ionode ) THEN
!
filename = TRIM( wfc_filename( dirname, 'evcm', ik, ispin ) )
!
ngw, filename, scalef_ )
END IF
!
IF( PRESENT( cm4 ) ) THEN
CALL read_wfc( iunout, ik, nk, kunit, ispin_ , nspin_ , &
cm4(:,:,ik,ispin), ngwt_ , nbnd_ , ig_l2g, &
ngw, filename )
ELSE IF( PRESENT( cm2 ) ) THEN
ib = iupdwn(ispin)
CALL read_wfc( iunout, ik, nk, kunit, ispin_ , nspin_ , &
cm2( :, ib: ), ngwt_ , nbnd_ , ig_l2g, &
ngw, filename )
END IF
IF( tread_cm ) THEN
!
IF ( ionode ) THEN
!
IF( nspin == 1 ) THEN
filename = TRIM( wfc_filename( dirname, 'evcm', ik ) )
ELSE
filename = TRIM( wfc_filename( dirname, 'evcm', ik, ispin ) )
END IF
!
END IF
!
IF( PRESENT( cm4 ) ) THEN
CALL read_wfc( iunout, ik, nk, kunit, ispin_ , nspin_ , &
cm4(:,:,ik,ispin), ngwt_ , nbnd_ , ig_l2g, &
ngw, filename, scalef_ )
ELSE IF( PRESENT( cm2 ) ) THEN
ib = iupdwn(ispin)
CALL read_wfc( iunout, ik, nk, kunit, ispin_ , nspin_ , &
cm2( :, ib: ), ngwt_ , nbnd_ , ig_l2g, &
ngw, filename, scalef_ )
END IF
!
ELSE
!
IF( PRESENT( cm4 ) ) THEN
cm4 = c04
ELSE IF( PRESENT( cm2 ) ) THEN
cm2 = c02
END IF
!
END IF
!
END DO
!
@ -828,9 +1022,8 @@ MODULE cp_restart
!
END DO k_points_loop
DO ispin = 1, nspin
IF( ionode .AND. tens ) THEN
IF( ionode .AND. PRESENT( mat_z ) ) THEN
call iotk_scan_dat( iunpun, "mat_z" // TRIM( iotk_index(ispin) ), mat_z(:,:,ispin) )
END IF
END DO
@ -889,7 +1082,7 @@ MODULE cp_restart
CALL mp_bcast(occ0( :, :, :), ionode_id)
CALL mp_bcast(occm( :, :, :), ionode_id)
!
IF( tens ) THEN
IF( PRESENT( mat_z ) ) THEN
CALL mp_bcast(mat_z( :, :, :), ionode_id)
END IF
!
@ -903,10 +1096,14 @@ MODULE cp_restart
!
filename = TRIM( dirname ) // '/lambda.dat'
IF( ionode ) THEN
OPEN( unit = 10, file = TRIM(filename), status = 'OLD', form = 'UNFORMATTED' )
READ( 10 ) lambda0
READ( 10 ) lambdam
CLOSE( unit = 10 )
INQUIRE( file = TRIM(filename), EXIST = found )
IF( found ) THEN
OPEN( unit = 10, file = TRIM(filename), status = 'OLD', form = 'UNFORMATTED' )
READ( 10 ) lambda0
READ( 10 ) lambdam
CLOSE( unit = 10 )
ELSE
END IF
END IF
CALL mp_bcast(lambda0, ionode_id)
CALL mp_bcast(lambdam, ionode_id)
@ -934,25 +1131,34 @@ MODULE cp_restart
!
CHARACTER(LEN=256) :: dirname, filename
INTEGER :: ib, kunit , ispin_ , nspin_ , ngwt_ , nbnd_
REAL(dbl) :: scalef
!
kunit = 1
!
dirname = restart_dir( scradir, ndr )
IF( tag /= 'm' ) THEN
filename = TRIM( wfc_filename( dirname, 'evc', ik, ispin ) )
IF( nspin == 1 ) THEN
filename = TRIM( wfc_filename( dirname, 'evc', ik ) )
ELSE
filename = TRIM( wfc_filename( dirname, 'evc', ik, ispin ) )
END IF
ELSE
filename = TRIM( wfc_filename( dirname, 'evcm', ik, ispin ) )
IF( nspin == 1 ) THEN
filename = TRIM( wfc_filename( dirname, 'evcm', ik ) )
ELSE
filename = TRIM( wfc_filename( dirname, 'evcm', ik, ispin ) )
END IF
END IF
!
IF( PRESENT( c4 ) ) THEN
CALL read_wfc( iunout, ik, nk, kunit, ispin_ , nspin_ , &
c4(:,:,ik,ispin), ngwt_ , nbnd_ , ig_l2g, &
ngw, filename )
ngw, filename, scalef )
ELSE IF( PRESENT( c2 ) ) THEN
ib = iupdwn(ispin)
CALL read_wfc( iunout, ik, nk, kunit, ispin_ , nspin_ , &
c2( :, ib: ), ngwt_ , nbnd_ , ig_l2g, &
ngw, filename )
ngw, filename, scalef )
END IF
!
RETURN
@ -990,9 +1196,12 @@ MODULE cp_restart
!
! Variables read for testing pourposes
!
INTEGER :: ibrav_
REAL(dbl) :: alat_
REAL(dbl) :: celldm_ ( 6 )
INTEGER :: ibrav_
REAL(dbl) :: alat_
REAL(dbl) :: celldm_ ( 6 )
REAL(dbl) :: a1_ ( 3 ), a2_ ( 3 ), a3_ ( 3 )
REAL(dbl) :: b1_ ( 3 ), b2_ ( 3 ), b3_ ( 3 )
CHARACTER(LEN=9) :: symm_type_
dirname = restart_dir( scradir, ndr )
filename = TRIM( dirname ) // '/' // 'restart.xml'
@ -1004,52 +1213,76 @@ MODULE cp_restart
CALL mp_bcast( ierr, ionode_id )
IF ( ierr /= 0 ) &
call errore(" cp_read_cell ", " cannot open restart file for reading ", ierr )
!
ierr = 0
!
IF( ionode ) THEN
!
call iotk_scan_begin(iunpun,"TIMESTEPS", attr)
call iotk_scan_attr (attr, "nt", nt_ )
call iotk_scan_begin(iunpun,"TIMESTEPS", attr, found = found )
IF( found ) THEN
call iotk_scan_attr (attr, "nt", nt_ )
!
IF( nt_ > 0 ) THEN
!
call iotk_scan_begin(iunpun,"STEP0")
!
call iotk_scan_begin(iunpun,"CELL_PARAMETERS")
call iotk_scan_dat (iunpun, "ht", ht)
call iotk_scan_dat (iunpun, "htvel", htvel)
call iotk_scan_dat (iunpun, "gvel", gvel, FOUND=found, IERR=ierr )
if ( .NOT. found ) gvel = 0.0d0
call iotk_scan_end(iunpun,"CELL_PARAMETERS")
!
call iotk_scan_begin(iunpun,"CELL_NOSE")
call iotk_scan_dat (iunpun, "xnhh", xnhh0)
call iotk_scan_dat (iunpun, "vnhh", vnhh)
call iotk_scan_end(iunpun,"CELL_NOSE")
!
call iotk_scan_end(iunpun,"STEP0")
!
ELSE
ierr = 40
GOTO 100
END IF
IF( nt_ > 1 ) THEN
!
call iotk_scan_begin(iunpun,"STEPM")
!
call iotk_scan_begin(iunpun,"CELL_PARAMETERS")
call iotk_scan_dat (iunpun, "ht", htm)
call iotk_scan_end(iunpun,"CELL_PARAMETERS")
!
call iotk_scan_begin(iunpun,"CELL_NOSE")
call iotk_scan_dat (iunpun, "xnhh", xnhhm)
call iotk_scan_end(iunpun,"CELL_NOSE")
!
call iotk_scan_end(iunpun,"STEPM")
!
END IF
!
call iotk_scan_end(iunpun,"TIMESTEPS")
ELSE
!
IF( nt_ > 0 ) THEN
!
call iotk_scan_begin(iunpun,"STEP0")
!
call iotk_scan_begin(iunpun,"CELL_PARAMETERS")
call iotk_scan_dat (iunpun, "ht", ht)
call iotk_scan_dat (iunpun, "htvel", htvel)
call iotk_scan_dat (iunpun, "gvel", gvel, FOUND=found, IERR=ierr )
if ( .NOT. found ) gvel = 0.0d0
call iotk_scan_end(iunpun,"CELL_PARAMETERS")
!
call iotk_scan_begin(iunpun,"CELL_NOSE")
call iotk_scan_dat (iunpun, "xnhh", xnhh0)
call iotk_scan_dat (iunpun, "vnhh", vnhh)
call iotk_scan_end(iunpun,"CELL_NOSE")
!
call iotk_scan_end(iunpun,"STEP0")
!
ELSE
ierr = 40
GOTO 100
END IF
IF( nt_ > 1 ) THEN
!
call iotk_scan_begin(iunpun,"STEPM")
!
call iotk_scan_begin(iunpun,"CELL_PARAMETERS")
call iotk_scan_dat (iunpun, "ht", htm)
call iotk_scan_end(iunpun,"CELL_PARAMETERS")
!
call iotk_scan_begin(iunpun,"CELL_NOSE")
call iotk_scan_dat (iunpun, "xnhh", xnhhm)
call iotk_scan_end(iunpun,"CELL_NOSE")
!
call iotk_scan_end(iunpun,"STEPM")
!
END IF
! MD steps have not been found, try to restart from cell data
!
call iotk_scan_end(iunpun,"TIMESTEPS")
CALL read_cell( ibrav_ , symm_type_ , celldm_ , alat_ , a1_ , a2_ , a3_ , b1_ , b2_ , b3_ )
!
ht( 1, : ) = a1_
ht( 2, : ) = a2_
ht( 3, : ) = a3_
!
htm = ht
htvel = 0.0d0
gvel = 0.0d0
xnhh0 = 0.0d0
vnhh = 0.0d0
xnhhm = 0.0d0
!
END IF
!
END IF
!
@ -1058,7 +1291,7 @@ MODULE cp_restart
CALL mp_bcast( ierr, ionode_id )
CALL mp_bcast( attr, ionode_id )
IF( ierr /= 0 ) THEN
CALL errore( " cp_readfile ", attr, ierr )
CALL errore( " cp_read_cell ", attr, ierr )
END IF
!
CALL mp_bcast( ht, ionode_id )
@ -1084,7 +1317,8 @@ MODULE cp_restart
! .. This subroutine write wavefunctions to the disk
!
SUBROUTINE write_wfc(iuni, ik, nk, kunit, ispin, nspin, wf, ngw, nbnd, igl, ngwl, filename )
SUBROUTINE write_wfc(iuni, ik, nk, kunit, ispin, nspin, wf, ngw, nbnd, igl, &
ngwl, filename, scalef )
!
USE kinds, ONLY: dbl
USE mp_wave
@ -1104,6 +1338,7 @@ MODULE cp_restart
INTEGER, INTENT(IN) :: ngwl
INTEGER, INTENT(IN) :: igl(:)
CHARACTER(LEN=256), INTENT(IN) :: filename
REAL(dbl), INTENT(IN) :: scalef
INTEGER :: i, j, ierr
INTEGER :: nkl, nkr, nkbl, iks, ike, nkt, ikt, igwx
@ -1209,6 +1444,7 @@ MODULE cp_restart
CALL iotk_write_attr( attr, "ispin", ispin )
CALL iotk_write_attr( attr, "nspin", nspin )
CALL iotk_write_attr( attr, "igwx", igwx )
CALL iotk_write_attr( attr, "scale_factor", scalef )
!
CALL iotk_write_empty( iuni, "INFO", attr )
!
@ -1252,7 +1488,7 @@ MODULE cp_restart
!=----------------------------------------------------------------------------=!
SUBROUTINE read_wfc( iuni, ik, nk, kunit, ispin, nspin, wf, ngw, nbnd, igl, &
ngwl, filename )
ngwl, filename, scalef )
!
USE kinds, ONLY: dbl
USE mp_wave
@ -1272,6 +1508,7 @@ MODULE cp_restart
INTEGER, INTENT(IN) :: ngwl
INTEGER, INTENT(IN) :: igl(:)
CHARACTER(LEN=256), INTENT(IN) :: filename
REAL(dbl), INTENT(OUT) :: scalef
INTEGER :: ierr_iotk
CHARACTER(LEN=iotk_attlenx) :: attr
@ -1381,6 +1618,7 @@ MODULE cp_restart
CALL iotk_scan_attr( attr, "ispin", ispin )
CALL iotk_scan_attr( attr, "nspin", nspin )
CALL iotk_scan_attr( attr, "igwx", igwx_ )
CALL iotk_scan_attr( attr, "scale_factor", scalef )
!
END IF
@ -1392,6 +1630,7 @@ MODULE cp_restart
CALL mp_bcast( ispin, ionode_id )
CALL mp_bcast( nspin, ionode_id )
CALL mp_bcast( igwx_ , ionode_id )
CALL mp_bcast( scalef , ionode_id )
ALLOCATE( wtmp( MAX(igwx_, igwx) ) )
@ -1514,6 +1753,67 @@ MODULE cp_restart
RETURN
END SUBROUTINE
!------------------------------------------------------------------------------!
SUBROUTINE read_ions( nsp, nat, atm, ityp, psfile, amass, tau, if_pos, pos_unit, ierr )
!
USE iotk_module
USE kinds, ONLY: dbl
USE io_files, ONLY: iunpun
!
INTEGER, INTENT(OUT) :: nsp, nat
CHARACTER(LEN=3), INTENT(OUT) :: atm( : )
INTEGER, INTENT(OUT) :: ityp( : )
CHARACTER(LEN=256), INTENT(OUT) :: psfile( : )
REAL(dbl), INTENT(OUT) :: amass( : )
REAL(dbl), INTENT(OUT) :: tau( :, : )
INTEGER, INTENT(OUT) :: if_pos( :, : )
INTEGER, INTENT(OUT) :: ierr
CHARACTER(LEN=*), INTENT(OUT) :: pos_unit
!
CHARACTER(iotk_attlenx) :: attr
LOGICAL :: found
INTEGER :: i
CHARACTER(LEN=3) :: lab
!
ierr = 0
!
CALL iotk_scan_begin( iunpun, "IONS", found = found )
IF( .NOT. found ) THEN
ierr = 1
GOTO 110
END IF
!
CALL iotk_scan_dat ( iunpun, "NUMBER_OF_ATOMS", nat )
CALL iotk_scan_dat ( iunpun, "NUMBER_OF_SPECIES", nsp )
IF( nsp > SIZE( atm ) .OR. nat > SIZE( ityp ) ) THEN
ierr = 10
GOTO 100
END IF
!
DO i = 1, nsp
CALL iotk_scan_dat ( iunpun, "ATOM_TYPE", atm( i ) )
CALL iotk_scan_dat( iunpun, TRIM( atm( i ) )//"_MASS", amass( i ), ATTR = attr )
END DO
!
CALL iotk_scan_empty( iunpun, "UNITS_FOR_ATOMIC_POSITIONS", attr )
CALL iotk_scan_attr( attr, "UNIT", pos_unit )
!
DO i = 1, nat
CALL iotk_scan_empty( iunpun, "ATOM" // TRIM( iotk_index(i) ), attr )
CALL iotk_scan_attr ( attr, "SPECIES", lab )
CALL iotk_scan_attr ( attr, "INDEX", ityp( i ) )
CALL iotk_scan_attr ( attr, "tau", tau( :, i ) )
CALL iotk_scan_attr ( attr, "if_pos", if_pos( :, i ) )
END DO
100 call iotk_scan_end( iunpun, "IONS" )
!
110 RETURN
!
END SUBROUTINE read_ions
!------------------------------------------------------------------------------!
LOGICAL FUNCTION check_restartfile( scradir, ndr )
@ -1541,7 +1841,6 @@ MODULE cp_restart
RETURN
END FUNCTION check_restartfile
!------------------------------------------------------------------------------!
END MODULE cp_restart
!------------------------------------------------------------------------------!

View File

@ -760,7 +760,7 @@ SUBROUTINE cprmain( tau, fion_out, etot_out )
vels, velsm, acc, lambda, lambdam, xnhe0, xnhem, &
vnhe, xnhp0, xnhpm, vnhp, nhpcl, ekincm, xnhh0, &
xnhhm, vnhh, velh, ecutp, ecutw, delt, pmass, ibrav,&
celldm, fion, tps, z0, f )
celldm, fion, tps, z0, f, rhor )
!
ELSE
!
@ -768,7 +768,7 @@ SUBROUTINE cprmain( tau, fion_out, etot_out )
tausm, vels, velsm, acc, lambda, lambdam, xnhe0, &
xnhem, vnhe, xnhp0, xnhpm, vnhp, nhpcl, ekincm, &
xnhh0, xnhhm, vnhh, velh, ecutp, ecutw, delt, pmass,&
ibrav, celldm, fion, tps, z0, f )
ibrav, celldm, fion, tps, z0, f, rhor )
!
END IF
!
@ -838,7 +838,7 @@ SUBROUTINE cprmain( tau, fion_out, etot_out )
velsm, acc, lambda, lambdam, xnhe0, xnhem, &
vnhe, xnhp0, xnhpm, vnhp, nhpcl, ekincm, &
xnhh0, xnhhm, vnhh, velh, ecutp, ecutw, &
delt, celldm, fion, tps, z0, f )
delt, celldm, fion, tps, z0, f, rhor )
!
IF ( ( nfi >= nomore ) .OR. tstop ) EXIT main_loop
!
@ -919,7 +919,7 @@ SUBROUTINE cprmain( tau, fion_out, etot_out )
velsm, acc, lambda, lambdam, xnhe0, xnhem, vnhe, xnhp0, &
xnhpm, vnhp, nhpcl, ekincm, xnhh0, xnhhm, vnhh, velh, &
ecutp, ecutw, delt, pmass, ibrav, celldm, fion, tps, &
z0, f )
z0, f, rhor )
!
ELSE
!
@ -927,7 +927,7 @@ SUBROUTINE cprmain( tau, fion_out, etot_out )
vels, velsm, acc, lambda, lambdam, xnhe0, xnhem, vnhe, &
xnhp0, xnhpm, vnhp, nhpcl, ekincm, xnhh0, xnhhm, vnhh, &
velh, ecutp, ecutw, delt, pmass, ibrav, celldm, fion, tps,&
z0, f )
z0, f, rhor )
!
END IF
!

View File

@ -543,9 +543,12 @@ CONTAINS
DO is = 1, nspnl
DO l = 1, nbeta( is )
DO l = 1, nbetax
CALL nullify_spline( wnl_sp( l, is ) )
CALL nullify_spline( wnla_sp( l, is ) )
END DO
DO l = 1, nbeta( is )
CALL allocate_spline( wnl_sp(l,is), pstab_size, xgmin, xgmax )
CALL allocate_spline( wnla_sp(l,is), pstab_size, xgmin, xgmax )
END DO

View File

@ -92,7 +92,7 @@
& ( ndw,h,hold,nfi,c0,cm,taus,tausm,vels,velsm,acc, &
& lambda,lambdam,xnhe0,xnhem,vnhe,xnhp0,xnhpm,vnhp,nhpcl,ekincm, &
& xnhh0,xnhhm,vnhh,velh,ecut,ecutw,delt,pmass,ibrav,celldm, &
& fion, tps, mat_z, occ_f )
& fion, tps, mat_z, occ_f, rho )
!-----------------------------------------------------------------------
!
! read from file and distribute data calculated in preceding iterations
@ -103,6 +103,7 @@
USE electrons_base, ONLY: nspin, nbnd, nbsp, iupdwn, nupdwn
USE electrons_module, ONLY: ei
USE io_files, ONLY: scradir
USE ensemble_dft, ONLY: tens
!
implicit none
integer, INTENT(IN) :: ndw, nfi
@ -119,8 +120,10 @@
real(kind=8), INTENT(in) :: pmass(:)
real(kind=8), INTENT(in) :: celldm(:)
real(kind=8), INTENT(in) :: tps
real(kind=8), INTENT(in) :: rho(:,:)
integer, INTENT(in) :: ibrav
real(kind=8), INTENT(in) :: mat_z(:,:,:), occ_f(:)
real(kind=8), INTENT(in) :: occ_f(:)
real(kind=8), INTENT(in) :: mat_z(:,:,:)
real(kind=8) :: ht(3,3), htm(3,3), htvel(3,3), gvel(3,3)
integer :: nk = 1, ispin, i, ib
@ -159,11 +162,19 @@
end do
end do
CALL cp_writefile( ndw, scradir, .TRUE., nfi, tps, acc, nk, xk, wk, &
ht, htm, htvel, gvel, xnhh0, xnhhm, vnhh, taui_ , cdmi_ , taus, &
vels, tausm, velsm, fion, vnhp, xnhp0, xnhpm, nhpcl, occ_ , &
occ_ , lambda, lambdam, xnhe0, xnhem, vnhe, ekincm, mat_z, ei, &
c02 = c0, cm2 = cm )
IF( tens ) THEN
CALL cp_writefile( ndw, scradir, .TRUE., nfi, tps, acc, nk, xk, wk, &
ht, htm, htvel, gvel, xnhh0, xnhhm, vnhh, taui_ , cdmi_ , taus, &
vels, tausm, velsm, fion, vnhp, xnhp0, xnhpm, nhpcl, occ_ , &
occ_ , lambda, lambdam, xnhe0, xnhem, vnhe, ekincm, ei, &
rho, c02 = c0, cm2 = cm, mat_z = mat_z )
ELSE
CALL cp_writefile( ndw, scradir, .TRUE., nfi, tps, acc, nk, xk, wk, &
ht, htm, htvel, gvel, xnhh0, xnhhm, vnhh, taui_ , cdmi_ , taus, &
vels, tausm, velsm, fion, vnhp, xnhp0, xnhpm, nhpcl, occ_ , &
occ_ , lambda, lambdam, xnhe0, xnhem, vnhe, ekincm, ei, &
rho, c02 = c0, cm2 = cm )
END IF
DEALLOCATE( taui_ )
DEALLOCATE( occ_ )
@ -237,11 +248,19 @@
ALLOCATE( taui_ ( 3, SIZE( taus, 2 ) ) )
ALLOCATE( occ_ ( nbnd, 1, nspin ) )
CALL cp_readfile( ndr, scradir, .TRUE., nfi, tps, acc, nk, xk, wk, &
ht, htm, htvel, gvel, xnhh0, xnhhm, vnhh, taui_ , cdmi_ , taus, &
vels, tausm, velsm, fion, vnhp, xnhp0, xnhpm, nhpcl , occ_ , &
occ_ , lambda, lambdam, b1, b2, b3, &
xnhe0, xnhem, vnhe, ekincm, mat_z, tens, c02 = c0, cm2 = cm )
IF( tens ) THEN
CALL cp_readfile( ndr, scradir, .TRUE., nfi, tps, acc, nk, xk, wk, &
ht, htm, htvel, gvel, xnhh0, xnhhm, vnhh, taui_ , cdmi_ , taus, &
vels, tausm, velsm, fion, vnhp, xnhp0, xnhpm, nhpcl , occ_ , &
occ_ , lambda, lambdam, b1, b2, b3, &
xnhe0, xnhem, vnhe, ekincm, c02 = c0, cm2 = cm, mat_z = mat_z )
ELSE
CALL cp_readfile( ndr, scradir, .TRUE., nfi, tps, acc, nk, xk, wk, &
ht, htm, htvel, gvel, xnhh0, xnhhm, vnhh, taui_ , cdmi_ , taus, &
vels, tausm, velsm, fion, vnhp, xnhp0, xnhpm, nhpcl , occ_ , &
occ_ , lambda, lambdam, b1, b2, b3, &
xnhe0, xnhem, vnhe, ekincm, c02 = c0, cm2 = cm )
END IF
! AutoPilot (Dynamic Rules) Implementation
@ -299,6 +318,7 @@
USE cp_restart, ONLY: cp_writefile
USE electrons_module, ONLY: ei
USE io_files, ONLY: scradir
USE grid_dimensions, ONLY: nr1, nr2, nr3, nr1x, nr2x, nr3x
IMPLICIT NONE
@ -317,9 +337,10 @@
REAL(dbl), INTENT(IN) :: trutime
REAL(dbl), ALLOCATABLE :: lambda(:,:)
REAL(dbl), ALLOCATABLE :: rhow(:,:)
REAL(dbl) S0, S1
REAL(dbl) :: ekincm
REAL(dbl) :: mat_z(1,1,nspin)
INTEGER :: i, j, k, iss, ir
s0 = cclock()
@ -330,16 +351,36 @@
! properties on the writefile subroutine
ALLOCATE( lambda(nbsp,nbsp) )
ALLOCATE( rhow( nr1x * nr2x * SIZE( rho, 3 ), nspin ) )
lambda = 0.0d0
ekincm = 0.0d0
mat_z = 0.0d0
!
!
IF( SIZE( rho, 1 ) /= nr1x .OR. SIZE( rho, 2 ) /= nr2x ) THEN
WRITE( stdout, * ) nr1x, nr2x
WRITE( stdout, * ) SIZE( rho, 1 ), SIZE( rho, 2 )
CALL errore( ' writefile_fpmd ', ' rho dimensions do not correspond ', 1 )
END IF
!
DO iss = 1, nspin
ir = 0
DO k = 1, SIZE( rho, 3 )
DO j = 1, nr2x
DO i = 1, nr1x
ir = ir + 1
rhow( ir, iss ) = rho( i, j, k, iss )
END DO
END DO
END DO
END DO
CALL cp_writefile( ndw, scradir, .TRUE., nfi, trutime, acc, kp%nkpt, kp%xk, kp%weight, &
ht_0%a, ht_m%a, ht_0%hvel, ht_0%gvel, xnhh0, xnhhm, vnhh, taui, cdmi, &
atoms_0%taus, atoms_0%vels, atoms_m%taus, atoms_m%vels, atoms_0%for, vnhp, &
xnhp0, xnhpm, nhpcl, occ, occ, lambda, lambda, &
xnhe0, xnhem, vnhe, ekincm, mat_z, ei, c04 = c0, cm4 = cm )
xnhe0, xnhem, vnhe, ekincm, ei, rhow, c04 = c0, cm4 = cm )
DEALLOCATE( rhow )
DEALLOCATE( lambda )
s1 = cclock()
@ -413,7 +454,6 @@
REAL(dbl) :: hm1_ (3,3)
REAL(dbl) :: gvel_ (3,3)
REAL(dbl) :: hvel_ (3,3)
REAL(dbl) :: mat_z_(1,1,nspin)
REAL(dbl) :: b1(3), b2(3), b3(3)
LOGICAL :: tens = .FALSE.
@ -427,7 +467,7 @@
hp0_ , hm1_ , hvel_ , gvel_ , xnhh0, xnhhm, vnhh, taui, cdmi, &
atoms_0%taus, atoms_0%vels, atoms_m%taus, atoms_m%vels, atoms_0%for, vnhp, &
xnhp0, xnhpm, nhpcl, occ, occ, lambda_ , lambda_ , b1, b2, &
b3, xnhe0, xnhem, vnhe, ekincm, mat_z_ , tens, c04 = c0, cm4 = cm )
b3, xnhe0, xnhem, vnhe, ekincm, c04 = c0, cm4 = cm )
DEALLOCATE( lambda_ )

View File

@ -1379,7 +1379,7 @@ SUBROUTINE smdmain( tau, fion_out, etot_out, nat_out )
& rep_el(sm_k)%lambda,rep_el(sm_k)%lambdam,xnhe0(sm_k),xnhem(sm_k), &
& vnhe(sm_k),xnhp0(:,sm_k),xnhpm(:,sm_k),vnhp(:,sm_k),nhpcl,ekincm(sm_k), &
& xnhh0,xnhhm,vnhh,velh,ecutp,ecutw,delt,pmass,ibrav,celldm, &
& rep(sm_k)%fion, tps, mat_z, f )
& rep(sm_k)%fion, tps, mat_z, f, rhor )
ENDIF
!
@ -1546,7 +1546,7 @@ SUBROUTINE smdmain( tau, fion_out, etot_out, nat_out )
& rep_el(sm_k)%lambda,rep_el(sm_k)%lambdam,xnhe0(sm_k),xnhem(sm_k), &
& vnhe(sm_k),xnhp0(:,sm_k),xnhpm(:,sm_k),vnhp(:,sm_k),nhpcl,ekincm(sm_k), &
& xnhh0,xnhhm,vnhh,velh,ecutp,ecutw,delt,pmass,ibrav,celldm, &
& rep(sm_k)%fion, tps, mat_z, f )
& rep(sm_k)%fion, tps, mat_z, f, rhor )
!
!

View File

@ -644,7 +644,7 @@ MODULE wannier_subroutines
velsm, acc, lambda, lambdam, xnhe0, xnhem, &
vnhe, xnhp0, xnhpm, vnhp, nhpcl, ekincm, &
xnhh0, xnhhm, vnhh, velh, ecut, ecutw, delt, &
celldm, fion, tps, mat_z, occ_f )
celldm, fion, tps, mat_z, occ_f, rho )
!--------------------------------------------------------------------------
!
USE efcalc, ONLY : wf_efield
@ -679,7 +679,7 @@ MODULE wannier_subroutines
REAL(KIND=dbl) :: xnhh0(:,:), xnhhm(:,:), vnhh(:,:)
REAL(KIND=dbl) :: ecut, ecutw, delt, celldm(:)
REAL(KIND=dbl) :: fion(:,:), tps
REAL(KIND=dbl) :: mat_z(:,:,:), occ_f(:)
REAL(KIND=dbl) :: mat_z(:,:,:), occ_f(:), rho(:,:)
!
!
! ... More Wannier Function Options
@ -706,7 +706,7 @@ MODULE wannier_subroutines
tausm, vels, velsm,acc, lambda, lambdam, xnhe0, xnhem, &
vnhe, xnhp0, xnhpm, vnhp, nhpcl, ekincm, xnhh0, xnhhm, &
vnhh, velh, ecut, ecutw, delt, pmass, ibrav, celldm, &
fion, tps, mat_z, occ_f )
fion, tps, mat_z, occ_f, rho )
!
CALL stop_run( .TRUE. )
!

View File

@ -260,7 +260,7 @@ MODULE xml_io_base
!
!------------------------------------------------------------------------
SUBROUTINE write_ions( nsp, nat, atm, ityp, psfile, &
pseudo_dir, amass, tau, if_pos, dirname )
pseudo_dir, amass, tau, if_pos, dirname, pos_unit )
!------------------------------------------------------------------------
!
INTEGER, INTENT(IN) :: nsp, nat
@ -272,6 +272,7 @@ MODULE xml_io_base
REAL(KIND=DP), INTENT(IN) :: amass(:)
REAL(KIND=DP), INTENT(IN) :: tau(:,:)
INTEGER, INTENT(IN) :: if_pos(:,:)
CHARACTER(LEN=*), INTENT(IN) :: pos_unit
!
INTEGER :: i, flen
CHARACTER(LEN=256) :: file_pseudo
@ -311,7 +312,7 @@ MODULE xml_io_base
!
END DO
!
CALL iotk_write_attr( attr, "UNIT", "Bohr", FIRST = .TRUE. )
CALL iotk_write_attr( attr, "UNIT", TRIM(pos_unit), FIRST = .TRUE. )
CALL iotk_write_empty( iunpun, "UNITS_FOR_ATOMIC_POSITIONS", attr )
!
DO i = 1, nat
@ -584,4 +585,113 @@ MODULE xml_io_base
!
END SUBROUTINE write_bz
!
!----------------------------------------------
SUBROUTINE write_rho_xml( rho_file, rho, nr1, nr2, nr3, nr1x, nr2x, ipp, npp )
!
! Writes charge density rho, one plane at a time.
! If ipp and npp are specified planes are collected one by one from
! all processors, avoiding an overall collect of the charge density
! on a single proc.
!
USE kinds, ONLY : dbl
USE io_files, ONLY : rhounit
USE io_global, ONLY : ionode, ionode_id
USE mp, ONLY : mp_sum, mp_get, mp_bcast, mp_max
USE mp_global, ONLY : mpime, nproc, root, me_pool, my_pool_id, &
nproc_pool, intra_pool_comm, root_pool, my_image_id
USE iotk_module
!
IMPLICIT NONE
!
INTEGER, INTENT(IN) :: nr1, nr2, nr3
INTEGER, INTENT(IN) :: nr1x, nr2x
CHARACTER(LEN=*), INTENT(IN) :: rho_file
REAL(dbl), INTENT(IN) :: rho( : )
INTEGER, OPTIONAL, INTENT(IN) :: ipp( : )
INTEGER, OPTIONAL, INTENT(IN) :: npp( : )
!
INTEGER :: ierr, i, j, k, kk, ldr, ip
INTEGER :: ierr_iotk
CHARACTER(LEN=iotk_attlenx) :: attr
REAL(dbl), ALLOCATABLE :: rho_plane( : )
INTEGER, ALLOCATABLE :: kowner( : )
IF( ionode ) THEN
!
CALL iotk_open_write( rhounit, FILE = TRIM( rho_file ), BINARY = .FALSE., ierr = ierr )
!
END IF
!
CALL mp_bcast( ierr, ionode_id )
IF( ierr /= 0 ) &
call errore(" write_rho_xml ", " cannot open rho_file file for writing ", ierr )
!
IF ( ionode ) THEN
!
CALL iotk_write_begin( rhounit, "CHARGE-DENSITY" )
!
CALL iotk_write_attr( attr, "nr1", nr1, FIRST = .TRUE. )
CALL iotk_write_attr( attr, "nr2", nr2 )
CALL iotk_write_attr( attr, "nr3", nr3 )
!
CALL iotk_write_empty( rhounit, "INFO", attr )
!
END IF
ALLOCATE( rho_plane( nr1 * nr2 ) )
ALLOCATE( kowner( nr3 ) )
! find out the owner of each "z" plane
!
IF( PRESENT( ipp ) .AND. PRESENT( npp ) ) THEN
DO ip = 1, nproc
kowner( ipp( ip ) + 1 : ipp( ip ) + npp( ip ) ) = ip - 1
END DO
ELSE
kowner = ionode_id
END IF
ldr = nr1x * nr2x
DO k = 1, nr3
!
IF( kowner( k ) == mpime ) THEN
!
kk = k
!
IF( PRESENT( ipp ) ) kk = k - ipp( mpime + 1 )
!
DO j = 1, nr2
DO i = 1, nr1
rho_plane( i + ( j - 1 ) * nr1 ) = rho( i + ( j - 1 ) * nr1x + ( kk - 1 ) * ldr )
END DO
END DO
END IF
!
IF( kowner( k ) /= ionode_id ) THEN
CALL mp_get( rho_plane, rho_plane, mpime, ionode_id, kowner( k ), k )
END IF
!
IF ( ionode ) &
CALL iotk_write_dat( rhounit, "z" // iotk_index( k ), rho_plane )
END DO
DEALLOCATE( rho_plane )
DEALLOCATE( kowner )
IF ( ionode ) THEN
!
CALL iotk_write_end( rhounit, "CHARGE-DENSITY" )
!
CALL iotk_close_write( rhounit )
!
END IF
RETURN
END SUBROUTINE write_rho_xml
!
!
END MODULE xml_io_base

View File

@ -86,7 +86,7 @@ SUBROUTINE punch()
!
iunpun = 999
!
CALL pw_writefile( '' )
CALL pw_writefile( 'all' )
!
#endif
!

View File

@ -229,7 +229,7 @@ MODULE pw_restart
! ... IONS
!
CALL write_ions( nsp, nat, atm, ityp, &
psfile, pseudo_dir, amass, tau, if_pos, dirname )
psfile, pseudo_dir, amass, tau, if_pos, dirname, "alat" )
!
! ... SYMMETRIES
!
@ -258,13 +258,16 @@ MODULE pw_restart
DEGAUSS = degauss, LTETRA = ltetra, NTETRA = ntetra, &
TFIXED_OCC = tfixed_occ, LSDA = lsda, NELUP = nbnd, &
NELDW = nbnd, F_INP = f_inp )
END IF
!
num_k_points = nkstot
!
IF ( nspin == 2 ) num_k_points = nkstot / 2
!
IF ( ionode ) THEN
!
! ... BRILLOUIN_ZONE
!
num_k_points = nkstot
!
IF ( nspin == 2 ) num_k_points = nkstot / 2
!
CALL write_bz( num_k_points, xk, wk )
!
! ... PARALLELISM
@ -419,7 +422,7 @@ MODULE pw_restart
!
CALL write_wfc( iunout, ik, nkstot, kunit, ispin, nspin, &
evc, npw_g, nbnd, igk_l2g(:,ik-iks+1), &
ngk(ik-iks+1), filename )
ngk(ik-iks+1), filename, 1.0d0 )
!
ik_eff = ik + num_k_points
!
@ -444,7 +447,7 @@ MODULE pw_restart
!
CALL write_wfc( iunout, ik_eff, nkstot, kunit, ispin, nspin, &
evc, npw_g, nbnd, igk_l2g(:,ik_eff-iks+1), &
ngk(ik_eff-iks+1), filename )
ngk(ik_eff-iks+1), filename, 1.0d0 )
!
ELSE
!
@ -480,7 +483,7 @@ MODULE pw_restart
CALL write_wfc( iunout, ik, nkstot, kunit, ispin, nspin, &
evc_nc(:,ipol,:), npw_g, nbnd, &
igk_l2g(:,ik-iks+1), ngk(ik-iks+1), &
filename )
filename, 1.0d0 )
!
END DO
!
@ -488,7 +491,7 @@ MODULE pw_restart
!
CALL write_wfc( iunout, ik, nkstot, kunit, ispin, nspin, &
evc, npw_g, nbnd, igk_l2g(:,ik-iks+1), &
ngk(ik-iks+1), filename )
ngk(ik-iks+1), filename, 1.0d0 )
!
END IF
!
@ -1619,7 +1622,7 @@ MODULE pw_restart
!
!------------------------------------------------------------------------
SUBROUTINE write_wfc( iuni, ik, nk, kunit, ispin, &
nspin, wf0, ngw, nbnd, igl, ngwl, filename )
nspin, wf0, ngw, nbnd, igl, ngwl, filename, scalef )
!------------------------------------------------------------------------
!
USE mp_wave
@ -1638,6 +1641,10 @@ MODULE pw_restart
INTEGER, INTENT(IN) :: igl(:)
CHARACTER(LEN=256), INTENT(IN) :: filename
!
! scale factor, usually 1.0 for pw and 1/SQRT( omega ) CP
!
REAL(KIND=DP), INTENT(IN) :: scalef
!
INTEGER :: i, j, ierr
INTEGER :: nkl, nkr, nkbl, iks, ike, nkt, ikt, igwx
INTEGER :: npool, ipmask(nproc), ipsour
@ -1751,6 +1758,7 @@ MODULE pw_restart
CALL iotk_write_attr( attr, "ispin", ispin )
CALL iotk_write_attr( attr, "nspin", nspin )
CALL iotk_write_attr( attr, "igwx", igwx )
CALL iotk_write_attr( attr, "scale_factor", scalef )
!
CALL iotk_write_empty( iuni, "INFO", attr )
!