TABLE OF CONTENTS


ABINIT/dfpt_v1zeeman [ Functions ]

[ Top ] [ Functions ]

NAME

  dfpt_v1zeeman

FUNCTION

  Calculate 1st order Zeeman potential = -vec{\sigma}.\vec{b}, where 
  sigma is the vector of Pauli matrices and \vec{b} is the unit 
  vector indicating the perturbing field direction.

COPYRIGHT

  Copyright (C) 2017-2018 ABINIT group (SPr)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

  nspden = number of density matrix components
  nfft   = numbder of fft grid points
  cplex  = complex or real density matrix
  idir   = direction of the perturbing field in Cartesian frame
           1: along x
           2: along y
           3: along z
           4: identity matrix at each fft point is returned (for density-density response)

OUTPUT

  v1zeeman(nfft*cplex,nspden)= 1st order Zeeman potential, or Identity matrix (electrostatic potential) for idir=4

SIDE EFFECTS

  None

NOTES

  The definition of components of the potential matrix differ depending on cplex
  for nspden=4:
  For cplex=1, the potential is defined as (V_upup,V_dndn,Re[V_updn],Im[V_updn])
  For cplex=2, the definition is (V_upup,V_dndn,V_updn,i.V_updn)

PARENTS

      dfpt_rhotov

CHILDREN

SOURCE

 48 #if defined HAVE_CONFIG_H
 49 #include "config.h"
 50 #endif
 51 
 52 #include "abi_common.h"
 53 
 54 
 55 subroutine dfpt_v1zeeman(nspden,nfft,cplex,idir,v1zeeman)
 56     
 57  use defs_basis
 58  use m_errors
 59  use m_profiling_abi
 60 
 61 !This section has been created automatically by the script Abilint (TD).
 62 !Do not modify the following lines by hand.
 63 #undef ABI_FUNC
 64 #define ABI_FUNC 'dfpt_v1zeeman'
 65 !End of the abilint section
 66 
 67  implicit none
 68 
 69 !Arguments ------------------------------------
 70  integer , intent(in)    :: idir,nfft,cplex,nspden
 71  real(dp), intent(inout) :: v1zeeman(cplex*nfft,nspden)                       
 72 
 73 !Local variables-------------------------------
 74  integer :: ifft                                
 75 !character(len=500) :: msg                   
 76  
 77 ! *************************************************************************
 78 
 79  DBG_ENTER("COLL")
 80  
 81 ! if (option/=1 .and. option/=2 ) then
 82 !   write(msg,'(3a,i0)')&
 83 !&   'The argument option should be 1 or 2,',ch10,&
 84 !&   'however, option=',option
 85 !   MSG_BUG(msg)
 86 ! end if
 87 !
 88 ! if (sizein<1) then
 89 !   write(msg,'(3a,i0)')&
 90 !&   'The argument sizein should be a positive number,',ch10,&
 91 !&   'however, sizein=',sizein
 92 !   MSG_ERROR(msg)
 93 ! end if
 94 
 95  DBG_EXIT("COLL")
 96 
 97  select case(cplex)
 98  case(1)
 99    if (nspden==4) then
100      if(idir==3)then       ! Zeeman field along the 3rd axis (z)   
101        v1zeeman(:,1)=-0.5d0
102        v1zeeman(:,2)=+0.5d0
103        v1zeeman(:,3)= 0.0d0
104        v1zeeman(:,4)= 0.0d0
105      else if(idir==2)then  ! Zeeman field along the 2nd axis (y)
106        v1zeeman(:,1)= 0.0d0
107        v1zeeman(:,2)= 0.0d0
108        v1zeeman(:,3)= 0.0d0
109        v1zeeman(:,4)=+0.5d0
110      else                  ! Zeeman field along the 1st axis (x)
111        v1zeeman(:,1)= 0.0d0
112        v1zeeman(:,2)= 0.0d0
113        v1zeeman(:,3)=-0.5d0
114        v1zeeman(:,4)= 0.0d0
115      end if
116    else if (nspden==2) then
117      v1zeeman(:,1)=-0.5e0
118      v1zeeman(:,2)= 0.5e0
119    else
120      v1zeeman(:,1)= 0.0e0
121    end if
122  case(2)
123    if (nspden==2) then
124      do ifft=1,nfft
125        v1zeeman(2*ifft-1,1)  =-0.5e0
126        v1zeeman(2*ifft  ,1)  = 0.0e0
127        v1zeeman(2*ifft-1,2)  = 0.5e0
128        v1zeeman(2*ifft  ,2)  = 0.0e0
129      end do
130    else if (nspden==4) then
131      select case(idir)
132      case(1) !along x, v1=-sigma_x
133        do ifft=1,nfft
134          v1zeeman(2*ifft-1,1)= 0.0e0 !Re[V^11]
135          v1zeeman(2*ifft  ,1)= 0.0e0 !Im[V^11]
136          v1zeeman(2*ifft-1,2)= 0.0e0 !Re[V^22]
137          v1zeeman(2*ifft  ,2)= 0.0e0 !Im[V^22]
138          v1zeeman(2*ifft-1,3)=-0.5e0 !Re[V^12]
139          v1zeeman(2*ifft  ,3)= 0.0e0 !Im[V^12]
140          v1zeeman(2*ifft-1,4)= 0.0e0 !Re[i.V^21]=Im[V^12]
141          v1zeeman(2*ifft  ,4)=-0.5e0 !Im[i.V^21]=Re[V^12]
142        end do
143      case(2) !along y, v1 = -sigma_y
144        do ifft=1,nfft
145          v1zeeman(2*ifft-1,1)= 0.0e0 !Re[V^11]
146          v1zeeman(2*ifft  ,1)= 0.0e0 !Im[V^11]
147          v1zeeman(2*ifft-1,2)= 0.0e0 !Re[V^22]
148          v1zeeman(2*ifft  ,2)= 0.0e0 !Im[V^22]
149          v1zeeman(2*ifft-1,3)= 0.0e0 !Re[V^12]
150          v1zeeman(2*ifft  ,3)=+0.5e0 !Im[V^12]
151          v1zeeman(2*ifft-1,4)=+0.5e0 !Re[i.V^21]=Im[V^12]
152          v1zeeman(2*ifft  ,4)= 0.0e0 !Im[i.V^21]=Re[V^12]
153        end do
154      case(3)
155        do ifft=1,nfft
156          v1zeeman(2*ifft-1,1)=-0.5e0 !Re[V^11]
157          v1zeeman(2*ifft  ,1)= 0.0e0 !Im[V^11]
158          v1zeeman(2*ifft-1,2)= 0.5e0 !Re[V^22]
159          v1zeeman(2*ifft  ,2)= 0.0e0 !Im[V^22]
160          v1zeeman(2*ifft-1,3)= 0.0e0 !Re[V^12]
161          v1zeeman(2*ifft  ,3)= 0.0e0 !Im[V^12]
162          v1zeeman(2*ifft-1,4)= 0.0e0 !Re[i.V^21]
163          v1zeeman(2*ifft  ,4)= 0.0e0 !Im[i.V^21]
164        end do
165      end select
166    end if
167  end select !cplex
168 
169 end subroutine dfpt_v1zeeman