! nsx = max number of different species ! nax = max number of atoms in one particular species module ion_parameters integer, parameter:: nsx=5, nax=180 end module ion_parameters module van_parameters ! nlx = combined angular momentum (for s,p,d states: nlx=9) ! lix = max angular momentum l+1 (lix=3 if s,p,d are included) ! lx = max 2*l+1 ! mx = 2*lx-1 integer, parameter:: lix=3, nlx=9, lx=2*lix-1, mx=2*lx-1 end module van_parameters module bhs ! analytical BHS pseudopotential parameters use ion_parameters real(kind=8) rc1(nsx), rc2(nsx), wrc1(nsx), wrc2(nsx), & rcl(3,nsx,3), al(3,nsx,3), bl(3,nsx,3) integer lloc(nsx) end module bhs module control ! iprsta = output verbosity (increasing from 0 to infinity) integer iprsta end module control module core ! rhoc = core charge (not yet used) real(kind=8), allocatable:: rhoc(:,:) end module core module cnst ! scmass = 1822.89d0 = mass of a proton, in a.u. real(kind=8), parameter:: pi=3.14159265358979d0, fpi=4.d0*pi real(kind=8) scmass end module cnst module cvan ! ionic pseudo-potential variables use ion_parameters use van_parameters ! ap = Clebsch-Gordan coefficients (?) ! lpx = max number of allowed Y_lm ! lpl = composite lm index of Y_lm real(kind=8) ap(25,nlx,nlx) integer lpx(nlx,nlx),lpl(nlx,nlx,mx) ! nvb = number of species with Vanderbilt PPs ! nh(is) = number of beta functions, including Y_lm, for species is ! ish(is)= used for indexing the nonlocal projectors betae ! with contiguous indices inl=ish(is)+(iv-1)*na(is)+1 ! where "is" is the species and iv=1,nh(is) ! nhx = max value of nh(np) ! nhsavb = total number of Vanderbilt nonlocal projectors ! nhsa = total number of nonlocal projectors for all atoms integer nvb, nhsavb, ish(nsx), nh(nsx), nhsa, nhx ! nhtol: nhtol(ind,is)=value of l for projector ind of species is ! indv : indv(ind,is) =beta function (without Y_lm) for projector ind ! indlm: indlm(ind,is)=Y_lm for projector ind integer, allocatable:: nhtol(:,:), indv(:,:), indlm(:,:) ! beta = nonlocal projectors in g space without e^(-ig.r) factor ! qq = ionic Q_ij for each species (Vanderbilt only) ! dvan = ionic D_ij for each species (Vanderbilt only) real(kind=8), allocatable:: beta(:,:,:), qq(:,:,:), dvan(:,:,:) end module cvan module dft_mod integer lda, blyp, becke, bp88, pw91, pbe parameter (lda=0, blyp=1, becke=2, bp88=3, pw91=4, pbe=5) integer dft end module dft_mod module elct ! f = occupation numbers ! qbac = background neutralizing charge real(kind=8), allocatable:: f(:) real(kind=8) qbac ! nel = number of electrons ! nspin = number of spins (1=no spin, 2=LSDA) ! nupdwn= number of states with spin up (1) and down (2) ! iupdwn= first state with spin (1) and down (2) ! n = total number of electronic states ! nx = if n is even, nx=n ; if it is odd, nx=n+1 ! nx is used only to dimension arrays ! ngw = number of plane waves for wavefunctions ! ngwl = number of G-vector shells up to ngw ! ng0 = first G-vector with nonzero modulus ! needed in the parallel case (G=0 is on one node only!) integer nel, nspin, nupdwn(2), iupdwn(2), n, nx integer ngw, ngwl, ng0 ! ispin = spin of each state integer, allocatable:: ispin(:) end module elct module ener real(kind=8) etot,ekin,eht,epseu,enl,exc,esr,eself end module ener module gvec ! G-vector quantities for the thick grid - see also doc in ggen ! g = G^2 in increasing order (in units of tpiba2=(2pi/a)^2) ! gl = shells of G^2 ( " " " " " ) ! gx = G-vectors ( " " " tpiba =(2pi/a) ) real(kind=8), allocatable:: gl(:), g(:), gx(:,:) ! tpiba = 2*pi/alat ! tpiba2 = (2*pi/alat)**2 real(kind=8) tpiba, tpiba2 ! np = fft index for G> ! nm = fft index for G< ! in1p,in2p,in3p = G components in crystal axis integer,allocatable:: np(:), nm(:), in1p(:),in2p(:),in3p(:), igl(:) ! ng = number of G vectors for density and potential ! ngl = number of shells of G integer ng, ngl end module gvec module gvecb ! As above, for the box grid real(kind=8), allocatable:: gb(:), gxb(:,:),gxnb(:,:), glb(:) integer, allocatable:: npb(:),nmb(:),iglb(:),in1pb(:),in2pb(:),in3pb(:) integer ngb, nglb end module gvecb module gvecs ! As above, for the smooth grid integer ngs, ngsl integer, allocatable:: nps(:), nms(:) end module gvecs module ions use ion_parameters ! nsp = number of species ! na(is) = number of atoms of species is ! nas = max number of atoms of a given species ! nat = total number of atoms of all species ! ipp(is) = PP type for species is (see INPUT) integer nat, nas, nsp, na(nsx), ipp(nsx) ! zv(is) = (pseudo-)atomic charge ! pmass(is) = mass (converted to a.u.) of ions ! rcmax(is) = Ewald radius (for ion-ion interactions) real(kind=8) zv(nsx), pmass(nsx), rcmax(nsx) end module ions module ncprm use ion_parameters use van_parameters ! ! lqx : maximum angular momentum of Q ! nqfx : maximum number of coefficients in Q smoothing ! nvalx: maximum number of pseudo wavefunctions ! nbrx : maximum number of distinct radial beta functions ! mmaxx: maximum number of points in the radial grid ! integer nqfx, nvalx, lqx, nbrx, mmaxx parameter (lqx=5, nqfx=8, nvalx=6, nbrx=7, mmaxx=921) ! ifpcor 1 if "partial core correction" of louie, froyen, ! & cohen to be used; 0 otherwise ! nbeta number of beta functions (sum over all l) ! kkbeta last radial mesh point used to describe functions ! which vanish outside core ! nqf coefficients in Q smoothing ! nqlc angular momenta present in Q smoothing ! lll lll(j) is l quantum number of j'th beta function integer ifpcor(nsx), nbeta(nsx), kkbeta(nsx), & nqf(nsx), nqlc(nsx), lll(nbrx,nsx) ! rscore partial core charge (if present) ! dion bare pseudopotential D_{\mu,\nu} parameters ! (ionic and screening parts subtracted out) ! betar the beta function on a r grid (actually, r*beta) ! qqq Q_ij matrix ! qfunc Q_ij(r) function (for r>rinner) ! rinner radius at which to cut off partial core or Q_ij ! ! qfcoef coefficients to pseudize qfunc for different total ! angular momentum (for r