1* 2* $Id$ 3* 4 subroutine lcao_init_dn(ispin,ne,n2ft3d,dn,phi) 5 implicit none 6#include "errquit.fh" 7 integer ispin,ne(2) 8 integer n2ft3d 9 real*8 dn(n2ft3d,2) 10 complex*16 phi(*) 11 12#include "bafdecls.fh" 13#include "global.fh" 14 15* **** local variables **** 16 logical value 17 integer i,k,nbasis,ms 18 integer nx,ny,nz 19 real*8 sum,scal,scal1,scal2,dv 20 real*8 dnscal(2) 21 integer tmp(2) 22 23 24* ***** external functions **** 25 integer aorbs_nbasis 26 real*8 aorbs_weight,lattice_omega,util_random 27 external aorbs_nbasis 28 external aorbs_weight,lattice_omega,util_random 29 30 dnscal(1) = dble(ne(1))/dble(ne(1)+ne(2)) 31 dnscal(2) = dble(ne(2))/dble(ne(1)+ne(2)) 32 33 call D3dB_nx(1,nx) 34 call D3dB_ny(1,ny) 35 call D3dB_nz(1,nz) 36 scal1 = 1.0d0/dble(nx*ny*nz) 37 scal2 = 1.0d0/lattice_omega() 38 dv = scal1*lattice_omega() 39 40 value = BA_push_get(mt_dbl,n2ft3d,'tmp',tmp(2),tmp(1)) 41 if (.not. value) 42 > call errquit('lcao_init_dn:out of stack memory',0, MA_ERR) 43 44 45 call dcopy(ispin*n2ft3d,0.0d0,0,dn,1) 46 nbasis = aorbs_nbasis() 47 do i=1,nbasis 48 49* **** get phi1 **** 50 call aorbs_aorb(i,phi) 51 call Pack_c_Copy(1,phi,dbl_mb(tmp(1))) 52 call Pack_c_unpack(1,dbl_mb(tmp(1))) 53 call D3dB_cr_pfft3b(1,1,dbl_mb(tmp(1))) 54 55 do ms=1,ispin 56 do k=1,n2ft3d 57 scal = aorbs_weight(i)*scal2*dnscal(ms) 58 dn(i,ms) = dn(i,ms) + scal*(dbl_mb(tmp(1)+k-1)**2) 59 end do 60 end do 61 62 end do 63 64 65* **** randomize dn(*,2) if ispin=2 and ne(1)==ne(2) **** 66 if ((ispin.eq.2).and.(ne(1).eq.ne(2))) then 67 scal = util_random(9) 68 do k=1,n2ft3d 69 scal = (0.5d0-util_random(0))/dsqrt(dble(n2ft3d)) 70 dn(i,2) = dn(i,2) + scal 71 end do 72 end if 73 74 75* **** normalize densities *** 76 do ms=1,ispin 77 call D3dB_r_Zero_Ends(1,dn(1,ms)) 78 call D3dB_r_dsum(1,dn(1,ms),sum) 79 sum = sum*dv 80 sum = dble(ne(ms))/sum 81c call D3dB_r_SMul(1,sum,dn(1,ms),dn(1,ms)) 82 call D3dB_r_SMul1(1,sum,dn(1,ms)) 83 end do 84 85 value = BA_pop_stack(tmp(2)) 86 if (.not. value) 87 > call errquit('lcao_init_dn:popping stack memory',0, MA_ERR) 88 return 89 end 90 91 92