Added to CP the possibility to remove the total torque acting on the system

(useful to simulate an isolated system). Documentation updated. C.S.


git-svn-id: http://qeforge.qe-forge.org/svn/q-e/trunk/espresso@3205 c92efa57-630b-4861-b058-cf58834340f0
This commit is contained in:
sbraccia 2006-06-25 23:15:20 +00:00
parent 6c06b0cbbc
commit eac4123a5b
3 changed files with 51 additions and 32 deletions

View File

@ -13,14 +13,13 @@ SUBROUTINE cprmain( tau_out, fion_out, etot_out )
!
USE kinds, ONLY : DP
USE constants, ONLY : bohr_radius_angs, amu_au
USE control_flags, ONLY : iprint, isave, thdyn, tpre, &
iprsta, tfor, tvlocw, &
taurdr, tprnfor, tsdc, lconstrain, lwf, &
lneb, lcoarsegrained, ndr, ndw, nomore, &
tsde, tortho, tnosee, tnosep, trane, &
tranp, tsdp, tcp, tcap, ampre, amprp, &
tnoseh, tolp, ortho_eps, ortho_max, &
printwfc
USE control_flags, ONLY : iprint, isave, thdyn, tpre, iprsta, &
tfor, tvlocw,remove_rigid_rot, taurdr, &
tprnfor, tsdc, lconstrain, lwf, lneb, &
lcoarsegrained, ndr, ndw, nomore, tsde, &
tortho, tnosee, tnosep, trane, tranp, &
tsdp, tcp, tcap, ampre, amprp, tnoseh, &
tolp, ortho_eps, ortho_max, printwfc
USE core, ONLY : nlcc_any, rhoc
USE uspp_param, ONLY : nhm, nh
USE cvan, ONLY : nvb, ish
@ -156,7 +155,7 @@ SUBROUTINE cprmain( tau_out, fion_out, etot_out )
!
! ... forces on ions
!
REAL(DP) :: maxfion
REAL(DP) :: maxfion, fion_tot(3)
!
! ... work variables
!
@ -168,11 +167,9 @@ SUBROUTINE cprmain( tau_out, fion_out, etot_out )
REAL(DP) :: delta_etot
REAL(DP) :: ftmp, enb, enbi
INTEGER :: is, nacc, ia, j, iter, i, isa, ipos
INTEGER :: k, ii, l, m, iss
!
INTEGER :: k, ii, l, m, iss
REAL(DP) :: hgamma(3,3), temphh(3,3)
REAL(DP) :: fcell(3,3)
!
REAL(DP) :: stress_gpa(3,3), thstress(3,3)
!
REAL(DP), ALLOCATABLE :: usrt_tau0(:,:), usrt_taup(:,:), usrt_fion(:,:)
@ -181,9 +178,10 @@ SUBROUTINE cprmain( tau_out, fion_out, etot_out )
REAL(DP), ALLOCATABLE :: tauw(:,:)
! temporary array used to printout positions
CHARACTER(LEN=3) :: labelw( nat )
! for force_pairing
! for force_pairing
INTEGER :: nspin_sub
!
!
dt2bye = dt2 / emass
etot_out = 0.D0
enow = 1.D9
@ -335,6 +333,13 @@ SUBROUTINE cprmain( tau_out, fion_out, etot_out )
!
IF ( lwf ) CALL ef_force( fion, na, nsp, zv )
!
fion_tot(:) = SUM( fion(:,:), DIM = 2 ) / DBLE( nat )
!
FORALL( ia = 1:nat ) fion(:,ia) = fion(:,ia) - fion_tot(:)
!
IF ( remove_rigid_rot ) &
CALL remove_tot_torque( nat, tau0, pmass(ityp(ind_srt(:))), fion )
!
IF ( lconstrain ) THEN
!
IF ( ionode ) THEN

View File

@ -261,6 +261,7 @@ MODULE input
etot_maxiter_ => etot_maxiter, &
forc_maxiter_ => forc_maxiter
USE control_flags, ONLY : force_pairing_ => force_pairing
USE control_flags, ONLY : remove_rigid_rot_ => remove_rigid_rot
!
! ... Other modules
!
@ -278,8 +279,9 @@ MODULE input
tefield2_ => tefield2, &
epol2_ => epol2, &
efield2_ => efield2
!
USE cvan, ONLY: nvb
!
USE input_parameters, ONLY: &
electron_dynamics, electron_damping, diis_rot, electron_temperature, &
ion_dynamics, ekin_conv_thr, etot_conv_thr, forc_conv_thr, ion_maxstep,&
@ -292,9 +294,7 @@ MODULE input
tdipole_card, toptical_card, tnewnfi_card, newnfi_card, &
ampre, nstep, restart_mode, ion_positions, startingwfc, printwfc, &
orthogonalization, electron_velocities, nat, if_pos, phase_space, &
tefield, epol, efield, tefield2, epol2, efield2
USE cvan, ONLY: nvb
tefield, epol, efield, tefield2, epol2, efield2, remove_rigid_rot
!
IMPLICIT NONE
!
@ -316,16 +316,15 @@ MODULE input
etot_conv_thr_ = etot_conv_thr
forc_conv_thr_ = forc_conv_thr
ekin_maxiter_ = electron_maxstep
tefield_ = tefield
epol_ = epol
efield_ = efield
tefield2_ = tefield2
epol2_ = epol2
efield2_ = efield2
!
remove_rigid_rot_ = remove_rigid_rot
!
tefield_ = tefield
epol_ = epol
efield_ = efield
tefield2_ = tefield2
epol2_ = epol2
efield2_ = efield2
!
! ... Set internal time step variables ( delt, twodelt, dt2 ... )
!
@ -354,8 +353,8 @@ MODULE input
trhor_ = ( TRIM( calculation ) == 'nscf' )
tvlocw_ = ( TRIM( disk_io ) == 'high' ) ! warning this is not working
! for now use reduce_io in CP to specify if the charge density is saved
! reduce_io_ = ( TRIM( disk_io ) == 'low' )
reduce_io_ = .not.( TRIM( disk_io ) == 'high' )
!
reduce_io_ = .NOT.( TRIM( disk_io ) == 'high' )
!
SELECT CASE( TRIM( verbosity ) )
CASE( 'minimal' )
@ -661,6 +660,8 @@ MODULE input
CALL errore(' control_flags ',' unknown ion_dynamics '//TRIM(ion_dynamics), 1 )
END SELECT
IF ( ANY( if_pos(:,1:nat) == 0 ) ) lfixatom = .TRUE.
! ... Ionic Temperature
@ -806,8 +807,6 @@ MODULE input
force_pairing_ = force_pairing
!
!
!
IF( ( nvb > 0 ) .and. ( program_name == 'FPMD' ) ) &
CALL errore(' iosys ',' USPP not yet implemented in FPMD ',1)

View File

@ -564,6 +564,21 @@ amprp(i) REAL ( default = 0.D0 )
greasp REAL ( default = 1.D0 )
same as "grease", for ionic damped dynamics.
remove_rigid_rot
LOGICAL ( default = .FALSE. )
this keyword is useful when simulating the dynamics and/or the
thermodynamics of an isolated system. If set to true the total
torque of the internal forces is set to zero by adding new forces
that compensate the spurious interaction with the periodic
images. This allowes for the use of smaller supercells.
BEWARE: since the potential energy is no longer consistent with
the forces (it still contains the spurious interaction with the
repeated images), the total energy is not conserved anymore.
However the dynamical and thermodynamical properties should be
in closer agreement with those of an isolated system.
Also the final energy of a structural relaxation will be higher,
but the relaxation itself should be faster.
!
! ... keywords used only in NEB calculations
!