TABLE OF CONTENTS


ABINIT/m_paw2wvl [ Modules ]

[ Top ] [ Modules ]

NAME

  paw2wvl

FUNCTION

COPYRIGHT

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

PARENTS

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_paw2wvl
28 
29  use defs_basis
30  use defs_wvltypes
31  use m_abicore
32  use m_errors
33 
34  use m_pawtab, only : pawtab_type
35  use m_paw_ij, only : paw_ij_type
36 
37  implicit none
38 
39  private

ABINIT/paw2wvl [ Functions ]

[ Top ] [ Functions ]

NAME

  paw2wvl

FUNCTION

  Points WVL objects to PAW objects

INPUTS

  argin(sizein)=description

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

      gstate

CHILDREN

SOURCE

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

ABINIT/paw2wvl_ij [ Functions ]

[ Top ] [ Functions ]

NAME

  paw2wvl_ij

FUNCTION

  FIXME: add description.

INPUTS

  argin(sizein)=description

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

      scfcv

CHILDREN

      nullify_paw_ij_objects

SOURCE

364 subroutine paw2wvl_ij(option,paw_ij,wvl)
365 
366 #if defined HAVE_BIGDFT
367  use BigDFT_API, only : nullify_paw_ij_objects
368 #endif
369 
370 !This section has been created automatically by the script Abilint (TD).
371 !Do not modify the following lines by hand.
372 #undef ABI_FUNC
373 #define ABI_FUNC 'paw2wvl_ij'
374 !End of the abilint section
375 
376  implicit none
377 
378 !Arguments ------------------------------------
379  integer,intent(in)::option
380  type(wvl_internal_type), intent(inout)::wvl
381  type(paw_ij_type),intent(in) :: paw_ij(:)
382 !Local variables-------------------------------
383 #if defined HAVE_BIGDFT
384  integer :: iatom,iaux,my_natom
385  character(len=500) :: message
386 #endif
387 
388 ! *************************************************************************
389 
390  DBG_ENTER("COLL")
391 
392 #if defined HAVE_BIGDFT
393  my_natom=size(paw_ij)
394 
395 !Option==1: allocate and copy
396  if(option==1) then
397    ABI_DATATYPE_ALLOCATE(wvl%paw%paw_ij,(my_natom))
398    do iatom=1,my_natom
399      call nullify_paw_ij_objects(wvl%paw%paw_ij(iatom))
400      wvl%paw%paw_ij(iatom)%cplex          =paw_ij(iatom)%cplex_rf
401      wvl%paw%paw_ij(iatom)%cplex_dij      =paw_ij(iatom)%cplex_dij
402      wvl%paw%paw_ij(iatom)%has_dij        =paw_ij(iatom)%has_dij
403      wvl%paw%paw_ij(iatom)%has_dijfr      =0
404      wvl%paw%paw_ij(iatom)%has_dijhartree =0
405      wvl%paw%paw_ij(iatom)%has_dijhat     =0
406      wvl%paw%paw_ij(iatom)%has_dijso      =0
407      wvl%paw%paw_ij(iatom)%has_dijU       =0
408      wvl%paw%paw_ij(iatom)%has_dijxc      =0
409      wvl%paw%paw_ij(iatom)%has_dijxc_val  =0
410      wvl%paw%paw_ij(iatom)%has_exexch_pot =0
411      wvl%paw%paw_ij(iatom)%has_pawu_occ   =0
412      wvl%paw%paw_ij(iatom)%lmn_size       =paw_ij(iatom)%lmn_size
413      wvl%paw%paw_ij(iatom)%lmn2_size      =paw_ij(iatom)%lmn2_size
414      wvl%paw%paw_ij(iatom)%ndij           =paw_ij(iatom)%ndij
415      wvl%paw%paw_ij(iatom)%nspden         =paw_ij(iatom)%nspden
416      wvl%paw%paw_ij(iatom)%nsppol         =paw_ij(iatom)%nsppol
417      if (paw_ij(iatom)%has_dij/=0) then
418        iaux=paw_ij(iatom)%cplex_dij*paw_ij(iatom)%lmn2_size
419        ABI_ALLOCATE(wvl%paw%paw_ij(iatom)%dij,(iaux,paw_ij(iatom)%ndij))
420        wvl%paw%paw_ij(iatom)%dij(:,:)=paw_ij(iatom)%dij(:,:)
421      end if
422    end do
423 
424 !  Option==2: deallocate
425  elseif(option==2) then
426    do iatom=1,my_natom
427      wvl%paw%paw_ij(iatom)%has_dij=0
428      if (associated(wvl%paw%paw_ij(iatom)%dij)) then
429        ABI_DEALLOCATE(wvl%paw%paw_ij(iatom)%dij)
430      end if
431    end do
432    ABI_DATATYPE_DEALLOCATE(wvl%paw%paw_ij)
433 
434 !  Option==3: only copy
435  elseif(option==3) then
436    do iatom=1,my_natom
437      wvl%paw%paw_ij(iatom)%cplex     =paw_ij(iatom)%cplex_rf
438      wvl%paw%paw_ij(iatom)%cplex_dij =paw_ij(iatom)%cplex_dij
439      wvl%paw%paw_ij(iatom)%lmn_size  =paw_ij(iatom)%lmn_size
440      wvl%paw%paw_ij(iatom)%lmn2_size =paw_ij(iatom)%lmn2_size
441      wvl%paw%paw_ij(iatom)%ndij      =paw_ij(iatom)%ndij
442      wvl%paw%paw_ij(iatom)%nspden    =paw_ij(iatom)%nspden
443      wvl%paw%paw_ij(iatom)%nsppol    =paw_ij(iatom)%nsppol
444      wvl%paw%paw_ij(iatom)%dij(:,:)  =paw_ij(iatom)%dij(:,:)
445    end do
446 
447  else
448    message = 'paw2wvl_ij: option should be equal to 1, 2 or 3'
449    MSG_ERROR(message)
450  end if
451 
452 #else
453  if (.false.) write(std_out,*) option,wvl%h(1),paw_ij(1)%ndij
454 #endif
455 
456  DBG_EXIT("COLL")
457 
458 end subroutine paw2wvl_ij

ABINIT/wvl_cprjreorder [ Functions ]

[ Top ] [ Functions ]

NAME

  wvl_cprjreorder

FUNCTION

 Change the order of a wvl-cprj datastructure
   From unsorted cprj to atom-sorted cprj (atm_indx=atindx)
   From atom-sorted cprj to unsorted cprj (atm_indx=atindx1)

INPUTS

  atm_indx(natom)=index table for atoms
   From unsorted wvl%paw%cprj to atom-sorted wvl%paw%cprj (atm_indx=atindx)
   From atom-sorted wvl%paw%cprj to unsorted wvl%paw%cprj (atm_indx=atindx1)

OUTPUT

PARENTS

      scfcv

CHILDREN

      cprj_clean,cprj_paw_alloc

SOURCE

569 subroutine wvl_cprjreorder(wvl,atm_indx)
570 
571 #if defined HAVE_BIGDFT
572  use BigDFT_API,only : cprj_objects,cprj_paw_alloc,cprj_clean
573  use dynamic_memory
574 #endif
575 
576 !This section has been created automatically by the script Abilint (TD).
577 !Do not modify the following lines by hand.
578 #undef ABI_FUNC
579 #define ABI_FUNC 'wvl_cprjreorder'
580 !End of the abilint section
581 
582  implicit none
583 
584 !Arguments ------------------------------------
585 !scalars
586 !arrays
587  integer,intent(in) :: atm_indx(:)
588  type(wvl_internal_type),intent(inout),target :: wvl
589 
590 !Local variables-------------------------------
591 #if defined HAVE_BIGDFT
592 !scalars
593  integer :: iexit,ii,jj,kk,n1atindx,n1cprj,n2cprj,ncpgr
594  character(len=100) :: msg
595 !arrays
596  integer,allocatable :: nlmn(:)
597  type(cprj_objects),pointer :: cprj(:,:)
598  type(cprj_objects),allocatable :: cprj_tmp(:,:)
599 #endif
600 
601 ! *************************************************************************
602 
603  DBG_ENTER("COLL")
604 
605 #if defined HAVE_BIGDFT
606  cprj => wvl%paw%cprj
607 
608  n1cprj=size(cprj,dim=1);n2cprj=size(cprj,dim=2)
609  n1atindx=size(atm_indx,dim=1)
610  if (n1cprj==0.or.n2cprj==0.or.n1atindx<=1) return
611  if (n1cprj/=n1atindx) then
612    msg='wrong sizes!'
613    MSG_BUG(msg)
614  end if
615 
616 !Nothing to do when the atoms are already sorted
617  iexit=1;ii=0
618  do while (iexit==1.and.ii<n1atindx)
619    ii=ii+1
620    if (atm_indx(ii)/=ii) iexit=0
621  end do
622  if (iexit==1) return
623 
624  ABI_ALLOCATE(nlmn,(n1cprj))
625  do ii=1,n1cprj
626    nlmn(ii)=cprj(ii,1)%nlmn
627  end do
628  ncpgr=cprj(1,1)%ncpgr
629 
630  ABI_DATATYPE_ALLOCATE(cprj_tmp,(n1cprj,n2cprj))
631  call cprj_paw_alloc(cprj_tmp,ncpgr,nlmn)
632  do jj=1,n2cprj
633    do ii=1,n1cprj
634      cprj_tmp(ii,jj)%nlmn=nlmn(ii)
635      cprj_tmp(ii,jj)%ncpgr=ncpgr
636      cprj_tmp(ii,jj)%cp(:,:)=cprj(ii,jj)%cp(:,:)
637      if (ncpgr>0) cprj_tmp(ii,jj)%dcp(:,:,:)=cprj(ii,jj)%dcp(:,:,:)
638    end do
639  end do
640 
641  call cprj_clean(cprj)
642 
643  do jj=1,n2cprj
644    do ii=1,n1cprj
645      kk=atm_indx(ii)
646      cprj(kk,jj)%nlmn=nlmn(ii)
647      cprj(kk,jj)%ncpgr=ncpgr
648      cprj(kk,jj)%cp=f_malloc_ptr((/2,nlmn(ii)/),id='cprj%cp')
649      cprj(kk,jj)%cp(:,:)=cprj_tmp(ii,jj)%cp(:,:)
650      if (ncpgr>0) then
651        cprj(kk,jj)%dcp=f_malloc_ptr((/2,ncpgr,nlmn(ii)/),id='cprj%dcp')
652        cprj(kk,jj)%dcp(:,:,:)=cprj_tmp(kk,jj)%dcp(:,:,:)
653      end if
654    end do
655  end do
656 
657  call cprj_clean(cprj_tmp)
658  ABI_DATATYPE_DEALLOCATE(cprj_tmp)
659  ABI_DEALLOCATE(nlmn)
660 
661 #else
662  if (.false.) write(std_out,*) atm_indx(1),wvl%h(1)
663 #endif
664 
665  DBG_EXIT("COLL")
666 
667 end subroutine wvl_cprjreorder

ABINIT/wvl_paw_free [ Functions ]

[ Top ] [ Functions ]

NAME

  wvl_paw_free

FUNCTION

  Frees memory for WVL+PAW implementation

INPUTS

  ntypat = number of atom types
  wvl= wvl type
  wvl_proj= wvl projector type

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

      gstate

CHILDREN

SOURCE

487 subroutine wvl_paw_free(wvl)
488 
489 
490 !This section has been created automatically by the script Abilint (TD).
491 !Do not modify the following lines by hand.
492 #undef ABI_FUNC
493 #define ABI_FUNC 'wvl_paw_free'
494 !End of the abilint section
495 
496  implicit none
497 
498 !Arguments ------------------------------------
499  type(wvl_internal_type),intent(inout) :: wvl
500 
501 !Local variables-------------------------------
502 
503 ! *************************************************************************
504 
505 #if defined HAVE_BIGDFT
506 
507 !PAW objects
508  if( associated(wvl%paw%spsi)) then
509    ABI_DEALLOCATE(wvl%paw%spsi)
510  end if
511  if( associated(wvl%paw%indlmn)) then
512    ABI_DEALLOCATE(wvl%paw%indlmn)
513  end if
514  if( associated(wvl%paw%sij)) then
515    ABI_DEALLOCATE(wvl%paw%sij)
516  end if
517  if( associated(wvl%paw%rpaw)) then
518    ABI_DEALLOCATE(wvl%paw%rpaw)
519  end if
520 
521 !rholoc
522  if( associated(wvl%rholoc%msz )) then
523    ABI_DEALLOCATE(wvl%rholoc%msz)
524  end if
525  if( associated(wvl%rholoc%d )) then
526    ABI_DEALLOCATE(wvl%rholoc%d)
527  end if
528  if( associated(wvl%rholoc%rad)) then
529    ABI_DEALLOCATE(wvl%rholoc%rad)
530  end if
531  if( associated(wvl%rholoc%radius)) then
532    ABI_DEALLOCATE(wvl%rholoc%radius)
533  end if
534 
535 #else
536  if (.false.) write(std_out,*) wvl%h(1)
537 #endif
538 
539 !paw%paw_ij and paw%cprj are allocated and deallocated inside vtorho
540 
541 end subroutine wvl_paw_free