LCOV - code coverage report
Current view: top level - hybrid - spmvec.F90 (source / functions) Hit Total Coverage
Test: combined.info Lines: 0 246 0.0 %
Date: 2019-09-08 04:53:50 Functions: 0 3 0.0 %

          Line data    Source code
       1             : MODULE m_spmvec
       2             : 
       3             :       CONTAINS
       4             :          !Note this module contains a real/complex version of spmvec
       5             : 
       6           0 :          SUBROUTINE spmvec_invs(&
       7             :         &           atoms, hybrid,&
       8             :         &           hybdat, ikpt, kpts, cell,&
       9           0 :         &           coulomb_mt1, coulomb_mt2, coulomb_mt3,&
      10           0 :         &           coulomb_mtir, vecin,&
      11           0 :         &           vecout)
      12             : 
      13             :             USE m_wrapper
      14             :             USE m_constants
      15             :             USE m_types
      16             :             IMPLICIT NONE
      17             :             TYPE(t_hybdat), INTENT(IN)   :: hybdat
      18             :             TYPE(t_hybrid), INTENT(IN)   :: hybrid
      19             :             TYPE(t_cell), INTENT(IN)     :: cell
      20             :             TYPE(t_kpts), INTENT(IN)     :: kpts
      21             :             TYPE(t_atoms), INTENT(IN)    :: atoms
      22             : 
      23             :             ! - scalars -
      24             :             INTEGER, INTENT(IN) ::  ikpt
      25             : 
      26             :             ! - arrays -
      27             :             REAL, INTENT(IN) ::  coulomb_mt1(hybrid%maxindxm1 - 1, hybrid%maxindxm1 - 1,&
      28             :            &                                    0:hybrid%maxlcutm1, atoms%ntype)
      29             :             REAL, INTENT(IN) ::  coulomb_mt2(hybrid%maxindxm1 - 1, -hybrid%maxlcutm1:hybrid%maxlcutm1,&
      30             :            &                                    0:hybrid%maxlcutm1 + 1, atoms%nat)
      31             :             REAL, INTENT(IN) ::  coulomb_mt3(hybrid%maxindxm1 - 1, atoms%nat, atoms%nat)
      32             : #ifdef CPP_IRCOULOMBAPPROX
      33             :             REAL, INTENT(IN) ::  coulomb_mtir(:, :)
      34             : #else
      35             :             REAL, INTENT(IN) ::  coulomb_mtir(:)
      36             : #endif
      37             :             REAL, INTENT(IN) ::  vecin(:)!(hybrid%nbasm)
      38             :             REAL, INTENT(INOUT)::  vecout(:)!(hybrid%nbasm)
      39             : 
      40             :             ! - local scalars -
      41             :             INTEGER             ::  itype, ieq, iatom, ishift
      42             :             INTEGER             ::  itype1, ieq1, iatom1, ishift1
      43             :             INTEGER             ::  indx0, indx1, indx2, indx3, indx4
      44             :             INTEGER             ::  i, ibasm, igptm
      45             :             INTEGER             ::  l
      46             :             INTEGER             ::  n, m
      47             : 
      48             :             REAL                ::  gnorm
      49             : 
      50             :             COMPLEX, PARAMETER  ::  img = (0.0, 1.0)
      51             :             ! - local arrays -
      52             : 
      53           0 :             REAL                ::  vecinhlp(hybrid%nbasm(ikpt))
      54             :             REAL, ALLOCATABLE ::  coulhlp(:, :)
      55             : 
      56           0 :             vecinhlp = vecin
      57             : 
      58           0 :             CALL reorder(hybrid%nbasm(ikpt), hybrid%nbasp, atoms, hybrid%lcutm1, hybrid%maxlcutm1, hybrid%nindxm1, 1, vecinhlp)
      59             : 
      60           0 :             ibasm = 0
      61           0 :             iatom = 0
      62           0 :             DO itype = 1, atoms%ntype
      63           0 :                DO ieq = 1, atoms%neq(itype)
      64           0 :                   iatom = iatom + 1
      65           0 :                   DO l = 0, hybrid%lcutm1(itype)
      66           0 :                      DO m = -l, l
      67           0 :                         ibasm = ibasm + hybrid%nindxm1(l, itype) - 1
      68             :                      END DO
      69             :                   END DO
      70             :                END DO
      71             :             END DO
      72             : 
      73             :             ! compute vecout for the indices from 0:ibasm
      74           0 :             iatom = 0
      75           0 :             indx1 = 0; indx2 = 0; indx3 = ibasm
      76           0 :             DO itype = 1, atoms%ntype
      77           0 :                DO ieq = 1, atoms%neq(itype)
      78           0 :                   iatom = iatom + 1
      79           0 :                   DO l = 0, hybrid%lcutm1(itype)
      80           0 :                      DO m = -l, l
      81           0 :                         indx1 = indx1 + 1
      82           0 :                         indx2 = indx2 + hybrid%nindxm1(l, itype) - 1
      83           0 :                         indx3 = indx3 + 1
      84             : 
      85             :                         vecout(indx1:indx2) = matmul(coulomb_mt1(:hybrid%nindxm1(l, itype) - 1, :hybrid%nindxm1(l, itype) - 1, l, itype),&
      86           0 :                &                vecinhlp(indx1:indx2))
      87             : 
      88           0 :                         vecout(indx1:indx2) = vecout(indx1:indx2) + coulomb_mt2(:hybrid%nindxm1(l, itype) - 1, m, l, iatom)*vecinhlp(indx3)
      89             : 
      90           0 :                         indx1 = indx2
      91             :                      END DO
      92             : 
      93             :                   END DO
      94             :                END DO
      95             :             END DO
      96             : 
      97           0 :             IF (indx2 /= ibasm) STOP 'spmvec: error counting basis functions'
      98             : 
      99           0 :             IF (ikpt == 1) THEN
     100             :                iatom = 0
     101             :                indx0 = 0
     102           0 :                DO itype = 1, atoms%ntype
     103           0 :                   ishift = sum((/((2*l + 1)*(hybrid%nindxm1(l, itype) - 1), l=0, hybrid%lcutm1(itype))/))
     104           0 :                   DO ieq = 1, atoms%neq(itype)
     105           0 :                      iatom = iatom + 1
     106           0 :                      l = 0
     107           0 :                      m = 0
     108             : 
     109           0 :                      indx1 = indx0 + 1
     110           0 :                      indx2 = indx1 + hybrid%nindxm1(l, itype) - 2
     111             : 
     112           0 :                      iatom1 = 0
     113           0 :                      indx3 = ibasm
     114           0 :                      DO itype1 = 1, atoms%ntype
     115           0 :                         ishift1 = (hybrid%lcutm1(itype1) + 1)**2
     116           0 :                         DO ieq1 = 1, atoms%neq(itype1)
     117           0 :                            iatom1 = iatom1 + 1
     118           0 :                            indx4 = indx3 + (ieq1 - 1)*ishift1 + 1
     119           0 :                            IF (iatom == iatom1) CYCLE
     120             : 
     121           0 :                            vecout(indx1:indx2) = vecout(indx1:indx2) + coulomb_mt3(:hybrid%nindxm1(l, itype) - 1, iatom1, iatom)*vecinhlp(indx4)
     122             : 
     123             :                         END DO
     124           0 :                         indx3 = indx3 + atoms%neq(itype1)*ishift1
     125             :                      END DO
     126             : 
     127           0 :                      IF (indx3 /= hybrid%nbasp) STOP 'spmvec: error counting index indx3'
     128             : 
     129           0 :                      vecout(indx1:indx2) = vecout(indx1:indx2) + coulomb_mt2(:hybrid%nindxm1(l, itype) - 1, 0, hybrid%maxlcutm1 + 1, iatom)*vecinhlp(indx3 + 1)
     130             : 
     131           0 :                      indx0 = indx0 + ishift
     132             :                   END DO
     133             : 
     134             :                END DO
     135             :             END IF
     136             : 
     137             :             ! compute vecout for the index-range from ibasm+1:nbasm
     138             : 
     139             : #ifdef CPP_IRCOULOMBAPPROX
     140             : 
     141             :             indx0 = sum((/(((2*l + 1)*atoms%neq(itype), l=0, hybrid%lcutm1(itype)), itype=1, atoms%ntype)/)) + hybrid%ngptm
     142             :             indx1 = sum((/(((2*l + 1)*atoms%neq(itype), l=0, hybrid%lcutm1(itype)), itype=1, atoms%ntype)/))
     143             : 
     144             :             CALL DGEMV('N', indx1, indx0, 1.0, coulomb_mtir, (hybrid%maxlcutm1 + 1)**2*atoms,&
     145             :            &          vecinhlp(ibasm + 1:), 1, 0.0, vecout(ibasm + 1:), 1)
     146             : 
     147             :             CALL DGEMV('T', indx1, hybrid, 1.0, coulomb_mtir(:indx1, indx1 + 1:),&
     148             :            &          indx1, vecinhlp(ibasm + 1:), 1, 0.0, vecout(ibasm + indx1 + 1:), 1)
     149             : 
     150             : !       vecout(ibasm+1:ibasm+indx1) = matmul( coulomb_mtir(:indx1,:indx0),vecinhlp(ibasm+1:ibasm+indx0) )
     151             : !       vecout(ibasm+indx1+1:ibasm+indx0) = matmul( conjg(transpose(coulomb_mtir(:indx1,indx1+1:indx0))),
     152             : !      &                                            vecinhlp(ibasm+1:ibasm+indx1) )
     153             : 
     154             :             indx0 = ibasm + indx1
     155             :             IF (indx0 /= hybrid%nbasp) STOP 'spmvec: error indx0'
     156             :             DO i = 1, hybrid%ngptm
     157             :                indx0 = indx0 + 1
     158             :                igptm = hybrid%pgptm(i)
     159             :                gnorm = sqrt(sum(matmul(kpts%bk(:) + hybrid%gptm(:, igptm), cell%bmat)**2))
     160             :                IF (gnorm == 0) CYCLE
     161             :                vecout(indx0) = vecout(indx0) + fpi*vecinhlp(indx0)/gnorm
     162             :             END DO
     163             : 
     164             : #else
     165             : 
     166             :             indx1 = sum((/(((2*l + 1)*atoms%neq(itype), l=0, hybrid%lcutm1(itype)),&
     167           0 :            &                                      itype=1, atoms%ntype)/)) + hybrid%ngptm(ikpt)
     168           0 :             CALL dspmv('U', indx1, 1.0, coulomb_mtir, vecinhlp(ibasm + 1:), 1, 0.0, vecout(ibasm + 1:), 1)
     169             : 
     170             : #endif
     171           0 :             iatom = 0
     172           0 :             indx1 = ibasm; indx2 = 0; indx3 = 0
     173           0 :             DO itype = 1, atoms%ntype
     174           0 :                DO ieq = 1, atoms%neq(itype)
     175           0 :                   iatom = iatom + 1
     176           0 :                   DO l = 0, hybrid%lcutm1(itype)
     177           0 :                      n = hybrid%nindxm1(l, itype)
     178           0 :                      DO m = -l, l
     179           0 :                         indx1 = indx1 + 1
     180           0 :                         indx2 = indx2 + 1
     181           0 :                         indx3 = indx3 + n - 1
     182             : 
     183           0 :                         vecout(indx1) = vecout(indx1) + dotprod(coulomb_mt2(:n - 1, m, l, iatom), vecinhlp(indx2:indx3))
     184           0 :                         indx2 = indx3
     185             :                      END DO
     186             : 
     187             :                   END DO
     188             :                END DO
     189             :             END DO
     190             : 
     191           0 :             IF (ikpt == 1) THEN
     192             :                iatom = 0
     193             :                indx0 = 0
     194           0 :                DO itype = 1, atoms%ntype
     195           0 :                   ishift = sum((/((2*l + 1)*(hybrid%nindxm1(l, itype) - 1), l=0, hybrid%lcutm1(itype))/))
     196           0 :                   DO ieq = 1, atoms%neq(itype)
     197           0 :                      iatom = iatom + 1
     198           0 :                      indx1 = indx0 + 1
     199           0 :                      indx2 = indx1 + hybrid%nindxm1(0, itype) - 2
     200           0 :                      vecout(hybrid%nbasp + 1) = vecout(hybrid%nbasp + 1) + dotprod(coulomb_mt2(:hybrid%nindxm1(0, itype) - 1, 0, hybrid%maxlcutm1 + 1, iatom), vecinhlp(indx1:indx2))
     201             : 
     202           0 :                      indx0 = indx0 + ishift
     203             :                   END DO
     204             :                END DO
     205             : 
     206             :                iatom = 0
     207             :                indx0 = ibasm
     208           0 :                DO itype = 1, atoms%ntype
     209           0 :                   ishift = (hybrid%lcutm1(itype) + 1)**2
     210           0 :                   DO ieq = 1, atoms%neq(itype)
     211           0 :                      iatom = iatom + 1
     212           0 :                      indx1 = indx0 + 1
     213             : 
     214           0 :                      iatom1 = 0
     215           0 :                      indx2 = 0
     216           0 :                      DO itype1 = 1, atoms%ntype
     217           0 :                         ishift1 = sum((/((2*l + 1)*(hybrid%nindxm1(l, itype1) - 1), l=0, hybrid%lcutm1(itype1))/))
     218           0 :                         DO ieq1 = 1, atoms%neq(itype1)
     219           0 :                            iatom1 = iatom1 + 1
     220           0 :                            IF (iatom1 == iatom) CYCLE
     221             : 
     222           0 :                            indx3 = indx2 + (ieq1 - 1)*ishift1 + 1
     223           0 :                            indx4 = indx3 + hybrid%nindxm1(0, itype1) - 2
     224             : 
     225           0 :                            vecout(indx1) = vecout(indx1) + dotprod(coulomb_mt3(:hybrid%nindxm1(0, itype1) - 1, iatom, iatom1), vecinhlp(indx3:indx4))
     226             : 
     227             :                         END DO
     228           0 :                         indx2 = indx2 + atoms%neq(itype1)*ishift1
     229             :                      END DO
     230           0 :                      indx0 = indx0 + ishift
     231             :                   END DO
     232             :                END DO
     233           0 :                IF (indx0 /= hybrid%nbasp) STOP 'spmvec: error index counting (indx0)'
     234             :             END IF
     235             : 
     236             :             CALL reorder(hybrid%nbasm(ikpt), hybrid%nbasp, atoms, hybrid%lcutm1, hybrid%maxlcutm1, hybrid%nindxm1,&
     237             :            &             2,&
     238           0 :            &             vecout)
     239             : 
     240           0 :          END SUBROUTINE spmvec_invs
     241             : 
     242           0 :          SUBROUTINE spmvec_noinvs(&
     243             :           &           atoms, hybrid,&
     244             :           &           hybdat, ikpt, kpts, cell,&
     245           0 :           &           coulomb_mt1, coulomb_mt2, coulomb_mt3,&
     246           0 :           &           coulomb_mtir, vecin,&
     247           0 :           &           vecout)
     248             : 
     249             :             USE m_wrapper
     250             :             USE m_constants
     251             :             USE m_types
     252             :             IMPLICIT NONE
     253             :             TYPE(t_hybdat), INTENT(IN)   :: hybdat
     254             :             TYPE(t_hybrid), INTENT(IN)   :: hybrid
     255             :             TYPE(t_cell), INTENT(IN)     :: cell
     256             :             TYPE(t_kpts), INTENT(IN)     :: kpts
     257             :             TYPE(t_atoms), INTENT(IN)    :: atoms
     258             : 
     259             :             ! - scalars -
     260             :             INTEGER, INTENT(IN) ::  ikpt
     261             : 
     262             :             ! - arrays -
     263             :             REAL, INTENT(IN) ::  coulomb_mt1(hybrid%maxindxm1 - 1, hybrid%maxindxm1 - 1,&
     264             :            &                                    0:hybrid%maxlcutm1, atoms%ntype)
     265             :             COMPLEX, INTENT(IN) ::  coulomb_mt2(hybrid%maxindxm1 - 1, -hybrid%maxlcutm1:hybrid%maxlcutm1,&
     266             :            &                                    0:hybrid%maxlcutm1 + 1, atoms%nat)
     267             :             COMPLEX, INTENT(IN) ::  coulomb_mt3(hybrid%maxindxm1 - 1, atoms%nat, atoms%nat)
     268             : #ifdef CPP_IRCOULOMBAPPROX
     269             :             COMPLEX, INTENT(IN) ::  coulomb_mtir(:, :)
     270             : #else
     271             :             COMPLEX, INTENT(IN) ::  coulomb_mtir(:)
     272             : #endif
     273             :             COMPLEX, INTENT(IN) ::  vecin(:)!(hybrid%nbasm)
     274             :             COMPLEX, INTENT(OUT)::  vecout(:)!(hybrid%nbasm)
     275             : 
     276             :             ! - local scalars -
     277             :             INTEGER             ::  itype, ieq, iatom, ishift
     278             :             INTEGER             ::  itype1, ieq1, iatom1, ishift1
     279             :             INTEGER             ::  indx0, indx1, indx2, indx3, indx4
     280             :             INTEGER             ::  i, ibasm, igptm
     281             :             INTEGER             ::  l
     282             :             INTEGER             ::  n, m
     283             : 
     284             :             REAL                ::  gnorm
     285             : 
     286             :             COMPLEX, PARAMETER  ::  img = (0.0, 1.0)
     287             :             ! - local arrays -
     288             : 
     289             :             REAL                ::  vecr(hybrid%maxindxm1 - 1), veci(hybrid%maxindxm1 - 1)
     290           0 :             COMPLEX             ::  vecinhlp(hybrid%nbasm(ikpt))
     291             :             COMPLEX, ALLOCATABLE ::  coulhlp(:, :)
     292             : 
     293           0 :             vecinhlp = vecin
     294             : 
     295           0 :             CALL reorder(hybrid%nbasm(ikpt), hybrid%nbasp, atoms, hybrid%lcutm1, hybrid%maxlcutm1, hybrid%nindxm1, 1, vec_c=vecinhlp)
     296             : 
     297           0 :             ibasm = 0
     298           0 :             iatom = 0
     299           0 :             DO itype = 1, atoms%ntype
     300           0 :                DO ieq = 1, atoms%neq(itype)
     301           0 :                   iatom = iatom + 1
     302           0 :                   DO l = 0, hybrid%lcutm1(itype)
     303           0 :                      DO m = -l, l
     304           0 :                         ibasm = ibasm + hybrid%nindxm1(l, itype) - 1
     305             :                      END DO
     306             :                   END DO
     307             :                END DO
     308             :             END DO
     309             : 
     310             :             ! compute vecout for the indices from 0:ibasm
     311           0 :             iatom = 0
     312           0 :             indx1 = 0; indx2 = 0; indx3 = ibasm
     313           0 :             DO itype = 1, atoms%ntype
     314           0 :                DO ieq = 1, atoms%neq(itype)
     315           0 :                   iatom = iatom + 1
     316           0 :                   DO l = 0, hybrid%lcutm1(itype)
     317           0 :                      DO m = -l, l
     318           0 :                         indx1 = indx1 + 1
     319           0 :                         indx2 = indx2 + hybrid%nindxm1(l, itype) - 1
     320           0 :                         indx3 = indx3 + 1
     321             : 
     322             :                         vecout(indx1:indx2) = matmul(coulomb_mt1(:hybrid%nindxm1(l, itype) - 1, :hybrid%nindxm1(l, itype) - 1, l, itype),&
     323           0 :                &                vecinhlp(indx1:indx2))
     324             : 
     325           0 :                         vecout(indx1:indx2) = vecout(indx1:indx2) + coulomb_mt2(:hybrid%nindxm1(l, itype) - 1, m, l, iatom)*vecinhlp(indx3)
     326             : 
     327           0 :                         indx1 = indx2
     328             :                      END DO
     329             : 
     330             :                   END DO
     331             :                END DO
     332             :             END DO
     333             : 
     334           0 :             IF (indx2 /= ibasm) &
     335           0 :            &  STOP 'spmvec: error counting basis functions'
     336             : 
     337           0 :             IF (ikpt == 1) THEN
     338             :                iatom = 0
     339             :                indx0 = 0
     340           0 :                DO itype = 1, atoms%ntype
     341           0 :                   ishift = sum((/((2*l + 1)*(hybrid%nindxm1(l, itype) - 1), l=0, hybrid%lcutm1(itype))/))
     342           0 :                   DO ieq = 1, atoms%neq(itype)
     343           0 :                      iatom = iatom + 1
     344           0 :                      l = 0
     345           0 :                      m = 0
     346             : 
     347           0 :                      indx1 = indx0 + 1
     348           0 :                      indx2 = indx1 + hybrid%nindxm1(l, itype) - 2
     349             : 
     350           0 :                      iatom1 = 0
     351           0 :                      indx3 = ibasm
     352           0 :                      DO itype1 = 1, atoms%ntype
     353           0 :                         ishift1 = (hybrid%lcutm1(itype1) + 1)**2
     354           0 :                         DO ieq1 = 1, atoms%neq(itype1)
     355           0 :                            iatom1 = iatom1 + 1
     356           0 :                            indx4 = indx3 + (ieq1 - 1)*ishift1 + 1
     357           0 :                            IF (iatom == iatom1) CYCLE
     358             : 
     359             :                            vecout(indx1:indx2) = vecout(indx1:indx2)&
     360           0 :                 &                              + coulomb_mt3(:hybrid%nindxm1(l, itype) - 1, iatom1, iatom)*vecinhlp(indx4)
     361             : 
     362             :                         END DO
     363           0 :                         indx3 = indx3 + atoms%neq(itype1)*ishift1
     364             :                      END DO
     365             : 
     366           0 :                      IF (indx3 /= hybrid%nbasp) STOP 'spmvec: error counting index indx3'
     367             : 
     368           0 :                      vecout(indx1:indx2) = vecout(indx1:indx2) + coulomb_mt2(:hybrid%nindxm1(l, itype) - 1, 0, hybrid%maxlcutm1 + 1, iatom)*vecinhlp(indx3 + 1)
     369             : 
     370           0 :                      indx0 = indx0 + ishift
     371             :                   END DO
     372             : 
     373             :                END DO
     374             :             END IF
     375             : 
     376             :             ! compute vecout for the index-range from ibasm+1:nbasm
     377             : 
     378             : #ifdef CPP_IRCOULOMBAPPROX
     379             : 
     380             :             indx0 = sum((/(((2*l + 1)*atoms%neq(itype), l=0, hybrid%lcutm1(itype)), itype=1, atoms%ntype)/)) + hybrid%ngptm
     381             :             indx1 = sum((/(((2*l + 1)*atoms%neq(itype), l=0, hybrid%lcutm1(itype)), itype=1, atoms%ntype)/))
     382             : 
     383             :             CALL ZGEMV('N', indx1, indx0, (1.0, 0.0), coulomb_mtir,&
     384             :            &          (hybrid%maxlcutm1 + 1)**2*atoms, vecinhlp(ibasm + 1:),&
     385             :            &          1, (0.0, 0.0), vecout(ibasm + 1:), 1)
     386             : 
     387             :             CALL ZGEMV('C', indx1, hybrid, (1.0, 0.0), coulomb_mtir(:indx1, indx1 + 1:)&
     388             :            &          , indx1, vecinhlp(ibasm + 1:), 1, (0.0, 0.0),&
     389             :            &          vecout(ibasm + indx1 + 1:), 1)
     390             : 
     391             : !       vecout(ibasm+1:ibasm+indx1) = matmul( coulomb_mtir(:indx1,:indx0),vecinhlp(ibasm+1:ibasm+indx0) )
     392             : !       vecout(ibasm+indx1+1:ibasm+indx0) = matmul( conjg(transpose(coulomb_mtir(:indx1,indx1+1:indx0))),
     393             : !      &                                            vecinhlp(ibasm+1:ibasm+indx1) )
     394             : 
     395             :             indx0 = ibasm + indx1
     396             :             IF (indx0 /= hybrid%nbasp) STOP 'spmvec: error indx0'
     397             :             DO i = 1, hybrid%ngptm
     398             :                indx0 = indx0 + 1
     399             :                igptm = hybrid%pgptm(i)
     400             :                gnorm = sqrt(sum(matmul(kpts%bk(:) + hybrid%gptm(:, igptm), cell%bmat)**2))
     401             :                IF (gnorm == 0) CYCLE
     402             :                vecout(indx0) = vecout(indx0) + fpi*vecinhlp(indx0)/gnorm
     403             :             END DO
     404             : 
     405             : #else
     406             : 
     407             :             indx1 = sum((/(((2*l + 1)*atoms%neq(itype), l=0, hybrid%lcutm1(itype)),&
     408           0 :            &                                      itype=1, atoms%ntype)/)) + hybrid%ngptm(ikpt)
     409           0 :             call zhpmv('U', indx1, (1.0, 0.0), coulomb_mtir, vecinhlp(ibasm + 1), 1, (0.0, 0.0), vecout(ibasm + 1), 1)
     410             : 
     411             : #endif
     412           0 :             iatom = 0
     413           0 :             indx1 = ibasm; indx2 = 0; indx3 = 0
     414           0 :             DO itype = 1, atoms%ntype
     415           0 :                DO ieq = 1, atoms%neq(itype)
     416           0 :                   iatom = iatom + 1
     417           0 :                   DO l = 0, hybrid%lcutm1(itype)
     418           0 :                      n = hybrid%nindxm1(l, itype)
     419           0 :                      DO m = -l, l
     420           0 :                         indx1 = indx1 + 1
     421           0 :                         indx2 = indx2 + 1
     422           0 :                         indx3 = indx3 + n - 1
     423             : 
     424           0 :                         vecout(indx1) = vecout(indx1) + dotprod(coulomb_mt2(:n - 1, m, l, iatom), vecinhlp(indx2:indx3))
     425           0 :                         indx2 = indx3
     426             :                      END DO
     427             : 
     428             :                   END DO
     429             :                END DO
     430             :             END DO
     431             : 
     432           0 :             IF (ikpt == 1) THEN
     433             :                iatom = 0
     434             :                indx0 = 0
     435           0 :                DO itype = 1, atoms%ntype
     436           0 :                   ishift = sum((/((2*l + 1)*(hybrid%nindxm1(l, itype) - 1), l=0, hybrid%lcutm1(itype))/))
     437           0 :                   DO ieq = 1, atoms%neq(itype)
     438           0 :                      iatom = iatom + 1
     439           0 :                      indx1 = indx0 + 1
     440           0 :                      indx2 = indx1 + hybrid%nindxm1(0, itype) - 2
     441           0 :                      vecout(hybrid%nbasp + 1) = vecout(hybrid%nbasp + 1) + dotprod(coulomb_mt2(:hybrid%nindxm1(0, itype) - 1, 0, hybrid%maxlcutm1 + 1, iatom), vecinhlp(indx1:indx2))
     442             : 
     443           0 :                      indx0 = indx0 + ishift
     444             :                   END DO
     445             :                END DO
     446             : 
     447             :                iatom = 0
     448             :                indx0 = ibasm
     449           0 :                DO itype = 1, atoms%ntype
     450           0 :                   ishift = (hybrid%lcutm1(itype) + 1)**2
     451           0 :                   DO ieq = 1, atoms%neq(itype)
     452           0 :                      iatom = iatom + 1
     453           0 :                      indx1 = indx0 + 1
     454             : 
     455           0 :                      iatom1 = 0
     456           0 :                      indx2 = 0
     457           0 :                      DO itype1 = 1, atoms%ntype
     458           0 :                         ishift1 = sum((/((2*l + 1)*(hybrid%nindxm1(l, itype1) - 1), l=0, hybrid%lcutm1(itype1))/))
     459           0 :                         DO ieq1 = 1, atoms%neq(itype1)
     460           0 :                            iatom1 = iatom1 + 1
     461           0 :                            IF (iatom1 == iatom) CYCLE
     462             : 
     463           0 :                            indx3 = indx2 + (ieq1 - 1)*ishift1 + 1
     464           0 :                            indx4 = indx3 + hybrid%nindxm1(0, itype1) - 2
     465             : 
     466           0 :                            vecout(indx1) = vecout(indx1) + dotprod(coulomb_mt3(:hybrid%nindxm1(0, itype1) - 1, iatom, iatom1), vecinhlp(indx3:indx4))
     467             : 
     468             :                         END DO
     469           0 :                         indx2 = indx2 + atoms%neq(itype1)*ishift1
     470             :                      END DO
     471           0 :                      indx0 = indx0 + ishift
     472             :                   END DO
     473             :                END DO
     474           0 :                IF (indx0 /= hybrid%nbasp) STOP 'spmvec: error index counting (indx0)'
     475             :             END IF
     476             : 
     477             :             CALL reorder(hybrid%nbasm(ikpt), hybrid%nbasp, atoms, hybrid%lcutm1, hybrid%maxlcutm1, hybrid%nindxm1,&
     478             :            &             2,&
     479           0 :            &             vec_c=vecout)
     480             : 
     481           0 :          END SUBROUTINE spmvec_noinvs
     482             : 
     483           0 :          SUBROUTINE reorder(nbasm, nbasp, atoms, lcutm, maxlcutm, nindxm, imode, vec_r, vec_c)
     484             :             USE m_types
     485             :             IMPLICIT NONE
     486             :             TYPE(t_atoms), INTENT(IN)   :: atoms
     487             : 
     488             :             ! - scalars -
     489             :             INTEGER, INTENT(IN)   ::  maxlcutm
     490             :             INTEGER, INTENT(IN)   ::  nbasm, nbasp
     491             :             INTEGER, INTENT(IN)   ::  imode
     492             : 
     493             :             ! - arrays -
     494             :             INTEGER, INTENT(IN)   ::  lcutm(atoms%ntype)
     495             :             INTEGER, INTENT(IN)   ::  nindxm(0:maxlcutm, atoms%ntype)
     496             :             REAL, INTENT(INOUT), OPTIONAL::  vec_r(nbasm)
     497             :             COMPLEX, INTENT(INOUT), OPTIONAL::  vec_c(nbasm)
     498             :             ! - local scalars -
     499             :             INTEGER               ::  itype, ieq
     500             :             INTEGER               ::  indx1, indx2
     501             :             INTEGER               ::  l
     502             :             INTEGER               ::  n, m
     503             :             LOGICAL               :: l_real
     504             :             ! - local arrays -
     505           0 :             REAL                  ::  vechlp_r(nbasm)
     506           0 :             COMPLEX               ::  vechlp_c(nbasm)
     507             : 
     508           0 :             l_real = PRESENT(vec_r)
     509             : 
     510           0 :             IF (imode /= 1 .and. imode /= 2) STOP 'reorder: imode equals neither 1 nor 2'
     511             : 
     512           0 :             if (l_real) THEN
     513           0 :                vechlp_r = vec_r
     514             :             else
     515           0 :                vechlp_c = vec_c
     516             :             end if
     517             : 
     518           0 :             IF (imode == 1) THEN
     519           0 :                indx1 = 0
     520           0 :                indx2 = 0
     521           0 :                DO itype = 1, atoms%ntype
     522           0 :                   DO ieq = 1, atoms%neq(itype)
     523           0 :                      DO l = 0, lcutm(itype)
     524           0 :                         DO m = -l, l
     525           0 :                            DO n = 1, nindxm(l, itype) - 1
     526           0 :                               indx1 = indx1 + 1
     527           0 :                               indx2 = indx2 + 1
     528           0 :                               if (l_real) THEN
     529           0 :                                  vec_r(indx1) = vechlp_r(indx2)
     530             :                               else
     531           0 :                                  vec_c(indx1) = vechlp_c(indx2)
     532             :                               endif
     533             :                            END DO
     534           0 :                            indx2 = indx2 + 1
     535             :                         END DO
     536             :                      END DO
     537             :                   END DO
     538             :                END DO
     539             : 
     540             :                indx2 = 0
     541           0 :                DO itype = 1, atoms%ntype
     542           0 :                   DO ieq = 1, atoms%neq(itype)
     543           0 :                      DO l = 0, lcutm(itype)
     544           0 :                         DO m = -l, l
     545           0 :                            indx1 = indx1 + 1
     546           0 :                            indx2 = indx2 + nindxm(l, itype)
     547           0 :                            if (l_real) THEN
     548           0 :                               vec_r(indx1) = vechlp_r(indx2)
     549             :                            else
     550           0 :                               vec_c(indx1) = vechlp_c(indx2)
     551             :                            endif
     552             : 
     553             :                         END DO
     554             :                      END DO
     555             :                   END DO
     556             :                END DO
     557           0 :             ELSE IF (imode == 2) THEN
     558           0 :                indx1 = 0
     559           0 :                indx2 = 0
     560           0 :                DO itype = 1, atoms%ntype
     561           0 :                   DO ieq = 1, atoms%neq(itype)
     562           0 :                      DO l = 0, lcutm(itype)
     563           0 :                         DO m = -l, l
     564           0 :                            DO n = 1, nindxm(l, itype) - 1
     565           0 :                               indx1 = indx1 + 1
     566           0 :                               indx2 = indx2 + 1
     567           0 :                               if (l_real) THEN
     568           0 :                                  vec_r(indx2) = vechlp_r(indx1)
     569             :                               else
     570           0 :                                  vec_c(indx2) = vechlp_c(indx1)
     571             :                               endif
     572             :                            END DO
     573           0 :                            indx2 = indx2 + 1
     574             :                         END DO
     575             :                      END DO
     576             :                   END DO
     577             :                END DO
     578             : 
     579             :                indx2 = 0
     580           0 :                DO itype = 1, atoms%ntype
     581           0 :                   DO ieq = 1, atoms%neq(itype)
     582           0 :                      DO l = 0, lcutm(itype)
     583           0 :                         DO m = -l, l
     584           0 :                            indx1 = indx1 + 1
     585           0 :                            indx2 = indx2 + nindxm(l, itype)
     586           0 :                            if (l_real) THEN
     587           0 :                               vec_r(indx2) = vechlp_r(indx1)
     588             :                            else
     589           0 :                               vec_c(indx2) = vechlp_c(indx1)
     590             :                            endif
     591             :                         END DO
     592             :                      END DO
     593             :                   END DO
     594             :                END DO
     595             :             END IF
     596             :             !IR must not be rearranged
     597             : 
     598           0 :          END SUBROUTINE reorder
     599             : 
     600             :       END MODULE

Generated by: LCOV version 1.13