scnc
/
quantum-espresso
mirror of https://gitlab.com/QEF/q-e.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

1) An experimental parallelization of NEB images has been implemented. 

This required a deep modification of the parallelism in PWscf:
there are two new communicators (intra_image_comm and inter_image_comm) and the
existing "pool" communicators (intra_pool_comm and inter_pool_comm) are now vectors
of length given by the number of parallel images. #ifdef __PARA is no longer
needed because all "parallel" variables are always initialized for a serial run
and all parallel routines are, in the case of a serial run, dummy routines.
The wrappers to MPI routines used only by PWscf are in the PW/para.f90 file.
The others (mp_***) are in the Modules/mp.f90. All explicit referencies to mpif.h
should be replaced by an "USE parallel_include" (in a serial run parallel_include
is simply a dummy module).

2) The extrapolation of both potential and wavefunctions has been rewritten in
order to be smarter than before: on the basis of the required extrapolation
order, on the basis of the history and on the basis of which files are really
present on the disk, the algorithm chooses the extrapolation order.
All the algorithms in which ions are moved can use the extrapolation.

These are both unstable features: I need the help of everybody to test them.

C.S.


git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@742 c92efa57-630b-4861-b058-cf58834340f0
This commit is contained in:
sbraccia 2004-03-24 09:36:50 +00:00
parent f94c73e766
commit 2a99b6fdac
46 changed files with 1810 additions and 1200 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

13
D3/d3_exc.f90
View File

@ -5,6 +5,7 @@
! in the root directory of the present distribution,
! or http://www.gnu.org/copyleft/gpl.txt .
!
#include "machine.h"
!
!-----------------------------------------------------------------------
subroutine d3_exc
@ -13,17 +14,18 @@ subroutine d3_exc
! Calculates the contribution to the derivative of the dynamical
! matrix due to the third derivative of the exchange and correlation
! energy
#include "machine.h"
USE kinds, only : DP
USE kinds, ONLY : DP
use pwcom
use phcom
use d3com
#ifdef __PARA
use para
USE io_global, ONLY : ionode_id
USE mp_global, ONLY : inter_pool_comm, my_image_id
USE mp, ONLY : mp_bcast
#endif
implicit none
integer :: errcode, ir, ipert, jpert, kpert, npert1, npert2
real (kind = dp) :: d2mxc, rhotot, xq0 (3)
real (kind = dp), allocatable :: d2muxc (:)
@ -74,7 +76,8 @@ subroutine d3_exc
enddo
#ifdef __PARA
100 continue
IF ( npool /= 1 ) CALL mp_bcast( d3dyn1, ionode_id, MPI_COMM_ROW )
IF ( npool /= 1 ) &
CALL mp_bcast( d3dyn1, ionode_id, inter_pool_comm(my_image_id) )
#endif
d3dyn = d3dyn + d3dyn1

15
D3/d3_init.f90
View File

@ -5,21 +5,20 @@
! in the root directory of the present distribution,
! or http://www.gnu.org/copyleft/gpl.txt .
!
#include"machine.h"
!
!-----------------------------------------------------------------------
subroutine d3_init
!-----------------------------------------------------------------------
#include"machine.h"
use pwcom
use phcom
use d3com
#ifdef __PARA
use para
#endif
USE mp, ONLY : mp_barrier
implicit none
#ifdef __PARA
include 'mpif.h'
#endif
integer :: nt, irr, irr1, ipert, imode0, errcode
! counter on atom types
#ifdef DEBUG
@ -130,9 +129,7 @@ subroutine d3_init
endif
#ifdef __PARA
100 continue
call MPI_barrier (MPI_COMM_WORLD, errcode)
call errore ('d3_init', 'at barrier', errcode)
call mp_barrier()
#endif
deallocate(drhoscf)

17
D3/davcio_drho2.f90
View File

@ -19,24 +19,20 @@ subroutine davcio_drho2 (drho, lrec, iunit, nrec, isw)
use pwcom
USE kinds, only : DP
use phcom
#ifdef __PARA
use para
USE io_global, ONLY : ionode_id
USE mp, ONLY : mp_bcast
#endif
USE mp_global, ONLY : intra_pool_comm, my_image_id
USE mp, ONLY : mp_bcast, mp_barrier
implicit none
#ifdef __PARA
include 'mpif.h'
#endif
integer :: iunit, lrec, nrec, isw
complex(kind=DP) :: drho (nrxx)
#ifdef __PARA
!
! local variables
!
integer :: root, errcode, itmp, proc
complex(kind=DP), allocatable :: ddrho (:)
allocate (ddrho( nrx1 * nrx2 * nrx3 ))
@ -47,8 +43,7 @@ subroutine davcio_drho2 (drho, lrec, iunit, nrec, isw)
!
call cgather_sym (drho, ddrho)
root = 0
call MPI_barrier (MPI_COMM_POOL, errcode)
call errore ('davcio_drho2', 'at barrier', errcode)
call mp_barrier()
if (me.eq.1) call davcio (ddrho, lrec, iunit, nrec, + 1)
elseif (isw < 0) then
!
@ -56,7 +51,7 @@ subroutine davcio_drho2 (drho, lrec, iunit, nrec, isw)
! processors of the pool
!
if (me == 1) call davcio (ddrho, lrec, iunit, nrec, - 1)
call mp_bcast( ddrho, ionode_id, MPI_COMM_POOL )
call mp_bcast( ddrho, ionode_id, intra_pool_comm(my_image_id) )
!
! Distributes ddrho between between the tasks of the pool
!

15
D3/drho_drc.f90
View File

@ -5,6 +5,7 @@
! in the root directory of the present distribution,
! or http://www.gnu.org/copyleft/gpl.txt .
!
#include "machine.h"
!
!-----------------------------------------------------------------------
subroutine drho_drc (iudrho_x, u_x, xq_x, drc_x, scale)
@ -13,18 +14,15 @@ subroutine drho_drc (iudrho_x, u_x, xq_x, drc_x, scale)
! it according to the variation of the core_charge
! It is used by drho_cc. Have a look there for more explanation
!
#include "machine.h"
USE kinds, only : DP
USE kinds, ONLY : DP
use pwcom
use phcom
use d3com
#ifdef __PARA
use para
#endif
USE mp, ONLY : mp_barrier
implicit none
#ifdef __PARA
include 'mpif.h'
#endif
integer :: iudrho_x
!input: the unit containing the charge variation
@ -79,8 +77,7 @@ subroutine drho_drc (iudrho_x, u_x, xq_x, drc_x, scale)
enddo
#ifdef __PARA
100 continue
call MPI_barrier (MPI_COMM_WORLD, errcode)
call errore ('drho_drc', 'at barrier', errcode)
call mp_barrier()
#endif
deallocate (drhoc)
deallocate (drhov)

17
D3/stop_d3.f90
View File

@ -16,16 +16,12 @@ subroutine stop_d3 (flag)
use pwcom
use phcom
use d3com
USE io_files, ONLY : iunigk
use mp, only: mp_end
#ifdef __PARA
USE io_files, ONLY : iunigk
use mp, ONLY : mp_end, mp_barrier
use para
#endif
implicit none
#ifdef __PARA
include 'mpif.h'
integer :: info
#endif
logical :: flag
close (unit = iuwfc, status = 'keep')
@ -56,10 +52,7 @@ subroutine stop_d3 (flag)
call print_clock_d3
call show_memory ()
#ifdef __PARA
call mpi_barrier (MPI_COMM_WORLD, info)
! call mpi_finalize (info)
#endif
call mp_barrier()
call mp_end()

235
Modules/check_stop.f90
View File

@ -1,116 +1,141 @@
!
! Copyright (C) 2002 FPMD group
! Copyright (C) 2002-2004 FPMD-PWSCF group
! This file is distributed under the terms of the
! GNU General Public License. See the file `License'
! in the root directory of the present distribution,
! or http://www.gnu.org/copyleft/gpl.txt .
!
! This module contains functions to check if the code should
! be stopped smootly.
! In particular th function check_stop_now return .TRUE. if
! either the user has created a given file or if the
! elapsed time is larger than max_seconds
!
! ... This module contains functions to check if the code should
! ... be stopped smootly.
! ... In particular th function check_stop_now return .TRUE. if
! ... either the user has created a given file or if the
! ... elapsed time is larger than max_seconds
!
!------------------------------------------------------------------------------!
MODULE check_stop
!------------------------------------------------------------------------------!
USE kinds
IMPLICIT NONE
SAVE
REAL(dbl) :: max_seconds = 1.d+7
LOGICAL, PRIVATE :: tinit = .FALSE.
MODULE check_stop
!------------------------------------------------------------------------------!
!
USE kinds
!
IMPLICIT NONE
!
SAVE
!
REAL(dbl) :: max_seconds = 1.D+7
LOGICAL, PRIVATE :: tinit = .FALSE.
!
!
CONTAINS
!
! ... internal procedures
!
!-----------------------------------------------------------------------
SUBROUTINE check_stop_init( val )
!-----------------------------------------------------------------------
!
USE io_global, ONLY : ionode, ionode_id, stdout
USE io_files, ONLY : exit_file, stopunit
!
IMPLICIT NONE
!
REAL(dbl), INTENT(IN) :: val
LOGICAL :: tex
REAL(dbl) :: seconds
REAL(dbl) :: elapsed_seconds
EXTERNAL elapsed_seconds
!
IF( tinit ) &
WRITE( UNIT = stdout, &
FMT = '("WARNING: check_stop already initialized *** ")' )
!
IF ( val > 0.D0 ) max_seconds = val
!
IF ( ionode ) THEN
!
INQUIRE( FILE = exit_file, EXIST = tex )
!
IF ( tex ) THEN
!
OPEN( stopunit, FILE = exit_file, STATUS = 'OLD' )
CLOSE( stopunit, STATUS = 'DELETE' )
!
END IF
!
END IF
!
seconds = elapsed_seconds()
tinit = .TRUE.
!
RETURN
!
END SUBROUTINE
!
!
!-----------------------------------------------------------------------
FUNCTION check_stop_now()
!-----------------------------------------------------------------------
!
USE mp, ONLY : mp_bcast
USE mp_global, ONLY : intra_image_comm
USE io_global, ONLY : ionode, ionode_id, stdout
USE io_files, ONLY : exit_file, stopunit, iunexit
!
IMPLICIT NONE
!
LOGICAL :: check_stop_now, tex
REAL(dbl) :: seconds
REAL(dbl) :: elapsed_seconds
EXTERNAL elapsed_seconds
!
!
! ... elapsed_seconds is a C function returning the elapsed solar
! ... time in seconds since the first call to the function itself
!
IF( .NOT. tinit ) &
CALL errore( 'check_stop_now', 'check_stop not initialized', 1 )
!
check_stop_now = .FALSE.
!
IF ( ionode ) THEN
!
INQUIRE( FILE = TRIM( exit_file ), EXIST = tex )
!
IF ( tex ) THEN
!
check_stop_now = .TRUE.
!
WRITE( UNIT = stdout, &
FMT = '(" *** Program stopped by user request *** ")' )
!
OPEN( UNIT = iunexit, FILE = TRIM( exit_file ) )
CLOSE( UNIT = iunexit, STATUS = 'DELETE' )
!
END IF
!
seconds = elapsed_seconds()
!
IF ( seconds > max_seconds ) THEN
!
check_stop_now = .TRUE.
!
WRITE( UNIT = stdout, &
FMT = '(" *** Maximum CPU time exceeded ***")' )
WRITE( UNIT = stdout, &
FMT = '(" *** max_seconds = ",D10.2," ***")' ) max_seconds
WRITE( UNIT = stdout, &
FMT = '(" *** elapsed seconds = ",D10.2," ***")' ) seconds
!
END IF
!
END IF
!
CALL mp_bcast( check_stop_now, ionode_id, intra_image_comm )
!
RETURN
!
END FUNCTION check_stop_now
!
!------------------------------------------------------------------------------!
SUBROUTINE check_stop_init( val )
USE io_global, ONLY: ionode, ionode_id, stdout
USE io_files, ONLY: exit_file, stopunit
IMPLICIT NONE
REAL(dbl), INTENT(IN) :: val
LOGICAL :: tex
REAL(dbl) :: seconds
REAL(dbl) :: elapsed_seconds
EXTERNAL :: elapsed_seconds
IF( tinit ) THEN
WRITE( stdout, fmt='("WARNING: check_stop already initialized *** ")' )
END IF
IF ( val > 0.0d0 ) max_seconds = val
IF ( ionode ) THEN
INQUIRE( FILE = exit_file, EXIST = tex )
IF ( tex ) THEN
OPEN( stopunit, FILE = exit_file, STATUS = 'OLD' )
CLOSE( stopunit, STATUS = 'DELETE' )
END IF
END IF
seconds = elapsed_seconds()
tinit = .TRUE.
RETURN
END SUBROUTINE
!------------------------------------------------------------------------------!
FUNCTION check_stop_now()
USE mp, ONLY: mp_bcast
USE io_global, ONLY: ionode, ionode_id, stdout
USE io_files, ONLY: exit_file, stopunit
IMPLICIT NONE
LOGICAL :: check_stop_now, tex
REAL(dbl) :: seconds
REAL(dbl) :: elapsed_seconds
EXTERNAL :: elapsed_seconds
! ... elapsed_seconds is a C function returning the elapsed solar time in
! ... seconds since the first call to the function itself
IF( .NOT. tinit ) THEN
CALL errore( ' check_stop_now ', ' check_stop not initialized ', 1 )
END IF
check_stop_now = .FALSE.
IF ( ionode ) THEN
INQUIRE( FILE = TRIM( exit_file ), EXIST = tex )
!!!! WRITE( *, * ) 'DEBUG from check_stop: ', exit_file, tex
IF ( tex ) THEN
check_stop_now = .TRUE.
WRITE( stdout, fmt='(" *** Program stopped by user request *** ")' )
END IF
seconds = elapsed_seconds()
IF( seconds > max_seconds ) THEN
check_stop_now = .TRUE.
WRITE( stdout, fmt='(" *** Maximum CPU time exceeded *** ")' )
WRITE( stdout, fmt='(" *** max_seconds = ", D10.2, " *** ")' ) max_seconds
WRITE( stdout, fmt='(" *** elapsed seconds = ", D10.2, " *** ")' ) seconds
END IF
END IF
CALL mp_bcast( check_stop_now, ionode_id )
RETURN
END FUNCTION check_stop_now
!------------------------------------------------------------------------------!
END MODULE check_stop
END MODULE check_stop
!------------------------------------------------------------------------------!

10
Modules/clocks.f90
View File

@ -179,7 +179,7 @@ SUBROUTINE print_this_clock( n )
USE kinds, ONLY : DP
USE mytime
USE mp, ONLY : mp_max, mp_min
USE mp_global, ONLY : group, inter_pool_comm
USE mp_global, ONLY : group, inter_pool_comm, my_image_id
!
IMPLICIT NONE
!
@ -212,7 +212,7 @@ SUBROUTINE print_this_clock( n )
!
#ifndef DEBUG
CALL mp_max( elapsed_cpu_time, group )
CALL mp_max( elapsed_cpu_time, inter_pool_comm )
CALL mp_max( elapsed_cpu_time, inter_pool_comm(my_image_id) )
#endif
#endif
IF ( n == 1 ) THEN
@ -256,7 +256,7 @@ FUNCTION get_clock( label )
USE kinds, ONLY : DP
USE mytime
USE mp, ONLY : mp_max, mp_min
USE mp_global, ONLY : group, inter_pool_comm
USE mp_global, ONLY : group, intra_image_comm
!
IMPLICIT NONE
!
@ -285,8 +285,8 @@ FUNCTION get_clock( label )
!
! ... In the parallel case, use the maximum over all nodes and pools
!
CALL mp_max( get_clock, group )
CALL mp_max( get_clock, inter_pool_comm )
CALL mp_max( get_clock, intra_image_comm )
! CALL mp_max( get_clock, inter_pool_comm(my_image_id) )
#endif
RETURN
END IF

1
Modules/control_flags.f90
View File

@ -203,6 +203,7 @@
max_cg_iter, &! maximum number of iterations in a CG di
diis_buff, &! dimension of the buffer in diis
diis_ndim, &! dimension of reduced basis in DIIS
history, &! number of old steps available for potential updating
order ! type of potential updating ( see update_pot )
!
LOGICAL, PUBLIC :: &

10
Modules/fft_base.f90
View File

@ -618,7 +618,7 @@ subroutine fft_scatter (f_in, nrx3, nxx_, f_aux, ncp_, npp_, sign)
! The output is overwritten on f_in ; f_aux is used as work space
!
#include "machine.h"
use mp_global, ONLY: nproc_pool, me_pool, intra_pool_comm, nproc
use mp_global, ONLY: nproc_pool, me_pool, intra_pool_comm, nproc, my_image_id
USE kinds, only : DP
implicit none
@ -689,9 +689,9 @@ subroutine fft_scatter (f_in, nrx3, nxx_, f_aux, ncp_, npp_, sign)
!
! step two: communication
!
call mpi_barrier (intra_pool_comm, ierr)
call mpi_barrier (intra_pool_comm(my_image_id), ierr)
call mpi_alltoallv (f_aux(1), sendcount, sdispls, MPI_DOUBLE_COMPLEX, f_in(1), &
recvcount, rdispls, MPI_DOUBLE_COMPLEX, intra_pool_comm, ierr)
recvcount, rdispls, MPI_DOUBLE_COMPLEX, intra_pool_comm(my_image_id), ierr)
if( ABS(ierr) /= 0 ) call errore ('fft_scatter', 'info<>0', ABS(ierr) )
!
else
@ -700,9 +700,9 @@ subroutine fft_scatter (f_in, nrx3, nxx_, f_aux, ncp_, npp_, sign)
!
! step two: communication
!
call mpi_barrier (intra_pool_comm, ierr)
call mpi_barrier (intra_pool_comm(my_image_id), ierr)
call mpi_alltoallv (f_in(1), recvcount, rdispls, MPI_DOUBLE_COMPLEX, f_aux(1), &
sendcount, sdispls, MPI_DOUBLE_COMPLEX, intra_pool_comm, ierr)
sendcount, sdispls, MPI_DOUBLE_COMPLEX, intra_pool_comm(my_image_id), ierr)
if( ABS(ierr) /= 0 ) call errore ('fft_scatter', 'info<>0', ABS(ierr) )
!
! step one: store contiguously the columns

4
Modules/formats.f90
View File

@ -28,8 +28,10 @@ MODULE formats
axsf_fmt = "(A2,6(2X,F14.10))"
!
CHARACTER (LEN=*), PARAMETER :: &
scf_fmt = "(5X,'tcpu = ',F10.2," // &
scf_fmt = "(5X,'tcpu = ',F8.2," // &
& "'; self-consistency for image ', I3)", &
scf_fmt_para = "(5X,'cpu = ',I2,'; tcpu = ',F8.2," // &
& "'; self-consistency for image ', I3)", &
run_output = "(/,5X,'iteration: ',I3,5X,'E activation ='," // &
& " F10.6,5X,'error =',F10.6,/)", &
run_output_T_const = "(/,5X,'iteration: ',I3,5X,'temperature ='," // &

18
Modules/io_base.f90
View File

@ -2646,7 +2646,7 @@
USE mp_wave
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
nproc_pool, intra_pool_comm, root_pool, my_image_id
USE io_global, ONLY: ionode, ionode_id
!
IMPLICIT NONE
@ -2729,7 +2729,7 @@
! get the maximum index within the pool
!
CALL mp_max( igwx, intra_pool_comm )
CALL mp_max( igwx, intra_pool_comm(my_image_id) )
! now notify all procs if an error has been found
!
@ -2745,7 +2745,7 @@
IF( ionode ) WRITE(iuni) ngw, nbnd, ik, nk, kunit, ispin, nspin, scal
IF( ionode ) WRITE(iuni) igwx
! write(200+mpime+ik*10,*) mpime, nproc, root, me_pool, my_pool_id, nproc_pool, intra_pool_comm, root_pool, npool
! write(200+mpime+ik*10,*) mpime, nproc, root, me_pool, my_pool_id, nproc_pool, intra_pool_comm(my_image_id), root_pool, npool
! write(200+mpime+ik*10,*) ngwl, nkbl, kunit, iks, ike, ngw, nbnd, ik, nk, kunit, ispin, nspin, scal, igwx, ierr
! close(200+mpime+ik*10)
@ -2758,7 +2758,7 @@
IF( t0 ) THEN
IF( npool > 1 ) THEN
IF( ( ikt >= iks ) .AND. ( ikt <= ike ) ) THEN
CALL mergewf(wf0(:,j), wtmp, ngwl, igl, me_pool, nproc_pool, root_pool, intra_pool_comm)
CALL mergewf(wf0(:,j), wtmp, ngwl, igl, me_pool, nproc_pool, root_pool, intra_pool_comm(my_image_id))
END IF
IF( ipsour /= ionode_id ) THEN
CALL mp_get( wtmp, wtmp, mpime, ionode_id, ipsour, j )
@ -2778,7 +2778,7 @@
IF( tm ) THEN
IF( npool > 1 ) THEN
IF( ( ikt >= iks ) .AND. ( ikt <= ike ) ) THEN
CALL mergewf(wfm(:,j), wtmp, ngwl, igl, me_pool, nproc_pool, root_pool, intra_pool_comm)
CALL mergewf(wfm(:,j), wtmp, ngwl, igl, me_pool, nproc_pool, root_pool, intra_pool_comm(my_image_id))
END IF
IF( ipsour /= ionode_id ) THEN
CALL mp_get( wtmp, wtmp, mpime, ionode_id, ipsour, j )
@ -2841,7 +2841,7 @@
USE mp_wave
USE mp, ONLY: mp_sum, mp_put, mp_bcast, mp_max, mp_get
USE mp_global, ONLY: mpime, nproc, root, me_pool, my_pool_id, &
nproc_pool, intra_pool_comm, root_pool
nproc_pool, intra_pool_comm, root_pool, my_image_id
USE io_global, ONLY: ionode, ionode_id
!
IMPLICIT NONE
@ -2940,7 +2940,7 @@
! get the maximum index within the pool
!
CALL mp_max( igwx, intra_pool_comm )
CALL mp_max( igwx, intra_pool_comm(my_image_id) )
! now notify all procs if an error has been found
!
@ -2978,7 +2978,7 @@
CALL mp_put( wtmp, wtmp, mpime, ionode_id, ipdest, j )
END IF
IF( ( ikt >= iks ) .AND. ( ikt <= ike ) ) THEN
CALL splitwf(wf0(:,j), wtmp, ngwl, igl, me_pool, nproc_pool, root_pool, intra_pool_comm)
CALL splitwf(wf0(:,j), wtmp, ngwl, igl, me_pool, nproc_pool, root_pool, intra_pool_comm(my_image_id))
END IF
ELSE
CALL splitwf(wf0(:,j), wtmp, ngwl, igl, mpime, nproc, ionode_id)
@ -3017,7 +3017,7 @@
CALL mp_put( wtmp, wtmp, mpime, ionode_id, ipdest, j )
END IF
IF( ( ik >= iks ) .AND. ( ik <= ike ) ) THEN
CALL splitwf(wfm(:,j), wtmp, ngwl, igl, me_pool, nproc_pool, root_pool, intra_pool_comm)
CALL splitwf(wfm(:,j), wtmp, ngwl, igl, me_pool, nproc_pool, root_pool, intra_pool_comm(my_image_id))
END IF
ELSE
CALL splitwf(wfm(:,j), wtmp, ngwl, igl, mpime, nproc, ionode_id)

8
Modules/io_files.f90
View File

@ -69,7 +69,7 @@
INTEGER :: iunwfc = 10 ! unit with wavefunctions
INTEGER :: iunat = 13 ! unit for saving orthogonal atomic wfcs
INTEGER :: iunocc = 14 ! unit for saving the atomic n_{ij}
INTEGER :: iunoldwfc = 11 ! unit with old wavefunctions (molecular dynamics)
INTEGER :: iunoldwfc = 11 ! unit with old wavefunctions
INTEGER :: iunoldwfc2 = 12 ! as above at step -2
INTEGER :: iunigk = 16 ! unit for saving indices
INTEGER :: iunres = 1 ! unit for the restart of the run
@ -78,6 +78,11 @@
INTEGER :: nwordwfc = 2 ! lenght of record in wavefunction file
INTEGER :: nwordatwfc = 2 ! lenght of record in atomic wfc file
!
INTEGER :: iunexit = 26 ! unit for a soft exit
INTEGER :: iunupdate = 27 ! unit for saving old positions (extrapolation)
INTEGER :: iunpara = 28 ! unit for parallelization among images
INTEGER :: iunblock = 29 ! as above (blocking file)
!
! ... NEB specific
!
INTEGER :: iunneb = 6 ! unit for NEB output ( stdout or what else )
@ -86,7 +91,6 @@
INTEGER :: iunint = 23 ! unit for saving the interpolated energy profile
INTEGER :: iunxyz = 24 ! unit for saving coordinates ( xyz format )
INTEGER :: iunaxsf = 25 ! unit for saving coordinates ( axsf format )
INTEGER :: iunexit = 26 ! unit for a soft exit
!=----------------------------------------------------------------------------=!
END MODULE io_files

28
Modules/io_routines.f90
View File

@ -28,8 +28,7 @@ MODULE io_routines
PES_gradient, suspended_image, &
Emax, Emin, Emax_index, &
lquick_min , ldamped_dyn, lmol_dyn
USE mp_global, ONLY : mpime, my_pool_id
USE io_global, ONLY : ionode_id
USE io_global, ONLY : ionode, ionode_id
USE mp, ONLY : mp_bcast
!
IMPLICIT NONE
@ -42,7 +41,7 @@ MODULE io_routines
! ... end of local variables
!
!
IF ( mpime == 0 .AND. my_pool_id == 0 ) THEN
IF ( ionode ) THEN
WRITE( UNIT = iunneb, &
FMT = '(/,5X,"reading file ", A,/)') TRIM( neb_file )
@ -147,20 +146,20 @@ MODULE io_routines
!
! ... broadcast to all nodes
!
CALL mp_bcast( istep, ionode_id )
CALL mp_bcast( nstep, ionode_id )
CALL mp_bcast( istep, ionode_id )
CALL mp_bcast( nstep, ionode_id )
CALL mp_bcast( suspended_image, ionode_id )
!
CALL mp_bcast( pos, ionode_id )
CALL mp_bcast( if_pos, ionode_id )
CALL mp_bcast( PES, ionode_id )
CALL mp_bcast( pos, ionode_id )
CALL mp_bcast( if_pos, ionode_id )
CALL mp_bcast( PES, ionode_id )
CALL mp_bcast( PES_gradient, ionode_id )
!
IF ( lquick_min .OR. ldamped_dyn .OR. lmol_dyn ) &
CALL mp_bcast( vel, ionode_id )
!
CALL mp_bcast( Emax, ionode_id )
CALL mp_bcast( Emin, ionode_id )
CALL mp_bcast( Emax, ionode_id )
CALL mp_bcast( Emin, ionode_id )
CALL mp_bcast( Emax_index, ionode_id )
!
END SUBROUTINE read_restart
@ -178,7 +177,7 @@ MODULE io_routines
lquick_min , ldamped_dyn, lmol_dyn
USE formats, ONLY : energy, restart_first, restart_others, &
velocities
USE mp_global, ONLY : mpime, my_pool_id
USE io_global, ONLY : ionode
!
IMPLICIT NONE
!
@ -189,7 +188,7 @@ MODULE io_routines
! ... end of local variables
!
!
IF ( mpime == 0 .AND. my_pool_id == 0 ) THEN
IF ( ionode ) THEN
OPEN( UNIT = iunrestart, FILE = neb_file, STATUS = "UNKNOWN", &
ACTION = "WRITE" )
@ -287,7 +286,7 @@ MODULE io_routines
USE io_files, ONLY : iundat, iunint, iunxyz, iunaxsf, &
dat_file, int_file, xyz_file, &
axsf_file
USE mp_global, ONLY : mpime, my_pool_id
USE io_global, ONLY : ionode
!
IMPLICIT NONE
!
@ -304,7 +303,8 @@ MODULE io_routines
! ... end of local variables
!
!
IF ( mpime == 0 .AND. my_pool_id == 0 ) THEN
IF ( ionode ) THEN
!
ALLOCATE( d_R( dim ) )
!
ALLOCATE( a( num_of_images - 1 ) )

128
Modules/mp_global.f90
View File

@ -1,55 +1,87 @@
!
! Copyright (C) 2002-2003 PWSCF-FPMD-CP90 group
! Copyright (C) 2002-2004 PWSCF-FPMD-CP90 group
! This file is distributed under the terms of the
! GNU General Public License. See the file `License'
! in the root directory of the present distribution,
! or http://www.gnu.org/copyleft/gpl.txt .
!
MODULE mp_global
USE shmem_include
USE parallel_include
USE kinds
!----------------------------------------------------------------------------
MODULE mp_global
!----------------------------------------------------------------------------
!
USE shmem_include
USE parallel_include
USE kinds
!
IMPLICIT NONE
!
SAVE
!
INTEGER :: mpime = 0 ! absolute processor index starting from 0
INTEGER :: root = 0 ! index of the absolute root processor
INTEGER :: nproc = 1 ! absolute number of processor
INTEGER :: group = 0 ! group communicator
INTEGER :: me_pool = 0 ! index of the processor within a pool (starting from 0 !!! )
INTEGER :: me_image = 0 ! index of the processor within an image (starting from 0 !!! )
INTEGER :: root_pool = 0 ! index of the root processor within a pool (starting from 0 !!! )
INTEGER :: root_image = 0 ! index of the root processor within an image (starting from 0 !!! )
INTEGER :: my_pool_id = 0 ! index of my pool (starting from 0 !!! )
INTEGER :: my_image_id = 0 ! index of my image (starting from 0 !!! )
INTEGER :: npool = 1 ! number of "k-points"-pools
INTEGER :: nimage = 1 ! number of "neb-images"-pools
INTEGER :: nproc_pool = 1 ! number of processor within a pool
INTEGER :: nproc_image = 1 ! number of processor within an image
INTEGER :: inter_image_comm = 0 ! inter image communicator
INTEGER :: intra_image_comm = 0 ! intra image communicator
!
INTEGER, ALLOCATABLE :: inter_pool_comm(:) ! inter pool communicator
INTEGER, ALLOCATABLE :: intra_pool_comm(:) ! intra pool communicator
!
CONTAINS
!
!-----------------------------------------------------------------------
SUBROUTINE mp_global_start( root_i, mpime_i, group_i, nproc_i )
!-----------------------------------------------------------------------
!
IMPLICIT NONE
SAVE
INTEGER :: mpime = 0 ! absolute processor index starting from 0
INTEGER :: root = 0 ! index of the absolute root processor
INTEGER :: nproc = 1 ! absolute number of processor
INTEGER :: group = 0 ! group communicator
INTEGER :: me_pool = 0 ! index of the processor within a pool (starting from 0 !!! )
INTEGER :: root_pool = 0 ! index of the root processor within a pool (starting from 0 !!! )
INTEGER :: my_pool_id = 0 ! index of my pool (starting from 0 !!! )
INTEGER :: nproc_pool = 1 ! number of processor within a pool
INTEGER :: inter_pool_comm = 0 ! inter pool communicator
INTEGER :: intra_pool_comm = 0 ! intra pool communicator
CONTAINS
SUBROUTINE mp_global_start(root_i, mpime_i, group_i, nproc_i )
INTEGER, INTENT(IN) :: root_i, mpime_i, group_i, nproc_i
root = root_i
mpime = mpime_i
group = group_i
nproc = nproc_i
nproc_pool = nproc_i
my_pool_id = 0
me_pool = mpime
root_pool = root
inter_pool_comm = group_i
intra_pool_comm = group_i
RETURN
END SUBROUTINE
SUBROUTINE mp_global_group_start( mep, myp, intrap, interp, nprocp )
INTEGER, INTENT(IN) :: mep, myp, intrap, interp, nprocp
me_pool = mep
my_pool_id = myp
intra_pool_comm = intrap
inter_pool_comm = interp
nproc_pool = nprocp
RETURN
END SUBROUTINE
END MODULE mp_global
!
INTEGER, INTENT(IN) :: root_i, mpime_i, group_i, nproc_i
!
root = root_i
mpime = mpime_i
group = group_i
nproc = nproc_i
nproc_pool = nproc_i
my_pool_id = 0
my_image_id = 0
me_pool = mpime
me_image = mpime
root_pool = root
root_image = root
inter_pool_comm = group_i
intra_pool_comm = group_i
inter_image_comm = group_i
intra_image_comm = group_i
!
RETURN
!
END SUBROUTINE
!
!-----------------------------------------------------------------------
SUBROUTINE mp_global_group_start( mep, myp, nprocp, num_of_pools )
!-----------------------------------------------------------------------
!
IMPLICIT NONE
!
INTEGER, INTENT(IN) :: mep, myp, nprocp, num_of_pools
!
me_pool = mep
my_pool_id = myp
nproc_pool = nprocp
npool = num_of_pools
!
RETURN
!
END SUBROUTINE
!
END MODULE mp_global

96
Modules/neb_base.f90
View File

@ -51,6 +51,7 @@ MODULE neb_base
USE parser, ONLY : int_to_char
USE io_routines, ONLY : read_restart
USE formats, ONLY : stringfmt
USE io_global, ONLY : ionode
!
IMPLICIT NONE
@ -172,31 +173,35 @@ MODULE neb_base
!
CALL compute_deg_of_freedom()
!
! ... details of the calculation are written on output
! ... details of the calculation are written on output (only by ionode)
!
WRITE( UNIT = iunneb, FMT = stringfmt ) &
"calculation", TRIM( calculation )
WRITE( UNIT = iunneb, FMT = stringfmt ) &
"restart_mode", TRIM( restart_mode )
WRITE( UNIT = iunneb, FMT = stringfmt ) &
"CI_scheme", TRIM( CI_scheme )
WRITE( UNIT = iunneb, FMT = stringfmt ) &
"VEC_scheme", TRIM( VEC_scheme )
WRITE( UNIT = iunneb, FMT = stringfmt ) &
"minimization_scheme", TRIM( minimization_scheme )
WRITE( UNIT = iunneb, &
FMT = '(5X,"optimization",T35," = ",L1))' ) optimization
WRITE( UNIT = iunneb, &
FMT = '(5X,"num_of_images",T35," = ",I3)' ) num_of_images
WRITE( UNIT = iunneb, &
FMT = '(5X,"ds",T35," = ",F6.4)' ) ds
WRITE( UNIT = iunneb, &
FMT = '(5X,"k_max",T35," = ",F6.4)' ) k_max
WRITE( UNIT = iunneb, &
FMT = '(5X,"k_min",T35," = ",F6.4)' ) k_min
WRITE( UNIT = iunneb, &
FMT = '(5X,"neb_thr",T35," = ",F6.4)' ) neb_thr
WRITE( UNIT = iunneb, FMT = '(/)' )
IF ( ionode ) THEN
!
WRITE( UNIT = iunneb, FMT = stringfmt ) &
"calculation", TRIM( calculation )
WRITE( UNIT = iunneb, FMT = stringfmt ) &
"restart_mode", TRIM( restart_mode )
WRITE( UNIT = iunneb, FMT = stringfmt ) &
"CI_scheme", TRIM( CI_scheme )
WRITE( UNIT = iunneb, FMT = stringfmt ) &
"VEC_scheme", TRIM( VEC_scheme )
WRITE( UNIT = iunneb, FMT = stringfmt ) &
"minimization_scheme", TRIM( minimization_scheme )
WRITE( UNIT = iunneb, &
FMT = '(5X,"optimization",T35," = ",L1))' ) optimization
WRITE( UNIT = iunneb, &
FMT = '(5X,"num_of_images",T35," = ",I3)' ) num_of_images
WRITE( UNIT = iunneb, &
FMT = '(5X,"ds",T35," = ",F6.4)' ) ds
WRITE( UNIT = iunneb, &
FMT = '(5X,"k_max",T35," = ",F6.4)' ) k_max
WRITE( UNIT = iunneb, &
FMT = '(5X,"k_min",T35," = ",F6.4)' ) k_min
WRITE( UNIT = iunneb, &
FMT = '(5X,"neb_thr",T35," = ",F6.4)' ) neb_thr
WRITE( UNIT = iunneb, FMT = '(/)' )
!
END IF
!
RETURN
!
@ -228,7 +233,6 @@ MODULE neb_base
END SUBROUTINE initialize_neb
!
!
!
!------------------------------------------------------------------------
SUBROUTINE compute_action( langevin_action )
!------------------------------------------------------------------------
@ -619,8 +623,7 @@ MODULE neb_base
!!! workaround for ifc8 compiler internal error
END SUBROUTINE path_tangent_
!
!
!-----------------------------------------------------------------------
SUBROUTINE born_oppenheimer_PES( flag, stat )
!-----------------------------------------------------------------------
@ -696,6 +699,7 @@ MODULE neb_base
ldamped_dyn, lmol_dyn, istep_neb, nstep_neb
USE io_routines, ONLY : write_restart, write_dat_files, write_output
USE check_stop, ONLY : check_stop_now
USE io_global, ONLY : ionode
#if defined (__LANGEVIN)
USE parser, ONLY : int_to_char
#endif
@ -811,6 +815,8 @@ MODULE neb_base
END IF
!
IF ( .NOT. stat ) THEN
!
conv_neb = .FALSE.
!
EXIT minimization
!
@ -885,28 +891,36 @@ MODULE neb_base
SUM( langevin_action ) * ( AU * BOHR_RADIUS_ANGS )
#endif
!
IF ( lmol_dyn ) THEN
IF ( ionode ) THEN
!
WRITE( UNIT = iunneb, FMT = run_output_T_const ) &
istep_neb, &
temp * AU * eV_to_kelvin, &
err * ( AU / BOHR_RADIUS_ANGS )
IF ( lmol_dyn ) THEN
!
WRITE( UNIT = iunneb, FMT = run_output_T_const ) &
istep_neb, &
temp * AU * eV_to_kelvin, &
err * ( AU / BOHR_RADIUS_ANGS )
!
ELSE
!
WRITE( UNIT = iunneb, FMT = run_output ) &
istep_neb, &
( Emax - PES(1) ) * AU, &
err * ( AU / BOHR_RADIUS_ANGS )
!
END IF
!
ELSE
!
WRITE( UNIT = iunneb, FMT = run_output ) &
istep_neb, &
( Emax - PES(1) ) * AU, &
err * ( AU / BOHR_RADIUS_ANGS )
CALL write_output()
!
END IF
!
CALL write_output()
!
!
! ... the program checks if the convergence has been achieved
!
IF ( ( err * AU / BOHR_RADIUS_ANGS ) <= neb_thr ) THEN
!
IF ( ionode ) &
WRITE( UNIT = iunneb, &
FMT = '(/,5X,"NEB convergence achieved")' )
!
conv_neb = .TRUE.
!
EXIT minimization

4
Modules/read_namelists.f90
View File

@ -1268,8 +1268,6 @@ MODULE read_namelists_module
!
! ... CONTROL namelist
!
CALL control_defaults( prog )
!
ios = 0
IF( ionode ) THEN
READ( 5, control, iostat = ios )
@ -1283,7 +1281,7 @@ MODULE read_namelists_module
CALL control_bcast( )
CALL control_checkin( prog )
!
! ... defaults values are changed according to the CONTROL namelist
! ... defaults values are changed here according to the CONTROL namelist
!
CALL fixval( prog )
!

8
Modules/stick_base.f90
View File

@ -36,7 +36,7 @@
SUBROUTINE sticks_maps( tk, ub, lb, b1, b2, b3, gcut, gcutw, gcuts, st, stw, sts )
USE mp, ONLY: mp_sum
USE mp_global, ONLY: me_pool, nproc_pool, intra_pool_comm
USE mp_global, ONLY: me_pool, nproc_pool, intra_pool_comm, my_image_id
LOGICAL, INTENT(IN) :: tk ! if true use the full space grid
INTEGER, INTENT(IN) :: ub(:) ! upper bounds for i-th grid dimension
@ -163,9 +163,9 @@
END IF
CALL mp_sum(st ,intra_pool_comm)
CALL mp_sum(stw ,intra_pool_comm)
CALL mp_sum(sts ,intra_pool_comm)
CALL mp_sum(st ,intra_pool_comm(my_image_id) )
CALL mp_sum(stw ,intra_pool_comm(my_image_id) )
CALL mp_sum(sts ,intra_pool_comm(my_image_id) )
! Test sticks
! WRITE( stdout,*) 'testtesttesttesttesttesttesttesttesttest'

134
PH/davcio_drho.f90
View File

@ -1,79 +1,99 @@
!
! Copyright (C) 2001 PWSCF group
! Copyright (C) 2001-2004 PWSCF group
! This file is distributed under the terms of the
! GNU General Public License. See the file `License'
! in the root directory of the present distribution,
! or http://www.gnu.org/copyleft/gpl.txt .
!
!
!-----------------------------------------------------------------------
subroutine davcio_drho (drho, lrec, iunit, nrec, isw)
!-----------------------------------------------------------------------
!
! reads/writes variation of the charge with respect to a perturbation
! on a file.
! isw = +1 : gathers data from the nodes and writes on a single file
! isw = -1 : reads data from a single file and distributes them
!
#include "machine.h"
!
!----------------------------------------------------------------------------
SUBROUTINE davcio_drho( drho, lrec, iunit, nrec, isw )
!----------------------------------------------------------------------------
!
! ... reads/writes variation of the charge with respect to a perturbation
! ... on a file.
! ... isw = +1 : gathers data from the nodes and writes on a single file
! ... isw = -1 : reads data from a single file and distributes them
!
use pwcom
USE kinds, only : DP
USE kinds, ONLY : DP
use phcom
#ifdef __PARA
use para
USE io_global, ONLY : ionode_id
USE mp, ONLY : mp_bcast
USE pfft, ONLY : npp, ncplane
USE mp_global, ONLY : intra_pool_comm, me_pool, root_pool, my_image_id
USE mp, ONLY : mp_bcast, mp_barrier
USE parallel_include
#endif
implicit none
integer :: iunit, lrec, nrec, isw
complex(kind=DP) :: drho (nrxx, nspin)
!
IMPLICIT NONE
!
INTEGER :: iunit, lrec, nrec, isw
COMPLEX(KIND=DP) :: drho(nrxx,nspin)
!
#ifdef __PARA
!
! local variables
! ... local variables
!
integer :: root, errcode, itmp, proc, is
complex(kind=DP), allocatable :: ddrho (:,:)
allocate (ddrho( nrx1 * nrx2 * nrx3 , nspin))
if (isw == 1) then
INTEGER :: itmp, proc, is, dim
COMPLEX(KIND=DP), ALLOCATABLE :: ddrho(:,:)
!
!
ALLOCATE( ddrho( nrx1 * nrx2 * nrx3 , nspin) )
!
IF ( isw == 1 ) THEN
!
! First task of the first pool is the only task allowed to write
! the file
! ... First task of each pool is the only task allowed to write
! ... the file
!
do is = 1, nspin
call cgather_sym (drho (1, is), ddrho (1, is) )
enddo
root = 0
call MPI_barrier (MPI_COMM_WORLD, errcode)
call errore ('davcio_drho', 'at barrier', errcode)
if (me == 1) call davcio (ddrho, lrec, iunit, nrec, + 1)
elseif (isw < 0) then
DO is = 1, nspin
!
CALL cgather_sym( drho(1,is), ddrho(1,is) )
!
END DO
!
! First task of the pool reads ddrho, and broadcasts to all the
! processors of the pool
call mp_barrier()
!
if (me == 1) call davcio (ddrho, lrec, iunit, nrec, - 1)
call mp_bcast( ddrho, ionode_id, MPI_COMM_POOL )
IF ( me_pool == root_pool ) CALL davcio( ddrho, lrec, iunit, nrec, + 1 )
!
! Distributes ddrho between between the tasks of the pool
ELSE IF ( isw < 0 ) THEN
!
! ... First task of the pool reads ddrho, and broadcasts to all the
! ... processors of the pool
!
IF ( me_pool == root_pool ) CALL davcio( ddrho, lrec, iunit, nrec, - 1 )
!
CALL mp_bcast( ddrho, root_pool, intra_pool_comm(my_image_id) )
!
! ... Distributes ddrho between between the tasks of the pool
!
itmp = 1
do proc = 1, me-1
itmp = itmp + ncplane * npp (proc)
enddo
do is = 1, nspin
drho (:, is) = (0.d0, 0.d0)
call ZCOPY (ncplane * npp (me), ddrho (itmp, is), 1, drho (1, is), 1)
enddo
endif
deallocate(ddrho)
!
DO proc = 1, me_pool
!
itmp = itmp + ncplane * npp(proc)
!
END DO
!
dim = ncplane * npp(me_pool+1)
!
DO is = 1, nspin
!
drho(:,is) = ( 0.D0, 0.D0 )
!
drho(1:dim,is) = ddrho(itmp:itmp+dim,is)
!CALL ZCOPY( ncplane*npp(me_pool+1), ddrho(itmp,is), 1, drho(1,is), 1 )
!
END DO
!
END IF
!
DEALLOCATE( ddrho )
!
#else
call davcio (drho, lrec, iunit, nrec, isw)
!
CALL davcio( drho, lrec, iunit, nrec, isw )
! !
#endif
return
end subroutine davcio_drho
!
RETURN
!
END SUBROUTINE davcio_drho

10
PP/pw2wan.f90
View File

@ -11,9 +11,9 @@ program wannier
!
USE io_global, ONLY : stdout, ionode
USE mp_global, ONLY : mpime
USE mp, ONLY : mp_bcast
USE mp, ONLY : mp_bcast
use pwcom
use para, only : kunit
use para, ONLY : kunit
use io_files
!
implicit none
@ -110,7 +110,7 @@ subroutine write_wannier (nk, s0, kunit, ispinw)
use io_base, only : write_restart_wfc
use io_global, only : ionode
use mp_global, only : nproc, nproc_pool, mpime
use mp_global, only : my_pool_id, intra_pool_comm, inter_pool_comm
use mp_global, only : my_pool_id, my_image_id, intra_pool_comm
use mp, only : mp_sum, mp_max
@ -171,7 +171,7 @@ subroutine write_wannier (nk, s0, kunit, ispinw)
! find out the global number of G vectors: ngm_g
ngm_g = ngm
call mp_sum( ngm_g , intra_pool_comm )
call mp_sum( ngm_g, intra_pool_comm(my_image_id) )
allocate ( ei_k ( nbnd, nkstot ) ) ! eigenvectors
allocate ( ei_kw( nbnd, nkstot/nspin ) ) ! eigenvectors
@ -238,7 +238,7 @@ subroutine write_wannier (nk, s0, kunit, ispinw)
itmp( 2, ig_l2g( ig ) ) = ig2( ig )
itmp( 3, ig_l2g( ig ) ) = ig3( ig )
end do
call mp_sum( itmp , intra_pool_comm )
call mp_sum( itmp, intra_pool_comm(my_image_id) )
! write G space parameters and vectors
if( ionode ) then

3
PW/add_efield.f90
View File

@ -47,6 +47,7 @@
USE io_global, ONLY : stdout
USE control_flags, ONLY: mixing_beta
#ifdef __PARA
USE mp_global, ONLY : intra_image_comm
use para
use mp
#endif
@ -116,7 +117,7 @@
dipold=dip
endif
#ifdef __PARA
call mp_bcast(dip,0)
call mp_bcast(dip,0,intra_image_comm)
#endif
if (.not.dipfield) then
etotefield=-2.d0*dipion*eamp*omega/fpi

21
PW/cdiagh.f90
View File

@ -15,11 +15,8 @@ SUBROUTINE cdiagh( n, h, ldh, e, v )
! ... hermitean matrix H . On output, the matrix is unchanged
!
USE kinds, ONLY : DP
#if defined (__PARA)
USE para, ONLY : me, mypool, npool, MPI_COMM_POOL
USE io_global, ONLY : ionode_id
USE mp_global, ONLY : npool, me_pool, root_pool, intra_pool_comm, my_image_id
USE mp, ONLY : mp_bcast
#endif
!
IMPLICIT NONE
!
@ -98,7 +95,7 @@ SUBROUTINE cdiagh( n, h, ldh, e, v )
!
! ... only the first processor diagonalize the matrix
!
IF ( me == 1 ) THEN
IF ( me_pool == root_pool ) THEN
!
# endif
!
@ -108,8 +105,8 @@ SUBROUTINE cdiagh( n, h, ldh, e, v )
!
END IF
!
CALL mp_bcast( e, ionode_id, MPI_COMM_POOL )
CALL mp_bcast( v, ionode_id, MPI_COMM_POOL )
CALL mp_bcast( e, root_pool, intra_pool_comm(my_image_id) )
CALL mp_bcast( v, root_pool, intra_pool_comm(my_image_id) )
!
# endif
!
@ -131,7 +128,7 @@ SUBROUTINE cdiagh( n, h, ldh, e, v )
!
! ... local variables (Cray Eispack/Scilib version)
!
INTEGER :: i, j, k, info
INTEGER :: i, j, k, info
REAL(KIND=DP) :: ar(ldh,n), ai(ldh,n), zr(ldh,n), zi(ldh,n)
! real and imaginary part of h(ldh,n) and of v(ldh,n)
! (used as auxiliary arrays)
@ -161,7 +158,7 @@ SUBROUTINE cdiagh( n, h, ldh, e, v )
!
IMPLICIT NONE
! workaround for Intel ifc8 bug:
COMPLEX(KIND=DP) :: v(ldh,n)
COMPLEX(KIND=DP) :: v(ldh,n)
!
! ... local variables (LAPACK version)
!
@ -207,7 +204,7 @@ SUBROUTINE cdiagh( n, h, ldh, e, v )
!
! ... else only the first processor diagonalize the matrix
!
IF ( me == 1 ) THEN
IF ( me_pool == root_pool ) THEN
# endif
!
! ... allocate workspace
@ -230,8 +227,8 @@ SUBROUTINE cdiagh( n, h, ldh, e, v )
!
END IF
!
CALL mp_bcast( e, ionode_id, MPI_COMM_POOL )
CALL mp_bcast( v, ionode_id, MPI_COMM_POOL )
CALL mp_bcast( e, root_pool, intra_pool_comm(my_image_id) )
CALL mp_bcast( v, root_pool, intra_pool_comm(my_image_id) )
!
# endif
!

10
PW/cdiaghg.f90
View File

@ -18,9 +18,7 @@ SUBROUTINE cdiaghg( n, m, h, s, ldh, e, v )
! ... LAPACK version - uses both ZHEGV and ZHEGVX
!
USE kinds, ONLY : DP
USE para, ONLY : me, npool
USE io_global, ONLY : ionode_id
USE mp_global, ONLY : intra_pool_comm
USE mp_global, ONLY : npool, me_pool, root_pool, intra_pool_comm, my_image_id
USE mp, ONLY : mp_bcast
!
IMPLICIT NONE
@ -115,7 +113,7 @@ SUBROUTINE cdiaghg( n, m, h, s, ldh, e, v )
!
! ... only the first processor diagonalize the matrix
!
IF ( me == 1 ) THEN
IF ( me_pool == root_pool ) THEN
!
IF ( all_eigenvalues ) THEN
!
@ -144,8 +142,8 @@ SUBROUTINE cdiaghg( n, m, h, s, ldh, e, v )
!
! ... broadcast the eigenvectors and the eigenvalues
!
CALL mp_bcast( e, ionode_id, intra_pool_comm )
CALL mp_bcast( v, ionode_id, intra_pool_comm )
CALL mp_bcast( e, root_pool, intra_pool_comm(my_image_id) )
CALL mp_bcast( v, root_pool, intra_pool_comm(my_image_id) )
!
! ... deallocate workspace
!

19
PW/close_files.f90
View File

@ -11,11 +11,10 @@ SUBROUTINE close_files()
!
! ... Close all files and synchronize processes for a new scf calculation.
!
USE control_flags, ONLY : order
USE io_files, ONLY : prefix, iunwfc, iunoldwfc, iunoldwfc2, iunigk
#ifdef __PARA
USE mp, ONLY : mp_barrier
#endif
USE control_flags, ONLY : order
USE io_files, ONLY : prefix, iunwfc, iunoldwfc, iunoldwfc2, iunigk
USE mp_global, ONLY : intra_image_comm
USE mp, ONLY : mp_barrier
!
IMPLICIT NONE
!
@ -27,19 +26,11 @@ SUBROUTINE close_files()
!
CLOSE( UNIT = iunwfc, STATUS = 'KEEP' )
!
IF ( order > 1 ) &
CLOSE( UNIT = iunoldwfc, STATUS = 'KEEP' )
!
IF ( order > 2 ) &
CLOSE( UNIT = iunoldwfc2, STATUS = 'KEEP' )
!
! ... iunigk is kept open during the execution - close and remove
!
CLOSE( UNIT = iunigk, STATUS = 'DELETE' )
!
#ifdef __PARA
CALL mp_barrier()
#endif
CALL mp_barrier( intra_image_comm )
!
#ifdef __T3E
!

249
PW/compute_scf.f90
View File

@ -9,10 +9,14 @@
SUBROUTINE compute_scf( N_in, N_fin, stat )
!----------------------------------------------------------------------------
!
! ... this subroutine is the main scf-driver for NEB
! ... ( called by Modules/neb_base.f90, born_oppenheimer() subroutine )
!
USE kinds, ONLY : DP
USE input_parameters, ONLY : if_pos, sp_pos, startingwfc, startingpot
USE input_parameters, ONLY : if_pos, sp_pos, startingwfc, startingpot, &
diago_thr_init
USE constants, ONLY : e2
USE control_flags, ONLY : conv_elec, istep, alpha0, beta0
USE control_flags, ONLY : conv_elec, istep, history, alpha0, beta0, ethr
USE check_stop, ONLY : check_stop_now
USE brilz, ONLY : alat
USE basis, ONLY : tau, ityp, nat, &
@ -23,16 +27,16 @@ SUBROUTINE compute_scf( N_in, N_fin, stat )
USE relax, ONLY : if_pos_ => if_pos
USE extfield, ONLY : tefield, forcefield
USE io_files, ONLY : prefix, tmp_dir, &
iunneb, iunexit, exit_file
iunneb, iunupdate
USE io_global, ONLY : stdout
USE formats, ONLY : scf_fmt
USE formats, ONLY : scf_fmt, scf_fmt_para
USE neb_variables, ONLY : pos, PES, PES_gradient, num_of_images, &
dim, suspended_image, istep_neb
USE parser, ONLY : int_to_char
USE para, ONLY : me, mypool
USE io_global, ONLY : ionode_id
USE mp_global, ONLY : inter_pool_comm, intra_pool_comm
USE mp, ONLY : mp_bcast, mp_barrier
USE io_global, ONLY : ionode
USE mp_global, ONLY : inter_image_comm, intra_image_comm, &
my_image_id, me_image, root_image, nimage
USE mp, ONLY : mp_bcast, mp_barrier, mp_sum
!
IMPLICIT NONE
!
@ -43,7 +47,7 @@ SUBROUTINE compute_scf( N_in, N_fin, stat )
!
! ... local variables definition
!
INTEGER :: image, ia
INTEGER :: image, ia, istat
REAL(KIND=DP), ALLOCATABLE :: tauold(:,:,:)
! previous positions of atoms
REAL (KIND=DP) :: tcpu
@ -59,36 +63,60 @@ SUBROUTINE compute_scf( N_in, N_fin, stat )
! ... end of external functions definition
!
!
istep = istep_neb
! ... all processes are syncronized (needed to have an ordered output)
!
stat = .TRUE.
CALL mp_barrier( intra_image_comm )
CALL mp_barrier( inter_image_comm )
!
istep = istep_neb
istat = 0
!
! ... only the first cpu on each image needs the tauold vector
!
IF ( me_image == root_image ) ALLOCATE( tauold( 3, nat, 3 ) )
!
tmp_dir_saved = tmp_dir
!
IF ( me == 1 .AND. mypool == 1 ) THEN
!
ALLOCATE( tauold( 3, nat, 3 ) )
!
END IF
! ... vectors PES and PES_gradient are initalized to zero for all images on
! ... all nodes: this is needed for the final mp_sum()
!
PES(N_in:N_fin) = 0.D0
PES_gradient(:,N_in:N_fin) = 0.D0
!
! ... only the first cpu initializes the file needed by parallelization
! ... among images
!
IF ( ionode ) CALL para_file_init()
!
image = N_in + my_image_id
!
DO image = N_in, N_fin
!
scf_loop: DO
!
suspended_image = image
!
IF( check_stop_now() ) THEN
IF ( check_stop_now() ) THEN
!
stat = .FALSE.
istat = 1
!
RETURN
EXIT scf_loop
!
END IF
!
tmp_dir = TRIM( tmp_dir_saved ) // TRIM( prefix ) // "_" // &
TRIM( int_to_char( image ) ) // "/"
!
tcpu = get_clock( 'PWSCF' )
!
WRITE( UNIT = iunneb, FMT = scf_fmt ) tcpu, image
IF ( nimage > 1 ) THEN
!
WRITE( UNIT = iunneb, FMT = scf_fmt_para ) my_image_id, tcpu, image
!
ELSE
!
WRITE( UNIT = iunneb, FMT = scf_fmt ) tcpu, image
!
END IF
!
CALL clean_pw()
!
@ -96,7 +124,7 @@ SUBROUTINE compute_scf( N_in, N_fin, stat )
!
! ... unit stdout is connected to the appropriate file
!
IF ( me == 1 .AND. mypool == 1 ) THEN
IF ( me_image == root_image ) THEN
!
INQUIRE( UNIT = stdout, OPENED = opnd )
IF ( opnd ) CLOSE( UNIT = stdout )
@ -122,41 +150,49 @@ SUBROUTINE compute_scf( N_in, N_fin, stat )
!
CALL init_run()
!
IF ( me == 1 .AND. mypool == 1 ) THEN
IF ( me_image == root_image ) THEN
!
! ... the file containing old positions is opened
! ... ( needed for extrapolation )
!
CALL seqopn( 4, TRIM( prefix ) // '.update', 'FORMATTED', file_exists )
CALL seqopn( iunupdate, TRIM( prefix ) // '.update', &
'FORMATTED', file_exists )
!
IF ( file_exists ) THEN
!
READ( UNIT = 4, FMT = * ) tauold
READ( UNIT = iunupdate, FMT = * ) history
READ( UNIT = iunupdate, FMT = * ) tauold
!
ELSE
!
tauold = 0.D0
history = 0
tauold = 0.D0
!
END IF
!
CLOSE( UNIT = 4, STATUS = 'KEEP' )
CLOSE( UNIT = iunupdate, STATUS = 'KEEP' )
!
! ... find the best coefficients for the extrapolation of the potential
!
! WRITE( UNIT = *, FMT = * )
! WRITE( UNIT = *, FMT = * ) tau(1,:)
! WRITE( UNIT = *, FMT = * )
! WRITE( UNIT = *, FMT = * ) tauold(1,:,:)
!
CALL find_alpha_and_beta( nat, tau, tauold, alpha0, beta0 )
!
END IF
!
IF ( me == 1 ) CALL mp_bcast( alpha0, ionode_id, inter_pool_comm )
IF ( me == 1 ) CALL mp_bcast( beta0, ionode_id, inter_pool_comm )
!
CALL mp_bcast( alpha0, ionode_id, intra_pool_comm )
CALL mp_bcast( beta0, ionode_id, intra_pool_comm )
CALL mp_bcast( alpha0, root_image, intra_image_comm )
CALL mp_bcast( beta0, root_image, intra_image_comm )
CALL mp_bcast( history, root_image, intra_image_comm )
!
! ... potential and wavefunctions are extrapolated
!
CALL update_pot()
!
! ... self-consistency loop
!
CALL electrons()
!
IF ( .NOT. conv_elec ) THEN
@ -164,12 +200,14 @@ SUBROUTINE compute_scf( N_in, N_fin, stat )
WRITE( iunneb, '(/,5X,"WARNING : scf convergence NOT achieved",/, &
& 5X,"stopping in compute_scf()...",/)' )
!
stat = .FALSE.
istat = 1
!
RETURN
EXIT scf_loop
!
END IF
!
! ... self-consistent forces
!
CALL forces()
!
! ... energy is converted from rydberg to hartree
@ -182,7 +220,7 @@ SUBROUTINE compute_scf( N_in, N_fin, stat )
PES_gradient(:,image) = - RESHAPE( SOURCE = force, &
SHAPE = (/ dim /) ) / e2
!
IF ( me == 1 .AND. mypool == 1 ) THEN
IF ( me_image == root_image ) THEN
!
! ... save the previous two steps ( a total of three steps is saved )
!
@ -190,33 +228,152 @@ SUBROUTINE compute_scf( N_in, N_fin, stat )
tauold(:,:,2) = tauold(:,:,1)
tauold(:,:,1) = tau(:,:)
!
CALL seqopn( 4, TRIM( prefix ) // '.update', 'FORMATTED', file_exists )
history = MIN( 3, ( history + 1 ) )
!
WRITE( UNIT = 4, FMT = * ) tauold
CALL seqopn( iunupdate, &
& TRIM( prefix ) // '.update', 'FORMATTED', file_exists )
!
CLOSE( UNIT = 4, STATUS = 'KEEP' )
WRITE( UNIT = iunupdate, FMT = * ) history
WRITE( UNIT = iunupdate, FMT = * ) tauold
!
CLOSE( UNIT = iunupdate, STATUS = 'KEEP' )
!
! ... the new image is obtained
!
CALL get_new_image( image )
!
END IF
!
CALL mp_bcast( image, root_image, intra_image_comm )
!
! ... input values are restored at the end of each iteration
!
startingpot_ = startingpot
startingwfc_ = startingwfc
!
ethr = diago_thr_init
!
CALL reset_k_points()
!
END DO
! ... exit if finished
!
IF ( image > N_fin ) EXIT scf_loop
!
END DO scf_loop
!
IF ( me == 1 .AND. mypool == 1 ) THEN
!
DEALLOCATE( tauold )
!
END IF
IF ( me_image == root_image ) DEALLOCATE( tauold )
!
tmp_dir = tmp_dir_saved
!
suspended_image = 0
!
CALL mp_barrier( intra_image_comm )
CALL mp_barrier( inter_image_comm )
!
! ... PES and PES_gradient are communicated among "image" pools
!
CALL mp_sum( PES(N_in:N_fin), inter_image_comm )
CALL mp_sum( PES_gradient(:,N_in:N_fin), inter_image_comm )
CALL mp_sum( istat, inter_image_comm )
!
! ... global status is computed here
!
IF ( istat == 0 ) THEN
!
stat = .TRUE.
!
ELSE
!
stat = .FALSE.
!
END IF
!
RETURN
!
CONTAINS
!
! ... internal procedures
!
SUBROUTINE para_file_init()
!
! ... this subroutine initializes the file needed for the
! ... parallelization among images
!
USE io_files, ONLY : iunpara
USE mp_global, ONLY : nimage
!
IMPLICIT NONE
!
!
OPEN( UNIT = iunpara, FILE = TRIM( tmp_dir_saved ) // &
& TRIM( prefix ) // '.para' , STATUS = 'UNKNOWN' )
!
WRITE( iunpara, * ) N_in + nimage
!
CLOSE( UNIT = iunpara, STATUS = 'KEEP' )
!
RETURN
!
END SUBROUTINE para_file_init
!
SUBROUTINE get_new_image( image )
!
! ... this subroutine is used to get the new image to work on
! ... the *.BLOCK file is other needed to avoid that other jobs
! ... try to read/write on file "para"
!
USE io_files, ONLY : iunpara, iunblock
USE mp_global, ONLY : my_image_id
!
IMPLICIT NONE
!
INTEGER, INTENT(OUT) :: image
INTEGER :: ioerr
LOGICAL :: opened, exists
!
!
open_loop: DO
!
OPEN( UNIT = iunblock, FILE = TRIM( tmp_dir_saved ) // &
& TRIM( prefix ) // '.BLOCK' , IOSTAT = ioerr, STATUS = 'NEW' )
!
IF ( ioerr > 0 ) CYCLE open_loop
!
INQUIRE( UNIT = iunpara, OPENED = opened )
!
IF ( .NOT. opened ) THEN
!
INQUIRE( FILE = TRIM( tmp_dir_saved ) // &
& TRIM( prefix ) // '.para', EXIST = exists )
!
IF ( exists ) THEN
!
OPEN( UNIT = iunpara, FILE = TRIM( tmp_dir_saved ) // &
& TRIM( prefix ) // '.para' , STATUS = 'OLD' )
!
READ( iunpara, * ) image
!
CLOSE( UNIT = iunpara, STATUS = 'DELETE' )
!
OPEN( UNIT = iunpara, FILE = TRIM( tmp_dir_saved ) // &
& TRIM( prefix ) // '.para' , STATUS = 'NEW' )
!
WRITE( iunpara, * ) image + 1
!
CLOSE( UNIT = iunpara, STATUS = 'KEEP' )
!
EXIT open_loop
!
END IF
!
END IF
!
END DO open_loop
!
CLOSE( UNIT = iunblock, STATUS = 'DELETE' )
!
RETURN
!
END SUBROUTINE get_new_image
!
END SUBROUTINE compute_scf

6
PW/data_structure.f90
View File

@ -28,7 +28,7 @@ subroutine data_structure( lgamma )
ncp0, ncp0s, nxx, nxxs, nct, ncts, ncp, ncps
#endif
use mp, only: mp_sum
use mp_global, only: intra_pool_comm, nproc_pool, me_pool
use mp_global, only: intra_pool_comm, nproc_pool, me_pool, my_image_id
use stick_base
use fft_scalar, only: good_fft_dimension
use fft_types, only: fft_dlay_allocate, fft_dlay_set, fft_dlay_scalar
@ -376,8 +376,8 @@ subroutine data_structure( lgamma )
ngms_l = ngms
ngm_g = ngm
ngms_g = ngms
call mp_sum( ngm_g , intra_pool_comm )
call mp_sum( ngms_g, intra_pool_comm )
call mp_sum( ngm_g , intra_pool_comm(my_image_id) )
call mp_sum( ngms_g, intra_pool_comm(my_image_id) )
return

25
PW/divide_et_impera.f90
View File

@ -7,18 +7,19 @@
!
#include "machine.h"
!
!-----------------------------------------------------------------------
!----------------------------------------------------------------------------
SUBROUTINE divide_et_impera( xk, wk, isk, lsda, nkstot, nks )
!-----------------------------------------------------------------------
!----------------------------------------------------------------------------
!
! ... This routine divides the k points (with granularity kunit) among no
! ... This routine divides the k points (with granularity kunit) among nodes
! ... and sets the variable nkstot equal to the total number of k-points
!
#if defined (__PARA)
!
USE kinds, ONLY : DP
USE para, ONLY : mypool, npool, kunit
USE io_global, only : stdout
USE kinds, ONLY : DP
USE mp_global, ONLY : my_pool_id, npool
USE para, ONLY : kunit
!
IMPLICIT NONE
!
@ -47,22 +48,24 @@ SUBROUTINE divide_et_impera( xk, wk, isk, lsda, nkstot, nks )
!
rest = ( nkstot - nks * npool ) / kunit
!
IF ( mypool <= rest ) nks = nks + kunit
IF ( ( my_pool_id + 1 ) <= rest ) nks = nks + kunit
!
! ... calculates nbase = the position in the list of the first point that
! ... belong to this npool - 1
!
nbase = nks * ( mypool - 1 )
nbase = nks * my_pool_id
!
IF ( mypool > rest ) nbase = nbase + rest * kunit
IF ( ( my_pool_id + 1 ) > rest ) nbase = nbase + rest * kunit
!
! ... displaces these points in the first positions of the list
!
IF ( nbase > 0 ) THEN
!
xk(:,1:nks) = xk(:,nbase:nbase+nks)
xk(:,1:nks) = xk(:,nbase+1:nbase+nks)
!
IF ( lsda ) isk(1:nks) = isk(nbase:nbase+nks)
wk(1:nks) = wk(nbase+1:nbase+nks)
!
IF ( lsda ) isk(1:nks) = isk(nbase+1:nbase+nks)
!
END IF
!

137
PW/electrons.f90
View File

@ -7,11 +7,6 @@
!
#include "machine.h"
!
! ... uncomment the following line to use the "old ethr"
! ... for the first iteration
!
!#define __OLDSTYLE
!
!----------------------------------------------------------------------------
SUBROUTINE electrons()
!----------------------------------------------------------------------------
@ -57,10 +52,8 @@ SUBROUTINE electrons()
USE extfield, ONLY : tefield, etotefield
USE bp, ONLY : lberry
USE wavefunctions_module, ONLY : evc
#if defined (__PARA)
USE para, ONLY : me, mypool, npp, ncplane
USE mp, ONLY : mp_barrier
#endif
USE mp_global, ONLY : me_image, root_image
USE para, ONLY : npp, ncplane
!
IMPLICIT NONE
!
@ -69,7 +62,7 @@ SUBROUTINE electrons()
#if defined (__PARA)
INTEGER :: &
ngkp(npk) ! number of plane waves summed on all nodes
#define NRXX ncplane*npp(me)
#define NRXX ncplane*npp(me_image+1)
! This is needed in mix_pot whenever nproc is not a divisor of nr3.
#else
#define NRXX nrxx
@ -83,8 +76,8 @@ SUBROUTINE electrons()
mag, &! local magnetization
magtot, &! total magnetization
absmag, &! total absolute magnetization
tcpu !
INTEGER :: &
tcpu ! cpu time
INTEGER :: &
i, &! counter on polarization
ir, &! counter on the mesh points
ig, &!
@ -102,11 +95,10 @@ SUBROUTINE electrons()
charge_new !
REAL (KIND=DP) :: &
ethr_min ! minimal threshold for diagonalization at the first scf
! iteration of a MD calculation
! iteration
REAL (KIND=DP), EXTERNAL :: ewald, get_clock
LOGICAL :: &
exst, &!
file_exists ! .TRUE. if a soft exit has been required
exst
!
!
CALL start_clock( 'electrons' )
@ -125,6 +117,7 @@ SUBROUTINE electrons()
! ...jump to the end
!
IF ( output_drho /= ' ' ) CALL remove_atomic_rho
!
CALL stop_clock( 'electrons' )
!
RETURN
@ -141,10 +134,15 @@ SUBROUTINE electrons()
g, gg, ngm, gcutm, gstart, gamma_only, strf )
!
IF ( reduce_io ) THEN
!
flmix = ' '
!
ELSE
!
flmix = 'flmix'
!
END IF
!
END IF
!
! ... Convergence threshold for iterative diagonalization
@ -152,15 +150,7 @@ SUBROUTINE electrons()
! ... for the first scf iteration of each ionic step (except than for the
! ... first) the threshold is fixed to a default value of 1.D-5
!
#if defined (__OLDSTYLE)
IF ( .FALSE. ) THEN
#else
IF ( istep > 1 ) THEN
#endif
!
ethr = 1.D-5
!
END IF
IF ( istep > 1 ) ethr = 1.D-5
!
!%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
!%%%%%%%%%%%%%%%%%%%% iterate ! %%%%%%%%%%%%%%%%%%%%%
@ -171,17 +161,24 @@ SUBROUTINE electrons()
DO idum = 1, niter
!
tcpu = get_clock( 'PWSCF' )
!
WRITE( stdout, 9000 ) tcpu
!
IF ( imix >= 0 ) CALL DCOPY( ( nspin * nrxx), rho, 1, rho_save, 1 )
!
IF ( imix >= 0 ) rho_save = rho
!IF ( imix >= 0 ) CALL DCOPY( ( nspin * nrxx), rho, 1, rho_save, 1 )
!
iter = iter + 1
!
IF ( lscf ) THEN
!
WRITE( stdout, 9010 ) iter, ecutwfc, mixing_beta
!
ELSE
!
WRITE( stdout, 9009 )
!
END IF
!
#if defined (FLUSH)
CALL flush( stdout )
#endif
@ -207,6 +204,8 @@ SUBROUTINE electrons()
!
END IF
!
! ... diagonalziation of the KS hamiltonian
!
CALL c_bands( iter, ik_, dr2 )
!
! ... skip all the rest if not lscf
@ -215,17 +214,20 @@ SUBROUTINE electrons()
!
conv_elec = .TRUE.
!
#if defined (__PARA)
CALL poolrecover( et, nbnd, nkstot, nks )
#endif
!
DO ik = 1, nkstot
!
IF ( lsda ) THEN
!
IF ( ik == 1 ) WRITE( stdout, 9015 )
IF ( ik == ( 1 + nkstot / 2 ) ) WRITE( stdout, 9016 )
!
END IF
!
WRITE( stdout, 9020 ) ( xk(i,ik), i = 1, 3 )
WRITE( stdout, 9030 ) ( et(ibnd,ik) * 13.6058, ibnd = 1, nbnd )
!
END DO
!
! ... do a Berry phase polarization calculation if required
@ -235,6 +237,7 @@ SUBROUTINE electrons()
! ... jump to the end
!
IF ( output_drho /= ' ' ) CALL remove_atomic_rho()
!
CALL stop_clock( 'electrons' )
!
RETURN
@ -254,7 +257,7 @@ SUBROUTINE electrons()
IF ( lda_plus_u ) CALL write_ns()
!
IF ( iter == 1 .AND. lda_plus_u .AND. &
startingpot=='atomic' .AND. istep == 1 ) CALL ns_adj()
startingpot == 'atomic' .AND. istep == 1 ) CALL ns_adj()
!
! ... calculate total and absolute magnetization
!
@ -272,15 +275,12 @@ SUBROUTINE electrons()
IF ( lda_plus_u .AND. iter <= niter_with_fixed_ns ) THEN
!
ldim2 = ( 2 * Hubbard_lmax + 1 )**2
CALL DCOPY( ( ldim2 * nspin * nat ), ns, 1, nsnew, 1 )
nsnew = ns
!CALL DCOPY( ( ldim2 * nspin * nat ), ns, 1, nsnew, 1 )
!
END IF
!
#if defined (__OLDSTYLE)
IF ( .FALSE. ) THEN
#else
IF ( iter == 1 ) THEN
#endif
!
! ... for the first scf iteration ethr_min is set for a check
! ... in mix_rho ( in mix_rho ethr_min = dr2 * ethr_min )
@ -302,11 +302,7 @@ SUBROUTINE electrons()
! ... for the first scf iteration it is controlled that the threshold
! ... is small enought for the diagonalization to be adequate
!
#if defined (__OLDSTYLE)
IF ( .FALSE. ) THEN
#else
IF ( iter == 1 .AND. ethr >= ethr_min ) THEN
#endif
!
! ... a new diagonalization is needed
!
@ -342,7 +338,8 @@ SUBROUTINE electrons()
IF ( lda_plus_u .AND. iter <= niter_with_fixed_ns ) THEN
!
ldim2 = ( 2 * Hubbard_lmax + 1 )**2
CALL DCOPY( ( ldim2 * nspin * nat ), ns, 1, nsnew, 1 )
nsnew = ns
!CALL DCOPY( ( ldim2 * nspin * nat ), ns, 1, nsnew, 1 )
!
END IF
!
@ -353,7 +350,8 @@ SUBROUTINE electrons()
!
END IF
!
CALL DAXPY( ( nspin * nrxx ), -1.D0, vr, 1, vnew, 1 )
vnew = vnew - vr
!CALL DAXPY( ( nspin * nrxx ), -1.D0, vr, 1, vnew, 1 )
!
CALL v_of_rho( rho_save, rho_core, nr1, nr2, nr3, nrx1, nrx2, nrx3, &
nrxx, nl, ngm, gstart, nspin, g, gg, alat, omega, &
@ -381,17 +379,19 @@ SUBROUTINE electrons()
!
ldim2 = ( 2 * Hubbard_lmax + 1 )**2
!
IF ( iter > niter_with_fixed_ns .AND. imix < 0 ) &
CALL DCOPY( ( ldim2 * nspin * nat ), nsnew, 1, ns, 1 )
#if defined (__PARA)
IF ( me == 1 .AND. mypool == 1 ) THEN
#endif
CALL seqopn( iunocc, TRIM( prefix )//'.occup', 'formatted', exst )
IF ( iter > niter_with_fixed_ns .AND. imix < 0 ) ns = nsnew
! CALL DCOPY( ( ldim2 * nspin * nat ), nsnew, 1, ns, 1 )
!
IF ( me_image == root_image ) THEN
!
CALL seqopn( iunocc, TRIM( prefix )//'.occup', 'FORMATTED', exst )
!
WRITE( iunocc, * ) ns
!
CLOSE( UNIT = iunocc, STATUS = 'KEEP' )
#if defined (__PARA)
!
END IF
#endif
!
END IF
!
! ... In the US case we need to recompute the self consistent term in
@ -401,8 +401,9 @@ SUBROUTINE electrons()
!
! ... write the potential (and rho) on file
!
IF ( imix >= 0 ) CALL io_pot( +1, TRIM( prefix )//'.rho', rho_save, nspin )
CALL io_pot( +1, TRIM( prefix )//'.pot', vr, nspin )
IF ( imix >= 0 ) CALL io_pot( 1, TRIM( prefix )//'.rho', rho_save, nspin )
!
CALL io_pot( 1, TRIM( prefix )//'.pot', vr, nspin )
!
! ... save converged wfc if they have not been written previously
!
@ -419,20 +420,22 @@ SUBROUTINE electrons()
!IF ( lda_plus_u ) CALL write_ns()
!
#if defined (__PARA)
DO ik = 1, nks
ngkp(ik) = ngk(ik)
END DO
!
ngkp(1:nks) = ngk(1:nks)
!
CALL ireduce( nks, ngkp )
CALL ipoolrecover( ngkp, 1, nkstot, nks )
CALL poolrecover( et, nbnd, nkstot, nks )
!
#endif
!
DO ik = 1, nkstot
!
IF ( lsda ) THEN
!
IF ( ik == 1 ) WRITE( stdout, 9015)
IF ( ik == ( 1 + nkstot / 2 ) ) WRITE( stdout, 9016)
!
END IF
!
IF ( conv_elec ) THEN
@ -466,9 +469,13 @@ SUBROUTINE electrons()
( iswitch <= 2 ) ) THEN
!
IF ( imix >= 0 ) THEN
!
WRITE( stdout, 9081 ) etot, dr2
!
ELSE
!
WRITE( stdout, 9086 ) etot, dr2
!
END IF
!
WRITE( stdout, 9060 ) &
@ -481,17 +488,25 @@ SUBROUTINE electrons()
ELSE IF ( conv_elec .AND. iswitch > 2 ) THEN
!
IF ( imix >= 0 ) THEN
!
WRITE( stdout, 9081 ) etot, dr2
!
ELSE
!
WRITE( stdout, 9086 ) etot, dr2
!
END IF
!
ELSE
!
IF ( imix >= 0 ) THEN
!
WRITE( stdout, 9080 ) etot, dr2
!
ELSE
!
WRITE( stdout, 9085 ) etot, dr2
!
END IF
!
END IF
@ -520,7 +535,8 @@ SUBROUTINE electrons()
!
!CALL forces()
!
IF ( imix >= 0 ) CALL DCOPY( ( nspin * nrxx), rho_save, 1, rho, 1 )
IF ( imix >= 0 ) rho = rho_save
!IF ( imix >= 0 ) CALL DCOPY( ( nspin * nrxx), rho_save, 1, rho, 1 )
!
END DO
!
@ -584,7 +600,8 @@ SUBROUTINE electrons()
!
DO ir = 1, nrxx
!
mag = rho(ir,1) - rho(ir,2)
mag = rho(ir,1) - rho(ir,2)
!
magtot = magtot + mag
absmag = absmag + ABS( mag )
!
@ -593,10 +610,8 @@ SUBROUTINE electrons()
magtot = magtot * omega / ( nr1 * nr2 * nr3 )
absmag = absmag * omega / ( nr1 * nr2 * nr3 )
!
#if defined (__PARA)
CALL reduce( 1, magtot )
CALL reduce( 1, absmag )
#endif
!
RETURN
!
@ -621,7 +636,9 @@ SUBROUTINE electrons()
tcpu, max_seconds
END IF
!
CALL stop_pw( .FALSE. )
conv_elec = .FALSE.
!
RETURN
!
END IF
!
@ -644,7 +661,9 @@ SUBROUTINE electrons()
& 5X,"stopping in electrons ...",/)' )
END IF
!
CALL stop_pw( .FALSE. )
conv_elec = .FALSE.
!
RETURN
!
END IF
!

18
PW/error.f90
View File

@ -25,16 +25,17 @@ subroutine errore (routin, messag, ierr)
!
USE io_global, ONLY : stdout
USE kinds
USE parallel_include
!
implicit none
#ifdef __PARA
include 'mpif.h'
#endif
!
character (len=*) :: routin, messag
! the name of the calling routine
! the output message
! the name of the calling routine
! the output message
integer :: ierr
! the error flag
! the error flag
!
!
if (ierr.eq.0) return
WRITE( * , * ) ' '
WRITE( * , '(1x,78("%"))' )
@ -56,7 +57,7 @@ subroutine errore (routin, messag, ierr)
call flush (6)
#endif
#ifdef __PARA
call mpi_abort (MPI_COMM_WORLD, ierr, ierr)
CALL MPI_ABORT( MPI_COMM_WORLD, ierr, ierr )
#endif
stop 2
else
@ -64,4 +65,3 @@ subroutine errore (routin, messag, ierr)
return
endif
end subroutine errore

156
PW/init_pool.f90
View File

@ -9,57 +9,153 @@
SUBROUTINE init_pool()
!----------------------------------------------------------------------------
!
! ... This routine initialize the pool
! ... This routine initialize the pool : MPI division in pools and images
!
#if defined (__PARA)
!
USE para, ONLY : me, mypool, npool, nprocp, MPI_COMM_POOL, MPI_COMM_ROW
USE mp, ONLY : mp_barrier
USE para, ONLY : me, mypool, npool, nprocp
USE mp, ONLY : mp_barrier, mp_bcast
USE mp_global, ONLY : mpime, me_image, my_image_id, nproc, &
nproc_image, nproc_pool, nimage, me_pool, my_pool_id, &
intra_image_comm, inter_image_comm, &
intra_pool_comm, inter_pool_comm
USE mp_global, ONLY : mp_global_group_start
USE parallel_include
!
IMPLICIT NONE
!
! ... MPI pool division in pools
!
INTEGER :: ierr, rank
INTEGER :: ierr = 0, numtask, taskid
!
!
! ... set "mypool" and reset "me"
!
rank = me
mypool = ( me - 1 ) / nprocp + 1
me = me - ( mypool - 1 ) * nprocp
! ... communicators are allocated
!
ALLOCATE( inter_pool_comm( 0 : nimage - 1 ) )
ALLOCATE( intra_pool_comm( 0 : nimage - 1 ) )
!
! ... and initialized
!
inter_pool_comm(:) = 0
intra_pool_comm(:) = 0
!
#if defined (__PARA)
!
! ... set "my_image_id", "mypool" and reset "me"
!
! ... my_image_id = 0 : nimage
! ... mypool = 1 : npool ==> nimage * npool * nprocp = nproc
! ... me = 1 : nprocp
!
! ... number of cpus per image
!
nproc_image = nproc / nimage
!
IF ( MOD( nproc, nimage ) /= 0 ) &
CALL errore( 'startup', 'nproc /= nproc_image * nimage', 1 )
!
my_image_id = INT( REAL( mpime ) / REAL( nproc_image ) )
me = MOD( mpime, nproc_image )
me_image = me
!
CALL mp_barrier()
!
CALL MPI_comm_split( MPI_COMM_WORLD, mypool, rank, MPI_COMM_POOL, ierr )
CALL MPI_COMM_SPLIT( MPI_COMM_WORLD, &
my_image_id, mpime, intra_image_comm, ierr )
!
CALL errore( 'init_pool', 'MPI_COMM_POOL is wrong', ierr )
CALL errore( 'init_pool', 'intra_image_comm is wrong', ierr )
!
IF ( npool > 1 ) THEN
!
CALL MP_barrier()
!
CALL MPI_comm_split( MPI_COMM_WORLD, me, rank, MPI_COMM_ROW, ierr )
!
call errore( 'init_pool', 'MPI_COMM_ROW is wrong', ierr )
!
END IF
CALL mp_barrier()
!
CALL MPI_COMM_SPLIT( MPI_COMM_WORLD, me, mpime, inter_image_comm, ierr )
!
CALL errore( 'init_pool', 'inter_image_comm is wrong', ierr )
!
! ... number of cpus per pool
!
nproc_pool = nproc_image / npool
!
IF ( MOD( nproc, npool ) /= 0 ) &
CALL errore( 'startup', 'nproc /= nproc_pool * npool', 1 )
!
mypool = INT( REAL( me ) / REAL( nproc_pool ) )
me = MOD( me, nproc_pool )
!
! ... This is added for compatibility with PVM notations
! ... parent process (source) will have me=1 - child process me=2,...,NPROC
!
me = me + 1
mypool = mypool + 1
!
CALL mp_barrier( intra_image_comm )
!
CALL MPI_COMM_SPLIT( intra_image_comm, mypool, &
me_image, intra_pool_comm(my_image_id), ierr )
!
CALL errore( 'init_pool', 'intra_pool_comm is wrong', ierr )
!
CALL mp_barrier( intra_image_comm )
!
CALL MPI_COMM_SPLIT( intra_image_comm, me, &
me_image, inter_pool_comm(my_image_id), ierr )
!
call errore( 'init_pool', 'inter_pool_comm is wrong', ierr )
!
! ... compatibility with old PWscf routines
!
nprocp = nproc_pool
!
! ... Initialize globally accessible pool variables
!
! ... me => me_pool + 1
! ... mypool => my_pool_id + 1
! ... MPI_COMM_POOL => intra_pool_comm
! ... MPI_COMM_ROW => inter_pool_comm
! ... nprocp => nproc_pool
! ... me => me_pool + 1
! ... mypool => my_pool_id + 1
! ... nprocp => nproc_pool
!
CALL mp_global_group_start( ( me - 1 ), ( mypool - 1 ), &
MPI_COMM_POOL, MPI_COMM_ROW, nprocp )
CALL mp_global_group_start( ( me - 1 ), ( mypool - 1 ), nprocp, npool )
!
#endif
!
#if defined (__NEW_PARALLEL_DEBUG)
PRINT *, ""
PRINT *, "MPIME = ", MPIME
PRINT *, "ME = ", ME
PRINT *, "MY_IMAGE_ID = ", MY_IMAGE_ID
PRINT *, "MYPOOL = ", MYPOOL
PRINT *, ""
CALL MPI_COMM_RANK( intra_image_comm, taskid, ierr )
CALL MPI_COMM_SIZE( intra_image_comm, numtask, ierr )
PRINT *, "intra_image_comm : ", taskid, numtask
CALL MPI_COMM_RANK( inter_image_comm, taskid, ierr )
CALL MPI_COMM_SIZE( inter_image_comm, numtask, ierr )
PRINT *, "inter_image_comm : ", taskid, numtask
CALL MPI_COMM_RANK( intra_pool_comm(my_image_id), taskid, ierr )
CALL MPI_COMM_SIZE( intra_pool_comm(my_image_id), numtask, ierr )
PRINT *, "intra_pool_comm : ", taskid, numtask
CALL MPI_COMM_RANK( inter_pool_comm(my_image_id), taskid, ierr )
CALL MPI_COMM_SIZE( inter_pool_comm(my_image_id), numtask, ierr )
PRINT *, "inter_pool_comm : ", taskid, numtask
PRINT *, ""
ierr = 666
IF ( me_image == 0 ) ierr = my_image_id
PRINT *, "1:MPIME = ", MPIME, " IERR = ", IERR
CALL mp_bcast( ierr, 0, inter_pool_comm(my_image_id) )
IF ( ierr /= 666 ) ierr = ierr + my_pool_id
PRINT *, "2:MPIME = ", MPIME, " IERR = ", IERR
CALL mp_bcast( ierr, 0, intra_pool_comm(my_image_id) )
IF ( ierr /= 666 ) ierr = ierr + me_pool
PRINT *, "3:MPIME = ", MPIME, " IERR = ", IERR
ierr = 999
IF ( me_image == 0 ) ierr = 100 + my_image_id
PRINT *, "4:MPIME = ", MPIME, " IERR = ", IERR
CALL mp_bcast( ierr, 0, intra_pool_comm(my_image_id) )
IF ( ierr /= 999 ) ierr = ierr + me_pool
PRINT *, "5:MPIME = ", MPIME, " IERR = ", IERR
CALL mp_bcast( ierr, 0, inter_pool_comm(my_image_id) )
IF ( ierr /= 999 ) ierr = ierr + my_pool_id
PRINT *, "6:MPIME = ", MPIME, " IERR = ", IERR
PRINT *, ""
!
CALL stop_pw( .FALSE. )
!
#endif
!
RETURN
!
END SUBROUTINE init_pool

50
PW/input.f90
View File

@ -22,7 +22,7 @@ SUBROUTINE iosys()
!
!
USE constants, ONLY : AU, eV_to_kelvin
USE para, ONLY : npool
USE mp_global, ONLY : npool
!
USE io_global, ONLY : stdout
USE bp, ONLY : nppstr_ => nppstr, &
@ -1227,10 +1227,11 @@ SUBROUTINE verify_tmpdir()
!-----------------------------------------------------------------------
!
USE input_parameters, ONLY : restart_mode
USE control_flags, ONLY : lneb
USE control_flags, ONLY : lneb
USE io_files, ONLY : prefix, tmp_dir, nd_nmbr
USE neb_variables, ONLY : num_of_images
USE para, ONLY : me, mypool
USE mp_global, ONLY : mpime
USE io_global, ONLY : ionode
USE mp, ONLY : mp_barrier
!
USE parser, ONLY : int_to_char, delete_if_present
@ -1238,7 +1239,7 @@ SUBROUTINE verify_tmpdir()
IMPLICIT NONE
!
INTEGER :: l, ios, image
CHARACTER (LEN=80) :: tmp_dir_saved
CHARACTER (LEN=80) :: file_path, tmp_dir_saved
INTEGER :: c_mkdir
EXTERNAL c_mkdir
!
@ -1256,8 +1257,10 @@ SUBROUTINE verify_tmpdir()
!
ios = 0
!
OPEN( UNIT = 4, FILE = TRIM( tmp_dir ) // 'pwscf' // nd_nmbr, &
STATUS = 'UNKNOWN', FORM = 'UNFORMATTED', IOSTAT = ios )
file_path = TRIM( tmp_dir ) // 'pwscf'
!
OPEN( UNIT = 4, FILE = TRIM( file_path ) // TRIM( int_to_char( mpime ) ), &
& STATUS = 'UNKNOWN', FORM = 'UNFORMATTED', IOSTAT = ios )
CLOSE( UNIT = 4, STATUS = 'DELETE' )
!
IF ( ios /= 0 ) CALL errore( 'outdir: ', TRIM( tmp_dir ) // &
@ -1266,21 +1269,23 @@ SUBROUTINE verify_tmpdir()
! ... if starting from scratch all temporary files are removed
! ... from tmp_dir ( only by the master node )
!
file_path = TRIM( tmp_dir ) // TRIM( prefix )
!
IF ( restart_mode == 'from_scratch' ) THEN
!
IF ( me == 1 .AND. mypool == 1 ) THEN
IF ( ionode ) THEN
!
! ... wfc-extrapolation file is removed
!
CALL delete_if_present( TRIM( tmp_dir ) // TRIM( prefix ) // '.update' )
CALL delete_if_present( TRIM( file_path ) // '.update' )
!
! ... MD restart file is removed
!
CALL delete_if_present( TRIM( tmp_dir ) // TRIM( prefix ) // '.md' )
CALL delete_if_present( TRIM( file_path ) // '.md' )
!
! ... BFGS rstart file is removed
!
CALL delete_if_present( TRIM( tmp_dir ) // TRIM( prefix ) // '.bfgs' )
CALL delete_if_present( TRIM( file_path ) // '.bfgs' )
!
END IF
!
@ -1291,6 +1296,15 @@ SUBROUTINE verify_tmpdir()
! ... created
!
IF ( lneb ) THEN
!
IF ( ionode ) THEN
!
! ... files needed by parallelization among images are removed
!
CALL delete_if_present( TRIM( file_path ) // '.BLOCK' )
CALL delete_if_present( TRIM( file_path ) // '.para' )
!
END IF
!
tmp_dir_saved = tmp_dir
!
@ -1300,7 +1314,7 @@ SUBROUTINE verify_tmpdir()
tmp_dir = TRIM( tmp_dir_saved ) // TRIM( prefix ) //"_" // &
TRIM( int_to_char( image ) ) // '/'
!
IF ( me == 1 .AND. mypool == 1 ) THEN
IF ( ionode ) THEN
!
! ... a scratch directory for this image of the elastic band is
! ... created ( only by the master node )
@ -1315,9 +1329,10 @@ SUBROUTINE verify_tmpdir()
!
! ... each job checks whether the scratch directory is accessible
! ... or not
!
OPEN( UNIT = 4, FILE = TRIM( tmp_dir ) // 'pwscf' // nd_nmbr, &
STATUS = 'UNKNOWN', FORM = 'UNFORMATTED', IOSTAT = ios )
!
OPEN( UNIT = 4, FILE = TRIM( tmp_dir ) // TRIM( prefix ) // &
& TRIM( int_to_char( mpime ) ), &
& STATUS = 'UNKNOWN', FORM = 'UNFORMATTED', IOSTAT = ios )
CLOSE( UNIT = 4, STATUS = 'DELETE' )
!
IF ( ios /= 0 ) &
@ -1329,7 +1344,12 @@ SUBROUTINE verify_tmpdir()
!
IF ( restart_mode == 'from_scratch' ) THEN
!
IF ( me == 1 .AND. mypool == 1 ) THEN
IF ( ionode ) THEN
!
! ... wfc-extrapolation file is removed
!
CALL delete_if_present( TRIM( tmp_dir ) // &
& TRIM( prefix ) // '.update' )
!
! ... standard output of the self-consistency is removed
!

34
PW/io_pot.f90
View File

@ -11,24 +11,19 @@
SUBROUTINE io_pot( iop, filename, pot, nc )
!----------------------------------------------------------------------------
!
! ... This routine reads (iop=-1) or write (iop=1) the potential
! ... in real space onto a file
!
! ... rewritten to use mp wrappers
! ... This routine reads ( iop = - 1 ) or write ( iop = + 1 ) the
! ... potential in real space onto a file
!
USE kinds, ONLY : DP
USE gvect, ONLY : nrxx, nrx1, nrx2, nrx3
#if defined (__PARA)
USE para, ONLY : me, mypool
USE mp_global, ONLY : intra_pool_comm, inter_pool_comm
USE io_global, ONLY : ionode_id
USE mp, ONLY : mp_bcast, mp_gather
#endif
USE mp_global, ONLY : inter_pool_comm, me_pool, &
root_pool, me_image, root_image, my_image_id
USE mp, ONLY : mp_bcast
!
IMPLICIT NONE
!
INTEGER :: iop, nc, ic
! option: write if +1, read if -1
! option: write if + 1, read if - 1
! number of components and index for them
CHARACTER (LEN=*) :: filename
REAL(KIND=DP) :: pot(nrxx,nc)
@ -42,13 +37,11 @@ SUBROUTINE io_pot( iop, filename, pot, nc )
!
! ... parallel case
!
IF ( me == 1 ) ALLOCATE( allv( nrx1*nrx2*nrx3, nc ) )
IF ( me_pool == root_pool ) ALLOCATE( allv( nrx1*nrx2*nrx3, nc ) )
!
! ... On writing: gather the potential on the first node of each pool
!
IF ( iop == 1 ) THEN
!
! ... CALL mp_gather( pot(:,:), allv(:,:), ionode_id, intra_pool_comm )
!
DO ic = 1, nc
!
@ -58,7 +51,7 @@ SUBROUTINE io_pot( iop, filename, pot, nc )
!
END IF
!
IF ( me == 1 .AND. mypool == 1 ) THEN
IF ( me_image == root_image ) THEN
!
! ... Only the first node of the first pool reads or writes the file
!
@ -81,9 +74,10 @@ SUBROUTINE io_pot( iop, filename, pot, nc )
! ... On reading: copy the potential on the first node of all pools
! ... scatter the potential on all nodes of each pool
!
IF ( iop == -1 ) THEN
IF ( iop == - 1 ) THEN
!
IF ( me == 1 ) CALL mp_bcast( allv, ionode_id, inter_pool_comm )
IF ( me_pool == root_pool ) &
CALL mp_bcast( allv, root_pool, inter_pool_comm(my_image_id) )
!
DO ic = 1, nc
!
@ -93,7 +87,7 @@ SUBROUTINE io_pot( iop, filename, pot, nc )
!
END IF
!
IF ( me == 1 ) DEALLOCATE( allv )
IF ( me_pool == root_pool ) DEALLOCATE( allv )
!
#else
!
@ -117,7 +111,7 @@ SUBROUTINE io_pot( iop, filename, pot, nc )
!
RETURN
!
10 CALL errore( 'io_pot', 'error writing '//filename, 1 )
20 CALL errore( 'io_pot', 'error reading '//filename, 2 )
10 CALL errore( 'io_pot', 'error writing ' // filename, 1 )
20 CALL errore( 'io_pot', 'error reading ' // filename, 2 )
!
END SUBROUTINE io_pot

4
PW/mix_rho.f90
View File

@ -17,13 +17,13 @@ SUBROUTINE mix_rho( rhout, rhoin, nsout, nsin, alphamix, dr2, ethr, ethr_min, &
! ... d.d. johnson prb 38, 12807 (1988)
! ... On output: the mixed density is in rhoin, rhout is UNCHANGED
!
USE kinds, ONLY : DP
USE kinds, ONLY : DP
USE basis, ONLY : nat
USE gvect, ONLY : nr1, nr2, nr3, nrx1, nrx2, nrx3, nrxx, &
nl, nlm
USE ldaU, ONLY : lda_plus_u, Hubbard_lmax
USE lsda_mod, ONLY : nspin
USE control_flags, ONLY : imix, ngm0, tr2
USE control_flags, ONLY : imix, ngm0, tr2
USE wvfct, ONLY : gamma_only
USE wavefunctions_module, ONLY : psic
!

92
PW/move_ions.f90
View File

@ -34,7 +34,7 @@ SUBROUTINE move_ions()
! ... also computed here
!
USE io_global, ONLY : stdout
USE io_files, ONLY : tmp_dir, prefix
USE io_files, ONLY : tmp_dir, prefix, iunupdate
USE bfgs_module, ONLY : lbfgs_ndim, new_bfgs => bfgs, lin_bfgs
USE kinds, ONLY : DP
USE brilz, ONLY : alat, at, bg
@ -45,12 +45,11 @@ SUBROUTINE move_ions()
USE ener, ONLY : etot
USE force_mod, ONLY : force
USE control_flags, ONLY : upscale, lbfgs, loldbfgs, lconstrain, &
lmd, conv_ions, alpha0, beta0, tr2
lmd, conv_ions, history, alpha0, beta0, tr2
USE relax, ONLY : epse, epsf, starting_scf_threshold
USE cellmd, ONLY : lmovecell, calc
USE mp_global, ONLY : inter_pool_comm, intra_pool_comm
USE para, ONLY : me, mypool, npool
USE io_global, ONLY : ionode_id
USE mp_global, ONLY : intra_image_comm
USE io_global, ONLY : ionode_id, ionode
USE mp, ONLY : mp_bcast
!
! ... external procedures
@ -74,7 +73,7 @@ SUBROUTINE move_ions()
!
! ... only one node does the calculation in the parallel case
!
IF ( me == 1 .AND. mypool == 1 ) THEN
IF ( ionode ) THEN
!
conv_ions = .FALSE.
!
@ -82,25 +81,26 @@ SUBROUTINE move_ions()
!
! ... constrains are imposed here
!
IF ( lconstrain ) &
CALL impose_constrains()
IF ( lconstrain ) CALL impose_constrains()
!
! ... the file containing old positions is opened
! ... ( needed for extrapolation )
!
CALL seqopn( 4, TRIM( prefix ) // '.update', 'FORMATTED', exst )
CALL seqopn( iunupdate, TRIM( prefix ) // '.update', 'FORMATTED', exst )
!
IF ( exst ) THEN
!
READ( UNIT = 4, FMT = * ) tauold
READ( UNIT = iunupdate, FMT = * ) history
READ( UNIT = iunupdate, FMT = * ) tauold
!
ELSE
!
tauold = 0.D0
history = 0
tauold = 0.D0
!
END IF
!
CLOSE( UNIT = 4, STATUS = 'KEEP' )
CLOSE( UNIT = iunupdate, STATUS = 'KEEP' )
!
! ... save the previous two steps ( a total of three steps is saved )
!
@ -123,7 +123,7 @@ SUBROUTINE move_ions()
ALLOCATE( pos( 3 * nat ) )
ALLOCATE( gradient( 3 * nat ) )
!
pos = RESHAPE( SOURCE = tau, SHAPE = (/ 3 * nat /) ) * alat
pos = RESHAPE( SOURCE = tau, SHAPE = (/ 3 * nat /) ) * alat
gradient = - RESHAPE( SOURCE = force, SHAPE = (/ 3 * nat /) )
!
IF ( lbfgs_ndim == 1 ) THEN
@ -157,8 +157,8 @@ SUBROUTINE move_ions()
!
END IF
!
tau = RESHAPE( SOURCE = pos, SHAPE = (/ 3 , nat /) ) / alat
force = - RESHAPE( SOURCE = gradient, SHAPE = (/ 3 , nat /) )
tau = RESHAPE( SOURCE = pos, SHAPE = (/ 3, nat /) ) / alat
force = - RESHAPE( SOURCE = gradient, SHAPE = (/ 3, nat /) )
!
CALL output_tau( conv_ions )
!
@ -196,11 +196,14 @@ SUBROUTINE move_ions()
!
CALL find_alpha_and_beta( nat, tau, tauold, alpha0, beta0 )
!
CALL seqopn( 4, TRIM( prefix ) // '.update', 'FORMATTED', exst )
history = MIN( 3, ( history + 1 ) )
!
WRITE( UNIT = 4, FMT = * ) tauold
CALL seqopn( iunupdate, TRIM( prefix ) // '.update', 'FORMATTED', exst )
!
CLOSE( UNIT = 4, STATUS = 'KEEP' )
WRITE( UNIT = iunupdate, FMT = * ) history
WRITE( UNIT = iunupdate, FMT = * ) tauold
!
CLOSE( UNIT = iunupdate, STATUS = 'KEEP' )
!
DEALLOCATE( tauold )
!
@ -208,17 +211,14 @@ SUBROUTINE move_ions()
!
! ... broadcast calculated quantities to all nodes
!
CALL mp_bcast( conv_ions, ionode_id )
CALL mp_bcast( tau, ionode_id )
CALL mp_bcast( force, ionode_id )
CALL mp_bcast( tr2, ionode_id )
CALL mp_bcast( conv_ions, ionode_id )
!
IF ( me == 1 ) CALL mp_bcast( alpha0, ionode_id, inter_pool_comm )
IF ( me == 1 ) CALL mp_bcast( beta0, ionode_id, inter_pool_comm )
!
CALL mp_bcast( alpha0, ionode_id, intra_pool_comm )
CALL mp_bcast( beta0, ionode_id, intra_pool_comm )
CALL mp_bcast( conv_ions, ionode_id, intra_image_comm )
CALL mp_bcast( tau, ionode_id, intra_image_comm )
CALL mp_bcast( force, ionode_id, intra_image_comm )
CALL mp_bcast( tr2, ionode_id, intra_image_comm )
CALL mp_bcast( conv_ions, ionode_id, intra_image_comm )
CALL mp_bcast( alpha0, ionode_id, intra_image_comm )
CALL mp_bcast( beta0, ionode_id, intra_image_comm )
CALL mp_bcast( history, ionode_id, intra_image_comm )
!
RETURN
!
@ -346,8 +346,8 @@ SUBROUTINE find_alpha_and_beta( nat, tau, tauold, alpha0, beta0 )
!
! ... | tau(t+dt) - tau' | is minimum, where
!
! ... tau' = alpha0 * ( tau(t) - tau(t-dt) ) +
! ... beta0 * ( tau(t-dt) - tau(t-2*dt) )
! ... tau' = tau(t) + alpha0 * ( tau(t) - tau(t-dt) )
! ... + beta0 * ( tau(t-dt) -tau(t-2*dt) )
!
USE constants, ONLY : eps8
USE kinds, ONLY : DP
@ -363,7 +363,7 @@ SUBROUTINE find_alpha_and_beta( nat, tau, tauold, alpha0, beta0 )
a11 = 0.D0
a12 = 0.D0
a21 = 0.D0
a22 = 0.D0 + eps8
a22 = 0.D0
b1 = 0.D0
b2 = 0.D0
c = 0.D0
@ -395,11 +395,20 @@ SUBROUTINE find_alpha_and_beta( nat, tau, tauold, alpha0, beta0 )
!
det = a11 * a22 - a12 * a21
!
IF ( det < 0.D0 ) CALL errore( 'find_alpha_and_beta', ' det < 0', 1 )
IF ( det < 0.D0 ) THEN
!
alpha0 = 0.D0
beta0 = 0.D0
!
PRINT *, "WARNING IN find_alpha_and_beta: det < 0"
!
!CALL errore( 'find_alpha_and_beta', ' det < 0', 1 )
!
END IF
!
! ... case det > 0: a well defined minimum exists
!
IF ( det > 0.D0 ) THEN
IF ( det > eps8 ) THEN
!
alpha0 = ( b1 * a22 - b2 * a12 ) / det
beta0 = ( a11 * b2 - a21 * b1 ) / det
@ -407,10 +416,10 @@ SUBROUTINE find_alpha_and_beta( nat, tau, tauold, alpha0, beta0 )
ELSE
!
! ... case det = 0 : the two increments are linearly dependent,
! ... chose solution with beta = 0
! ... chose solution with alpha = 0 and beta = 0
! ... ( discard oldest configuration )
!
alpha0 = 1.D0
alpha0 = 0.D0
beta0 = 0.D0
!
IF ( a11 > 0.D0 ) alpha0 = b1 / a11
@ -431,7 +440,16 @@ SUBROUTINE find_alpha_and_beta( nat, tau, tauold, alpha0, beta0 )
!
END DO
!
!WRITE( stdout, * ) chi, alpha0, beta0
#if defined (__DEBUG_EXTR)
PRINT *, ""
PRINT *, "chi = ", chi, " det = ", det
PRINT *, "alpha = ", alpha0, " beta = ", beta0
PRINT *, ""
PRINT *, "PREDICTED POSITIONS:"
PRINT *, tauold(1,:,1) + alpha0 * ( tauold(1,:,1) - tauold(1,:,2) ) + &
beta0 * ( tauold(1,:,2) - tauold(1,:,3) )
PRINT *, ""
#endif
!
RETURN
!

22
PW/openfil.f90
View File

@ -23,9 +23,6 @@ SUBROUTINE openfil()
iunat, iunocc, iunwfc, iunoldwfc, iunoldwfc2, &
iunigk, nwordwfc, nwordatwfc, iunneb
USE restart_module, ONLY : readfile_new
!#ifdef __PARA
USE para, ONLY :
!#endif
!
IMPLICIT NONE
!
@ -34,15 +31,6 @@ SUBROUTINE openfil()
REAL(KIND=DP) :: edum(1,1), wdum(1,1)
!
!
! ... iunwfc contains the wavefunctions
!
iunwfc = 10
!
! ... iunoldwfc contains the old wavefunctions, used in molecular dynamics
!
iunoldwfc = 11
iunoldwfc2= 12
!
! ... nwordwfc is the record length for the direct-access file
! ... containing wavefunctions
!
@ -88,15 +76,6 @@ SUBROUTINE openfil()
!
iunocc = 14
!
! ... if extrapolation of wfc's is requested (order=2)
! ... another file is needed to store the "old" wfc's
!
IF ( order > 1 ) &
CALL diropn( iunoldwfc, TRIM( prefix )//'.oldwfc', nwordwfc, exst )
!
IF ( order > 2 ) &
CALL diropn( iunoldwfc2, TRIM( prefix )//'.oldwfc2', nwordwfc, exst )
!
! ... iunigk contains the number of PW and the indices igk
! ... Note that unit 15 is reserved for error messages
!
@ -106,4 +85,3 @@ SUBROUTINE openfil()
RETURN
!
END SUBROUTINE openfil

240
PW/para.f90
View File

@ -90,9 +90,9 @@ MODULE para
!
! ... number of processors = # of tasks
!
INTEGER :: &
MPI_COMM_POOL = 0, &! comunicator handle intra-pool
MPI_COMM_ROW = 0 ! " " inter-pool
! INTEGER :: &
! MPI_COMM_POOL = 0, &! comunicator handle intra-pool
! MPI_COMM_ROW = 0 ! " " inter-pool
!
! ... general parallel information
!
@ -120,7 +120,7 @@ SUBROUTINE reduce( dim, ps )
!
#if defined (__PARA)
!
USE mp_global, ONLY : intra_pool_comm, my_pool_id, nproc_pool
USE mp_global, ONLY : intra_pool_comm, my_pool_id, nproc_pool, my_image_id
USE mp, ONLY : mp_barrier
USE kinds, ONLY : DP
USE parallel_include
@ -155,7 +155,7 @@ SUBROUTINE reduce( dim, ps )
!
! ... syncronize processes - maybe unneeded on T3D but necessary on T3E !!!
!
CALL mp_barrier( intra_pool_comm )
CALL mp_barrier( intra_pool_comm(my_image_id) )
!
nbuf = dim / maxb
!
@ -174,8 +174,8 @@ SUBROUTINE reduce( dim, ps )
!
# else
!
CALL MPI_allreduce( ps(1+(n-1)*maxb), buff, maxb, &
MPI_REAL8, MPI_SUM, intra_pool_comm, info )
CALL MPI_ALLREDUCE( ps(1+(n-1)*maxb), buff, maxb, MPI_REAL8, &
MPI_SUM, intra_pool_comm(my_image_id), info )
!
CALL errore( 'reduce', 'error in allreduce1', info )
!
@ -196,8 +196,8 @@ SUBROUTINE reduce( dim, ps )
!
# else
!
CALL MPI_allreduce( ps(1+nbuf*maxb), buff, (dim-nbuf*maxb), &
MPI_REAL8, MPI_SUM, intra_pool_comm, info )
CALL MPI_ALLREDUCE( ps(1+nbuf*maxb), buff, (dim-nbuf*maxb), MPI_REAL8, &
MPI_SUM, intra_pool_comm(my_image_id), info )
!
CALL errore( 'reduce', 'error in allreduce2', info )
!
@ -223,7 +223,7 @@ SUBROUTINE ireduce( dim, is )
!
#if defined (__PARA)
!
USE mp_global, ONLY : intra_pool_comm, nproc_pool
USE mp_global, ONLY : intra_pool_comm, nproc_pool, my_image_id
USE mp, ONLY : mp_barrier
USE parallel_include
!
@ -239,14 +239,14 @@ SUBROUTINE ireduce( dim, is )
!
! ... syncronize processes
!
CALL mp_barrier( intra_pool_comm )
CALL mp_barrier( intra_pool_comm(my_image_id) )
!
nbuf = dim / maxi
!
DO n = 1, nbuf
!
CALL MPI_allreduce( is(1+(n-1)*maxi), buff, maxi, &
MPI_INTEGER, MPI_SUM, intra_pool_comm, info )
CALL MPI_ALLREDUCE( is(1+(n-1)*maxi), buff, maxi, MPI_INTEGER, &
MPI_SUM, intra_pool_comm(my_image_id), info )
!
CALL errore( 'ireduce', 'error in allreduce 1', info )
!
@ -258,8 +258,8 @@ SUBROUTINE ireduce( dim, is )
!
IF ( ( dim - nbuf * maxi ) > 0 ) THEN
!
CALL MPI_allreduce( is(1+nbuf*maxi), buff, (dim-nbuf*maxi), &
MPI_INTEGER, MPI_SUM, intra_pool_comm, info )
CALL MPI_ALLREDUCE( is(1+nbuf*maxi), buff, (dim-nbuf*maxi), MPI_INTEGER, &
MPI_SUM, intra_pool_comm(my_image_id), info )
!
CALL errore( 'reduce', 'error in allreduce 2', info )
!
@ -282,7 +282,8 @@ SUBROUTINE poolreduce( dim, ps )
!
#if defined (__PARA)
!
USE mp_global, ONLY : inter_pool_comm, my_pool_id, nproc_pool
USE mp_global, ONLY : inter_pool_comm, intra_image_comm, &
my_pool_id, nproc_pool, my_image_id
USE mp, ONLY : mp_barrier
USE kinds, ONLY : DP
USE parallel_include
@ -302,14 +303,14 @@ SUBROUTINE poolreduce( dim, ps )
!
! ... MPI syncronize processes
!
CALL mp_barrier()
CALL mp_barrier( intra_image_comm )
!
nbuf = dim / maxb
!
DO n = 1, nbuf
!
CALL MPI_allreduce( ps(1+(n-1)*maxb), buff, maxb, &
MPI_REAL8, MPI_SUM, inter_pool_comm, info )
CALL MPI_ALLREDUCE( ps(1+(n-1)*maxb), buff, maxb, MPI_REAL8, &
MPI_SUM, inter_pool_comm(my_image_id), info )
!
CALL errore( 'poolreduce', 'info<>0 at allreduce1', info )
!
@ -320,8 +321,8 @@ SUBROUTINE poolreduce( dim, ps )
!
IF ( ( dim - nbuf * maxb ) > 0 ) THEN
!
CALL MPI_allreduce( ps(1+nbuf*maxb), buff, (dim-nbuf*maxb), &
MPI_REAL8, MPI_SUM, inter_pool_comm, info )
CALL MPI_ALLREDUCE( ps(1+nbuf*maxb), buff, (dim-nbuf*maxb), MPI_REAL8, &
MPI_SUM, inter_pool_comm(my_image_id), info )
!
CALL errore( 'poolreduce', 'info<>0 at allreduce2', info )
!
@ -352,9 +353,8 @@ SUBROUTINE gather( f_in, f_out )
#if defined (__PARA)
!
USE pfft, ONLY : ncplane, npp, nxx
USE mp_global, ONLY : intra_pool_comm
USE io_global, ONLY : ionode_id
USE para, ONLY : me, nprocp
USE mp_global, ONLY : intra_pool_comm, nproc_pool, me_pool, &
root_pool, my_image_id
USE mp, ONLY : mp_barrier
USE kinds, ONLY : DP
USE parallel_include
@ -362,16 +362,17 @@ SUBROUTINE gather( f_in, f_out )
IMPLICIT NONE
!
REAL (KIND=DP) :: f_in(nxx), f_out(*)
INTEGER :: proc, info, displs(nprocp), recvcount(nprocp)
INTEGER :: proc, info
INTEGER :: displs(0:nproc_pool-1), recvcount(0:nproc_pool-1)
!
!
CALL start_clock( 'gather' )
!
DO proc = 1, nprocp
DO proc = 0, ( nproc_pool - 1 )
!
recvcount(proc) = ncplane * npp(proc)
recvcount(proc) = ncplane * npp(proc+1)
!
IF ( proc == 1 ) THEN
IF ( proc == 0 ) THEN
!
displs(proc) = 0
!
@ -383,10 +384,11 @@ SUBROUTINE gather( f_in, f_out )
!
END DO
!
CALL mp_barrier( intra_pool_comm )
CALL mp_barrier( intra_pool_comm(my_image_id) )
!
CALL MPI_gatherv( f_in, recvcount(me), MPI_REAL8, f_out, recvcount, &
displs, MPI_REAL8, ionode_id, intra_pool_comm, info )
CALL MPI_GATHERV( f_in, recvcount(me_pool), MPI_REAL8, f_out, &
recvcount, displs, MPI_REAL8, root_pool, &
intra_pool_comm(my_image_id), info )
!
CALL errore( 'gather', 'info<>0', info )
!
@ -409,23 +411,25 @@ SUBROUTINE cgather_sym( f_in, f_out )
#if defined (__PARA)
!
USE pfft, ONLY : ncplane, npp, nxx
USE mp_global, ONLY : intra_pool_comm
USE para, ONLY : me, mypool, nprocp
USE mp_global, ONLY : intra_pool_comm, intra_image_comm, &
nproc_pool, me_pool, my_image_id
USE mp, ONLY : mp_barrier
USE parallel_include
!
IMPLICIT NONE
!
COMPLEX(KIND=DP) :: f_in(nxx), f_out(*)
INTEGER :: proc, info, displs(nprocp), recvcount(nprocp)
INTEGER :: proc, info
INTEGER :: displs(0:nproc_pool-1), recvcount(0:nproc_pool-1)
!
!
CALL start_clock( 'cgather' )
!
DO proc = 1, nprocp
DO proc = 0, ( nproc_pool - 1 )
!
recvcount(proc) = 2 * ncplane * npp(proc)
recvcount(proc) = 2 * ncplane * npp(proc+1)
!
IF ( proc == 1 ) THEN
IF ( proc == 0 ) THEN
!
displs(proc) = 0
!
@ -437,14 +441,15 @@ SUBROUTINE cgather_sym( f_in, f_out )
!
END DO
!
CALL mp_barrier( intra_pool_comm )
CALL mp_barrier( intra_pool_comm(my_image_id) )
!
CALL MPI_allgatherv( f_in, recvcount(me), MPI_REAL8, f_out, &
recvcount, displs, MPI_REAL8, intra_pool_comm, info )
CALL MPI_ALLGATHERV( f_in, recvcount(me_pool), MPI_REAL8, &
f_out, recvcount, displs, MPI_REAL8, &
intra_pool_comm(my_image_id), info )
!
CALL errore( 'cgather_sym', 'info<>0', info )
!
CALL mp_barrier()
CALL mp_barrier( intra_image_comm )
!
CALL stop_clock( 'cgather' )
!
@ -468,9 +473,8 @@ SUBROUTINE scatter( f_in, f_out )
#if defined (__PARA)
!
USE pfft, ONLY : ncplane, npp, nxx
USE mp_global, ONLY : intra_pool_comm
USE io_global, ONLY : ionode_id
USE para, ONLY : me, nprocp
USE mp_global, ONLY : intra_pool_comm, nproc_pool, &
me_pool, root_pool, my_image_id
USE mp, ONLY : mp_barrier
USE kinds, ONLY : DP
USE parallel_include
@ -478,16 +482,17 @@ SUBROUTINE scatter( f_in, f_out )
IMPLICIT NONE
!
REAL (KIND=DP) :: f_in(*), f_out(nxx)
INTEGER :: proc, info, displs(nprocp), sendcount(nprocp)
INTEGER :: proc, info
INTEGER :: displs(0:nproc_pool-1), sendcount(0:nproc_pool-1)
!
!
CALL start_clock( 'scatter' )
!
DO proc = 1, nprocp
DO proc = 0, ( nproc_pool - 1 )
!
sendcount(proc) = ncplane * npp(proc)
sendcount(proc) = ncplane * npp(proc+1)
!
IF ( proc == 1 ) THEN
IF ( proc == 0 ) THEN
!
displs(proc) = 0
!
@ -499,14 +504,15 @@ SUBROUTINE scatter( f_in, f_out )
!
END DO
!
CALL mp_barrier( intra_pool_comm )
CALL mp_barrier( intra_pool_comm(my_image_id) )
!
CALL MPI_scatterv( f_in, sendcount, displs, MPI_REAL8, f_out, sendcount(me), &
MPI_REAL8, ionode_id, intra_pool_comm, info )
CALL MPI_SCATTERV( f_in, sendcount, displs, MPI_REAL8, &
f_out, sendcount(me_pool), MPI_REAL8, &
root_pool, intra_pool_comm(my_image_id), info )
!
CALL errore( 'scatter', 'info<>0', info )
!
IF ( sendcount(me) /= nxx ) f_out(sendcount(me)+1:nxx) = 0.D0
IF ( sendcount(me_pool) /= nxx ) f_out(sendcount(me_pool)+1:nxx) = 0.D0
!
CALL stop_clock( 'scatter' )
!
@ -529,9 +535,9 @@ SUBROUTINE poolscatter( nsize, nkstot, f_in, nks, f_out )
#if defined (__PARA)
!
USE kinds, ONLY : DP
USE mp_global, ONLY : intra_pool_comm, inter_pool_comm, my_pool_id
USE para, ONLY : me, npool, kunit
USE io_global, ONLY : ionode_id
USE mp_global, ONLY : intra_pool_comm, inter_pool_comm, &
my_pool_id, npool, me_pool, root_pool, my_image_id
USE para, ONLY : kunit
USE mp, ONLY : mp_bcast
!
IMPLICIT NONE
@ -551,7 +557,8 @@ SUBROUTINE poolscatter( nsize, nkstot, f_in, nks, f_out )
! ... copy from the first node of the first pool
! ... to the first node of all the other pools
!
IF ( me == 1 ) CALL mp_bcast( f_in, ionode_id, inter_pool_comm )
IF ( me_pool == root_pool ) &
CALL mp_bcast( f_in, root_pool, inter_pool_comm(my_image_id) )
!
! ... distribute the vector on the first node of each pool
!
@ -564,10 +571,10 @@ SUBROUTINE poolscatter( nsize, nkstot, f_in, nks, f_out )
f_out(:,1:nks) = f_in(:,(nbase+1):(nbase+nks))
!CALL DCOPY( nsize * nks, f_in(1,nbase+1), 1, f_out, 1 )
!
! ... copy from the first node of every pool
! ... copy from the first node of every pool
! ... to the other nodes of every pool
!
CALL mp_bcast( f_out, ionode_id, intra_pool_comm )
CALL mp_bcast( f_out, root_pool, intra_pool_comm(my_image_id) )
!
#endif
!
@ -611,8 +618,7 @@ SUBROUTINE fft_scatter1( f_in, nrx3, nxx_, f_aux, ncp_, npp_, sign )
#if defined (__PARA)
!
USE para_const, ONLY : maxproc
USE mp_global, ONLY : intra_pool_comm, nproc_pool
USE para, ONLY : me, mypool, nprocp
USE mp_global, ONLY : intra_pool_comm, nproc_pool, me_pool, my_image_id
USE mp, ONLY : mp_barrier
USE parallel_include
!
@ -620,9 +626,10 @@ SUBROUTINE fft_scatter1( f_in, nrx3, nxx_, f_aux, ncp_, npp_, sign )
!
INTEGER :: nrx3, nxx_, sign, ncp_(maxproc), npp_(maxproc)
REAL (KIND=DP) :: f_in(2*nxx_), f_aux(2*nxx_)
INTEGER :: dest, from, k, offset1(maxproc), sendcount(maxproc), &
sdispls(maxproc), recvcount(maxproc), rdispls(maxproc), &
proc, ierr
INTEGER :: dest, from, k, proc, ierr
INTEGER :: offset1(0:maxproc-1), sendcount(0:maxproc-1), &
sdispls(0:maxproc-1), recvcount(0:maxproc-1), &
rdispls(0:maxproc-1)
!
!
IF ( nproc_pool == 1 ) RETURN
@ -632,10 +639,10 @@ SUBROUTINE fft_scatter1( f_in, nrx3, nxx_, f_aux, ncp_, npp_, sign )
! ... sendcount(proc): amount of data processor "me" must send to processor
! ... recvcount(proc): amount of data processor "me" must receive from
!
DO proc = 1, nprocp
DO proc = 0, ( nproc_pool - 1 )
!
sendcount(proc) = 2 * npp_(proc) * ncp_(me)
recvcount(proc) = 2 * npp_(me) * ncp_(proc)
sendcount(proc) = 2 * npp_(proc+1) * ncp_(me_pool+1)
recvcount(proc) = 2 * npp_(me_pool+1) * ncp_(proc+1)
!
END DO
!
@ -643,13 +650,13 @@ SUBROUTINE fft_scatter1( f_in, nrx3, nxx_, f_aux, ncp_, npp_, sign )
! ... sdispls(proc)+1 is the beginning of data that must be sent to proc
! ... rdispls(proc)+1 is the beginning of data that must be received from proc
!
offset1(1) = 1
sdispls(1) = 0
rdispls(1) = 0
offset1(0) = 1
sdispls(0) = 0
rdispls(0) = 0
!
DO proc = 2, nprocp
DO proc = 1, ( nproc_pool - 1 )
!
offset1(proc) = offset1(proc-1) + 2 * npp_(proc-1)
offset1(proc) = offset1(proc-1) + 2 * npp_(proc)
sdispls(proc) = sdispls(proc-1) + sendcount(proc-1)
rdispls(proc) = rdispls(proc-1) + recvcount(proc-1)
!
@ -663,15 +670,15 @@ SUBROUTINE fft_scatter1( f_in, nrx3, nxx_, f_aux, ncp_, npp_, sign )
!
! ... step one: store contiguously the slices
!
DO proc = 1, nprocp
DO proc = 0, ( nproc_pool - 1 )
!
from = offset1(proc)
dest = 1 + sdispls(proc)
!
DO k = 1, ncp_(me)
DO k = 1, ncp_(me_pool+1)
!
CALL DCOPY( 2*npp_(proc), f_in(from+2*(k-1)*nrx3), &
1, f_aux(dest+2*(k-1)*npp_(proc)), 1 )
CALL DCOPY( 2*npp_(proc+1), f_in(from+2*(k-1)*nrx3), &
1, f_aux(dest+2*(k-1)*npp_(proc+1)), 1 )
!
END DO
!
@ -683,10 +690,11 @@ SUBROUTINE fft_scatter1( f_in, nrx3, nxx_, f_aux, ncp_, npp_, sign )
!
! ... step two: communication
!
CALL mp_barrier( intra_pool_comm )
CALL mp_barrier( intra_pool_comm(my_image_id) )
!
CALL MPI_alltoallv( f_aux, sendcount, sdispls, MPI_REAL8, f_in, &
recvcount, rdispls, MPI_REAL8, intra_pool_comm, ierr )
CALL MPI_ALLTOALLV( f_aux, sendcount, sdispls, MPI_REAL8, &
f_in, recvcount, rdispls, MPI_REAL8, &
intra_pool_comm(my_image_id), ierr )
!
CALL errore( 'fft_scatter', 'info<>0', ierr )
!
@ -696,10 +704,11 @@ SUBROUTINE fft_scatter1( f_in, nrx3, nxx_, f_aux, ncp_, npp_, sign )
!
! ... step two: communication
!
CALL mp_barrier( intra_pool_comm )
CALL mp_barrier( intra_pool_comm(my_image_id) )
!
CALL MPI_alltoallv( f_in, recvcount, rdispls, MPI_REAL8, f_aux, &
sendcount, sdispls, MPI_REAL8, intra_pool_comm, ierr )
CALL MPI_ALLTOALLV( f_in, recvcount, rdispls, MPI_REAL8, &
f_aux, sendcount, sdispls, MPI_REAL8, &
intra_pool_comm(my_image_id), ierr )
!
CALL errore( 'fft_scatter', 'info<>0', ierr )
!
@ -707,15 +716,15 @@ SUBROUTINE fft_scatter1( f_in, nrx3, nxx_, f_aux, ncp_, npp_, sign )
!
f_in(:) = ( 0.D0, 0.D0 )
!
DO proc = 1, nprocp
DO proc = 0,( nproc_pool - 1 )
!
from = 1 + sdispls(proc)
!
dest = offset1(proc)
!
DO k = 1, ncp_(me)
DO k = 1, ncp_(me_pool+1)
!
CALL DCOPY( 2*npp_(proc), f_aux(from+2*(k-1)*npp_(proc)), &
CALL DCOPY( 2*npp_(proc+1), f_aux(from+2*(k-1)*npp_(proc+1)), &
1, f_in (dest+2*(k-1)*nrx3), 1 )
!
END DO
@ -744,8 +753,7 @@ SUBROUTINE poolextreme( ps, iflag )
!
#if defined (__PARA)
!
USE mp_global, ONLY : inter_pool_comm
USE para, ONLY : npool
USE mp_global, ONLY : inter_pool_comm, intra_image_comm, npool, my_image_id
USE mp, ONLY : mp_barrier
USE parallel_include
!
@ -757,19 +765,19 @@ SUBROUTINE poolextreme( ps, iflag )
!
IF ( npool <= 1 ) RETURN
!
CALL mp_barrier()
CALL mp_barrier( intra_image_comm )
!
IF ( iflag > 0 ) THEN
!
CALL MPI_allreduce( ps, psr, 1, MPI_REAL8, MPI_MAX, &
inter_pool_comm, info )
CALL MPI_ALLREDUCE( ps, psr, 1, MPI_REAL8, MPI_MAX, &
inter_pool_comm(my_image_id), info )
!
CALL errore( 'poolextreme', 'info<>0 in allreduce1', info )
!
ELSE
!
CALL MPI_allreduce( ps, psr, 1, MPI_REAL8, MPI_MIN, &
inter_pool_comm, info )
CALL MPI_ALLREDUCE( ps, psr, 1, MPI_REAL8, MPI_MIN, &
inter_pool_comm(my_image_id), info )
!
CALL errore( 'poolextreme', 'info<>0 in allreduce2', info )
!
@ -792,8 +800,9 @@ SUBROUTINE poolrecover( vec, length, nkstot, nks )
!
#if defined (__PARA)
!
USE mp_global, ONLY : inter_pool_comm
USE para, ONLY : me, mypool, kunit, npool
USE mp_global, ONLY : inter_pool_comm, intra_image_comm, &
npool, me_pool, root_pool, my_pool_id, my_image_id
USE para, ONLY : kunit
USE mp, ONLY : mp_barrier
USE parallel_include
!
@ -813,12 +822,12 @@ SUBROUTINE poolrecover( vec, length, nkstot, nks )
!
rest = ( nkstot - nks1 * npool ) / kunit
!
CALL mp_barrier()
CALL mp_barrier( intra_image_comm )
!
IF ( me == 1 .AND. mypool /= 1 ) THEN
IF ( me_pool == root_pool .AND. my_pool_id > 0 ) THEN
!
CALL MPI_send( vec, (length*nks), MPI_REAL8, 0, 17, &
inter_pool_comm, info )
CALL MPI_SEND( vec, (length*nks), MPI_REAL8, 0, 17, &
inter_pool_comm(my_image_id), info )
!
CALL errore( 'poolrecover', 'info<>0 in send', info )
!
@ -840,10 +849,10 @@ SUBROUTINE poolrecover( vec, length, nkstot, nks )
!
END IF
!
IF ( me == 1 .AND. mypool == 1 ) THEN
IF ( me_pool == root_pool .AND. my_pool_id == 0 ) THEN
!
CALL MPI_recv( vec(1,nbase+1), (length*fine), MPI_REAL8, &
(i-1), 17, inter_pool_comm, status, info )
CALL MPI_RECV( vec(1,nbase+1), (length*fine), MPI_REAL8, &
(i-1), 17, inter_pool_comm(my_image_id), status, info )
!
CALL errore( 'poolrecover', 'info<>0 in recv', info )
!
@ -865,8 +874,9 @@ SUBROUTINE ipoolrecover( ivec, length, nkstot, nks )
!
#if defined (__PARA)
!
USE mp_global, ONLY : inter_pool_comm
USE para, ONLY : me, mypool, kunit, npool
USE mp_global, ONLY : inter_pool_comm, intra_image_comm, &
npool, me_pool, root_pool, my_pool_id, my_image_id
USE para, ONLY : kunit
USE mp, ONLY : mp_barrier
USE parallel_include
!
@ -886,12 +896,12 @@ SUBROUTINE ipoolrecover( ivec, length, nkstot, nks )
!
rest = ( nkstot - nks1 * npool ) / kunit
!
CALL mp_barrier()
CALL mp_barrier( intra_image_comm )
!
IF ( me == 1 .AND. mypool /= 1 ) THEN
IF ( me_pool == root_pool .AND. my_pool_id > 0 ) THEN
!
CALL MPI_send( ivec, (length*nks), MPI_INTEGER, 0, 17, &
inter_pool_comm, info )
CALL MPI_SEND( ivec, (length*nks), MPI_INTEGER, 0, 17, &
inter_pool_comm(my_image_id), info )
!
CALL errore( 'ipoolrecover', 'info<>0 in send', info )
!
@ -913,10 +923,10 @@ SUBROUTINE ipoolrecover( ivec, length, nkstot, nks )
!
END IF
!
IF ( me == 1 .AND. mypool == 1 ) THEN
IF ( me_pool == root_pool .AND. my_pool_id == 0 ) THEN
!
CALL MPI_recv( ivec(1,nbase+1), (length*fine), MPI_INTEGER, &
(i-1), 17, inter_pool_comm, status, info )
CALL MPI_RECV( ivec(1,nbase+1), (length*fine), MPI_INTEGER, &
(i-1), 17, inter_pool_comm(my_image_id), status, info )
!
CALL errore( 'ipoolrecover', 'info<>0 in recv', info )
!
@ -939,7 +949,7 @@ SUBROUTINE extreme( ps, iflag )
!
#if defined (__PARA)
!
USE mp_global, ONLY : intra_pool_comm, inter_pool_comm
USE mp_global, ONLY : intra_image_comm
USE mp, ONLY : mp_barrier
USE parallel_include
!
@ -949,17 +959,17 @@ SUBROUTINE extreme( ps, iflag )
INTEGER :: iflag, info
!
!
CALL mp_barrier()
CALL mp_barrier( intra_image_comm )
!
IF ( iflag > 0 ) THEN
!
CALL MPI_allreduce( ps, psr, 1, MPI_REAL8, MPI_MAX, &
MPI_COMM_WORLD, info )
CALL MPI_ALLREDUCE( ps, psr, 1, MPI_REAL8, MPI_MAX, &
intra_image_comm, info )
!
ELSE
!
CALL MPI_allreduce( ps, psr, 1, MPI_REAL8, MPI_MIN, &
MPI_COMM_WORLD, info )
CALL MPI_ALLREDUCE( ps, psr, 1, MPI_REAL8, MPI_MIN, &
intra_image_comm, info )
!
END IF
!

99
PW/potinit.f90
View File

@ -11,35 +11,32 @@
SUBROUTINE potinit()
!----------------------------------------------------------------------------
!
! This routine initializes the self consistent potential in the array
! vr. There are three possible cases:
! ... This routine initializes the self consistent potential in the array
! ... vr. There are three possible cases:
!
! a) In this run the code is restarting from a broken run
! b) The potential (or rho) is read from file
! c) if a and b are both false, the total charge is computed
! as a sum of atomic charges, and the corresponding potential
! is saved in vr
! ... a) In this run the code is restarting from a broken run
! ... b) The potential (or rho) is read from file
! ... c) if a and b are both false, the total charge is computed
! ... as a sum of atomic charges, and the corresponding potential
! ... is saved in vr
!
USE kinds, ONLY : DP
USE io_global, ONLY : stdout
USE brilz, ONLY : alat, omega
USE basis, ONLY : nat, startingpot
USE klist, ONLY : nelec
USE lsda_mod, ONLY : lsda, nspin
USE gvect, ONLY : ngm, gstart, nr1, nr2, nr3, nrx1, nrx2, nrx3, &
nrxx, nl, g, gg
USE gsmooth, ONLY : doublegrid
USE control_flags, ONLY : imix, lscf
USE scf, ONLY : rho, rho_core, vltot, vr, vrs
USE ener, ONLY : ehart, etxc, vtxc
USE ldaU, ONLY : niter_with_fixed_ns
USE ldaU, ONLY : lda_plus_u, Hubbard_lmax, ns, nsnew
USE io_files, ONLY : prefix, iunocc, input_drho
#ifdef __PARA
USE para, ONLY : me, mypool
USE mp, ONLY : mp_bcast
USE io_global, ONLY : ionode_id
#endif
USE kinds, ONLY : DP
USE io_global, ONLY : stdout
USE brilz, ONLY : alat, omega
USE basis, ONLY : nat, startingpot
USE klist, ONLY : nelec
USE lsda_mod, ONLY : lsda, nspin
USE gvect, ONLY: ngm, gstart, nr1, nr2, nr3, nrx1, nrx2, nrx3, &
nrxx, nl, g, gg
USE gsmooth, ONLY : doublegrid
USE control_flags, ONLY : imix, lscf
USE scf, ONLY : rho, rho_core, vltot, vr, vrs
USE ener, ONLY : ehart, etxc, vtxc
USE ldaU, ONLY : niter_with_fixed_ns
USE ldaU, ONLY : lda_plus_u, Hubbard_lmax, ns, nsnew
USE io_files, ONLY : prefix, iunocc, input_drho
USE mp, ONLY : mp_bcast
USE mp_global, ONLY : intra_image_comm, me_image, root_image
!
IMPLICIT NONE
!
@ -53,9 +50,8 @@ SUBROUTINE potinit()
! ... end of local variables
!
!
#ifdef __PARA
IF ( me == 1 .AND. mypool == 1 ) THEN
#endif
IF ( me_image == root_image ) THEN
!
IF ( imix >= 0 .AND. lscf ) THEN
CALL seqopn( 4, TRIM( prefix )//'.rho', 'UNFORMATTED', exst )
ELSE
@ -66,20 +62,21 @@ SUBROUTINE potinit()
ELSE
CLOSE( UNIT = 4, STATUS = 'DELETE' )
END IF
#ifdef __PARA
!
END IF
CALL mp_bcast( exst, ionode_id )
#endif
!
CALL mp_bcast( exst, root_image, intra_image_comm )
!
IF ( startingpot == 'file' .AND. exst ) THEN
!
!
! ... First case, the potential is read from file
! ... NB: this case applies also for a restarting run, in which case
! ... potential and rho files have been read from the restart file
!
IF ( imix >= 0 .AND. lscf ) THEN
!
!
CALL io_pot( -1, TRIM( prefix )//'.rho', rho, nspin )
!
WRITE( stdout, '(/5X,"The initial density is read from file ", A14)' ) &
TRIM( prefix ) // '.rho'
!
@ -89,38 +86,47 @@ SUBROUTINE potinit()
CALL v_of_rho( rho, rho_core, nr1, nr2, nr3, nrx1, nrx2, nrx3, &
nrxx, nl, ngm, gstart, nspin, g, gg, alat, omega, &
ehart, etxc, vtxc, charge, vr )
!
!
IF ( ABS( charge - nelec ) / charge > 1.0D-4 ) &
WRITE( stdout, '(/5X,"starting charge =",F10.5)') charge
WRITE( stdout, '(/5X,"starting charge =",F10.5)') charge
!
ELSE
!
CALL io_pot( -1, TRIM( prefix )//'.pot', vr, nspin )
!
WRITE( stdout, '(/5X,"The initial potential is read from file ", A14)' ) &
TRIM( prefix ) // '.pot'
!
END IF
!
! ... The occupations ns also need to be read in order to build up
! ... the poten
!
IF ( lda_plus_u ) THEN
!
ldim = 2 * Hubbard_lmax + 1
#ifdef __PARA
IF ( me == 1 .AND. mypool == 1 ) THEN
#endif
!
IF ( me_image == root_image ) THEN
!
CALL seqopn( iunocc, TRIM( prefix )//'.occup', 'FORMATTED', exst )
READ( UNIT = iunocc, FMT = * ) ns
CLOSE( UNIT = iunocc, STATUS = 'KEEP' )
#ifdef __PARA
ELSE
ns(:,:,:,:) = 0.d0
!
ELSE
!
ns(:,:,:,:) = 0.D0
!
END IF
!
CALL reduce( ( ldim * ldim * nspin * nat ), ns )
CALL poolreduce( ( ldim * ldim * nspin * nat ), ns )
#endif
!
CALL DCOPY( ( ldim * ldim * nspin * nat ), ns, 1, nsnew, 1 )
!
END IF
!
ELSE
!
!
! ... Second case, the potential is built from a superposition
! ... of atomic charges contained in the array rho_at and already
! ... set in readin-readva
@ -158,7 +164,7 @@ SUBROUTINE potinit()
CALL v_of_rho( rho, rho_core, nr1, nr2, nr3, nrx1, nrx2, nrx3, &
nrxx, nl, ngm, gstart, nspin, g, gg, alat, omega, &
ehart, etxc, vtxc, charge, vr )
!
!
IF ( ABS( charge - nelec ) / charge > 1.0D-4 ) &
WRITE( stdout, '(/5X,"starting charge =",F10.5)') charge
!
@ -181,6 +187,7 @@ SUBROUTINE potinit()
END IF
!
IF ( imix >= 0 ) CALL io_pot( +1, TRIM( prefix )//'.rho', rho, nspin )
!
CALL io_pot( +1, TRIM( prefix )//'.pot', vr, nspin )
!
RETURN

209
PW/print_clock_pw.f90
View File

@ -1,111 +1,144 @@
!
! Copyright (C) 2001 PWSCF group
! Copyright (C) 2001-2004 PWSCF group
! This file is distributed under the terms of the
! GNU General Public License. See the file `License'
! in the root directory of the present distribution,
! or http://www.gnu.org/copyleft/gpl.txt .
!
!
subroutine print_clock_pw
!----------------------------------------------------------------------------
SUBROUTINE print_clock_pw()
!---------------------------------------------------------------------------
!
! ... this routine prints out the clocks at the end of the run
! ... it tries to construct the calling tree of the program.
!
USE io_global, ONLY : stdout
USE control_flags, ONLY : isolve, imix
USE force_mod, ONLY : lforce, lstres
USE mp_global, ONLY : mpime, root
!
IMPLICIT NONE
!
!
IF ( mpime /= root ) &
OPEN( UNIT = stdout, FILE = '/dev/null', STATUS = 'UNKNOWN' )
!
! this routine prints out the clocks at the end of the run
! it tries to construct the calling tree of the program.
USE io_global, ONLY : stdout
USE control_flags, ONLY: isolve, imix
USE force_mod, ONLY: lforce, lstres
implicit none
WRITE( stdout, * )
call print_clock ('PWSCF')
call print_clock ('init_run')
call print_clock ('electrons')
if (lforce) call print_clock ('forces')
if (lstres) call print_clock ('stress')
!
CALL print_clock( 'PWSCF' )
CALL print_clock( 'init_run' )
CALL print_clock( 'electrons' )
!
IF ( lforce ) CALL print_clock( 'forces' )
IF ( lstres ) CALL print_clock( 'stress' )
!
WRITE( stdout, * )
call print_clock ('electrons')
call print_clock ('c_bands')
call print_clock ('sum_band')
call print_clock ('v_of_rho')
call print_clock ('newd')
!
CALL print_clock( 'electrons' )
CALL print_clock( 'c_bands' )
CALL print_clock( 'sum_band' )
CALL print_clock( 'v_of_rho' )
CALL print_clock( 'newd' )
!
#ifdef DEBUG_NEWD
WRITE( stdout,*) "nhm*(nhm+1)/2 = ", nhm*(nhm+1)/2, nhm
WRITE( stdout,*) "nbrx*(nbrx+1)/2*lqx = ", nbrx*(nbrx+1)/2*lqx, nbrx,lqx
call print_clock ('newd:fftvg')
call print_clock ('newd:qvan2')
call print_clock ('newd:int1')
call print_clock ('newd:int2')
!
CALL print_clock( 'newd:fftvg' )
CALL print_clock( 'newd:qvan2' )
CALL print_clock( 'newd:int1' )
CALL print_clock( 'newd:int2' )
#endif
if (imix.ge.0) then
call print_clock ('mix_rho')
else
call print_clock ('mix_pot')
endif
!
IF ( imix >= 0 ) THEN
CALL print_clock( 'mix_rho' )
ELSE
CALL print_clock( 'mix_pot' )
END IF
!
WRITE( stdout, * )
call print_clock ('c_bands')
call print_clock ('init_us_2')
call print_clock ('cegterg')
call print_clock ('ccgdiagg')
call print_clock ('diis')
!
CALL print_clock( 'c_bands' )
CALL print_clock( 'init_us_2' )
CALL print_clock( 'cegterg' )
CALL print_clock( 'ccgdiagg' )
CALL print_clock( 'diis' )
!
WRITE( stdout, * )
call print_clock ('sum_band')
call print_clock ('sumbec')
call print_clock ('addusdens')
!
CALL print_clock( 'sum_band' )
CALL print_clock( 'sumbec' )
!
CALL print_clock( 'addusdens' )
!
#ifdef DEBUG_ADDUSDENS
call print_clock ('addus:qvan2')
call print_clock ('addus:strf')
call print_clock ('addus:aux2')
call print_clock ('addus:aux')
CALL print_clock( 'addus:qvan2' )
CALL print_clock( 'addus:strf' )
CALL print_clock( 'addus:aux2' )
CALL print_clock( 'addus:aux' )
#endif
!
WRITE( stdout, * )
call print_clock ('wfcrot')
call print_clock ('wfcrot1')
call print_clock ('cegterg')
call print_clock ('ccdiagg')
call print_clock ('cdiisg')
if (isolve == 0) then
call print_clock ('h_psi')
call print_clock ('g_psi')
call print_clock ('overlap')
call print_clock ('cdiaghg')
call print_clock ('update')
call print_clock ('last')
!
CALL print_clock( 'wfcrot' )
CALL print_clock( 'wfcrot1' )
CALL print_clock( 'cegterg' )
CALL print_clock( 'ccdiagg' )
CALL print_clock( 'cdiisg' )
!
IF ( isolve == 0 ) THEN
!
CALL print_clock( 'h_psi' )
CALL print_clock( 'g_psi' )
CALL print_clock( 'overlap' )
CALL print_clock( 'cdiaghg' )
CALL print_clock( 'update' )
CALL print_clock( 'last' )
!
WRITE( stdout, * )
call print_clock ('h_psi')
call print_clock ('init')
call print_clock ('firstfft')
call print_clock ('secondfft')
call print_clock ('add_vuspsi')
call print_clock ('s_psi')
else
call print_clock ('h_1psi')
call print_clock ('s_1psi')
call print_clock ('cdiaghg')
!
CALL print_clock( 'h_psi' )
CALL print_clock( 'init' )
CALL print_clock( 'firstfft' )
CALL print_clock( 'secondfft' )
CALL print_clock( 'add_vuspsi' )
CALL print_clock( 's_psi' )
!
ELSE
!
CALL print_clock( 'h_1psi' )
CALL print_clock( 's_1psi' )
CALL print_clock( 'cdiaghg' )
!
WRITE( stdout, * )
call print_clock ('h_1psi')
call print_clock ('init')
call print_clock ('firstfft')
call print_clock ('secondfft')
call print_clock ('add_vuspsi')
endif
!
CALL print_clock( 'h_1psi' )
CALL print_clock( 'init' )
CALL print_clock( 'firstfft' )
CALL print_clock( 'secondfft' )
CALL print_clock( 'add_vuspsi' )
!
END IF
!
WRITE( stdout, * )
WRITE( stdout, * ) ' General routines'
call print_clock ('ccalbec')
call print_clock ('cft3')
call print_clock ('cft3s')
call print_clock ('interpolate')
call print_clock ('davcio')
WRITE( stdout, '(5X,"General routines")' )
!
CALL print_clock( 'ccalbec' )
CALL print_clock( 'cft3' )
CALL print_clock( 'cft3s' )
CALL print_clock( 'interpolate' )
CALL print_clock( 'davcio' )
!
WRITE( stdout, * )
#ifdef __PARA
WRITE( stdout, * ) ' Parallel routines'
call print_clock ('reduce')
call print_clock ('fft_scatter')
! call print_clock('poolreduce')
!
#if defined (__PARA)
WRITE( stdout, '(5X,"Parallel routines")' )
!
CALL print_clock( 'reduce' )
CALL print_clock( 'fft_scatter' )
CALL print_clock( 'poolreduce' )
#endif
return
end subroutine print_clock_pw
!
RETURN
!
END SUBROUTINE print_clock_pw

38
PW/pwscf.f90
View File

@ -1,5 +1,5 @@
!
! Copyright (C) 2001-2003 PWSCF group
! Copyright (C) 2001-2004 PWSCF group
! This file is distributed under the terms of the
! GNU General Public License. See the file `License'
! in the root directory of the present distribution,
@ -23,9 +23,8 @@ PROGRAM pwscf
USE neb_variables, ONLY : neb_deallocation
USE input_parameters, ONLY : deallocate_input_parameters
USE neb_routines, ONLY : initialize_neb, search_mep
#ifdef __PARA
USE para, ONLY : me, mypool
#endif
USE mp_global, ONLY : me_image, root_image
USE io_global, ONLY : ionode
!
IMPLICIT NONE
!
@ -40,17 +39,22 @@ PROGRAM pwscf
!
CALL startup( nd_nmbr, code, version_number )
!
WRITE( UNIT = stdout, &
FMT = '(/5X,"Ultrasoft (Vanderbilt) Pseudopotentials")')
!
WRITE( unit = stdout, FMT = 9010 ) ntypx, npk, lmaxx, nchix, ndm, nbrx, nqfm
IF ( ionode) THEN
!
WRITE( UNIT = stdout, &
FMT = '(/5X,"Ultrasoft (Vanderbilt) Pseudopotentials")')
!
WRITE( unit = stdout, FMT = 9010 ) &
ntypx, npk, lmaxx, nchix, ndm, nbrx, nqfm
!
END IF
!
CALL iosys()
!
IF ( noncolin ) &
IF ( ionode .AND. noncolin ) &
WRITE( UNIT = stdout, &
& FMT = '(/,5X,"non-colinear magnetization allowed",/)' )
IF ( gamma_only ) &
IF ( ionode .AND. gamma_only ) &
WRITE( UNIT = stdout, &
& FMT = '(/,5X,"gamma-point specific algorithms are used",/)' )
!
@ -59,17 +63,9 @@ PROGRAM pwscf
IF ( lneb ) THEN
!
! ... stdout is connected to a file ( specific for each image )
! ... via unit 17
! ... via unit 17 ( only root_image performes I/O )
!
#ifdef __PARA
IF ( me == 1 .AND. mypool == 1 ) THEN
#endif
!
stdout = 17
!
#ifdef __PARA
END IF
#endif
IF ( me_image == root_image ) stdout = 17
!
CALL initialize_neb( 'PW' )
!
@ -79,7 +75,7 @@ PROGRAM pwscf
!
! ... stdout is reconnected to standard output
!
stdout = 6
stdout = 6
!
CALL stop_pw( conv_neb )
!

10
PW/rdiaghg.f90
View File

@ -17,10 +17,8 @@ SUBROUTINE rdiaghg( n, m, h, s, ldh, e, v )
! ... Uses LAPACK routines
!
USE kinds, ONLY : DP
USE para, ONLY : me
USE mp, ONLY : mp_bcast
USE io_global, ONLY : ionode_id
USE mp_global, ONLY : intra_pool_comm
USE mp_global, ONLY : npool, me_pool, root_pool, intra_pool_comm, my_image_id
!
IMPLICIT NONE
!
@ -90,7 +88,7 @@ SUBROUTINE rdiaghg( n, m, h, s, ldh, e, v )
!
! ... only the first processor diagonalize the matrix
!
IF ( me == 1 ) THEN
IF ( me_pool == root_pool ) THEN
!
IF ( all_eigenvalues ) THEN
!
@ -129,8 +127,8 @@ SUBROUTINE rdiaghg( n, m, h, s, ldh, e, v )
!
! ... broadcast eigenvectors and eigenvalues to all other processors
!
CALL mp_bcast( e, ionode_id, intra_pool_comm )
CALL mp_bcast( v, ionode_id, intra_pool_comm )
CALL mp_bcast( e, root_pool, intra_pool_comm(my_image_id) )
CALL mp_bcast( v, root_pool, intra_pool_comm(my_image_id) )
!
! ... deallocate workspace
!

20
PW/restart.f90
View File

@ -59,8 +59,7 @@ subroutine writefile_new( what, ndw, et_g, wg_g, kunit )
use io_global, only: ionode
use mp, only: mp_sum, mp_max, mp_end
use mp_global, only: mpime, nproc, root, me_pool, my_pool_id, &
nproc_pool, intra_pool_comm, root_pool, inter_pool_comm
nproc_pool, intra_pool_comm, root_pool, inter_pool_comm, my_image_id
USE io_base, only: write_restart_header, write_restart_ions, &
write_restart_cell, write_restart_electrons, &
@ -215,13 +214,13 @@ subroutine writefile_new( what, ndw, et_g, wg_g, kunit )
ike = iks + nkl - 1
ngk_g(iks:ike) = ngk(1:nkl)
CALL mp_sum( ngk_g, intra_pool_comm )
CALL mp_sum( ngk_g, intra_pool_comm(my_image_id) )
!write(400+mpime,*) what, ngk(1:nks), nkstot
!write(400+mpime,*) what, ngk_g(1:nkstot)
IF( npool > 1 ) THEN
CALL mp_sum( ngk_g, inter_pool_comm )
CALL mp_sum( ngk_g, inter_pool_comm(my_image_id) )
!write(400+mpime,*) what, ngk_g(1:nkstot)
END IF
@ -529,7 +528,7 @@ subroutine readfile_new( what, ndr, et_g, wg_g, kunit, nsizwfc, iunitwfc, ierr )
USE pseudo_types, ONLY: pseudo_upf
use mp, only: mp_sum, mp_bcast, mp_max, mp_end
use mp_global, only: mpime, nproc, root, me_pool, my_pool_id, &
nproc_pool, intra_pool_comm, root_pool
nproc_pool, intra_pool_comm, root_pool, intra_image_comm, my_image_id
use io_global, only: ionode, ionode_id
USE io_base, only: read_restart_header, read_restart_ions, &
@ -606,7 +605,7 @@ subroutine readfile_new( what, ndr, et_g, wg_g, kunit, nsizwfc, iunitwfc, ierr )
endif
rewind ndr
end if
call mp_bcast( ierr, ionode_id )
call mp_bcast( ierr, ionode_id, intra_image_comm )
if( ierr /= 0 ) then
return
end if
@ -984,7 +983,7 @@ subroutine readfile_new( what, ndr, et_g, wg_g, kunit, nsizwfc, iunitwfc, ierr )
IF( ( ik >= iks ) .AND. ( ik <= ike ) ) THEN
IF( me_pool == root_pool ) ipmask( mpime + 1 ) = 1
END IF
CALL mp_sum( ipmask )
CALL mp_sum( ipmask, intra_image_comm )
DO i = 1, nproc
IF( ipmask(i) == 1 ) ipdest = ( i - 1 )
END DO
@ -994,10 +993,10 @@ subroutine readfile_new( what, ndr, et_g, wg_g, kunit, nsizwfc, iunitwfc, ierr )
CALL gk_sort (xk (1, ik), ngm, g, ecutwfc / tpiba2, npw, igk, g2kin)
CALL gk_l2gmap (ngm, ig_l2g(1), npw, igk, igk_l2g(1,ik-iks+1))
npw_g = MAXVAL( igk_l2g(:,ik-iks+1) )
CALL mp_max( npw_g, intra_pool_comm )
CALL mp_max( npw_g, intra_pool_comm(my_image_id) )
END IF
CALL mp_bcast( npw_g, ipdest )
CALL mp_bcast( npw_g, ipdest, intra_image_comm )
CALL read_restart_wfc(ndr, ik, nkstot, kunit, ispin_, nspin_, &
wfc_scal, evc, twf0, evc, twfm, npw_g, nbnd_, igk_l2g(:,ik-iks+1), npw )
@ -1046,6 +1045,7 @@ subroutine readfile_config( ndr, ibrav, nat, alat, at, tau, ierr )
USE parameters, only: npk, nacx, nsx
use io_files, only: prefix, tmp_dir
use io_global, only: ionode, ionode_id
USE mp_global, ONLY : intra_image_comm
use mp, only: mp_bcast
USE io_base, only: read_restart_header, read_restart_ions, &
@ -1109,7 +1109,7 @@ subroutine readfile_config( ndr, ibrav, nat, alat, at, tau, ierr )
endif
rewind ndr
end if
call mp_bcast( ierr, ionode_id )
call mp_bcast( ierr, ionode_id, intra_image_comm )
!
! if the file is not present or unreadable
! return immediately

111
PW/restart_from_file.f90
View File

@ -1,52 +1,81 @@
!
! Copyright (C) 2001,2003 PWSCF group
! Copyright (C) 2001-2004 PWSCF group
! This file is distributed under the terms of the
! GNU General Public License. See the file `License'
! in the root directory of the present distribution,
! or http://www.gnu.org/copyleft/gpl.txt .
!
!-----------------------------------------------------------------------
subroutine restart_from_file
!-----------------------------------------------------------------------
USE io_global, ONLY : stdout
USE io_files, ONLY : iunres
USE control_flags, ONLY: restart
implicit none
character :: where * 20 ! parameter indicating from where to restart
logical :: exst
!----------------------------------------------------------------------------
SUBROUTINE restart_from_file
!----------------------------------------------------------------------------
!
! check if restart file is present
USE io_global, ONLY : stdout
USE io_files, ONLY : iunres, tmp_dir
USE control_flags, ONLY : restart
USE parser, ONLY : delete_if_present
USE mp_global, ONLY : mpime
!
iunres = 1
if (.not.restart) then
! WRITE( stdout, '(/5x,"RECOVER from restart file has been switched off on input")')
call seqopn (iunres, 'restart', 'unformatted', exst)
! if (exst) WRITE( stdout,'(/5x,"Existing restart file has been removed")')
close (unit = iunres, status = 'delete')
return
endif
call seqopn (iunres, 'restart', 'unformatted', restart)
if (.not.restart) then
WRITE( stdout, '(/5x,"RECOVER from restart file failed: file not found")')
close (unit = iunres, status = 'delete')
return
endif
IMPLICIT NONE
!
WRITE( stdout, '(/5x,"read information from restart file")')
read (iunres, err = 10, end = 10) where
WRITE( stdout, '(5x,"Restarting in ",a)') where
if (where.ne.'ELECTRONS'.and.where.ne.'IONS') then
WRITE( stdout,*) where, '......?'
call errore ('readin', ' wrong recover file ', 1)
endif
CHARACTER(LEN=20) :: where_restart
! parameter indicating from where to restart
LOGICAL :: exst
!
! close the file for later use
!
close (unit = iunres, status = 'keep')
return
10 call errore ('readin', 'problems in reading recover file', 1)
end subroutine restart_from_file
! ... check if restart file is present
!
IF ( mpime == 0 ) THEN
!
iunres = 1
!
IF ( .NOT. restart ) THEN
!
!WRITE( UNIT = stdout, &
! & FMT = '(/5X,"RECOVER from restart file has been", &
! & " switched off on input")' )
!
CALL delete_if_present( TRIM( tmp_dir ) // 'restart' )
!
RETURN
!
END IF
!
CALL seqopn( iunres, 'restart', 'UNFORMATTED', restart )
!
IF ( .NOT. restart ) THEN
!
WRITE( UNIT = stdout, &
& FMT = '(/5X,"RECOVER from restart file failed:", &
& " file not found")')
!
CLOSE( UNIT = iunres, STATUS = 'DELETE' )
!
RETURN
!
END IF
!
WRITE( UNIT = stdout, FMT = '(/5X,"read information from restart file")' )
!
READ( iunres, ERR = 10, END = 10 ) where_restart
!
WRITE( UNIT = stdout, FMT = '(5X,"Restarting in ",A)' ) where_restart
!
IF ( where_restart /= 'ELECTRONS' .AND. where_restart /= 'IONS' ) THEN
!
WRITE( UNIT = stdout, FMT = * ) where_restart , '......?'
!
CALL errore( 'restart_from_file', ' wrong recover file ', 1 )
!
END IF
!
! ... close the file for later use
!
CLOSE( UNIT = iunres, STATUS = 'KEEP' )
!
END IF
!
RETURN
!
10 CALL errore( 'restart_from_file', 'problems in reading recover file', 1 )
!
END SUBROUTINE restart_from_file

65
PW/setup.f90
View File

@ -38,38 +38,36 @@ SUBROUTINE setup()
! invsym if true the system has inversion symmetry
! + LDA+U-related quantities.
!
!
USE kinds, ONLY : DP
USE parameters, ONLY : npsx, nchix, npk
USE io_global, ONLY : stdout
USE constants, ONLY : pi
USE brilz, ONLY : at, bg, alat, tpiba, tpiba2, ibrav, symm_type
USE basis, ONLY : nat, tau, ntyp, ityp, startingwfc, startingpot, &
natomwfc
USE gvect, ONLY : gcutm, ecutwfc, dual, nr1, nr2, nr3
USE gsmooth, ONLY : doublegrid, gcutms
USE klist, ONLY : xk, wk, xqq, nks, nelec, degauss, lgauss, lxkcry, &
nkstot
USE lsda_mod, ONLY : lsda, nspin, current_spin, isk
USE ktetra, ONLY : nk1, nk2, nk3, k1, k2, k3, tetra, ntetra, ltetra
USE symme, ONLY : s, irt, ftau, nsym, invsym
USE atom, ONLY : r, oc, nchi, lchi, mesh, msh
USE pseud, ONLY : zv, zp, nlc, nnl, bhstype, alps, aps, lmax
USE wvfct, ONLY : nbnd, nbndx
USE control_flags, ONLY : tr2, ethr, alpha0, beta0, iswitch, lscf, lmd, &
lphonon, david, isolve, imix, niter, noinv, &
restart, nosym, modenum
USE relax, ONLY : dtau_ref, starting_diag_threshold
USE cellmd, ONLY : calc
USE us, ONLY : tvanp, okvan, newpseudo
USE ldaU, ONLY : d1, d2, d3, lda_plus_u, Hubbard_U, Hubbard_l, &
Hubbard_alpha, Hubbard_lmax
USE bp, ONLY : gdir, lberry, nppstr
USE fixed_occ, ONLY : f_inp, tfixed_occ
USE char, ONLY : sname, psd
#ifdef __PARA
USE para
#endif
USE kinds, ONLY : DP
USE parameters, ONLY : npsx, nchix, npk
USE io_global, ONLY : stdout
USE constants, ONLY : pi
USE brilz, ONLY : at, bg, alat, tpiba, tpiba2, ibrav, symm_type
USE basis, ONLY : nat, tau, ntyp, ityp, startingwfc, startingpot, &
natomwfc
USE gvect, ONLY : gcutm, ecutwfc, dual, nr1, nr2, nr3
USE gsmooth, ONLY : doublegrid, gcutms
USE klist, ONLY : xk, wk, xqq, nks, nelec, degauss, lgauss, lxkcry, &
nkstot
USE lsda_mod, ONLY : lsda, nspin, current_spin, isk
USE ktetra, ONLY : nk1, nk2, nk3, k1, k2, k3, tetra, ntetra, ltetra
USE symme, ONLY : s, irt, ftau, nsym, invsym
USE atom, ONLY : r, oc, nchi, lchi, mesh, msh
USE pseud, ONLY : zv, zp, nlc, nnl, bhstype, alps, aps, lmax
USE wvfct, ONLY : nbnd, nbndx
USE control_flags, ONLY : tr2, ethr, alpha0, beta0, iswitch, lscf, lmd, &
lneb, lphonon, david, isolve, imix, niter, noinv, &
restart, nosym, modenum
USE relax, ONLY : dtau_ref, starting_diag_threshold
USE cellmd, ONLY : calc
USE us, ONLY : tvanp, okvan, newpseudo
USE ldaU, ONLY : d1, d2, d3, lda_plus_u, Hubbard_U, Hubbard_l, &
Hubbard_alpha, Hubbard_lmax
USE bp, ONLY : gdir, lberry, nppstr
USE fixed_occ, ONLY : f_inp, tfixed_occ
USE char, ONLY : sname, psd
USE para, ONLY : kunit
USE mp_global, ONLY : nimage
!
IMPLICIT NONE
!
@ -107,6 +105,9 @@ SUBROUTINE setup()
! ... end of local variables
!
!
IF ( nimage > 1 .AND. .NOT. lneb ) &
CALL errore( 'setup', 'images parallelization not permitted', 1 )
!
DO nt = 1, ntyp
DO ir = 1, mesh(nt)
IF ( r(ir,nt) > rcut) THEN

125
PW/startup.f90
View File

@ -44,24 +44,27 @@ SUBROUTINE startup( nd_nmbr, code, version )
! ... The following two modules hold global information about processors
! ... number, IDs and communicators
!
USE io_global, ONLY : stdout, io_global_start, ionode_id
USE mp_global, ONLY : mp_global_start, nproc
USE io_global, ONLY : stdout, io_global_start, ionode, ionode_id
USE mp_global, ONLY : nproc, nimage, mpime, me_image, root, root_image
USE mp_global, ONLY : mp_global_start
USE mp, ONLY : mp_start, mp_env, mp_barrier, mp_bcast
USE para_const, ONLY : maxproc
USE para, ONLY : me, mypool, npool, nprocp
USE para, ONLY : me, npool, nprocp
!
IMPLICIT NONE
!
CHARACTER :: nd_nmbr*3, code*9, version*6
INTEGER :: gid
CHARACTER :: np*80, cdate*9, ctime*9
EXTERNAL date_and_tim
INTEGER :: ierr, ilen, iargc, nargs, iiarg
CHARACTER (LEN=3) :: nd_nmbr
CHARACTER (LEN=6) :: version
CHARACTER (LEN=9) :: code, cdate, ctime
CHARACTER (LEN=80) :: np
INTEGER :: gid
EXTERNAL date_and_tim
INTEGER :: ierr = 0, ilen, iargc, nargs, iiarg
!
!
#if defined (__PARA)
!
! ... prallel case
! ... prallel case setup : MPI environment is initialized
!
# if defined (__T3E)
!
@ -83,12 +86,7 @@ SUBROUTINE startup( nd_nmbr, code, version )
!
CALL mp_global_start( 0, me, gid, nproc )
!
! ... This is added for compatibility with PVM notations
! ... parent process (source) will have me=1 - child process me=2,...,NPROC
!
me = me + 1
!
IF ( me == 1 ) THEN
IF ( ionode ) THEN
!
! ... How many pools ?
!
@ -111,101 +109,120 @@ SUBROUTINE startup( nd_nmbr, code, version )
npool = MAX( npool, 1 )
npool = MIN( npool, nproc )
!
! ... set number of processes per pool ( must be equal for all pools )
! ... How many parallel images ?
!
nprocp = nproc / npool
nargs = iargc()
!
IF ( nproc /= ( nprocp * npool ) ) &
CALL errore( 'startup', 'nproc /= nprocp*npool', 1 )
DO iiarg = 1, ( nargs - 1 )
!
CALL getarg( iiarg, np )
!
IF ( TRIM( np ) == '-nimage' .OR. TRIM( np ) == '-nimages' ) THEN
!
CALL getarg( ( iiarg + 1 ), np )
READ( np, * ) nimage
!
END IF
!
END DO
!
nimage = MAX( nimage, 1 )
nimage = MIN( nimage, nproc )
!
END IF
!
CALL mp_barrier( gid )
!
! ... transmit nprocp and npool
! ... transmit npool and nimage
!
CALL mp_bcast( nprocp, ionode_id, gid )
CALL mp_bcast( npool, ionode_id, gid )
!
! ... set the processor label for files
CALL mp_bcast( npool, ionode_id, gid )
CALL mp_bcast( nimage, ionode_id, gid )
!
IF ( nproc > maxproc ) &
CALL errore( 'startup', ' too many processors', nproc )
!
!
! ... all pools are initialized here
!
CALL init_pool()
!
! ... set the processor label for files
!
nd_nmbr = ' '
!
IF ( nproc < 10 ) THEN
!
WRITE( nd_nmbr(1:1) , '(I1)' ) me
WRITE( nd_nmbr(1:1) , '(I1)' ) ( me_image + 1 )
!
ELSE IF ( nproc < 100 ) THEN
!
IF ( me < 10 ) THEN
nd_nmbr = '0'
WRITE( nd_nmbr(2:2) , '(I1)' ) me
WRITE( nd_nmbr(2:2) , '(I1)' ) ( me_image + 1 )
ELSE
WRITE( nd_nmbr(1:2) , '(I2)' ) me
WRITE( nd_nmbr(1:2) , '(I2)' ) ( me_image + 1 )
END IF
!
ELSE
!
IF ( me < 10 ) THEN
nd_nmbr = '00'
WRITE( nd_nmbr(3:3) , '(I1)' ) me
WRITE( nd_nmbr(3:3) , '(I1)' ) ( me_image + 1 )
ELSE IF ( me < 100 ) THEN
nd_nmbr = '0'
WRITE( nd_nmbr(2:3) , '(I2)' ) me
WRITE( nd_nmbr(2:3) , '(I2)' ) ( me_image + 1 )
ELSE
WRITE( nd_nmbr, '(I3)' ) me
WRITE( nd_nmbr, '(I3)' ) ( me_image + 1 )
END IF
!
END IF
END IF
!
! ... pools are initialized here
!
CALL init_pool()
!
! ... stdout is printed only by the first cpu ( me == 1, mypool == 1 )
! ... stdout is printed only by the root_image (set in init_pool())
!
# if defined (DEBUG)
!
IF ( me /= 1 .OR. mypool /= 1 ) &
IF ( me_image /= root_image ) &
OPEN( UNIT = stdout, FILE = './out_'//nd_nmbr, STATUS = 'UNKNOWN' )
!
# else
!
IF ( me /= 1 .OR. mypool /= 1 ) &
IF ( me_image /= root_image ) &
OPEN( UNIT = stdout, FILE = '/dev/null', STATUS = 'UNKNOWN' )
!
# endif
!
! ... information printout
!
CALL date_and_tim( cdate, ctime )
!
WRITE( stdout, 9000 ) code, version, cdate, ctime
WRITE( stdout, '(/5X,"Parallel version (MPI)")' )
WRITE( stdout, '(5X,"Number of processors in use: ",I4)' ) nproc
IF ( npool /= 1 ) &
WRITE( stdout, '(5X,"K-points division: npool = ",i4)' ) npool
IF ( nprocp /= 1 ) &
WRITE( stdout, '(5X,"R & G space division: nprocp = ",i4/)' ) nprocp
!
IF ( mpime == root ) THEN
!
CALL date_and_tim( cdate, ctime )
!
WRITE( stdout, '(/5X,"Program ",A9," v.",A6," starts ...",&
&/5X,"Today is ",A9," at ",A9)' ) &
code, version, cdate, ctime
WRITE( stdout, '(/5X,"Parallel version (MPI)")' )
WRITE( stdout, '(5X,"Number of processors in use: ",I4)' ) nproc
IF ( nimage > 1 ) &
WRITE( stdout, '(5X,"NEB images division: nimage = ",i4)' ) nimage
IF ( npool > 1 ) &
WRITE( stdout, '(5X,"K-points division: npool = ",i4)' ) npool
IF ( nprocp > 1 ) &
WRITE( stdout, '(5X,"R & G space division: nprocp = ",i4/)' ) nprocp
!
END IF
!
#else
!
! ... serial case : only information printout
! ... serial case setup : only information printout
!
nd_nmbr = ' '
!
CALL date_and_tim( cdate, ctime )
!
WRITE( stdout, 9000 ) code, version, cdate, ctime
WRITE( stdout, '(/5X,"Program ",A9," v.",A6," starts ...",&
&/5X,"Today is ",A9," at ",A9)' ) code, version, cdate, ctime
!
#endif
!
RETURN
!
9000 FORMAT( /5X,'Program ',A9,' v.',A6,' starts ...',/5X, &
& 'Today is ',A9,' at ',A9)
!
END SUBROUTINE startup

37
PW/stop_pw.f90
View File

@ -1,5 +1,5 @@
!
! Copyright (C) 2001-2003 PWSCF group
! Copyright (C) 2001-2004 PWSCF group
! This file is distributed under the terms of the
! GNU General Public License. See the file `License'
! in the root directory of the present distribution,
@ -14,16 +14,13 @@ SUBROUTINE stop_pw( flag )
! ... Called at the end of the run with flag = .TRUE. (removes 'restart')
! ... or during execution with flag = .FALSE. (does not remove 'restart')
!
USE io_global, ONLY : stdout
USE control_flags, ONLY : order, lneb
USE io_files, ONLY : prefix, &
iunwfc, iunoldwfc, iunoldwfc2, iunigk, iunres
USE io_global, ONLY : stdout, ionode
USE control_flags, ONLY : lneb
USE io_files, ONLY : prefix, iunwfc, iunigk, iunres
USE input_parameters, ONLY : deallocate_input_parameters
USE io_routines, ONLY : write_restart
USE neb_variables, ONLY : neb_deallocation
#ifdef __PARA
USE mp, ONLY : mp_barrier, mp_end
#endif
!
IMPLICIT NONE
!
@ -39,24 +36,24 @@ SUBROUTINE stop_pw( flag )
! ... the execution - close and save the file
!
CLOSE( UNIT = iunwfc, STATUS = 'KEEP' )
!
IF ( order > 1 ) &
CLOSE( UNIT = iunoldwfc, STATUS = 'KEEP' )
!
IF ( order > 2 ) &
CLOSE( UNIT = iunoldwfc2, STATUS = 'KEEP' )
!
IF ( flag ) THEN
IF ( flag .AND. ionode ) THEN
!
! ... all other files must be reopened and removed
!
CALL seqopn( iunres, 'restart', 'UNFORMATTED', exst )
CLOSE( UNIT = iunres, STATUS = 'DELETE' )
!
CALL seqopn( 4, TRIM( prefix )//'.bfgs', 'UNFORMATTED', exst )
CALL seqopn( 4, TRIM( prefix ) // '.bfgs', 'UNFORMATTED', exst )
CLOSE( UNIT = 4, STATUS = 'DELETE' )
!
CALL seqopn( 4, TRIM( prefix )//'.md', 'FORMATTED', exst )
CALL seqopn( 4, TRIM( prefix ) // '.md', 'FORMATTED', exst )
CLOSE( UNIT = 4, STATUS = 'DELETE' )
!
CALL seqopn( 4, TRIM( prefix ) // '.para', 'FORMATTED', exst )
CLOSE( UNIT = 4, STATUS = 'DELETE' )
!
CALL seqopn( 4, TRIM( prefix ) // '.BLOCK', 'FORMATTED', exst )
CLOSE( UNIT = 4, STATUS = 'DELETE' )
!
END IF
@ -64,18 +61,18 @@ SUBROUTINE stop_pw( flag )
! ... iunigk is kept open during the execution - close and remove
!
CLOSE( UNIT = iunigk, STATUS = 'DELETE' )
CALL print_clock_pw
!
CALL print_clock_pw()
!
! ... NEB specific
!
IF ( lneb ) CALL write_restart()
!
CALL show_memory ()
CALL show_memory()
!
#ifdef __PARA
CALL mp_barrier()
!
CALL mp_end()
#endif
!
#ifdef __T3E
!

20
PW/sum_band.f90
View File

@ -35,12 +35,8 @@ SUBROUTINE sum_band()
USE us, ONLY : okvan, tvanp, becsum, nh, nkb, vkb
USE wavefunctions_module, ONLY : evc, psic
USE wvfct, ONLY : nbnd, npwx, npw, igk, wg, et
#if defined (__PARA)
USE para, ONLY : me, mypool
USE mp_global, ONLY : intra_pool_comm, inter_pool_comm
USE io_global, ONLY : ionode_id
USE mp_global, ONLY : intra_image_comm, me_image, root_image
USE mp, ONLY : mp_bcast
#endif
!
IMPLICIT NONE
!
@ -73,24 +69,18 @@ SUBROUTINE sum_band()
!
ELSE IF ( ltetra ) THEN
!
#if defined (__PARA)
CALL poolrecover( et, nbnd, nkstot, nks )
!
IF ( me == 1 .AND. mypool == 1 ) THEN
IF ( me_image == root_image ) THEN
!
#endif
CALL tweights( nkstot, nspin, nbnd, nelec, ntetra, tetra, et, ef, wg )
#if defined (__PARA)
!
END IF
!
CALL poolscatter( nbnd, nkstot, wg, nks, wg )
!
IF ( me == 1 ) CALL mp_bcast( ef, ionode_id, inter_pool_comm )
CALL mp_bcast( ef, root_image, intra_image_comm )
!
CALL mp_bcast( ef, ionode_id, intra_pool_comm )
!
#endif
ELSE IF ( lgauss ) THEN
!
CALL gweights( nks, wk, nbnd, nelec, degauss, ngauss, et, ef, demet, wg )
@ -146,13 +136,9 @@ SUBROUTINE sum_band()
!
IF ( okvan ) CALL addusdens()
!
#if defined (__PARA)
!
CALL poolreduce( 1, eband )
CALL poolreduce( 1, demet )
!
#endif
!
! ... symmetrization of the charge density (and local magnetization)
!
#if defined (__PARA)

419
PW/update_pot.f90
View File

@ -13,35 +13,43 @@ SUBROUTINE update_pot()
!
! ... update potential, use the integer variable order to decide the way
!
! order = 0 copy the old potential (nothing is done)
! ... order = 0 copy the old potential (nothing is done)
!
! order = 1 subtract old atomic charge density and sum the new
! if dynamics is done the routine extrapolates also
! the difference between the the scf charge and the
! atomic one,
! ... order = 1 subtract old atomic charge density and sum the new
! ... if dynamics is done the routine extrapolates also
! ... the difference between the the scf charge and the
! ... atomic one,
!
! order = 2 extrapolate the wavefunctions:
! |psi(t+dt)> = 2*|psi(t)> - |psi(t-dt)>
! ... order = 2 extrapolate the wavefunctions:
! ... |psi(t+dt)> = 2*|psi(t)> - |psi(t-dt)>
!
! order = 3 extrapolate the wavefunctions with the second-order
! formula:
! |psi(t+dt)> = |psi(t) +
! + alpha0*(|psi(t)> - |psi(t-dt)>
! + beta0* (|psi(t-dt)> - |psi(t-2*dt)>
! ... order = 3 extrapolate the wavefunctions with the second-order
! ... formula:
! ... |psi(t+dt)> = |psi(t) +
! ... + alpha0*(|psi(t)> - |psi(t-dt)>
! ... + beta0* (|psi(t-dt)> - |psi(t-2*dt)>
!
! where alpha0 and beta0 are calculated in "dynamics" so
! that |tau'-tau(t+dt)| is minimum; tau' and tau(t+dt)
! are respectively the atomic positions at time t+dt
! and the extrapolated one:
! tau(t+dt) = tau(t) +
! + alpha0*( tau(t) - tau(t-dt) )
! + beta0*( tau(t-dt) -tau(t-2*dt) )
! ... where alpha0 and beta0 are calculated in "move_ions" so
! ... that |tau'-tau(t+dt)| is minimum; tau' and tau(t+dt)
! ... are respectively the atomic positions at time t+dt
! ... and the extrapolated one:
!
! ... tau(t+dt) = tau(t) +
! ... + alpha0*( tau(t) - tau(t-dt) )
! ... + beta0*( tau(t-dt) -tau(t-2*dt) )
!
!
USE control_flags, ONLY : order
USE control_flags, ONLY : order, history
USE io_files, ONLY : prefix, tmp_dir
! USE mp_global, ONLY : mpime
!
IMPLICIT NONE
!
! ... local variables
!
INTEGER :: rho_order, wfc_order
LOGICAL :: exists
!
!
CALL start_clock( 'update_pot' )
!
@ -53,9 +61,59 @@ SUBROUTINE update_pot()
!
END IF
!
CALL extrapolate_charge()
! ... determines the maximum effective order of the extrapolation on the
! ... basis of the files that are really available
!
IF ( order >= 2 ) CALL extrapolate_wfcs()
rho_order = MIN( 1, history )
!
INQUIRE( FILE = TRIM( tmp_dir ) // &
& TRIM( prefix ) // '.oldrho', EXIST = exists )
!
IF ( exists ) THEN
!
rho_order = MIN( 2, history )
!
INQUIRE( FILE = TRIM( tmp_dir ) // &
& TRIM( prefix ) // '.oldrho2', EXIST = exists )
!
IF ( exists ) THEN
!
rho_order = MIN( 3, history )
!
END IF
!
END IF
!
wfc_order = MIN( 1, history, order )
!
INQUIRE( FILE = TRIM( tmp_dir ) // &
& TRIM( prefix ) // '.oldwfc', EXIST = exists )
!
IF ( exists ) THEN
!
wfc_order = MIN( 2, history, order )
!
INQUIRE( FILE = TRIM( tmp_dir ) // &
& TRIM( prefix ) // '.oldwfc2', EXIST = exists )
!
IF ( exists ) THEN
!
wfc_order = MIN( 3, history, order )
!
END IF
!
END IF
!
CALL extrapolate_charge( rho_order )
!
! PRINT *, "HISTORY = ", HISTORY
! PRINT *, "ORDER = ", ORDER
! PRINT *, "RHO_ORDER = ", RHO_ORDER
! PRINT *, "WFC_ORDER = ", WFC_ORDER
!
IF ( order >= 2 ) CALL extrapolate_wfcs( wfc_order )
!
! PRINT *, "mpime = ", mpime, "EXTRAPOLATION COMPLETED"
!
CALL stop_clock( 'update_pot' )
!
@ -65,25 +123,29 @@ END SUBROUTINE update_pot
!
!
!----------------------------------------------------------------------------
SUBROUTINE extrapolate_charge()
SUBROUTINE extrapolate_charge( rho_order )
!----------------------------------------------------------------------------
!
USE io_global, ONLY : stdout
USE kinds, ONLY : DP
USE brilz, ONLY : omega, bg, alat
USE basis, ONLY : nat, tau, ntyp, ityp
USE gvect, ONLY : nrxx, ngm, g, gg, gstart, nr1, nr2, nr3, nl, &
eigts1, eigts2, eigts3, nrx1, nrx2, nrx3
USE lsda_mod, ONLY : lsda, nspin
USE scf, ONLY : rho, rho_core, vr
USE control_flags, ONLY : lbfgs, istep, alpha0, beta0, imix
USE ener, ONLY : ehart, etxc, vtxc
USE cellmd, ONLY : lmovecell, omega_old
USE vlocal, ONLY : strf
USE io_files, ONLY : prefix
USE io_global, ONLY : stdout
USE kinds, ONLY : DP
USE brilz, ONLY : omega, bg, alat
USE basis, ONLY : nat, tau, ntyp, ityp
USE gvect, ONLY : nrxx, ngm, g, gg, gstart, nr1, nr2, nr3, nl, &
eigts1, eigts2, eigts3, nrx1, nrx2, nrx3
USE lsda_mod, ONLY : lsda, nspin
USE scf, ONLY : rho, rho_core, vr
USE control_flags, ONLY : alpha0, beta0, imix
USE ener, ONLY : ehart, etxc, vtxc
USE cellmd, ONLY : lmovecell, omega_old
USE vlocal, ONLY : strf
USE io_files, ONLY : prefix
!
IMPLICIT NONE
!
INTEGER, INTENT(IN) :: rho_order
!
! ... local variables
!
REAL(KIND=DP), ALLOCATABLE :: work(:), work1(:)
! work is the difference between charge density and atomic charge
! at time t
@ -91,9 +153,10 @@ SUBROUTINE extrapolate_charge()
REAL(KIND=DP) :: charge
!
!
IF ( istep == 0 ) RETURN
IF ( rho_order == 0 ) RETURN
!
ALLOCATE( work(nrxx) )
!
work(:) = 0.D0
!
! ... if order = 1 update the potential subtracting to the charge density
@ -116,39 +179,53 @@ SUBROUTINE extrapolate_charge()
!
IF ( lmovecell ) rho(:,1) = rho(:,1) * omega_old
!
! ... if dynamics extrapolate the difference between the atomic charge a
! ... extrapolate the difference between the atomic charge a
! ... the self-consistent one
!
IF ( .NOT. lbfgs ) THEN
!
IF ( istep == 1 ) THEN
!
CALL io_pot( + 1, TRIM( prefix )//'.oldrho', rho, 1 )
!
ELSE
!
ALLOCATE( work1(nrxx) )
!
work1(:) = 0.D0
!
CALL io_pot( - 1, TRIM( prefix )//'.oldrho', work, 1 )
CALL io_pot( + 1, TRIM( prefix )//'.oldrho', rho, 1 )
!
IF ( istep == 2 ) &
CALL io_pot( + 1, TRIM( prefix )//'.oldrho2', work, 1 )
!
CALL io_pot( - 1, TRIM( prefix )//'.oldrho2', work1, 1 )
CALL io_pot( + 1, TRIM( prefix )//'.oldrho2', work, 1 )
!
! ... alpha0 and beta0 have been calculated in dynamics
! ... or in vcsmd subs.
!
rho(:,1) = rho(:,1) + alpha0 * ( rho(:,1) - work(:) ) + &
beta0 * ( work(:) - work1(:) )
!
DEALLOCATE( work1 )
!
END IF
IF ( rho_order == 1 ) THEN
!
CALL io_pot( + 1, TRIM( prefix )//'.oldrho', rho, 1 )
!
ELSE IF ( rho_order == 2 ) THEN
!
! ... oldrho -> work
!
CALL io_pot( - 1, TRIM( prefix )//'.oldrho', work, 1 )
!
! ... rho -> oldrho
! ... work -> oldrho2
!
CALL io_pot( + 1, TRIM( prefix )//'.oldrho', rho, 1 )
CALL io_pot( + 1, TRIM( prefix )//'.oldrho2', work, 1 )
!
! ... alpha0 has been calculated in move_ions
!
rho(:,1) = rho(:,1) + alpha0 * ( rho(:,1) - work(:) )
!
ELSE IF ( rho_order == 3 ) THEN
!
ALLOCATE( work1(nrxx) )
!
work1(:) = 0.D0
!
! ... oldrho2 -> work1
! ... oldrho -> work
!
CALL io_pot( - 1, TRIM( prefix )//'.oldrho2', work1, 1 )
CALL io_pot( - 1, TRIM( prefix )//'.oldrho', work, 1 )
!
! ... rho -> oldrho
! ... work -> oldrho2
!
CALL io_pot( + 1, TRIM( prefix )//'.oldrho', rho, 1 )
CALL io_pot( + 1, TRIM( prefix )//'.oldrho2', work, 1 )
!
! ... alpha0 and beta0 have been calculated in move_ions
!
rho(:,1) = rho(:,1) + alpha0 * ( rho(:,1) - work(:) ) + &
beta0 * ( work(:) - work1(:) )
!
DEALLOCATE( work1 )
!
END IF
!
@ -164,7 +241,9 @@ SUBROUTINE extrapolate_charge()
! ... add atomic charges in the new positions
!
CALL atomic_rho( work, 1 )
!
rho(:,1) = rho(:,1) + work(:)
!
CALL set_rhoc()
!
! ... reset up and down charge densities in the LSDA case
@ -189,7 +268,7 @@ END SUBROUTINE extrapolate_charge
!
!
!-----------------------------------------------------------------------
SUBROUTINE extrapolate_wfcs()
SUBROUTINE extrapolate_wfcs( wfc_order )
!-----------------------------------------------------------------------
!
! ... This routine extrapolate the wfc's after a "parallel alignment"
@ -202,15 +281,16 @@ SUBROUTINE extrapolate_wfcs()
USE io_global, ONLY : stdout
USE kinds, ONLY : DP
USE klist, ONLY : nks
USE control_flags, ONLY : isolve, istep, order, alpha0, beta0
USE basis, ONLY : startingwfc
USE control_flags, ONLY : isolve, alpha0, beta0, order
USE wvfct, ONLY : nbnd, npw, npwx, igk
USE io_files, ONLY : nwordwfc, iunigk, iunwfc, iunoldwfc, &
iunoldwfc2
iunoldwfc2, prefix
USE wavefunctions_module, ONLY : evc
!
IMPLICIT NONE
!
INTEGER, INTENT(IN) :: wfc_order
!
! ... local variables
!
INTEGER :: j, i, ik, zero_ew, lwork, info
@ -219,8 +299,8 @@ SUBROUTINE extrapolate_wfcs()
! number of zero 'eigenvalues' of the s_m matrix
! used by singular value decomposition (ZGESVD)
! flag returned by ZGESVD
COMPLEX(KIND=DP), ALLOCATABLE :: s_m(:,:), sp_m(:,:), u_m(:,:), w_m(:,:), &
work(:)
COMPLEX(KIND=DP), ALLOCATABLE :: s_m(:,:), sp_m(:,:), &
u_m(:,:), w_m(:,:), work(:)
! the overlap matrix s (eq. 3.24)
! its dagger
! left unitary matrix in the SVD of sp_m
@ -231,33 +311,41 @@ SUBROUTINE extrapolate_wfcs()
REAL(KIND=DP), ALLOCATABLE :: ew(:), rwork(:)
! the eigenvalues of s_m
! workspace for ZGESVD
! ... istep = 0 when extrapolate_wfcs() is called in neb
LOGICAL :: exst
!
IF ( istep == 0 ) THEN
!
IF ( wfc_order == 0 ) THEN
!
RETURN
!
ELSE IF ( istep == 1 ) THEN
ELSE IF ( wfc_order == 1 ) THEN
!
CALL diropn( iunoldwfc, TRIM( prefix ) // '.oldwfc', nwordwfc, exst )
!
DO ik = 1, nks
!
CALL davcio( evc, nwordwfc, iunwfc, ik, - 1 )
CALL davcio( evc, nwordwfc, iunoldwfc, ik, 1 )
! ... "now" -> "old"
!
CALL davcio( evc, nwordwfc, iunwfc, ik, - 1 )
CALL davcio( evc, nwordwfc, iunoldwfc, ik, + 1 )
!
END DO
!
ELSE
CLOSE( UNIT = iunoldwfc, STATUS = 'KEEP' )
!
ELSE IF ( wfc_order == 2 ) THEN
!
CALL diropn( iunoldwfc, TRIM( prefix ) // '.oldwfc', nwordwfc, exst )
IF ( order > 2 ) &
CALL diropn( iunoldwfc2, TRIM( prefix ) // '.oldwfc2', nwordwfc, exst )
!
ALLOCATE( evcold(npwx,nbnd) )
!
IF ( order == 2 ) THEN
WRITE( stdout, '(5X,"Extrapolating wave-functions (first order) ...")' )
ELSE
WRITE( stdout, '(5X,"Extrapolating wave-functions (second order) ...")' )
END IF
WRITE( UNIT = stdout, &
FMT = '(5X,"Extrapolating wave-functions (first order) ...")' )
!
lwork = 5 * nbnd
!
lwork = 5*nbnd
ALLOCATE( s_m(nbnd,nbnd), sp_m(nbnd,nbnd), u_m(nbnd,nbnd), &
w_m(nbnd,nbnd), work(lwork), ew(nbnd), rwork(lwork) )
!
@ -270,15 +358,14 @@ SUBROUTINE extrapolate_wfcs()
IF ( nks > 1 ) READ( iunigk ) npw, igk
!
CALL davcio( evcold, nwordwfc, iunoldwfc, ik, - 1 )
CALL davcio( evc, nwordwfc, iunwfc, ik, - 1 )
CALL davcio( evc, nwordwfc, iunwfc, ik, - 1 )
!
! ... construct s_m = <evcold|evc>
!
CALL ZGEMM( 'C', 'N', nbnd, nbnd, npw, ONE, evcold, npwx, evc, &
npwx, ZERO, s_m, nbnd )
#ifdef __PARA
!
CALL reduce( 2 * nbnd * nbnd, s_m )
#endif
!
! ... construct sp_m
!
@ -318,49 +405,27 @@ SUBROUTINE extrapolate_wfcs()
!
! ... save on file the aligned wavefcts
!
CALL davcio( evcold, nwordwfc, iunwfc, ik, 1 )
CALL davcio( evcold, nwordwfc, iunwfc, ik, + 1 )
!
! ... re-read from file the wavefcts at (t-dt)
!
CALL davcio( evc, nwordwfc, iunoldwfc, ik, - 1 )
!
! ... extrapolate the wfc's,
! ... if order=3 use the second order extrapolation formula
! ... alpha0 and beta0 are calculated in "move_ions"
!
IF ( order > 2 ) THEN
!
IF ( istep == 2 ) THEN
!
evc = ( 1 + alpha0 ) * evcold - alpha0 * evc
!
ELSE
!
evc = ( 1 + alpha0 ) * evcold + ( beta0 - alpha0 ) * evc
!
CALL davcio( evcold, nwordwfc, iunoldwfc2, ik, - 1 )
!
evc = evc - beta0 * evcold
!
END IF
!
ELSE
!
evc = 2.D0 * evcold - evc
!
END IF
evc = 2.D0 * evcold - evc
!
! ... move the files: "old" -> "old1" and "now" -> "old"
!
IF ( order > 2 ) THEN
!
CALL davcio( evcold, nwordwfc, iunoldwfc, ik, - 1 )
CALL davcio( evcold, nwordwfc, iunoldwfc2, ik, 1 )
CALL davcio( evcold, nwordwfc, iunoldwfc, ik, - 1 )
CALL davcio( evcold, nwordwfc, iunoldwfc2, ik, + 1 )
!
END IF
!
CALL davcio( evcold, nwordwfc, iunwfc, ik, - 1 )
CALL davcio( evcold, nwordwfc, iunoldwfc, ik, 1 )
CALL davcio( evcold, nwordwfc, iunwfc, ik, - 1 )
CALL davcio( evcold, nwordwfc, iunoldwfc, ik, + 1 )
!
! ... save evc on file iunwfc
!
@ -377,6 +442,124 @@ SUBROUTINE extrapolate_wfcs()
!
DEALLOCATE( evcold )
!
CLOSE( UNIT = iunoldwfc, STATUS = 'KEEP' )
IF ( order > 2 ) &
CLOSE( UNIT = iunoldwfc2, STATUS = 'KEEP' )
!
ELSE
!
! ... case : wfc_order = 3
!
CALL diropn( iunoldwfc, TRIM( prefix ) // '.oldwfc', nwordwfc, exst )
CALL diropn( iunoldwfc2, TRIM( prefix ) // '.oldwfc2', nwordwfc, exst )
!
ALLOCATE( evcold(npwx,nbnd) )
!
WRITE( UNIT = stdout, &
FMT = '(5X,"Extrapolating wave-functions (second order) ...")' )
!
lwork = 5 * nbnd
!
ALLOCATE( s_m(nbnd,nbnd), sp_m(nbnd,nbnd), u_m(nbnd,nbnd), &
w_m(nbnd,nbnd), work(lwork), ew(nbnd), rwork(lwork) )
!
IF ( nks > 1 ) REWIND( iunigk )
!
zero_ew = 0
!
DO ik = 1, nks
!
IF ( nks > 1 ) READ( iunigk ) npw, igk
!
CALL davcio( evcold, nwordwfc, iunoldwfc, ik, - 1 )
CALL davcio( evc, nwordwfc, iunwfc, ik, - 1 )
!
! ... construct s_m = <evcold|evc>
!
CALL ZGEMM( 'C', 'N', nbnd, nbnd, npw, ONE, evcold, npwx, evc, &
npwx, ZERO, s_m, nbnd )
!
CALL reduce( 2 * nbnd * nbnd, s_m )
!
! ... construct sp_m
!
DO i = 1, nbnd
!
sp_m(:,i) = CONJG( s_m (i,:) )
!
END DO
!
! ... the unitary matrix [sp_m*s_m]^(-1/2)*sp_m (eq. 3.29)
! ... by means the singular value decomposition (SVD) of
! ... sp_m = u_m * diag(ew) * w_m
! ... becomes u_m * w_m
!
CALL ZGESVD( 'A', 'A', nbnd, nbnd, sp_m, nbnd, ew, u_m, nbnd, &
w_m, nbnd, work, lwork, rwork, info )
!
! ... check on eigenvalues
!
DO i = 1, nbnd
!
IF ( ew(i) < 0.1D0 ) zero_ew = zero_ew + 1
!
END DO
!
! ... use sp_m to store u_m * w_m
!
CALL ZGEMM( 'N', 'N', nbnd, nbnd, nbnd, ONE, u_m, nbnd, w_m, &
nbnd, ZERO, sp_m, nbnd )
!
! ... now use evcold as workspace to calculate "aligned" wavefcts:
!
! ... evcold_i = sum_j evc_j*sp_m_ji (eq.3.21)
!
CALL ZGEMM( 'N', 'N', npw, nbnd, nbnd, ONE, evc, npwx, sp_m, &
nbnd, ZERO, evcold, npwx )
!
! ... save on file the aligned wavefcts
!
CALL davcio( evcold, nwordwfc, iunwfc, ik, + 1 )
!
! ... re-read from file the wavefcts at (t-dt)
!
CALL davcio( evc, nwordwfc, iunoldwfc, ik, - 1 )
!
! ... extrapolate the wfc's,
! ... if wfc_order == 3 use the second order extrapolation formula
! ... alpha0 and beta0 are calculated in "move_ions"
!
evc = ( 1 + alpha0 ) * evcold + ( beta0 - alpha0 ) * evc
!
CALL davcio( evcold, nwordwfc, iunoldwfc2, ik, - 1 )
!
evc = evc - beta0 * evcold
!
! ... move the files: "old" -> "old1" and "now" -> "old"
!
CALL davcio( evcold, nwordwfc, iunoldwfc, ik, - 1 )
CALL davcio( evcold, nwordwfc, iunoldwfc2, ik, + 1 )
CALL davcio( evcold, nwordwfc, iunwfc, ik, - 1 )
CALL davcio( evcold, nwordwfc, iunoldwfc, ik, + 1 )
!
! ... save evc on file iunwfc
!
CALL davcio( evc, nwordwfc, iunwfc, ik, 1 )
!
END DO
!
IF ( zero_ew > 0 ) &
WRITE( stdout, '(/,5X,"Message from extrapolate_wfcs: ",/, &
& 5X,"the matrix <psi(t-dt)|psi(t)> has ",I2, &
& " zero eigenvalues")' ) zero_ew
!
DEALLOCATE( s_m, sp_m, u_m, w_m, work, ew, rwork )
!
DEALLOCATE( evcold )
!
CLOSE( UNIT = iunoldwfc, STATUS = 'KEEP' )
CLOSE( UNIT = iunoldwfc2, STATUS = 'KEEP' )
!
END IF
!
RETURN