Merge branch 'develop' of ifflinux.iff.kfa-juelich.de:fleur into develop

.gitlab-ci.yml

+stages:
+   - build
+   - test
+   - deploy
+   - build-pgi
+   - test-pgi
+
 build-gfortran:
-  image: ubuntu:16.10
+  image: iffregistry.fz-juelich.de/fleur/fleur:latest
   stage: build
   cache: 
      paths:
        - build
   script:
-    - apt-get update -qq && apt-get install -y -qq build-essential gfortran gfortran-6 g++-6 gcc-6 binutils-gold git cmake libxml2-dev liblapack-dev libblas-dev doxygen libopenmpi-dev libscalapack-mpi-dev
     - cd /builds/fleur/fleur; ./configure.sh GITLAB; cd build; make 
 test-gfortran:
-  image: ubuntu:16.10
+  image: iffregistry.fz-juelich.de/fleur/fleur:latest
   stage: test
   cache:
      paths:
        - build
   script:
-    - apt-get update -qq && apt-get install -y -qq build-essential gfortran gfortran-6 g++-6 gcc-6 binutils-gold git cmake libxml2-dev liblapack-dev libblas-dev libopenmpi-dev libscalapack-mpi-dev
     - ulimit -s unlimited ;export juDFT_MPI="mpirun -n 2 --allow-run-as-root ";cd /builds/fleur/fleur/build;ctest
 
+  
+pages:
+  image: iffregistry.fz-juelich.de/fleur/fleur:latest
+  stage: deploy
+  cache:
+     paths:
+       - build
+  script:
+    - cd /builds/fleur/fleur/build ; make doc
+    - mv docs/html/ ../public/
+  artifacts:
+    paths:
+      - public
+      - build/fleur
+      - build/fleur_MPI
+      - build/inpgen
+ 
 
-#build-pgi:
-#  image: hopobcn/pgi-ce
-#  stage: build
-#  cache: 
-#     paths:
-#       - build.debug
-#  script:
-#    - apt-get update -qq && apt-get install -y -qq cmake libxml2-dev
-#    - cd /builds/fleur/fleur; ./configure.sh AUTO debug; cd build.debug; make 
-#test-pgi:
-#  image: hopobcn/pgi-ce
-#  stage: test
-#  cache:
-#     paths:
-#       - build.debug
-#  script:
-#    - apt-get update -qq && apt-get install -y -qq cmake libxml2-dev 
-#    - cd /builds/fleur/fleur/build.debug;ctest
+build-pgi:
+  image: iffregistry.fz-juelich.de/fleur/fleur:pgi
+  stage: build-pgi
+  cache: 
+     paths:
+       - build.debug
+  script:
+    - apt-get update -qq && apt-get install -y -qq cmake libxml2-dev
+    - cd /builds/fleur/fleur; ./configure.sh AUTO debug; cd build.debug; make 
+  only:
+    - schedules
+    - triggers
+    - web
+ 
+    
+test-pgi:
+  image: iffregistry.fz-juelich.de/fleur/fleur:pgi
+  stage: test-pgi
+  cache:
+     paths:
+       - build.debug
+  script:
+    - apt-get update -qq && apt-get install -y -qq cmake libxml2-dev 
+    - cd /builds/fleur/fleur/build.debug;ctest
+  only:
+     - schedules
+     - web
+     - triggers

docs/ContributorsGuide.md

+Contributors guide
+======================================
+
+Everyone is very welcome to contribute to the enhancement of FLEUR.
+Please use the [gitlab service] (https://iffgit.fz-juelich.de/fleur/fleur) to obtain the
+latest development version of FLEUR.
+
+##Coding rules for FLEUR:
+In no particular order we try to collect items to consider when writing code for FLEUR
+
+- Instead of 'stop' use calls to judft_error, judft_warn, judft_end
+
+- Do not read and write any files. Files are neither replacements for common-blocks nor a
+storage for status variables.
+ Only exceptions: 
++ you create nice IO subroutines in the subdirectory io
++ you write to the typical FLEUR output files
+

docs/Main.md

+\mainpage Welcome to the developers documentation of FLEUR
+
+These pages are complementing the <A HREF="http://www.flapw.de/pm/index.php?n=User-Documentation.Documentation"> User documentation </A> in which the actual use of FLEUR is documented.
+
+Here we collect information useful for developers or advanced users that actually will modify the source code.
+
+## GitLab
+The development process is performed using gitlab. You can access the  <A HREF="https://iffgit.fz-juelich.de/fleur/fleur"> main gitlab page here.</A>
+
+If you checkout the code please be aware that there are several branches. The release branch contains the code of the last release published on  <A HREF="http://www.flapw.de/pm/index.php?n=FLEUR.Downloads"> the FLEUR webpage</A>. You can *not* push to this branch directly.
+You probably want to use the development branch to insert your changes. If your changes are large, it might be a good idea to create your own branch first.
+
+## More info
+In the menu to the left you will find information on various topics as well as the links to the source-code documentation.
+
+

docs/Tips+Tricks.md

+Useful info for developers
+==============================================
+
+## Using fleur with the HDF5 library and debugging it with valgrind
+
+HDF5 has to be built with the same compiler that is also used to compile fleur. If adapted the following commands can be used to compile a HDF5 library for fleur:
+
++ 'FC=/usr/local/intel/impi/4.0.3.008/intel64/bin/mpiifort CC=/usr/local/intel/impi/4.0.3.008/intel64/bin/mpiicc CXX=/usr/local/intel/impi/4.0.3.008/intel64/bin/mpiicc ./configure --enable-fortran --enable-fortran2003 --enable-parallel --enable-using-memchecker --enable-clear-file-buffers'
++ 'make'
++ 'make install'
++ 'make check' (optional)
+
+Note:
+
++ The paths have to be adjusted such that that compiler is used which is also used to compile fleur.
++ The parallel switch is not needed for every calculation: Only for parallel calculations in which HDF5 is also used for the eigenvector IO.
++ The last two command line switches in the configure command turn on initializations of irrelevant array parts in HDF5. If valgrind is not needed it is probably the better choice to leave them away. If left away valgrind will complain about missing initializations in the HDF5 library.
++ valgrind gives partially strange behavior if used together with the intel compiler. It would be better to use it together with gfortran.
++ At the moment HDF5 is needed in version 1.8.*. Usage of version 1.10.* yields some problems.
+
+Furthermore to configure and start fleur with HDF5 the following has to be done:
+
++ In your .bashrc the HDF5 library has to be added to the LD_LIBRARY_PATH. This implies a line like 'export LD_LIBRARY_PATH=$LD_LIBRARY_PATH:~/hdf5/current/hdf5/lib'
++ configure fleur with some line like 'CMAKE_Fortran_FLAGS="-I~/hdf5/current/hdf5/include" FLEUR_LIBRARIES="-L~/hdf5/current/hdf5/lib;-lhdf5_fortran;-lhdf5" ./fleur/configure.sh IFF'
\ No newline at end of file

eigen/eigen.F90

@@ -17,7 +17,7 @@ CONTAINS
     ! nmat, nbasfcn                                   including LO's
     !        g. bihlmayer '96
     !**********************************************************************
-    USE m_constants, ONLY : pi_const,sfp_const
+    USE m_constants
     USE m_types
     USE m_lodpot
     USE m_tlmplm
@@ -75,7 +75,7 @@ CONTAINS
     !+odim
     !     ..
     !     .. Local Scalars ..
-    INTEGER jsp,nk,nred,ne_all,n_u_in,ne_found
+    INTEGER jsp,nk,nred,ne_all,ne_found
     INTEGER ne,matsize  ,nrec,lh0
     INTEGER nspins,isp,i,j,err
     INTEGER mlotot,mlolotot
@@ -83,13 +83,11 @@ CONTAINS
     
     !     ..
     !     .. Local Arrays ..
-    INTEGER, PARAMETER :: lmaxb=3
     INTEGER, ALLOCATABLE :: matind(:,:),kveclo(:)
     INTEGER, ALLOCATABLE :: nv2(:)
     REAL,    ALLOCATABLE :: bkpt(:)
     REAL,    ALLOCATABLE :: eig(:)
 
-    COMPLEX, ALLOCATABLE :: vs_mmp(:,:,:,:)
     TYPE(t_tlmplm)  :: td
     TYPE(t_usdus)   :: ud
     TYPE(t_lapw)    :: lapw
@@ -203,7 +201,10 @@ CONTAINS
           matsize = ((DIMENSION%nbasfcn+1)/2)*DIMENSION%nbasfcn
        ENDIF
     ENDIF
-    IF (matsize<2) CALL judft_error("Wrong size of matrix",calledby="eigen",hint="Your basis might be too large or the parallelization fail or ??")
+    IF (matsize<2) THEN
+       CALL judft_error("Wrong size of matrix",calledby="eigen",&
+                        hint="Your basis might be too large or the parallelization fail or ??")
+    END IF
     ne = DIMENSION%neigd
 
     eig_id=open_eig(&
@@ -247,12 +248,11 @@ CONTAINS
     !
     !  ..
     !  LDA+U
-    n_u_in=atoms%n_u
     IF ((atoms%n_u.GT.0)) THEN
-       ALLOCATE( vs_mmp(-lmaxb:lmaxb,-lmaxb:lmaxb,atoms%n_u,input%jspins) )
-       CALL u_setup(sym,atoms,lmaxb,sphhar,input, enpara%el0(0:,:,:),v%mt,mpi, vs_mmp,results)
+       ALLOCATE( v%mmpMat(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,atoms%n_u,input%jspins) )
+       CALL u_setup(sym,atoms,sphhar,input, enpara%el0(0:,:,:),v,mpi,results)
     ELSE
-       ALLOCATE( vs_mmp(-lmaxb:-lmaxb,-lmaxb:-lmaxb,1,2) )
+       ALLOCATE( v%mmpMat(-lmaxU_const:-lmaxU_const,-lmaxU_const:-lmaxU_const,1,2) )
     ENDIF
     !
     !--->    loop over k-points: each can be a separate task
@@ -282,7 +282,7 @@ CONTAINS
        ENDIF
        lh0=1
        CALL tlmplm(sphhar,atoms,DIMENSION,enpara, jsp,1,mpi, v%mt(1,0,1,jsp),lh0,input, td,ud)
-       IF (input%l_f) CALL write_tlmplm(td,vs_mmp,atoms%n_u>0,1,jsp,input%jspins)
+       IF (input%l_f) CALL write_tlmplm(td,v%mmpMat,atoms%n_u>0,1,jsp,input%jspins)
        CALL timestop("tlmplm")
 
        !---> pk non-collinear
@@ -293,7 +293,7 @@ CONTAINS
           isp = 2
           CALL timestart("tlmplm")
           CALL tlmplm(sphhar,atoms,DIMENSION,enpara,isp,isp,mpi, v%mt(1,0,1,isp),lh0,input, td,ud)
-          IF (input%l_f) CALL write_tlmplm(td,vs_mmp,atoms%n_u>0,2,2,input%jspins)
+          IF (input%l_f) CALL write_tlmplm(td,v%mmpMat,atoms%n_u>0,2,2,input%jspins)
           CALL timestop("tlmplm")
        ENDIF
        !
@@ -332,7 +332,7 @@ CONTAINS
           IF (.NOT.l_wu) THEN
              CALL timestart("MT Hamiltonian&Overlap")
              CALL hsmt(DIMENSION,atoms,sphhar,sym,enpara, mpi%SUB_COMM,mpi%n_size,mpi%n_rank,jsp,input,mpi,&
-                  lmaxb, noco,cell, lapw, bkpt,v%mt, vs_mmp, oneD,ud, kveclo,td,l_real,hamOvlp)
+                  lmaxU_const, noco,cell, lapw, bkpt,v%mt,v%mmpMat, oneD,ud, kveclo,td,l_real,hamOvlp)
              CALL timestop("MT Hamiltonian&Overlap")
           ENDIF
           !
@@ -423,7 +423,7 @@ ENDIF
        DEALLOCATE (td%ind,td%tuulo,td%tdulo)
        DEALLOCATE (td%tuloulo)
     END DO ! spin loop ends
-    DEALLOCATE( vs_mmp )
+    DEALLOCATE(v%mmpMat)
     DEALLOCATE (matind)
     IF (l_real) THEN
        DEALLOCATE(hamOvlp%a_r,hamOvlp%b_r)
@@ -439,8 +439,6 @@ ENDIF
     CALL MPI_BARRIER(mpi%MPI_COMM,ierr)
 #endif
     !IF (hybrid%l_hybrid.OR.hybrid%l_calhf) CALL close_eig(eig_id)
-    atoms%n_u=n_u_in
-
 
     IF( input%jspins .EQ. 1 .AND. hybrid%l_hybrid ) THEN
        results%te_hfex%valence = 2*results%te_hfex%valence

force/force_a21.F90

@@ -54,7 +54,6 @@ CONTAINS
     COMPLEX, INTENT (IN) :: cveccof(3,-atoms%llod:atoms%llod,nobd,atoms%nlod,atoms%nat)
     !     ..
     !     .. Local Scalars ..
-    INTEGER, PARAMETER :: lmaxb=3
     COMPLEX dtd,dtu,utd,utu
     INTEGER lo, mlotot, mlolotot, mlot_d, mlolot_d
     INTEGER i,ie,im,in,l1,l2,ll1,ll2,lm1,lm2,m1,m2,n,natom,m,i_u
@@ -88,7 +87,7 @@ CONTAINS
          a21(3,atoms%nat),b4(3,atoms%nat),tlmplm%ind(0:DIMENSION%lmd,0:DIMENSION%lmd,atoms%ntype,1) )
     !
     IF(atoms%n_u.GT.0) THEN
-       ALLOCATE(v_mmp(-lmaxb:lmaxb,-lmaxb:lmaxb,atoms%n_u))
+       ALLOCATE(v_mmp(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,atoms%n_u))
        v_mmp = CMPLX(0.0,0.0)
        CALL read_tlmplm_vs_mmp(jsp,atoms%n_u,v_mmp)
     END IF
@@ -197,7 +196,7 @@ CONTAINS
                cveccof, tlmplm,usdus, a21)
 
           IF ((atoms%n_u.GT.0).AND.(i_u.LE.atoms%n_u)) THEN
-             CALL force_a21_U(nobd,atoms,lmaxb,i_u,n,jsp,we,ne,&
+             CALL force_a21_U(nobd,atoms,i_u,n,jsp,we,ne,&
                               usdus,v_mmp,acof,bcof,ccof,&
                               aveccof,bveccof,cveccof, a21)
           END IF

force/force_a21_U.f90

 MODULE m_forcea21U
 CONTAINS
-  SUBROUTINE force_a21_U(nobd,atoms,lmaxb,i_u,itype,isp,we,ne,&
+  SUBROUTINE force_a21_U(nobd,atoms,i_u,itype,isp,we,ne,&
        usdus,v_mmp, acof,bcof,ccof,aveccof,bveccof,cveccof, a21)
     !
     !***********************************************************************
@@ -9,6 +9,7 @@ CONTAINS
     ! Comp.Phys.Comm. 179 (2008) 784-790
     !***********************************************************************
     !
+    USE m_constants
     USE m_types
     IMPLICIT NONE
 
@@ -17,13 +18,13 @@ CONTAINS
     !     ..
     !     .. Scalar Arguments ..
     INTEGER, INTENT (IN)    :: nobd   
-    INTEGER, INTENT (IN)    :: itype,isp,ne,lmaxb
+    INTEGER, INTENT (IN)    :: itype,isp,ne
     INTEGER, INTENT (INOUT) :: i_u ! on input: index for the first U for atom type "itype or higher"
                                    ! on exit: index for the first U for atom type "itype+1 or higher"
     !     ..
     !     .. Array Arguments ..
     REAL,    INTENT (IN) :: we(nobd) 
-    COMPLEX, INTENT (IN) :: v_mmp(-lmaxb:lmaxb,-lmaxb:lmaxb,atoms%n_u)
+    COMPLEX, INTENT (IN) :: v_mmp(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,atoms%n_u)
     COMPLEX, INTENT (IN) :: acof(:,0:,:)!(nobd,0:dimension%lmd,atoms%nat)
     COMPLEX, INTENT (IN) :: bcof(:,0:,:)!(nobd,0:dimension%lmd,atoms%nat)
     COMPLEX, INTENT (IN) :: ccof(-atoms%llod:atoms%llod,nobd,atoms%nlod,atoms%nat)

global/constants.f

@@ -10,6 +10,7 @@
       INTEGER,          PARAMETER :: noState_const = 0
       INTEGER,          PARAMETER :: coreState_const = 1
       INTEGER,          PARAMETER :: valenceState_const = 2
+      INTEGER,          PARAMETER :: lmaxU_const = 3
       REAL,             PARAMETER :: pi_const=3.1415926535897932
       REAL,             PARAMETER :: tpi_const=2.*3.1415926535897932
       REAL,             PARAMETER :: fpi_const=4.*3.1415926535897932

global/types.F90

@@ -262,7 +262,7 @@ MODULE m_types
      INTEGER,          ALLOCATABLE :: igrd(:)
      INTEGER,          ALLOCATABLE :: krla(:)
      LOGICAL,          ALLOCATABLE :: relcor(:)
-     !lda_u information(ntype)
+     !lda_u information(size: atoms%n_u)
      TYPE(t_utype),ALLOCATABLE::lda_u(:)
      INTEGER,ALLOCATABLE :: relax(:,:) !<(3,ntype)
      INTEGER, ALLOCATABLE :: nflip(:) !<flip magnetisation of this atom
@@ -893,11 +893,21 @@ MODULE m_types
 
 
   TYPE t_potden
+     ! General variables and arrays
      INTEGER             :: iter
      COMPLEX,ALLOCATABLE :: pw(:,:)
      REAL,ALLOCATABLE    :: mt(:,:,:,:)
      REAL,ALLOCATABLE    :: vacz(:,:,:)
      COMPLEX,ALLOCATABLE :: vacxy(:,:,:,:)
+
+     ! For density only (noco case)
+     COMPLEX, ALLOCATABLE :: cdom(:)
+     COMPLEX, ALLOCATABLE :: cdomvz(:,:)
+     COMPLEX, ALLOCATABLE :: cdomvxy(:,:,:)
+
+     ! For density matrix and associated potential matrix
+     COMPLEX, ALLOCATABLE :: mmpMat(:,:,:,:)
+
      !this type contains two init routines that should be used to allocate
      !memory. You can either specify the datatypes or give the dimensions as integers
      !See implementation below!

io/r_inpXML.F90

@@ -119,7 +119,7 @@ SUBROUTINE r_inpXML(&
   REAL               :: tauTemp(3,48)
   REAL               :: bk(3)
   LOGICAL            :: flipSpin, l_eV, invSym, l_qfix, relaxX, relaxY, relaxZ
-  LOGICAL            :: l_vca, coreConfigPresent, l_enpara, l_orbcomp
+  LOGICAL            :: l_vca, coreConfigPresent, l_enpara, l_orbcomp, tempBool
   REAL               :: magMom, radius, logIncrement, qsc(3), latticeScale, dr
   REAL               :: aTemp, zp, rmtmax, sumWeight, ldau_u(4), ldau_j(4), tempReal
   REAL               :: weightScale, eParamUp, eParamDown
@@ -1794,12 +1794,15 @@ SUBROUTINE r_inpXML(&
      END IF
 
      IF (numberNodes.EQ.1) THEN
-        coreSpecInput%verb = evaluateFirstBoolOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@verbose'))
+        coreSpecInput%verb = 0
+        tempBool = evaluateFirstBoolOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@verbose'))
+        IF(tempBool) coreSpecInput%verb = 1
         coreSpecInput%ek0 = evaluateFirstOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@eKin'))
         coreSpecInput%atomType = evaluateFirstIntOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@atomType'))
         coreSpecInput%lx = evaluateFirstIntOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@lmax'))
         coreSpecInput%edge = TRIM(ADJUSTL(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@edgeType')))
-        coreSpecInput%edgeidx = evaluateFirstIntOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@edgeIndex'))
+        coreSpecInput%edgeidx(:) = 0
+        coreSpecInput%edgeidx(1) = evaluateFirstIntOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@edgeIndex'))
         coreSpecInput%emn = evaluateFirstOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@eMin'))
         coreSpecInput%emx = evaluateFirstOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@eMax'))
         tempInt = evaluateFirstIntOnly(xmlGetAttributeValue(TRIM(ADJUSTL(xPathA))//'/@numPoints'))

ldau/u_setup.f90

@@ -17,12 +17,13 @@ MODULE m_usetup
   !     Extension to multiple U per atom type  G.M. 2017              |
   !-------------------------------------------------------------------+
 CONTAINS
-  SUBROUTINE u_setup(sym,atoms,lmaxb,sphhar, input,el,vr,mpi, vs_mmp,results,number)
+  SUBROUTINE u_setup(sym,atoms,sphhar, input,el,pot,mpi,results)
     USE m_umtx
     USE m_uj2f
     USE m_nmat_rot
     USE m_vmmp
     USE m_types
+    USE m_constants
     IMPLICIT NONE
     TYPE(t_sym),INTENT(IN)          :: sym
     TYPE(t_results),INTENT(INOUT)   :: results
@@ -30,14 +31,9 @@ CONTAINS
     TYPE(t_input),INTENT(IN)        :: input
     TYPE(t_sphhar),INTENT(IN)       :: sphhar
     TYPE(t_atoms),INTENT(INOUT)     :: atoms !n_u might be modified if no density matrix is found
+    TYPE(t_potden),INTENT(INOUT)    :: pot
 
-    ! ... Arguments ...
-    INTEGER, INTENT (IN) :: lmaxb  
-
-    REAL,    INTENT (IN) :: el(0:,:,:) !(0:atoms%lmaxd,ntype,jspd)
-    REAL,    INTENT (IN) :: vr(:,0:,:,:) !(atoms%jmtd,0:sphhar%nlhd,ntype,jspd)
-    COMPLEX,INTENT (OUT)::vs_mmp(-lmaxb:lmaxb,-lmaxb:lmaxb,atoms%n_u,input%jspins)
-    INTEGER,INTENT(IN),OPTIONAL::number
+    REAL,    INTENT(IN)           :: el(0:,:,:) !(0:atoms%lmaxd,ntype,jspd)
     ! ... Local Variables ...
     INTEGER itype,ispin,j,k,l,jspin,urec,i_u
     INTEGER noded,nodeu,ios,lty(atoms%n_u)
@@ -53,9 +49,6 @@ CONTAINS
     ! look, whether density matrix exists already:
     !
     filename="n_mmp_mat"
-    IF (PRESENT(number)) THEN
-       WRITE(filename,"('n_mmp_mat.',i0)") number
-    ENDIF
     INQUIRE (file=filename,exist=n_mmp_exist)
     IF (n_mmp_exist.AND.atoms%n_u>0) THEN
        !
@@ -66,21 +59,21 @@ CONTAINS
        !
        ! set up e-e- interaction matrix
        !
-       ALLOCATE ( u(-lmaxb:lmaxb,-lmaxb:lmaxb,&
-            -lmaxb:lmaxb,-lmaxb:lmaxb,atoms%n_u,input%jspins) )
+       ALLOCATE (u(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,&
+                   -lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,atoms%n_u,input%jspins))
        DO ispin = 1, 1 ! input%jspins
           f0(:,1) = (f0(:,1) + f0(:,input%jspins) ) / 2
           f2(:,1) = (f2(:,1) + f2(:,input%jspins) ) / 2
           f4(:,1) = (f4(:,1) + f4(:,input%jspins) ) / 2
           f6(:,1) = (f6(:,1) + f6(:,input%jspins) ) / 2
-          CALL umtx(atoms,lmaxb,f0(1,ispin),&
+          CALL umtx(atoms,f0(1,ispin),&
                     f2(1,ispin),f4(1,ispin),f6(1,ispin),&
-                    u(-lmaxb,-lmaxb,-lmaxb,-lmaxb,1,ispin) )
+                    u(-lmaxU_const,-lmaxU_const,-lmaxU_const,-lmaxU_const,1,ispin) )
        END DO
        !
        ! read density matrix
        !
-       ALLOCATE (ns_mmp(-lmaxb:lmaxb,-lmaxb:lmaxb,atoms%n_u,input%jspins))
+       ALLOCATE (ns_mmp(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,atoms%n_u,input%jspins))
        OPEN (69,file=filename,status='unknown',form='formatted')
        READ (69,9000) ns_mmp
 9000   FORMAT(7f20.13)
@@ -110,8 +103,8 @@ CONTAINS
        !
        ! calculate potential matrix and total energy correction
        !
-       CALL v_mmp(sym,atoms,input%jspins,lmaxb,ns_mmp,u,f0,f2,&
-                  vs_mmp,results)
+       CALL v_mmp(sym,atoms,input%jspins,ns_mmp,u,f0,f2,&
+                  pot%mmpMat,results)
        IF (mpi%irank.EQ.0) THEN
           DO jspin = 1,input%jspins
              WRITE (6,'(a7,i3)') 'spin #',jspin
@@ -128,7 +121,7 @@ CONTAINS
                 ELSE
                    WRITE (6,*) 'using the atomic limit of LDA+U '
                 ENDIF
-                WRITE (6,l_form) ((vs_mmp(k,j,i_u,jspin),k=-l,l),j=-l,l)
+                WRITE (6,l_form) ((pot%mmpMat(k,j,i_u,jspin),k=-l,l),j=-l,l)
              END DO
           END DO
           WRITE (6,*) results%e_ldau
@@ -139,9 +132,8 @@ CONTAINS
           WRITE (*,*) 'no density matrix found ... skipping LDA+U'
           WRITE(*,*) "File:",filename
        ENDIF
-       vs_mmp(:,:,:,:) = CMPLX(0.0,0.0)
+       pot%mmpMat(:,:,:,:) = CMPLX(0.0,0.0)
        results%e_ldau = 0.0
-       atoms%n_u = 0
     ENDIF
 
     RETURN

ldau/umtx.f90

@@ -6,7 +6,7 @@ MODULE m_umtx
   !* Extension to multiple U per atom type by G.M. 2017                *
   !*********************************************************************
 CONTAINS
-  SUBROUTINE umtx(atoms,lmaxb,f0,f2,f4,f6,&
+  SUBROUTINE umtx(atoms,f0,f2,f4,f6,&
                   u)
 
     USE m_constants
@@ -14,16 +14,17 @@ CONTAINS
     USE m_types
     IMPLICIT NONE
 
-    INTEGER, PARAMETER   :: lmaxw=3,lmmaxw1=(2*lmaxw+2)**2
-    TYPE(t_atoms),  INTENT(IN)  :: atoms
-    INTEGER, INTENT (IN) :: lmaxb
-    REAL,    INTENT (IN) :: f0(atoms%n_u),f2(atoms%n_u),f4(atoms%n_u),f6(atoms%n_u)
-    REAL,    INTENT (OUT) :: u(-lmaxb:lmaxb,-lmaxb:lmaxb,-lmaxb:lmaxb,-lmaxb:lmaxb,atoms%n_u)
+    TYPE(t_atoms), INTENT(IN)  :: atoms
+    REAL,          INTENT(IN)  :: f0(atoms%n_u),f2(atoms%n_u),f4(atoms%n_u),f6(atoms%n_u)
+    REAL,          INTENT(OUT) :: u(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,&
+                                    -lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,atoms%n_u)
+
+    INTEGER, PARAMETER         :: lmaxw=3,lmmaxw1=(2*lmaxw+2)**2
 
     INTEGER i,j,k,l,m,mk,nfk,itype,i_u
     INTEGER m1,m2,m3,m4,lm1,lm2,lm3,lm4,kf
     REAL    uk,uq,avu,avj,cgk1,cgk2,tol
-    REAL    fk(lmaxb+1,atoms%n_u)
+    REAL    fk(lmaxU_const+1,atoms%n_u)
     REAL,   ALLOCATABLE :: c(:,:,:)
     !
     tol = 1.0e-14
@@ -54,7 +55,7 @@ CONTAINS
        END DO
     END DO
 
-    CALL sgaunt(lmaxw,lmmaxw1,lmaxb,c)
+    CALL sgaunt(lmaxw,lmmaxw1,lmaxU_const,c)
 
     DO i_u = 1, atoms%n_u                     !!! over U parameters
        l = atoms%lda_u(i_u)%l

ldau/v_mmp.F90

@@ -19,7 +19,7 @@ MODULE m_vmmp
   !     Extension to multiple U per atom type  G.M. 2017
   !     ************************************************************
 CONTAINS
-  SUBROUTINE v_mmp(sym,atoms,jspins,lmaxb,ns_mmp,u,f0,f2, vs_mmp,results)
+  SUBROUTINE v_mmp(sym,atoms,jspins,ns_mmp,u,f0,f2, vs_mmp,results)
 
     USE m_types
     USE m_constants
@@ -29,13 +29,14 @@ CONTAINS
     TYPE(t_atoms),INTENT(IN)        :: atoms
     !
     ! ..  Arguments ..
-    INTEGER, INTENT (IN) :: lmaxb,jspins 
-    REAL,    INTENT (IN) :: u(-lmaxb:lmaxb,-lmaxb:lmaxb, -lmaxb:lmaxb,-lmaxb:lmaxb,atoms%n_u)
-    REAL,    INTENT (IN) :: f0(atoms%n_u),f2(atoms%n_u)
+    INTEGER, INTENT(IN)  :: jspins 
+    REAL,    INTENT(IN)  :: u(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,&
+                              -lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,atoms%n_u)
+    REAL,    INTENT(IN)  :: f0(atoms%n_u),f2(atoms%n_u)
+    COMPLEX, INTENT(OUT) :: vs_mmp(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,atoms%n_u,jspins)
+
+    COMPLEX              :: ns_mmp(-lmaxU_const:lmaxU_const,-lmaxU_const:lmaxU_const,atoms%n_u,jspins)
 
-    COMPLEX           :: ns_mmp(-lmaxb:lmaxb,-lmaxb:lmaxb,atoms%n_u,jspins)
-    COMPLEX,INTENT(OUT)::vs_mmp(-lmaxb:lmaxb,-lmaxb:lmaxb,atoms%n_u,jspins)
-    !
     ! ..  Local Variables ..
     INTEGER ispin,jspin,l ,mp,p,q,itype,m,i_u
     REAL rho_tot,u_htr,j_htr,e_ee,ns_sum,spin_deg,e_dc,e_dcc

tests/CMakeLists.txt

@@ -4,7 +4,7 @@ set(SerialParallelTests CuBulk CuBulkXML Fe_1l Fe_1lXML Fe-Atom CuBand CuBandXML
 Fe_bct Fe_bctXML PTO PTOXML Fe_1l_SOC Fe_1l_SOCXML PTO-SOC PTO-SOCXML Fe_bct_SOC Fe_bct_SOCXML
 GaAsMultiUForceXML SiFilmSlicePlotXML)
 
-set(SerialOnlyTests Fe_bct_LO Fe_bct_LOXML Fe_fcc Fe_fccXML TiO2eels)
+set(SerialOnlyTests Fe_bct_LO Fe_bct_LOXML Fe_fcc Fe_fccXML TiO2eels TiO2eelsXML)
 
 set(Testdirs ${SerialParallelTests} ${SerialOnlyTests})
 set(ParTestdirs ${SerialParallelTests})

tests/tests/TiO2eelsXML/files/fort.37

+    0    1 430.098      4.3808E-03      0.0000E+00      5.3804E-03      0.0000E+00
+    0    1 431.098      1.1145E-03      0.0000E+00      2.2368E-03      0.0000E+00
+    0    1 432.098      5.4398E-07      0.0000E+00      2.4874E-06      0.0000E+00
+    0    1 433.098      5.7937E-15      0.0000E+00      3.4247E-14      0.0000E+00
+    0    1 434.098      2.8353E-27      0.0000E+00      1.4684E-26      0.0000E+00
+    0    1 435.098      8.8574E-14      0.0000E+00      4.1369E-13      0.0000E+00
+    0    1 436.098      2.1025E-05      0.0000E+00      9.0098E-05      0.0000E+00
+    0    1 437.098      4.6458E-02      0.0000E+00      1.4574E-01      0.0000E+00
+    0    1 438.098      1.9176E-01      0.0000E+00      2.2600E-01      0.0000E+00
+    0    1 439.098      7.0534E-02      0.0000E+00      1.1324E-01      0.0000E+00
+    0    1 440.098      7.8333E-02      0.0000E+00      3.2395E-01      0.0000E+00
+    0    1 441.098      7.8703E-02      0.0000E+00      1.5333E-01      0.0000E+00
+    0    1 442.098      2.8849E-02      0.0000E+00      1.2373E-01      0.0000E+00
+    0    1 443.098      8.9866E-04      0.0000E+00      5.0424E-03      0.0000E+00
+    0    1 444.098      4.6961E-03      0.0000E+00      1.9480E-02      0.0000E+00
+    0    1 445.098      1.7447E-03      0.0000E+00      1.9236E-03      0.0000E+00
+    0    1 446.098      7.4332E-03      0.0000E+00      7.6095E-03      0.0000E+00
+ 
+ 
+    1    1 430.098      6.3784E-02      0.0000E+00      3.5411E-02      0.0000E+00
+    1    1 431.098      1.4520E-02      0.0000E+00      1.4683E-02      0.0000E+00
+    1    1 432.098      3.7277E-06      0.0000E+00      1.6287E-05      0.0000E+00
+    1    1 433.098      2.0842E-14      0.0000E+00      2.2391E-13      0.0000E+00
+    1    1 434.098      1.7334E-26      0.0000E+00      9.6432E-26      0.0000E+00
+    1    1 435.098      6.6101E-13      0.0000E+00      2.7179E-12      0.0000E+00
+    1    1 436.098      1.7587E-04      0.0000E+00      5.9129E-04      0.0000E+00
+    1    1 437.098      5.1795E-01      0.0000E+00      9.5547E-01      0.0000E+00
+    1    1 438.098      3.0585E+00      0.0000E+00      1.4812E+00      0.0000E+00
+    1    1 439.098      1.0505E+00      0.0000E+00      7.4181E-01      0.0000E+00
+    1    1 440.098      6.7892E-01      0.0000E+00      2.1229E+00      0.0000E+00
+    1    1 441.098      1.1043E+00      0.0000E+00      1.0042E+00      0.0000E+00
+    1    1 442.098      2.3789E-01      0.0000E+00      8.1016E-01      0.0000E+00
+    1    1 443.098      4.3919E-03      0.0000E+00      3.2996E-02      0.0000E+00
+    1    1 444.098      4.0029E-02      0.0000E+00      1.2750E-01      0.0000E+00
+    1    1 445.098      2.6406E-02      0.0000E+00      1.2632E-02      0.0000E+00
+    1    1 446.098      1.1761E-01      0.0000E+00      4.9857E-02      0.0000E+00
+ 
+ 

tests/tests/TiO2eelsXML/files/inp-2.xml

+<?xml version="1.0" encoding="UTF-8" standalone="no"?>
+<fleurInput fleurInputVersion="0.28">
+   <comment>
+      TiO2 Rutile                                                                     
+   </comment>
+   <calculationSetup>
+      <cutoffs Kmax="3.80000000" Gmax="9.0" GmaxXC="8.00000000" numbands="40"/>
+      <scfLoop itmax="1" minDistance=".00000000" maxIterBroyd="99" imix="Anderson" alpha=".05000000" spinf="2.00000000"/>
+      <coreElectrons ctail="T" frcor="F" kcrel="0" coretail_lmax="99"/>
+      <magnetism jspins="1" l_noco="F" l_J="F" swsp="F" lflip="F"/>
+      <soc theta=".00000000" phi=".00000000" l_soc="F" spav="F" off="F"/>
+      <nocoParams l_ss="F" l_mperp="F" l_constr="F" l_disp="F" sso_opt="FFF" mix_b=".00000000" thetaJ=".00000000" nsh="0">
+         <qss>.0000000000 .0000000000 .0000000000</qss>
+      </nocoParams>
+      <expertModes gw="0" pot8="F" isec1="99" secvar="F"/>
+      <geometryOptimization l_f="F" xa="2.00000000" thetad="330.00000000" epsdisp=".00001000" epsforce=".00001000"/>
+      <bzIntegration valenceElectrons="44.00000000" mode="hist" fermiSmearingEnergy=".00100000">
+         <kPointList posScale="1.00000000" weightScale="1.00000000" count="2">
+            <kPoint weight="    0.500000">    0.250000     0.250000     0.375000</kPoint>
+            <kPoint weight="    0.500000">    0.250000     0.250000     0.125000</kPoint>
+         </kPointList>
+      </bzIntegration>
+      <energyParameterLimits ellow="-1.80000000" elup="1.00000000"/>
+   </calculationSetup>
+   <cell>
+      <symmetryOperations>
+         <symOp>
+            <row-1>1 0 0 .0000000000</row-1>
+            <row-2>0 1 0 .0000000000</row-2>
+            <row-3>0 0 1 .0000000000</row-3>
+         </symOp>
+         <symOp>
+            <row-1>-1 0 0 .0000000000</row-1>
+            <row-2>0 -1 0 .0000000000</row-2>
+            <row-3>0 0 -1 .0000000000</row-3>
+         </symOp>
+         <symOp>
+            <row-1>-1 0 0 .5000000000</row-1>
+            <row-2>0 1 0 .5000000000</row-2>
+            <row-3>0 0 -1 .5000000000</row-3>
+         </symOp>
+         <symOp>
+            <row-1>1 0 0 .5000000000</row-1>
+            <row-2>0 -1 0 .5000000000</row-2>
+            <row-3>0 0 1 .5000000000</row-3>
+         </symOp>
+         <symOp>
+            <row-1>0 -1 0 .5000000000</row-1>
+            <row-2>1 0 0 .5000000000</row-2>
+            <row-3>0 0 1 .5000000000</row-3>
+         </symOp>
+         <symOp>
+            <row-1>-1 0 0 .0000000000</row-1>
+            <row-2>0 -1 0 .0000000000</row-2>
+            <row-3>0 0 1 .0000000000</row-3>