1! 2! Copyright (C) 2002-2011 Quantum ESPRESSO group 3! This file is distributed under the terms of the 4! GNU General Public License. See the file `License' 5! in the root directory of the present distribution, 6! or http://www.gnu.org/copyleft/gpl.txt . 7! 8!------------------------------------------------------------------------------! 9 MODULE cell_nose 10!------------------------------------------------------------------------------! 11 12 USE kinds, ONLY : DP 13! 14 IMPLICIT NONE 15 SAVE 16 17 REAL(DP) :: xnhh0(3,3) = 0.0_DP 18 REAL(DP) :: xnhhm(3,3) = 0.0_DP 19 REAL(DP) :: xnhhp(3,3) = 0.0_DP 20 REAL(DP) :: vnhh(3,3) = 0.0_DP 21 REAL(DP) :: temph = 0.0_DP ! Thermostat temperature (from input) 22 REAL(DP) :: fnoseh = 0.0_DP ! Thermostat frequency (from input) 23 REAL(DP) :: qnh = 0.0_DP ! Thermostat mass (computed) 24 25CONTAINS 26 27 subroutine cell_nose_init( temph_init, fnoseh_init ) 28 USE constants, ONLY: pi, au_terahertz, k_boltzmann_au 29 REAL(DP), INTENT(IN) :: temph_init, fnoseh_init 30 ! set thermostat parameter for cell 31 qnh = 0.0_DP 32 temph = temph_init 33 fnoseh = fnoseh_init 34 if( fnoseh > 0.0_DP ) qnh = 2.0_DP * ( 3 * 3 ) * temph * k_boltzmann_au / & 35 (fnoseh*(2.0_DP*pi)*au_terahertz)**2 36 return 37 end subroutine cell_nose_init 38 39 subroutine cell_nosezero( vnhh, xnhh0, xnhhm ) 40 real(DP), intent(out) :: vnhh(3,3), xnhh0(3,3), xnhhm(3,3) 41 xnhh0=0.0_DP 42 xnhhm=0.0_DP 43 vnhh =0.0_DP 44 return 45 end subroutine cell_nosezero 46 47 subroutine cell_nosevel( vnhh, xnhh0, xnhhm, delt ) 48 implicit none 49 REAL(DP), intent(inout) :: vnhh(3,3) 50 REAL(DP), intent(in) :: xnhh0(3,3), xnhhm(3,3), delt 51 vnhh(:,:)=2.0_DP*(xnhh0(:,:)-xnhhm(:,:))/delt-vnhh(:,:) 52 return 53 end subroutine cell_nosevel 54 55 subroutine cell_noseupd( xnhhp, xnhh0, xnhhm, delt, qnh, temphh, temph, vnhh ) 56 use constants, only: k_boltzmann_au 57 implicit none 58 REAL(DP), intent(out) :: xnhhp(3,3), vnhh(3,3) 59 REAL(DP), intent(in) :: xnhh0(3,3), xnhhm(3,3), delt, qnh, temphh(3,3), temph 60 integer :: i, j 61 do j=1,3 62 do i=1,3 63 xnhhp(i,j) = 2.0_DP*xnhh0(i,j)-xnhhm(i,j) + & 64 (delt**2/qnh)* k_boltzmann_au * (temphh(i,j)-temph) 65 vnhh(i,j) =(xnhhp(i,j)-xnhhm(i,j))/( 2.0_DP * delt ) 66 end do 67 end do 68 return 69 end subroutine cell_noseupd 70 71 72 REAL(DP) function cell_nose_nrg( qnh, xnhh0, vnhh, temph, iforceh ) 73 use constants, only: k_boltzmann_au 74 implicit none 75 REAL(DP) :: qnh, vnhh( 3, 3 ), temph, xnhh0( 3, 3 ) 76 integer :: iforceh( 3, 3 ) 77 integer :: i, j 78 REAL(DP) :: enij 79 cell_nose_nrg = 0.0_DP 80 do i=1,3 81 do j=1,3 82 enij = 0.5_DP*qnh*vnhh(i,j)*vnhh(i,j)+temph*k_boltzmann_au*xnhh0(i,j) 83 cell_nose_nrg = cell_nose_nrg + iforceh( i, j ) * enij 84 enddo 85 enddo 86 return 87 end function cell_nose_nrg 88 89 subroutine cell_nose_shiftvar( xnhhp, xnhh0, xnhhm ) 90 ! shift values of nose variables to start a new step 91 implicit none 92 REAL(DP), intent(out) :: xnhhm(3,3) 93 REAL(DP), intent(inout) :: xnhh0(3,3) 94 REAL(DP), intent(in) :: xnhhp(3,3) 95 xnhhm = xnhh0 96 xnhh0 = xnhhp 97 return 98 end subroutine cell_nose_shiftvar 99 100 101 SUBROUTINE cell_nose_info ( delt ) 102 103 use constants, only: au_terahertz, pi 104 USE io_global, ONLY: stdout 105 USE control_flags, ONLY: tnoseh 106 107 IMPLICIT NONE 108 109 REAL(DP), INTENT (IN) :: delt 110 111 INTEGER :: nsvar 112 REAL(DP) :: wnoseh 113 114 IF( tnoseh ) THEN 115 ! 116 IF( fnoseh <= 0.0_DP) & 117 CALL errore(' cell_nose_info ', ' fnoseh less than zero ', 1) 118 IF( delt <= 0.0_DP) & 119 CALL errore(' cell_nose_info ', ' delt less than zero ', 1) 120 121 wnoseh = fnoseh * ( 2.0_DP * pi ) * au_terahertz 122 nsvar = ( 2.0_DP * pi ) / ( wnoseh * delt ) 123 124 WRITE( stdout,563) temph, nsvar, fnoseh, qnh 125 END IF 126 127 563 format( //, & 128 & 3X,'cell dynamics with nose` temperature control:', /, & 129 & 3X,'Kinetic energy required = ', f10.5, ' (Kelvin) ', /, & 130 & 3X,'time steps per nose osc. = ', i5, /, & 131 & 3X,'nose` frequency = ', f10.3, ' (THz) ', /, & 132 & 3X,'nose` mass(es) = ', 20(1X,f10.3),//) 133 134 RETURN 135 END SUBROUTINE cell_nose_info 136 137 138! 139!------------------------------------------------------------------------------! 140 END MODULE cell_nose 141!------------------------------------------------------------------------------! 142 143