TABLE OF CONTENTS
ABINIT/etheta [ Functions ]
NAME
etheta
FUNCTION
Computes the energy per unit cell and its first derivative for a given angle theta. More precisely, computes only the part of the energy that changes with theta.
COPYRIGHT
Copyright (C) 2000-2017 ABINIT group (MVeithen,ISouza,JIniguez) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
bcut(ifor,idir) = branch cut of the ellipse associated with (ifor,idir) chc = <C|H_0|C> where |C> is the wavefunction of the current band detovc = determinant of the overlap matrix S detovd = determinant of the overlap matrix where for the band that is being updated <C| is replaced by <D| (search direction) dhc = Re[<D|H_0|C>] dhd = <D|H_0|D> efield_dot = reciprocal lattice coordinates of the electric field hel(ifor,idir) = helicity of the ellipse associated with (ifor,idir) nkpt = number of k-points nsppol = 1 for unpolarized, 2 for spin-polarized nstr(idir) = number of strings along the idir-th direction sdeg = spin degeneracy theta = value of the angle for which the energy (e0) and its derivative (e1) are computed
OUTPUT
e0 = energy for the given value of theta e1 = derivative of the energy with respect to theta
PARENTS
cgwf,linemin
CHILDREN
rhophi
SOURCE
47 #if defined HAVE_CONFIG_H 48 #include "config.h" 49 #endif 50 51 #include "abi_common.h" 52 53 54 subroutine etheta(bcut,chc,detovc,detovd,dhc,dhd,efield_dot,e0,e1,& 55 & hel,nkpt,nstr,sdeg,theta) 56 57 use defs_basis 58 use m_profiling_abi 59 60 use m_numeric_tools, only : rhophi 61 62 !This section has been created automatically by the script Abilint (TD). 63 !Do not modify the following lines by hand. 64 #undef ABI_FUNC 65 #define ABI_FUNC 'etheta' 66 !End of the abilint section 67 68 implicit none 69 70 !Arguments ------------------------------------ 71 !scalars 72 integer,intent(in) :: nkpt 73 real(dp),intent(in) :: chc,dhc,dhd,sdeg,theta 74 real(dp),intent(out) :: e0,e1 75 !arrays 76 integer,intent(in) :: hel(2,3),nstr(3) 77 real(dp),intent(in) :: bcut(2,3),detovc(2,2,3),detovd(2,2,3),efield_dot(3) 78 79 !Local variables ------------------------- 80 !scalars 81 integer :: idir,ifor 82 real(dp) :: c2theta,ctheta,dphase,gnorm,phase,rho,s2theta,sgn,stheta 83 !arrays 84 real(dp) :: dg_theta(2),g_theta(2) 85 86 ! *********************************************************************** 87 88 e0 = zero ; e1 = zero 89 90 ctheta = cos(theta) 91 stheta = sin(theta) 92 c2theta = ctheta*ctheta - stheta*stheta ! cos(2*theta) 93 s2theta = two*ctheta*stheta ! sin(2*theta) 94 95 e0 = chc*ctheta*ctheta + dhd*stheta*stheta + dhc*s2theta 96 e0 = e0*sdeg/nkpt 97 98 !DEBUG 99 !e0 = zero 100 !ENDDEBUG 101 102 e1 = (dhd - chc)*s2theta + two*dhc*c2theta 103 e1 = e1*sdeg/nkpt 104 105 sgn = -1_dp 106 do idir = 1, 3 107 108 if (abs(efield_dot(idir)) < tol12) cycle 109 110 do ifor = 1, 2 111 112 g_theta(:) = ctheta*detovc(:,ifor,idir) + & 113 & stheta*detovd(:,ifor,idir) 114 dg_theta(:) = -1_dp*stheta*detovc(:,ifor,idir) + & 115 & ctheta*detovd(:,ifor,idir) 116 117 ! Compute E(theta) 118 119 call rhophi(g_theta,phase,rho) 120 if (theta >= bcut(ifor,idir)) phase = phase + hel(ifor,idir)*two_pi 121 122 ! DEBUG 123 ! unit = 100 + 10*idir + ifor 124 ! write(unit,'(4(f16.9))')theta,g_theta(:),phase 125 ! ENDDEBUG 126 127 e0 = e0 + sgn*sdeg*efield_dot(idir)*phase/(two_pi*nstr(idir)) 128 129 130 ! Compute dE/dtheta 131 132 ! imaginary part of the derivative of ln(g_theta) 133 gnorm = g_theta(1)*g_theta(1) + g_theta(2)*g_theta(2) 134 dphase = (dg_theta(2)*g_theta(1) - dg_theta(1)*g_theta(2))/gnorm 135 136 e1 = e1 + sgn*sdeg*efield_dot(idir)*dphase/(two_pi*nstr(idir)) 137 138 sgn = -1_dp*sgn 139 140 end do 141 end do 142 143 end subroutine etheta