TABLE OF CONTENTS


ABINIT/m_paw_lmn [ Modules ]

[ Top ] [ Modules ]

NAME

  m_paw_lmn

FUNCTION

  This module provides tools to calculate tables commonly used to iterate
  over the the (l,m,n) channels of the PAW partial waves.

COPYRIGHT

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

SOURCE

17 #if defined HAVE_CONFIG_H
18 #include "config.h"
19 #endif
20 
21 #include "abi_common.h"
22 
23 MODULE m_paw_lmn
24 
25  use defs_basis
26  use m_abicore
27  use m_errors
28 
29  implicit none
30 
31  private
32 
33 ! Public procedures.
34  public :: ilm2lm          ! Returns the value of l and m in from the index ilm
35  public :: make_indlmn     ! Calculates the indlmn(6,lmn_size) table giving l,m,n,lm,ln,spin for i=lmn.
36  public :: make_indklmn    ! Calculates the indklmn(8,lmn2_size) table giving
37                            !   klm, kln, abs(il-jl), (il+jl), ilm and jlm, ilmn and jlmn for each symmetric klmn=(ilmn,jlmn)
38  public :: make_kln2ln     ! Calculates the kln2ln(6,ln2_size) table giving
39                            !   il, jl ,in, jn, iln, jln for each symmetric kln=(iln,jln)
40  public :: make_klm2lm     ! Calculates the klm2lm(6,lm2_size) table giving
41                            !   il, jl ,im, jm, ilm, jlm for each symmetric klm=(ilm,jlm)
42  public :: klmn2ijlmn      ! Calculates ilmn and jlmn from klmn.
43  public :: make_indln      ! Calculates indln(2,ln_size) giving l and n for i=ln
44  public :: uppert_index    ! The sequential index of an element in the upper triangle of a matrix 
45 
46 CONTAINS  !========================================================================================

m_paw_lmn/ilm2lm [ Functions ]

[ Top ] [ m_paw_lmn ] [ Functions ]

NAME

  ilm2lm

FUNCTION

  Returns the value of l and m in from the index ilm

INPUTS

  ilm=The contracted index for (l,m).

OUTPUT

  ll=The angular momentun defined in [0,1,2,3,4].
  mm=The magnetic number in the interval [-l,....+l].

SOURCE

65 subroutine ilm2lm(ilm,ll,mm)
66 
67  implicit none
68 
69 !Arguments ------------------------------------
70 !scalars
71  integer,intent(in) :: ilm
72  integer,intent(out) :: ll,mm
73 
74 !Local variables-------------------------------
75  integer :: ii
76 ! *********************************************************************
77 
78  if (ilm<1) then 
79    ABI_ERROR("Wrong ilm")
80  end if
81 
82  ll = -1
83  do ii=0,100
84   if ( (ii+1)**2 >= ilm) then
85    ll = ii
86    EXIT
87   end if
88  end do
89 
90  mm = ilm - ll**2 -ll-1
91 
92  if (ll==-1) then
93    ABI_ERROR("l>100 not programmed!")
94  end if
95 
96 end subroutine ilm2lm

m_paw_lmn/klmn2ijlmn [ Functions ]

[ Top ] [ m_paw_lmn ] [ Functions ]

NAME

  klmn2ijlmn

FUNCTION

  Find ilmn and jlmn from klmn and lmn_size.

INPUTS

  lmn_size=Number of (l,m,n) elements for the PAW basis set.
  klmn=The index corresponding to (ilmn,jlmn) in packed form.

OUTPUT

  jlmn, ilmn=The two symmetrix indices corresponding to klmn. NB: jlmn >= ilmn

SOURCE

430 subroutine klmn2ijlmn(klmn,lmn_size,ilmn,jlmn)
431 
432  implicit none
433 
434 !Arguments ------------------------------------
435 !scalars
436  integer,intent(in) :: klmn,lmn_size
437  integer,intent(out) :: ilmn,jlmn
438 
439 !Local variables ------------------------------
440 !scalars
441  integer :: ii,jj,k0
442 !************************************************************************
443 
444  ilmn=-1; jlmn=-1
445 
446  do jj=1,lmn_size
447    k0=jj*(jj-1)/2
448    do ii=1,jj
449      if (klmn==ii+k0) then
450        ilmn = ii
451        jlmn = jj; RETURN
452      end if
453    end do
454  end do
455 
456  ABI_BUG("Not able to found ilmn and jlmn")
457 
458 end subroutine klmn2ijlmn

m_paw_lmn/make_indklmn [ Functions ]

[ Top ] [ m_paw_lmn ] [ Functions ]

NAME

  make_indklmn

FUNCTION

  Performs the setup of the indklmn table for PAW calculations.
  Compute the indklmn indexes giving klm, kln, abs(il-jl) and (il+jl), ilm and jlm, ilmn and jlmn
  for each klmn=(ilmn,jlmn) with jlmn >= ilmn

INPUTS

  lcutdens=Maximum l for densities/potentials moments computations
  lmn_size=Number of (l,m,n) elements for the PAW basis set.
  lmn2_size=Number of elements in the symmetric basis set: lmn2_size=lmn_size*(lmn_size+1)/2
  indlmn(6,lmn_size)=Array giving l,m,n,lm,ln,spin for i=lmn.

OUTPUT

  indklmn(6,lmn2_size)=Array giving klm, kln, abs(il-jl), (il+jl), ilm and jlm, ilmn and jlmn
    for each klmn=(ilmn,jlmn). Note: ilmn=(il,im,in) and ilmn<=jlmn
  klm_diag(lmn2_size)=1 il==jl and im==jm, 0 otherwise.

SOURCE

184 subroutine make_indklmn(lcutdens,lmn_size,lmn2_size,indlmn,indklmn,klm_diag)
185 
186  implicit none
187 
188 !Arguments ------------------------------------
189 !scalars
190  integer,intent(in) :: lmn2_size,lmn_size
191  integer,intent(in) ::  lcutdens
192 !scalars
193  integer,intent(in) :: indlmn(6,lmn_size)
194  integer,intent(out) :: indklmn(8,lmn2_size)
195  integer,intent(out) :: klm_diag(lmn2_size)
196 
197 !Local variables ------------------------------
198 !scalars
199  integer :: i0lm,i0ln,il,ilm,ilmn,iln
200  integer :: j0lm,j0lmn,j0ln,jl,jlm,jlmn,jln,klmn
201 
202 !************************************************************************
203 
204  klm_diag=0
205 
206  do jlmn=1,lmn_size
207   jl= indlmn(1,jlmn); jlm=indlmn(4,jlmn); jln=indlmn(5,jlmn)
208   j0lmn=jlmn*(jlmn-1)/2
209   j0lm =jlm *(jlm -1)/2
210   j0ln =jln *(jln -1)/2
211   do ilmn=1,jlmn
212    il=indlmn(1,ilmn); ilm=indlmn(4,ilmn); iln=indlmn(5,ilmn)
213    klmn=j0lmn+ilmn
214    if (ilm<=jlm) then
215     indklmn(1,klmn)=j0lm+ilm     !klm
216    else
217     i0lm=ilm*(ilm-1)/2
218     indklmn(1,klmn)=i0lm+jlm
219    end if
220    if (iln<=jln) then
221     indklmn(2,klmn)=j0ln+iln     !kln
222    else
223     i0ln=iln*(iln-1)/2
224     indklmn(2,klmn)=i0ln+jln
225    end if
226    !MG This is not safe, what happens if lcutdens < |il-jl|?
227    indklmn(3,klmn)=MIN(ABS(il-jl),lcutdens)  ! abs(li-lj) NB this is a l-value, not an index >=1.
228    indklmn(4,klmn)=MIN(il+jl,lcutdens)       ! abs(li+lj) NB this is a l-value, not an index >=1.
229 
230    indklmn(5,klmn)=ilm                       ! ilm
231    indklmn(6,klmn)=jlm                       ! jlm
232 
233    indklmn(7,klmn)=ilmn                      ! ilmn
234    indklmn(8,klmn)=jlmn                      ! jlmn
235 
236 
237    if (ilm==jlm) klm_diag(klmn)=1
238   end do
239  end do
240 
241 end subroutine make_indklmn

m_paw_lmn/make_indlmn [ Functions ]

[ Top ] [ m_paw_lmn ] [ Functions ]

NAME

  make_indlmn

FUNCTION

  Performs the setup of the indlmn table for PAW calculations (indices for (l,m,n) basis)

INPUTS

  ln_size= Total number of nl components
  lmn_size= Second dimension in indlmn. Total number of (l,m,n) components for this atom.
  orbitals(ln_size)=Give the value of l for each element of the augmented basis set.

OUTPUT

  indlmn(6,lmn_size)=array giving l,m,n,lm,ln,s for i=lmn

SOURCE

118 subroutine make_indlmn(ln_size,lmn_size,orbitals,indlmn)
119 
120  implicit none
121 
122 !Arguments ------------------------------------
123 !scalars
124  integer,intent(in) :: ln_size,lmn_size
125  integer,intent(in) :: orbitals(ln_size)
126 !scalars
127  integer,intent(out) :: indlmn(6,lmn_size)
128 
129 !Local variables ------------------------------
130 !scalars
131  integer :: ilmn,ib,il,iln,ilm
132 !arrays
133  integer,allocatable :: nprj(:)
134 
135 !************************************************************************
136 
137  ABI_MALLOC(nprj,(0:MAXVAL(orbitals)))
138 
139  ilmn=0; iln=0; nprj=0
140  do ib=1,ln_size
141    il=orbitals(ib)
142    nprj(il)=nprj(il)+1
143    iln=iln+1
144    do ilm=1,2*il+1
145      indlmn(1,ilmn+ilm)=il           ! l
146      indlmn(2,ilmn+ilm)=ilm-(il+1)   ! m
147      indlmn(3,ilmn+ilm)=nprj(il)     ! n
148      indlmn(4,ilmn+ilm)=il*il+ilm    ! lm index
149      indlmn(5,ilmn+ilm)=iln          ! ln index
150      indlmn(6,ilmn+ilm)=1            ! spin (not yet used!)
151    end do
152    ilmn=ilmn+2*il+1
153  end do
154 
155  ABI_FREE(nprj)
156 
157 end subroutine make_indlmn

m_paw_lmn/make_indln [ Functions ]

[ Top ] [ m_paw_lmn ] [ Functions ]

NAME

  make_indln

FUNCTION

  Performs the setup of the indln table for PAW calculations.
  Compute the indln indexes giving ilmn=(l,m,n)

INPUTS

  lmn_size=Number of (l,m,n) elements for the PAW basis set.
  ln_size=Number of (l,n) elements
  indlmn(6,lmn_size)=Array giving l,m,n,lm,ln,spin for i=lmn.

OUTPUT

  indln(2,ln_size)=Array giving l and n for i=ln

SOURCE

481 subroutine make_indln(lmn_size,ln_size,indlmn,indln)
482 
483  implicit none
484 
485 !Arguments ------------------------------------
486 !scalars
487  integer,intent(in) :: lmn_size,ln_size
488 !arrays
489  integer,intent(in) :: indlmn(6,lmn_size)
490  integer,intent(out) :: indln(2,ln_size)
491 
492 !Local variables ------------------------------
493 !scalars
494  integer :: ilmn,ll,nn,ll_old,nn_old,ii
495 !************************************************************************
496 
497  ii=0; ll_old=-1; nn_old=-1
498  do ilmn=1,lmn_size
499   ll=indlmn(1,ilmn)
500   nn=indlmn(3,ilmn)
501   if (ll/=ll_old.or.nn/=nn_old) then
502    ll_old=ll
503    nn_old=nn
504    ii=ii+1
505    indln(1,ii)=ll
506    indln(2,ii)=nn
507   end if
508  end do
509 
510  ABI_CHECK(ii==ln_size,"ii/=ln_size")
511 
512 end subroutine make_indln

m_paw_lmn/make_klm2lm [ Functions ]

[ Top ] [ m_paw_lmn ] [ Functions ]

NAME

  make_klm2lm

FUNCTION

  Performs the setup of the klm2lm table for PAW calculations.

INPUTS

  lmn_size=Number of (l,m,n) elements for the PAW basis set.
  lmn2_size=Number of elements in the symmetric basis set: lmn2_size=lmn_size*(lmn_size+1)/2
  lm2_size)=Number of (l.m) elements in the symmetric basis set.
  indlmn(6,lmn_size)=Array giving l,m,n,lm,ln,spin for i=lmn.
  indklmn(8,lmn2_size)=Array giving klm, kln, abs(il-jl), (il+jl), ilm and jlm
   for each klmn=(ilmn,jlmn). Note: ilmn=(il,im,in) and ilmn<=jlmn

OUTPUT

  klm2lm(6,lm2_size)=Table giving il, jl ,im, jm, ilm, jlm for each klm=(ilm,jlm)
  where ilm=(il,im) and ilm<=jlm. NB: klm2lm is an application and not a bijection.

NOTES

  klm2lm can be calculated easily if we assume that all (l,m) channels
  are ordered by increasing l and m. This is the standard convention
  used in most of the PAW datasets. This routines, howevever, works
  works also in the unlikely case in with (l,m) are not ordered.

SOURCE

347 subroutine make_klm2lm(lmn_size,lmn2_size,lm2_size,indlmn,indklmn,klm2lm)
348 
349  implicit none
350 
351 !Arguments ------------------------------------
352 !scalars
353  integer,intent(in) :: lmn2_size,lmn_size,lm2_size
354 !arrays
355  integer,intent(in) :: indlmn(6,lmn_size)
356  integer,intent(in) :: indklmn(8,lmn2_size)
357  integer,intent(out) :: klm2lm(6,lm2_size)
358 
359 !Local variables ------------------------------
360 !scalars
361  integer :: il,ilm,ilmn,jl,jlm,jlmn,im,jm
362  integer :: klmn,klm,klm_old !,iklm
363 
364 !************************************************************************
365 
366  klm2lm = 0
367  klm_old = -1
368  do klmn=1,lmn2_size
369    klm = indklmn(1,klmn)
370    if (klm /= klm_old) then
371      klm_old = klm
372 
373      ilm = indklmn(5,klmn)
374      jlm = indklmn(6,klmn)
375 
376      call klmn2ijlmn(klmn,lmn_size,ilmn,jlmn)
377 
378      il=indlmn(1,ilmn); im=indlmn(2,ilmn)
379      jl=indlmn(1,jlmn); jm=indlmn(2,jlmn)
380 
381      !shift to have the index instead of l or m.
382      klm2lm(1,klm) = il +1       ! il
383      klm2lm(2,klm) = jl +1       ! jl
384      klm2lm(3,klm) = im + il +1  ! im
385      klm2lm(4,klm) = jm + jl +1  ! jm
386      klm2lm(5,klm) = ilm         ! ilm
387      klm2lm(6,klm) = jlm         ! jlm
388    end if
389  end do !klmn
390 
391  if (ANY(klm2lm==0)) then
392    ABI_BUG("check klm2lm")
393  end if
394 
395 !DEBUG
396 #if 0
397  write(std_out,*)"Debugging make_klm2lm:"
398  do iklm=1,lm2_size
399   il  = klm2lm(1,iklm) -1      ! li
400   jl  = klm2lm(2,iklm) -1      ! lj
401   im  = klm2lm(3,iklm) -il -1  ! mi
402   jm  = klm2lm(4,iklm) -jl -1  ! mj
403   ilm = klm2lm(5,iklm)         ! ilm
404   jlm = klm2lm(6,iklm)         ! jlm
405   write(std_out,'(i3,2(a,2i3,a))')"iklm ",iklm," l m (",il,im,")"," l m (",jl,jm,")"
406  end do
407 #endif
408 
409 end subroutine make_klm2lm

m_paw_lmn/make_kln2ln [ Functions ]

[ Top ] [ m_paw_lmn ] [ Functions ]

NAME

  make_kln2ln

FUNCTION

  Performs the setup of the kln2ln table for PAW calculations.

INPUTS

  lmn_size=Number of (l,m,n) elements for the PAW basis set.
  lmn2_size=Number of elements in the symmetric basis set: lmn2_size=lmn_size*(lmn_size+1)/2
  ln2_size=Number of symmetric (l,n) channels i.e. ln_size*(ln_size+1)/2
  indlmn(6,lmn_size)=Array giving l,m,n,lm,ln,spin for i=lmn.
  indklmn(8,lmn2_size)=Array giving klm, kln, abs(il-jl), (il+jl), ilm and jlm, ilmn and jlmn
   for each klmn=(ilmn,jlmn). Note: ilmn=(il,im,in) and ilmn<=jlmn

OUTPUT

  kln2ln(6,ln2_size)=Table giving il, jl ,in, jn, iln, jln for each kln=(iln,jln)
  where iln=(il,in) and iln<=jln. NB: kln2ln is an application and not a bijection

SOURCE

267 subroutine make_kln2ln(lmn_size,lmn2_size,ln2_size,indlmn,indklmn,kln2ln)
268 
269  implicit none
270 
271 !Arguments ------------------------------------
272 !scalars
273  integer,intent(in) :: lmn2_size,lmn_size,ln2_size
274 !arrays
275  integer,intent(in) :: indlmn(6,lmn_size)
276  integer,intent(in) :: indklmn(8,lmn2_size)
277  integer,intent(out) :: kln2ln(6,ln2_size)
278 
279 !Local variables ------------------------------
280 !scalars
281  integer :: il,in,ilmn,iln
282  integer :: jl,jn,jlmn,jln
283  integer :: klmn,kln,kln_old
284 
285 !************************************************************************
286 
287  kln2ln = 0
288 
289  kln_old = -1
290  do klmn=1,lmn2_size
291    kln = indklmn(2,klmn)
292    if (kln /= kln_old) then
293      kln_old = kln
294 
295      call klmn2ijlmn(klmn,lmn_size,ilmn,jlmn)
296 
297      il= indlmn(1,ilmn); iln=indlmn(5,ilmn)
298      jl= indlmn(1,jlmn); jln=indlmn(5,jlmn)
299      !$in = indklmn(9 ,klmn) !=indlmn(3,ilmn)
300      !$jn = indklmn(10,klmn) !=indlmn(3,jlmn)
301 
302      in = indlmn(3,ilmn)
303      jn = indlmn(3,jlmn)
304 
305      ! Shift il to have the index instead of the value of l.
306      kln2ln(1,kln) = il +1
307      kln2ln(2,kln) = jl +1
308      kln2ln(3,kln) = in
309      kln2ln(4,kln) = jn
310      kln2ln(5,kln) = iln
311      kln2ln(6,kln) = jln
312    end if
313  end do !klmn
314 
315 end subroutine make_kln2ln

m_paw_lmn/uppert_index [ Functions ]

[ Top ] [ m_paw_lmn ] [ Functions ]

NAME

   uppert_index

FUNCTION

  Helper function returning the sequential index of an element in the upper triangle of a matrix 
  given the row-index ii and the column-index jj. If ii>jj the index of the element a_{jj,ii} is returned.

INPUTS

  ii=Row index
  jj=column index

SOURCE

531 function uppert_index(ii,jj)
532 
533  implicit none
534 
535 !Arguments ------------------------------------
536 !scalars
537  integer,intent(in) :: ii,jj
538  integer :: uppert_index
539 !scalars
540 
541 !************************************************************************
542 
543  if (jj>=jj) then
544    uppert_index = ii + jj*(jj-1)/2
545  else
546    uppert_index = jj + ii*(ii-1)/2
547  end if
548 
549 end function uppert_index