1#define TCGMSG
2
3
4*     ***************************
5*     *                         *
6*     *    semicore_xc_F        *
7*     *                         *
8*     ***************************
9
10      subroutine semicore_xc_F(ispin,xcp,fion)
11      implicit none
12#include "errquit.fh"
13      integer ispin
14      real*8  xcp(*)
15
16      real*8 fion(3,*)
17
18#include "bafdecls.fh"
19
20*     **** semicore common block ****
21c     real*8  ncore(nfft3d,nkatmx),rcore(nkatmx)
22c     logocal semicore(0:nkatmx)
23      integer ncore(2),rcore(2)
24      integer semicore(2)
25      common / ccore / ncore,rcore,semicore
26
27*     **** local variables ****
28      logical value
29      integer npack0,nfft3d,n2ft3d
30      integer ii,ia,nx,ny,nz
31      real*8  sumx,sumy,sumz
32      real*8  scal1,scal2
33      integer exi(2),vxcG(2)
34      integer Gx(2),Gy(2),Gz(2)
35      integer xtmp(2),dng(2)
36
37*     **** external functions ****
38      integer  ion_nion,ion_katm
39      real*8   lattice_omega
40      external ion_nion,ion_katm
41      external lattice_omega
42
43
44      call D3dB_nx(1,nx)
45      call D3dB_ny(1,ny)
46      call D3dB_nz(1,nz)
47      scal1 = 1.0d0/dble(nx*ny*nz)
48      scal2 = 1.0d0/lattice_omega()
49
50      call D3dB_nfft3d(1,nfft3d)
51      call D3dB_n2ft3d(1,n2ft3d)
52      call Pack_npack(0,npack0)
53
54      value = BA_push_get(mt_dcpl,nfft3d,'exi', exi(2), exi(1))
55      value = value.and.
56     >        BA_push_get(mt_dcpl,nfft3d,'vxcG',vxcG(2),vxcG(1))
57
58      value = value.and.
59     >        BA_push_get(mt_dbl, nfft3d,'Gx',Gx(2),Gx(1))
60      value = value.and.
61     >        BA_push_get(mt_dbl, nfft3d,'Gy',Gy(2),Gy(1))
62      value = value.and.
63     >        BA_push_get(mt_dbl, nfft3d,'Gz',Gz(2),Gz(1))
64
65      value = value.and.
66     >        BA_push_get(mt_dbl, npack0,'xtmp',xtmp(2),xtmp(1))
67
68      value = value.and.
69     >        BA_push_get(mt_dcpl, nfft3d,'dng',dng(2),dng(1))
70      if (.not. value) call errquit('out of stack memory',0, MA_ERR)
71
72
73      call dcopy(n2ft3d,0.0d0,0,dcpl_mb(vxcG(1)),1)
74      call D3dB_rr_Sum(1,xcp(1),
75     >                   xcp(1+(ispin-1)*n2ft3d),
76     >                   dcpl_mb(vxcG(1)))
77      call D3dB_rc_fft3f(1,dcpl_mb(vxcG(1)))
78      call Pack_c_pack(0,dcpl_mb(vxcG(1)))
79
80*     **** define Gx, Gy, and Gz in packed space ****
81      call D3dB_t_Copy(1,dbl_mb(G_indx(1)),dbl_mb(Gx(1)))
82      call D3dB_t_Copy(1,dbl_mb(G_indx(2)),dbl_mb(Gy(1)))
83      call D3dB_t_Copy(1,dbl_mb(G_indx(3)),dbl_mb(Gz(1)))
84      call Pack_t_pack(0,dbl_mb(Gx(1)))
85      call Pack_t_pack(0,dbl_mb(Gy(1)))
86      call Pack_t_pack(0,dbl_mb(Gz(1)))
87
88
89      do ii=1,ion_nion()
90         ia = ion_katm(ii)
91
92         if (log_mb(semicore(1)+ia)) then
93
94*          **** structure factor and local pseudopotential ****
95           call strfac(ii,dcpl_mb(exi(1)))
96           call Pack_c_pack(0,dcpl_mb(exi(1)))
97
98*          **** put sqrt(core-density) at atom position ****
99           call Pack_tc_iMul(0,
100     >               dbl_mb(ncore(1)+(ia-1)*5*npack0),
101     >              dcpl_mb(exi(1)),
102     >              dcpl_mb(dng(1)))
103
104*          **** put dng in real space and square it ****
105           call Pack_c_unpack(0,dcpl_mb(dng(1)))
106           call D3dB_cr_fft3b(1,dcpl_mb(dng(1)))
107           call D3dB_rr_Sqr(1,dcpl_mb(dng(1)),dcpl_mb(dng(1)))
108
109*          **** put dng back in complex space ****
110           call D3dB_r_Zero_Ends(1,dcpl_mb(dng(1)))
111           call D3dB_rc_fft3f(1,dcpl_mb(dng(1)))
112           call Pack_c_pack(0,dcpl_mb(dng(1)))
113
114           do k=1,npack0
115              dbl_mb(xtmp(1)+k-1) =
116     >         dimag(dcpl_mb(VxcG(1)+k-1))* dble(dcpl_mb(dng(1)+k-1))
117     >        - dble(dcpl_mb(VxcG(1)+k-1))*dimag(dcpl_mb(dng(1)+k-1))
118           end do
119
120           call Pack_tt_dot(0,dbl_mb(Gx(1)),dbl_mb(xtmp(1)),sumx)
121           call Pack_tt_dot(0,dbl_mb(Gy(1)),dbl_mb(xtmp(1)),sumy)
122           call Pack_tt_dot(0,dbl_mb(Gz(1)),dbl_mb(xtmp(1)),sumz)
123
124           fion(1,ii) = fion(1,ii) + sumx*scal2*scal2*dsqrt(scal1)
125           fion(2,ii) = fion(2,ii) + sumy*scal2*scal2*dsqrt(scal1)
126           fion(3,ii) = fion(3,ii) + sumz*scal2*scal2*dsqrt(scal1)
127         end if
128
129      end do
130
131      value = BA_pop_stack(dng(2))
132      value = BA_pop_stack(xtmp(2))
133      value = BA_pop_stack(Gz(2))
134      value = BA_pop_stack(Gy(2))
135      value = BA_pop_stack(Gx(2))
136      value = BA_pop_stack(vxcG(2))
137      value = BA_pop_stack(exi(2))
138
139
140      return
141      end
142
143
144c $Id$
145