LCOV - code coverage report
Current view: top level - shared/libpaw/src - m_paw_lmn.F90 (source / functions) Coverage Total Hit
Test: coverage.info Lines: 88.1 % 151 133
Test Date: 2026-09-19 17:42:43 Functions: 75.0 % 8 6

            Line data    Source code
       1              : !!****m* ABINIT/m_paw_lmn
       2              : !! NAME
       3              : !!  m_paw_lmn
       4              : !!
       5              : !! FUNCTION
       6              : !!  This module provides tools to calculate tables commonly used to iterate
       7              : !!  over the the (l,m,n) channels of the PAW partial waves.
       8              : !!
       9              : !! COPYRIGHT
      10              : !! Copyright (C) 2008-2026 ABINIT group (MG)
      11              : !! This file is distributed under the terms of the
      12              : !! GNU General Public License, see ~abinit/COPYING
      13              : !! or http://www.gnu.org/copyleft/gpl.txt .
      14              : !!
      15              : !! SOURCE
      16              : 
      17              : #include "libpaw.h"
      18              : 
      19              : MODULE m_paw_lmn
      20              : 
      21              :  USE_DEFS
      22              :  USE_MSG_HANDLING
      23              :  USE_MPI_WRAPPERS
      24              :  USE_MEMORY_PROFILING
      25              : 
      26              :  implicit none
      27              : 
      28              :  private
      29              : 
      30              : ! Public procedures.
      31              :  public :: ilm2lm          ! Returns the value of l and m in from the index ilm
      32              :  public :: make_indlmn     ! Calculates the indlmn(6,lmn_size) table giving l,m,n,lm,ln,spin for i=lmn.
      33              :  public :: make_indklmn    ! Calculates the indklmn(8,lmn2_size) table giving
      34              :                            !   klm, kln, abs(il-jl), (il+jl), ilm and jlm, ilmn and jlmn for each symmetric klmn=(ilmn,jlmn)
      35              :  public :: make_kln2ln     ! Calculates the kln2ln(6,ln2_size) table giving
      36              :                            !   il, jl ,in, jn, iln, jln for each symmetric kln=(iln,jln)
      37              :  public :: make_klm2lm     ! Calculates the klm2lm(6,lm2_size) table giving
      38              :                            !   il, jl ,im, jm, ilm, jlm for each symmetric klm=(ilm,jlm)
      39              :  public :: klmn2ijlmn      ! Calculates ilmn and jlmn from klmn.
      40              :  public :: make_indln      ! Calculates indln(2,ln_size) giving l and n for i=ln
      41              :  public :: uppert_index    ! The sequential index of an element in the upper triangle of a matrix
      42              : 
      43              : CONTAINS  !========================================================================================
      44              : !!***
      45              : 
      46              : !!****f* m_paw_lmn/ilm2lm
      47              : !! NAME
      48              : !!  ilm2lm
      49              : !!
      50              : !! FUNCTION
      51              : !!  Returns the value of l and m in from the index ilm
      52              : !!
      53              : !! INPUTS
      54              : !!  ilm=The contracted index for (l,m).
      55              : !!
      56              : !! OUTPUT
      57              : !!  ll=The angular momentun defined in [0,1,2,3,4].
      58              : !!  mm=The magnetic number in the interval [-l,....+l].
      59              : !!
      60              : !! SOURCE
      61              : 
      62            0 : subroutine ilm2lm(ilm,ll,mm)
      63              : 
      64              : !Arguments ------------------------------------
      65              : !scalars
      66              :  integer,intent(in) :: ilm
      67              :  integer,intent(out) :: ll,mm
      68              : 
      69              : !Local variables-------------------------------
      70              :  integer :: ii
      71              : ! *********************************************************************
      72              : 
      73            0 :  if (ilm<1) then
      74            0 :    LIBPAW_ERROR("Wrong ilm")
      75              :  end if
      76              : 
      77            0 :  ll = -1
      78            0 :  do ii=0,100
      79            0 :   if ( (ii+1)**2 >= ilm) then
      80            0 :    ll = ii
      81            0 :    EXIT
      82              :   end if
      83              :  end do
      84              : 
      85            0 :  mm = ilm - ll**2 -ll-1
      86              : 
      87            0 :  if (ll==-1) then
      88            0 :    LIBPAW_ERROR("l>100 not programmed!")
      89              :  end if
      90              : 
      91            0 : end subroutine ilm2lm
      92              : !!***
      93              : 
      94              : !----------------------------------------------------------------------
      95              : 
      96              : !!****f* m_paw_lmn/make_indlmn
      97              : !! NAME
      98              : !!  make_indlmn
      99              : !!
     100              : !! FUNCTION
     101              : !!  Performs the setup of the indlmn table for PAW calculations (indices for (l,m,n) basis)
     102              : !!
     103              : !! INPUTS
     104              : !!  ln_size= Total number of nl components
     105              : !!  lmn_size= Second dimension in indlmn. Total number of (l,m,n) components for this atom.
     106              : !!  orbitals(ln_size)=Give the value of l for each element of the augmented basis set.
     107              : !!
     108              : !! OUTPUT
     109              : !!  indlmn(6,lmn_size)=array giving l,m,n,lm,ln,s for i=lmn
     110              : !!  indlmn(8,lmn_size)=array giving l,m,sign of kappa,ln,spinor,2*j,2*m_j in case of dirac relativism
     111              : !!
     112              : !! SOURCE
     113              : 
     114           10 : subroutine make_indlmn(ln_size,lmn_size,orbitals,indlmn,kappa)
     115              : 
     116              : !Arguments ------------------------------------
     117              : !scalars
     118              :  integer,intent(in) :: ln_size,lmn_size
     119              :  integer,intent(in) :: orbitals(ln_size)
     120              : !scalars
     121              :  integer,allocatable,intent(inout) :: indlmn(:,:)
     122              :  integer, intent(in), optional :: kappa(ln_size)
     123              : 
     124              : !Local variables ------------------------------
     125              : !scalars
     126              :  integer :: ilmn,ib,il,iln,ilm,kappa_sign,spinor,i2j,i2mj,im,ilmn_ws
     127              : !arrays
     128           10 :  integer,allocatable :: nprj(:)
     129              : 
     130              : !************************************************************************
     131              : 
     132           10 :  if(.not.present(kappa)) then
     133           46 :    LIBPAW_ALLOCATE(nprj,(0:MAXVAL(orbitals)))
     134           24 :    LIBPAW_ALLOCATE(indlmn,(6,lmn_size))
     135          260 :    indlmn=0
     136           23 :    ilmn=0; iln=0; nprj=0
     137           30 :    do ib=1,ln_size
     138           22 :      il=orbitals(ib)
     139           22 :      nprj(il)=nprj(il)+1
     140           22 :      iln=iln+1
     141           58 :      do ilm=1,2*il+1
     142           36 :        indlmn(1,ilmn+ilm)=il           ! l
     143           36 :        indlmn(2,ilmn+ilm)=ilm-(il+1)   ! m
     144           36 :        indlmn(3,ilmn+ilm)=nprj(il)     ! n
     145           36 :        indlmn(4,ilmn+ilm)=il*il+ilm    ! lm index
     146           36 :        indlmn(5,ilmn+ilm)=iln          ! ln index
     147           58 :        indlmn(6,ilmn+ilm)=1            ! spin (not yet used!)
     148              :      end do
     149           30 :      ilmn=ilmn+2*il+1
     150              :    end do
     151            8 :    LIBPAW_DEALLOCATE(nprj)
     152              :  else
     153            6 :    LIBPAW_ALLOCATE(indlmn,(9,lmn_size))
     154          402 :    indlmn=0;ilmn=0;iln=0;ilmn_ws=0
     155           18 :    do ib=1,2*ln_size
     156           16 :      iln=iln+modulo(ib,2)
     157           16 :      il=orbitals(iln)
     158           16 :      kappa_sign=sign(1,kappa(iln)) ! sgn(kappa)=+1 or -1
     159           16 :      spinor=2-modulo(ib,2)       ! spinor= 1 or 2
     160           16 :      i2j=2*il-kappa_sign              ! j=l-sgn(kappa)/2 = l-1/2 or l+1/2
     161           56 :      do ilm=1,i2j+1
     162           40 :        ilmn_ws=ilmn_ws+1
     163              :        !mj= -j,...,j
     164           40 :        i2mj=-i2j+2*(ilm-1)       ! 2m_j= -jc ... +jc
     165           40 :        im=(i2mj-3+2*spinor)/2    ! m=m_j-1/2 (spinor=1) or m_j+1/2 (spinor=2)
     166           56 :        if(abs(im)<=il) then
     167              :          !Valid value for sph. harm., i.e. abs(m)<=l
     168           28 :          indlmn(1,ilmn+ilm)=il !l
     169           28 :          indlmn(2,ilmn+ilm)=im !m
     170           28 :          indlmn(3,ilmn+ilm)=kappa_sign !sign of kappa
     171           28 :          indlmn(4,ilmn+ilm)=il*il+im+il+1 !lm
     172           28 :          indlmn(5,ilmn+ilm)=iln !ln also includes the two spinor values here
     173           28 :          indlmn(6,ilmn+ilm)=spinor !spinor index (1 up, 2 down)
     174           28 :          indlmn(7,ilmn+ilm)=i2j !2*j (times 2 to make it an integer)
     175           28 :          indlmn(8,ilmn+ilm)=i2mj !2*m_j (times 2 to make it an integer)
     176           28 :          indlmn(9,ilmn+ilm)=ilmn_ws ! ilmn without spinor
     177              :        else
     178              :          !Invalid value for sph. harm. ; will be multiplied by zero later
     179           12 :          indlmn(1,ilmn+ilm)=-1 !Invalid value that should be checked later
     180          108 :          indlmn(2:9,ilmn+ilm)=-1 ; indlmn(3,ilmn+ilm)=0
     181              :        endif
     182              :      end do
     183           16 :      ilmn=ilmn+i2j+1
     184           18 :      if(modulo(ib,2)==1) ilmn_ws=ilmn_ws-i2j-1
     185              :    end do
     186              :  endif
     187              : 
     188           10 : end subroutine make_indlmn
     189              : !!***
     190              : 
     191              : !----------------------------------------------------------------------
     192              : 
     193              : !!****f* m_paw_lmn/make_indklmn
     194              : !! NAME
     195              : !!  make_indklmn
     196              : !!
     197              : !! FUNCTION
     198              : !!  Performs the setup of the indklmn table for PAW calculations.
     199              : !!  Compute the indklmn indexes giving klm, kln, abs(il-jl) and (il+jl), ilm and jlm, ilmn and jlmn
     200              : !!  for each klmn=(ilmn,jlmn) with jlmn >= ilmn
     201              : !!
     202              : !! INPUTS
     203              : !!  lcutdens=Maximum l for densities/potentials moments computations
     204              : !!  lmn_size=Number of (l,m,n) elements for the PAW basis set.
     205              : !!  lmn2_size=Number of elements in the symmetric basis set: lmn2_size=lmn_size*(lmn_size+1)/2
     206              : !!  indlmn(6,lmn_size)=Array giving l,m,n,lm,ln,spin for i=lmn.
     207              : !!
     208              : !! OUTPUT
     209              : !!  indklmn(6,lmn2_size)=Array giving klm, kln, abs(il-jl), (il+jl), ilm and jlm, ilmn and jlmn
     210              : !!    for each klmn=(ilmn,jlmn). Note: ilmn=(il,im,in) and ilmn<=jlmn
     211              : !!  klm_diag(lmn2_size)=1 il==jl and im==jm, 0 otherwise.
     212              : !!
     213              : !! SOURCE
     214              : 
     215           10 : subroutine make_indklmn(lcutdens,lmn_size,lmn2_size,indlmn,indklmn,klm_diag)
     216              : 
     217              : !Arguments ------------------------------------
     218              : !scalars
     219              :  integer,intent(in) :: lmn2_size,lmn_size
     220              :  integer,intent(in) ::  lcutdens
     221              : !scalars
     222              :  integer,intent(in) :: indlmn(6,lmn_size)
     223              :  integer,intent(out) :: indklmn(8,lmn2_size)
     224              :  integer,intent(out) :: klm_diag(lmn2_size)
     225              : 
     226              : !Local variables ------------------------------
     227              : !scalars
     228              :  integer :: i0lm,i0ln,il,ilm,ilmn,iln
     229              :  integer :: j0lm,j0lmn,j0ln,jl,jlm,jlmn,jln,klmn
     230              : 
     231              : !************************************************************************
     232              : 
     233          536 :  klm_diag=0
     234              : 
     235           86 :  do jlmn=1,lmn_size
     236           76 :   jl= indlmn(1,jlmn); jlm=indlmn(4,jlmn); jln=indlmn(5,jlmn)
     237           76 :   j0lmn=jlmn*(jlmn-1)/2
     238           76 :   j0lm =jlm *(jlm -1)/2
     239           76 :   j0ln =jln *(jln -1)/2
     240          612 :   do ilmn=1,jlmn
     241          526 :    il=indlmn(1,ilmn); ilm=indlmn(4,ilmn); iln=indlmn(5,ilmn)
     242          526 :    klmn=j0lmn+ilmn
     243          526 :    if (ilm<=jlm) then
     244          440 :     indklmn(1,klmn)=j0lm+ilm     !klm
     245              :    else
     246           86 :     i0lm=ilm*(ilm-1)/2
     247           86 :     indklmn(1,klmn)=i0lm+jlm
     248              :    end if
     249          526 :    if (iln<=jln) then
     250          406 :     indklmn(2,klmn)=j0ln+iln     !kln
     251              :    else
     252          120 :     i0ln=iln*(iln-1)/2
     253          120 :     indklmn(2,klmn)=i0ln+jln
     254              :    end if
     255              :    !MG This is not safe, what happens if lcutdens < |il-jl|?
     256          526 :    indklmn(3,klmn)=MIN(ABS(il-jl),lcutdens)  ! abs(li-lj) NB this is a l-value, not an index >=1.
     257          526 :    indklmn(4,klmn)=MIN(il+jl,lcutdens)       ! abs(li+lj) NB this is a l-value, not an index >=1.
     258              : 
     259          526 :    indklmn(5,klmn)=ilm                       ! ilm
     260          526 :    indklmn(6,klmn)=jlm                       ! jlm
     261              : 
     262          526 :    indklmn(7,klmn)=ilmn                      ! ilmn
     263          526 :    indklmn(8,klmn)=jlmn                      ! jlmn
     264              : 
     265              : 
     266          602 :    if (ilm==jlm) klm_diag(klmn)=1
     267              :   end do
     268              :  end do
     269              : 
     270           10 : end subroutine make_indklmn
     271              : !!***
     272              : 
     273              : !----------------------------------------------------------------------
     274              : 
     275              : !!****f* m_paw_lmn/make_kln2ln
     276              : !! NAME
     277              : !!  make_kln2ln
     278              : !!
     279              : !! FUNCTION
     280              : !!  Performs the setup of the kln2ln table for PAW calculations.
     281              : !!
     282              : !! INPUTS
     283              : !!  lmn_size=Number of (l,m,n) elements for the PAW basis set.
     284              : !!  lmn2_size=Number of elements in the symmetric basis set: lmn2_size=lmn_size*(lmn_size+1)/2
     285              : !!  ln2_size=Number of symmetric (l,n) channels i.e. ln_size*(ln_size+1)/2
     286              : !!  indlmn(6,lmn_size)=Array giving l,m,n,lm,ln,spin for i=lmn.
     287              : !!  indklmn(8,lmn2_size)=Array giving klm, kln, abs(il-jl), (il+jl), ilm and jlm, ilmn and jlmn
     288              : !!   for each klmn=(ilmn,jlmn). Note: ilmn=(il,im,in) and ilmn<=jlmn
     289              : !!
     290              : !! OUTPUT
     291              : !!  kln2ln(6,ln2_size)=Table giving il, jl ,in, jn, iln, jln for each kln=(iln,jln)
     292              : !!  where iln=(il,in) and iln<=jln. NB: kln2ln is an application and not a bijection
     293              : !!
     294              : !! SOURCE
     295              : 
     296            2 : subroutine make_kln2ln(lmn_size,lmn2_size,ln2_size,indlmn,indklmn,kln2ln)
     297              : 
     298              : !Arguments ------------------------------------
     299              : !scalars
     300              :  integer,intent(in) :: lmn2_size,lmn_size,ln2_size
     301              : !arrays
     302              :  integer,intent(in) :: indlmn(6,lmn_size)
     303              :  integer,intent(in) :: indklmn(8,lmn2_size)
     304              :  integer,intent(out) :: kln2ln(6,ln2_size)
     305              : 
     306              : !Local variables ------------------------------
     307              : !scalars
     308              :  integer :: il,in,ilmn,iln
     309              :  integer :: jl,jn,jlmn,jln
     310              :  integer :: klmn,kln,kln_old
     311              : 
     312              : !************************************************************************
     313              : 
     314          142 :  kln2ln = 0
     315              : 
     316            2 :  kln_old = -1
     317           74 :  do klmn=1,lmn2_size
     318           72 :    kln = indklmn(2,klmn)
     319           74 :    if (kln /= kln_old) then
     320           48 :      kln_old = kln
     321              : 
     322           48 :      call klmn2ijlmn(klmn,lmn_size,ilmn,jlmn)
     323              : 
     324           48 :      il= indlmn(1,ilmn); iln=indlmn(5,ilmn)
     325           48 :      jl= indlmn(1,jlmn); jln=indlmn(5,jlmn)
     326              :      !$in = indklmn(9 ,klmn) !=indlmn(3,ilmn)
     327              :      !$jn = indklmn(10,klmn) !=indlmn(3,jlmn)
     328              : 
     329           48 :      in = indlmn(3,ilmn)
     330           48 :      jn = indlmn(3,jlmn)
     331              : 
     332              :      ! Shift il to have the index instead of the value of l.
     333           48 :      kln2ln(1,kln) = il +1
     334           48 :      kln2ln(2,kln) = jl +1
     335           48 :      kln2ln(3,kln) = in
     336           48 :      kln2ln(4,kln) = jn
     337           48 :      kln2ln(5,kln) = iln
     338           48 :      kln2ln(6,kln) = jln
     339              :    end if
     340              :  end do !klmn
     341              : 
     342            2 : end subroutine make_kln2ln
     343              : !!***
     344              : 
     345              : !----------------------------------------------------------------------
     346              : 
     347              : !!****f* m_paw_lmn/make_klm2lm
     348              : !! NAME
     349              : !!  make_klm2lm
     350              : !!
     351              : !! FUNCTION
     352              : !!  Performs the setup of the klm2lm table for PAW calculations.
     353              : !!
     354              : !! INPUTS
     355              : !!  lmn_size=Number of (l,m,n) elements for the PAW basis set.
     356              : !!  lmn2_size=Number of elements in the symmetric basis set: lmn2_size=lmn_size*(lmn_size+1)/2
     357              : !!  lm2_size)=Number of (l.m) elements in the symmetric basis set.
     358              : !!  indlmn(6,lmn_size)=Array giving l,m,n,lm,ln,spin for i=lmn.
     359              : !!  indklmn(8,lmn2_size)=Array giving klm, kln, abs(il-jl), (il+jl), ilm and jlm
     360              : !!   for each klmn=(ilmn,jlmn). Note: ilmn=(il,im,in) and ilmn<=jlmn
     361              : !!
     362              : !! OUTPUT
     363              : !!  klm2lm(6,lm2_size)=Table giving il, jl ,im, jm, ilm, jlm for each klm=(ilm,jlm)
     364              : !!  where ilm=(il,im) and ilm<=jlm. NB: klm2lm is an application and not a bijection.
     365              : !!
     366              : !! NOTES
     367              : !!  klm2lm can be calculated easily if we assume that all (l,m) channels
     368              : !!  are ordered by increasing l and m. This is the standard convention
     369              : !!  used in most of the PAW datasets. This routines, howevever, works
     370              : !!  works also in the unlikely case in with (l,m) are not ordered.
     371              : !!
     372              : !! SOURCE
     373              : 
     374            1 : subroutine make_klm2lm(lmn_size,lmn2_size,lm2_size,indlmn,indklmn,klm2lm)
     375              : 
     376              : !Arguments ------------------------------------
     377              : !scalars
     378              :  integer,intent(in) :: lmn2_size,lmn_size,lm2_size
     379              : !arrays
     380              :  integer,intent(in) :: indlmn(6,lmn_size)
     381              :  integer,intent(in) :: indklmn(8,lmn2_size)
     382              :  integer,intent(out) :: klm2lm(6,lm2_size)
     383              : 
     384              : !Local variables ------------------------------
     385              : !scalars
     386              :  integer :: il,ilm,ilmn,jl,jlm,jlmn,im,jm
     387              :  integer :: klmn,klm,klm_old !,iklm
     388              : 
     389              : !************************************************************************
     390              : 
     391           71 :  klm2lm = 0
     392            1 :  klm_old = -1
     393           37 :  do klmn=1,lmn2_size
     394           36 :    klm = indklmn(1,klmn)
     395           37 :    if (klm /= klm_old) then
     396           28 :      klm_old = klm
     397              : 
     398           28 :      ilm = indklmn(5,klmn)
     399           28 :      jlm = indklmn(6,klmn)
     400              : 
     401           28 :      call klmn2ijlmn(klmn,lmn_size,ilmn,jlmn)
     402              : 
     403           28 :      il=indlmn(1,ilmn); im=indlmn(2,ilmn)
     404           28 :      jl=indlmn(1,jlmn); jm=indlmn(2,jlmn)
     405              : 
     406              :      !shift to have the index instead of l or m.
     407           28 :      klm2lm(1,klm) = il +1       ! il
     408           28 :      klm2lm(2,klm) = jl +1       ! jl
     409           28 :      klm2lm(3,klm) = im + il +1  ! im
     410           28 :      klm2lm(4,klm) = jm + jl +1  ! jm
     411           28 :      klm2lm(5,klm) = ilm         ! ilm
     412           28 :      klm2lm(6,klm) = jlm         ! jlm
     413              :    end if
     414              :  end do !klmn
     415              : 
     416           71 :  if (ANY(klm2lm==0)) then
     417            0 :    LIBPAW_BUG("check klm2lm")
     418              :  end if
     419              : 
     420              : !DEBUG
     421              : #if 0
     422              :  write(std_out,*)"Debugging make_klm2lm:"
     423              :  do iklm=1,lm2_size
     424              :   il  = klm2lm(1,iklm) -1      ! li
     425              :   jl  = klm2lm(2,iklm) -1      ! lj
     426              :   im  = klm2lm(3,iklm) -il -1  ! mi
     427              :   jm  = klm2lm(4,iklm) -jl -1  ! mj
     428              :   ilm = klm2lm(5,iklm)         ! ilm
     429              :   jlm = klm2lm(6,iklm)         ! jlm
     430              :   write(std_out,'(i3,2(a,2i3,a))')"iklm ",iklm," l m (",il,im,")"," l m (",jl,jm,")"
     431              :  end do
     432              : #endif
     433              : 
     434            1 : end subroutine make_klm2lm
     435              : !!***
     436              : 
     437              : !----------------------------------------------------------------------
     438              : 
     439              : !!****f* m_paw_lmn/klmn2ijlmn
     440              : !! NAME
     441              : !!  klmn2ijlmn
     442              : !!
     443              : !! FUNCTION
     444              : !!  Find ilmn and jlmn from klmn and lmn_size.
     445              : !!
     446              : !! INPUTS
     447              : !!  lmn_size=Number of (l,m,n) elements for the PAW basis set.
     448              : !!  klmn=The index corresponding to (ilmn,jlmn) in packed form.
     449              : !!
     450              : !! OUTPUT
     451              : !!  jlmn, ilmn=The two symmetrix indices corresponding to klmn. NB: jlmn >= ilmn
     452              : !!
     453              : !! SOURCE
     454              : 
     455           76 : subroutine klmn2ijlmn(klmn,lmn_size,ilmn,jlmn)
     456              : 
     457              : !Arguments ------------------------------------
     458              : !scalars
     459              :  integer,intent(in) :: klmn,lmn_size
     460              :  integer,intent(out) :: ilmn,jlmn
     461              : 
     462              : !Local variables ------------------------------
     463              : !scalars
     464              :  integer :: ii,jj,k0
     465              : !************************************************************************
     466              : 
     467           76 :  ilmn=-1; jlmn=-1
     468              : 
     469          417 :  do jj=1,lmn_size
     470          417 :    k0=jj*(jj-1)/2
     471         1653 :    do ii=1,jj
     472         1653 :      if (klmn==ii+k0) then
     473           76 :        ilmn = ii
     474           76 :        jlmn = jj; RETURN
     475              :      end if
     476              :    end do
     477              :  end do
     478              : 
     479            0 :  LIBPAW_BUG("Not able to found ilmn and jlmn")
     480              : 
     481              : end subroutine klmn2ijlmn
     482              : !!***
     483              : 
     484              : !----------------------------------------------------------------------
     485              : 
     486              : !!****f* m_paw_lmn/make_indln
     487              : !! NAME
     488              : !!  make_indln
     489              : !!
     490              : !! FUNCTION
     491              : !!  Performs the setup of the indln table for PAW calculations.
     492              : !!  Compute the indln indexes giving ilmn=(l,m,n)
     493              : !!
     494              : !! INPUTS
     495              : !!  lmn_size=Number of (l,m,n) elements for the PAW basis set.
     496              : !!  ln_size=Number of (l,n) elements
     497              : !!  indlmn(6,lmn_size)=Array giving l,m,n,lm,ln,spin for i=lmn.
     498              : !!
     499              : !! OUTPUT
     500              : !!  indln(2,ln_size)=Array giving l and n for i=ln
     501              : !!
     502              : !! SOURCE
     503              : 
     504            1 : subroutine make_indln(lmn_size,ln_size,indlmn,indln)
     505              : 
     506              : !Arguments ------------------------------------
     507              : !scalars
     508              :  integer,intent(in) :: lmn_size,ln_size
     509              : !arrays
     510              :  integer,intent(in) :: indlmn(6,lmn_size)
     511              :  integer,intent(out) :: indln(2,ln_size)
     512              : 
     513              : !Local variables ------------------------------
     514              : !scalars
     515              :  integer :: ilmn,ll,nn,ll_old,nn_old,ii
     516              : !************************************************************************
     517              : 
     518            1 :  ii=0; ll_old=-1; nn_old=-1
     519            9 :  do ilmn=1,lmn_size
     520            8 :   ll=indlmn(1,ilmn)
     521            8 :   nn=indlmn(3,ilmn)
     522            9 :   if (ll/=ll_old.or.nn/=nn_old) then
     523            4 :    ll_old=ll
     524            4 :    nn_old=nn
     525            4 :    ii=ii+1
     526            4 :    indln(1,ii)=ll
     527            4 :    indln(2,ii)=nn
     528              :   end if
     529              :  end do
     530              : 
     531            1 :  if(ii/=ln_size) LIBPAW_ERROR("ii/=ln_size")
     532              : 
     533            1 : end subroutine make_indln
     534              : !!***
     535              : 
     536              : !----------------------------------------------------------------------
     537              : 
     538              : !!****f* m_paw_lmn/uppert_index
     539              : !! NAME
     540              : !!   uppert_index
     541              : !!
     542              : !! FUNCTION
     543              : !!  Helper function returning the sequential index of an element in the upper triangle of a matrix
     544              : !!  given the row-index ii and the column-index jj. If ii>jj the index of the element a_{jj,ii} is returned.
     545              : !!
     546              : !! INPUTS
     547              : !!  ii=Row index
     548              : !!  jj=column index
     549              : !!
     550              : !! SOURCE
     551              : 
     552            0 : function uppert_index(ii,jj)
     553              : 
     554              : !Arguments ------------------------------------
     555              : !scalars
     556              :  integer,intent(in) :: ii,jj
     557              :  integer :: uppert_index
     558              : !scalars
     559              : 
     560              : !************************************************************************
     561              : 
     562            0 :  if (jj>=jj) then
     563            0 :    uppert_index = ii + jj*(jj-1)/2
     564              :  else
     565              :    uppert_index = jj + ii*(ii-1)/2
     566              :  end if
     567              : 
     568            0 : end function uppert_index
     569              : !!***
     570              : 
     571              : END MODULE m_paw_lmn
     572              : !!***
        

Generated by: LCOV version 2.3-1