TABLE OF CONTENTS
ABINIT/invars1 [ Functions ]
NAME
invars1
FUNCTION
Initialize the dimensions needed to allocate the input arrays for one dataset characterized by jdtset, by taking from string the necessary data. Perform some preliminary checks and echo these dimensions.
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR, AR, MKV) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt.
INPUTS
iout=unit number of output file jdtset=number of the dataset looked for lenstr=actual length of string msym=default maximal number of symmetries npsp1= number of pseudopotential files zionpsp(npsp1)= valence charge over all psps
OUTPUT
mband_upper=estimation of the maximum number of bands for any k-point
SIDE EFFECTS
Input/Output (the default value is given in the calling routine) dtset=<type datafiles_type>contains all input variables, some of which are initialized here, while other were already initialized, while some others will still be initialized later. The list of records of dtset initialized in the present routine is: acell_orig,densty,iatfix,kptopt,kptrlatt, mkmem,mkqmem,mk1mem,natsph,natvshift,nconeq,nkpt,nkptgw,nkpthf, nqptdm,nshiftk,nucdipmom,nzchempot,optdriver, rprim_orig,rprimd_orig,shiftk, spgroup,spinat,typat,vel_orig,vel_cell_orig,xred_orig bravais(11)=characteristics of Bravais lattice (see symlatt.F90) symafm(1:msym)=(anti)ferromagnetic part of symmetry operations symrel(3,3,1:msym)=symmetry operations in real space in terms of primitive translations tnons(3,1:msym)=nonsymmorphic translations for symmetry operations string*(*)=string of characters containing all input variables and data
NOTES
Must set up the geometry of the system, needed to compute k point grids in an automatic fashion. Treat separately mband_upper, since fband, charge and zionpsp must be known for being able to initialize it. Defaults are provided in the calling routine. Defaults are also provided here for the following variables : mband_upper, occopt, fband, charge They should be kept consistent with defaults of the same variables provided to the invars routines.
PARENTS
invars1m
CHILDREN
atomdata_from_znucl,chkint_ge,ingeo,inkpts,inqpt,intagm,inupper invacuum,mkrdim,wrtout
SOURCE
69 #if defined HAVE_CONFIG_H 70 #include "config.h" 71 #endif 72 73 #include "abi_common.h" 74 75 subroutine invars1(bravais,dtset,iout,jdtset,lenstr,mband_upper,msym,npsp1,& 76 & string,symafm,symrel,tnons,zionpsp) 77 78 use defs_basis 79 use defs_abitypes 80 use m_errors 81 use m_profiling_abi 82 use m_xmpi 83 use m_atomdata 84 85 use m_fstrings, only : inupper 86 use m_geometry, only : mkrdim 87 use m_parser, only : intagm, chkint_ge 88 89 !This section has been created automatically by the script Abilint (TD). 90 !Do not modify the following lines by hand. 91 #undef ABI_FUNC 92 #define ABI_FUNC 'invars1' 93 use interfaces_14_hidewrite 94 use interfaces_57_iovars, except_this_one => invars1 95 !End of the abilint section 96 97 implicit none 98 99 !Arguments ------------------------------------ 100 !scalars 101 integer,intent(in) :: iout,jdtset,lenstr,msym,npsp1 102 integer,intent(out) :: mband_upper 103 character(len=*),intent(inout) :: string 104 type(dataset_type),intent(inout) :: dtset 105 !arrays 106 integer,intent(inout) :: bravais(11),symafm(msym),symrel(3,3,msym) 107 real(dp),intent(inout) :: tnons(3,msym) 108 real(dp),intent(in) :: zionpsp(npsp1) 109 110 !Local variables------------------------------- 111 character :: blank=' ',string1 112 !scalars 113 integer :: chksymbreak,found,ierr,iatom,ii,ikpt,iimage,index_blank,index_lower 114 integer :: index_typsymb,index_upper,ipsp,iscf,intimage,itypat,leave,marr 115 integer :: natom,nkpt,nkpthf,npsp,npspalch 116 integer :: nqpt,nspinor,nsppol,ntypat,ntypalch,ntyppure,occopt,response 117 integer :: rfddk,rfelfd,rfphon,rfstrs,rfuser,rf2_dkdk,rf2_dkde,rfmagn 118 integer :: tfband,tnband,tread,tread_alt 119 real(dp) :: charge,fband,kptnrm,kptrlen,zelect,zval 120 character(len=2) :: string2,symbol 121 character(len=500) :: message 122 type(atomdata_t) :: atom 123 !arrays 124 integer :: cond_values(4),vacuum(3) 125 integer,allocatable :: iatfix(:,:),intarr(:),istwfk(:),nband(:),typat(:) 126 real(dp) :: acell(3),rprim(3,3) 127 !real(dp) :: field(3) 128 real(dp),allocatable :: amu(:),dprarr(:),kpt(:,:),kpthf(:,:),mixalch(:,:),nucdipmom(:,:) 129 real(dp),allocatable :: ratsph(:),reaalloc(:),spinat(:,:) 130 real(dp),allocatable :: vel(:,:),vel_cell(:,:),wtk(:),xred(:,:),znucl(:) 131 character(len=32) :: cond_string(4) 132 133 134 !************************************************************************ 135 136 !Some initialisations 137 ierr=0 138 cond_string(1:4)=' ' 139 cond_values(1:4)=(/0,0,0,0/) 140 141 !Read parameters 142 marr=dtset%npsp;if (dtset%npsp<3) marr=3 143 marr=max(marr,dtset%nimage) 144 ABI_ALLOCATE(intarr,(marr)) 145 ABI_ALLOCATE(dprarr,(marr)) 146 147 !--------------------------------------------------------------------------- 148 149 rfddk=0; rfelfd=0; rfphon=0; rfmagn=0; rfstrs=0; rfuser=0; rf2_dkdk=0; rf2_dkde=0 150 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfddk',tread,'INT') 151 if(tread==1) rfddk=intarr(1) 152 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfelfd',tread,'INT') 153 if(tread==1) rfelfd=intarr(1) 154 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfmagn',tread,'INT') 155 if(tread==1) rfmagn=intarr(1) 156 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfphon',tread,'INT') 157 if(tread==1) rfphon=intarr(1) 158 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfstrs',tread,'INT') 159 if(tread==1) rfstrs=intarr(1) 160 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfuser',tread,'INT') 161 if(tread==1) rfuser=intarr(1) 162 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rf2_dkdk',tread,'INT') 163 if(tread==1) rf2_dkdk=intarr(1) 164 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rf2_dkde',tread,'INT') 165 if(tread==1) rf2_dkde=intarr(1) 166 167 response=0 168 if(rfddk/=0.or.rf2_dkdk/=0.or.rf2_dkde/=0.or.rfelfd/=0.or.rfphon/=0.or.rfstrs/=0.or.rfuser/=0.or.rfmagn/=0)response=1 169 170 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'optdriver',tread,'INT') 171 if (tread==1) then 172 dtset%optdriver=intarr(1) 173 else 174 ! If optdriver was not read, while response=1, set optdriver to 1 175 if(response==1)dtset%optdriver=1 176 end if 177 178 !--------------------------------------------------------------------------- 179 !For now, waiting express parallelisation for recursion 180 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tfkinfunc',tread,'INT') 181 if(tread==1) dtset%tfkinfunc=intarr(1) 182 183 !--------------------------------------------------------------------------- 184 ! wvl_bigdft_comp, done here since default values of nline, nwfshist and iscf 185 ! depend on its value (see indefo) 186 if(dtset%usewvl==1) then 187 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'wvl_bigdft_comp',tread,'INT') 188 if(tread==1) dtset%wvl_bigdft_comp=intarr(1) 189 end if 190 191 !--------------------------------------------------------------------------- 192 193 natom=dtset%natom 194 npsp=dtset%npsp 195 ntypat=dtset%ntypat 196 197 !No default value for znucl 198 call intagm(dprarr,intarr,jdtset,marr,dtset%npsp,string(1:lenstr),'znucl',tread,'DPR') 199 if(tread==1)then 200 dtset%znucl(1:dtset%npsp)=dprarr(1:dtset%npsp) 201 end if 202 if(tread/=1)then 203 write(message, '(3a)' )& 204 & 'The array znucl MUST be initialized in the input file while this is not done.',ch10,& 205 & 'Action: initialize znucl in your input file.' 206 MSG_ERROR(message) 207 end if 208 209 !The default for ratsph has already been initialized 210 call intagm(dprarr,intarr,jdtset,marr,dtset%ntypat,string(1:lenstr),'ratsph',tread,'LEN') 211 if(tread==1)then 212 do ii=1,dtset%ntypat 213 dtset%ratsph(ii)=dprarr(ii) 214 end do 215 end if 216 ABI_ALLOCATE(ratsph,(dtset%ntypat)) 217 do ii=1,dtset%ntypat 218 ratsph(ii)=dtset%ratsph(ii) 219 end do 220 221 !Special treatment of _TYPAX (from a XYZ file), taking into account 222 !the fact that znucl does NOT depend on the dataset 223 !Examine all occurences of '_TYPAX' 224 225 do 226 index_typsymb=index(string(1:lenstr),'_TYPAX') 227 if(index_typsymb==0)exit 228 ! Replace '_TYPAX' by '_TYPAT' 229 string(index_typsymb:index_typsymb+5)='_TYPAT' 230 index_upper=index_typsymb+5 231 ! Must start from the first blank after the tag (including possible dtset_char) 232 index_upper=index(string(index_upper:lenstr),blank)+index_upper-1 233 index_lower=index_upper 234 235 ! Examine all atoms (the end of the symbol string is delimited by a XX ) 236 do 237 index_blank=index(string(index_upper:lenstr),blank)+index_upper-1 238 string2=string(index_blank+1:index_blank+2) 239 if(string2=="XX")exit 240 found=0 241 ! Find the matching symbol 242 do ipsp=1,dtset%npsp 243 call atomdata_from_znucl(atom,dtset%znucl(ipsp)) 244 symbol = atom%symbol 245 call inupper(symbol) 246 call inupper(string2) 247 248 ! DEBUG 249 ! write(std_out,'(a)')' invars1 : before test, trim(adjustl(symbol)),trim(adjustl(string2))' 250 ! write(std_out,'(5a)' )'"',trim(adjustl(symbol)),'","',trim(adjustl(string2)),'"' 251 ! ENDDEBUG 252 253 if(trim(adjustl(symbol))==trim(adjustl(string2)))then 254 found=1 255 index_upper=index_blank+1 256 ! Cannot deal properly with more that 9 psps 257 if(ipsp>=10)then 258 message = 'Need to use a pseudopotential with number larger than 9. Not allowed yet.' 259 MSG_ERROR(message) 260 end if 261 262 ! DEBUG 263 ! write(std_out,*)' invars1 : found ipsp=',ipsp 264 ! ENDDEBUG 265 266 write(string1,'(i1)')ipsp 267 string(index_lower:index_lower+1)=blank//string1 268 index_lower=index_lower+2 269 end if 270 end do ! ipsp 271 ! if not found ... 272 if(found==0)then 273 write(message,'(6a)' )& 274 & 'Did not find matching pseudopotential for XYZ atomic symbol,',ch10,& 275 & 'with value ',string2,ch10,& 276 & 'Action: check that the atoms required by the XYZ file correspond to one psp file.' 277 MSG_ERROR(message) 278 end if 279 end do ! Loop on atoms 280 ! One should find blanks after the last significant type value 281 string(index_lower:index_blank+2)=blank 282 end do ! loop to identify _TYPAX 283 284 !--------------------------------------------------------------------------- 285 286 !Here, set up quantities that are related to geometrical description 287 !of the system (acell,rprim,xred), as well as 288 !initial velocity(vel), and spin of atoms (spinat), nuclear dipole moments 289 ! of atoms (nucdipmom), 290 !the symmetries (symrel,symafm, and tnons) 291 !and the list of fixed atoms (iatfix,iatfixx,iatfixy,iatfixz). 292 !Arrays have already been 293 !dimensioned thanks to the knowledge of msym and mxnatom 294 295 !ji: We need to read the electric field before calling ingeo 296 !****** Temporary ****** 297 298 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'berryopt',tread,'INT') 299 if(tread==1) dtset%berryopt=intarr(1) 300 301 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'berrysav',tread,'INT') 302 if(tread==1) dtset%berrysav=intarr(1) 303 304 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'bfield',tread,'DPR') 305 if (tread==1) dtset%bfield(1:3) = dprarr(1:3) 306 307 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'dfield',tread,'DPR') 308 if (tread==1) dtset%dfield(1:3) = dprarr(1:3) 309 310 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'efield',tread,'DPR') 311 if (tread==1) dtset%efield(1:3) = dprarr(1:3) 312 313 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'red_dfield',tread,'DPR') 314 if (tread==1) dtset%red_dfield(1:3) = dprarr(1:3) 315 316 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'red_efield',tread,'DPR') 317 if (tread==1) dtset%red_efield(1:3) = dprarr(1:3) 318 319 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'red_efieldbar',tread,'DPR') 320 if (tread==1) dtset%red_efieldbar(1:3) = dprarr(1:3) 321 322 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'jfielddir',tread,'INT') 323 if(tread==1) dtset%jfielddir(1:3)=intarr(1:3) 324 325 !We need to know nsppol/nspinor/nspden before calling ingeo 326 nsppol=dtset%nsppol 327 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nsppol',tread,'INT') 328 if(tread==1) nsppol=intarr(1) 329 330 !Alternate SIESTA definition of nsppol 331 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'SpinPolarized',tread_alt,'LOG') 332 if(tread_alt==1)then 333 if(tread==1)then 334 write(message, '(a,a,a,a,a,a,a,a)' ) ch10,& 335 & ' invars1: ERROR -',ch10,& 336 & ' nsppol and SpinPolarized cannot be specified simultaneously',ch10,& 337 & ' for the same dataset.',ch10,& 338 & ' Action : check the input file.' 339 call wrtout(std_out, message,'COLL') 340 leave=1 341 else 342 ! Note that SpinPolarized is a logical input variable 343 nsppol=1 344 if(intarr(1)==1)nsppol=2 345 tread=1 346 end if 347 end if 348 dtset%nsppol=nsppol 349 350 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nspinor',tread,'INT') 351 if(tread==1) dtset%nspinor=intarr(1) 352 353 !Has to read pawspnorb now, in order to adjust nspinor 354 !Also, if nspinor=1, turn on spin-orbit coupling by default, here for the PAW case. NC case is treated elsewhere. 355 if (dtset%usepaw>0)then 356 ! Change the default value 357 if(dtset%nspinor==2)dtset%pawspnorb=1 358 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawspnorb',tread,'INT') 359 if(tread==1)then 360 dtset%pawspnorb=intarr(1) 361 if(dtset%pawspnorb>0) dtset%nspinor=2 362 else 363 if(dtset%nspinor==2)then 364 write(message, '(4a)' ) ch10,& 365 & ' invars1: COMMENT -',ch10,& 366 & ' With nspinor=2 and usepaw=1, pawspnorb=1 has been switched on by default.' 367 call wrtout(iout, message,'COLL') 368 end if 369 end if 370 end if 371 nspinor=dtset%nspinor 372 373 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nspden',tread,'INT') 374 if(tread==1) then 375 dtset%nspden=intarr(1) 376 else 377 dtset%nspden=dtset%nsppol 378 end if 379 380 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ntypalch',tread,'INT') 381 if(tread==1) dtset%ntypalch=intarr(1) 382 383 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nzchempot',tread,'INT') 384 if(tread==1) dtset%nzchempot=intarr(1) 385 386 ntypalch=dtset%ntypalch 387 if(ntypalch>ntypat)then 388 write(message, '(3a,i0,a,i0,a,a)' )& 389 & 'The input variable ntypalch must be smaller than ntypat, while it is',ch10,& 390 & 'ntypalch=',dtset%ntypalch,', and ntypat=',ntypat,ch10,& 391 & 'Action: check ntypalch vs ntypat in your input file.' 392 MSG_ERROR(message) 393 end if 394 395 ntyppure=ntypat-ntypalch 396 dtset%ntyppure=ntyppure 397 npspalch=npsp-ntyppure 398 dtset%npspalch=npspalch 399 if(npspalch<0)then 400 write(message, '(a,i0,2a,i0,a,a)' )& 401 & 'The number of available pseudopotentials, npsp=',npsp,ch10,& 402 & 'is smaller than the requested number of types of pure atoms, ntyppure=',ntyppure,ch10,& 403 & 'Action: check ntypalch versus ntypat and npsp in your input file.' 404 MSG_ERROR(message) 405 end if 406 407 if(ntypalch>0)then 408 call intagm(dprarr,intarr,jdtset,marr,ntypalch,string(1:lenstr),'algalch',tread,'INT') 409 if(tread==1) dtset%algalch(1:ntypalch)=intarr(1:ntypalch) 410 end if 411 412 !Read the Zeeman field 413 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'zeemanfield',tread,'BFI') 414 if(tread==1) then 415 if(dtset%nspden == 2)then 416 write(message,'(7a)')& 417 & 'A Zeeman field has been specified without noncollinear spins.',ch10,& 418 & 'Only the z-component of the magnetic field will be used.' 419 MSG_WARNING(message) 420 else if (dtset%nspden == 1)then 421 write(message, '(a,a,a)' )& 422 & 'A Zeeman field has been specified for a non-spin-polarized calculation.',ch10,& 423 & 'Action: check the input file.' 424 MSG_ERROR(message) 425 end if 426 427 dtset%zeemanfield(1:3) = dprarr(1:3) 428 end if 429 430 431 ABI_ALLOCATE(amu,(ntypat)) 432 ABI_ALLOCATE(mixalch,(npspalch,ntypalch)) 433 ABI_ALLOCATE(vel,(3,natom)) 434 ABI_ALLOCATE(vel_cell,(3,3)) 435 ABI_ALLOCATE(xred,(3,natom)) 436 intimage=2 ; if(dtset%nimage==1)intimage=1 437 do ii=1,dtset%nimage+1 438 iimage=ii 439 if(dtset%nimage==1 .and. ii==2)exit 440 if(dtset%nimage==2 .and. ii==3)exit 441 if(dtset%nimage> 2 .and. ii==intimage)cycle ! Will do the intermediate reference image at the last reading 442 if(dtset%nimage>=2 .and. ii==dtset%nimage+1)iimage=intimage 443 444 write(message,'(a,i0)')' invars1 : treat image number: ',iimage 445 call wrtout(std_out,message,'COLL') 446 447 ! Need to reset nsym to default value for each image 448 dtset%nsym=0 449 450 ! Call ingeo for each image in turn, with the possible default values 451 acell=dtset%acell_orig(1:3,iimage) 452 amu=dtset%amu_orig(1:ntypat,iimage) 453 mixalch=dtset%mixalch_orig(1:npspalch,1:ntypalch,iimage) 454 rprim=dtset%rprim_orig(1:3,1:3,iimage) 455 vel=dtset%vel_orig(1:3,1:natom,iimage) 456 vel_cell=dtset%vel_cell_orig(1:3,1:3,iimage) 457 xred=dtset%xred_orig(1:3,1:natom,iimage) 458 ABI_ALLOCATE(iatfix,(3,natom)) 459 ABI_ALLOCATE(nucdipmom,(3,natom)) 460 ABI_ALLOCATE(spinat,(3,natom)) 461 ABI_ALLOCATE(typat,(natom)) 462 ABI_ALLOCATE(znucl,(dtset%npsp)) 463 nucdipmom(1:3,1:natom)=dtset%nucdipmom(1:3,1:natom) 464 spinat(1:3,1:natom)=dtset%spinat(1:3,1:natom) 465 znucl(1:dtset%npsp)=dtset%znucl(1:dtset%npsp) 466 467 call ingeo(acell,amu,dtset,bravais,dtset%genafm(1:3),iatfix,& 468 & dtset%icoulomb,iimage,iout,jdtset,dtset%jellslab,lenstr,mixalch,& 469 & msym,natom,dtset%nimage,dtset%npsp,npspalch,dtset%nspden,dtset%nsppol,& 470 & dtset%nsym,ntypalch,dtset%ntypat,nucdipmom,dtset%nzchempot,& 471 & dtset%pawspnorb,dtset%ptgroupma,ratsph,& 472 & rprim,dtset%slabzbeg,dtset%slabzend,dtset%spgroup,spinat,& 473 & string,symafm,dtset%symmorphi,symrel,tnons,dtset%tolsym,typat,vel,vel_cell,xred,znucl) 474 dtset%iatfix(1:3,1:natom)=iatfix(1:3,1:natom) 475 dtset%nucdipmom(1:3,1:natom)=nucdipmom(1:3,1:natom) 476 dtset%spinat(1:3,1:natom)=spinat(1:3,1:natom) 477 dtset%typat(1:natom)=typat(1:natom) 478 ABI_DEALLOCATE(iatfix) 479 ABI_DEALLOCATE(nucdipmom) 480 ABI_DEALLOCATE(spinat) 481 ABI_DEALLOCATE(typat) 482 ABI_DEALLOCATE(znucl) 483 dtset%acell_orig(1:3,iimage)=acell 484 dtset%amu_orig(1:ntypat,iimage)=amu 485 dtset%mixalch_orig(1:npspalch,1:ntypalch,iimage)=mixalch 486 dtset%rprim_orig(1:3,1:3,iimage)=rprim 487 dtset%vel_orig(1:3,1:natom,iimage)=vel 488 dtset%vel_cell_orig(1:3,1:3,iimage)=vel_cell 489 dtset%xred_orig(1:3,1:natom,iimage)=xred 490 call mkrdim(dtset%acell_orig(1:3,iimage),dtset%rprim_orig(1:3,1:3,iimage),dtset%rprimd_orig(1:3,1:3,iimage)) 491 end do 492 ABI_DEALLOCATE(amu) 493 ABI_DEALLOCATE(mixalch) 494 ABI_DEALLOCATE(vel) 495 ABI_DEALLOCATE(vel_cell) 496 ABI_DEALLOCATE(xred) 497 498 !Examine whether there is some vacuum space in the unit cell 499 call invacuum(jdtset,lenstr,natom,dtset%rprimd_orig(1:3,1:3,intimage),string,vacuum,& 500 & dtset%xred_orig(1:3,1:natom,intimage)) 501 502 !DEBUG 503 !write(std_out,*)' invars1: before inkpts, dtset%mixalch_orig(1:npspalch,1:ntypalch,:)=',& 504 !dtset%mixalch_orig(1:npspalch,1:ntypalch,1:dtset%nimage) 505 !stop 506 !ENDDEBUG 507 508 !--------------------------------------------------------------------------- 509 510 !Set up k point grid number 511 !First, get additional information 512 dtset%kptopt=1 513 if(dtset%nspden==4)dtset%kptopt=4 514 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'kptopt',tread,'INT') 515 if(tread==1) dtset%kptopt=intarr(1) 516 517 dtset%qptopt=1 518 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'qptopt',tread,'INT') 519 if(tread==1) dtset%qptopt=intarr(1) 520 521 iscf=5 522 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'iscf',tread,'INT') 523 if(tread==1) iscf=intarr(1) 524 525 526 dtset%natsph=dtset%natom 527 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natsph',tread,'INT') 528 if(tread==1) dtset%natsph=intarr(1) 529 530 dtset%natsph_extra=0 531 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natsph_extra',tread,'INT') 532 if(tread==1) dtset%natsph_extra=intarr(1) 533 534 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natvshift',tread,'INT') 535 if(tread==1) dtset%natvshift=intarr(1) 536 537 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nconeq',tread,'INT') 538 if(tread==1) dtset%nconeq=intarr(1) 539 540 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nkptgw',tread,'INT') 541 if(tread==1) dtset%nkptgw=intarr(1) 542 if (dtset%nkptgw<0) then 543 write(message, '(a,i0,a,a,a,a)' )& 544 & 'Input nkptgw must be >= 0, but was ',dtset%nkptgw,ch10,& 545 & 'This is not allowed.',ch10,& 546 & 'Action: check the input file.' 547 MSG_ERROR(message) 548 end if 549 550 !Number of points for long wavelength limit. 551 !Default is dtset%gw_nqlwl=0 552 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_nqlwl',tread,'INT') 553 if(tread==1) dtset%gw_nqlwl=intarr(1) 554 if (dtset%gw_nqlwl<0) then 555 write(message, '(a,i12,a,a,a,a)' )& 556 & 'Input gw_nqlwl must be > 0, but was ',dtset%gw_nqlwl,ch10,& 557 & 'This is not allowed.',ch10,& 558 & 'Action: check the input file.' 559 MSG_ERROR(message) 560 end if 561 562 !Read number of k-points (if specified) 563 nkpt=0 564 if(dtset%kptopt==0)nkpt=1 565 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nkpt',tread,'INT') 566 if(tread==1) nkpt=intarr(1) 567 dtset%nkpt=nkpt 568 569 call chkint_ge(0,0,cond_string,cond_values,ierr,'nkpt',nkpt,0,iout) 570 if(dtset%kptopt==0)then 571 cond_string(1)='kptopt' ; cond_values(1)=0 572 call chkint_ge(1,1,cond_string,cond_values,ierr,'nkpt',nkpt,1,iout) 573 end if 574 575 nkpthf=nkpt 576 dtset%nkpthf=nkpt 577 578 !Will compute the actual value of nkpt, if needed. Otherwise, 579 !test that the value of nkpt is OK, if kptopt/=0 580 !Set up dummy arrays istwfk, kpt, wtk 581 582 if(nkpt/=0 .or. dtset%kptopt/=0)then 583 ABI_ALLOCATE(istwfk,(nkpt)) 584 ABI_ALLOCATE(kpt,(3,nkpt)) 585 ABI_ALLOCATE(kpthf,(3,nkpthf)) 586 ABI_ALLOCATE(wtk,(nkpt)) 587 ! Here, occopt is also a dummy argument 588 occopt=1 ; dtset%nshiftk=1 ; dtset%kptrlatt(:,:)=0 589 590 kptrlen=20.0_dp ; wtk(:)=1.0_dp 591 dtset%shiftk(:,:)=half 592 593 nqpt=0 594 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nqpt',tread,'INT') 595 if(tread==1) nqpt=intarr(1) 596 597 chksymbreak=1 598 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'chksymbreak',tread,'INT') 599 if(tread==1) chksymbreak=intarr(1) 600 601 602 ! Use the first image to predict k and/or q points, 603 ! except if an intermediate image is available 604 intimage=1 ; if(dtset%nimage>2)intimage=(1+dtset%nimage)/2 605 606 ! Find the q-point, if any. 607 if(nqpt==1)then 608 call inqpt(chksymbreak,std_out,jdtset,lenstr,msym,natom,dtset%qptn,dtset%wtq,& 609 & dtset%rprimd_orig(1:3,1:3,intimage),dtset%spinat,string,dtset%typat,& 610 & vacuum,dtset%xred_orig(1:3,1:natom,intimage),dtset%qptrlatt) 611 end if 612 613 ! Find the k point grid 614 call inkpts(bravais,chksymbreak,dtset%fockdownsampling,iout,iscf,istwfk,jdtset,& 615 & kpt,kpthf,dtset%kptopt,kptnrm,dtset%kptrlatt_orig,dtset%kptrlatt,kptrlen,lenstr,msym,& 616 & nkpt,nkpthf,nqpt,dtset%ngkpt,dtset%nshiftk,dtset%nshiftk_orig,dtset%shiftk_orig,dtset%nsym,& 617 & occopt,dtset%qptn,response,dtset%rprimd_orig(1:3,1:3,intimage),dtset%shiftk,& 618 & string,symafm,symrel,vacuum,wtk) 619 620 ABI_DEALLOCATE(istwfk) 621 ABI_DEALLOCATE(kpt) 622 ABI_DEALLOCATE(kpthf) 623 ABI_DEALLOCATE(wtk) 624 ! nkpt and nkpthf have been computed, as well as the k point grid, if needed 625 dtset%nkpt=nkpt 626 dtset%nkpthf=nkpthf 627 end if 628 629 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nqptdm',tread,'INT') 630 if(tread==1) dtset%nqptdm=intarr(1) 631 632 if (dtset%nqptdm<-1) then 633 write(message, '(a,i12,a,a,a,a)' )& 634 & 'Input nqptdm must be >= 0, but was ',dtset%nqptdm,ch10,& 635 & 'This is not allowed.',ch10,& 636 & 'Action: check the input file.' 637 MSG_ERROR(message) 638 end if 639 640 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nzchempot',tread,'INT') 641 if(tread==1) dtset%nzchempot=intarr(1) 642 643 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'cd_customnimfrqs',tread,'INT') 644 if(tread==1) dtset%cd_customnimfrqs=intarr(1) 645 646 if (dtset%cd_customnimfrqs<0) then 647 write(message, '(a,i0,a,a,a,a)' )& 648 & 'Input cd_customnimfrqs must be >= 0, but was ',dtset%cd_customnimfrqs,ch10,& 649 & 'This is not allowed.',ch10,& 650 & 'Action: check the input file.' 651 MSG_ERROR(message) 652 end if 653 654 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_customnfreqsp',tread,'INT') 655 if(tread==1) dtset%gw_customnfreqsp=intarr(1) 656 657 if (dtset%gw_customnfreqsp<0) then 658 write(message, '(a,i0,a,a,a,a)' )& 659 & 'Input gw_customnfreqsp must be >= 0, but was ',dtset%gw_customnfreqsp,ch10,& 660 & 'This is not allowed.',ch10,& 661 & 'Action: check the input file.' 662 MSG_ERROR(message) 663 end if 664 665 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gwls_n_proj_freq',tread,'INT') 666 if(tread==1) dtset%gwls_n_proj_freq=intarr(1) 667 668 if (dtset%gwls_n_proj_freq<0) then 669 write(message, '(a,i0,a,a,a,a)' )& 670 & 'Input gwls_n_proj_freq must be >= 0, but was ',dtset%gwls_n_proj_freq,ch10,& 671 & 'This is not allowed.',ch10,& 672 & 'Action : check the input file.' 673 MSG_ERROR(message) 674 end if 675 676 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'efmas_calc_dirs',tread,'INT') 677 if(tread==1) dtset%efmas_calc_dirs=intarr(1) 678 679 if (ABS(dtset%efmas_calc_dirs)>3) then 680 write(message, '(a,i0,a,a,a,a)' )& 681 & 'Input efmas_calc_dirs must be between -3 and 3, but was ',dtset%efmas_calc_dirs,ch10,& 682 & 'This is not allowed.',ch10,& 683 & 'Action: check the input file.' 684 MSG_ERROR(message) 685 end if 686 687 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'efmas_n_dirs',tread,'INT') 688 if(tread==1) dtset%efmas_n_dirs=intarr(1) 689 690 if (dtset%efmas_n_dirs<0) then 691 write(message, '(a,i0,a,a,a,a)' )& 692 & 'Input efmas_n_dirs must be >= 0, but was ',dtset%efmas_n_dirs,ch10,& 693 & 'This is not allowed.',ch10,& 694 & 'Action: check the input file.' 695 MSG_ERROR(message) 696 end if 697 !--------------------------------------------------------------------------- 698 699 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nnos',tread,'INT') 700 if(tread==1) dtset%nnos=intarr(1) 701 702 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ga_n_rules',tread,'INT') 703 if(tread==1) dtset%ga_n_rules=intarr(1) 704 705 !Perform the first checks 706 707 leave=0 708 709 !Check that nkpt is greater than 0 710 if (nkpt<=0) then 711 write(message, '(a,a,a,a,i12,a,a,a,a)' ) ch10,& 712 & ' invars1: ERROR -',ch10,& 713 & ' After inkpts, nkpt must be > 0, but was ',nkpt,ch10,& 714 & ' This is not allowed.',ch10,& 715 & ' Action : check the input file.' 716 call wrtout(std_out, message,'COLL') 717 leave=1 718 end if 719 720 !Check that nsppol is 1 or 2 721 if (nsppol/=1 .and. nsppol/=2) then 722 write(message, '(a,a,a,a,i12,a,a,a,a)' ) ch10,& 723 & ' invars1: ERROR -',ch10,& 724 & ' Input nsppol must be 1 or 2, but was ',nsppol,ch10,& 725 & ' This is not allowed.',ch10,& 726 & ' Action : check the input file.' 727 call wrtout(std_out,message,'COLL') 728 leave=1 729 end if 730 731 !Check that nspinor is 1 or 2 732 if (nspinor/=1 .and. nspinor/=2) then 733 write(message, '(a,a,a,a,i12,a,a,a,a)' ) ch10,& 734 & ' invars1: ERROR -',ch10,& 735 & ' Input nspinor must be 1 or 2, but was ',nspinor,ch10,& 736 & ' This is not allowed.',ch10,& 737 & ' Action : check the input file.' 738 call wrtout(std_out,message,'COLL') 739 leave=1 740 end if 741 742 !Check that nspinor and nsppol are not 2 together 743 if (nsppol==2 .and. nspinor==2) then 744 write(message, '(8a)' ) ch10,& 745 & ' invars1: ERROR -',ch10,& 746 & ' nspinor and nsappol cannot be 2 together !',ch10,& 747 & ' This is not allowed.',ch10,& 748 & ' Action : check the input file.' 749 call wrtout(std_out,message,'COLL') 750 leave=1 751 end if 752 753 !Here, leave if an error has been detected earlier 754 if(leave==1) then 755 message = ' Other errors might be present in the input file. ' 756 MSG_ERROR(message) 757 end if 758 759 760 !Now, take care of mband_upper 761 762 mband_upper=1 763 occopt=1 764 fband=0.5_dp 765 766 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'occopt',tread,'INT') 767 if(tread==1) occopt=intarr(1) 768 769 !Also read fband, that is an alternative to nband. The default 770 !is different for occopt==1 and for metallic occupations. 771 if(occopt==1)fband=0.125_dp 772 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fband',tfband,'DPR') 773 if(tfband==1)fband=dprarr(1) 774 775 !fband cannot be used when occopt==0 or occopt==2 776 if(tfband==1 .and. (occopt==0 .or. occopt==2) )then 777 write(message, '(3a)' )& 778 & 'fband cannot be used if occopt==0 or occopt==2 ',ch10,& 779 & 'Action: correct your input file, suppress fband, or change occopt.' 780 MSG_ERROR(message) 781 end if 782 783 ABI_ALLOCATE(nband,(nkpt*nsppol)) 784 tnband=0 785 786 !Compute ziontypat 787 !When the pseudo-atom is pure, simple copy 788 if(ntyppure>0)then 789 do itypat=1,ntyppure 790 dtset%ziontypat(itypat)=zionpsp(itypat) 791 end do 792 end if 793 !When the pseudo-atom is alchemical, must make mixing 794 if(ntypalch>0)then 795 do itypat=ntyppure+1,ntypat 796 dtset%ziontypat(itypat)=zero 797 do ipsp=ntyppure+1,npsp 798 dtset%ziontypat(itypat)=dtset%ziontypat(itypat) & 799 & +dtset%mixalch_orig(ipsp-ntyppure,itypat-ntyppure,1)*zionpsp(ipsp) 800 end do 801 end do 802 end if 803 804 805 806 if (occopt==0 .or. occopt==1 .or. (occopt>=3 .and. occopt<=8) ) then 807 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nband',tnband,'INT') 808 ! Note : mband_upper is initialized, not nband 809 if(tnband==1) mband_upper=intarr(1) 810 811 if(tfband==1 .and. tnband==1)then 812 write(message, '(3a)' )& 813 & 'fband and nband cannot be used together. ',ch10,& 814 & 'Action: correct your input file, suppress either fband or nband.' 815 MSG_ERROR(message) 816 end if 817 818 ! In case nband was not read, use fband, either read, or the default, 819 ! to provide an upper limit for mband_upper 820 if(tnband==0)then 821 822 charge=0.0_dp 823 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'charge',tread,'DPR') 824 if(tread==1) charge=dprarr(1) 825 826 ! Only take into account negative charge, to compute maximum number of bands 827 if(charge > 0.0_dp)charge=0.0_dp 828 829 ! mband_upper=nspinor*((nint(zion_max)*natom+1)/2 - floor(charge/2.0_dp)& 830 !& + ceiling(fband*natom-1.0d-10)) 831 zval=0.0_dp 832 do iatom=1,natom 833 zval=zval+dtset%ziontypat(dtset%typat(iatom)) 834 end do 835 zelect=zval-charge 836 mband_upper=nspinor * ((ceiling(zelect-1.0d-10)+1)/2 + ceiling( fband*natom - 1.0d-10 )) 837 838 nband(:)=mband_upper 839 840 ! DEBUG 841 ! write(std_out,*)' invars1 : zion_max,natom,fband,mband_upper ' 842 ! write(std_out,*)zion_max,natom,fband,mband_upper 843 ! ENDDEBUG 844 845 end if 846 847 nband(:)=mband_upper 848 849 else if (occopt==2) then 850 ABI_ALLOCATE(reaalloc,(nkpt*nsppol)) 851 call intagm(reaalloc,nband,jdtset,nkpt*nsppol,nkpt*nsppol,string(1:lenstr),'nband',tnband,'INT') 852 if(tnband==1)then 853 do ikpt=1,nkpt*nsppol 854 if (nband(ikpt)>mband_upper) mband_upper=nband(ikpt) 855 end do 856 end if 857 ABI_DEALLOCATE(reaalloc) 858 else 859 write(message, '(a,i0,3a)' )& 860 & 'occopt=',occopt,' is not an allowed value.',ch10,& 861 & 'Action: correct your input file.' 862 MSG_ERROR(message) 863 end if 864 865 !Check that mband_upper is greater than 0 866 if (mband_upper<=0) then 867 write(message, '(a,i0,4a)' )& 868 & 'Maximal nband must be > 0, but was ',mband_upper,ch10,& 869 & 'This is not allowed.',ch10,& 870 & 'Action: check the input file.' 871 MSG_ERROR(message) 872 end if 873 874 !The following 3 values are needed to dimension the parallelism over images 875 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'imgmov',tread,'INT') 876 if(tread==1) dtset%imgmov=intarr(1) 877 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ntimimage',tread,'INT') 878 if(tread==1) dtset%ntimimage=intarr(1) 879 call intagm(dprarr,intarr,jdtset,marr,dtset%nimage,string(1:lenstr),'dynimage',tread,'INT') 880 if(tread==1)then 881 dtset%dynimage(1:dtset%nimage)=intarr(1:dtset%nimage) 882 else if (dtset%imgmov==2.or.dtset%imgmov==5) then 883 dtset%dynimage(1)=0;dtset%dynimage(dtset%nimage)=0 884 end if 885 dtset%ndynimage=count(dtset%dynimage(1:dtset%nimage)/=0) 886 887 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'wfoptalg',tread,'INT') 888 if(tread==1) then 889 dtset%wfoptalg=intarr(1) 890 else 891 if (dtset%usepaw==0) dtset%wfoptalg=0 892 if (dtset%usepaw/=0) dtset%wfoptalg=10 893 if (dtset%optdriver==RUNL_GSTATE) then 894 if (dtset%paral_kgb/=0) dtset%wfoptalg=14 895 end if 896 end if 897 898 !--------------------------------------------------------------------------- 899 !Some PAW+DMFT keywords 900 dtset%usedmft=0 901 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usedmft',tread,'INT') 902 if(tread==1) dtset%usedmft=intarr(1) 903 904 !Some ucrpa keywords 905 dtset%ucrpa=0 906 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ucrpa',tread,'INT') 907 if(tread==1) dtset%ucrpa=intarr(1) 908 909 if (dtset%ucrpa > 0 .and. dtset%usedmft > 0) then 910 write(message, '(7a)' )& 911 & 'usedmft and ucrpa are both activated in the input file ',ch10,& 912 & 'In the following, abinit assume you are doing a ucrpa calculation and ',ch10,& 913 & 'you define Wannier functions as in DFT+DMFT calculation',ch10,& 914 & 'If instead, you want to do a full dft+dmft calculation and not only the Wannier construction, use ucrpa=0' 915 MSG_WARNING(message) 916 end if 917 918 !Some PAW+U keywords 919 dtset%usepawu=0 920 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usepawu',tread,'INT') 921 if(tread==1) dtset%usepawu=intarr(1) 922 if ( dtset%usedmft > 0 .and. dtset%usepawu >= 0 ) dtset%usepawu = 1 923 924 if (dtset%usepawu > 0 ) then 925 write(message, '(7a)' )& 926 & 'usedmft and usepawu are both activated ',ch10,& 927 & 'This is not an usual calculation:',ch10,& 928 & 'usepawu will be put to a value >= 10:',ch10,& 929 & 'LDA+U potential and energy will be put to zero' 930 MSG_WARNING(message) 931 end if 932 933 dtset%usedmatpu=0 934 dtset%lpawu(1:dtset%ntypat)=-1 935 if (dtset%usepawu>0.or.dtset%usedmft>0) then 936 937 call intagm(dprarr,intarr,jdtset,marr,dtset%ntypat,string(1:lenstr),'lpawu',tread,'INT') 938 if(tread==1) dtset%lpawu(1:dtset%ntypat)=intarr(1:dtset%ntypat) 939 940 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usedmatpu',tread,'INT') 941 if(tread==1) dtset%usedmatpu=intarr(1) 942 943 end if 944 945 !Some PAW+Exact exchange keywords 946 dtset%useexexch=0 947 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'useexexch',tread,'INT') 948 if(tread==1) dtset%useexexch=intarr(1) 949 950 951 dtset%lexexch(1:dtset%ntypat)=-1 952 953 if (dtset%useexexch>0) then 954 call intagm(dprarr,intarr,jdtset,marr,dtset%ntypat,string(1:lenstr),'lexexch',tread,'INT') 955 if(tread==1) dtset%lexexch(1:dtset%ntypat)=intarr(1:dtset%ntypat) 956 end if 957 ! LDA minus half keyword 958 call intagm(dprarr,intarr,jdtset,marr,dtset%ntypat,string(1:lenstr),'ldaminushalf',tread,'INT') 959 if(tread==1) dtset%ldaminushalf(1:dtset%ntypat)=intarr(1:dtset%ntypat) 960 !Some plowan data 961 dtset%plowan_natom=0 962 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'plowan_natom',tread,'INT') 963 if(tread==1) dtset%plowan_natom=intarr(1) 964 965 dtset%plowan_nt=0 966 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'plowan_nt',tread,'INT') 967 if(tread==1) dtset%plowan_natom=intarr(1) 968 969 970 !PAW potential zero keyword 971 dtset%usepotzero=0 972 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usepotzero',tread,'INT') 973 if(tread==1) dtset%usepotzero=intarr(1) 974 975 !Macro_uj (determination of U in PAW+U), governs also allocation of atvshift 976 dtset%macro_uj = 0 977 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'macro_uj',tread,'INT') 978 if(tread==1) dtset%macro_uj=intarr(1) 979 980 ABI_DEALLOCATE(nband) 981 ABI_DEALLOCATE(ratsph) 982 ABI_DEALLOCATE(intarr) 983 ABI_DEALLOCATE(dprarr) 984 985 end subroutine invars1