TABLE OF CONTENTS
ABINIT/dfpt_mkvxcstr [ 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 non_magnetic_xc= if true, handle density/potential as non-magnetic (even if it is) 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)
SOURCE
87 subroutine dfpt_mkvxcstr(cplex,idir,ipert,kxc,mpi_enreg,natom,nfft,ngfft,nhat,nhat1,& 88 & nkxc,non_magnetic_xc,nspden,n3xccc,option,qphon,& 89 & rhor,rhor1,rprimd,usepaw,usexcnhat,vxc1,xccc3d1) 90 91 !Arguments ------------------------------------ 92 !scalars 93 integer,intent(in) :: cplex,idir,ipert,n3xccc,natom,nfft,nkxc,nspden,option 94 integer,intent(in) :: usepaw,usexcnhat 95 logical,intent(in) :: non_magnetic_xc 96 type(MPI_type),intent(in) :: mpi_enreg 97 !arrays 98 integer,intent(in) :: ngfft(18) 99 real(dp),target,intent(in) :: nhat(nfft,nspden) 100 real(dp),target,intent(in) :: nhat1(cplex*nfft,nspden) 101 real(dp),intent(in) :: kxc(nfft,nkxc),qphon(3) 102 real(dp),target,intent(in) :: rhor(nfft,nspden),rhor1(cplex*nfft,nspden) 103 real(dp),intent(in) :: rprimd(3,3) 104 real(dp),intent(in) :: xccc3d1(cplex*n3xccc) 105 real(dp),intent(out) :: vxc1(cplex*nfft,nspden) 106 107 !Local variables------------------------------- 108 !scalars 109 integer :: ii,ir,istr 110 real(dp) :: rho1_dn,rho1_up,spin_scale,str_scale 111 character(len=500) :: message 112 !arrays 113 real(dp) :: gprimd(3,3),tsec(2) 114 real(dp),allocatable :: rhor1tmp(:,:),rhowk1(:,:) 115 real(dp),pointer :: rhor_(:,:),rhor1_(:,:) 116 117 ! ************************************************************************* 118 119 call timab(181,1,tsec) 120 121 if(nspden/=1 .and. nspden/=2) then 122 message = ' dfpt_mkvxc, Only for nspden==1 and 2.' 123 ABI_BUG(message) 124 end if 125 126 if (usepaw==1.and.usexcnhat==0) then 127 ABI_MALLOC(rhor_,(nfft,nspden)) 128 rhor_(:,:)=rhor(:,:)-nhat(:,:) 129 else 130 rhor_ => rhor 131 end if 132 133 if (usepaw==1.and.usexcnhat==0.and.option==1) then 134 ABI_MALLOC(rhor1_,(nfft,nspden)) 135 rhor1_(:,:)=rhor1(:,:)-nhat1(:,:) 136 else 137 rhor1_ => rhor1 138 end if 139 140 !Inhomogeneous term for diagonal strain 141 ABI_MALLOC(rhowk1,(nfft,nspden)) 142 if(option==0 .or. option==2) then 143 if(ipert==natom+3) then 144 rhowk1(:,:)=-rhor_(:,:) 145 else 146 rhowk1(:,:)=zero 147 end if 148 else if(option==1) then 149 if(ipert==natom+3) then 150 rhowk1(:,:)=rhor1_(:,:)-rhor_(:,:) 151 else 152 rhowk1(:,:)=rhor1_(:,:) 153 end if 154 end if 155 156 if (non_magnetic_xc) then 157 if(nspden==2) rhowk1(:,2)=rhowk1(:,1)*half 158 if(nspden==4) rhowk1(:,2:4)=zero 159 end if 160 161 !Treat first LDA 162 if(nkxc==1.or.nkxc==3)then 163 164 ! Case without non-linear core correction 165 if(n3xccc==0)then 166 167 ! Non-spin-polarized 168 if(nspden==1)then 169 do ir=1,nfft 170 vxc1(ir,1)=kxc(ir,1)*rhowk1(ir,1) 171 end do 172 173 ! Spin-polarized 174 else 175 do ir=1,nfft 176 rho1_dn=rhowk1(ir,1)-rhowk1(ir,2) 177 vxc1(ir,1)=kxc(ir,1)*rhowk1(ir,2)+kxc(ir,2)*rho1_dn 178 vxc1(ir,2)=kxc(ir,2)*rhowk1(ir,2)+kxc(ir,3)*rho1_dn 179 end do 180 end if ! nspden==1 181 182 ! Treat case with non-linear core correction 183 else 184 if(nspden==1)then 185 do ir=1,nfft 186 vxc1(ir,1)=kxc(ir,1)*(rhowk1(ir,1)+xccc3d1(ir)) 187 end do 188 else 189 do ir=1,nfft 190 rho1_dn=rhowk1(ir,1)-rhowk1(ir,2) + xccc3d1(ir)*half 191 rho1_up=rhowk1(ir,2) + xccc3d1(ir)*half 192 vxc1(ir,1)=kxc(ir,1)*rho1_up+kxc(ir,2)*rho1_dn 193 vxc1(ir,2)=kxc(ir,2)*rho1_up+kxc(ir,3)*rho1_dn 194 end do 195 end if ! nspden==1 196 197 end if ! n3xccc==0 198 199 ! Treat GGA 200 else if (nkxc==7.or.nkxc==19) then 201 202 ! Generates gprimd and its strain derivative 203 ! Note that unlike the implicitly symmetric metric tensor strain 204 ! derivatives, we must explicltly symmetrize the strain derivative 205 ! here. 206 call matr3inv(rprimd,gprimd) 207 istr=idir + 3*(ipert-natom-3) 208 if(istr<1 .or. istr>6)then 209 write(message, '(a,i10,a,a,a)' )& 210 & 'Input dir gives istr=',istr,' not allowed.',ch10,& 211 & 'Possible values are 1,2,3,4,5,6 only.' 212 ABI_BUG(message) 213 end if 214 215 ! Rescalling needed for use in dfpt_eltfrxc for elastic tensor (not internal strain tensor). 216 str_scale=one;if(option==2) str_scale=two 217 218 ! Transfer the data to spin-polarized storage 219 ABI_MALLOC(rhor1tmp,(cplex*nfft,nspden)) 220 if(nspden==1)then 221 do ir=1,cplex*nfft 222 rhor1tmp(ir,1)=rhowk1(ir,1) 223 end do 224 else 225 do ir=1,cplex*nfft 226 rho1_dn=rhowk1(ir,1)-rhowk1(ir,2) 227 rhor1tmp(ir,1)=rhowk1(ir,2) 228 rhor1tmp(ir,2)=rho1_dn 229 end do 230 end if ! nspden==1 231 if(n3xccc/=0)then 232 spin_scale=one;if (nspden==2) spin_scale=half 233 do ii=1,nspden 234 do ir=1,cplex*nfft 235 rhor1tmp(ir,ii)=rhor1tmp(ir,ii)+xccc3d1(ir)*spin_scale 236 end do 237 end do 238 end if 239 240 call dfpt_mkvxcstrgga(cplex,gprimd,istr,kxc,mpi_enreg,nfft,ngfft,nkxc,& 241 & nspden,qphon,rhor1tmp,str_scale,vxc1) 242 ABI_FREE(rhor1tmp) 243 244 else 245 ABI_BUG('Invalid nkxc!') 246 247 end if ! LDA or GGA 248 249 if (usepaw==1.and.usexcnhat==0) then 250 ABI_FREE(rhor_) 251 end if 252 if (usepaw==1.and.usexcnhat==0.and.option==1) then 253 ABI_FREE(rhor1_) 254 end if 255 256 ABI_FREE(rhowk1) 257 258 call timab(181,2,tsec) 259 260 end subroutine dfpt_mkvxcstr
ABINIT/dfpt_mkvxcstrgga [ 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)
SOURCE
327 subroutine dfpt_mkvxcstrgga(cplex,gprimd,istr,kxc,mpi_enreg,nfft,ngfft,& 328 & nkxc,nspden,qphon,rhor1tmp,str_scale,vxc1) 329 330 !Arguments ------------------------------------ 331 !scalars 332 integer,intent(in) :: cplex,istr,nfft,nkxc,nspden 333 real(dp) :: str_scale 334 type(MPI_type),intent(in) :: mpi_enreg 335 !arrays 336 integer,intent(in) :: ngfft(18) 337 real(dp),intent(in) :: gprimd(3,3),kxc(nfft,nkxc) 338 real(dp),intent(in) :: qphon(3),rhor1tmp(cplex*nfft,2) 339 real(dp),intent(out) :: vxc1(cplex*nfft,nspden) 340 !Local variables------------------------------- 341 !scalars 342 integer :: ii,ir,ishift,ispden,mgga,ngrad,nspgrad 343 real(dp) :: coeff_grho,coeff_grho_corr,coeff_grho_dn,coeff_grho_up 344 real(dp) :: gradrho_gradrho1,gradrho_gradrho1_dn,gradrho_gradrho1_up 345 character(len=500) :: msg 346 !arrays 347 real(dp) :: r0(3),r0_dn(3),r0_up(3),r1(3),r1_dn(3),r1_up(3) 348 real(dp),allocatable :: dnexcdn(:,:),rho1now(:,:,:),rhodgnow(:,:,:) 349 350 ! ************************************************************************* 351 352 DBG_ENTER("COLL") 353 354 if (nkxc/=12*min(nspden,2)-5) then 355 msg='Wrong nkxc value for GGA!' 356 ABI_BUG(msg) 357 end if 358 if (nspden>2) then 359 msg='Not compatible with non-collinear magnetism!' 360 ABI_ERROR(msg) 361 end if 362 363 !metaGGA contributions are not taken into account here 364 mgga=0 365 366 !if you uncomment the following line, you will have to modify 367 !the original function call to pass in gmet and gsqcut 368 !call filterpot(cplex,gmet,gsqcut,nfft,ngfft,2,qphon,rhor1tmp) 369 370 !Compute the gradients of the first-order density 371 !rho1now(:,:,1) contains the first-order density, and 372 !rho1now(:,:,2:4) contains the gradients of the first-order density 373 ishift=0 ; ngrad=2 374 ABI_MALLOC(rho1now,(cplex*nfft,nspden,ngrad*ngrad)) 375 call xcden(cplex,gprimd,ishift,mpi_enreg,nfft,ngfft,ngrad,nspden,qphon,rhor1tmp,rho1now) 376 377 !Calculate the 1st-order contribution to grad(n) from the strain derivative 378 ! acting on the gradient operator acting on the GS charge density, 379 !Simply use the following formula: 380 ! (dGprim/ds_alpha_beta)(i,j) = -half.( delta_alpha,i Gprim(beta,j) + delta_beta,i Gprim(alpha,j) ) 381 !To finally get: 382 ! (nabla)^(alpha,beta)_i[n] = -half ( delta_alpha,i nabla_beta[n] + delta_beta,i nabla_alpha[n] ) 383 ABI_MALLOC(rhodgnow,(cplex*nfft,nspden,3)) 384 rhodgnow(1:nfft,1:nspden,1:3)=zero 385 if (nspden==1) then 386 if (istr==1) rhodgnow(1:nfft,1,1)=- kxc(1:nfft,5) 387 if (istr==2) rhodgnow(1:nfft,1,2)=- kxc(1:nfft,6) 388 if (istr==3) rhodgnow(1:nfft,1,3)=- kxc(1:nfft,7) 389 if (istr==4) rhodgnow(1:nfft,1,2)=-half*kxc(1:nfft,7) 390 if (istr==4) rhodgnow(1:nfft,1,3)=-half*kxc(1:nfft,6) 391 if (istr==5) rhodgnow(1:nfft,1,1)=-half*kxc(1:nfft,7) 392 if (istr==5) rhodgnow(1:nfft,1,3)=-half*kxc(1:nfft,5) 393 if (istr==6) rhodgnow(1:nfft,1,1)=-half*kxc(1:nfft,6) 394 if (istr==6) rhodgnow(1:nfft,1,2)=-half*kxc(1:nfft,5) 395 else 396 if (istr==1) rhodgnow(1:nfft,1,1)=- kxc(1:nfft,15) 397 if (istr==2) rhodgnow(1:nfft,1,2)=- kxc(1:nfft,17) 398 if (istr==3) rhodgnow(1:nfft,1,3)=- kxc(1:nfft,19) 399 if (istr==4) rhodgnow(1:nfft,1,2)=-half*kxc(1:nfft,19) 400 if (istr==4) rhodgnow(1:nfft,1,3)=-half*kxc(1:nfft,17) 401 if (istr==5) rhodgnow(1:nfft,1,1)=-half*kxc(1:nfft,19) 402 if (istr==5) rhodgnow(1:nfft,1,3)=-half*kxc(1:nfft,15) 403 if (istr==6) rhodgnow(1:nfft,1,1)=-half*kxc(1:nfft,17) 404 if (istr==6) rhodgnow(1:nfft,1,2)=-half*kxc(1:nfft,15) 405 if (istr==1) rhodgnow(1:nfft,2,1)=- (kxc(1:nfft,14)-kxc(1:nfft,15)) 406 if (istr==2) rhodgnow(1:nfft,2,2)=- (kxc(1:nfft,16)-kxc(1:nfft,17)) 407 if (istr==3) rhodgnow(1:nfft,2,3)=- (kxc(1:nfft,18)-kxc(1:nfft,19)) 408 if (istr==4) rhodgnow(1:nfft,2,2)=-half*(kxc(1:nfft,18)-kxc(1:nfft,19)) 409 if (istr==4) rhodgnow(1:nfft,2,3)=-half*(kxc(1:nfft,16)-kxc(1:nfft,17)) 410 if (istr==5) rhodgnow(1:nfft,2,1)=-half*(kxc(1:nfft,18)-kxc(1:nfft,19)) 411 if (istr==5) rhodgnow(1:nfft,2,3)=-half*(kxc(1:nfft,14)-kxc(1:nfft,15)) 412 if (istr==6) rhodgnow(1:nfft,2,1)=-half*(kxc(1:nfft,16)-kxc(1:nfft,17)) 413 if (istr==6) rhodgnow(1:nfft,2,2)=-half*(kxc(1:nfft,14)-kxc(1:nfft,15)) 414 end if 415 416 !Add to the gradients of the first-order density 417 do ii=1,3 418 do ispden=1,nspden 419 rhodgnow(1:nfft,ispden,ii)=str_scale*rhodgnow(1:nfft,ispden,ii) 420 rho1now(1:nfft,ispden,1+ii)=rho1now(1:nfft,ispden,1+ii)+rhodgnow(1:nfft,ispden,ii) 421 end do 422 end do 423 424 !rho1now(:,:,1) contains the 1st-order density, and rho1now(:,:,2:4) contains the grads of the 1st-order density 425 426 !Apply the XC kernel 427 nspgrad=2; if (nspden==2) nspgrad=5 428 ABI_MALLOC(dnexcdn,(cplex*nfft,nspgrad)) 429 430 !== Non polarized 431 if (nspden==1) then 432 do ir=1,nfft 433 r0(:)=kxc(ir,5:7) ; r1(:)=rho1now(ir,1,2:4) 434 gradrho_gradrho1=dot_product(r0,r1) 435 dnexcdn(ir,1)=kxc(ir,1)*rho1now(ir,1,1) + kxc(ir,3)*gradrho_gradrho1 436 coeff_grho=kxc(ir,3)*rho1now(ir,1,1) + kxc(ir,4)*gradrho_gradrho1 437 ! Grad strain derivative contribution enters the following term with a 438 ! factor of two compared to above terms, so add it again. 439 r1(:)=r1(:)+rhodgnow(ir,1,1:3) 440 ! Reuse the storage in rho1now 441 rho1now(ir,1,2:4)=r1(:)*kxc(ir,2)+r0(:)*coeff_grho 442 end do 443 444 !== Spin-polarized 445 else ! nspden==2 446 do ir=1,nfft 447 do ii=1,3 ! grad of spin-up ans spin_dn GS rho 448 r0_up(ii)=kxc(ir,13+2*ii);r0_dn(ii)=kxc(ir,12+2*ii)-kxc(ir,13+2*ii) 449 end do 450 r0(:)=r0_up(:)+r0_dn(:) ! grad of GS rho 451 r1_up(:)=rho1now(ir,1,2:4) ! grad of spin-up rho1 452 r1_dn(:)=rho1now(ir,2,2:4) ! grad of spin-down rho1 453 r1(:)=r1_up(:)+r1_dn(:) ! grad of GS rho1 454 gradrho_gradrho1_up=dot_product(r0_up,r1_up) 455 gradrho_gradrho1_dn=dot_product(r0_dn,r1_dn) 456 gradrho_gradrho1 =dot_product(r0,r1) 457 458 dnexcdn(ir,1)=kxc(ir, 1)*rho1now(ir,1,1) & 459 & +kxc(ir, 2)*rho1now(ir,2,1) & 460 & +kxc(ir, 6)*gradrho_gradrho1_up & 461 & +kxc(ir,11)*gradrho_gradrho1 462 dnexcdn(ir,2)=kxc(ir, 3)*rho1now(ir,2,1) & 463 & +kxc(ir, 2)*rho1now(ir,1,1) & 464 & +kxc(ir, 7)*gradrho_gradrho1_dn & 465 & +kxc(ir,12)*gradrho_gradrho1 466 coeff_grho_corr=kxc(ir,11)*rho1now(ir,1,1) & 467 & +kxc(ir,12)*rho1now(ir,2,1) & 468 & +kxc(ir,13)*gradrho_gradrho1 469 coeff_grho_up=kxc(ir,6)*rho1now(ir,1,1)+kxc(ir,8)*gradrho_gradrho1_up 470 coeff_grho_dn=kxc(ir,7)*rho1now(ir,2,1)+kxc(ir,9)*gradrho_gradrho1_dn 471 472 ! Grad strain derivative contribution enters the following term with a 473 ! factor of two compared to above terms, so add it again. 474 r1_up(:)=r1_up(:)+rhodgnow(ir,1,1:3) 475 r1_dn(:)=r1_dn(:)+rhodgnow(ir,2,1:3) 476 477 ! Reuse the storage in rho1now 478 rho1now(ir,1,2:4)=(kxc(ir,4)+kxc(ir,10))*r1_up(:) & 479 & +kxc(ir,10) *r1_dn(:) & 480 & +coeff_grho_up *r0_up(:) & 481 & +coeff_grho_corr *r0(:) 482 rho1now(ir,2,2:4)=(kxc(ir,5)+kxc(ir,10))*r1_dn(:) & 483 & +kxc(ir,10) *r1_up(:) & 484 & +coeff_grho_dn *r0_dn(:) & 485 & +coeff_grho_corr *r0(:) 486 end do 487 488 end if ! nspden 489 ABI_FREE(rhodgnow) 490 491 vxc1(:,:)=zero 492 call xcpot(cplex,gprimd,ishift,mgga,mpi_enreg,nfft,ngfft,ngrad,nspden,& 493 & nspgrad,qphon,depsxc=dnexcdn,rhonow=rho1now,vxc=vxc1) 494 495 !if you uncomment the following line, you will have to modify 496 !the original function call to pass in gmet and gsqcut 497 !call filterpot(cplex,gmet,gsqcut,nfft,ngfft,nspden,qphon,vxc1) 498 499 ABI_FREE(dnexcdn) 500 ABI_FREE(rho1now) 501 502 end subroutine dfpt_mkvxcstrgga
ABINIT/m_dfpt_mkvxcstr [ Modules ]
NAME
m_dfpt_mkvxcstr
FUNCTION
COPYRIGHT
Copyright (C) 2001-2024 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 .
SOURCE
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_dfpt_mkvxcstr 23 24 use defs_basis 25 use m_errors 26 use m_abicore 27 28 use defs_abitypes, only : MPI_type 29 use m_time, only : timab 30 use m_symtk, only : matr3inv 31 use m_xctk, only : xcden, xcpot 32 33 implicit none 34 35 private