Merge lp:~nickpapior/siesta/4.0-trs into lp:siesta/4.0

[Nick Papior]

Removed TRS in the diag2k code.

Also made it default to not apply TRS for the kpoint generation for non-colinear spin.

I ran the fe_noncol_kp test. And only saw that it converges faster using the new code.
The results are not exactly the same (even when mpirun -np 1), so I can't comment on this.

Probably this code needs revision from Ramon.

[lp:~nickpapior/siesta/4.0-trs]

556. By Nick Papior on 2018-01-31

Further cleaned interface to clear out an inconsistency in the code

[Alberto Garcia (albertog)]

It looks good to me, but I will wait for Ramón's opinion.

[Nick Papior (nickpapior)]

Super, I agree! :)

[Ramon Cuadrado (ramon-cuadrado)]

Hi Nick, Alberto,

   It seems that diag2k.F is ok. However, do you mind to send me the output that you got for fe_noncoll.fdf? The one that you comment that it does not give the same results? Thanks.

Cheers,
Ramón

review: Needs Information

[Nick Papior (nickpapior)]

I cannot attach the output, but I send it offsite.

[lp:~nickpapior/siesta/4.0-trs]

557. By Nick Papior on 2018-04-11

Merged 4.0

558. By Nick Papior on 2018-04-11

Stream-lined NC mulliken output with 4.1

559. By Nick Papior on 2018-04-12

Merged 4.0 into trs

560. By Nick Papior on 2018-04-12

Fixed documentation for TRS

561. By Nick Papior on 2018-04-12

Patched diag2[gk] to conform to sign-conventions in 4.1

This also lead to the discovery of a sign error in the phases.
So now everything in 4.0 is the same as 4.1 (for non-collinear).

I have checked fe_clust_noncollinear and FeO.
Both yield, exactly the same output before and after this commit.

[Nick Papior (nickpapior)]

Amended!

review: Needs Resubmitting

=== modified file 'Docs/siesta.tex'
--- Docs/siesta.tex	2018-04-09 13:45:49 +0000
+++ Docs/siesta.tex	2018-04-12 11:54:08 +0000
@@ -2779,13 +2779,14 @@
 \item[{\bf TimeReversalSymmetryForKpoints}] ({\it boolean}):
 \index{TimeReversalSymmetryForKpoints@{\bf TimeReversalSymmetryForKpoints}}
 
-If {\tt .true.}, the k-points in the BZ generated by the methods above
-are paired as (k,-k) and only one member of the pair is retained. This
+If {\tt .true.}, the $k$-points in the BZ generated by the methods above
+are paired as ($k$,$-k$) and only one member of the pair is retained. This
 symmetry is valid in the absence of external magnetic fields or
-spin-orbit interaction.
+non-collinear/spin-orbit interaction.
 
 {\it Default value:} {\tt .true.} unless the option {\bf SpinSpiral} 
-is used. In this case time-reversal-symmetry is broken explicitly.
+is used. In this case time-reversal-symmetry is broken explicitly. It is also default {\tt
+    .false.} for non-collinear calculations.
 
 \end{description}
 
 
=== modified file 'Src/Makefile'
--- Src/Makefile	2018-04-11 09:23:27 +0000
+++ Src/Makefile	2018-04-12 11:54:08 +0000
@@ -558,10 +558,10 @@
 kinefsm.o: alloc.o atmfuncs.o mneighb.o new_matel.o parallel.o parallelsubs.o
 kinefsm.o: precision.o
 kpoint_convert.o: precision.o sys.o units.o
-kpoint_grid.o: find_kgrid.o minvec.o parallel.o precision.o siesta_cml.o
-kpoint_grid.o: siesta_options.o sys.o units.o
-kpoint_pdos.o: find_kgrid.o minvec.o parallel.o precision.o siesta_options.o
-kpoint_pdos.o: sys.o units.o
+kpoint_grid.o: find_kgrid.o m_spin.o minvec.o parallel.o precision.o
+kpoint_grid.o: siesta_cml.o siesta_options.o sys.o units.o
+kpoint_pdos.o: find_kgrid.o m_spin.o minvec.o parallel.o precision.o
+kpoint_pdos.o: siesta_options.o sys.o units.o
 ksv.o: alloc.o atmfuncs.o densematrix.o ksvinit.o parallel.o precision.o sys.o
 ksvinit.o: alloc.o parallel.o precision.o sys.o
 listsc.o: alloc.o
 
=== modified file 'Src/diag2g.F'
--- Src/diag2g.F	2016-01-25 16:00:16 +0000
+++ Src/diag2g.F	2018-04-12 11:54:08 +0000
@@ -99,7 +99,7 @@
       logical               getD
       
 !     Internal variables .............................................
-      real(dp)                       psi(2,2,nuotot,2*nuo)
+      complex(dp) :: psi(2,nuotot,2*nuo)
       complex(dp),dimension(:,:,:,:), allocatable :: Haux, Saux
       complex(dp),dimension(:,:),     allocatable :: caux
       
@@ -196,8 +196,8 @@
         if (abs(qe).gt.occtol) then
           if (Node.eq.BNode) then
             do j = 1,nuotot
-              caux(1,j)=dcmplx( psi(1,1,j,iie), psi(2,1,j,iie)) ! c_{i,up}
-              caux(2,j)=dcmplx( psi(1,2,j,iie), psi(2,2,j,iie)) ! c_{i,dn}
+              caux(1,j) = psi(1,j,iie) ! c_{i,up}
+              caux(2,j) = psi(2,j,iie) ! c_{i,dn}
             enddo
           endif
 #ifdef MPI
@@ -211,45 +211,42 @@
               ind = listdptr(io) + j
               jo = listd(ind)
       
-!                 | ------- 1,1 -------     ------- 2,1 ------- |
-!                 | c_{j,up} c_{i,up}^*     c_{j,dn} c_{i,up)^* |
-!     D_{j,i} =   |                                             |
-!                 | ------- 1,2 -------     ------- 2,2 ------- |
-!                 | c_{j,up} c_{i,dn}^*     c_{j,dn} c_{i,dn)^* |
+!                | ------- 1,1 -------     ------- 1,2 ------- |
+!                | c_{j,up}^* c_{i,up}     c_{j,dn}^* c_{i,up) |
+!     D_{j,i} =  |                                             |
+!                | ------- 2,1 -------     ------- 2,2 ------- |
+!                | c_{j,up}^* c_{i,dn}     c_{j,dn}^* c_{i,dn) |
 !
 !------- 1,1 -----------------------------------------------------------
-              D11  = dconjg(caux(1,iio)) * caux(1,jo)
+              D11  = dconjg(caux(1,jo)) * caux(1,iio)
 !------- 2,2 -----------------------------------------------------------
-              D22  = dconjg(caux(2,iio)) * caux(2,jo)
+              D22  = dconjg(caux(2,jo)) * caux(2,iio)
+!------- 1,2 -----------------------------------------------------------
+              D12  = dconjg(caux(2,jo)) * caux(1,iio)
 !------- 2,1 -----------------------------------------------------------
-              D12  = dconjg(caux(1,iio)) * caux(2,jo)
-!------- 1,2 -----------------------------------------------------------
-              D21  = dconjg(caux(2,iio)) * caux(1,jo)
-      
-!------------ Density matrix has to be hermitian ----------------------
-              D11 = 0.5_dp * (D11 + dconjg(D11))
-              D22 = 0.5_dp * (D22 + dconjg(D22))
+              D21  = dconjg(caux(1,jo)) * caux(2,iio)
+
+              ! Make D12/D21 Hermitian
               D12 = 0.5_dp * (D12 + dconjg(D21))
-              D21 = dconjg(D12)
       
 !       Add contribution to density matrices of unit-cell orbitals
 !       ----------------------------------------------------------------
-!       | D11 = D_{j,i}(1)               D21 = D_{j,i}(3)-i D_{j,i}(4) |
-!       | D12 = D_{i,j}(3)+i D_{i,j}(4)  D22 = D_{j,i}(2)              |
-!       ----------------------------------------------------------------
-!       ----------------------------------------------------------------
-!       | D11 = Dnew(1)                  D21 = Dnew(3)-i Dnew(4)       |
-!       | D12 = Dnew(3)+i Dnew(4)        D22 = Dnew(2)                 |
+!       | D11 = D_{j,i}(1)               D12 = D_{j,i}(3)+i D_{j,i}(4) |
+!       | D12 = D_{i,j}(3)-i D_{i,j}(4)  D22 = D_{j,i}(2)              |
+!       ----------------------------------------------------------------
+!       ----------------------------------------------------------------
+!       | D11 = Dnew(1)                  D12 = Dnew(3)+i Dnew(4)       |
+!       | D21 = Dnew(3)-i Dnew(4)        D22 = Dnew(2)                 |
 !       ----------------------------------------------------------------
               Dnew(ind,1) = Dnew(ind,1) + dreal(D11) * qe
               Dnew(ind,2) = Dnew(ind,2) + dreal(D22) * qe
               Dnew(ind,3) = Dnew(ind,3) + dreal(D12) * qe
-              Dnew(ind,4) = Dnew(ind,4) + dimag(D12) * qe
+              Dnew(ind,4) = Dnew(ind,4) - dimag(D12) * qe
       
               Enew(ind,1) = Enew(ind,1) + dreal(D11) * ee
               Enew(ind,2) = Enew(ind,2) + dreal(D22) * ee
               Enew(ind,3) = Enew(ind,3) + dreal(D12) * ee
-              Enew(ind,4) = Enew(ind,4) + dimag(D12) * ee
+              Enew(ind,4) = Enew(ind,4) - dimag(D12) * ee
       
             enddo
           enddo
 
=== modified file 'Src/diag2k.F'
--- Src/diag2k.F	2016-01-25 16:00:16 +0000
+++ Src/diag2k.F	2018-04-12 11:54:08 +0000
@@ -109,13 +109,12 @@
       real(dp) kpoint(3,nk), wk(nk), xij(3,maxnh)
       real(dp) e1, e2, ef, eo(maxo*4,nk), qo(maxo*4,nk)
       real(dp) occtol, qtot, temp, Entropy
-      
+      complex(dp) :: psi(2,nuotot,nuo*2)
+
       logical               getD
       
 !     TEMPOS, INTERNAL VARIABLES etc.
       
-      real(dp) psi(2,2,nuotot,nuo*2)
-      
       integer  BNode, BTest, ie, ierror, iie, ik, ind, io, iio
       integer  iuo, j, jo, juo, neigneeded
       real(dp) ee, qe, t
@@ -123,9 +122,9 @@
 !     Haux(js,juo,is,iuo) = <js,juo|H|is,iuo>
 !     Indices is and js are for spin components
 !     Indices iuo and juo are for orbital components
-      complex(dp),dimension(:,:,:,:), allocatable :: Haux, Saux
+      complex(dp),dimension(:,:,:,:), allocatable, target :: Haux, Saux
       
-      complex(dp),dimension(:,:,:,:), allocatable :: Dkc, Ekc
+      complex(dp),dimension(:,:,:,:), pointer :: Dkc, Ekc
       complex(dp)                                 :: cicj
       complex(dp)                                 :: D11, D22, D12, D21
       complex(dp),dimension(:,:),     allocatable :: caux
@@ -144,12 +143,10 @@
       call memory('A','Z',MColl*nuotot*MColl*nuo,'diag2k')
       allocate(Saux(MColl,nuotot,MColl,nuo))
       call memory('A','Z',MColl*nuotot*MColl*nuo,'diag2k')
-      allocate(Dkc(MColl,nuotot,MColl,nuo))
-      call memory('A','Z',MColl*nuotot*MColl*nuo,'diag2k')
-      allocate(Ekc(MColl,nuotot,MColl,nuo))
-      call memory('A','Z',MColl*nuotot*MColl*nuo,'diag2k')
-      allocate(caux(MColl,MColl*nuotot))
-      call memory('A','Z',MColl*MColl*nuotot,'diag2k')
+      Dkc => Haux
+      Ekc => Saux
+      allocate(caux(MColl,nuotot))
+      call memory('A','Z',MColl*nuotot,'diag2k')
       
       eo = 0.0_dp
       
@@ -174,7 +171,7 @@
             kxij = kpoint(1,ik) * xij(1,ind) +
      .             kpoint(2,ik) * xij(2,ind) +
      .             kpoint(3,ik) * xij(3,ind)
-            kphs = cdexp(dcmplx(0.0_dp, 1.0_dp)*kxij)
+            kphs = cdexp(dcmplx(0.0_dp, -1.0_dp)*kxij)
       
             Saux(1,juo,1,iuo) = Saux(1,juo,1,iuo) + S(ind)   * kphs
             Saux(2,juo,2,iuo) = Saux(2,juo,2,iuo) + S(ind)   * kphs
@@ -243,7 +240,7 @@
             kxij = kpoint(1,ik) * xij(1,ind) +
      .             kpoint(2,ik) * xij(2,ind) +
      .             kpoint(3,ik) * xij(3,ind)
-            kphs = cdexp(dcmplx(0.0_dp, 1.0_dp)*kxij)
+            kphs = cdexp(dcmplx(0.0_dp, -1.0_dp)*kxij)
       
             Saux(1,juo,1,iuo) = Saux(1,juo,1,iuo) + S(ind)   * kphs
             Saux(2,juo,2,iuo) = Saux(2,juo,2,iuo) + S(ind)   * kphs
@@ -274,7 +271,7 @@
               kxij = kpoint(1,ik) * xij(1,ind) +
      .               kpoint(2,ik) * xij(2,ind) +
      .               kpoint(3,ik) * xij(3,ind)
-              kphs = cdexp(dcmplx(0.0_dp, 1.0_dp)*kxij)
+              kphs = cdexp(dcmplx(0.0_dp, -1.0_dp)*kxij)
       
               Saux(1,juo,1,iuo) = Saux(1,juo,1,iuo) + S(ind)   * kphs
               Saux(2,juo,2,iuo) = Saux(2,juo,2,iuo) + S(ind)   * kphs
@@ -290,14 +287,13 @@
      .               neigneeded,iscf,ierror)
         endif
       
-!-----------------------------------------------------------------------	     
+!-----------------------------------------------------------------------
 !     Store the products of eigenvectors in matrices Dk and Ek
 !     WARNING: Dk and Ek may be EQUIVALENCE'd to Haux and Saux
-!-----------------------------------------------------------------------	     
+!-----------------------------------------------------------------------
       
         Dkc  = dcmplx(0.0_dp,0.0_dp)
         Ekc  = dcmplx(0.0_dp,0.0_dp)
-        caux = dcmplx(0.0_dp,0.0_dp)
       
         BNode = 0
         iie = 0
@@ -311,8 +307,8 @@
           if (qe.gt.occtol) then
             if (Node.eq.BNode) then
               do j = 1,nuotot
-                caux(1,j)=dcmplx(psi(1,1,j,iie),psi(2,1,j,iie)) ! c_{i,up}
-                caux(2,j)=dcmplx(psi(1,2,j,iie),psi(2,2,j,iie)) ! c_{i,dn}
+                caux(1,j) = psi(1,j,iie) ! c_{i,up}
+                caux(2,j) = psi(2,j,iie) ! c_{i,dn}
               enddo
             endif
 #ifdef MPI
@@ -324,37 +320,33 @@
               call LocalToGlobalOrb(iuo,Node,Nodes,iio)
               do juo = 1,nuotot
       
-!     Relation:      c_{i,up} = caux(1,i)
-!                    c_{i,dn} = caux(2,i)
-!
-!                 | ----- 1,1 -------     ----- 2,1 ------- |
-!                 | c_{j,u} c_{i,u}^*     c_{j,u} c_{i,d)^* |
-!     D_{j,i} =   |                                         |
-!                 | ----- 1,2 -------     ----- 2,2 ------- |
-!                 | c_{j,d} c_{i,u}^*     c_{j,d} c_{i,d)^* |
-!
-!
-!                 | ----- 1,1 -------          ----- 2,1 -------       |
-!                 | D_{j,i}(1)                 D_{j,i}(3)+i D_{j,i}(4) |
-!             =   |                                                    |
-!                 | ----- 1,2 -------          ----- 2,2 -------       |
-!                 | D_{i,j}(3)-i D_{i,j}(4)    D_{j,i}(2)              |
-!
+!                | ------- 1,1 -------     ------- 1,2 ------- |
+!                | c_{j,up}^* c_{i,up}     c_{j,dn}^* c_{i,up) |
+!     D_{j,i} =  |                                             |
+!                | ------- 2,1 -------     ------- 2,2 ------- |
+!                | c_{j,up}^* c_{i,dn}     c_{j,dn}^* c_{i,dn) |
+!
+!
+!                | ------- 1,1 -------        ------- 1,2 -------     |
+!                | D_{j,i}(1)                 D_{j,i}(3)+i D_{j,i}(4) |
+!             =  |                                                    |
+!                | ------- 2,1 -------        ------- 2,2 -------     |
+!                | D_{i,j}(3)-i D_{i,j}(4)    D_{j,i}(2)              |
       
 !------- 1,1 -----------------------------------------------------------
-               cicj = dconjg(caux(1,iio)) * caux(1,juo)
+               cicj = dconjg(caux(1,juo)) * caux(1,iio)
                Dkc(1,juo,1,iuo) = Dkc(1,juo,1,iuo) + qe * cicj
                Ekc(1,juo,1,iuo) = Ekc(1,juo,1,iuo) + ee * cicj
 !------- 2,2 -----------------------------------------------------------
-               cicj = dconjg(caux(2,iio)) * caux(2,juo)
+               cicj = dconjg(caux(2,juo)) * caux(2,iio)
                Dkc(2,juo,2,iuo) = Dkc(2,juo,2,iuo) + qe * cicj
                Ekc(2,juo,2,iuo) = Ekc(2,juo,2,iuo) + ee * cicj
-!------- 2,1 -----------------------------------------------------------
-               cicj = dconjg(caux(1,iio)) * caux(2,juo)
+!------- 1,2 -----------------------------------------------------------
+               cicj = dconjg(caux(2,juo)) * caux(1,iio)
                Dkc(1,juo,2,iuo) = Dkc(1,juo,2,iuo) + qe * cicj
                Ekc(1,juo,2,iuo) = Ekc(1,juo,2,iuo) + ee * cicj
-!------- 1,2 -----------------------------------------------------------
-               cicj = dconjg(caux(2,iio)) * caux(1,juo)
+!------- 2,1 -----------------------------------------------------------
+               cicj = dconjg(caux(1,juo)) * caux(2,iio)
                Dkc(2,juo,1,iuo) = Dkc(2,juo,1,iuo) + qe * cicj
                Ekc(2,juo,1,iuo) = Ekc(2,juo,1,iuo) + ee * cicj
               enddo
@@ -383,47 +375,31 @@
      .           kpoint(3,ik) * xij(3,ind)
           kphs = cdexp(dcmplx(0.0_dp,-1.0_dp)*kxij)
       
-!         Average k and -k solutions because time-reversal symetry
-          D11 = 0.5_dp * (Dkc(1,juo,1,iuo) * kphs
-     .        +    dconjg(Dkc(1,juo,1,iuo) * kphs))
-      
-          D22 = 0.5_dp * (Dkc(2,juo,2,iuo) * kphs
-     .        +    dconjg(Dkc(2,juo,2,iuo) * kphs))
-          
-          D12 = 0.5_dp * (Dkc(1,juo,2,iuo) * kphs
-     .        +    dconjg(Dkc(2,juo,1,iuo) * kphs))
-      
-          D21 = 0.5_dp * (Dkc(2,juo,1,iuo) * kphs
-     .        +    dconjg(Dkc(1,juo,2,iuo) * kphs))
-      
+          D11 = Dkc(1,juo,1,iuo) * kphs
+          D22 = Dkc(2,juo,2,iuo) * kphs
+          D12 = Dkc(1,juo,2,iuo) * kphs
+          D21 = Dkc(2,juo,1,iuo) * kphs
+
+          ! Make D12/D21 Hermitian
           D12 = 0.5_dp * (D12 + dconjg(D21))
-          D21 = dconjg(D12)
-      
+
           Dnew(ind,1) = Dnew(ind,1) + dreal(D11)
           Dnew(ind,2) = Dnew(ind,2) + dreal(D22)
           Dnew(ind,3) = Dnew(ind,3) + dreal(D12)
-          Dnew(ind,4) = Dnew(ind,4) + dimag(D12)
-      
-!         Average k and -k solutions because time-reversal symetry
-          D11 = 0.5_dp * (Ekc(1,juo,1,iuo) * kphs
-     .        +    dconjg(Ekc(1,juo,1,iuo) * kphs))
-      
-          D22 = 0.5_dp * (Ekc(2,juo,2,iuo) * kphs
-     .        +    dconjg(Ekc(2,juo,2,iuo) * kphs))
-          
-          D12 = 0.5_dp * (Ekc(1,juo,2,iuo) * kphs
-     .        +    dconjg(Ekc(2,juo,1,iuo) * kphs))
-      
-          D21 = 0.5_dp * (Ekc(2,juo,1,iuo) * kphs
-     .        +    dconjg(Ekc(1,juo,2,iuo) * kphs))
-      
+          Dnew(ind,4) = Dnew(ind,4) - dimag(D12)
+
+          D11 = Ekc(1,juo,1,iuo) * kphs
+          D22 = Ekc(2,juo,2,iuo) * kphs
+          D12 = Ekc(1,juo,2,iuo) * kphs
+          D21 = Ekc(2,juo,1,iuo) * kphs
+
+          ! Make D12/D21 Hermitian
           D12 = 0.5_dp * (D12 + dconjg(D21))
-          D21 = dconjg(D12)
       
           Enew(ind,1) = Enew(ind,1) + dreal(D11)
           Enew(ind,2) = Enew(ind,2) + dreal(D22)
           Enew(ind,3) = Enew(ind,3) + dreal(D12) 
-          Enew(ind,4) = Enew(ind,4) + dimag(D12) 
+          Enew(ind,4) = Enew(ind,4) - dimag(D12) 
       
         enddo
         enddo
@@ -436,30 +412,9 @@
       deallocate(Haux)
       call memory('D','Z',MColl*nuotot*MColl*nuo,'diag2k')
       deallocate(Saux)
-      call memory('D','Z',MColl*nuotot*MColl*nuo,'diag2k')
-      deallocate(Dkc)
-      call memory('D','Z',MColl*nuotot*MColl*nuo,'diag2k')
-      deallocate(Ekc)
-      call memory('D','Z',MColl*MColl*nuotot,'diag2k')
+      call memory('D','Z',MColl*nuotot,'diag2k')
       deallocate(caux)
       
-!        write(160+node,*)'ITERATION'
-!        write(170+node,*)'ITERATION'
-!        do iuo = 1,nuo
-!          do j = 1,numd(iuo)
-!            ind = listdptr(iuo) + j
-!            jo = listd(ind)
-!            juo = indxuo(jo)
-!             write(160+node,'(3i5,4f15.9)')iuo,juo,jo,Dnew(ind,1:4)
-!             write(170+node,'(3i5,4f15.9)')iuo,juo,jo,
-!     .    (Dnew(ind,1)+Dnew(ind,2))*0.5,
-!     .    (Dnew(ind,1)-Dnew(ind,2))*0.5,Dnew(ind,3:4)
-!          enddo
-!        enddo
-!      
-!      call die('maslyuk')
-!***********************************************************************      
       return
       end subroutine diag2k
-!***********************************************************************      
       
 
=== modified file 'Src/kpoint_grid.F90'
--- Src/kpoint_grid.F90	2017-10-05 10:01:03 +0000
+++ Src/kpoint_grid.F90	2018-04-12 11:54:08 +0000
@@ -44,6 +44,7 @@
   USE parallel, only  : Node
   USE fdf, only       : fdf_defined, fdf_get
   USE m_find_kgrid, only : find_kgrid
+  use m_spin, only: MColl
 
     implicit none
     real(dp) :: ucell(3,3)
@@ -52,11 +53,12 @@
     if (scf_kgrid_first_time) then
        nullify(kweight,kpoint)
        spiral = fdf_defined('SpinSpiral')
-          ! Allow the user to control the use of time-reversal-symmetry
-          ! By default, it is on, except for "spin-spiral" calculations
-       time_reversal_symmetry = fdf_get(             &
-            "TimeReversalSymmetryForKpoints",   &
-            (.not. spiral))
+       ! Allow the user to control the use of time-reversal-symmetry
+       ! By default, it is on, except for "spin-spiral" calculations
+       ! and/or non-colinear calculations
+       time_reversal_symmetry = fdf_get( &
+            "TimeReversalSymmetryForKpoints", &
+            (.not. spiral) .and. MColl == 1)
        call setup_scf_kscell(ucell, firm_displ)
 
        scf_kgrid_first_time = .false.
 
=== modified file 'Src/kpoint_pdos.F90'
--- Src/kpoint_pdos.F90	2017-10-05 10:01:03 +0000
+++ Src/kpoint_pdos.F90	2018-04-12 11:54:08 +0000
@@ -43,6 +43,7 @@
   USE parallel, only  : Node
   USE fdf, only       : fdf_defined, fdf_get
   USE m_find_kgrid, only : find_kgrid
+  use m_spin, only: MColl
 
   implicit none
   real(dp), intent(in)  :: ucell(3,3)
@@ -52,11 +53,12 @@
   if (pdos_kgrid_first_time) then
     nullify(kweight_pdos,kpoints_pdos)
     spiral = fdf_defined('SpinSpiral')
-      ! Allow the user to control the use of time-reversal-symmetry
-      ! By default, it is on, except for "spin-spiral" calculations
-      time_reversal_symmetry = fdf_get(             &
-           "TimeReversalSymmetryForKpoints",   &
-           (.not. spiral))
+    ! Allow the user to control the use of time-reversal-symmetry
+    ! By default, it is on, except for "spin-spiral" calculations
+    ! and/or non-colinear calculations
+    time_reversal_symmetry = fdf_get( &
+         "TimeReversalSymmetryForKpoints", &
+         (.not. spiral) .and. MColl == 1)
 
     call setup_pdos_kscell(ucell, firm_displ)
 
 
=== modified file 'Src/mulliken.F'
--- Src/mulliken.F	2018-04-11 23:32:54 +0000
+++ Src/mulliken.F	2018-04-12 11:54:08 +0000
@@ -513,7 +513,6 @@
         end if ! nspin .ge. 4
       end if
 
-C ...................
 
 #ifdef DEBUG
       call write_debug( '    POS mulliken' )
 
=== modified file 'version.info'
--- version.info	2018-04-11 23:32:54 +0000
+++ version.info	2018-04-12 11:54:08 +0000
@@ -1,1 +1,1 @@
-siesta-4.0--563
+siesta-4.0--563--trs-7