TABLE OF CONTENTS


ABINIT/dfpt_rhotov [ Functions ]

[ Top ] [ Functions ]

NAME

 dfpt_rhotov

FUNCTION

 This routine is called to compute, from a given 1st-order total density
   - the trial (local) 1st-order potential and/or the residual potential,
   - some contributions to the 2nd-order energy

INPUTS

  cplex: if 1, real space 1-order WF on FFT grid are REAL; if 2, COMPLEX
  gsqcut=cutoff on (k+G)^2 (bohr^-2)
  idir=direction of atomic displacement (=1,2 or 3 : displacement of atom ipert along the 1st, 2nd or 3rd axis).
  ipert=type of the perturbation
  ixc= choice of exchange-correlation scheme
  kxc(nfft,nkxc)=exchange-correlation kernel
  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
  nhat1gr(cplex*nfft,nspden,3*nhat1grdim)= -PAW only- gradients of 1st-order compensation density
  nhat1grdim= -PAW only- 1 if nhat1gr array is used ; 0 otherwise
  nkxc=second dimension of the array kxc, see rhotoxc.f for a description
  nspden=number of spin-density components
  n3xccc=dimension of xccc3d1 ; 0 if no XC core correction is used
  optene=0: the contributions to the 2nd order energy are not computed
         1: the contributions to the 2nd order energy are computed
  optres=0: the trial potential residual is computed ; the input potential value is kept
         1: the new value of the trial potential is computed in place of the input value
  paral_kgb=flag controlling (k,g,bands) parallelization
  qphon(3)=reduced coordinates for the phonon wavelength
  rhog(2,nfft)=array for Fourier transform of GS electron density
  rhog1(2,nfft)=RF electron density in reciprocal space
  rhor(nfft,nspden)=array for GS electron density in electrons/bohr**3.
  rhor1(cplex*nfft,nspden)=RF electron density in real space (electrons/bohr**3).
  rprimd(3,3)=dimensional primitive translations in real space (bohr)
  ucvol=unit cell volume in ($\textrm{bohr}^{3}$)
  usepaw= 0 for non paw calculation; =1 for paw calculation
  usexcnhat= -PAW only- flag controling use of compensation density in Vxc
  vpsp1(cplex*nfft)=first-order derivative of the ionic potential
  xccc3d1(cplex*n3xccc)=3D change in core charge density, see n3xccc

OUTPUT

  vhartr1(cplex*nfft)=1-order Hartree potential (not output if size=0)
  vxc1(cplex*nfft,nspden)= 1st-order XC potential (not output if size=0)
  ==== if optene==1
    ehart01=inhomogeneous 1st-order Hartree part of 2nd-order total energy
    ehart1=1st-order Hartree part of 2nd-order total energy
    exc1=1st-order exchange-correlation part of 2nd-order total energy
    elpsp1=1st-order local pseudopot. part of 2nd-order total energy.
  ==== if optres==0
    vresid1(cplex*nfft,nspden)=potential residual
    vres2=square of the norm of the residual

SIDE EFFECTS

  ==== if optres==1
    vtrial1(cplex*nfft,nspden)= new value of 1st-order trial potential

PARENTS

      dfpt_scfcv

CHILDREN

      dfpt_mkvxc,dfpt_mkvxc_noncoll,dfpt_mkvxcstr,dfpt_v1zeeman,dotprod_vn
      hartre,hartrestr,sqnorm_v,timab

SOURCE

120  subroutine dfpt_rhotov(cplex,ehart01,ehart1,elpsp1,exc1,elmag1,gsqcut,idir,ipert,&
121 &           ixc,kxc,mpi_enreg,natom,nfft,ngfft,nhat,nhat1,nhat1gr,nhat1grdim,nkxc,nspden,n3xccc,&
122 &           optene,optres,paral_kgb,qphon,rhog,rhog1,rhor,rhor1,rprimd,ucvol,&
123 &           usepaw,usexcnhat,vhartr1,vpsp1,vresid1,vres2,vtrial1,vxc,vxc1,xccc3d1,ixcrot)
124 
125 
126 !This section has been created automatically by the script Abilint (TD).
127 !Do not modify the following lines by hand.
128 #undef ABI_FUNC
129 #define ABI_FUNC 'dfpt_rhotov'
130 !End of the abilint section
131 
132  implicit none
133 
134 !Arguments ------------------------------------
135 !scalars
136  integer,intent(in) :: cplex,idir,ipert,ixc,n3xccc,natom,nfft,nhat1grdim,nkxc,nspden
137  integer,intent(in) :: optene,optres,paral_kgb,usepaw,usexcnhat,ixcrot
138  real(dp),intent(in) :: gsqcut,ucvol
139  real(dp),intent(inout) :: ehart01 !vz_i
140  real(dp),intent(out) :: vres2
141  type(MPI_type),intent(in) :: mpi_enreg
142 !arrays
143  real(dp),intent(in) :: kxc(nfft,nkxc)
144  real(dp),intent(in) :: vxc(nfft,nspden)
145  real(dp),intent(in) :: nhat(nfft,nspden)
146  real(dp),intent(in) :: nhat1(cplex*nfft,nspden)  !vz_d
147  real(dp),intent(in) :: nhat1gr(cplex*nfft,nspden,3*nhat1grdim)
148  real(dp),intent(in) :: qphon(3),rhog(2,nfft)
149  real(dp),intent(in) :: rhog1(2,nfft)
150  real(dp),target,intent(in) :: rhor(nfft,nspden),rhor1(cplex*nfft,nspden)
151  real(dp),intent(in) :: rprimd(3,3),vpsp1(cplex*nfft)
152  real(dp),intent(in) :: xccc3d1(cplex*n3xccc)
153  real(dp),intent(inout) :: vtrial1(cplex*nfft,nspden),elpsp1,ehart1,exc1,elmag1
154  real(dp),intent(out) :: vresid1(cplex*nfft,nspden)
155  real(dp),target,intent(out) :: vhartr1(:),vxc1(:,:)
156 
157 !Local variables-------------------------------
158 !scalars
159  integer :: ifft,ispden,nfftot,option
160  integer :: optnc,nkxc_cur
161  logical :: vhartr1_allocated,vxc1_allocated
162  real(dp) :: doti,elpsp10
163 !arrays
164  integer,intent(in)   :: ngfft(18)
165  real(dp)             :: tsec(20)
166  real(dp),allocatable :: rhor1_nohat(:,:),vhartr01(:),vxc1val(:,:)
167  real(dp),pointer     :: rhor1_(:,:),vhartr1_(:),vxc1_(:,:),v1zeeman(:,:)
168 
169 ! *********************************************************************
170 
171  call timab(157,1,tsec)
172 
173  !FR EB SPr
174  if (nspden==4) then
175    if(usepaw==1) then
176      MSG_ERROR('DFPT with nspden=4 works only for norm-conserving psp!')
177    end if
178  end if
179 
180 !Get size of FFT grid
181  nfftot=ngfft(1)*ngfft(2)*ngfft(3)
182 
183 !Eventually allocate temporary memory space
184  vhartr1_allocated=(size(vhartr1)>0)
185  if (vhartr1_allocated) then
186    vhartr1_ => vhartr1
187  else
188    ABI_ALLOCATE(vhartr1_,(cplex*nfft))
189  end if
190  vxc1_allocated=(size(vxc1)>0)
191  if (vxc1_allocated) then
192    vxc1_ => vxc1
193  else
194    ABI_ALLOCATE(vxc1_,(cplex*nfft,nspden))
195  end if
196 
197 !If needed, store pseudo density without charge compensation
198  if (usepaw==1.and.usexcnhat==0) then
199    ABI_ALLOCATE(rhor1_,(cplex*nfft,nspden))
200    rhor1_(:,:)=rhor1(:,:)-nhat1(:,:)
201  else
202    rhor1_ => rhor1
203  end if
204 
205 
206  if(ipert==natom+5)then
207    ABI_ALLOCATE(v1zeeman,(cplex*nfft,nspden))
208    call dfpt_v1zeeman(nspden,nfft,cplex,idir,v1zeeman)
209  end if
210 
211 !------ Compute 1st-order Hartree potential (and energy) ----------------------
212 
213  call hartre(cplex,gsqcut,0,mpi_enreg,nfft,ngfft,paral_kgb,rhog1,rprimd,vhartr1_,qpt=qphon)
214 
215  if (optene>0) then
216    call dotprod_vn(cplex,rhor1,ehart1,doti,nfft,nfftot,1,1,vhartr1_,ucvol)
217  end if
218 
219  if (optene>0) ehart01=zero
220  if(ipert==natom+3 .or. ipert==natom+4) then
221    ABI_ALLOCATE(vhartr01,(cplex*nfft))
222    call hartrestr(gsqcut,idir,ipert,mpi_enreg,natom,nfft,ngfft,paral_kgb,rhog,rprimd,vhartr01)
223    if (optene>0) then
224      call dotprod_vn(cplex,rhor1,ehart01,doti,nfft,nfftot,1,1,vhartr01,ucvol)
225      ehart01=two*ehart01
226      ehart1=ehart1+ehart01
227    end if
228 !  Note that there is a factor 2.0_dp difference with the similar GS formula
229    vhartr1_(:)=vhartr1_(:)+vhartr01(:)
230 
231    ABI_DEALLOCATE(vhartr01)
232  end if
233 
234 !------ Compute 1st-order XC potential (and energy) ----------------------
235 !(including the XC core correction)
236 
237 !Compute Vxc^(1) (with or without valence contribution according to options)
238  option=0;if (optene==0) option=1
239  if(ipert==natom+3.or.ipert==natom+4) then
240    call dfpt_mkvxcstr(cplex,idir,ipert,kxc,mpi_enreg,natom,nfft,ngfft,nhat,&
241 &   nhat1,nkxc,nspden,n3xccc,option,paral_kgb,qphon,rhor,rhor1,rprimd,&
242 &   usepaw,usexcnhat,vxc1_,xccc3d1)
243  else
244 ! FR EB non-collinear magnetism
245 ! the second nkxc should be nkxc_cur (see 67_common/nres2vres.F90)
246    if (nspden==4) then
247      optnc=1
248      nkxc_cur=nkxc ! TODO: remove nkxc_cur?
249 
250      call dfpt_mkvxc_noncoll(cplex,ixc,kxc,mpi_enreg,nfft,ngfft,nhat,usepaw,nhat1,usepaw,nhat1gr,nhat1grdim,nkxc,&
251 &     nspden,n3xccc,optnc,option,paral_kgb,qphon,rhor,rhor1,rprimd,usexcnhat,vxc,vxc1_,xccc3d1,ixcrot=ixcrot)
252 
253    else
254      call dfpt_mkvxc(cplex,ixc,kxc,mpi_enreg,nfft,ngfft,nhat1,usepaw,nhat1gr,nhat1grdim,nkxc,&
255 &     nspden,n3xccc,option,paral_kgb,qphon,rhor1,rprimd,usexcnhat,vxc1_,xccc3d1)
256    end if !nspden==4
257  end if
258 
259 !Compute local contribution to 2nd-order energy (includes Vxc and Vpsp and Vmag)
260  if (optene>0) then
261    if (usepaw==0) then
262      call dotprod_vn(cplex,rhor1,elpsp10,doti,nfft,nfftot,nspden,1,vxc1_,ucvol)
263      call dotprod_vn(cplex,rhor1,elpsp1 ,doti,nfft,nfftot,1     ,1,vpsp1,ucvol)
264      if (ipert==natom+5) then
265        call dotprod_vn(cplex,rhor1,elmag1 ,doti,nfft,nfftot,nspden,1,v1zeeman,ucvol)
266      end if
267    else
268      if (usexcnhat/=0) then
269        ABI_ALLOCATE(rhor1_nohat,(cplex*nfft,1))
270        rhor1_nohat(:,1)=rhor1(:,1)-nhat1(:,1)
271        call dotprod_vn(cplex,rhor1      ,elpsp10,doti,nfft,nfftot,nspden,1,vxc1_,ucvol)
272        call dotprod_vn(cplex,rhor1_nohat,elpsp1 ,doti,nfft,nfftot,1     ,1,vpsp1,ucvol)
273        ABI_DEALLOCATE(rhor1_nohat)
274      else
275        call dotprod_vn(cplex,rhor1_,elpsp10,doti,nfft,nfftot,nspden,1,vxc1_,ucvol)
276        call dotprod_vn(cplex,rhor1_,elpsp1 ,doti,nfft,nfftot,1     ,1,vpsp1,ucvol)
277      end if
278    end if
279 
280 !  Note that there is a factor 2 difference with the similar GS formula
281    elpsp1=two*(elpsp1+elpsp10)
282  end if
283 
284 
285 !Compute XC valence contribution exc1 and complete eventually Vxc^(1)
286  if (optene>0) then
287    ABI_ALLOCATE(vxc1val,(cplex*nfft,nspden))
288    vxc1val=zero
289    option=2
290 !FR SPr EB non-collinear magnetism
291    if (nspden==4) then
292      optnc=1
293      nkxc_cur=nkxc
294      call dfpt_mkvxc_noncoll(cplex,ixc,kxc,mpi_enreg,nfft,ngfft,nhat,usepaw,nhat1,usepaw,nhat1gr,nhat1grdim,nkxc,&
295 &     nspden,n3xccc,optnc,option,paral_kgb,qphon,rhor,rhor1,rprimd,usexcnhat,vxc,vxc1val,xccc3d1,ixcrot=ixcrot)
296    else
297      call dfpt_mkvxc(cplex,ixc,kxc,mpi_enreg,nfft,ngfft,nhat1,usepaw,nhat1gr,nhat1grdim,nkxc,&
298 &     nspden,n3xccc,option,paral_kgb,qphon,rhor1,rprimd,usexcnhat,vxc1val,xccc3d1)
299    end if !nspden==4
300    vxc1_(:,:)=vxc1_(:,:)+vxc1val(:,:)
301    call dotprod_vn(cplex,rhor1_,exc1,doti,nfft,nfftot,nspden,1,vxc1val,ucvol)
302    ABI_DEALLOCATE(vxc1val)
303  end if
304 
305  if (usepaw==1.and.usexcnhat==0) then
306    ABI_DEALLOCATE(rhor1_)
307  end if
308 
309 !DEBUG (do not take away)
310 !Compute NSC energy ensc1 associated with rhor1 in vtrial1, for debugging purposes
311 !call dotprod_vn(cplex,rhor1,ensc1,doti,nfft,nfftot,nspden,1,vtrial1,ucvol)
312 !write(std_out,*)' ek0+eeig0+eloc0=',ek0+eeig0+eloc0
313 !write(std_out,*)' ensc1=',ensc1
314 !Compute NSC energy associated with vtrial1, for debugging purposes
315 !call dotprod_vn(cplex,rhor1,ensc1,doti,mpi_enreg,nfft,nfftot,nspden,1,vtrial1,ucvol)
316 !ensc1=ensc1+half*enl1
317 !write(std_out,*)' dfpt_rhotov : check NSC energy, diff=',&
318 !&  ek0+edocc+eeig0+eloc0+enl0+ensc1
319 !write(std_out,*)' evarNSC=',ek0+edocc+eeig0+eloc0+enl0
320 !write(std_out,*)' ensc1,exc1=',ensc1,exc1
321 !ENDDEBUG
322 
323 !Here, vhartr1 contains Hartree potential, vpsp1 contains local psp,
324 !while vxc1 contain xc potential
325 
326 !------ Produce residual vector and square of norm of it -------------
327 !(only if requested ; if optres==0)
328  if (optres==0) then
329 !$OMP PARALLEL DO COLLAPSE(2)
330    do ispden=1,min(nspden,2)
331      do ifft=1,cplex*nfft
332        vresid1(ifft,ispden)=vhartr1_(ifft)+vxc1_(ifft,ispden)+vpsp1(ifft)-vtrial1(ifft,ispden)
333      end do
334    end do
335    if(nspden==4)then
336 !$OMP PARALLEL DO COLLAPSE(2)
337      do ispden=3,4
338        do ifft=1,cplex*nfft
339          vresid1(ifft,ispden)=vxc1_(ifft,ispden)-vtrial1(ifft,ispden)
340        end do
341      end do
342    end if
343 
344    if (ipert==natom+5) then
345      vresid1 = vresid1 + v1zeeman
346    end if
347 !  Compute square norm vres2 of potential residual vresid
348    call sqnorm_v(cplex,nfft,vres2,nspden,optres,vresid1)
349 
350  else
351 
352 !  ------ Produce new value of trial potential-------------
353 !  (only if requested ; if optres==1)
354 
355 !$OMP PARALLEL DO COLLAPSE(2)
356    do ispden=1,min(nspden,2)
357      do ifft=1,cplex*nfft
358        vtrial1(ifft,ispden)=vhartr1_(ifft)+vxc1_(ifft,ispden)+vpsp1(ifft)
359      end do
360    end do
361    if(nspden==4)then
362 !$OMP PARALLEL DO COLLAPSE(2)
363      do ispden=3,4
364        do ifft=1,cplex*nfft
365          vtrial1(ifft,ispden)=vxc1_(ifft,ispden)
366        end do
367      end do
368    end if
369 
370    if (ipert==natom+5) then
371      vtrial1 = vtrial1 + v1zeeman
372    end if
373 
374  end if
375 
376 !Release temporary memory space
377  if (.not.vhartr1_allocated) then
378    ABI_DEALLOCATE(vhartr1_)
379  end if
380  if (.not.vxc1_allocated) then
381    ABI_DEALLOCATE(vxc1_)
382  end if
383 
384  if (ipert==natom+5) then
385    ABI_DEALLOCATE(v1zeeman)
386  end if
387 
388  call timab(157,2,tsec)
389 
390 end subroutine dfpt_rhotov

ABINIT/dfpt_v1zeeman [ Functions ]

[ Top ] [ Functions ]

NAME

  dfpt_v1zeeman

FUNCTION

  Calculate 1st order Zeeman potential = -vec{\sigma}.\vec{b}, where
  sigma is the vector of Pauli matrices and \vec{b} is the unit
  vector indicating the perturbing field direction.

INPUTS

  nspden = number of density matrix components
  nfft   = numbder of fft grid points
  cplex  = complex or real density matrix
  idir   = direction of the perturbing field in Cartesian frame
           1: along x
           2: along y
           3: along z
           4: identity matrix at each fft point is returned (for density-density response)

OUTPUT

  v1zeeman(nfft*cplex,nspden)= 1st order Zeeman potential, or Identity matrix (electrostatic potential) for idir=4

SIDE EFFECTS

NOTES

  The definition of components of the potential matrix differ depending on cplex
  for nspden=4:
  For cplex=1, the potential is defined as (V_upup,V_dndn,Re[V_updn],Im[V_updn])
  For cplex=2, the definition is (V_upup,V_dndn,V_updn,i.V_updn)

PARENTS

      dfpt_rhotov

CHILDREN

SOURCE

430 subroutine dfpt_v1zeeman(nspden,nfft,cplex,idir,v1zeeman)
431 
432 
433 !This section has been created automatically by the script Abilint (TD).
434 !Do not modify the following lines by hand.
435 #undef ABI_FUNC
436 #define ABI_FUNC 'dfpt_v1zeeman'
437 !End of the abilint section
438 
439  implicit none
440 
441 !Arguments ------------------------------------
442  integer , intent(in)    :: idir,nfft,cplex,nspden
443  real(dp), intent(inout) :: v1zeeman(cplex*nfft,nspden)
444 
445 !Local variables-------------------------------
446  integer :: ifft
447 !character(len=500) :: msg
448 
449 ! *************************************************************************
450 
451  DBG_ENTER("COLL")
452 
453 ! if (option/=1 .and. option/=2 ) then
454 !   write(msg,'(3a,i0)')&
455 !&   'The argument option should be 1 or 2,',ch10,&
456 !&   'however, option=',option
457 !   MSG_BUG(msg)
458 ! end if
459 !
460 ! if (sizein<1) then
461 !   write(msg,'(3a,i0)')&
462 !&   'The argument sizein should be a positive number,',ch10,&
463 !&   'however, sizein=',sizein
464 !   MSG_ERROR(msg)
465 ! end if
466 
467  DBG_EXIT("COLL")
468 
469  select case(cplex)
470  case(1)
471    if (nspden==4) then
472      if(idir==3)then       ! Zeeman field along the 3rd axis (z)
473        v1zeeman(:,1)=-0.5d0
474        v1zeeman(:,2)=+0.5d0
475        v1zeeman(:,3)= 0.0d0
476        v1zeeman(:,4)= 0.0d0
477      else if(idir==2)then  ! Zeeman field along the 2nd axis (y)
478        v1zeeman(:,1)= 0.0d0
479        v1zeeman(:,2)= 0.0d0
480        v1zeeman(:,3)= 0.0d0
481        v1zeeman(:,4)=+0.5d0
482      else                  ! Zeeman field along the 1st axis (x)
483        v1zeeman(:,1)= 0.0d0
484        v1zeeman(:,2)= 0.0d0
485        v1zeeman(:,3)=-0.5d0
486        v1zeeman(:,4)= 0.0d0
487      end if
488    else if (nspden==2) then
489      v1zeeman(:,1)=-0.5e0
490      v1zeeman(:,2)= 0.5e0
491    else
492      v1zeeman(:,1)= 0.0e0
493    end if
494  case(2)
495    if (nspden==2) then
496      do ifft=1,nfft
497        v1zeeman(2*ifft-1,1)  =-0.5e0
498        v1zeeman(2*ifft  ,1)  = 0.0e0
499        v1zeeman(2*ifft-1,2)  = 0.5e0
500        v1zeeman(2*ifft  ,2)  = 0.0e0
501      end do
502    else if (nspden==4) then
503      select case(idir)
504      case(1) !along x, v1=-sigma_x
505        do ifft=1,nfft
506          v1zeeman(2*ifft-1,1)= 0.0e0 !Re[V^11]
507          v1zeeman(2*ifft  ,1)= 0.0e0 !Im[V^11]
508          v1zeeman(2*ifft-1,2)= 0.0e0 !Re[V^22]
509          v1zeeman(2*ifft  ,2)= 0.0e0 !Im[V^22]
510          v1zeeman(2*ifft-1,3)=-0.5e0 !Re[V^12]
511          v1zeeman(2*ifft  ,3)= 0.0e0 !Im[V^12]
512          v1zeeman(2*ifft-1,4)= 0.0e0 !Re[i.V^21]=Im[V^12]
513          v1zeeman(2*ifft  ,4)=-0.5e0 !Im[i.V^21]=Re[V^12]
514        end do
515      case(2) !along y, v1 = -sigma_y
516        do ifft=1,nfft
517          v1zeeman(2*ifft-1,1)= 0.0e0 !Re[V^11]
518          v1zeeman(2*ifft  ,1)= 0.0e0 !Im[V^11]
519          v1zeeman(2*ifft-1,2)= 0.0e0 !Re[V^22]
520          v1zeeman(2*ifft  ,2)= 0.0e0 !Im[V^22]
521          v1zeeman(2*ifft-1,3)= 0.0e0 !Re[V^12]
522          v1zeeman(2*ifft  ,3)=+0.5e0 !Im[V^12]
523          v1zeeman(2*ifft-1,4)=+0.5e0 !Re[i.V^21]=Im[V^12]
524          v1zeeman(2*ifft  ,4)= 0.0e0 !Im[i.V^21]=Re[V^12]
525        end do
526      case(3)
527        do ifft=1,nfft
528          v1zeeman(2*ifft-1,1)=-0.5e0 !Re[V^11]
529          v1zeeman(2*ifft  ,1)= 0.0e0 !Im[V^11]
530          v1zeeman(2*ifft-1,2)= 0.5e0 !Re[V^22]
531          v1zeeman(2*ifft  ,2)= 0.0e0 !Im[V^22]
532          v1zeeman(2*ifft-1,3)= 0.0e0 !Re[V^12]
533          v1zeeman(2*ifft  ,3)= 0.0e0 !Im[V^12]
534          v1zeeman(2*ifft-1,4)= 0.0e0 !Re[i.V^21]
535          v1zeeman(2*ifft  ,4)= 0.0e0 !Im[i.V^21]
536        end do
537      end select
538    end if
539  end select !cplex
540 
541 end subroutine dfpt_v1zeeman

ABINIT/m_dfpt_rhotov [ Modules ]

[ Top ] [ Modules ]

NAME

  m_dfpt_rhotov

FUNCTION

COPYRIGHT

  Copyright (C) 1999-2018 ABINIT group (XG, DRH, MT, SPr)
  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_rhotov
28 
29  use defs_basis
30  use defs_abitypes
31  use m_abicore
32  use m_errors
33  use m_cgtools
34 
35  use m_time,        only : timab
36  use m_spacepar,    only : hartrestr, hartre
37  use m_dfpt_mkvxc,    only : dfpt_mkvxc, dfpt_mkvxc_noncoll
38  use m_dfpt_mkvxcstr, only : dfpt_mkvxcstr
39 
40  implicit none
41 
42  private