TABLE OF CONTENTS


ABINIT/etheta [ Functions ]

[ Top ] [ 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