create add_l0_fun

hybrid/hybrid.F90

@@ -85,13 +85,12 @@ CONTAINS
          if(allocated(hybrid%ne_eig)) deallocate(hybrid%ne_eig)
          allocate(hybrid%ne_eig(kpts%nkpt), source=0)
 
-         write (*,*) "allocated(hybrid%nbands): ",allocated(hybrid%nbands)
-         write (*,*) "shape(hybrid%nbands)", shape(hybrid%nbands)
-         if(allocated(hybrid%nbands)) deallocate(hybrid%nbands, stat=err, errmsg=msg)
-
-         if(err /= 0) THEN
-            write (*,*) "errorcode", err
-            write (*,*) "errormessage", msg
+         if(allocated(hybrid%nbands)) then
+            deallocate(hybrid%nbands, stat=err, errmsg=msg)
+            if(err /= 0) THEN
+               write (*,*) "errorcode", err
+               write (*,*) "errormessage", msg
+            endif
          endif
 
          allocate(hybrid%nbands(kpts%nkpt), source=0)
@@ -102,7 +101,7 @@ CONTAINS
          if(allocated(hybrid%nbasm)) deallocate(hybrid%nbasm)
          allocate(hybrid%nbasm(kpts%nkptf), source=0)
 
-         if(allocated(hybrid%nbasm)) deallocate(hybrid%nbasm)
+         if(allocated(hybrid%div_vv)) deallocate(hybrid%div_vv)
          allocate(hybrid%div_vv(DIMENSION%neigd, kpts%nkpt, input%jspins), source=0.0)
          init_vex = .FALSE.
       END IF

hybrid/mixedbasis.F90

@@ -279,36 +279,7 @@ CONTAINS
             nn = mpbasis%num_rad_bas_fun(l, itype)
 
             ! Add constant function to l=0 basis and then do a Gram-Schmidt orthonormalization
-            IF (l == 0) THEN
-
-               ! Check if radbasfn_mt must be reallocated
-               IF (nn + 1 > SIZE(mpbasis%radbasfn_mt, 2)) THEN
-                  allocate(basmhlp(atoms%jmtd, nn + 1, 0:maxval(hybrid%lcutm1), atoms%ntype))
-                  basmhlp(:,1:nn, :,:) = mpbasis%radbasfn_mt
-                  deallocate(mpbasis%radbasfn_mt)
-                  allocate(mpbasis%radbasfn_mt(atoms%jmtd, nn + 1, 0:maxval(hybrid%lcutm1), atoms%ntype))
-                  mpbasis%radbasfn_mt(:,1:nn, :,:) = basmhlp(:,1:nn, :,:)
-                  deallocate(basmhlp)
-               END IF
-
-               mpbasis%radbasfn_mt(:n_grid_pt, nn + 1, 0, itype) = atoms%rmsh(:n_grid_pt, itype)/SQRT(atoms%rmsh(n_grid_pt, itype)**3/3)
-               DO i = nn, 1, -1
-                  DO j = i + 1, nn + 1
-                     mpbasis%radbasfn_mt(:n_grid_pt, i, 0, itype) = mpbasis%radbasfn_mt(:n_grid_pt, i, 0, itype) &
-                                                      - intgrf(mpbasis%radbasfn_mt(:n_grid_pt, i, 0, itype)*mpbasis%radbasfn_mt(:n_grid_pt, j, 0, itype), &
-                                                               atoms%jri, atoms%jmtd, atoms%rmsh, atoms%dx, atoms%ntype, itype, gridf) &
-                                                      *mpbasis%radbasfn_mt(:n_grid_pt, j, 0, itype)
-
-                  END DO
-                  mpbasis%radbasfn_mt(:n_grid_pt, i, 0, itype) = mpbasis%radbasfn_mt(:n_grid_pt, i, 0, itype) &
-                                                   /SQRT(intgrf(mpbasis%radbasfn_mt(:n_grid_pt, i, 0, itype)**2, &
-                                                                atoms%jri, atoms%jmtd, atoms%rmsh, atoms%dx, atoms%ntype, itype, gridf))
-               END DO
-               nn = nn + 1
-               deallocate(olap)
-               allocate(olap(nn, nn))
-               mpbasis%num_rad_bas_fun(l, itype) = nn
-            END IF
+            call add_l0_fun(atoms, hybrid, n_grid_pt, l, itype, gridf, nn, mpbasis)
 
             ! Check orthonormality of product basis
             rdum = 0
@@ -479,6 +450,7 @@ CONTAINS
       END DO
       hybrid%maxbasm1 = hybrid%nbasp + maxval(mpbasis%ngptm)
       DO nk = 1, kpts%nkptf
+         if(.not. allocated(hybrid%nbasm)) call judft_error("hybrid%nbasm not alloced")
          hybrid%nbasm(nk) = hybrid%nbasp + mpbasis%ngptm(nk)
       END DO
 
@@ -683,4 +655,48 @@ CONTAINS
             = MATMUL(mpbasis%radbasfn_mt(i, 1:n_radbasfn, l, itype), eigv(:,1:nn))/SQRT(eig(:nn))
       END DO
    end subroutine trafo_to_orthonorm_bas
+
+   subroutine add_l0_fun(atoms, hybrid, n_grid_pt, l, itype, gridf, nn, mpbasis)
+      use m_types
+      USE m_intgrf, ONLY: intgrf
+      implicit none
+      type(t_atoms), intent(in)       :: atoms
+      type(t_hybrid), intent(in)      :: hybrid
+      integer, intent(in)             :: n_grid_pt, l, itype
+      real, intent(in)                :: gridf(:,:)
+      integer, intent(inout)          :: nn
+      type(t_mpbasis), intent(inout)  :: mpbasis
+
+      integer                         :: i, j
+      REAL, ALLOCATABLE               :: basmhlp(:,:,:,:)
+
+      IF (l == 0) THEN
+
+         ! Check if radbasfn_mt must be reallocated
+         IF (nn + 1 > SIZE(mpbasis%radbasfn_mt, 2)) THEN
+            allocate(basmhlp(atoms%jmtd, nn + 1, 0:maxval(hybrid%lcutm1), atoms%ntype))
+            basmhlp(:,1:nn, :,:) = mpbasis%radbasfn_mt
+            deallocate(mpbasis%radbasfn_mt)
+            allocate(mpbasis%radbasfn_mt(atoms%jmtd, nn + 1, 0:maxval(hybrid%lcutm1), atoms%ntype))
+            mpbasis%radbasfn_mt(:,1:nn, :,:) = basmhlp(:,1:nn, :,:)
+            deallocate(basmhlp)
+         END IF
+
+         mpbasis%radbasfn_mt(:n_grid_pt, nn + 1, 0, itype) = atoms%rmsh(:n_grid_pt, itype)/SQRT(atoms%rmsh(n_grid_pt, itype)**3/3)
+         DO i = nn, 1, -1
+            DO j = i + 1, nn + 1
+               mpbasis%radbasfn_mt(:n_grid_pt, i, 0, itype) = mpbasis%radbasfn_mt(:n_grid_pt, i, 0, itype) &
+                                                - intgrf(mpbasis%radbasfn_mt(:n_grid_pt, i, 0, itype)*mpbasis%radbasfn_mt(:n_grid_pt, j, 0, itype), &
+                                                         atoms%jri, atoms%jmtd, atoms%rmsh, atoms%dx, atoms%ntype, itype, gridf) &
+                                                *mpbasis%radbasfn_mt(:n_grid_pt, j, 0, itype)
+
+            END DO
+            mpbasis%radbasfn_mt(:n_grid_pt, i, 0, itype) = mpbasis%radbasfn_mt(:n_grid_pt, i, 0, itype) &
+                                             /SQRT(intgrf(mpbasis%radbasfn_mt(:n_grid_pt, i, 0, itype)**2, &
+                                                          atoms%jri, atoms%jmtd, atoms%rmsh, atoms%dx, atoms%ntype, itype, gridf))
+         END DO
+         nn = nn + 1
+         mpbasis%num_rad_bas_fun(l, itype) = nn
+      END IF
+   end subroutine add_l0_fun
 END MODULE m_mixedbasis