TABLE OF CONTENTS


ABINIT/paw2wvl [ Functions ]

[ Top ] [ Functions ]

NAME

  paw2wvl

FUNCTION

  Points WVL objects to PAW objects

COPYRIGHT

  Copyright (C) 2011-2018 ABINIT group (T. Rangel)
  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

  argin(sizein)=description

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

      gstate

CHILDREN

SOURCE

 32 #if defined HAVE_CONFIG_H
 33 #include "config.h"
 34 #endif
 35 
 36 #include "abi_common.h"
 37 
 38 subroutine paw2wvl(pawtab,proj,wvl)
 39 
 40  use m_profiling_abi
 41  use m_errors
 42  use defs_basis
 43  use m_pawtab, only : pawtab_type
 44  use defs_wvltypes
 45 
 46 !This section has been created automatically by the script Abilint (TD).
 47 !Do not modify the following lines by hand.
 48 #undef ABI_FUNC
 49 #define ABI_FUNC 'paw2wvl'
 50 !End of the abilint section
 51 
 52  implicit none
 53 
 54 !Arguments ------------------------------------
 55  type(pawtab_type),intent(in)::pawtab(:)
 56  type(wvl_internal_type), intent(inout) :: wvl
 57  type(wvl_projectors_type),intent(inout)::proj
 58 
 59 !Local variables-------------------------------
 60 #if defined HAVE_BIGDFT
 61  integer :: ib,ig
 62  integer :: itypat,jb,ll,lmnmax,lmnsz,nn,ntypat,ll_,nn_
 63  integer ::maxmsz,msz1,ptotgau,max_lmn2_size
 64  logical :: test_wvl
 65  real(dp) :: a1
 66  integer,allocatable :: msz(:)
 67 #endif
 68 
 69 !extra variables, use to debug
 70 !  integer::nr,unitp,i_shell,ii,ir,ng
 71 !  real(dp)::step,rmax
 72 !  real(dp),allocatable::r(:), y(:)
 73 !  complex::fac,arg
 74 !  complex(dp),allocatable::f(:),g(:,:)
 75 
 76 ! *********************************************************************
 77 
 78 !DEBUG
 79 !write (std_out,*) ' paw2wvl : enter'
 80 !ENDDEBUG
 81 
 82 #if defined HAVE_BIGDFT
 83 
 84  ntypat=size(pawtab)
 85 
 86 !wvl object must be allocated in pawtab
 87  if (ntypat>0) then
 88    test_wvl=.true.
 89    do itypat=1,ntypat
 90      if (pawtab(itypat)%has_wvl==0.or.(.not.associated(pawtab(itypat)%wvl))) test_wvl=.false.
 91    end do
 92    if (.not.test_wvl) then
 93      MSG_BUG('pawtab%wvl must be allocated!')
 94    end if
 95  end if
 96 
 97 !Find max mesh size
 98  ABI_ALLOCATE(msz,(ntypat))
 99  do itypat=1,ntypat
100    msz(itypat)=pawtab(itypat)%wvl%rholoc%msz
101  end do
102  maxmsz=maxval(msz(1:ntypat))
103 
104  ABI_ALLOCATE(wvl%rholoc%msz,(ntypat))
105  ABI_ALLOCATE(wvl%rholoc%d,(maxmsz,4,ntypat))
106  ABI_ALLOCATE(wvl%rholoc%rad,(maxmsz,ntypat))
107  ABI_ALLOCATE(wvl%rholoc%radius,(ntypat))
108 
109  do itypat=1,ntypat
110    msz1=pawtab(itypat)%wvl%rholoc%msz
111    msz(itypat)=msz1
112    wvl%rholoc%msz(itypat)=msz1
113    wvl%rholoc%d(1:msz1,:,itypat)=pawtab(itypat)%wvl%rholoc%d(1:msz1,:)
114    wvl%rholoc%rad(1:msz1,itypat)=pawtab(itypat)%wvl%rholoc%rad(1:msz1)
115    wvl%rholoc%radius(itypat)=pawtab(itypat)%rpaw
116  end do
117  ABI_DEALLOCATE(msz)
118 !
119 !Now fill projectors type:
120 !
121  ABI_DATATYPE_ALLOCATE(proj%G,(ntypat))
122 !
123 !nullify all for security:
124 !do itypat=1,ntypat
125 !call nullify_gaussian_basis(proj%G(itypat))
126 !end do
127  proj%G(:)%nat=1 !not used
128  proj%G(:)%ncplx=2 !Complex gaussians
129 !
130 !Obtain dimensions:
131 !jb=0
132 !ptotgau=0
133  do itypat=1,ntypat
134 !  do ib=1,pawtab(itypat)%basis_size
135 !  jb=jb+1
136 !  end do
137 !  ptotgau=ptotgau+pawtab(itypat)%wvl%ptotgau
138    ptotgau=pawtab(itypat)%wvl%ptotgau
139    proj%G(itypat)%nexpo=ptotgau
140    proj%G(itypat)%nshltot=pawtab(itypat)%basis_size
141  end do
142 !proj%G%nexpo=ptotgau
143 !proj%G%nshltot=jb
144 !
145 !Allocations
146  do itypat=1,ntypat
147    ABI_ALLOCATE(proj%G(itypat)%ndoc ,(proj%G(itypat)%nshltot))
148    ABI_ALLOCATE(proj%G(itypat)%nam  ,(proj%G(itypat)%nshltot))
149    ABI_ALLOCATE(proj%G(itypat)%xp   ,(proj%G(itypat)%ncplx,proj%G(itypat)%nexpo))
150    ABI_ALLOCATE(proj%G(itypat)%psiat,(proj%G(itypat)%ncplx,proj%G(itypat)%nexpo))
151  end do
152 
153 !jb=0
154  do itypat=1,ntypat
155    proj%G(itypat)%ndoc(:)=pawtab(itypat)%wvl%pngau(:)
156  end do
157 !
158  do itypat=1,ntypat
159    proj%G(itypat)%xp(:,:)=pawtab(itypat)%wvl%parg(:,:)
160    proj%G(itypat)%psiat(:,:)=pawtab(itypat)%wvl%pfac(:,:)
161  end do
162 
163 !Change the real part of psiat to the form adopted in BigDFT:
164 !Here we use exp{(a+ib)x^2}, where a is a negative number.
165 !In BigDFT: exp{-0.5(x/c)^2}exp{i(bx)^2}, and c is positive
166 !Hence c=sqrt(0.5/abs(a))
167  do itypat=1,ntypat
168    do ig=1,proj%G(itypat)%nexpo
169      a1=proj%G(itypat)%xp(1,ig)
170      a1=(sqrt(0.5/abs(a1)))
171      proj%G(itypat)%xp(1,ig)=a1
172    end do
173  end do
174 
175 !debug
176 !write(*,*)'paw2wvl 178: erase me set gaussians real equal to hgh (for Li)'
177 !ABI_DEALLOCATE(proj%G(1)%xp)
178 !ABI_DEALLOCATE(proj%G(1)%psiat)
179 !ABI_ALLOCATE(proj%G(1)%psiat,(2,2))
180 !ABI_ALLOCATE(proj%G(1)%xp,(2,2))
181 !proj%G(1)%ndoc(:)=1 !1 gaussian per shell
182 !proj%G(1)%nexpo=2 ! two gaussians in total
183 !proj%G(1)%xp=zero
184 !proj%G(1)%xp(1,1)=0.666375d0
185 !proj%G(1)%xp(1,2)=1.079306d0
186 !proj%G(1)%psiat(1,:)=1.d0
187 !proj%G(1)%psiat(2,:)=zero
188 
189 !
190 !begin debug
191 !
192 !write(*,*)'paw2wvl, comment me'
193 !and comment out variables 
194 !rmax=2.d0
195 !nr=rmax/(0.0001d0)
196 !ABI_ALLOCATE(r,(nr))
197 !ABI_ALLOCATE(f,(nr))
198 !ABI_ALLOCATE(y,(nr))
199 !step=rmax/real(nr-1,dp)
200 !do ir=1,nr
201 !r(ir)=real(ir-1,dp)*step
202 !end do
203 !!
204 !unitp=400
205 !do itypat=1,ntypat
206 !ig=0
207 !do i_shell=1,proj%G(itypat)%nshltot
208 !unitp=unitp+1
209 !f(:)=czero
210 !ng=proj%G(itypat)%ndoc(i_shell)
211 !ABI_ALLOCATE(g,(nr,ng))
212 !do ii=1,ng
213 !ig=ig+1
214 !fac=cmplx(proj%G(itypat)%psiat(1,ig),proj%G(itypat)%psiat(2,ig))
215 !a1=-0.5d0/(proj%G(itypat)%xp(1,ig)**2)
216 !arg=cmplx(a1,proj%G(itypat)%xp(2,ig))
217 !g(:,ii)=fac*exp(arg*r(:)**2)
218 !f(:)=f(:)+g(:,ii)
219 !end do
220 !do ir=1,nr
221 !write(unitp,'(9999f16.7)')r(ir),real(f(ir)),(real(g(ir,ii)),&
222 !&    ii=1,ng)
223 !end do
224 !ABI_DEALLOCATE(g)
225 !end do
226 !end do
227 !ABI_DEALLOCATE(r)
228 !ABI_DEALLOCATE(f)
229 !ABI_DEALLOCATE(y)
230 !stop 
231 !end debug
232 
233 
234 !jb=0
235  do itypat=1,ntypat
236    jb=0
237    ll_ = 0
238    nn_ = 0
239    do ib=1,pawtab(itypat)%lmn_size 
240      ll=pawtab(itypat)%indlmn(1,ib)
241      nn=pawtab(itypat)%indlmn(3,ib)
242 !    write(*,*)ll,pawtab(itypat)%indlmn(2,ib),nn
243      if(ib>1 .and. ll == ll_ .and. nn==nn_) cycle
244      jb=jb+1
245 !    proj%G%nam(jb)= pawtab(itypat)%indlmn(1,ib)+1!l quantum number
246 !    write(*,*)jb,proj%G%nshltot,proj%G%nam(jb)
247      proj%G(itypat)%nam(jb)= pawtab(itypat)%indlmn(1,ib)+1 !l quantum number
248 !    1 is added due to BigDFT convention
249 !    write(*,*)jb,pawtab(itypat)%indlmn(1:3,ib)
250      ll_=ll
251      nn_=nn
252    end do
253  end do
254 !Nullify remaining objects
255  do itypat=1,ntypat
256    nullify(proj%G(itypat)%rxyz)
257    nullify(proj%G(itypat)%nshell)
258  end do
259 !nullify(proj%G%rxyz)
260 
261 !now the index l,m,n objects:
262  lmnmax=maxval(pawtab(:)%lmn_size)
263  ABI_ALLOCATE(wvl%paw%indlmn,(6,lmnmax,ntypat))
264  wvl%paw%lmnmax=lmnmax
265  wvl%paw%ntypes=ntypat
266  wvl%paw%indlmn=0
267  do itypat=1,ntypat
268    lmnsz=pawtab(itypat)%lmn_size
269    wvl%paw%indlmn(1:6,1:lmnsz,itypat)=pawtab(itypat)%indlmn(1:6,1:lmnsz)
270  end do
271  
272 !allocate and copy sij
273 !max_lmn2_size=max(pawtab(:)%lmn2_size,1)
274  max_lmn2_size=lmnmax*(lmnmax+1)/2
275  ABI_ALLOCATE(wvl%paw%sij,(max_lmn2_size,ntypat))
276 !sij is not yet calculated here.
277 !We copy this in gstate after call to pawinit
278 !do itypat=1,ntypat
279 !wvl%paw%sij(1:pawtab(itypat)%lmn2_size,itypat)=pawtab(itypat)%sij(:)
280 !end do
281 
282  ABI_ALLOCATE(wvl%paw%rpaw,(ntypat))
283  do itypat=1,ntypat
284    wvl%paw%rpaw(itypat)= pawtab(itypat)%rpaw
285  end do
286 
287  if(allocated(wvl%npspcode_paw_init_guess)) then
288    ABI_DEALLOCATE(wvl%npspcode_paw_init_guess)
289  end if
290  ABI_ALLOCATE(wvl%npspcode_paw_init_guess,(ntypat))
291  do itypat=1,ntypat
292    wvl%npspcode_paw_init_guess(itypat)=pawtab(itypat)%wvl%npspcode_init_guess
293  end do
294 
295 #else
296  if (.false.) write(std_out) pawtab(1)%mesh_size,wvl%h(1),proj%nlpsp
297 #endif
298 
299 !DEBUG
300 !write (std_out,*) ' paw2wvl : exit'
301 !stop
302 !ENDDEBUG
303 
304 end subroutine paw2wvl