1!========================================================================== 2! 3! Modules: 4! 5! (1) typedefs Originally By GMR Last Modified 7/8/2008 (JRD) 6! 7!> Derived types that are used throughout the code. 8! 9!========================================================================== 10 11#include "f_defs.h" 12 13module typedefs_m 14 15 use nrtype_m 16 17 implicit none 18 19 public ! only types in this module 20 21!--------------------------- 22 23 type crystal 24 real(DP) :: celvol !< cell volume in real space (a.u.) 25 real(DP) :: recvol !< cell volume in reciprocal space (a.u.) 26 real(DP) :: alat !< lattice constant in real space (a.u.) 27 real(DP) :: blat !< lattice constant in reciprocal space (a.u.) 28 real(DP) :: avec(3,3) !< lattice vectors in real space (alat) 29 real(DP) :: bvec(3,3) !< lattice vectors in reciprocal space (blat) 30 real(DP) :: adot(3,3) !< metric tensor in real space (a.u.) 31 real(DP) :: bdot(3,3) !< metric tensor in reciprocal space (a.u.) 32 integer :: nat !< number of atoms 33 integer, pointer :: atyp(:) !< atomic species, atyp(1:nat) 34 real(DP), pointer :: apos(:,:) !< atomic positions, apos(1:3,1:nat) (alat) 35 end type crystal 36 37!--------------------------- 38 39 type kpoints 40 integer :: nspinor = 1 !< nspinor = 2 if doing two-component spinor calculation; 1 is default 41 integer :: nspin !< nspin = 1 or 2; nspin = 1 when npsinor = 2 42 integer :: nrk !< number of k-points 43 integer :: mnband !< max number of bands 44 integer :: nvband !< number of valence bands 45 integer :: ncband !< number of conduction bands 46 integer :: kgrid(3) !< Monkhorst-Pack number of k-points in each direction 47 real(DP) :: shift(3) !< Monkhorst-Pack shift of grid 48 real(DP) :: ecutwfc !< wave-function cutoff, in Ry 49 integer, pointer :: ngk(:) !< number of g-vectors for each k-point 50 integer :: ngkmax !< max(ngk(:)) 51 integer, pointer :: ifmin(:,:) !< lowest occupied band (kpoint,spin) 52 integer, pointer :: ifmax(:,:) !< highest occupied band (kpoint,spin) 53 real(DP), pointer :: w(:) !< weights (kpoint) (between 0 and 1) 54 real(DP), pointer :: rk(:,:) !< k-vector (3, kpoint) in crystal coords 55 real(DP), pointer :: el(:,:,:) !< band energies (band, kpoint, spin) 56 real(DP), pointer :: elda(:,:,:) !< band energies before eqp correction 57 real(DP), pointer :: occ(:,:,:) !< occupations (between 0 and 1) 58 integer, pointer :: degeneracy(:,:,:) !< size of deg. subspace for (band, kpoint, spin) 59 end type kpoints 60 61!--------------------------- 62 63 type symmetry 64 integer :: ntran !< number of operations in full group 65 integer :: ntranq !< number of operations in small group of q 66 real(DP) :: rq(3) !< The q-point this ntranq belongs to 67 integer :: mtrx(3,3,48) !< symmetry matrix 68 real(DP) :: tnp(3,48) !< fractional translations 69 integer :: indsub(48) !< symmetry operations in subgroup of q 70 integer :: kgzero(3,48) !< Umklapp vectors for subgroup symmetry operations 71 integer :: cell_symmetry !< 0 = cubic, 1 = hexagonal 72 end type symmetry 73 74!--------------------------- 75 76 type grid 77 integer :: nr !< number in reduced zone 78 integer :: nf !< number in full zone 79 real(DP) :: sz !< radius of a spherical subzone equivalent to 80 !! one point in the set f 81 integer, pointer :: itran(:) !< sym op to go from irrbz to fullbz 82 integer, pointer :: indr(:) !< Index of the irrbz k-point to which each 83 !! k-point of the fullbz relates 84 integer, pointer :: kg0(:,:) !< Umklapp vectors (for Wigner-Seitz cell) 85 real(DP), pointer :: r(:,:) !< k/q-points in reduced zone 86 real(DP), pointer :: f(:,:) !< k/q-points in full zone 87 end type grid 88 89!----------------------------------- 90 91 type gspace 92 integer :: ng !< number of G-vectors 93 integer :: nFFTgridpts !< number in FFT grid = product(FFTgrid(1:3)) 94 real(DP) :: ecutrho !< charge-density cutoff, in Ry 95 integer, pointer :: components(:,:) !< the G-vectors, in units of 2 pi / a 96 integer :: FFTgrid(3) !< gsm: FFTgrid is the size of the FFT grid, not the maximum G-vector 97 integer, pointer :: index_vec(:) ! mapping to FFT grid 98 real(DP), pointer :: ekin(:) !< kinetic energy of G-vectors 99 end type gspace 100 101!--------------------------- 102!> Parameters for scissors operators 103!> e_cor = e_in + es + edel * (e_in - e0) 104!> es and e0 are in eV. edel is a dimensionless slope. 105 type sub_scissors_t 106 real(DP) :: es 107 real(DP) :: edel 108 real(DP) :: e0 109 end type sub_scissors_t 110 111 type scissors_t 112 type(sub_scissors_t) :: val 113 type(sub_scissors_t) :: cond 114 end type scissors_t 115 116!--------------------------- 117 118 type wavefunction 119 integer :: ng !< Number of planewaves 120 integer :: nband !< Number of bands 121 integer :: nspin !< Number of spins 122 integer :: nspinor = 1 !< nspinor = 2 if doing two-component spinor calculation; 1 is default 123 integer, pointer :: isort(:) !< (ng) 124 SCALAR, pointer :: cg(:,:,:) !< Planewave coefficients (ng,nband,ns*nspinor) 125 end type wavefunction 126 127!--------------------------- 128 129!> For Epsilon: this is the wavefunction before unfolding the irr. BZ. 130!! For BSE: ?? 131 type int_wavefunction 132 integer :: nspin 133 integer :: nspinor = 1 !< nspinor = 2 if doing two-component spinor calculation; 1 is default 134 integer, pointer :: ng(:) !< (nk) 135 integer, pointer :: isort(:,:) !< (ngmax, nk) 136 integer, pointer :: cbi(:) 137 !> I think this can be decommissioned if we use kp%rk instead 138 real(DP), pointer :: qk(:,:) 139 SCALAR, pointer :: cg(:,:,:) 140 SCALAR, pointer :: cgk(:,:,:,:) !< (ngkmax, nband, nspin*nspinor,ikt) 141 !! where ikt is the number of k-points owned locally 142 end type int_wavefunction 143 144!------------------------------------ 145 146 !> FHJ: valence WFNs for 1 particular kpt and 1 band 147 !! It stores all bands in the case of real-space WFNs 148 type valence_wfns 149 integer :: nband !< This is indeed the number of valence bands 150 integer :: ncore_excl 151 !>This is the number of core states that are not included in the valence in polarizability calculation. 152 integer :: ngv !< Number of G-vectors 153 integer :: idx_kp !< Idx of current kpt in kp/kpq structure 154 integer, pointer :: isort(:) 155 !> (nband+ncrit,spin). Note: the nband+ncrit index is actually useless! 156 real(DP), pointer :: ev(:,:) 157 SCALAR, pointer :: zv(:,:) !< (ngv,spin) 158 !> real-space wavefunction for all "local" val. bands (fft1,fft2,fft3,band) 159 complex(DPC), pointer :: wfn_fft(:,:,:,:) 160 end type valence_wfns 161 162!------------------------------- 163 164 !> FHJ: conduction WFNs for 1 particular kpt and all bands (!) the processor owns 165 type conduction_wfns 166 integer :: nband !< This is actually the number of valence+conduction bands! 167 integer :: ngc !< Number of G-vectors 168 integer :: idx_kp !< Idx of current kpt in kp structure 169 integer, pointer :: isort(:) 170 real(DP), pointer :: ec(:,:) !< (nband,nspin) 171 SCALAR, pointer :: zc(:,:) !< (ngc*ncownactual,spin) 172 !> real-space wavefunction for all "local" cond. bands (fft1,fft2,fft3,band) 173 complex(DPC), pointer :: wfn_fft(:,:,:,:) 174 end type conduction_wfns 175 176!---------------------------- 177 178 !> splines knots and coefficients 179 type spline_tck 180 integer :: n !< number of knots 181 integer :: k !< degree of spline (1=linear, etc.) 182 real(DP), pointer :: t(:) !< position of knots 183 real(DP), pointer :: c(:) !< splines coefficient 184 end type spline_tck 185 186!------------------------------------- 187 188 type coulomb_modifier_t 189 190 real(DP) :: short_range_frac_fock !< Short range exchange fraction 191 real(DP) :: long_range_frac_fock !< Long range exchange fraction 192 real(DP) :: screening_length !< Screening length 193 !< The above 3 parameters are used 194 !< only for TDDFT and sigma calculations. 195 196 end type coulomb_modifier_t 197 198!---------------------------- 199 200 type eps_sub_info 201 !> arrays necessary for applying subspace in sigma (only FF-CD) 202 integer :: actual_nm 203 integer :: nb_sub ! will not use the block-cyclic distr, defined here just to be sure will be 1 always 204 integer :: neps 205 integer :: Nbas_own, Nbas_own_max 206 integer :: ngpown_sub, ngpown_sub_max 207 integer, pointer :: eps_sub_info(:,:,:) ! (3,2,1:(nq+1)) indx_start/end/size for each q point (Neig / nmtx) 208 integer, pointer :: wing_pos(:) ! (1:(nq+1)) position of the wings for each q point 209 complex(DPC), pointer :: eigenvec_sub(:,:,:) ! (1:gvec%ng,1:Nbas_own,1:(nq+1)) 210 complex(DPC), pointer :: eps_sub(:,:,:,:) ! (1:Nbas,1:Nbas_own,1:sig%nFreq,1:(nq+1)) 211 complex(DPC), pointer :: eps_wings_rows(:,:,:) ! (1:gvec%ng,1:sig%nFreq,1:(nq+1)) 212 complex(DPC), pointer :: eps_wings_cols(:,:,:) ! (1:gvec%ng,1:sig%nFreq,1:(nq+1)) 213 complex(DPC), pointer :: eps_wings_correction_rows(:,:) ! (1:gvec%ng,1:sig%nFreq) 214 complex(DPC), pointer :: eps_wings_correction_cols(:,:) ! (1:gvec%ng,1:sig%nFreq) 215 real(DP), pointer :: vcoul_sub(:,:) 216 end type eps_sub_info 217 218!---------------------------- 219 220 type siginfo 221 integer :: freq_dep !< frequency dependence of the inverse dielectric matrix 222 integer :: freq_dep_method !< frequency dependence method of the inverse dielectric matrix 223 integer :: nFreq !< number of frequencies used in full frequency calculation 224 real(DP) :: dDeltaFreq !< frequency increment (eV) for polarizability energy denominator 225 real(DP) :: dBrdning !< Lorentzian broadening (eV) for polarizability energy denominator 226 real(DP), pointer :: dFreqGrid(:) !< Grid of Frequencies for Full Frequency 227 complex(DPC), pointer :: dFreqBrd(:) !< Corresponding Broadenings for Full Frequency 228 229 integer :: nSFreq !< number of frequencies used in spectral functions 230 real(DP) :: dDeltaSFreq !< frequency increment (eV) for spectral functions 231 real(DP), pointer :: dSFreqGrid(:) !< Grid of Frequencies for spectral functions 232 233 integer :: exact_ch !< compute the exact static CH 234 integer :: fullConvLog !< logging CH convergence 235 integer :: iwritecoul !< flag to write vcoul 236 real(DP) :: tol !< energy tolerance for degeneracy 237 logical :: use_hdf5 !< with -DHDF5, whether or not we actually use hdf5 238 logical :: wfn_hdf5 !< with -DHDF5, use HDF5 for WFN files? 239 logical :: use_xdat !< use saved exchange matrix elements from file x.dat 240 logical :: use_vxcdat !< use saved exchange-correlation matrix elements from file vxc.dat 241 logical :: use_vxc2dat !< use saved exchange-correlation matrix elements from file vxc2.dat 242 logical :: is_EXX !< is this calculation only of bare exchange for EXX calculation? 243 integer :: nkn !< number of k-points on which to calculate Sigma (from sigma.inp) 244 integer :: nphonq !< ZL: number of phonon_q points. Current implementation supports only nphonq = 1 245 integer :: nfreq_imag 246 logical :: need_advanced 247 integer :: cd_int_method !< Integration method for CD calculations 248 integer :: invalid_gpp_mode !< How to treat invalid GPP modes? See sigma.inp 249 integer :: nq0, nq1, nq !< Number of q->0 points, q/=0 points, and total number of q-points 250 logical :: subsample !< whether we perform a subsampling of the voronoi cell of q=0 251 real(DP), pointer :: subweights(:) !< (nq0) weight for each subsampling q-point 252 253 ! ZL: add for electron phonon calculation 254 logical :: elph !< whether we do electron phonon calculation for 255 ! electron-phonon matrix elements under GW approximation 256 integer :: ep_bracket ! evaluate elph dSig at <bra| or |ket> band energy 257 258 integer :: nvband !< number of bands in bare exchange 259 integer :: ntband !< number of bands in dynamical sigma summation 260 integer :: ncore_excl 261 !> number of core states that are not included in both the bare exchange and dynamical sigma summations 262 integer :: igamma !< nonzero if Gamma is the only q-point, 0 otherwise 263 integer :: nspin 264 integer :: spin_index(2) 265 integer :: icutv !< icutv encodes presence and type of truncation 266 integer :: iuseremainder 267 integer :: qgrid(3) 268 integer :: iscreen !< what type of screening is present. 0 = semiconductor, 1 = graphene, 2 = metal 269 integer :: fdf !< finite difference form for numerical derivative of Sigma 270 integer, pointer :: indkn(:) !< mapping of k-points from sigma.inp to those in kp%rk from WFN files 271 integer, pointer :: indk_phonq(:) !< ZL: mapping of k+phonq points from sigma.inp to those in kp%rk 272 integer, pointer :: indk_phonq_g0(:,:) !< ZL: mapping of k+phonq points from sigma.inp to those in kp%rk, needs umklapp 273 integer, pointer :: diag(:) !< energy bands for which Sigma diagonal matrix elements are calculated 274 integer :: ndiag !< number of bands contained in diag(:) 275 integer :: noffdiag !< offdiag 276 integer :: loff 277 integer :: toff 278 integer :: bmin 279 integer :: bmax 280 integer, pointer :: off1(:) !< offdiag <bra| 281 integer, pointer :: off2(:) !< offdiag |ket> 282 integer, pointer :: off3(:) !< offdiag energy at which to evaluate 283 integer, pointer :: offmap(:,:) !< sig%off1(ii) = sig%diag(sig%offmap(ii,1)) 284 real(DP) :: dw !< finite difference spacing for numerical derivative of Sigma in eV 285 real(DP) :: ecutb !< energy cutoff of bare coulomb interaction in Ry 286 real(DP) :: ecuts !< energy cutoff of screened coulomb interaction in Ry 287 real(DP) :: xfrac !< fraction of bare exchange 288 real(DP) :: gamma !< GPP broadening 289 real(DP) :: sexcut !< GPP SX cutoff 290 real(DP) :: q0vec(3) 291 real(DP) :: phonq(3) !< ZL: phonon_q point. Now only support one point 292 integer :: freq_grid_shift !< How to shift the requency grid. See sigma.inp for more info. 293 integer :: nfreqeval 294 real(DP) :: freqevalmin 295 real(DP) :: freqevalstep 296 logical :: eqp_corrections !< are we using eqp.dat 297 logical :: eqp_outer_corrections !< are we using eqp_outer.dat 298 type(scissors_t) :: scis 299 type(scissors_t) :: scis_outer 300 logical :: wfn_outer_present 301 type(spline_tck) :: spl_tck !< scissors b-spline coefficients 302 type(spline_tck) :: spl_tck_outer !< scissors b-spline coeffs for WFN_outer 303 real(DP) :: avgpot 304 real(DP) :: avgpot_outer 305 real(DP) :: truncval(3) !< in Bohr (au) 306 real(DP) :: avgcut !< Cut in which we do cell averages on 307 real(DP), pointer :: kpt(:,:) 308 309 real(DP), pointer :: k_phonq(:,:) !< ZL: for LR, (3,nkn*nphonq), we consider nphonq=1 only for now, 310 ! so k_phonq will have the same size as kpt, it stores (kpoint(3,nkn) + phonq(3)) 311 312 real(DP), pointer :: qpt(:,:) !< (3,nq) q-points in eps0mat+epsmat files, or provided in sigma.inp 313 integer :: ncrit !< number of partially occupied bands 314 real(DP) :: efermi !< Fermi level 315 real(DP) :: efermi_input !< The value to set E_Fermi in the input file, in eV 316 logical :: rfermi !< Measure the new Fermi level relative to that of the neutral system 317 SCALAR, pointer :: vxc(:,:) !< Vxc(G) 318 SCALAR, pointer :: vxc2(:,:) !< Vxc(G) for hybrid functional type calculations 319 SCALAR, pointer :: dvxc(:,:) !< dVxc(G), ZL adds for EP 320 SCALAR :: wcoul0 321 logical :: freplacebz 322 logical :: fwritebz 323 logical :: degeneracy_check_override 324 logical :: offdiagsym 325 logical :: qgridsym 326 logical :: symmetrize 327 logical :: die_outside_sphere 328 logical :: averagew 329 type(coulomb_modifier_t) :: coulomb_mod 330 logical :: coul_mod_flag !< Flag which tells if the coulomb interaction has been 331 !! modified (for hybrid functional like calculation in sigma) 332 logical :: sigma_correction !< if .true., compute only a correction to the QP self energy, ie, 333 !! don`t subtract vxc and don`t compute the bare exchange. 334 ! variables for subspace truncation 335 logical :: subspace_q0, matrix_in_subspace_basis_q0, keep_full_eps_static_q0 336 logical :: subspace, matrix_in_subspace_basis, keep_full_eps_static 337 integer :: neig_sub_max 338 ! activate the computation of Sigma directly in subspace basis 339 logical :: do_sigma_subspace 340 type(eps_sub_info) :: epssub 341 logical :: sub_collective_redistr 342 end type siginfo 343 344!--------------------------- 345 !> Dielectric matrix info using comm_mpi 346 !! In BLACS terms, we distribute the columns of the epsinv matrix in a 347 !! 1D-block-cyclic fashion. Currently, the block size nb is 1 348 !! The following old quantities were removed, as they can be calculated on the fly: 349 !! epsmpi%igp_owner(igp) = INDXG2P(igp, epsmpi%nb, peinf%pool_rank, 0, peinf%npes_pool) 350 !! epsmpi%igp_index(igp) = INDXG2L(igp, epsmpi%nb, peinf%pool_rank, 0, peinf%npes_pool) 351 !! The following array should be removed soon: 352 !! epsmpi%inv_igp_index(igp_loc) = INDXL2G(igp_loc, epsmpi%nb, peinf%pool_rank, 0, peinf%npes_pool) 353 type epsmpiinfo 354 integer :: nb !< block size. Currently set to 1 (round robin) 355 integer :: ngpown !< number of columns I own. Same as numroc(neps, nb, pool_rank, 0, npes_pool) 356 integer :: ngpown_max !< number of columns owned by pool_rank==0. 357 integer, pointer :: isrtq(:,:) !< These 3 arrays have a dimension of (1:gvec%ng,1:(nq+1)) where 358 integer, pointer :: isrtqi(:,:) !! (nq+1) is the total # of q`s including q0. 359 360 integer, pointer :: inv_igp_index(:) 361 integer, pointer :: nmtx(:) !< dimension of eps(q) for each q 362 real(DP), pointer :: qk(:,:) 363 SCALAR, pointer :: eps(:,:,:) !< eps(1:gvec%ng,1:ngpown,1:(nq+1)) 364 365 !> dimension of epsR and epsA (1:gvec%ng,1:ngpown,1:sig%nFreq,1:(nq+1)) 366 complex(DPC), pointer :: epsR(:,:,:,:) 367 complex(DPC), pointer :: epsA(:,:,:,:) 368 end type epsmpiinfo 369 370!--------------------------- 371 372 type wfnkqmpiinfo 373 integer, pointer :: nkptotal(:) ! ngk 374 integer, pointer :: isort(:,:) 375 integer, pointer :: band_index(:,:) 376 real(DP), pointer :: qk(:,:) 377 real(DP), pointer :: el(:,:,:) 378 SCALAR, pointer :: cg(:,:,:,:) 379 end type wfnkqmpiinfo 380 381!--------------------------- 382 383 type wfnkmpiinfo 384 integer, pointer :: nkptotal(:) 385 integer, pointer :: isort(:,:) 386 real(DP), pointer :: qk(:,:) 387 real(DP), pointer :: el(:,:,:) 388 real(DP), pointer :: elda(:,:,:) 389 SCALAR, pointer :: cg(:,:,:) 390 end type wfnkmpiinfo 391 392!--------------------------- 393 394 type wpgen 395 real(DP) :: wpsq(2) !< square of free el plasma freq for each spin 396 real(DP) :: nelec(2) !< number of electrons for each spin per cell 397 SCALAR, pointer :: rho(:,:) !< density, (ig, ispin) 398 end type wpgen 399 400!--------------------------- 401 402 type polarizability 403 integer :: freq_dep !< frequency dependence of the inverse dielectric matrix 404 ! 0: static calculation 2: full frequency 3: two imaginary frequencies 405 integer :: freq_dep_method !< full frequency calculation. 0: Adler-Wiser; 1: Shishkin and Kresse 2006 406 integer :: nFreq !< number of frequencies used in full frequency calculation 407 integer :: nfreq_imag !< number of imaginary freqs for CD (also 1 for GN GPP) 408 real(DP) :: dInitFreq !< initial frequency (eV) for polarizability energy denominator 409 real(DP) :: dDeltaFreq !< frequency increment (eV) for polarizability energy denominator 410 real(DP) :: dBrdning !< Lorentzian broadening (eV) for polarizability energy denominator 411 real(DP), pointer :: dFreqGrid(:) !< Grid of Frequencies for Full Frequency 412 real(DP) :: dFreqStepIncrease 413 real(DP) :: dFreqCutoff1 414 real(DP) :: dFreqCutoff2 415 416 integer :: nSFreq !< number of frequencies used in spectral function 417 real(DP) :: dInitSFreq !< initial frequency (eV) for polarizability spectral function 418 real(DP) :: dDeltaSFreq !< frequency increment (eV) for polarizability spectral function 419 real(DP), pointer :: dSFreqGrid(:) !< Grid of Frequencies for spectral function 420 real(DP) :: dSFreqStepIncrease 421 real(DP) :: dSFreqCutoff1 422 real(DP) :: dSFreqCutoff2 423 424 logical :: has_advanced !< Do we store eps_A or just eps_R? 425 integer :: matrix_type !< 0 to write epsilon^{-1}, 1 for epsilon, 2 for chi0. 426 integer :: nmatrix !< has_advanced+1. Multiply by nspin if matrix_type==2 427 integer :: matrix_flavor !< 2 (=CMPLX), unless we have freq_dep==0 and SCALARSIZE==1. 428 429 type(scissors_t) :: scis 430 logical :: eqp_corrections !< are we using eqp.dat and eqp_q.dat files 431 complex(DPC), pointer :: dFreqBrd(:) !< Corresponding Broadenings for Full Frequency 432 integer :: fullConvLog !< logging pol matrix head & tail convergence 433 integer :: iwritecoul !< flag to write vcoul 434 integer :: nmtx 435 integer, pointer :: nmtx_of_q(:) 436 integer :: qgrid(3) 437 integer, pointer :: qflags(:) 438 integer :: nq0, nq1, nq !< Number of q->0 points, q/=0 points, and total number of q-points 439 logical :: subsample !< whether we have more than one q0 point (used in subsampled calculation) 440 logical :: non_uniform !< do non-uniform sampling using Voronoi decomposition of BZ 441 integer :: gcomm 442 logical :: min_fftgrid !< use the smallest possible fftbox 443 ! FHJ: These flags control some experimental optimizations 444 integer :: os_opt_ffts !< optimizes calculation/reuse of FFTs (real-space WFNs) 445 integer :: nfreq_group !< num. of frequencies to calculate in parallel 446 integer :: nfreq_in_group !< num. of epsilon frequencies held by any processor 447 integer :: os_nsfreq_para !< num. of spectral frequencies held by any processor 448 logical :: os_hdf5 !< use parallel IO? 449 logical :: restart !< are we restarting the calculation? Only ok with HDF5 450 integer :: stop_after_qpt !< pretend the calculation was prematurely killed after this qpt (-1=don`t kill) 451 integer :: protection_window(2) !< the two band indices that define the protection window (see epsilon.inp) 452 integer :: intraband_flag !< 0=regular calculation, 1=only include intraband, 2=only interband 453 real(DP) :: intraband_overlap_min !< a transition is intraband if |<uvk|uvk+q>|^2 is larger than this 454 integer :: num_cond_bands_ignore !< num. of cond bands to ignore. default is 0. 455 logical :: patched_sampling !< Do we have only a patch in the BZ? 456 ! 457 integer :: WFN_FFTgrid(3)!< max. size FFTgrid that holds all WFNs 458 integer :: FFTgrid(3) !< FFT grid to use (RHO or economical one) 459 !! 460 logical :: skip_epsilon 461 logical :: skip_chi 462 logical :: use_hdf5 !< with -DHDF5, whether or not we actually use hdf5 463 logical :: need_WFNq !< will we need the WFNq file? (nq0>0.and.valueq0==1.and.iqexactlyzero==0) 464 integer :: iqexactlyzero !< 1 if the q->0 point is *exactly* zero and will be read from WFN; 0 otherwise 465 integer :: valueq0 !< 1=semiconductor (read from WFNq); 2=metal (read from WFN) 466 integer, pointer :: irow(:) 467 integer, pointer :: isrtx(:) 468 integer, pointer :: isrtxi(:) 469 integer :: icutv !< icutv encodes presence and type of truncation 470 real(DP) :: truncval(3) !< in Bohr (au) 471 real(DP), pointer :: qpt(:,:) 472 !> FHJ: gme = <c,k|e^(-i(q+G).r)|v,k+q>, and the indices are: 473 !! (nmtx, ncownactual, nvownactual, nspin, nrk, nfreq_group) 474 SCALAR, pointer :: gme(:,:,:,:,:,:) 475 SCALAR, pointer :: chi(:,:,:) 476 integer :: ncrit 477 real(DP) :: efermi 478 real(DP) :: efermi_input 479 logical :: rfermi 480 real(DP) :: ecuts !< energy cutoff of screened coulomb interaction in Ry 481 real(DP) :: ecutsExtra 482!> Reference regarding retarded/advanced functions: Catalin`s thesis, Eq. (1.44) 483 complex(DPC), pointer :: chiRDyn(:,:,:,:) !< Retarded polarizability 484 complex(DPC), pointer :: chiTDyn(:,:,:,:) !< Spectral function of polarizability 485 486 real(DP), pointer :: edenDyn(:,:,:,:,:) !< Dynamic energy denominator 487 logical :: freplacebz 488 logical :: fwritebz 489 logical :: degeneracy_check_override 490 real(DP) :: lin_denominator !< energy threshold below which to activate lin_denominator 491 type(cvpair_info), pointer :: lin_edenDyn(:,:,:,:,:) !< energies and 492 real(DP) :: de_min, de_max 493 ! velocities for calculating linearized denominator in dynamic case 494 real(DP) :: imaginary_frequency !< purely imaginary frequency used in Godby-Needs GPP 495 ! variables for subspace truncation method in epsilon 496 logical :: subspace 497 real(DP) :: chi_eigenvalue_cutoff 498 integer :: neig_sub_input 499 logical :: use_elpa 500 logical :: need_full_chi 501 logical :: keep_full_eps_static 502 logical :: matrix_in_subspace_basis 503 integer :: nrow_local_sub, ncol_local_sub, neig_sub 504 complex(DPC), allocatable :: chiRDyn_sym_omega0(:,:) 505 complex(DPC), allocatable :: eigenvect_omega0(:,:) 506 real(DP), allocatable :: eigenval_omega0(:) 507 real(DP), allocatable :: vcoul_sub(:) 508 ! keep trak of the number of eigenvalue for each q point (used in hdf5 I/O) 509 integer, allocatable :: neigen_of_q(:) 510 logical :: do_rpa ! Calculate RPA correlation energy in subspace calculation 511 real(DP), allocatable :: rpa_freq_grid(:) 512 real(DP), allocatable :: E_rpa_qp(:) ! RPA energy contribution from each q-point 513 real(DP), allocatable :: qw_rpa(:) ! q-point weight; needed for RPA calculations 514 ! variables for nonblocking scheme 515 logical :: nonblocking_cyclic 516 logical :: dont_keep_fix_buffers 517 logical :: sub_collective_eigen_redistr 518 end type polarizability 519 520 521!-------------------------------- 522 523 type cvpair_info 524 real(DP) :: vc(2) !< conduction band velocity 525 real(DP) :: vv(2) !< valence band velocity 526 real(DP) :: ec !< conduction band energy 527 real(DP) :: ev !< valence band energy 528 integer :: idx_kp !< kpoint index 529 logical :: vltc !< ev - ec < TOL_Degeneracy 530 end type cvpair_info 531 532!-------------------------------- 533 534 type wfnkstates 535 integer :: nkpt 536 integer :: ndv 537 integer, pointer :: isrtk(:) 538 real(DP) :: k(3) 539 real(DP), pointer :: ek(:,:) 540 real(DP), pointer :: elda(:,:) 541 SCALAR, pointer :: zk(:,:) 542 end type wfnkstates 543 544!--------------------------- 545 546 type wfnkqstates 547 integer :: nkpt 548 integer, pointer :: isrtkq(:) 549 real(DP), pointer :: ekq(:,:) 550 SCALAR, pointer :: zkq(:,:) 551 end type wfnkqstates 552 553!--------------------------------- 554 555!> Used in haydock/diag only (see epsdiag.f90) 556 type epsinfo 557 integer :: nq !< number of q-vectors stored 558 real(DP) :: emax !< maximum length of the stored q-vectors 559 real(DP), pointer :: q(:,:) !< (3, nq) coordinates of q-vectors 560 real(DP), pointer :: eps(:) !< (nq) head of dielectric matrix at each q-vector 561 SCALAR :: epshead !< head of dielectric matrix at q->0 562 real(DP) :: q0vec(3) !< coordinates of the q->0 vector 563 end type epsinfo 564 565!------------------------------------ 566 567!> Used in haydock/diag only 568 type eqpinfo 569 type(scissors_t) :: scis 570 type(spline_tck) :: spl_tck !< scissors spline coefficients 571 real(DP), pointer :: evqp(:,:,:) 572 real(DP), pointer :: ecqp(:,:,:) 573 real(DP), pointer :: evqp_co(:,:,:) 574 real(DP), pointer :: ecqp_co(:,:,:) 575 real(DP), pointer :: evqp_co_q(:,:,:) 576 real(DP), pointer :: ecqp_co_q(:,:,:) 577 real(DP), pointer :: evlda(:,:,:) 578 real(DP), pointer :: eclda(:,:,:) 579 real(DP), pointer :: evlda_co(:,:,:) 580 real(DP), pointer :: eclda_co(:,:,:) 581 real(DP), pointer :: evlda_co_q(:,:,:) 582 real(DP), pointer :: eclda_co_q(:,:,:) 583 real(DP), pointer :: evshift(:,:,:) 584 real(DP), pointer :: ecshift(:,:,:) 585 real(DP), pointer :: evshift_co(:,:,:) 586 real(DP), pointer :: ecshift_co(:,:,:) 587 real(DP), pointer :: evshift_co_q(:,:,:) 588 real(DP), pointer :: ecshift_co_q(:,:,:) 589 590 integer :: band_ordering=0 !< 0 = counting bands from fermi levels, with 591 !! conduction bands going up in energy 592 !! valence bands down in energy. 593 !! i.e., eqpv(1,:,:) is the VBM, 594 !! eqpv(2,:,:) is VMB-1, etc. 595 !! 1 = counting all bands from bottom-up, 596 !! starting at the first valence band 597 !! and going up in energy. 598 end type eqpinfo 599 600!------------------------------------ 601 602!> moments for Haydock 603 type mmtsinfo 604 integer :: nmax 605 integer :: nmaxp 606 real(DP) :: norm 607 real(DP) :: vol 608 real(DP), pointer :: an(:) !< FIXME document this 609 real(DP), pointer :: bn(:) !< FIXME document this 610 end type mmtsinfo 611 612!------------------------------------ 613 614 type xctinfo 615 logical :: is_absorption !< whether we are running the absorption code 616 integer :: algo !< algorithm to use to solve BSE. See Common/nrtype.f90 617 logical :: inteqp !< whether we are interpolating 618 logical :: is_periodic(3) !< which dimensions are periodic 619 integer :: idimensions !< how many total periodic dimensions 620 integer :: nkpt_co !< number of kpts in the coarse grid 621 integer :: nkptq_co !< number of kpts in the q-shifted coarse grid 622 integer :: nvb_co !< number of valence bands in the coarse grid 623 integer :: ncb_co !< number of conduction bands in the coarse grid 624 integer :: n1b_co !< nvb_co for TDA calculations, nvb_co + ncb_co for non-TDA 625 integer :: n2b_co !< ncb_co for TDA calculations, nvb_co + ncb_co for non-TDA 626 integer :: nspin 627 integer :: nspinor = 1 !< nspinor = 2 if doing two-component spinor calculation; 1 is default 628 integer :: qflag !< =0 for finite Q calculation with arbitrary Q (deprecated) 629 !< =1 for Q=0 calculation (default) 630 !< =2 use Q commensurate with WFN_co kgrid (under construction) 631 logical :: no_mtxel !< if .true. we will not use the dipole matrix elements in an 632 !! absorption calculation. Useful for finite-q calculations 633 logical :: read_kpoints 634 integer :: ipar 635 integer :: iscreen !< what type of screening is present. 0 = semiconductor, 1 = graphene, 2 = metal 636 logical :: renorm_transf !< renormalize the dcc/dvv interpolation transformation? 637 !> Calculate kernel blocks other than (vc)->(v'c') transitions? This will 638 !! even include transitions such as (e,e)->(e',e'). In principle, we should 639 !! always include these blocks if we are interpolating later on, but they 640 !! are typically not important for semiconductors within TDA. 641 logical :: extended_kernel 642 !> If true, we extend the co/fi transf. to WFNs of different characters: 643 !! |v> = \sum_n` d_vn`|n`>, where |n`> can be |v`> or |c`> 644 !! If false, we restrict the character of the expansion coefficients: 645 !! |v> = \sum_v` d_vv`|v`> 646 logical :: unrestricted_transf 647 !> Zero out dvc/dcv coefficients 648 logical :: zero_unrestricted_contrib 649 logical :: patched_sampling !< simplest case of non-uniform sampling. See absorption.inp. 650 logical :: patched_sampling_co !< Use non-uniform sampling for coarse grid for Kernel. See absorption.inp. 651 integer :: zero_q0_element !< Zero q=0 matrix element of BSE Hamiltonian? See absorption.inp 652 logical :: tda !< use Tamm-Dancoff approximation? (Absorption only) 653 logical :: zero_coupling_block !< If true, zero hbse_b before calling p*bseig 654 integer :: iabsorp0 !< 1 means noeh_only, 0 otherwise 655 integer :: iwriteint = 1 !< = 0 for comm_disk, = 1 for comm_mpi 656 logical :: eqp_corrections !< do we use eqp.dat and eqp_q.dat 657 logical :: eqp_co_corrections !< do we use eqp_co.dat 658 logical :: eqp_co_q_corrections !< do we use eqp_co_q.dat 659 660!> For Coulomb interaction truncation 661 integer :: iwritecoul 662 integer :: icutv !< icutv encodes presence and type of truncation 663 real(DP) :: truncval(3) !< in Bohr (au) 664 integer :: nint !< number of intervals used in 665 !! double integral truncated_factor 666 logical :: use_hdf5 !< with -DHDF5, whether or not we actually use hdf5 667 logical :: bLowComm !< If this is true, each processor will store the entire epsilon matrix 668 logical :: delaunay_interp!< use Delaunay interpolation? 669 integer :: neps !< Number of G vectors to capture the dielectric cutoff 670 integer :: ilowmem 671 logical :: skipinterp 672 integer :: ivpar, icpar 673 integer :: nn !< PlotXct restrict_kpoints 674 integer :: ng 675 integer :: nktotal !< total number of unit cells 676 !> Number of vertices in co k-grid that are used to expand each k-point in 677 !! the fine grid for the **kernel** interpolation. This is 1 for the 678 !! greedy interpolation (previous behaviour of the code), and ndims+1 679 !! if we are performing Delaunay interpolation. 680 integer :: npts_intp_kernel 681 real(DP) :: eta !< energy resolution 682 real(DP) :: sigma !< (used to calculate the optical spectrum) 683 real(DP) :: gamma !< (used to calculate the optical spectrum) 684 real(DP) :: qshift 685 real(DP) :: shift(3) !< shift vector (this is the small shift, 686 !< used to generate WFNq_fi, referenced only if xct%read_kpoints) 687 real(DP) :: finiteq(3) !< center-of-mass momentum of exciton 688 logical :: energy_loss !< calculate energy loss spectrum 689 real(DP) :: lpol !< norm of pol 690 real(DP) :: pol(3) !< light polarization for transition matrix elements 691 integer :: npol !< number of polarizations we have. Either 1 or 3. 692 integer :: nmtxmax !< max. number of columns in epsmat or eps0mat 693 integer :: theory !< theory level in kernel calculation 694 !< 0 - GW-BSE, 1 - TDDFT 695 integer :: qgrid(3) 696 real(DP) :: q0vec(3) ! This is a hack for passing q0vec for 697 ! TDDFT calculations (never used otherwise) 698 ! when there is no epsilon 699 type(coulomb_modifier_t) :: coulomb_mod 700 logical :: coul_mod_flag !< Flag which tells if the coulomb interaction has been 701 !< modified (for TDDFT calculations) 702 integer, pointer :: indexq(:), indexq_fi(:) !< When exciton has finite center-of-mass momentum, 703 !< maps between valence states at k+Q and conduction states at k 704 integer, pointer :: isrtqi(:,:), nmtxa(:) 705 integer, pointer :: ifmax(:,:), ifmaxq(:,:) 706 real(DP) :: ecute !< energy cutoff used in dielectric matrix 707 real(DP) :: scaling !< multiply kernel by arbitrary factor 708 real(DP) :: ecutg !< energy cutoff used in wavefunctions and interaction 709 !< kernel, see Rohlfing & Louie, PRB 62(8),p. 4938 710 !< (must be slightly longer than xct%ecute because of umklapp vectors) 711 real(DP) :: efermi !< computed efermi 712 real(DP) :: efermi_input !< as set in input file 713 logical :: rfermi !< relative or absolute Fermi level 714 SCALAR, pointer :: epsdiag(:,:) !< (nmtxmax, nq+1) 715 type (wpgen) :: wpg !< charge density/fxc for TDDFT 716 ! FHJ: TODO - move the following quantities to a separate derived type, 717 ! ir reuse epsmpi 718 !> Regular comm: (nmtxmax, ngpown_max, nq+1). The local processor stores row ig 719 !! and a "local column" igp_l from epsilon in epscol(ig, igp_l, ik). 720 !! Low comm: (nmtxmax, nmtxmax, nq+1). Each processor stores all eps(0)mat. 721 !! local column = global epsilon column. 722 SCALAR, pointer :: epscol(:,:,:) 723 integer :: ngpown_max !< max. number of eps columns a PE can have 724 integer :: ngpown !< number of columns of epsinv I own, for largest matrix 725 integer :: nb !< block size for column distribution 726 ! The arrays epsown and epsowni were removed, as they can be calculated on the fly: 727 ! iowner = INDXG2P(igp, xct%nb, peinf%inode, 0, peinf%npes) 728 ! xct%epsown(igp) = INDXG2P(igp, xct%nb, peinf%inode, 0, peinf%npes) 729 ! xct%epsowni(igp_loc, iowner+1) = INDXG2L(igp, xct%nb, peinf%inode, 0, peinf%npes) 730!> Used for screened_exchange 731 logical :: screen_exchange !< add background screening to exchange term 732 SCALAR, pointer :: epscol_bg(:,:,:) !< stores columns of the substrate dielectric matrix 733 SCALAR, pointer :: epsdiag_bg(:,:) !< same as epsdiag for substrate dielectric matrix 734 integer :: nmtxmax_bg !< max. number of columns in epsmat or eps0mat 735 integer, pointer :: isrtqi_bg(:,:), nmtxa_bg(:) 736 real(DP) :: ecute_bg !< energy cutoff used for background dielectric matrix 737 integer :: ngpown_max_bg !< max. number of columns of background epsmat a PE can have 738 integer :: ngpown_bg !< number of columns of background epsinv I own, for largest matrix 739 740!> Used in haydock/diag only 741 integer :: nkpt_fi !< number of kpts in the fine grid 742 integer :: nkptq_fi !< number of kpts in the q-shifted fine grid 743 integer :: nvb_fi !< number of valence bands in the fine grid 744 integer :: ncb_fi !< number of conduction bands in the fine grid 745 real(DP) :: avgcut !< Cut in which we do cell averages on 746 real(DP) :: wplasmon 747 SCALAR :: wcoul0 748 integer :: vmin,vmax !< Indices of the lowest and highest 749 !! occupied bands, if specified in input file. 750 integer :: rgrid(3) !< regular grid used to calculate qpt_averages (default is kgrid) 751 logical :: freplacebz 752 logical :: fwritebz 753 logical :: degeneracy_check_override 754 logical :: die_outside_sphere 755 logical :: averagew 756 logical :: subsample_line !< during kernel interpolation, replace interpolated matrix 757 !< elements with matrix elements from a precalculated bsemat file 758 !< on a subsampled grid 759 real(DP) :: subsample_cutoff !< use subsampled BSE matrix elements when |q| is less than cutoff 760 real(DP) :: exchange_fact !< multiplies the BSE exchange term by this factor 761 real(DP) :: direct_fact !< multiplies the BSE direct term by this factor 762 real(DP) :: delta_frequency !< Frequency step for absorption spectrum 763 end type xctinfo 764 765!------------------------------------ 766 767 type flags 768 769!> 770!> Used in haydock, diag, nonlinearoptics 771!> 772!> Defined flags: 773!> 774!> bz0 = 0 --> use symmetries to unfold the Brillouin zone in WFN_fi file 775!> 1 --> do not unfold the BZ in WFN_fi file (default) 776!> bzq = 0 --> use symmetries to unfold the BZ in WFNq_fi file 777!> 1 --> do not unfold the BZ in WFNq_fi file (default) 778!> bzc = 0 --> use symmetries to unfold the BZ in WFN_co file 779!> 1 --> do not unfold the BZ in WFN_co file (default) 780!> bzcq = 0 --> use symmetries to unfold the BZ in WFNq_co file 781!> 1 --> do not unfold the BZ in WFN_co file (default) 782!> 783!> read_dtmat = false --> calculate dcc,dvv matrices (default) 784!> true --> take dcc,dvv matrices from file dtmat 785!> 786!> eig = 0 --> do not write eigenvectors (default) 787!> < 0 --> write all eigenvectors 788!> > 0 --> write the first flag%eig eigenvectors 789!> 790!> read_epsdiag = false --> read files 'eps0mat'/'epsmat' (default) 791!> true --> read file 'epsdiag.dat' 792!> 793!> krnl = 0 --> spin triplet kernel, direct kernel only (only allow for nspin = 1) 794!> 1 --> spin singlet kernel (default) 795!> 2 --> local-fields + RPA, exchange kernel only 796!> 3 --> spinor kernel 797!> 798!> opr = 0 --> use velocity operator 799!> 1 --> use momentum operator 800!> 2 --> use JDOS operator (Haydock only) 801!> 802!> lor = 0 --> use Lorentzian broadening 803!> 1 --> use Gaussian broadening 804!> 2 --> use Voigt broadening 805!> 806!> spec = 0 --> go through the whole exciton calculation (default) 807!> 1 --> calculate only absorption spectrum (this option skips 808!> all calculation and goes right to the end of the code) 809!> 810!> vm = 0 --> calculate velocity/momentum matrix elements (default) 811!> 1 --> read velocity/momentum matrix elements from file vmtxel 812!> 2 --> use vectors from previous iteration (Haydock only!) 813!> 814!> job = 0 --> ultrafast calculation 815!> job = 1 --> two-photon calculation 816!> 817 818 integer :: bz0 819 integer :: lor 820 integer :: bzq 821 integer :: bzc 822 integer :: bzcq 823 logical :: read_dtmat 824 logical :: read_dtmat_sub 825 integer :: eig 826 logical :: read_epsdiag 827 integer :: krnl 828 integer :: opr 829 integer :: spec 830 integer :: vm 831 integer :: job 832 !> Use averaged Gauss quadrature in Lanczos algorithm? Default is true. 833 logical :: lanczos_gauss_quad 834 logical :: debug_lanczos 835 end type flags 836 837!--------------------------------- 838 839 type otherinfo 840 integer :: ithreeD 841 integer :: knx 842 integer :: kny 843 integer :: knz 844 real(DP) :: keta 845 end type otherinfo 846 847!--------------------------------- 848 849 type windowinfo 850 real(DP) :: evalue 851 real(DP) :: emax 852 real(DP), pointer :: cstates(:) 853 real(DP), pointer :: estates(:) 854 integer, pointer :: istates(:) 855 integer :: nstates 856 end type windowinfo 857 858!----------------------------- 859 !> coarse-grid wavefunctions for diag/haydock 860 type tdistgwf 861 862 integer :: block_sz !< block size for BLACS distribution = DIVUP(ngm, npes) 863 integer :: ngm !< maximum number of G-vectors = kp*%ngkmax 864 integer :: ngl !< local number of G-vectors = NUNROC(...). At most block_sz. 865 integer :: tgl !< local to global translation = block_sz * peinf%inode 866 867!> local to global index translation : ig_g = ig_l + tgl 868!> ig_g = 1 ... ng(ik) is the global G-index 869!> ig_l = 1 ... ngl is the local G-index 870 871 integer :: nk !< number of k-points 872 integer :: ns !< number of spin components 873 integer :: nspinor = 1 !< nspinor = 2 if doing two-component spinor calculation; 1 is default 874 integer :: nv !< number of valence bands 875 integer :: nc !< number of conduction bands 876 877 integer, pointer :: ng(:) !< (nk) 878 integer, pointer :: isort(:,:) !< (ngl,nk) 879 SCALAR, pointer :: zv(:,:,:,:) !< (ngl,nv,ns*nspinor,nk) 880 SCALAR, pointer :: zc(:,:,:,:) !< (ngl,nc,ns*nspinor,nk) 881 882 end type tdistgwf 883 884!----------------------------- 885!> MJ: work arrays - getting rid of save statements 886!! Mostly used in genwf 887 type work_genwf 888 integer :: ikold = 0 889 integer :: nb !< Number of bands 890 integer :: ng !< Number of planewaves 891 integer :: ns !< Number of spins 892 integer :: nspinor = 1 !< nspinor = 2 if doing two-component spinor calculation; 1 is default 893 SCALAR, pointer :: cg(:,:,:) !< Planewave coefficients (ng,nb,ns*nspinor) 894 SCALAR, pointer :: ph(:) !< (ng) 895 integer, pointer :: ind(:) !< (ng) 896 integer, pointer :: isort(:) !< (ng) 897 end type work_genwf 898 899!----------------------------- 900!> (gsm) work arrays - getting rid of save statements 901 902 type twork_scell 903 integer :: dNfft(3) 904 integer :: Nplane 905 integer :: Nrod 906 complex(DPC), pointer :: fftbox_1D(:,:) 907 end type twork_scell 908 909 !> FHJ: mean-field header 910 type mf_header_t 911 integer :: version 912 character(len=3) :: sheader 913 character(len=32) :: sdate 914 character(len=32) :: stime 915 integer :: iflavor 916 type(crystal) :: crys 917 type(kpoints) :: kp 918 type(symmetry) :: syms 919 type(gspace):: gvec 920 end type mf_header_t 921 922 !> FHJ: header information for kernel files (bsedmat, bsexmat, bsemat.h5) 923 type kernel_header_t 924 925 ! Mean-field and other general information 926 type(mf_header_t) :: mf !< mf header containing number of k-points, WFN cutoff, etc. 927 928 integer :: version 929 integer :: iflavor 930 931 integer :: iscreen !< screening flag 932 integer :: icutv !< truncation flag 933 real(DP) :: ecuts !< epsilon cutoff 934 real(DP) :: ecutg !< WFN cutoff 935 real(DP) :: efermi !< Fermi energy found by the code after any shift 936 integer :: theory !< 0 for GW-BSE, 1 for TD-HF, 2 for TD-DFT 937 !> How many transitions blocks are there in the kernel matrix? 938 !! 1 for restricted TDA kernel: vc -> v`c` 939 !! 2 for restricted non-TDA kernel: {vc,cv} -> {v`c`,c`v`} [not implemented] 940 !! 4 for extended kernel: {n1,n2} -> {n1`,n2`} 941 integer :: nblocks 942 integer :: storage !< 0 if storing full matrix (only option supported now) 943 integer :: nmat !< number of matrices in the file (1 for bsexmat, 3 for bsedmat) 944 logical :: energy_loss !< is this an energy-loss calculation? 945 946 integer :: nvb !< number of valence bands in the coarse grid 947 integer :: ncb !< number of conduction bands in the coarse grid 948 integer :: n1b !< nvb_co if kernel_sz==1; nvb_co + ncb_co if kernel_sz=4 949 integer :: n2b !< ncb_co if kernel_sz==1; nvb_co + ncb_co if kernel_sz=4 950 integer :: ns !< number of spins 951 integer :: nspinor = 1 !< nspinor = 2 if doing two-component spinor calculation; 1 is default 952 logical :: patched_sampling !< are we doing a calculation on a patch? 953 954 ! Variables specific to kernel files 955 integer :: nk !< number of k-points 956 real(DP), pointer :: kpts(:,:) 957 integer :: kgrid(3) 958 !> 0 for finite Q calculation with arbitrary Q (deprecated) 959 !! 1 for Q=0 calculation (default) 960 !! 2 use Q commensurate with WFN_co kgrid (under construction) 961 integer :: qflag 962 real(DP) :: center_mass_q(3) 963 964 end type kernel_header_t 965 966end module typedefs_m 967