TABLE OF CONTENTS


ABINIT/dfpt_mkvxcstr [ Functions ]

[ Top ] [ Functions ]

NAME

 dfpt_mkvxcstr

FUNCTION

 Compute the first-order change of exchange-correlation potential
 due to strain: assemble the first-order density change with the
 frozen-core density change, then use the exchange-correlation kernel.

INPUTS

  cplex= if 1, real space 1-order functions on FFT grid are REAL, if 2, COMPLEX
  idir=direction of the current perturbation
  ipert=type of the perturbation
  kxc(nfft,nkxc)=exchange and correlation kernel (see rhotoxc.f)
  mpi_enreg=information about MPI parallelization
  natom=number of atoms in cell.
  nfft=(effective) number of FFT grid points (for this processor)
  ngfft(18)=contain all needed information about 3D FFT,
     see ~abinit/doc/variables/vargs.htm#ngfft
  nhat(nfft,nspden*nhatdim)= -PAW only- compensation density
  nhat1(cplex*nfft,2nspden*usepaw)= -PAW only- 1st-order compensation density
  nkxc=second dimension of the kxc array
  nspden=number of spin-density components
  n3xccc=dimension of xccc3d1 ; 0 if no XC core correction is used, otherwise, nfft
  option=if 0, work only with strain-derivative frozen-wavefunction
    charge and the XC core-correction,
   if 1, treat both density change and XC core correction
   if 2, like 0 but multiply gradient strain derivative term by 2.0 for GGA.
  qphon(3)=reduced coordinates for the phonon wavelength (needed if cplex==2).
  rhor(nfft,nspden)=array for GS electron density in electrons/bohr**3.
  rhor1(cplex*nfft,nspden)=array for electron density in electrons/bohr**3.
  rprimd(3,3)=dimensional primitive translations in real space (bohr)
  usepaw= 0 for non paw calculation; =1 for paw calculation
  usexcnhat= -PAW only- flag controling use of compensation density in Vxc
  xccc3d1(cplex*n3xccc)=3D change in core charge density, see n3xccc

OUTPUT

  vxc1(cplex*nfft,nspden)=change in exchange-correlation potential (including
   core-correction, if applicable)

PARENTS

      dfpt_eltfrxc,dfpt_nselt,dfpt_nstpaw,dfpt_rhotov

CHILDREN

      dfpt_mkvxcstrgga,matr3inv,timab

SOURCE

 97 subroutine dfpt_mkvxcstr(cplex,idir,ipert,kxc,mpi_enreg,natom,nfft,ngfft,nhat,nhat1,&
 98 & nkxc,nspden,n3xccc,option,paral_kgb,qphon,rhor,rhor1,rprimd,usepaw,usexcnhat,vxc1,xccc3d1)
 99 
100 
101 !This section has been created automatically by the script Abilint (TD).
102 !Do not modify the following lines by hand.
103 #undef ABI_FUNC
104 #define ABI_FUNC 'dfpt_mkvxcstr'
105 !End of the abilint section
106 
107  implicit none
108 
109 !Arguments ------------------------------------
110 !scalars
111  integer,intent(in) :: cplex,idir,ipert,n3xccc,natom,nfft,nkxc,nspden,option
112  integer,intent(in) :: paral_kgb,usepaw,usexcnhat
113  type(MPI_type),intent(in) :: mpi_enreg
114 !arrays
115  integer,intent(in) :: ngfft(18)
116  real(dp),target,intent(in) :: nhat(nfft,nspden)
117  real(dp),target,intent(in) :: nhat1(cplex*nfft,nspden)
118  real(dp),intent(in) :: kxc(nfft,nkxc),qphon(3)
119  real(dp),target,intent(in) :: rhor(nfft,nspden),rhor1(cplex*nfft,nspden)
120  real(dp),intent(in) :: rprimd(3,3)
121  real(dp),intent(in) :: xccc3d1(cplex*n3xccc)
122  real(dp),intent(out) :: vxc1(cplex*nfft,nspden)
123 
124 !Local variables-------------------------------
125 !scalars
126  integer :: ii,ir,istr
127  real(dp) :: rho1_dn,rho1_up,spin_scale,str_scale
128  character(len=500) :: message
129 !arrays
130  real(dp) :: gprimd(3,3),tsec(2)
131  real(dp),allocatable :: rhor1tmp(:,:),rhowk1(:,:)
132  real(dp),pointer :: rhor_(:,:),rhor1_(:,:)
133 
134 ! *************************************************************************
135 
136  call timab(181,1,tsec)
137 
138  if(nspden/=1 .and. nspden/=2) then
139    message = ' dfpt_mkvxc, Only for nspden==1 and 2.'
140    MSG_BUG(message)
141  end if
142 
143  if (usepaw==1.and.usexcnhat==0) then
144    ABI_ALLOCATE(rhor_,(nfft,nspden))
145    rhor_(:,:)=rhor(:,:)-nhat(:,:)
146  else
147    rhor_ => rhor
148  end if
149 
150  if (usepaw==1.and.usexcnhat==0.and.option==1) then
151    ABI_ALLOCATE(rhor1_,(nfft,nspden))
152    rhor1_(:,:)=rhor1(:,:)-nhat1(:,:)
153  else
154    rhor1_ => rhor1
155  end if
156 
157 
158 !Inhomogeneous term for diagonal strain
159  ABI_ALLOCATE(rhowk1,(nfft,nspden))
160  if(option==0 .or. option==2) then
161    if(ipert==natom+3) then
162      rhowk1(:,:)=-rhor_(:,:)
163    else
164      rhowk1(:,:)=zero
165    end if
166  else if(option==1) then
167    if(ipert==natom+3) then
168      rhowk1(:,:)=rhor1_(:,:)-rhor_(:,:)
169    else
170      rhowk1(:,:)=rhor1_(:,:)
171    end if
172  end if
173 
174 !Treat first LDA
175  if(nkxc==1.or.nkxc==3)then
176 
177 !  Case without non-linear core correction
178    if(n3xccc==0)then
179 
180 !    Non-spin-polarized
181      if(nspden==1)then
182        do ir=1,nfft
183          vxc1(ir,1)=kxc(ir,1)*rhowk1(ir,1)
184        end do
185 
186 !      Spin-polarized
187      else
188        do ir=1,nfft
189          rho1_dn=rhowk1(ir,1)-rhowk1(ir,2)
190          vxc1(ir,1)=kxc(ir,1)*rhowk1(ir,2)+kxc(ir,2)*rho1_dn
191          vxc1(ir,2)=kxc(ir,2)*rhowk1(ir,2)+kxc(ir,3)*rho1_dn
192        end do
193      end if ! nspden==1
194 
195 !    Treat case with non-linear core correction
196    else
197      if(nspden==1)then
198        do ir=1,nfft
199          vxc1(ir,1)=kxc(ir,1)*(rhowk1(ir,1)+xccc3d1(ir))
200        end do
201      else
202        do ir=1,nfft
203          rho1_dn=rhowk1(ir,1)-rhowk1(ir,2) + xccc3d1(ir)*half
204          rho1_up=rhowk1(ir,2)              + xccc3d1(ir)*half
205          vxc1(ir,1)=kxc(ir,1)*rho1_up+kxc(ir,2)*rho1_dn
206          vxc1(ir,2)=kxc(ir,2)*rho1_up+kxc(ir,3)*rho1_dn
207        end do
208      end if ! nspden==1
209 
210    end if ! n3xccc==0
211 
212 !  Treat GGA
213  else if (nkxc==7.or.nkxc==19) then
214 
215 !  Generates gprimd and its strain derivative
216 !  Note that unlike the implicitly symmetric metric tensor strain
217 !  derivatives, we must explicltly symmetrize the strain derivative
218 !  here.
219    call matr3inv(rprimd,gprimd)
220    istr=idir + 3*(ipert-natom-3)
221    if(istr<1 .or. istr>6)then
222      write(message, '(a,i10,a,a,a)' )&
223 &     'Input dir gives istr=',istr,' not allowed.',ch10,&
224 &     'Possible values are 1,2,3,4,5,6 only.'
225      MSG_BUG(message)
226    end if
227 
228 !  Rescalling needed for use in dfpt_eltfrxc for elastic tensor (not internal strain tensor).
229    str_scale=one;if(option==2) str_scale=two
230 
231 !  FTransfer the data to spin-polarized storage
232    ABI_ALLOCATE(rhor1tmp,(cplex*nfft,nspden))
233    if(nspden==1)then
234      do ir=1,cplex*nfft
235        rhor1tmp(ir,1)=rhowk1(ir,1)
236      end do
237    else
238      do ir=1,cplex*nfft
239        rho1_dn=rhowk1(ir,1)-rhowk1(ir,2)
240        rhor1tmp(ir,1)=rhowk1(ir,2)
241        rhor1tmp(ir,2)=rho1_dn
242      end do
243    end if ! nspden==1
244    if(n3xccc/=0)then
245      spin_scale=one;if (nspden==2) spin_scale=half
246      do ii=1,nspden
247        do ir=1,cplex*nfft
248          rhor1tmp(ir,ii)=rhor1tmp(ir,ii)+xccc3d1(ir)*spin_scale
249        end do
250      end do
251    end if
252 
253    call dfpt_mkvxcstrgga(cplex,gprimd,istr,kxc,mpi_enreg,nfft,ngfft,nkxc,&
254 &   nspden,paral_kgb,qphon,rhor1tmp,str_scale,vxc1)
255    ABI_DEALLOCATE(rhor1tmp)
256 
257  else
258    MSG_BUG('Invalid nkxc!')
259 
260  end if ! LDA or GGA
261 
262  if (usepaw==1.and.usexcnhat==0) then
263    ABI_DEALLOCATE(rhor_)
264  end if
265  if (usepaw==1.and.usexcnhat==0.and.option==1) then
266    ABI_DEALLOCATE(rhor1_)
267  end if
268 
269  ABI_DEALLOCATE(rhowk1)
270 
271  call timab(181,2,tsec)
272 
273 end subroutine dfpt_mkvxcstr

ABINIT/dfpt_mkvxcstrgga [ Functions ]

[ Top ] [ Functions ]

NAME

 dfpt_mkvxcstrgga

FUNCTION

 Compute the first-order change of exchange-correlation potential
 for the strain perturbation in case of GGA functionals
 Use the exchange-correlation kernel.

INPUTS

  cplex= if 1, real space 1-order functions on FFT grid are REAL,
    if 2, COMPLEX
  gprimd(3,3)=dimensional primitive translations in reciprocal space (bohr^-1)
  istr=index of the strain perturbation (1..6)
  kxc(nfft,nkxc)=exchange and correlation kernel (see rhotoxc.f)
  mpi_enreg=information about MPI parallelization
  nfft=(effective) number of FFT grid points (for this processor)
  ngfft(18)=contain all needed information about 3D FFT,
    see ~abinit/doc/variables/vargs.htm#ngfft
  nkxc=second dimension of the kxc array
  nspden=number of spin-density components
  qphon(3)=reduced coordinates for the phonon wavelength (needed if cplex==2).
  rhor1tmp(cplex*nfft,2)=array for first-order electron spin-density
   in electrons/bohr**3 (second index corresponds to spin-up and spin-down)
  str_scale=scaling factor for gradient operator strain-derivative (1. or 2.)

OUTPUT

  vxc1(cplex*nfft,nspden)=change in exchange-correlation potential

NOTES

  Closely related to dfpt_mkvxcgga.
  Content of Kxc array:
   ===== if GGA
    if nspden==1:
       kxc(:,1)= d2Exc/drho2
       kxc(:,2)= 1/|grad(rho)| dExc/d|grad(rho)|
       kxc(:,3)= 1/|grad(rho)| d2Exc/d|grad(rho)| drho
       kxc(:,4)= 1/|grad(rho)| * d/d|grad(rho)| ( 1/|grad(rho)| dExc/d|grad(rho)| )
       kxc(:,5)= gradx(rho)
       kxc(:,6)= grady(rho)
       kxc(:,7)= gradz(rho)
    if nspden>=2:
       kxc(:,1)= d2Exc/drho_up drho_up
       kxc(:,2)= d2Exc/drho_up drho_dn
       kxc(:,3)= d2Exc/drho_dn drho_dn
       kxc(:,4)= 1/|grad(rho_up)| dEx/d|grad(rho_up)|
       kxc(:,5)= 1/|grad(rho_dn)| dEx/d|grad(rho_dn)|
       kxc(:,6)= 1/|grad(rho_up)| d2Ex/d|grad(rho_up)| drho_up
       kxc(:,7)= 1/|grad(rho_dn)| d2Ex/d|grad(rho_dn)| drho_dn
       kxc(:,8)= 1/|grad(rho_up)| * d/d|grad(rho_up)| ( 1/|grad(rho_up)| dEx/d|grad(rho_up)| )
       kxc(:,9)= 1/|grad(rho_dn)| * d/d|grad(rho_dn)| ( 1/|grad(rho_dn)| dEx/d|grad(rho_dn)| )
       kxc(:,10)=1/|grad(rho)| dEc/d|grad(rho)|
       kxc(:,11)=1/|grad(rho)| d2Ec/d|grad(rho)| drho_up
       kxc(:,12)=1/|grad(rho)| d2Ec/d|grad(rho)| drho_dn
       kxc(:,13)=1/|grad(rho)| * d/d|grad(rho)| ( 1/|grad(rho)| dEc/d|grad(rho)| )
       kxc(:,14)=gradx(rho_up)
       kxc(:,15)=gradx(rho_dn)
       kxc(:,16)=grady(rho_up)
       kxc(:,17)=grady(rho_dn)
       kxc(:,18)=gradz(rho_up)
       kxc(:,19)=gradz(rho_dn)

PARENTS

      dfpt_mkvxcstr

CHILDREN

      xcden,xcpot

SOURCE

346 subroutine dfpt_mkvxcstrgga(cplex,gprimd,istr,kxc,mpi_enreg,nfft,ngfft,&
347 & nkxc,nspden,paral_kgb,qphon,rhor1tmp,str_scale,vxc1)
348 
349 
350 !This section has been created automatically by the script Abilint (TD).
351 !Do not modify the following lines by hand.
352 #undef ABI_FUNC
353 #define ABI_FUNC 'dfpt_mkvxcstrgga'
354 !End of the abilint section
355 
356  implicit none
357 
358 !Arguments ------------------------------------
359 !scalars
360  integer,intent(in) :: cplex,istr,nfft,nkxc,nspden,paral_kgb
361  real(dp) :: str_scale
362  type(MPI_type),intent(in) :: mpi_enreg
363 !arrays
364  integer,intent(in) :: ngfft(18)
365  real(dp),intent(in) :: gprimd(3,3),kxc(nfft,nkxc)
366  real(dp),intent(in) :: qphon(3),rhor1tmp(cplex*nfft,2)
367  real(dp),intent(out) :: vxc1(cplex*nfft,nspden)
368 !Local variables-------------------------------
369 !scalars
370  integer :: ii,ir,ishift,ispden,mgga,ngrad,nspgrad
371  real(dp) :: coeff_grho,coeff_grho_corr,coeff_grho_dn,coeff_grho_up
372  real(dp) :: gradrho_gradrho1,gradrho_gradrho1_dn,gradrho_gradrho1_up
373  character(len=500) :: msg
374 !arrays
375  real(dp) :: r0(3),r0_dn(3),r0_up(3),r1(3),r1_dn(3),r1_up(3)
376  real(dp),allocatable :: dnexcdn(:,:),rho1now(:,:,:),rhodgnow(:,:,:)
377 
378 ! *************************************************************************
379 
380  DBG_ENTER("COLL")
381 
382  if (nkxc/=12*min(nspden,2)-5) then
383    msg='Wrong nkxc value for GGA!'
384    MSG_BUG(msg)
385  end if
386  if (nspden>2) then
387    msg='Not compatible with non-collinear magnetism!'
388    MSG_ERROR(msg)
389  end if
390 
391 !metaGGA contributions are not taken into account here
392  mgga=0
393 
394 !if you uncomment the following line, you will have to modify
395 !the original function call to pass in gmet and gsqcut
396 !call filterpot(cplex,gmet,gsqcut,nfft,ngfft,2,qphon,rhor1tmp)
397 
398 !Compute the gradients of the first-order density
399 !rho1now(:,:,1) contains the first-order density, and
400 !rho1now(:,:,2:4) contains the gradients of the first-order density
401  ishift=0 ; ngrad=2
402  ABI_ALLOCATE(rho1now,(cplex*nfft,nspden,ngrad*ngrad))
403  call xcden(cplex,gprimd,ishift,mpi_enreg,nfft,ngfft,ngrad,nspden,paral_kgb,qphon,rhor1tmp,rho1now)
404 
405 !Calculate the 1st-order contribution to grad(n) from the strain derivative
406 !  acting on the gradient operator acting on the GS charge density,
407 !Simply use the following formula:
408 !   (dGprim/ds_alpha_beta)(i,j) = -half.( delta_alpha,i Gprim(beta,j) + delta_beta,i Gprim(alpha,j) )
409 !To finally get:
410 !   (nabla)^(alpha,beta)_i[n] = -half ( delta_alpha,i nabla_beta[n] + delta_beta,i nabla_alpha[n] )
411  ABI_ALLOCATE(rhodgnow,(cplex*nfft,nspden,3))
412  rhodgnow(1:nfft,1:nspden,1:3)=zero
413  if (nspden==1) then
414    if (istr==1) rhodgnow(1:nfft,1,1)=-     kxc(1:nfft,5)
415    if (istr==2) rhodgnow(1:nfft,1,2)=-     kxc(1:nfft,6)
416    if (istr==3) rhodgnow(1:nfft,1,3)=-     kxc(1:nfft,7)
417    if (istr==4) rhodgnow(1:nfft,1,2)=-half*kxc(1:nfft,7)
418    if (istr==4) rhodgnow(1:nfft,1,3)=-half*kxc(1:nfft,6)
419    if (istr==5) rhodgnow(1:nfft,1,1)=-half*kxc(1:nfft,7)
420    if (istr==5) rhodgnow(1:nfft,1,3)=-half*kxc(1:nfft,5)
421    if (istr==6) rhodgnow(1:nfft,1,1)=-half*kxc(1:nfft,6)
422    if (istr==6) rhodgnow(1:nfft,1,2)=-half*kxc(1:nfft,5)
423  else
424    if (istr==1) rhodgnow(1:nfft,1,1)=-     kxc(1:nfft,15)
425    if (istr==2) rhodgnow(1:nfft,1,2)=-     kxc(1:nfft,17)
426    if (istr==3) rhodgnow(1:nfft,1,3)=-     kxc(1:nfft,19)
427    if (istr==4) rhodgnow(1:nfft,1,2)=-half*kxc(1:nfft,19)
428    if (istr==4) rhodgnow(1:nfft,1,3)=-half*kxc(1:nfft,17)
429    if (istr==5) rhodgnow(1:nfft,1,1)=-half*kxc(1:nfft,19)
430    if (istr==5) rhodgnow(1:nfft,1,3)=-half*kxc(1:nfft,15)
431    if (istr==6) rhodgnow(1:nfft,1,1)=-half*kxc(1:nfft,17)
432    if (istr==6) rhodgnow(1:nfft,1,2)=-half*kxc(1:nfft,15)
433    if (istr==1) rhodgnow(1:nfft,2,1)=-     (kxc(1:nfft,14)-kxc(1:nfft,15))
434    if (istr==2) rhodgnow(1:nfft,2,2)=-     (kxc(1:nfft,16)-kxc(1:nfft,17))
435    if (istr==3) rhodgnow(1:nfft,2,3)=-     (kxc(1:nfft,18)-kxc(1:nfft,19))
436    if (istr==4) rhodgnow(1:nfft,2,2)=-half*(kxc(1:nfft,18)-kxc(1:nfft,19))
437    if (istr==4) rhodgnow(1:nfft,2,3)=-half*(kxc(1:nfft,16)-kxc(1:nfft,17))
438    if (istr==5) rhodgnow(1:nfft,2,1)=-half*(kxc(1:nfft,18)-kxc(1:nfft,19))
439    if (istr==5) rhodgnow(1:nfft,2,3)=-half*(kxc(1:nfft,14)-kxc(1:nfft,15))
440    if (istr==6) rhodgnow(1:nfft,2,1)=-half*(kxc(1:nfft,16)-kxc(1:nfft,17))
441    if (istr==6) rhodgnow(1:nfft,2,2)=-half*(kxc(1:nfft,14)-kxc(1:nfft,15))
442  end if
443 
444 !Add to the gradients of the first-order density
445  do ii=1,3
446    do ispden=1,nspden
447      rhodgnow(1:nfft,ispden,ii)=str_scale*rhodgnow(1:nfft,ispden,ii)
448      rho1now(1:nfft,ispden,1+ii)=rho1now(1:nfft,ispden,1+ii)+rhodgnow(1:nfft,ispden,ii)
449    end do
450  end do
451 
452 !rho1now(:,:,1) contains the 1st-order density, and rho1now(:,:,2:4) contains the grads of the 1st-order density
453 
454 !Apply the XC kernel
455  nspgrad=2; if (nspden==2) nspgrad=5
456  ABI_ALLOCATE(dnexcdn,(cplex*nfft,nspgrad))
457 
458 !== Non polarized
459  if (nspden==1) then
460    do ir=1,nfft
461      r0(:)=kxc(ir,5:7) ; r1(:)=rho1now(ir,1,2:4)
462      gradrho_gradrho1=dot_product(r0,r1)
463      dnexcdn(ir,1)=kxc(ir,1)*rho1now(ir,1,1) + kxc(ir,3)*gradrho_gradrho1
464      coeff_grho=kxc(ir,3)*rho1now(ir,1,1) + kxc(ir,4)*gradrho_gradrho1
465 !    Grad strain derivative contribution enters the following term with a
466 !    factor of two compared to above terms, so add it again.
467      r1(:)=r1(:)+rhodgnow(ir,1,1:3)
468 !    Reuse the storage in rho1now
469      rho1now(ir,1,2:4)=r1(:)*kxc(ir,2)+r0(:)*coeff_grho
470    end do
471 
472 !== Spin-polarized
473  else ! nspden==2
474    do ir=1,nfft
475      do ii=1,3  ! grad of spin-up ans spin_dn GS rho
476        r0_up(ii)=kxc(ir,13+2*ii);r0_dn(ii)=kxc(ir,12+2*ii)-kxc(ir,13+2*ii)
477      end do
478      r0(:)=r0_up(:)+r0_dn(:)      ! grad of GS rho
479      r1_up(:)=rho1now(ir,1,2:4)   ! grad of spin-up rho1
480      r1_dn(:)=rho1now(ir,2,2:4)   ! grad of spin-down rho1
481      r1(:)=r1_up(:)+r1_dn(:)      ! grad of GS rho1
482      gradrho_gradrho1_up=dot_product(r0_up,r1_up)
483      gradrho_gradrho1_dn=dot_product(r0_dn,r1_dn)
484      gradrho_gradrho1   =dot_product(r0,r1)
485 
486      dnexcdn(ir,1)=kxc(ir, 1)*rho1now(ir,1,1)     &
487 &     +kxc(ir, 2)*rho1now(ir,2,1)     &
488 &     +kxc(ir, 6)*gradrho_gradrho1_up &
489 &     +kxc(ir,11)*gradrho_gradrho1
490      dnexcdn(ir,2)=kxc(ir, 3)*rho1now(ir,2,1)     &
491 &     +kxc(ir, 2)*rho1now(ir,1,1)     &
492 &     +kxc(ir, 7)*gradrho_gradrho1_dn &
493 &     +kxc(ir,12)*gradrho_gradrho1
494      coeff_grho_corr=kxc(ir,11)*rho1now(ir,1,1) &
495 &     +kxc(ir,12)*rho1now(ir,2,1) &
496 &     +kxc(ir,13)*gradrho_gradrho1
497      coeff_grho_up=kxc(ir,6)*rho1now(ir,1,1)+kxc(ir,8)*gradrho_gradrho1_up
498      coeff_grho_dn=kxc(ir,7)*rho1now(ir,2,1)+kxc(ir,9)*gradrho_gradrho1_dn
499 
500 !    Grad strain derivative contribution enters the following term with a
501 !    factor of two compared to above terms, so add it again.
502      r1_up(:)=r1_up(:)+rhodgnow(ir,1,1:3)
503      r1_dn(:)=r1_dn(:)+rhodgnow(ir,2,1:3)
504 
505 !    Reuse the storage in rho1now
506      rho1now(ir,1,2:4)=(kxc(ir,4)+kxc(ir,10))*r1_up(:) &
507 &     +kxc(ir,10)            *r1_dn(:) &
508 &     +coeff_grho_up         *r0_up(:) &
509 &     +coeff_grho_corr       *r0(:)
510      rho1now(ir,2,2:4)=(kxc(ir,5)+kxc(ir,10))*r1_dn(:) &
511 &     +kxc(ir,10)            *r1_up(:) &
512 &     +coeff_grho_dn         *r0_dn(:) &
513 &     +coeff_grho_corr       *r0(:)
514    end do
515 
516  end if ! nspden
517  ABI_DEALLOCATE(rhodgnow)
518 
519  vxc1(:,:)=zero
520  call xcpot(cplex,dnexcdn,gprimd,ishift,mgga,mpi_enreg,nfft,ngfft,ngrad,nspden,&
521 & nspgrad,paral_kgb,qphon,rho1now,vxc1)
522 
523 !if you uncomment the following line, you will have to modify
524 !the original function call to pass in gmet and gsqcut
525 !call filterpot(cplex,gmet,gsqcut,nfft,ngfft,nspden,qphon,vxc1)
526 
527  ABI_DEALLOCATE(dnexcdn)
528  ABI_DEALLOCATE(rho1now)
529 
530 end subroutine dfpt_mkvxcstrgga

ABINIT/m_dfpt_mkvxcstr [ Modules ]

[ Top ] [ Modules ]

NAME

  m_dfpt_mkvxcstr

FUNCTION

COPYRIGHT

  Copyright (C) 2001-2018 ABINIT group (DRH,XG)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_dfpt_mkvxcstr
28 
29  use defs_basis
30  use defs_abitypes
31  use m_errors
32  use m_abicore
33 
34  use m_time,      only : timab
35  use m_symtk,     only : matr3inv
36  use m_xctk,      only : xcden, xcpot
37 
38  implicit none
39 
40  private