TABLE OF CONTENTS
- ABINIT/m_extfpmd
- ABINIT/m_extfpmd/compute_e_kinetic
- ABINIT/m_extfpmd/compute_entropy
- ABINIT/m_extfpmd/compute_nelect
- ABINIT/m_extfpmd/compute_shiftfactor
- ABINIT/m_extfpmd/destroy
- ABINIT/m_extfpmd/extfpmd_dos
- ABINIT/m_extfpmd/extfpmd_e_fg
- ABINIT/m_extfpmd/init
- m_extfpmd/extfpmd_type
ABINIT/m_extfpmd [ Modules ]
NAME
m_extfpmd
FUNCTION
This module provides routines to run computations at very high temperature with reduced number of bands. High energy orbitals are represented as pure plane waves.
COPYRIGHT
Copyright (C) 2018-2024 ABINIT group (A. Blanchet) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
TODO
1) Add contribution to conductivity. 2) Smooth the contributions.
SOURCE
22 #if defined HAVE_CONFIG_H 23 #include "config.h" 24 #endif 25 26 #include "abi_common.h" 27 28 module m_extfpmd 29 use defs_basis 30 use defs_abitypes 31 use m_io_tools 32 use m_errors 33 use m_geometry 34 use m_special_funcs 35 use m_specialmsg 36 use m_xmpi 37 use m_energies, only : energies_type 38 use m_gsphere, only : getkpgnorm 39 use m_kg, only : mkkin,kpgio 40 use m_mpinfo, only : ptabs_fourdp,proc_distrb_cycle 41 use m_numeric_tools, only : simpson,simpson_int 42 use m_spacepar, only : meanvalue_g 43 44 implicit none 45 public :: extfpmd_dos,extfpmd_e_fg
ABINIT/m_extfpmd/compute_e_kinetic [ Functions ]
NAME
compute_e_kinetic
FUNCTION
Computes the value of the integral corresponding to the missing kinetic energy contribution of free electrons after band cut, with an order 3/2 incomplete Fermi-Dirac integral.
INPUTS
this=extfpmd_type object concerned fermie=chemical potential (Hartree) tsmear=smearing width (or temperature)
OUTPUT
this=extfpmd_type object concerned
SOURCE
401 subroutine compute_e_kinetic(this,fermie,tsmear) 402 ! Arguments ------------------------------- 403 ! Scalars 404 class(extfpmd_type),intent(inout) :: this 405 real(dp),intent(in) :: fermie,tsmear 406 407 ! Local variables ------------------------- 408 ! Scalars 409 integer :: ifft,ispden 410 real(dp) :: factor,gamma,xcut 411 real(dp) :: e_kinetic_hybrid_tf 412 character(len=500) :: msg 413 ! Arrays 414 real(dp),allocatable :: gamma_hybrid_tf(:,:) 415 real(dp),allocatable :: xcut_hybrid_tf(:,:) 416 417 ! ********************************************************************* 418 419 this%e_kinetic=zero 420 factor=sqrt(2.)/(PI*PI)*this%ucvol*tsmear**(2.5) 421 422 ! Computes extfpmd contribution to kinetic energy integrating 423 ! over accessible states from bcut to infinity with 424 ! order 3/2 incomplete Fermi-Dirac integral. 425 if(this%version==4.or.this%version==2) then 426 gamma=(fermie-this%shiftfactor)/tsmear 427 xcut=extfpmd_e_fg(dble(this%bcut),this%ucvol)/tsmear 428 this%e_kinetic=this%e_kinetic+factor*djp32(xcut,gamma) 429 end if 430 431 ! Computes extfpmd contribution to kinetic energy integrating 432 ! over energy from e_bcut to infinity with order 3/2 433 ! incomplete Fermi-Dirac integral. 434 if(this%version==1.or.this%version==3) then 435 gamma=(fermie-this%shiftfactor)/tsmear 436 xcut=(this%e_bcut-this%shiftfactor)/tsmear 437 this%e_kinetic=this%e_kinetic+factor*djp32(xcut,gamma) 438 end if 439 440 ! Computes extfpmd contribution to kinetic energy using a sum 441 ! of Fermi gas contributions for each point of the fftf grid. 442 ! Warning: This is not yet operational. Work in progress. 443 if(this%version==10) then 444 ABI_MALLOC(gamma_hybrid_tf,(this%nfft,this%nspden)) 445 ABI_MALLOC(xcut_hybrid_tf,(this%nfft,this%nspden)) 446 gamma_hybrid_tf(:,:)=(fermie-this%vtrial(:,:))/tsmear 447 xcut_hybrid_tf(:,:)=(this%e_bcut-this%vtrial(:,:))/tsmear 448 e_kinetic_hybrid_tf=zero 449 450 !$OMP PARALLEL DO REDUCTION (+:e_kinetic_hybrid_tf) 451 do ifft=1,this%nfft 452 do ispden=1,this%nspden 453 e_kinetic_hybrid_tf=e_kinetic_hybrid_tf+factor*djp32(xcut_hybrid_tf(ifft,ispden),gamma_hybrid_tf(ifft,ispden))/& 454 & (this%nfft*this%nspden) 455 end do 456 end do 457 !$OMP END PARALLEL DO 458 459 this%e_kinetic=e_kinetic_hybrid_tf 460 gamma_hybrid_tf(:,:)=zero 461 xcut_hybrid_tf(:,:)=zero 462 ABI_FREE(gamma_hybrid_tf) 463 ABI_FREE(xcut_hybrid_tf) 464 end if 465 466 467 ! Computes the double counting term from the shiftfactor, and 468 ! from the contributions to the kinetic energy and 469 ! the number of electrons 470 if(this%version==1.or.this%version==2.or.this%version==3.or.this%version==4) then 471 this%edc_kinetic=this%e_kinetic+this%nelect*this%shiftfactor 472 else if(this%version==10) then 473 this%edc_kinetic=this%e_kinetic+sum(this%nelectarr(:,:)*this%vtrial(:,:)/(this%nfft*this%nspden)) 474 end if 475 476 if((this%e_bcut.lt.this%shiftfactor).and.& 477 & (this%version==1.or.this%version==3.or.this%version==10)) then 478 write(msg,'(11a)')& 479 & 'Extended FPMD could not properly compute the contribution to the energy.',ch10,& 480 & 'This can be due to a too low number of bands in the calculation.',ch10,& 481 & 'This can also happen when restarting from a previous calculation.',ch10,& 482 & 'Poor prediction of the electron density based on forces may results in this error.',ch10,& 483 & 'Action: slightly increase nband if the electron density is supposed to be converged.',ch10,& 484 & 'Otherwise: wait for the density to be converged.' 485 ABI_WARNING(msg) 486 end if 487 end subroutine compute_e_kinetic
ABINIT/m_extfpmd/compute_entropy [ Functions ]
NAME
compute_entropy
FUNCTION
Computes the value of the integral corresponding to the missing entropy contribution of free electrons after band cut using incomplete Fermi-Dirac integrals.
INPUTS
this=extfpmd_type object concerned fermie=chemical potential (Hartree) tsmear=smearing width (or temperature)
OUTPUT
this=extfpmd_type object concerned
SOURCE
508 subroutine compute_entropy(this,fermie,tsmear) 509 ! Arguments ------------------------------- 510 ! Scalars 511 class(extfpmd_type),intent(inout) :: this 512 real(dp),intent(in) :: fermie,tsmear 513 514 ! Local variables ------------------------- 515 ! Scalars 516 integer :: ii,ifft,ispden 517 real(dp) :: ix,step,factor,fn,gamma 518 ! Arrays 519 real(dp),dimension(:),allocatable :: valuesent 520 real(dp),dimension(:,:),allocatable :: gamma_hybrid_tf 521 real(dp),dimension(:,:),allocatable :: step_hybrid_tf 522 523 ! ********************************************************************* 524 525 this%entropy=zero 526 527 ! Computes extfpmd contribution to the entropy integrating 528 ! over accessible states with Fermi-Dirac complete integrals and 529 ! substracting 0 to bcut contribution with numeric integration. 530 if(this%version==2.or.this%version==4) then 531 factor=sqrt(2.)/(PI*PI)*this%ucvol*tsmear**(2.5) 532 gamma=(fermie-this%shiftfactor)/tsmear 533 ABI_MALLOC(valuesent,(this%bcut+1)) 534 535 !$OMP PARALLEL DO PRIVATE(fn,ix) SHARED(valuesent) 536 do ii=1,this%bcut+1 537 ix=dble(ii)-one 538 fn=fermi_dirac(extfpmd_e_fg(ix,this%ucvol)+this%shiftfactor,fermie,tsmear) 539 if(fn>tol16.and.(one-fn)>tol16) then 540 valuesent(ii)=-two*(fn*log(fn)+(one-fn)*log(one-fn)) 541 else 542 valuesent(ii)=zero 543 end if 544 end do 545 !$OMP END PARALLEL DO 546 547 ! We need at least 6 elements in valuesent to call simpson function. 548 if(size(valuesent)>=6) then 549 this%entropy=5./3.*factor*dip32(gamma)/tsmear-& 550 & gamma*factor*dip12(gamma)/tsmear-& 551 simpson(one,valuesent) 552 end if 553 ABI_FREE(valuesent) 554 end if 555 556 ! Computes extfpmd contribution to the entropy integrating 557 ! over energy with Fermi-Dirac complete integrals and 558 ! substracting 0 to bcut contribution with numeric integration. 559 if(this%version==1.or.this%version==3) then 560 factor=sqrt(2.)/(PI*PI)*this%ucvol*tsmear**(2.5) 561 gamma=(fermie-this%shiftfactor)/tsmear 562 ABI_MALLOC(valuesent,(this%bcut+1)) 563 564 step=(this%e_bcut-this%shiftfactor)/(this%bcut) 565 !$OMP PARALLEL DO PRIVATE(fn,ix) SHARED(valuesent) 566 do ii=1,this%bcut+1 567 ix=this%shiftfactor+(dble(ii)-one)*step 568 fn=fermi_dirac(ix,fermie,tsmear) 569 if(fn>tol16.and.(one-fn)>tol16) then 570 valuesent(ii)=-(fn*log(fn)+(one-fn)*log(one-fn))*& 571 & extfpmd_dos(ix,this%shiftfactor,this%ucvol) 572 else 573 valuesent(ii)=zero 574 end if 575 end do 576 !$OMP END PARALLEL DO 577 578 ! We need at least 6 elements in valuesent to call simpson function. 579 if(size(valuesent)>=6) then 580 this%entropy=5./3.*factor*dip32(gamma)/tsmear-& 581 & gamma*factor*dip12(gamma)/tsmear-simpson(step,valuesent) 582 end if 583 ABI_FREE(valuesent) 584 end if 585 586 ! Computes extfpmd contribution to the entropy using a sum 587 ! of Fermi gas contributions for each point of the fftf grid, 588 ! as we do for version=1 and version=2. 589 ! Warning: This is not yet operational. Work in progress. 590 if(this%version==10) then 591 ABI_MALLOC(valuesent,(this%bcut+1)) 592 ABI_MALLOC(gamma_hybrid_tf,(this%nfft,this%nspden)) 593 ABI_MALLOC(step_hybrid_tf,(this%nfft,this%nspden)) 594 gamma_hybrid_tf(:,:)=(fermie-this%vtrial(:,:))/tsmear 595 step_hybrid_tf(:,:)=(this%e_bcut-this%vtrial(:,:))/(this%bcut) 596 this%entropy=zero 597 factor=sqrt(2.)/(PI*PI)*this%ucvol*tsmear**(2.5) 598 599 do ifft=1,this%nfft 600 do ispden=1,this%nspden 601 !$OMP PARALLEL DO PRIVATE(ix,fn) SHARED(valuesent) 602 do ii=1,this%bcut+1 603 ix=this%vtrial(ifft,ispden)+(dble(ii)-one)*step_hybrid_tf(ifft,ispden) 604 fn=fermi_dirac(ix,fermie,tsmear) 605 if(fn>tol16.and.(one-fn)>tol16) then 606 valuesent(ii)=-(fn*log(fn)+(one-fn)*log(one-fn))*& 607 & extfpmd_dos(ix,this%vtrial(ifft,ispden),this%ucvol) 608 else 609 valuesent(ii)=zero 610 end if 611 end do 612 !$OMP END PARALLEL DO 613 614 ! We need at least 6 elements in valuesent to call simpson function. 615 if(size(valuesent)>=6) then 616 this%entropy=this%entropy+(5./3.*factor*dip32(gamma_hybrid_tf(ifft,ispden))/tsmear-& 617 & gamma_hybrid_tf(ifft,ispden)*factor*dip12(gamma_hybrid_tf(ifft,ispden))/tsmear-& 618 & simpson(step_hybrid_tf(ifft,ispden),valuesent))/(this%nfft*this%nspden) 619 end if 620 end do 621 end do 622 623 gamma_hybrid_tf(:,:)=zero 624 step_hybrid_tf(:,:)=zero 625 ABI_FREE(step_hybrid_tf) 626 ABI_FREE(gamma_hybrid_tf) 627 ABI_FREE(valuesent) 628 end if 629 end subroutine compute_entropy
ABINIT/m_extfpmd/compute_nelect [ Functions ]
NAME
compute_nelect
FUNCTION
Computes the value of the integral corresponding to the missing free electrons contribution after band cut, with an order 1/2 incomplete Fermi-Dirac integral.
INPUTS
this=extfpmd_type object concerned fermie=chemical potential (Hartree) nelect=number of electrons per unit cell tsmear=smearing width (or temperature)
OUTPUT
this=extfpmd_type object concerned
SOURCE
318 subroutine compute_nelect(this,fermie,nelect,tsmear) 319 ! Arguments ------------------------------- 320 ! Scalars 321 real(dp),intent(in) :: fermie,tsmear 322 real(dp),intent(inout) :: nelect 323 class(extfpmd_type),intent(inout) :: this 324 325 ! Local variables ------------------------- 326 ! Scalars 327 integer :: ifft,ispden 328 real(dp) :: factor,gamma,xcut 329 ! Arrays 330 real(dp),allocatable :: gamma_hybrid_tf(:,:) 331 real(dp),allocatable :: xcut_hybrid_tf(:,:) 332 333 ! ********************************************************************* 334 335 factor=sqrt(2.)/(PI*PI)*this%ucvol*tsmear**(1.5) 336 337 ! Computes extfpmd contribution to nelect integrating 338 ! over accessible states from bcut to infinity with 339 ! order 1/2 incomplete Fermi-Dirac integral. 340 if(this%version==4.or.this%version==2) then 341 gamma=(fermie-this%shiftfactor)/tsmear 342 xcut=extfpmd_e_fg(dble(this%bcut),this%ucvol)/tsmear 343 nelect=nelect+factor*djp12(xcut,gamma) 344 end if 345 346 ! Computes extfpmd contribution to nelect integrating 347 ! over energy from e_bcut to infinity with order 1/2 348 ! incomplete Fermi-Dirac integral. 349 if(this%version==1.or.this%version==3) then 350 gamma=(fermie-this%shiftfactor)/tsmear 351 xcut=(this%e_bcut-this%shiftfactor)/tsmear 352 if(this%e_bcut.lt.this%shiftfactor) xcut=zero 353 nelect=nelect+factor*djp12(xcut,gamma) 354 end if 355 356 ! Computes extfpmd contribution to nelect using a sum 357 ! of Fermi gas contributions for each point of the fftf grid. 358 ! Warning: This is not yet operational. Work in progress. 359 if(this%version==10) then 360 ABI_MALLOC(gamma_hybrid_tf,(this%nfft,this%nspden)) 361 ABI_MALLOC(xcut_hybrid_tf,(this%nfft,this%nspden)) 362 gamma_hybrid_tf(:,:)=(fermie-this%vtrial(:,:))/tsmear 363 xcut_hybrid_tf(:,:)=(this%e_bcut-this%vtrial(:,:))/tsmear 364 if(ANY(this%e_bcut.lt.this%vtrial(:,:))) xcut_hybrid_tf(:,:)=zero 365 366 !$OMP PARALLEL DO 367 do ifft=1,this%nfft 368 do ispden=1,this%nspden 369 this%nelectarr(ifft,ispden)=factor*djp12(xcut_hybrid_tf(ifft,ispden),gamma_hybrid_tf(ifft,ispden)) 370 end do 371 end do 372 !$OMP END PARALLEL DO 373 374 nelect=nelect+sum(this%nelectarr)/(this%nfft*this%nspden) 375 gamma_hybrid_tf(:,:)=zero 376 xcut_hybrid_tf(:,:)=zero 377 ABI_FREE(gamma_hybrid_tf) 378 ABI_FREE(xcut_hybrid_tf) 379 end if 380 end subroutine compute_nelect
ABINIT/m_extfpmd/compute_shiftfactor [ Functions ]
NAME
compute_shiftfactor
FUNCTION
Computes the energy shift factor $U_0$ corresponding to constant potential contribution.
INPUTS
this=extfpmd_type object concerned eigen(mband*nkpt*nsppol)=eigenvalues (hartree) eknk(mband*nkpt*nsppol)=kinetic energies (hartree) mband=maximum number of bands nband(nkpt*nsppol)=desired number of bands at each k point nkpt=number of k points nsppol=1 for unpolarized, 2 for spin-polarized wtk(nkpt)=k point weights
OUTPUT
this=extfpmd_type object concerned
SOURCE
193 subroutine compute_shiftfactor(this,eigen,eknk,mband,me,nband,nkpt,nsppol,wtk) 194 ! Arguments ------------------------------- 195 ! Scalars 196 class(extfpmd_type),intent(inout) :: this 197 integer,intent(in) :: mband,me,nkpt,nsppol 198 ! Arrays 199 integer,intent(in) :: nband(nkpt*nsppol) 200 real(dp),intent(in) :: eigen(mband*nkpt*nsppol) 201 real(dp),intent(in) :: eknk(mband*nkpt*nsppol) 202 real(dp),intent(in) :: wtk(nkpt) 203 204 ! Local variables ------------------------- 205 ! Scalars 206 integer :: band_index,ii,ikpt,isppol,nband_k 207 real(dp) :: abs_err,rel_err 208 character(len=500) :: msg 209 210 ! ********************************************************************* 211 212 ! Computes U_0 from the sum of local 213 ! potentials (vtrial), averaging over all space. 214 ! Simplest and most precise way to evaluate U_0. 215 if(this%version==1.or.this%version==4.or.this%version==10) then 216 this%shiftfactor=sum(this%vtrial)/(this%nfft*this%nspden) 217 218 ! Computes the relative error of the model vs last eigenvalues 219 if(me==0.and.this%version==4) then 220 band_index=0 221 rel_err=zero 222 abs_err=zero 223 do isppol=1,nsppol 224 do ikpt=1,nkpt 225 nband_k=nband(ikpt+(isppol-1)*nkpt) 226 rel_err=rel_err+wtk(ikpt)*abs((eigen(band_index+nband_k-this%nbdbuf)-& 227 & extfpmd_e_fg(dble(nband_k-this%nbdbuf),this%ucvol)-this%shiftfactor)/& 228 & eigen(band_index+nband_k-this%nbdbuf)) 229 abs_err=abs_err+wtk(ikpt)*abs(eigen(band_index+nband_k-this%nbdbuf)-& 230 & extfpmd_e_fg(dble(nband_k-this%nbdbuf),this%ucvol)-this%shiftfactor) 231 band_index=band_index+nband_k 232 end do 233 end do 234 if(rel_err.gt.tol1) then 235 write(msg,'(a,es8.2,3a,es8.2,a,es8.2,7a)')& 236 & 'Relative difference between eigenvalues and Fermi gas energy (',rel_err,')',ch10,& 237 & 'is over ',tol1,' at band cut. Absolute difference is ',abs_err,' Ha.',ch10,& 238 & 'Execution will continue as the code will still add contributions in the right',ch10,& 239 & 'direction, but you should likely increase nband to make sure electrons of last',ch10,& 240 & 'band can be considered as free fermions in a constant potential.' 241 ABI_WARNING(msg) 242 end if 243 end if 244 end if 245 246 ! Computes U_0^{HEG} from the difference between 247 ! eigenvalues and Fermi gas energies, averaged 248 ! over lasts nbcut bands. 249 if(this%version==2) then 250 this%shiftfactor=zero 251 band_index=0 252 do isppol=1,nsppol 253 do ikpt=1,nkpt 254 nband_k=nband(ikpt+(isppol-1)*nkpt) 255 do ii=nband_k-this%nbdbuf-this%nbcut+1,nband_k-this%nbdbuf 256 this%shiftfactor=this%shiftfactor+& 257 & wtk(ikpt)*(eigen(band_index+ii)-extfpmd_e_fg(dble(ii),this%ucvol)) 258 end do 259 band_index=band_index+nband_k 260 end do 261 end do 262 this%shiftfactor=this%shiftfactor/this%nbcut 263 end if 264 265 ! Computes U_0^K from the difference between 266 ! eigenvalues and kinetic energies, averaged 267 ! over lasts nbcut bands. 268 if(this%version==3) then 269 this%shiftfactor=zero 270 band_index=0 271 do isppol=1,nsppol 272 do ikpt=1,nkpt 273 nband_k=nband(ikpt+(isppol-1)*nkpt) 274 do ii=nband_k-this%nbdbuf-this%nbcut+1,nband_k-this%nbdbuf 275 this%shiftfactor=this%shiftfactor+& 276 & wtk(ikpt)*(eigen(band_index+ii)-eknk(band_index+ii)) 277 end do 278 band_index=band_index+nband_k 279 end do 280 end do 281 this%shiftfactor=this%shiftfactor/this%nbcut 282 end if 283 284 ! Get extended FPMD band energy cutoff for version 1, 3 and 10. 285 if(this%version==1.or.this%version==3.or.this%version==10) then 286 this%e_bcut=0 287 band_index=0 288 do isppol=1,nsppol 289 do ikpt=1,nkpt 290 nband_k=nband(ikpt+(isppol-1)*nkpt) 291 this%e_bcut=this%e_bcut+wtk(ikpt)*eigen(band_index+nband_k-this%nbdbuf) 292 band_index=band_index+nband_k 293 end do 294 end do 295 end if 296 end subroutine compute_shiftfactor
ABINIT/m_extfpmd/destroy [ Functions ]
NAME
destroy
FUNCTION
Destroy extfpmd_type object, memory deallocation of arrays...
INPUTS
this=extfpmd_type object concerned
OUTPUT
this=extfpmd_type object concerned
SOURCE
143 subroutine destroy(this) 144 145 ! Arguments ------------------------------- 146 ! Scalars 147 class(extfpmd_type),intent(inout) :: this 148 149 ! ********************************************************************* 150 151 this%vtrial(:,:)=zero 152 ABI_FREE(this%vtrial) 153 this%nelectarr(:,:)=zero 154 ABI_FREE(this%nelectarr) 155 this%nfft=0 156 this%nspden=0 157 this%bcut=0 158 this%nbcut=0 159 this%nbdbuf=0 160 this%version=1 161 this%e_bcut=zero 162 this%edc_kinetic=zero 163 this%e_kinetic=zero 164 this%entropy=zero 165 this%nelect=zero 166 this%shiftfactor=zero 167 this%ucvol=zero 168 end subroutine destroy
ABINIT/m_extfpmd/extfpmd_dos [ Functions ]
NAME
extfpmd_dos
FUNCTION
Returns the free particle density of states for a given energy.
INPUTS
energy=get the value of the free particle density of states at this energy shiftfactor=energy shift factor ucvol=unit cell volume (bohr^3)
OUTPUT
extfpmd_dos=value of free particle density of states at given energy
SOURCE
650 function extfpmd_dos(energy,shiftfactor,ucvol) 651 ! Arguments ------------------------------- 652 ! Scalars 653 real(dp),intent(in) :: energy,shiftfactor,ucvol 654 real(dp) :: extfpmd_dos 655 656 ! ********************************************************************* 657 658 extfpmd_dos=sqrt(2.)*ucvol*sqrt(energy-shiftfactor)/(PI*PI) 659 end function extfpmd_dos
ABINIT/m_extfpmd/extfpmd_e_fg [ Functions ]
NAME
extfpmd_e_fg
FUNCTION
Returns the energy of the Fermi gas for a given number of accessible states.
INPUTS
iband=number of accessible states ucvol=unit cell volume (bohr^3)
OUTPUT
extfpmd_e_fg=energy of homogeneous electron gas for a given number of accessible states
SOURCE
678 function extfpmd_e_fg(iband,ucvol) 679 ! Arguments ------------------------------- 680 ! Scalars 681 real(dp),intent(in) :: iband,ucvol 682 real(dp) :: extfpmd_e_fg 683 684 ! ********************************************************************* 685 686 extfpmd_e_fg=.5*(iband*6*PI*PI/ucvol)**(2./3.) 687 end function extfpmd_e_fg
ABINIT/m_extfpmd/init [ Functions ]
NAME
init
FUNCTION
Initialize extfpmd_type object, memory allocation of arrays...
INPUTS
this=extfpmd_type object concerned mband=maximum number of bands nbcut=number of states used to average the constant potential value nbdbuf=Number of bands in the buffer to converge scf cycle with extfpmd models rprimd(3,3)=dimensional primitive translations in real space (bohr) version=extfpmd implementation version
OUTPUT
this=extfpmd_type object concerned
SOURCE
95 subroutine init(this,mband,nbcut,nbdbuf,nfft,nspden,rprimd,version) 96 ! Arguments ------------------------------- 97 ! Scalars 98 class(extfpmd_type),intent(inout) :: this 99 integer,intent(in) :: mband,nbcut,nbdbuf,nfft,nspden,version 100 ! Arrays 101 real(dp),intent(in) :: rprimd(3,3) 102 103 ! Local variables ------------------------- 104 ! Arrays 105 real(dp) :: gprimd(3,3),rmet(3,3), gmet(3,3) 106 107 ! ********************************************************************* 108 109 this%bcut=mband-nbdbuf 110 this%nbcut=nbcut 111 this%nbdbuf=nbdbuf 112 this%version=version 113 ABI_MALLOC(this%vtrial,(nfft,nspden)) 114 this%vtrial(:,:)=zero 115 this%nfft=nfft 116 this%nspden=nspden 117 this%e_bcut=zero 118 this%edc_kinetic=zero 119 this%e_kinetic=zero 120 this%entropy=zero 121 this%nelect=zero 122 ABI_MALLOC(this%nelectarr,(nfft,nspden)) 123 this%nelectarr(:,:)=zero 124 this%shiftfactor=zero 125 call metric(gmet,gprimd,-1,rmet,rprimd,this%ucvol) 126 end subroutine init
m_extfpmd/extfpmd_type [ Types ]
[ Top ] [ m_extfpmd ] [ Types ]
NAME
extfpmd_type
FUNCTION
Store extfpmd functions and parameters.
SOURCE
58 type,public :: extfpmd_type 59 integer :: bcut,nbcut,nbdbuf,nfft,nspden,version 60 real(dp) :: e_bcut,edc_kinetic,e_kinetic,entropy 61 real(dp) :: nelect,shiftfactor,ucvol 62 real(dp),allocatable :: vtrial(:,:) 63 real(dp),allocatable :: nelectarr(:,:) 64 contains 65 procedure :: compute_e_kinetic 66 procedure :: compute_entropy 67 procedure :: compute_nelect 68 procedure :: compute_shiftfactor 69 procedure :: init 70 procedure :: destroy 71 end type extfpmd_type