The Elk Code
readinput.f90
Go to the documentation of this file.
1 
2 ! Copyright (C) 2002-2008 J. K. Dewhurst, S. Sharma and C. Ambrosch-Draxl.
3 ! This file is distributed under the terms of the GNU General Public License.
4 ! See the file COPYING for license details.
5 
6 !BOP
7 ! !ROUTINE: readinput
8 ! !INTERFACE:
9 subroutine readinput
10 ! !USES:
11 use modmain
12 use moddftu
13 use modrdm
14 use modphonon
15 use modtest
16 use modrandom
17 use modpw
18 use modtddft
19 use modulr
20 use modvars
21 use modgw
22 use modbog
23 use modw90
24 use modtdhfc
25 use modmpi
26 use modomp
27 use modramdisk
28 ! !DESCRIPTION:
29 ! Reads in the input parameters from the file {\tt elk.in}. Also sets default
30 ! values for the input parameters.
31 !
32 ! !REVISION HISTORY:
33 ! Created September 2002 (JKD)
34 !EOP
35 !BOC
36 implicit none
37 ! local variables
38 logical lv
39 integer is,ia,ias,ios
40 integer i,j,k,l
41 real(8) sc,sc1,sc2,sc3
42 real(8) scx,scy,scz
43 real(8) scu,scu1,scu2,scu3
44 real(8) solscf,zn
45 real(8) axang(4),rot(3,3)
46 real(8) rndavec,v(3),t1
47 character(256) block,symb,str
48 
49 !------------------------!
50 ! default values !
51 !------------------------!
52 ntasks=0
53 avec(:,:)=0.d0
54 avec(1,1)=1.d0
55 avec(2,2)=1.d0
56 avec(3,3)=1.d0
57 davec(:,:)=0.d0
58 sc=1.d0
59 sc1=1.d0
60 sc2=1.d0
61 sc3=1.d0
62 scx=1.d0
63 scy=1.d0
64 scz=1.d0
65 epslat=1.d-6
66 primcell=.false.
67 tshift=.true.
68 ngridk(:)=1
69 dngridk(:)=0
70 vkloff(:)=0.d0
71 autokpt=.false.
72 radkpt=40.d0
73 reducek=1
74 ngridq(:)=-1
75 reduceq=1
76 rgkmax=7.d0
77 drgkmax=0.d0
78 gmaxvr=12.d0
79 dgmaxvr=0.d0
80 lmaxapw=8
81 dlmaxapw=0
82 lmaxo=6
83 dlmaxo=0
84 lmaxi=1
85 fracinr=0.01d0
86 trhonorm=.true.
87 xctype(1)=3
88 xctype(2:3)=0
89 xctsp(1)=3
90 xctsp(2:3)=0
91 ktype(1)=52
92 ktype(2:3)=0
93 stype=3
94 swidth=0.001d0
95 autoswidth=.false.
96 mstar=10.d0
97 epsocc=1.d-10
98 epschg=1.d-3
99 nempty0=4.d0
100 dnempty0=0.d0
101 maxscl=200
102 mixtype=3
103 mixsave=.false.
104 amixpm(1)=0.05d0
105 amixpm(2)=1.d0
106 ! Broyden parameters recommended by M. Meinert
107 mixsdb=5
108 broydpm(1)=0.4d0
109 broydpm(2)=0.15d0
110 mixrho=.false.
111 epspot=1.d-6
112 epsengy=1.d-4
113 epsforce=5.d-3
114 epsstress=2.d-3
115 molecule=.false.
116 nspecies=0
117 natoms(:)=0
118 atposl(:,:,:)=0.d0
119 datposl(:,:,:)=0.d0
120 atposc(:,:,:)=0.d0
121 bfcmt0(:,:,:)=0.d0
122 sppath=''
123 scrpath=''
124 nvp1d=2
125 if (allocated(vvlp1d)) deallocate(vvlp1d)
126 allocate(vvlp1d(3,nvp1d))
127 vvlp1d(:,1)=0.d0
128 vvlp1d(:,2)=1.d0
129 npp1d=200
130 ip01d=1
131 vclp2d(:,:)=0.d0
132 vclp2d(1,1)=1.d0
133 vclp2d(2,2)=1.d0
134 np2d(:)=40
135 vclp3d(:,:)=0.d0
136 vclp3d(1,1)=1.d0
137 vclp3d(2,2)=1.d0
138 vclp3d(3,3)=1.d0
139 np3d(:)=20
140 nwplot=500
141 ngrkf=100
142 nswplot=1
143 wplot(1)=-0.5d0
144 wplot(2)=0.5d0
145 dosocc=.false.
146 tpdos=.true.
147 dosmsum=.false.
148 dosssum=.false.
149 lmirep=.true.
150 spinpol=.false.
151 spinorb=.false.
152 socscf=1.d0
153 bforb=.false.
154 bfdmag=.false.
155 atpopt=1
156 maxatpstp=200
157 tau0atp=0.2d0
158 deltast=0.005d0
159 avecref(:,:)=0.d0
160 latvopt=0
161 maxlatvstp=30
162 tau0latv=0.2d0
163 lradstp=4
164 chgexs=0.d0
165 dchgexs=0.d0
166 scissor=0.d0
167 noptcomp=1
168 ! list of all optical tensor components
169 do k=1,3; do j=1,3; do i=1,3
170  l=(k-1)*9+(j-1)*3+i
171  optcomp(:,l)=[i,j,k]
172 end do; end do; end do
173 optcomp(:,1)=1
174 intraband=.false.
175 epsband=1.d-12
176 demaxbnd=2.5d0
177 autolinengy=.false.
178 dlefe=-0.1d0
179 autodlefe=.true.
180 deapw=0.2d0
181 delorb=0.05d0
182 bfieldc0(:)=0.d0
183 dbfieldc0(:)=0.d0
184 efieldc(:)=0.d0
185 dmaxefc=1.d6
186 afieldc(:)=0.d0
187 dafieldc(:)=0.d0
188 afspc(:,:)=0.d0
189 dafspc(:,:)=0.d0
190 fsmtype=0
191 momfix(:)=0.d0
192 momfixm=0.d0
193 dmomfix(:)=0.d0
194 mommtfix(:,:,:)=1.d6
195 mommtfixm(:,:)=-1.d0
196 taufsm=0.01d0
197 rmtdelta=0.05d0
198 isgkmax=-1
199 symtype=1
200 deltaph=0.01d0
201 nphwrt=1
202 if (allocated(vqlwrt)) deallocate(vqlwrt)
203 allocate(vqlwrt(3,nphwrt))
204 vqlwrt(:,:)=0.d0
205 notelns=0
206 tforce=.false.
207 maxitoep=400
208 tau0oep=0.1d0
209 nkstlist=1
210 kstlist(:,1)=1
211 vklem(:)=0.d0
212 deltaem=0.025d0
213 ndspem=1
214 nosource=.false.
215 spinsprl=.false.
216 ssdph=.true.
217 vqlss(:)=0.d0
218 dvqlss(:)=0.d0
219 nwrite=0
220 dftu=0
221 inpdftu=1
222 ndftu=0
223 ujdu(:,:)=0.d0
224 fdu(:,:)=0.d0
225 edu(:,:)=0.d0
226 lamdu(:)=0.d0
227 udufix(:)=0.d0
228 dudufix(:)=0.d0
229 tmwrite=.false.
230 rdmxctype=2
231 rdmmaxscl=2
232 maxitn=200
233 maxitc=0
234 taurdmn=0.5d0
235 taurdmc=0.25d0
236 rdmalpha=0.656d0
237 rdmtemp=0.d0
238 reducebf=1.d0
239 ptnucl=.true.
240 tefvr=.true.
241 tefvs=.false.
242 mefvs=-2
243 vecql(:)=0.d0
244 mustar=0.15d0
245 sqaxis(1:2)=0.d0
246 sqaxis(3)=1.d0
247 test=.false.
248 spincore=.false.
249 solscf=1.d0
250 emaxelnes=-1.2d0
251 wsfac(1)=-1.1d6; wsfac(2)=1.1d6
252 vhmat(:,:)=0.d0
253 vhmat(1,1)=1.d0
254 vhmat(2,2)=1.d0
255 vhmat(3,3)=1.d0
256 reduceh=.true.
257 hybrid0=.false.
258 hybridc=1.d0
259 ecvcut=-3.5d0
260 esccut=-0.4d0
261 gmaxrf=3.d0
262 mbwgrf=-1
263 emaxrf=1.d6
264 ntemp=40
265 nvbse0=2
266 ncbse0=3
267 nvxbse=0
268 ncxbse=0
269 bsefull=.false.
270 hxbse=.true.
271 hdbse=.true.
272 fxctype=-1
273 fxclrc(1)=0.d0
274 fxclrc(2)=0.d0
275 rndatposc=0.d0
276 rndbfcmt=0.d0
277 rndavec=0.d0
278 c_tb09=0.d0
279 tc_tb09=.false.
280 hmaxvr=20.d0
281 hkmax=12.d0
282 lorbcnd=.false.
283 lorbordc=3
284 nrmtscf=1.d0
285 dnrmtscf=0.d0
286 lmaxdb=3
287 epsdev=0.0025d0
288 npmae0=-1
289 wrtvars=.false.
290 ftmtype=0
291 ntmfix=0
292 tauftm=0.1d0
293 cmagz=.false.
294 axang(:)=0.d0
295 dncgga=1.d-8
296 tstime=1000.d0
297 dtimes=0.1d0
298 npulse=0
299 nramp=0
300 nstep=0
301 ntswrite(1)=500
302 ntswrite(2)=1
303 nxoapwlo=0
304 nxlo=0
305 tdrho1d=.false.
306 tdrho2d=.false.
307 tdrho3d=.false.
308 tdmag1d=.false.
309 tdmag2d=.false.
310 tdmag3d=.false.
311 tdjr1d=.false.
312 tdjr2d=.false.
313 tdjr3d=.false.
314 tddos=.false.
315 tdlsj=.false.
316 tdjtk=.false.
317 tdxrmk=.false.
318 rndevt0=0.d0
319 sxcscf=1.d0
320 dsxcscf=0.d0
321 avecu(:,:)=0.d0
322 avecu(1,1)=1.d0
323 avecu(2,2)=1.d0
324 avecu(3,3)=1.d0
325 scu=1.d0
326 scu1=1.d0
327 scu2=1.d0
328 scu3=1.d0
329 q0cut=0.d0
330 ngridkpa(:)=-1
331 rndbfcu=0.d0
332 bfieldcu(:)=0.d0
333 efieldcu(:)=0.d0
334 tplotq0=.true.
335 trdvclr=.false.
336 trdbfcr=.false.
337 wmaxgw=-10.d0
338 tsediag=.false.
339 actype=10
340 npole=3
341 nspade=100
342 tfav0=.true.
343 rmtscf=1.d0
344 mrmtav=0
345 rmtall=-1.d0
346 maxthd=0
347 maxthd1=0
348 maxthdmkl=0
349 maxlvl=4
350 tdphi=0.d0
351 thetamld=45.d0*pi/180.d0
352 ntsbackup=0
353 ! Wannier90 variables
354 seedname='wannier'
355 num_wann=0
356 num_bands=0
357 projw90=.false.
358 lprojw90(:)=3
359 num_iter=500
360 dis_num_iter=500
361 trial_step=1.d-3
362 nxlwin=0
363 wrtunk=.false.
364 tbdip=.false.
365 tjr=.false.
366 tauefm=0.01d0
367 epsefm=1.d-6
368 ehfb=1.d0
369 tafindt=.false.
370 afindpm(:)=0.d0
371 afindpm(2)=1.d0
372 nkspolar=4
373 ntsforce=100
374 wphcut=1.d-6
375 ephscf(1)=8.d0
376 ephscf(2)=0.02d0
377 anomalous=.false.
378 tephde=.false.
379 bdiag=.false.
380 ecutb=0.001d0
381 ediag=.false.
382 pwxpsn=2
383 ramdisk=.true.
384 wrtdisk=.true.
385 epsdmat=1.d-8
386 tm3vdl=.false.
387 batch=.false.
388 tafspt=.false.
389 tbaspat=.false.
390 trdatdv=.false.
391 atdfc=0.d0
392 maxforce=-1.d0
393 msmgmt=0
394 ntsorth=1000
395 deltabf=0.5d0
396 jtconst0=.false.
397 trmt0=.true.
398 ksgwrho=.false.
399 npfftg=4
400 npfftgc=4
401 npfftq=4
402 npfftw=4
403 tphnat=.false.
404 ecutthc=0.01d0
405 tbdipu=.false.
406 bdipscf=1.d0
407 
408 !--------------------------!
409 ! read from elk.in !
410 !--------------------------!
411 open(50,file='elk.in',status='OLD',form='FORMATTED',iostat=ios)
412 if (ios /= 0) then
413  write(*,*)
414  write(*,'("Error(readinput): error opening elk.in")')
415  write(*,*)
416  stop
417 end if
418 10 continue
419 read(50,*,end=30) block
420 ! check for a comment
421 if ((block(1:1) == '!').or.(block(1:1) == '#')) goto 10
422 select case(trim(block))
423 case('tasks')
424  do i=1,maxtasks
425  read(50,'(A)',err=20) str
426  if (trim(str) == '') then
427  if (i == 1) then
428  write(*,*)
429  write(*,'("Error(readinput): no tasks to perform")')
430  write(*,*)
431  stop
432  end if
433  ntasks=i-1
434  goto 10
435  end if
436  read(str,*,iostat=ios) tasks(i)
437  if (ios /= 0) then
438  write(*,*)
439  write(*,'("Error(readinput): error reading tasks")')
440  write(*,'("(blank line required after tasks block)")')
441  write(*,*)
442  stop
443  end if
444  end do
445  write(*,*)
446  write(*,'("Error(readinput): too many tasks")')
447  write(*,'("Adjust maxtasks in modmain and recompile code")')
448  write(*,*)
449  stop
450 case('species')
451 ! generate a species file
452  call genspecies(50)
453 case('fspecies')
454 ! generate fractional species files
455  do is=1,maxspecies
456  read(50,'(A)',err=20) str
457  if (trim(str) == '') goto 10
458  read(str,*,iostat=ios) zn,symb
459  if (ios /= 0) then
460  write(*,*)
461  write(*,'("Error(readinput): error reading fractional species")')
462  write(*,'("(blank line required after fspecies block)")')
463  write(*,*)
464  stop
465  end if
466  if (zn > 0.d0) then
467  write(*,*)
468  write(*,'("Error(readinput): fractional nuclear Z > 0 : ",G18.10)') zn
469  write(*,*)
470  stop
471  end if
472  call genfspecies(zn,symb)
473  end do
474  write(*,*)
475  write(*,'("Error(readinput): too many fractional nucleus species")')
476  write(*,*)
477  stop
478 case('avec')
479  do i=1,3
480  read(50,'(A)',err=20) str
481  read(str,*,err=20) avec(:,i)
482  read(str,*,iostat=ios) avec(:,i),davec(:,i)
483  end do
484 case('scale')
485  read(50,*,err=20) sc
486 case('scale1')
487  read(50,*,err=20) sc1
488 case('scale2')
489  read(50,*,err=20) sc2
490 case('scale3')
491  read(50,*,err=20) sc3
492 case('scalex')
493  read(50,*,err=20) scx
494 case('scaley')
495  read(50,*,err=20) scy
496 case('scalez')
497  read(50,*,err=20) scz
498 case('epslat')
499  read(50,*,err=20) epslat
500  if (epslat <= 0.d0) then
501  write(*,*)
502  write(*,'("Error(readinput): epslat <= 0 : ",G18.10)') epslat
503  write(*,*)
504  stop
505  end if
506 case('primcell')
507  read(50,*,err=20) primcell
508 case('tshift')
509  read(50,*,err=20) tshift
510 case('autokpt')
511  read(50,*,err=20) autokpt
512 case('radkpt')
513  read(50,*,err=20) radkpt
514  if (radkpt <= 0.d0) then
515  write(*,*)
516  write(*,'("Error(readinput): radkpt <= 0 : ",G18.10)') radkpt
517  write(*,*)
518  stop
519  end if
520 case('ngridk')
521  read(50,'(A)',err=20) str
522  read(str,*,err=20) ngridk(:)
523  read(str,*,iostat=ios) ngridk(:),dngridk(:)
524  if (any(ngridk(:) < 1)) then
525  write(*,*)
526  write(*,'("Error(readinput): invalid ngridk :",3(X,I0))') ngridk
527  write(*,*)
528  stop
529  end if
530  autokpt=.false.
531 case('vkloff')
532  read(50,*,err=20) vkloff(:)
533  if (any(vkloff(:) < 0.d0).or.any(vkloff(:) >= 1.d0)) then
534  write(*,*)
535  write(*,'("Error(readinput): vkloff components should be in [0,1) : ",&
536  &3G18.10)') vkloff
537  write(*,*)
538  stop
539  end if
540 case('reducek')
541  read(50,*,err=20) reducek
542 case('ngridq')
543  read(50,*,err=20) ngridq(:)
544  if (any(ngridq(:) < 1)) then
545  write(*,*)
546  write(*,'("Error(readinput): invalid ngridq :",3(X,I0))') ngridq
547  write(*,*)
548  stop
549  end if
550 case('reduceq')
551  read(50,*,err=20) reduceq
552 case('rgkmax')
553  read(50,'(A)',err=20) str
554  read(str,*,err=20) rgkmax
555  read(str,*,iostat=ios) rgkmax,drgkmax
556  if (rgkmax <= 0.d0) then
557  write(*,*)
558  write(*,'("Error(readinput): rgkmax <= 0 : ",G18.10)') rgkmax
559  write(*,*)
560  stop
561  end if
562 case('gmaxvr')
563  read(50,'(A)',err=20) str
564  read(str,*,err=20) gmaxvr
565  read(str,*,iostat=ios) gmaxvr,dgmaxvr
566 case('lmaxapw')
567  read(50,'(A)',err=20) str
568  read(str,*,err=20) lmaxapw
569  read(str,*,iostat=ios) lmaxapw,dlmaxapw
570  if (lmaxapw < 0) then
571  write(*,*)
572  write(*,'("Error(readinput): lmaxapw < 0 : ",I0)') lmaxapw
573  write(*,*)
574  stop
575  end if
576  if (lmaxapw >= maxlapw) then
577  write(*,*)
578  write(*,'("Error(readinput): lmaxapw too large : ",I0)') lmaxapw
579  write(*,'("Adjust maxlapw in modmain and recompile code")')
580  write(*,*)
581  stop
582  end if
583 case('lmaxo','lmaxvr')
584  read(50,'(A)',err=20) str
585  read(str,*,err=20) lmaxo
586  read(str,*,iostat=ios) lmaxo,dlmaxo
587  if (lmaxo < 3) then
588  write(*,*)
589  write(*,'("Error(readinput): lmaxo < 3 : ",I0)') lmaxo
590  write(*,*)
591  stop
592  end if
593 case('lmaxi','lmaxinr')
594  read(50,*,err=20) lmaxi
595  if (lmaxi < 1) then
596  write(*,*)
597  write(*,'("Error(readinput): lmaxi < 1 : ",I0)') lmaxi
598  write(*,*)
599  stop
600  end if
601 case('lmaxmat')
602  read(50,*,err=20)
603  write(*,'("Info(readinput): variable ''lmaxmat'' is no longer used")')
604 case('fracinr')
605  read(50,*,err=20) fracinr
606 case('trhonorm')
607  read(50,*,err=20) trhonorm
608 case('spinpol')
609  read(50,*,err=20) spinpol
610 case('spinorb')
611  read(50,*,err=20) spinorb
612 case('socscf')
613  read(50,*,err=20) socscf
614  if (socscf < 0.d0) then
615  write(*,*)
616  write(*,'("Error(readinput): socscf < 0 : ",G18.10)') socscf
617  write(*,*)
618  stop
619  end if
620 case('bforb')
621  read(50,*,err=20) bforb
622 case('bfdmag')
623  read(50,*,err=20) bfdmag
624 case('xctype')
625  read(50,'(A)',err=20) str
626  str=trim(str)//' 0 0'
627  read(str,*,err=20) xctype(:)
628 case('xctsp')
629  read(50,'(A)',err=20) str
630  str=trim(str)//' 0 0'
631  read(str,*,err=20) xctsp(:)
632 case('ktype')
633  read(50,'(A)',err=20) str
634  str=trim(str)//' 0 0'
635  read(str,*,err=20) ktype(:)
636  if (ktype(3) /= 0) then
637  write(*,*)
638  write(*,'("Error(readinput): ktype(3) should be zero : ",I0)') ktype(3)
639  write(*,*)
640  stop
641  end if
642 case('stype')
643  read(50,*,err=20) stype
644 case('swidth')
645  read(50,*,err=20) swidth
646  if (swidth < 1.d-9) then
647  write(*,*)
648  write(*,'("Error(readinput): swidth too small or negative : ",G18.10)') &
649  swidth
650  write(*,*)
651  stop
652  end if
653 case('autoswidth')
654  read(50,*,err=20) autoswidth
655 case('mstar')
656  read(50,*,err=20) mstar
657  if (mstar <= 0.d0) then
658  write(*,*)
659  write(*,'("Error(readinput): mstar <= 0 : ",G18.10)') mstar
660  write(*,*)
661  stop
662  end if
663 case('epsocc')
664  read(50,*,err=20) epsocc
665  if (epsocc <= 0.d0) then
666  write(*,*)
667  write(*,'("Error(readinput): epsocc <= 0 : ",G18.10)') epsocc
668  write(*,*)
669  stop
670  end if
671 case('epschg')
672  read(50,*,err=20) epschg
673  if (epschg <= 0.d0) then
674  write(*,*)
675  write(*,'("Error(readinput): epschg <= 0 : ",G18.10)') epschg
676  write(*,*)
677  stop
678  end if
679 case('nempty','nempty0')
680  read(50,'(A)',err=20) str
681  read(str,*,err=20) nempty0
682  read(str,*,iostat=ios) nempty0,dnempty0
683  if (nempty0 <= 0.d0) then
684  write(*,*)
685  write(*,'("Error(readinput): nempty <= 0 : ",G18.10)') nempty0
686  write(*,*)
687  stop
688  end if
689 case('mixtype')
690  read(50,*,err=20) mixtype
691 case('mixsave')
692  read(50,*,err=20) mixsave
693 case('amixpm','beta0','betamax')
694  if (trim(block) == 'amixpm') then
695  read(50,*,err=20) amixpm(:)
696  else if (trim(block) == 'beta0') then
697  read(50,*,err=20) amixpm(1)
698  else
699  read(50,*,err=20) amixpm(2)
700  end if
701  if (amixpm(1) < 0.d0) then
702  write(*,*)
703  write(*,'("Error(readinput): beta0 [amixpm(1)] < 0 : ",G18.10)') amixpm(1)
704  write(*,*)
705  stop
706  end if
707  if ((amixpm(2) < 0.d0).or.(amixpm(2) > 1.d0)) then
708  write(*,*)
709  write(*,'("Error(readinput): betamax [amixpm(2)] not in [0,1] : ",G18.10)')&
710  amixpm(2)
711  write(*,*)
712  stop
713  end if
714 case('mixsdb')
715  read(50,*,err=20) mixsdb
716  if (mixsdb < 2) then
717  write(*,*)
718  write(*,'("Error(readinput): mixsdb < 2 : ",I0)') mixsdb
719  write(*,*)
720  stop
721  end if
722 case('broydpm')
723  read(50,*,err=20) broydpm(:)
724  if ((broydpm(1) < 0.d0).or.(broydpm(1) > 1.d0).or. &
725  (broydpm(2) < 0.d0).or.(broydpm(2) > 1.d0)) then
726  write(*,*)
727  write(*,'("Error(readinput): invalid Broyden mixing parameters : ",&
728  &2G18.10)') broydpm
729  write(*,*)
730  stop
731  end if
732 case('mixrho')
733  read(50,*,err=20) mixrho
734 case('maxscl')
735  read(50,*,err=20) maxscl
736  if (maxscl < 0) then
737  write(*,*)
738  write(*,'("Error(readinput): maxscl < 0 : ",I0)') maxscl
739  write(*,*)
740  stop
741  end if
742 case('epspot')
743  read(50,*,err=20) epspot
744 case('epsengy')
745  read(50,*,err=20) epsengy
746 case('epsforce')
747  read(50,*,err=20) epsforce
748 case('epsstress')
749  read(50,*,err=20) epsstress
750 case('sppath')
751  read(50,*,err=20) sppath
752  sppath=adjustl(sppath)
753 case('scrpath')
754  read(50,*,err=20) scrpath
755 case('molecule')
756  read(50,*,err=20) molecule
757 case('atoms')
758  read(50,*,err=20) nspecies
759  if (nspecies < 1) then
760  write(*,*)
761  write(*,'("Error(readinput): nspecies < 1 : ",I0)') nspecies
762  write(*,*)
763  stop
764  end if
765  if (nspecies > maxspecies) then
766  write(*,*)
767  write(*,'("Error(readinput): nspecies too large : ",I0)') nspecies
768  write(*,'("Adjust maxspecies in modmain and recompile code")')
769  write(*,*)
770  stop
771  end if
772  do is=1,nspecies
773  read(50,*,err=20) spfname(is)
774  spfname(is)=adjustl(spfname(is))
775  read(50,*,err=20) natoms(is)
776  if (natoms(is) < 1) then
777  write(*,*)
778  write(*,'("Error(readinput): natoms < 1 : ",I0)') natoms(is)
779  write(*,'(" for species ",I0)') is
780  write(*,*)
781  stop
782  end if
783  if (natoms(is) > maxatoms) then
784  write(*,*)
785  write(*,'("Error(readinput): natoms too large : ",I0)') natoms(is)
786  write(*,'(" for species ",I0)') is
787  write(*,'("Adjust maxatoms in modmain and recompile code")')
788  write(*,*)
789  stop
790  end if
791  do ia=1,natoms(is)
792  read(50,'(A)',err=20) str
793  read(str,*,err=20) atposl(:,ia,is)
794  read(str,*,iostat=ios) atposl(:,ia,is),bfcmt0(:,ia,is),datposl(:,ia,is)
795  end do
796  end do
797 case('plot1d')
798  read(50,*,err=20) nvp1d,npp1d
799  if (nvp1d < 1) then
800  write(*,*)
801  write(*,'("Error(readinput): nvp1d < 1 : ",I0)') nvp1d
802  write(*,*)
803  stop
804  end if
805  if (npp1d < nvp1d) then
806  write(*,*)
807  write(*,'("Error(readinput): npp1d < nvp1d :",2(X,I0))') npp1d,nvp1d
808  write(*,*)
809  stop
810  end if
811  if (allocated(vvlp1d)) deallocate(vvlp1d)
812  allocate(vvlp1d(3,nvp1d))
813  do i=1,nvp1d
814  read(50,*,err=20) vvlp1d(:,i)
815  end do
816 case('ip01d')
817  read(50,*,err=20) ip01d
818  if (ip01d < 1) then
819  write(*,*)
820  write(*,'("Error(readinput): ip01d < 1 : ",I0)') ip01d
821  write(*,*)
822  stop
823  end if
824 case('plot2d')
825  read(50,*,err=20) vclp2d(:,0)
826  read(50,*,err=20) vclp2d(:,1)
827  read(50,*,err=20) vclp2d(:,2)
828  read(50,*,err=20) np2d(:)
829  if ((np2d(1) < 1).or.(np2d(2) < 1)) then
830  write(*,*)
831  write(*,'("Error(readinput): np2d < 1 :",2(X,I0))') np2d
832  write(*,*)
833  stop
834  end if
835 case('plot3d')
836  read(50,*,err=20) vclp3d(:,0)
837  read(50,*,err=20) vclp3d(:,1)
838  read(50,*,err=20) vclp3d(:,2)
839  read(50,*,err=20) vclp3d(:,3)
840  read(50,*,err=20) np3d(:)
841  if ((np3d(1) < 1).or.(np3d(2) < 1).or.(np3d(3) < 1)) then
842  write(*,*)
843  write(*,'("Error(readinput): np3d < 1 :",3(X,I0))') np3d
844  write(*,*)
845  stop
846  end if
847 case('wplot','dos')
848  read(50,*,err=20) nwplot,ngrkf,nswplot
849  if (nwplot < 2) then
850  write(*,*)
851  write(*,'("Error(readinput): nwplot < 2 : ",I0)') nwplot
852  write(*,*)
853  stop
854  end if
855  if (ngrkf < 1) then
856  write(*,*)
857  write(*,'("Error(readinput): ngrkf < 1 : ",I0)') ngrkf
858  write(*,*)
859  stop
860  end if
861  if (nswplot < 0) then
862  write(*,*)
863  write(*,'("Error(readinput): nswplot < 0 : ",I0)') nswplot
864  write(*,*)
865  stop
866  end if
867  read(50,*,err=20) wplot(:)
868  if (wplot(1) > wplot(2)) then
869  write(*,*)
870  write(*,'("Error(readinput): wplot(1) > wplot(2) : ",2G18.10)') wplot
871  write(*,*)
872  stop
873  end if
874 case('dosocc')
875  read(50,*,err=20) dosocc
876 case('tpdos')
877  read(50,*,err=20) tpdos
878 case('dosmsum')
879  read(50,*,err=20) dosmsum
880 case('dosssum')
881  read(50,*,err=20) dosssum
882 case('lmirep')
883  read(50,*,err=20) lmirep
884 case('atpopt')
885  read(50,*,err=20) atpopt
886 case('maxatpstp','maxatmstp')
887  read(50,*,err=20) maxatpstp
888  if (maxatpstp < 1) then
889  write(*,*)
890  write(*,'("Error(readinput): maxatpstp < 1 : ",I0)') maxatpstp
891  write(*,*)
892  stop
893  end if
894 case('tau0atp','tau0atm')
895  read(50,*,err=20) tau0atp
896 case('deltast')
897  read(50,*,err=20) deltast
898  if (deltast <= 0.d0) then
899  write(*,*)
900  write(*,'("Error(readinput): deltast <= 0 : ",G18.10)') deltast
901  write(*,*)
902  stop
903  end if
904 case('avecref')
905  read(50,*,err=20) avecref(:,1)
906  read(50,*,err=20) avecref(:,2)
907  read(50,*,err=20) avecref(:,3)
908 case('latvopt')
909  read(50,*,err=20) latvopt
910 case('maxlatvstp')
911  read(50,*,err=20) maxlatvstp
912  if (maxlatvstp < 1) then
913  write(*,*)
914  write(*,'("Error(readinput): maxlatvstp < 1 : ",I0)') maxlatvstp
915  write(*,*)
916  stop
917  end if
918 case('tau0latv')
919  read(50,*,err=20) tau0latv
920 case('nstfsp')
921  read(50,*,err=20)
922  write(*,'("Info(readinput): variable ''nstfsp'' is no longer used")')
923 case('lradstp')
924  read(50,*,err=20) lradstp
925  if (lradstp < 1) then
926  write(*,*)
927  write(*,'("Error(readinput): lradstp < 1 : ",I0)') lradstp
928  write(*,*)
929  stop
930  end if
931 case('chgexs')
932  read(50,'(A)',err=20) str
933  read(str,*,err=20) chgexs
934  read(str,*,iostat=ios) chgexs,dchgexs
935 case('nprad')
936  read(50,*,err=20)
937  write(*,'("Info(readinput): variable ''nprad'' is no longer used")')
938 case('scissor')
939  read(50,*,err=20) scissor
940 case('noptcomp')
941  read(50,*,err=20) noptcomp
942  if ((noptcomp < 1).or.(noptcomp > 27)) then
943  write(*,*)
944  write(*,'("Error(readinput): noptcomp should be from 1 to 27 : ",I0)') &
945  noptcomp
946  write(*,*)
947  stop
948  end if
949 case('optcomp')
950  do i=1,27
951  read(50,'(A)',err=20) str
952  if (trim(str) == '') then
953  if (i == 1) then
954  write(*,*)
955  write(*,'("Error(readinput): empty optical component list")')
956  write(*,*)
957  stop
958  end if
959  noptcomp=i-1
960  goto 10
961  end if
962  str=trim(str)//' 1 1'
963  read(str,*,iostat=ios) optcomp(:,i)
964  if (ios /= 0) then
965  write(*,*)
966  write(*,'("Error(readinput): error reading optical component list")')
967  write(*,'("(blank line required after optcomp block)")')
968  write(*,*)
969  stop
970  end if
971  if (any(optcomp(:,i) < 1).or.any(optcomp(:,i) > 3)) then
972  write(*,*)
973  write(*,'("Error(readinput): invalid optcomp :",3(X,I0))') optcomp(:,i)
974  write(*,*)
975  stop
976  end if
977  end do
978  write(*,*)
979  write(*,'("Error(readinput): optical component list too long")')
980  write(*,*)
981  stop
982 case('intraband')
983  read(50,*,err=20) intraband
984 case('evaltol')
985  read(50,*,err=20)
986  write(*,'("Info(readinput): variable ''evaltol'' is no longer used")')
987 case('deband')
988  read(50,*,err=20)
989  write(*,'("Info(readinput): variable ''deband'' is no longer used")')
990 case('epsband')
991  read(50,*,err=20) epsband
992  if (epsband <= 0.d0) then
993  write(*,*)
994  write(*,'("Error(readinput): epsband <= 0 : ",G18.10)') epsband
995  write(*,*)
996  stop
997  end if
998 case('demaxbnd')
999  read(50,*,err=20) demaxbnd
1000  if (demaxbnd <= 0.d0) then
1001  write(*,*)
1002  write(*,'("Error(readinput): demaxbnd <= 0 : ",G18.10)') demaxbnd
1003  write(*,*)
1004  stop
1005  end if
1006 case('autolinengy')
1007  read(50,*,err=20) autolinengy
1008 case('dlefe')
1009  read(50,*,err=20) dlefe
1010 case('autodlefe')
1011  read(50,*,err=20) autodlefe
1012 case('deapw')
1013  read(50,*,err=20) deapw
1014  if (abs(deapw) < 1.d-8) then
1015  write(*,*)
1016  write(*,'("Error(readinput): invalid deapw : ",G18.10)') deapw
1017  write(*,*)
1018  stop
1019  end if
1020 case('delorb')
1021  read(50,*,err=20) delorb
1022  if (abs(delorb) < 1.d-8) then
1023  write(*,*)
1024  write(*,'("Error(readinput): invalid delorb : ",G18.10)') delorb
1025  write(*,*)
1026  stop
1027  end if
1028 case('bfieldc')
1029  read(50,'(A)',err=20) str
1030  read(str,*,err=20) bfieldc0(:)
1031  read(str,*,iostat=ios) bfieldc0(:),dbfieldc0(:)
1032 case('efieldc')
1033  read(50,*,err=20) efieldc(:)
1034 case('dmaxefc')
1035  read(50,*,err=20) dmaxefc
1036  if (dmaxefc < 0) then
1037  write(*,*)
1038  write(*,'("Error(readinput): dmaxefc < 0 : ",G18.10)') dmaxefc
1039  write(*,*)
1040  stop
1041  end if
1042 case('afieldc')
1043  read(50,'(A)',err=20) str
1044  read(str,*,err=20) afieldc(:)
1045  read(str,*,iostat=ios) afieldc(:),dafieldc(:)
1046 case('afspc')
1047  do i=1,3
1048  read(50,'(A)',err=20) str
1049  read(str,*,err=20) afspc(i,:)
1050  read(str,*,iostat=ios) afspc(i,:),dafspc(i,:)
1051  end do
1052 case('fsmtype','fixspin')
1053  read(50,*,err=20) fsmtype
1054 case('momfix')
1055  read(50,'(A)',err=20) str
1056  read(str,*,err=20) momfix(:)
1057  read(str,*,iostat=ios) momfix(:),dmomfix(:)
1058 case('momfixm')
1059  read(50,*,err=20) momfixm
1060  if (momfixm < 0.d0) then
1061  write(*,*)
1062  write(*,'("Error(readinput): momfixm < 0 : ",G18.10)') momfixm
1063  write(*,*)
1064  stop
1065  end if
1066 case('mommtfix')
1067  do ias=1,maxspecies*maxatoms
1068  read(50,'(A)',err=20) str
1069  if (trim(str) == '') goto 10
1070  read(str,*,iostat=ios) is,ia,mommtfix(:,ia,is)
1071  if (ios /= 0) then
1072  write(*,*)
1073  write(*,'("Error(readinput): error reading muffin-tin fixed spin &
1074  &moments")')
1075  write(*,'("(blank line required after mommtfix block)")')
1076  write(*,*)
1077  stop
1078  end if
1079  end do
1080 case('mommtfixm')
1081  do ias=1,maxspecies*maxatoms
1082  read(50,'(A)',err=20) str
1083  if (trim(str) == '') goto 10
1084  read(str,*,iostat=ios) is,ia,mommtfixm(ia,is)
1085  if (ios /= 0) then
1086  write(*,*)
1087  write(*,'("Error(readinput): error reading muffin-tin fixed spin &
1088  &moment magnitudes")')
1089  write(*,'("(blank line required after mommtfixm block)")')
1090  write(*,*)
1091  stop
1092  end if
1093  end do
1094 case('taufsm')
1095  read(50,*,err=20) taufsm
1096  if (taufsm < 0.d0) then
1097  write(*,*)
1098  write(*,'("Error(readinput): taufsm < 0 : ",G18.10)') taufsm
1099  write(*,*)
1100  stop
1101  end if
1102 case('autormt')
1103  read(50,*,err=20)
1104  write(*,'("Info(readinput): variable ''autormt'' is no longer used")')
1105 case('rmtdelta')
1106  read(50,*,err=20) rmtdelta
1107  if (rmtdelta < 0.d0) then
1108  write(*,*)
1109  write(*,'("Warning(readinput): rmtdelta < 0 : ",G18.10)') rmtdelta
1110  end if
1111 case('isgkmax')
1112  read(50,*,err=20) isgkmax
1113 case('nosym')
1114  read(50,*,err=20) lv
1115  if (lv) symtype=0
1116 case('symtype')
1117  read(50,*,err=20) symtype
1118  if ((symtype < 0).or.(symtype > 2)) then
1119  write(*,*)
1120  write(*,'("Error(readinput): symtype not defined : ",I0)') symtype
1121  write(*,*)
1122  stop
1123  end if
1124 case('deltaph')
1125  read(50,*,err=20) deltaph
1126  if (deltaph <= 0.d0) then
1127  write(*,*)
1128  write(*,'("Error(readinput): deltaph <= 0 : ",G18.10)') deltaph
1129  write(*,*)
1130  stop
1131  end if
1132 case('phwrite')
1133  read(50,*,err=20) nphwrt
1134  if (nphwrt < 1) then
1135  write(*,*)
1136  write(*,'("Error(readinput): nphwrt < 1 : ",I0)') nphwrt
1137  write(*,*)
1138  stop
1139  end if
1140  if (allocated(vqlwrt)) deallocate(vqlwrt)
1141  allocate(vqlwrt(3,nphwrt))
1142  do i=1,nphwrt
1143  read(50,*,err=20) vqlwrt(:,i)
1144  end do
1145 case('notes')
1146  if (allocated(notes)) deallocate(notes)
1147  allocate(notes(0))
1148  notelns=0
1149  do
1150  read(50,'(A)') str
1151  if (trim(str) == '') goto 10
1152  notes=[notes(1:notelns),str]
1153  notelns=notelns+1
1154  end do
1155 case('tforce')
1156  read(50,*,err=20) tforce
1157 case('tfibs')
1158  read(50,*,err=20)
1159  write(*,'("Info(readinput): variable ''tfibs'' is no longer used")')
1160 case('maxitoep')
1161  read(50,*,err=20) maxitoep
1162  if (maxitoep < 1) then
1163  write(*,*)
1164  write(*,'("Error(readinput): maxitoep < 1 : ",I0)') maxitoep
1165  write(*,*)
1166  stop
1167  end if
1168 case('tauoep')
1169  read(50,*,err=20)
1170  write(*,'("Info(readinput): variable ''tauoep'' is no longer used")')
1171 case('tau0oep')
1172  read(50,*,err=20) tau0oep
1173  if (tau0oep < 0.d0) then
1174  write(*,*)
1175  write(*,'("Error(readinput): tau0oep < 0 : ",G18.10)') tau0oep
1176  write(*,*)
1177  stop
1178  end if
1179 case('kstlist')
1180  do i=1,maxkst
1181  read(50,'(A)',err=20) str
1182  if (trim(str) == '') then
1183  if (i == 1) then
1184  write(*,*)
1185  write(*,'("Error(readinput): empty k-point and state list")')
1186  write(*,*)
1187  stop
1188  end if
1189  nkstlist=i-1
1190  goto 10
1191  end if
1192  str=trim(str)//' 1'
1193  read(str,*,iostat=ios) kstlist(:,i)
1194  if (ios /= 0) then
1195  write(*,*)
1196  write(*,'("Error(readinput): error reading k-point and state list")')
1197  write(*,'("(blank line required after kstlist block)")')
1198  write(*,*)
1199  stop
1200  end if
1201  end do
1202  write(*,*)
1203  write(*,'("Error(readinput): k-point and state list too long")')
1204  write(*,*)
1205  stop
1206 case('vklem')
1207  read(50,*,err=20) vklem
1208 case('deltaem')
1209  read(50,*,err=20) deltaem
1210  if (deltaem <= 0.d0) then
1211  write(*,*)
1212  write(*,'("Error(readinput): deltaem <= 0 : ",G18.10)') deltaem
1213  write(*,*)
1214  stop
1215  end if
1216 case('ndspem')
1217  read(50,*,err=20) ndspem
1218  if ((ndspem < 1).or.(ndspem > 4)) then
1219  write(*,*)
1220  write(*,'("Error(readinput): ndspem out of range : ",I0)') ndspem
1221  write(*,*)
1222  stop
1223  end if
1224 case('nosource')
1225  read(50,*,err=20) nosource
1226 case('spinsprl')
1227  read(50,*,err=20) spinsprl
1228 case('ssdph')
1229  read(50,*,err=20) ssdph
1230 case('vqlss')
1231  read(50,'(A)',err=20) str
1232  read(str,*,err=20) vqlss
1233  read(str,*,iostat=ios) vqlss,dvqlss
1234 case('nwrite')
1235  read(50,*,err=20) nwrite
1236 case('DFT+U','dft+u','lda+u')
1237  read(50,*,err=20) dftu,inpdftu
1238  do i=1,maxdftu
1239  read(50,'(A)',err=20) str
1240  if (trim(str) == '') then
1241  ndftu=i-1
1242  goto 10
1243  end if
1244  select case(inpdftu)
1245  case(1)
1246  read(str,*,iostat=ios) is,l,ujdu(1:2,i)
1247  case(2)
1248  read(str,*,iostat=ios) is,l,(fdu(k,i),k=0,2*l,2)
1249  case(3)
1250  read(str,*,iostat=ios) is,l,(edu(k,i),k=0,l)
1251  case(4)
1252  read(str,*,iostat=ios) is,l,lamdu(i)
1253  case(5)
1254  read(str,*,iostat=ios) is,l,udufix(i),dudufix(i)
1255  read(str,*,iostat=ios) is,l,udufix(i)
1256  case default
1257  write(*,*)
1258  write(*,'("Error(readinput): invalid inpdftu : ",I0)') inpdftu
1259  write(*,*)
1260  stop
1261  end select
1262  if (ios /= 0) then
1263  write(*,*)
1264  write(*,'("Error(readinput): error reading DFT+U parameters")')
1265  write(*,'("(blank line required after dft+u block)")')
1266  write(*,*)
1267  stop
1268  end if
1269  if ((is < 1).or.(is >= maxspecies)) then
1270  write(*,*)
1271  write(*,'("Error(readinput): invalid species number in dft+u block : ", &
1272  &I0)') is
1273  write(*,*)
1274  stop
1275  end if
1276  if (l < 0) then
1277  write(*,*)
1278  write(*,'("Error(readinput): l < 0 in dft+u block : ",I0)') l
1279  write(*,*)
1280  stop
1281  end if
1282  if (l > lmaxdm) then
1283  write(*,*)
1284  write(*,'("Error(readinput): l > lmaxdm in dft+u block :",2(X,I0))') l, &
1285  lmaxdm
1286  write(*,*)
1287  stop
1288  end if
1289 ! check for repeated entries
1290  do j=1,i-1
1291  if ((is == isldu(1,j)).and.(l == isldu(2,j))) then
1292  write(*,*)
1293  write(*,'("Error(readinput): repeated entry in DFT+U block")')
1294  write(*,*)
1295  stop
1296  end if
1297  end do
1298  isldu(1,i)=is
1299  isldu(2,i)=l
1300  end do
1301  write(*,*)
1302  write(*,'("Error(readinput): too many DFT+U entries")')
1303  write(*,'("Adjust maxdftu in modmain and recompile code")')
1304  write(*,*)
1305  stop
1306 case('tmwrite','tmomlu')
1307  read(50,*,err=20) tmwrite
1308 case('readadu','readalu')
1309  read(50,*,err=20)
1310  write(*,'("Info(readinput): variable ''readadu'' is no longer used")')
1311 case('rdmxctype')
1312  read(50,*,err=20) rdmxctype
1313 case('rdmmaxscl')
1314  read(50,*,err=20) rdmmaxscl
1315  if (rdmmaxscl < 0) then
1316  write(*,*)
1317  write(*,'("Error(readinput): rdmmaxscl < 0 : ",I0)') rdmmaxscl
1318  write(*,*)
1319  end if
1320 case('maxitn')
1321  read(50,*,err=20) maxitn
1322 case('maxitc')
1323  read(50,*,err=20) maxitc
1324 case('taurdmn')
1325  read(50,*,err=20) taurdmn
1326  if (taurdmn < 0.d0) then
1327  write(*,*)
1328  write(*,'("Error(readinput): taurdmn < 0 : ",G18.10)') taurdmn
1329  write(*,*)
1330  stop
1331  end if
1332 case('taurdmc')
1333  read(50,*,err=20) taurdmc
1334  if (taurdmc < 0.d0) then
1335  write(*,*)
1336  write(*,'("Error(readinput): taurdmc < 0 : ",G18.10)') taurdmc
1337  write(*,*)
1338  stop
1339  end if
1340 case('rdmalpha')
1341  read(50,*,err=20) rdmalpha
1342  if ((rdmalpha <= 0.d0).or.(rdmalpha >= 1.d0)) then
1343  write(*,*)
1344  write(*,'("Error(readinput): rdmalpha not in (0,1) : ",G18.10)') rdmalpha
1345  write(*,*)
1346  stop
1347  end if
1348 case('rdmtemp')
1349  read(50,*,err=20) rdmtemp
1350  if (rdmtemp < 0.d0) then
1351  write(*,*)
1352  write(*,'("Error(readinput): rdmtemp < 0 : ",G18.10)') rdmtemp
1353  write(*,*)
1354  stop
1355  end if
1356 case('reducebf')
1357  read(50,*,err=20) reducebf
1358  if ((reducebf < 0.5d0).or.(reducebf > 1.d0)) then
1359  write(*,*)
1360  write(*,'("Error(readinput): reducebf not in [0.5,1] : ",G18.10)') reducebf
1361  write(*,*)
1362  stop
1363  end if
1364 case('ptnucl')
1365  read(50,*,err=20) ptnucl
1366 case('tefvr','tseqr')
1367  read(50,*,err=20) tefvr
1368 case('tefvs')
1369  read(50,*,err=20) tefvs
1370 case('mefvs')
1371  read(50,*,err=20) mefvs
1372 case('tefvit','tseqit')
1373  read(50,*,err=20)
1374  write(*,'("Info(readinput): variable ''tefvit'' is no longer used")')
1375 case('minitefv','minseqit')
1376  read(50,*,err=20)
1377  write(*,'("Info(readinput): variable ''minitefv'' is no longer used")')
1378 case('nefvit','maxitefv','maxseqit','nseqit')
1379  read(50,*,err=20)
1380  write(*,'("Info(readinput): variable ''nefvit'' is no longer used")')
1381 case('befvit','bseqit')
1382  read(50,*,err=20)
1383  write(*,'("Info(readinput): variable ''befvit'' is no longer used")')
1384 case('epsefvit','epsseqit')
1385  read(50,*,err=20)
1386  write(*,'("Info(readinput): variable ''epsefvit'' is no longer used")')
1387 case('tauseq')
1388  read(50,*,err=20)
1389  write(*,'("Info(readinput): variable ''tauseq'' is no longer used")')
1390 case('vecql')
1391  read(50,*,err=20) vecql(:)
1392 case('mustar')
1393  read(50,*,err=20) mustar
1394 case('sqaxis','sqados')
1395  read(50,*,err=20) sqaxis(:)
1396 case('test')
1397  read(50,*,err=20) test
1398 case('frozencr')
1399  read(50,*,err=20)
1400  write(*,'("Info(readinput): variable ''frozencr'' is no longer used")')
1401 case('spincore')
1402  read(50,*,err=20) spincore
1403 case('solscf')
1404  read(50,*,err=20) solscf
1405  if (solscf < 0.d0) then
1406  write(*,*)
1407  write(*,'("Error(readinput): solscf < 0 : ",G18.10)') solscf
1408  write(*,*)
1409  stop
1410  end if
1411 case('emaxelnes')
1412  read(50,*,err=20) emaxelnes
1413 case('wsfac')
1414  read(50,*,err=20) wsfac(:)
1415 case('vhmat')
1416  read(50,*,err=20) vhmat(1,:)
1417  read(50,*,err=20) vhmat(2,:)
1418  read(50,*,err=20) vhmat(3,:)
1419 case('reduceh')
1420  read(50,*,err=20) reduceh
1421 case('hybrid')
1422  read(50,*,err=20) hybrid0
1423 case('hybridc','hybmix')
1424  read(50,*,err=20) hybridc
1425  if ((hybridc < 0.d0).or.(hybridc > 1.d0)) then
1426  write(*,*)
1427  write(*,'("Error(readinput): invalid hybridc : ",G18.10)') hybridc
1428  write(*,*)
1429  stop
1430  end if
1431 case('ecvcut')
1432  read(50,*,err=20) ecvcut
1433 case('esccut')
1434  read(50,*,err=20) esccut
1435 case('nvbse')
1436  read(50,*,err=20) nvbse0
1437  if (nvbse0 < 0) then
1438  write(*,*)
1439  write(*,'("Error(readinput): nvbse < 0 : ",I0)') nvbse0
1440  write(*,*)
1441  stop
1442  end if
1443 case('ncbse')
1444  read(50,*,err=20) ncbse0
1445  if (ncbse0 < 0) then
1446  write(*,*)
1447  write(*,'("Error(readinput): ncbse < 0 : ",I0)') ncbse0
1448  write(*,*)
1449  stop
1450  end if
1451 case('istxbse')
1452  do i=1,maxxbse
1453  read(50,'(A)',err=20) str
1454  if (trim(str) == '') then
1455  if (i == 1) then
1456  write(*,*)
1457  write(*,'("Error(readinput): empty BSE extra valence state list")')
1458  write(*,*)
1459  stop
1460  end if
1461  nvxbse=i-1
1462  goto 10
1463  end if
1464  read(str,*,iostat=ios) istxbse(i)
1465  if (ios /= 0) then
1466  write(*,*)
1467  write(*,'("Error(readinput): error reading BSE valence state list")')
1468  write(*,'("(blank line required after istxbse block)")')
1469  write(*,*)
1470  stop
1471  end if
1472  end do
1473  write(*,*)
1474  write(*,'("Error(readinput): BSE extra valence state list too long")')
1475  write(*,*)
1476  stop
1477 case('jstxbse')
1478  do i=1,maxxbse
1479  read(50,'(A)',err=20) str
1480  if (trim(str) == '') then
1481  if (i == 1) then
1482  write(*,*)
1483  write(*,'("Error(readinput): empty BSE extra conduction state list")')
1484  write(*,*)
1485  stop
1486  end if
1487  ncxbse=i-1
1488  goto 10
1489  end if
1490  read(str,*,iostat=ios) jstxbse(i)
1491  if (ios /= 0) then
1492  write(*,*)
1493  write(*,'("Error(readinput): error reading BSE conduction state list")')
1494  write(*,'("(blank line required after jstxbse block)")')
1495  write(*,*)
1496  stop
1497  end if
1498  end do
1499  write(*,*)
1500  write(*,'("Error(readinput): BSE extra conduction state list too long")')
1501  write(*,*)
1502  stop
1503 case('bsefull')
1504  read(50,*,err=20) bsefull
1505 case('hxbse')
1506  read(50,*,err=20) hxbse
1507 case('hdbse')
1508  read(50,*,err=20) hdbse
1509 case('gmaxrf','gmaxrpa')
1510  read(50,*,err=20) gmaxrf
1511  if (gmaxrf < 0.d0) then
1512  write(*,*)
1513  write(*,'("Error(readinput): gmaxrf < 0 : ",G18.10)') gmaxrf
1514  write(*,*)
1515  stop
1516  end if
1517 case('mbwgrf')
1518  read(50,*,err=20) mbwgrf
1519 case('emaxrf')
1520  read(50,*,err=20) emaxrf
1521  if (emaxrf < 0.d0) then
1522  write(*,*)
1523  write(*,'("Error(readinput): emaxrf < 0 : ",G18.10)') emaxrf
1524  write(*,*)
1525  stop
1526  end if
1527 case('fxctype')
1528  read(50,'(A)',err=20) str
1529  str=trim(str)//' 0 0'
1530  read(str,*,err=20) fxctype
1531 case('fxclrc')
1532  read(50,'(A)',err=20) str
1533  str=trim(str)//' 0.0'
1534  read(str,*,err=20) fxclrc(:)
1535 case('ntemp')
1536  read(50,*,err=20) ntemp
1537  if (ntemp < 1) then
1538  write(*,*)
1539  write(*,'("Error(readinput): ntemp < 1 : ",I0)') ntemp
1540  write(*,*)
1541  stop
1542  end if
1543 case('trimvg')
1544  write(*,'("Info(readinput): variable ''trimvg'' is no longer used")')
1545  read(50,*,err=20)
1546 case('rndstate','rndseed')
1547  read(50,*,err=20) rndstate(0)
1548  rndstate(0)=abs(rndstate(0))
1549 case('rndatposc')
1550  read(50,*,err=20) rndatposc
1551 case('rndbfcmt')
1552  read(50,*,err=20) rndbfcmt
1553 case('rndavec')
1554  read(50,*,err=20) rndavec
1555 case('c_tb09')
1556  read(50,*,err=20) c_tb09
1557 ! set flag to indicate Tran-Blaha constant has been read in
1558  tc_tb09=.true.
1559 case('lowq','highq','vhighq','uhighq')
1560  read(50,*,err=20) lv
1561  if (lv) then
1562  if (trim(block) == 'lowq') then
1563  rgkmax=6.5d0
1564  gmaxvr=10.d0
1565  lmaxapw=7
1566  lmaxo=5
1567  nxlo=2
1568  lorbcnd=.true.
1569  radkpt=25.d0
1570  autokpt=.true.
1571  vkloff(:)=0.5d0
1572  nempty0=4.d0
1573  epspot=1.d-5
1574  epsengy=5.d-4
1575  epsforce=1.d-2
1576  epsstress=3.d-3
1577  autolinengy=.true.
1578  gmaxrf=2.5d0
1579  lradstp=6
1580  else if (trim(block) == 'highq') then
1581 ! parameter set for high-quality calculation
1582  rgkmax=max(rgkmax,8.d0)
1583  gmaxvr=max(gmaxvr,16.d0)
1584  lmaxapw=max(lmaxapw,9)
1585  lmaxo=max(lmaxo,7)
1586  nrmtscf=max(nrmtscf,1.5d0)
1587  nxlo=max(nxlo,2)
1588  lorbcnd=.true.
1589  radkpt=max(radkpt,50.d0)
1590  autokpt=.true.
1591  vkloff(:)=0.d0
1592  nempty0=max(nempty0,10.d0)
1593  epspot=min(epspot,1.d-7)
1594  epsengy=min(epsengy,1.d-5)
1595  epsforce=min(epsforce,5.d-4)
1596  epsstress=min(epsstress,1.d-3)
1597  autolinengy=.true.
1598  gmaxrf=max(gmaxrf,4.d0)
1599  else if (trim(block) == 'vhighq') then
1600 ! parameter set for very high-quality calculation
1601  rgkmax=max(rgkmax,9.d0)
1602  gmaxvr=max(gmaxvr,18.d0)
1603  lmaxapw=max(lmaxapw,11)
1604  lmaxo=max(lmaxo,9)
1605  nrmtscf=max(nrmtscf,2.d0)
1606  nxlo=max(nxlo,3)
1607  lorbcnd=.true.
1608  radkpt=max(radkpt,90.d0)
1609  autokpt=.true.
1610  vkloff(:)=0.d0
1611  nempty0=max(nempty0,20.d0)
1612  epspot=min(epspot,1.d-7)
1613  epsengy=min(epsengy,1.d-6)
1614  epsforce=min(epsforce,2.d-4)
1615  epsstress=min(epsstress,5.d-4)
1616  autolinengy=.true.
1617  gmaxrf=max(gmaxrf,5.d0)
1618  else
1619 ! parameter set for ultra high-quality calculation
1620  rgkmax=max(rgkmax,10.d0)
1621  gmaxvr=max(gmaxvr,20.d0)
1622  lmaxapw=max(lmaxapw,12)
1623  lmaxo=max(lmaxo,9)
1624  nrmtscf=max(nrmtscf,4.d0)
1625  nxlo=max(nxlo,3)
1626  lorbcnd=.true.
1627  radkpt=max(radkpt,120.d0)
1628  autokpt=.true.
1629  vkloff(:)=0.d0
1630  nempty0=max(nempty0,40.d0)
1631  epspot=min(epspot,1.d-7)
1632  epsengy=min(epsengy,1.d-6)
1633  epsforce=min(epsforce,1.d-4)
1634  epsstress=min(epsstress,2.d-4)
1635  autolinengy=.true.
1636  gmaxrf=max(gmaxrf,6.d0)
1637  end if
1638  if (mp_mpi) then
1639  write(*,*)
1640  write(*,'("Info(readinput): parameters set by ",A," option")') trim(block)
1641  write(*,'(" rgkmax : ",G18.10)') rgkmax
1642  write(*,'(" gmaxvr : ",G18.10)') gmaxvr
1643  write(*,'(" lmaxapw : ",I0)') lmaxapw
1644  write(*,'(" lmaxo : ",I0)') lmaxo
1645  write(*,'(" nrmtscf : ",G18.10)') nrmtscf
1646  write(*,'(" nxlo : ",I0)') nxlo
1647  write(*,'(" lorbcnd : ",L1)') lorbcnd
1648  write(*,'(" radkpt : ",G18.10)') radkpt
1649  write(*,'(" autokpt : ",L1)') autokpt
1650  write(*,'(" vkloff : ",3G18.10)') vkloff
1651  write(*,'(" nempty0 : ",G18.10)') nempty0
1652  write(*,'(" epspot : ",G18.10)') epspot
1653  write(*,'(" epsengy : ",G18.10)') epsengy
1654  write(*,'(" epsforce : ",G18.10)') epsforce
1655  write(*,'(" epsstress : ",G18.10)') epsstress
1656  write(*,'(" autolinengy : ",L1)') autolinengy
1657  write(*,'(" gmaxrf : ",G18.10)') gmaxrf
1658  if (trim(block) == 'lowq') then
1659  write(*,'(" lradstp : ",I0)') lradstp
1660  end if
1661  end if
1662  end if
1663 case('hmaxvr')
1664  read(50,*,err=20) hmaxvr
1665  if (hmaxvr < 0.d0) then
1666  write(*,*)
1667  write(*,'("Error(readinput): hmaxvr < 0 : ",G18.10)') hmaxvr
1668  write(*,*)
1669  stop
1670  end if
1671 case('hkmax')
1672  read(50,*,err=20) hkmax
1673  if (hkmax <= 0.d0) then
1674  write(*,*)
1675  write(*,'("Error(readinput): hkmax <= 0 : ",G18.10)') hkmax
1676  write(*,*)
1677  stop
1678  end if
1679 case('lorbcnd')
1680  read(50,*,err=20) lorbcnd
1681 case('lorbordc')
1682  read(50,*,err=20) lorbordc
1683  if (lorbordc < 2) then
1684  write(*,*)
1685  write(*,'("Error(readinput): lorbordc < 2 : ",I0)') lorbordc
1686  write(*,*)
1687  stop
1688  end if
1689  if (lorbordc > maxlorbord) then
1690  write(*,*)
1691  write(*,'("Error(readinput): lorbordc too large : ",I0)') lorbordc
1692  write(*,'("Adjust maxlorbord in modmain and recompile code")')
1693  write(*,*)
1694  stop
1695  end if
1696 case('nrmtscf')
1697  read(50,'(A)',err=20) str
1698  read(str,*,err=20) nrmtscf
1699  read(str,*,iostat=ios) nrmtscf,dnrmtscf
1700  if (nrmtscf < 0.5d0) then
1701  write(*,*)
1702  write(*,'("Error(readinput): nrmtscf < 0.5 : ",G18.10)') nrmtscf
1703  write(*,*)
1704  stop
1705  end if
1706 case('lmaxdb','lmaxdos')
1707  read(50,*,err=20) lmaxdb
1708  if (lmaxdb < 0) then
1709  write(*,*)
1710  write(*,'("Error(readinput): lmaxdb < 0 : ",I0)') lmaxdb
1711  write(*,*)
1712  stop
1713  end if
1714 case('epsdev')
1715  read(50,*,err=20) epsdev
1716  if (epsdev <= 0.d0) then
1717  write(*,*)
1718  write(*,'("Error(readinput): epsdev <= 0 : ",G18.10)') epsdev
1719  write(*,*)
1720  stop
1721  end if
1722 case('msmooth')
1723  read(50,*,err=20)
1724  write(*,'("Info(readinput): variable ''msmooth'' is no longer used")')
1725 case('npmae')
1726  read(50,*,err=20) npmae0
1727 case('wrtvars')
1728  read(50,*,err=20) wrtvars
1729 case('ftmtype')
1730  read(50,*,err=20) ftmtype
1731 case('tmomfix')
1732  write(*,*)
1733  write(*,'("Error(readinput): variable ''tmomfix'' is no longer used")')
1734  write(*,'(" use tm3fix instead")')
1735  write(*,*)
1736  stop
1737 case('tm3fix')
1738  read(50,*,err=20) ntmfix
1739  if (ntmfix < 1) then
1740  write(*,*)
1741  write(*,'("Error(readinput): ntmfix < 1 : ",I0)') ntmfix
1742  write(*,*)
1743  stop
1744  end if
1745  if (allocated(itmfix)) deallocate(itmfix)
1746  allocate(itmfix(7,ntmfix))
1747  if (allocated(wkprfix)) deallocate(wkprfix)
1748  allocate(wkprfix(ntmfix))
1749  do i=1,ntmfix
1750  read(50,*,err=20) is,ia,l
1751  if ((is < 1).or.(ia < 1).or.(l < 0)) then
1752  write(*,*)
1753  write(*,'("Error(readinput): invalid is, ia or l in tm3fix block :",&
1754  &3(X,I0))') is,ia,l
1755  write(*,*)
1756  stop
1757  end if
1758  itmfix(1,i)=is
1759  itmfix(2,i)=ia
1760  itmfix(3,i)=l
1761 ! read k, p, r, t for the 3-index tensor
1762  read(50,*,err=20) itmfix(4:7,i)
1763 ! read 3-index tensor component with conventional normalisation
1764  read(50,*,err=20) wkprfix(i)
1765  end do
1766 case('tauftm')
1767  read(50,*,err=20) tauftm
1768  if (tauftm < 0.d0) then
1769  write(*,*)
1770  write(*,'("Error(readinput): tauftm < 0 : ",G18.10)') tauftm
1771  write(*,*)
1772  stop
1773  end if
1774 case('ftmstep')
1775  read(50,*,err=20)
1776  write(*,'("Info(readinput): variable ''ftmstep'' is no longer used")')
1777 case('cmagz','forcecmag')
1778  read(50,*,err=20) cmagz
1779 case('rotavec')
1780  read(50,*,err=20) axang(:)
1781 case('tstime')
1782  read(50,*,err=20) tstime
1783  if (tstime <= 0.d0) then
1784  write(*,*)
1785  write(*,'("Error(readinput): tstime <= 0 : ",G18.10)') tstime
1786  write(*,*)
1787  stop
1788  end if
1789 case('dtimes')
1790  read(50,*,err=20) dtimes
1791  if (dtimes <= 0.d0) then
1792  write(*,*)
1793  write(*,'("Error(readinput): dtimes <= 0 : ",G18.10)') dtimes
1794  write(*,*)
1795  stop
1796  end if
1797 case('pulse')
1798  read(50,*,err=20) npulse
1799  if (npulse < 1) then
1800  write(*,*)
1801  write(*,'("Error(readinput): npulse < 1 : ",I0)') npulse
1802  write(*,*)
1803  stop
1804  end if
1805  if (allocated(pulse)) deallocate(pulse)
1806  allocate(pulse(12,npulse))
1807  do i=1,npulse
1808  read(50,'(A)',err=20) str
1809  str=trim(str)//' 1.0 0.0 0.0 0.0'
1810  read(str,*,err=20) pulse(:,i)
1811  end do
1812 case('ramp')
1813  read(50,*,err=20) nramp
1814  if (nramp < 1) then
1815  write(*,*)
1816  write(*,'("Error(readinput): nramp < 1 : ",I0)') nramp
1817  write(*,*)
1818  stop
1819  end if
1820  if (allocated(ramp)) deallocate(ramp)
1821  allocate(ramp(12,nramp))
1822  do i=1,nramp
1823  read(50,'(A)',err=20) str
1824  str=trim(str)//' 1.0 0.0 0.0 0.0'
1825  read(str,*,err=20) ramp(:,i)
1826  end do
1827 case('step')
1828  read(50,*,err=20) nstep
1829  if (nstep < 1) then
1830  write(*,*)
1831  write(*,'("Error(readinput): nstep < 1 : ",I0)') nstep
1832  write(*,*)
1833  stop
1834  end if
1835  if (allocated(step)) deallocate(step)
1836  allocate(step(9,nstep))
1837  do i=1,nstep
1838  read(50,'(A)',err=20) str
1839  str=trim(str)//' 1.0 0.0 0.0 0.0'
1840  read(str,*,err=20) step(:,i)
1841  end do
1842 case('ncgga')
1843  read(50,*,err=20)
1844  write(*,'("Info(readinput): variable ''ncgga'' is no longer used")')
1845 case('dncgga')
1846  read(50,*,err=20) dncgga
1847  if (dncgga < 0.d0) then
1848  write(*,*)
1849  write(*,'("Error(readinput): dncgga < 0 : ",G18.10)') dncgga
1850  write(*,*)
1851  stop
1852  end if
1853 case('ntswrite')
1854  read(50,'(A)',err=20) str
1855  str=trim(str)//' 1'
1856  read(str,*,err=20) ntswrite(:)
1857 case('nxoapwlo','nxapwlo')
1858  read(50,*,err=20) nxoapwlo
1859  if (nxoapwlo < 0) then
1860  write(*,*)
1861  write(*,'("Error(readinput): nxoapwlo < 0 : ",I0)') nxoapwlo
1862  write(*,*)
1863  stop
1864  end if
1865 case('nxlo')
1866  read(50,*,err=20) nxlo
1867  if (nxlo < 0) then
1868  write(*,*)
1869  write(*,'("Error(readinput): nxlo < 0 : ",I0)') nxlo
1870  write(*,*)
1871  stop
1872  end if
1873 case('tdrho1d')
1874  read(50,*,err=20) tdrho1d
1875 case('tdrho2d')
1876  read(50,*,err=20) tdrho2d
1877 case('tdrho3d')
1878  read(50,*,err=20) tdrho3d
1879 case('tdmag1d')
1880  read(50,*,err=20) tdmag1d
1881 case('tdmag2d')
1882  read(50,*,err=20) tdmag2d
1883 case('tdmag3d')
1884  read(50,*,err=20) tdmag3d
1885 case('tdjr1d','tdcd1d')
1886  read(50,*,err=20) tdjr1d
1887 case('tdjr2d','tdcd2d')
1888  read(50,*,err=20) tdjr2d
1889 case('tdjr3d','tdcd3d')
1890  read(50,*,err=20) tdjr3d
1891 case('tddos')
1892  read(50,*,err=20) tddos
1893 case('tdlsj')
1894  read(50,*,err=20) tdlsj
1895 case('tdjtk')
1896  read(50,*,err=20) tdjtk
1897 case('tdxrmk')
1898  read(50,*,err=20) tdxrmk
1899 case('epseph')
1900  read(50,*,err=20)
1901  write(*,'("Info(readinput): variable ''epseph'' is no longer used")')
1902 case('rndevt0')
1903  read(50,*,err=20) rndevt0
1904 case('sxcscf','ssxc','rstsf')
1905  read(50,'(A)',err=20) str
1906  read(str,*,err=20) sxcscf
1907  read(str,*,iostat=ios) sxcscf,dsxcscf
1908 case('tempk')
1909  read(50,*,err=20) tempk
1910  if (tempk <= 0.d0) then
1911  write(*,*)
1912  write(*,'("Error(readinput): tempk <= 0 : ",G18.10)') tempk
1913  write(*,*)
1914  stop
1915  end if
1916 ! set Fermi-Dirac smearing
1917  stype=3
1918 ! set the smearing width
1920 case('avecu')
1921  read(50,*,err=20) avecu(:,1)
1922  read(50,*,err=20) avecu(:,2)
1923  read(50,*,err=20) avecu(:,3)
1924 case('scaleu')
1925  read(50,*,err=20) scu
1926 case('scaleu1')
1927  read(50,*,err=20) scu1
1928 case('scaleu2')
1929  read(50,*,err=20) scu2
1930 case('scaleu3')
1931  read(50,*,err=20) scu3
1932 case('q0cut')
1933  read(50,*,err=20) q0cut
1934 case('ngridkpa')
1935  read(50,*,err=20) ngridkpa
1936 case('rndbfcu')
1937  read(50,*,err=20) rndbfcu
1938 case('bfieldcu','bfielduc')
1939  read(50,*,err=20) bfieldcu
1940 case('efieldcu','efielduc')
1941  read(50,*,err=20) efieldcu
1942 case('tplotq0')
1943  read(50,*,err=20) tplotq0
1944 case('trdvclr')
1945  read(50,*,err=20) trdvclr
1946 case('trdbfcr')
1947  read(50,*,err=20) trdbfcr
1948 case('evtype')
1949  read(50,*,err=20)
1950  write(*,'("Info(readinput): variable ''evtype'' is no longer used")')
1951 case('wmaxgw')
1952  read(50,*,err=20) wmaxgw
1953 case('twdiag')
1954  read(50,*,err=20)
1955  write(*,'("Info(readinput): variable ''twdiag'' is no longer used")')
1956 case('tsediag')
1957  read(50,*,err=20) tsediag
1958 case('actype')
1959  read(50,*,err=20) actype
1960 case('npole')
1961  read(50,*,err=20) npole
1962  if (npole < 1) then
1963  write(*,*)
1964  write(*,'("Error(readinput): npole < 1 : ",I0)') npole
1965  write(*,*)
1966  stop
1967  end if
1968 case('nspade')
1969  read(50,*,err=20) nspade
1970  if (nspade < 1) then
1971  write(*,*)
1972  write(*,'("Error(readinput): nspade < 1 : ",I0)') nspade
1973  write(*,*)
1974  stop
1975  end if
1976 case('tfav0')
1977  read(50,*,err=20) tfav0
1978 case('rmtscf')
1979  read(50,*,err=20) rmtscf
1980  if (rmtscf <= 0.d0) then
1981  write(*,*)
1982  write(*,'("Error(readinput): rmtscf <= 0 : ",G18.10)') rmtscf
1983  write(*,*)
1984  stop
1985  end if
1986 case('mrmtav')
1987  read(50,*,err=20) mrmtav
1988 case('rmtall')
1989  read(50,*,err=20) rmtall
1990 case('maxthd','omp_num_threads','OMP_NUM_THREADS')
1991  read(50,*,err=20) maxthd
1992 case('maxthd1')
1993  read(50,*,err=20) maxthd1
1994 case('maxthdmkl')
1995  read(50,*,err=20) maxthdmkl
1996 case('maxlvl','omp_max_active_levels','OMP_MAX_ACTIVE_LEVELS')
1997  read(50,*,err=20) maxlvl
1998  if (maxlvl < 1) then
1999  write(*,*)
2000  write(*,'("Error(readinput): maxlvl < 1 : ",I0)') maxlvl
2001  write(*,*)
2002  stop
2003  end if
2004 case('stable')
2005  read(50,*,err=20) lv
2006  if (lv) then
2007  autolinengy=.true.
2008  mrmtav=1
2009  lmaxapw=max(lmaxapw,10)
2010  gmaxvr=max(gmaxvr,24.d0)
2011  msmgmt=max(msmgmt,1)
2012  if (mp_mpi) then
2013  write(*,*)
2014  write(*,'("Info(readinput): parameters set by stable option")')
2015  write(*,'(" autolinengy : ",L1)') autolinengy
2016  write(*,'(" mrmtav : ",I0)') mrmtav
2017  write(*,'(" lmaxapw : ",I0)') lmaxapw
2018  write(*,'(" gmaxvr : ",G18.10)') gmaxvr
2019  write(*,'(" msmgmt : ",I0)') msmgmt
2020  end if
2021  end if
2022 case('metagga')
2023  read(50,*,err=20) lv
2024  if (lv) then
2025  lmaxi=max(lmaxi,2)
2026  gmaxvr=max(gmaxvr,16.d0)
2027  nrmtscf=max(nrmtscf,3.d0)
2028  msmgmt=max(msmgmt,4)
2029  epspot=1.d6
2030  epsengy=min(epsengy,1.d-6)
2031  if (mp_mpi) then
2032  write(*,*)
2033  write(*,'("Info(readinput): parameters set by metagga option")')
2034  write(*,'(" lmaxi : ",I0)') lmaxi
2035  write(*,'(" gmaxvr : ",G18.10)') gmaxvr
2036  write(*,'(" nrmtscf : ",G18.10)') nrmtscf
2037  write(*,'(" msmgmt : ",I0)') msmgmt
2038  write(*,'(" epspot : ",G18.10)') epspot
2039  write(*,'(" epsengy : ",G18.10)') epsengy
2040  end if
2041  end if
2042 case('t0tdlr')
2043  read(50,*,err=20)
2044  write(*,'("Info(readinput): variable ''t0tdlr'' is no longer used")')
2045 case('tdphi')
2046  read(50,*,err=20) tdphi
2047 ! convert phase from degrees to radians
2048  tdphi=tdphi*pi/180.d0
2049 case('thetamld')
2050  read(50,*,err=20) thetamld
2051 ! convert MLD angle from degrees to radians
2052  thetamld=thetamld*pi/180.d0
2053 case('ntsbackup')
2054  read(50,*,err=20) ntsbackup
2055 case('seedname')
2056  read(50,*,err=20) seedname
2057  seedname=adjustl(seedname)
2058 case('num_wann')
2059  read(50,*,err=20) num_wann
2060 case('idxw90','wann_bands')
2061  read(50,'(A)',err=20) str
2062  num_bands=1024
2063  if (allocated(idxw90)) deallocate(idxw90)
2064  allocate(idxw90(num_bands))
2065  call numlist(str,num_bands,idxw90)
2066 case('projw90')
2067  read(50,*,err=20) projw90
2068 case('lprojw90')
2069  do
2070  read(50,'(A)',err=20) str
2071  if (trim(str) == '') goto 10
2072  read(str,*) is,l
2073  if ((is < 1).or.(is > maxspecies)) then
2074  write(*,*)
2075  write(*,'("Error(readinput): invalid projection species : ",I0)') is
2076  write(*,*)
2077  stop
2078  end if
2079  if ((l < 0).or.(l > 3)) then
2080  write(*,*)
2081  write(*,'("Error(readinput): invalid projection l : ",I0)') l
2082  write(*,*)
2083  stop
2084  end if
2085  lprojw90(is)=l
2086  end do
2087 case('num_iter')
2088  read(50,*,err=20) num_iter
2089 case('dis_num_iter')
2090  read(50,*,err=20) dis_num_iter
2091 case('trial_step')
2092  read(50,*,err=20) trial_step
2093 case('xlwin','wannierExtra')
2094  if (allocated(xlwin)) deallocate(xlwin)
2095  allocate(xlwin(0))
2096  nxlwin=0
2097  do
2098  read(50,'(A)',err=20) str
2099  if (trim(str) == '') goto 10
2100  xlwin=[xlwin(1:nxlwin),str]
2101  nxlwin=nxlwin+1
2102  end do
2103 case('wrtunk')
2104  read(50,*,err=20) wrtunk
2105 case('tbdip')
2106  read(50,*,err=20) tbdip
2107 case('tjr','tcden')
2108  read(50,*,err=20) tjr
2109 case('tauefm')
2110  read(50,*,err=20) tauefm
2111 case('epsefm')
2112  read(50,*,err=20) epsefm
2113 case('ehfb')
2114  read(50,*,err=20) ehfb
2115 case('t0gclq0')
2116  read(50,*,err=20)
2117  write(*,'("Info(readinput): variable ''t0gclq0'' is no longer used")')
2118 case('tafindt')
2119  read(50,*,err=20) tafindt
2120 case('afindscf')
2121  read(50,*,err=20)
2122  write(*,'("Info(readinput): variable ''afindscf'' is no longer used")')
2123 case('afindpm')
2124  read(50,*,err=20) afindpm(:)
2125  if (afindpm(2) == 0.d0) then
2126  write(*,*)
2127  write(*,'("Error(readinput): afindpm(2) = 0")')
2128  write(*,*)
2129  stop
2130  end if
2131 case('nkspolar')
2132  read(50,*,err=20) nkspolar
2133  if (nkspolar < 1) then
2134  write(*,*)
2135  write(*,'("Error(readinput): nkspolar < 1 : ",I0)') nkspolar
2136  write(*,*)
2137  stop
2138  end if
2139 case('ntsforce')
2140  read(50,*,err=20) ntsforce
2141  if (ntsforce < 1) then
2142  write(*,*)
2143  write(*,'("Error(readinput): ntsforce < 1 : ",I0)') ntsforce
2144  write(*,*)
2145  stop
2146  end if
2147 case('wphcut')
2148  read(50,*,err=20) wphcut
2149  if (wphcut <= 0.d0) then
2150  write(*,*)
2151  write(*,'("Error(readinput): wphcut <= 0 : ",G18.10)') wphcut
2152  write(*,*)
2153  stop
2154  end if
2155 case('ephscf')
2156  read(50,*,err=20) ephscf(:)
2157 case('anomalous')
2158  read(50,*,err=20) anomalous
2159 case('tephde')
2160  read(50,*,err=20) tephde
2161 case('bdiag')
2162  read(50,*,err=20) bdiag
2163 case('ecutb')
2164  read(50,*,err=20) ecutb
2165  if (ecutb <= 0.d0) then
2166  write(*,*)
2167  write(*,'("Error(readinput): ecutb <= 0 : ",G18.10)') ecutb
2168  write(*,*)
2169  stop
2170  end if
2171 case('ediag')
2172  read(50,*,err=20) ediag
2173 case('pwxpsn')
2174  read(50,*,err=20) pwxpsn
2175  if (pwxpsn < 1) then
2176  write(*,*)
2177  write(*,'("Error(readinput): pwxpsn < 1 : ",I0)') pwxpsn
2178  write(*,*)
2179  stop
2180  end if
2181 case('ramdisk')
2182  read(50,*,err=20) ramdisk
2183 case('wrtdisk','wrtdsk')
2184  read(50,*,err=20) wrtdisk
2185 case('epsdmat')
2186  read(50,*,err=20) epsdmat
2187 case('tm3vdl','tm3old')
2188  read(50,*,err=20) tm3vdl
2189 case('batch')
2190  read(50,*,err=20) batch
2191 case('tafspt')
2192  read(50,*,err=20) tafspt
2193 case('tbaspat')
2194  read(50,*,err=20) tbaspat
2195 case('trdatdv')
2196  read(50,*,err=20) trdatdv
2197 case('atdfc')
2198  read(50,*,err=20) atdfc
2199  if (atdfc < 0.d0) then
2200  write(*,*)
2201  write(*,'("Error(readinput): atdfc < 0 : ",G18.10)') atdfc
2202  write(*,*)
2203  stop
2204  end if
2205 case('maxforce')
2206  read(50,*,err=20) maxforce
2207 case('msmgmt','msmg2mt')
2208  read(50,*,err=20) msmgmt
2209 case('ntsorth')
2210  read(50,*,err=20) ntsorth
2211 case('deltabf')
2212  read(50,*,err=20) deltabf
2213  if (deltabf <= 0.d0) then
2214  write(*,*)
2215  write(*,'("Error(readinput): deltabf <= 0 : ",G18.10)') deltabf
2216  write(*,*)
2217  stop
2218  end if
2219 case('jtconst0')
2220  read(50,*,err=20) jtconst0
2221 case('trmt0')
2222  read(50,*,err=20) trmt0
2223 case('ksgwrho')
2224  read(50,*,err=20) ksgwrho
2225 case('npfftg')
2226  read(50,*,err=20) npfftg
2227 case('npfftgc')
2228  read(50,*,err=20) npfftgc
2229 case('npfftq')
2230  read(50,*,err=20) npfftq
2231 case('npfftw')
2232  read(50,*,err=20) npfftw
2233 case('tphnat')
2234  read(50,*,err=20) tphnat
2235 case('ecutthc')
2236  read(50,*,err=20) ecutthc
2237  if (ecutthc <= 0.d0) then
2238  write(*,*)
2239  write(*,'("Error(readinput): ecutthc <= 0 : ",G18.10)') ecutthc
2240  write(*,*)
2241  stop
2242  end if
2243 case('tbdipu')
2244  read(50,*,err=20) tbdipu
2245 case('bdipscf')
2246  read(50,*,err=20) bdipscf
2247 case('')
2248  goto 10
2249 case default
2250  write(*,*)
2251  write(*,'("Error(readinput): invalid block name : ",A)') trim(block)
2252  write(*,*)
2253  stop
2254 end select
2255 goto 10
2256 20 continue
2257 write(*,*)
2258 write(*,'("Error(readinput): error reading from elk.in")')
2259 write(*,'("Problem occurred in ''",A,"'' block")') trim(block)
2260 write(*,'("Check input convention in manual")')
2261 write(*,*)
2262 stop
2263 30 continue
2264 close(50)
2265 ! scale the speed of light
2266 solsc=sol*solscf
2267 ! scale and rotate the lattice vectors (not referenced again in code)
2268 avec(:,:)=sc*avec(:,:)
2269 avec(:,1)=sc1*avec(:,1)
2270 avec(:,2)=sc2*avec(:,2)
2271 avec(:,3)=sc3*avec(:,3)
2272 avec(1,:)=scx*avec(1,:)
2273 avec(2,:)=scy*avec(2,:)
2274 avec(3,:)=scz*avec(3,:)
2275 t1=axang(4)
2276 if (t1 /= 0.d0) then
2277  t1=t1*pi/180.d0
2278  call axangrot(axang(:),t1,rot)
2279  do i=1,3
2280  v(:)=avec(:,i)
2281  call r3mv(rot,v,avec(:,i))
2282  end do
2283 end if
2284 ! randomise lattice vectors if required
2285 if (rndavec > 0.d0) then
2286  do i=1,3
2287  do j=1,3
2288  t1=rndavec*(randomu()-0.5d0)
2289  avec(i,j)=avec(i,j)+t1
2290  end do
2291  end do
2292 end if
2293 ! check if reference lattice vectors should be used
2294 tavref=(any(abs(avecref(:,:)) > epslat))
2295 ! case of isolated molecule
2296 if (molecule) then
2297 ! convert atomic positions from Cartesian to lattice coordinates
2298  call r3minv(avec,ainv)
2299  do is=1,nspecies
2300  do ia=1,natoms(is)
2301  call r3mv(ainv,atposl(:,ia,is),v)
2302  atposl(:,ia,is)=v(:)
2303  end do
2304  end do
2305 end if
2306 ! randomise atomic positions if required
2307 if (rndatposc > 0.d0) then
2308  call r3minv(avec,ainv)
2309  do is=1,nspecies
2310  do ia=1,natoms(is)
2311  call r3mv(avec,atposl(:,ia,is),v)
2312  do i=1,3
2313  t1=rndatposc*(randomu()-0.5d0)
2314  v(i)=v(i)+t1
2315  end do
2316  call r3mv(ainv,v,atposl(:,ia,is))
2317  end do
2318  end do
2319 end if
2320 ! randomise the muffin-tin magnetic fields if required
2321 if (rndbfcmt > 0.d0) then
2322  do is=1,nspecies
2323  do ia=1,natoms(is)
2324  do i=1,3
2325  t1=rndbfcmt*(randomu()-0.5d0)
2326  bfcmt0(i,ia,is)=bfcmt0(i,ia,is)+t1
2327  end do
2328  end do
2329  end do
2330 end if
2331 ! set fxctype to fxctype if required
2332 if (fxctype(1) == -1) fxctype(:)=xctype(:)
2333 ! find primitive cell if required
2334 if (primcell) call findprimcell
2335 ! scale the ultracell vectors if required
2336 avecu(:,1)=scu1*avecu(:,1)
2337 avecu(:,2)=scu2*avecu(:,2)
2338 avecu(:,3)=scu3*avecu(:,3)
2339 avecu(:,:)=scu*avecu(:,:)
2340 ! read in atomic species data
2341 call readspecies
2342 return
2343 
2344 end subroutine
2345 !EOC
2346 
real(8) socscf
Definition: modmain.f90:232
integer ngrkf
Definition: modmain.f90:1067
real(8) ecvcut
Definition: modmain.f90:117
integer maxlatvstp
Definition: modmain.f90:1031
real(8), dimension(3, 3) dafspc
Definition: modmain.f90:331
real(8) epsband
Definition: modmain.f90:818
subroutine genfspecies(zn, symb)
Definition: genfspecies.f90:7
character(256) scrpath
Definition: modmain.f90:1295
real(8) dchgexs
Definition: modmain.f90:722
real(8), dimension(3, maxatoms, maxspecies) bfcmt0
Definition: modmain.f90:275
integer, parameter maxspecies
Definition: modmain.f90:30
integer ncxbse
Definition: modmain.f90:1177
real(8) scissor
Definition: modmain.f90:905
real(8), dimension(3, 3) afspc
Definition: modmain.f90:331
subroutine findprimcell
real(8) dmaxefc
Definition: modmain.f90:318
real(8) epsforce
Definition: modmain.f90:1057
logical tjr
Definition: modmain.f90:618
integer, dimension(maxxbse) jstxbse
Definition: modmain.f90:1179
integer lorbordc
Definition: modmain.f90:835
integer maxthd1
Definition: modomp.f90:13
integer mixtype
Definition: modmain.f90:693
integer, dimension(3) ktype
Definition: modmain.f90:604
real(8), dimension(3, 3) davec
Definition: modmain.f90:12
real(8), dimension(3, maxatoms, maxspecies) datposl
Definition: modmain.f90:52
integer maxthd
Definition: modomp.f90:11
integer, dimension(3) xctype
Definition: modmain.f90:586
logical spinpol
Definition: modmain.f90:228
integer, parameter maxlapw
Definition: modmain.f90:195
integer npfftg
Definition: modmain.f90:404
integer mbwgrf
Definition: modmain.f90:1158
integer, dimension(3) xctsp
Definition: modmain.f90:142
logical tshift
Definition: modmain.f90:352
character(256) sppath
Definition: modmain.f90:72
real(8) rmtall
Definition: modmain.f90:158
real(8) taufsm
Definition: modmain.f90:265
real(8) reducebf
Definition: modmain.f90:279
logical autokpt
Definition: modmain.f90:444
integer msmgmt
Definition: modmain.f90:222
logical tfav0
Definition: modmain.f90:995
real(8) dnrmtscf
Definition: modmain.f90:148
real(8) sxcscf
Definition: modmain.f90:666
real(8), dimension(3, 3) ainv
Definition: modmain.f90:14
real(8), dimension(2) amixpm
Definition: modmain.f90:700
logical lorbcnd
Definition: modmain.f90:833
real(8) rmtscf
Definition: modmain.f90:154
subroutine r3minv(a, b)
Definition: r3minv.f90:10
logical spinsprl
Definition: modmain.f90:283
Definition: modomp.f90:6
real(8), parameter kboltz
Definition: modmain.f90:1255
real(8) swidth
Definition: modmain.f90:889
logical hybrid0
Definition: modmain.f90:1144
real(8) epsocc
Definition: modmain.f90:897
real(8) dncgga
Definition: modmain.f90:602
character(256), dimension(maxspecies) spfname
Definition: modmain.f90:74
real(8) dnempty0
Definition: modmain.f90:877
real(8), dimension(3, 0:3) vclp3d
Definition: modmain.f90:1124
integer, dimension(maxxbse) istxbse
Definition: modmain.f90:1179
integer nvxbse
Definition: modmain.f90:1177
integer notelns
Definition: modmain.f90:1297
integer, parameter maxxbse
Definition: modmain.f90:1175
subroutine numlist(str, n, list)
Definition: numlist.f90:10
logical mixrho
Definition: modmain.f90:685
logical nosource
Definition: modmain.f90:662
real(8), dimension(3) vkloff
Definition: modmain.f90:450
integer maxthdmkl
Definition: modomp.f90:15
integer lmaxo
Definition: modmain.f90:201
real(8), dimension(3) vqlss
Definition: modmain.f90:293
logical trhonorm
Definition: modmain.f90:616
integer, parameter maxkst
Definition: modmain.f90:917
real(8) c_tb09
Definition: modmain.f90:676
real(8), parameter pi
Definition: modmain.f90:1224
real(8), dimension(3) vecql
Definition: modmain.f90:1096
real(8) nempty0
Definition: modmain.f90:877
integer nkstlist
Definition: modmain.f90:919
logical dosmsum
Definition: modmain.f90:1079
integer nxlo
Definition: modmain.f90:839
logical tc_tb09
Definition: modmain.f90:678
integer dlmaxo
Definition: modmain.f90:201
integer ntasks
Definition: modmain.f90:1285
integer lmaxapw
Definition: modmain.f90:197
real(8) nrmtscf
Definition: modmain.f90:148
integer symtype
Definition: modmain.f90:340
real(8) epschg
Definition: modmain.f90:710
real(8), dimension(2) broydpm
Definition: modmain.f90:704
logical tefvr
Definition: modmain.f90:866
real(8), dimension(3) momfix
Definition: modmain.f90:253
integer mixsdb
Definition: modmain.f90:702
logical bforb
Definition: modmain.f90:234
logical tforce
Definition: modmain.f90:981
real(8) rmtdelta
Definition: modmain.f90:160
real(8), dimension(3, 3) avecref
Definition: modmain.f90:1022
Definition: modrdm.f90:6
integer ncbse0
Definition: modmain.f90:1173
integer maxlvl
Definition: modomp.f90:17
real(8), dimension(2) wplot
Definition: modmain.f90:1071
real(8), dimension(3, maxatoms, maxspecies) atposl
Definition: modmain.f90:51
real(8) emaxrf
Definition: modmain.f90:1154
real(8) gmaxrf
Definition: modmain.f90:1152
real(8) delorb
Definition: modmain.f90:800
logical lmirep
Definition: modmain.f90:1090
real(8) esccut
Definition: modmain.f90:119
Definition: modw90.f90:6
real(8), dimension(3) sqaxis
Definition: modmain.f90:1093
logical intraband
Definition: modmain.f90:1087
real(8) rndatposc
Definition: modmain.f90:56
real(8) radkpt
Definition: modmain.f90:446
integer, parameter maxlorbord
Definition: modmain.f90:780
logical tefvs
Definition: modmain.f90:869
integer nvbse0
Definition: modmain.f90:1173
Definition: modbog.f90:6
real(8) deltaem
Definition: modmain.f90:481
integer lmaxdb
Definition: modmain.f90:1073
real(8), dimension(3) dbfieldc0
Definition: modmain.f90:271
pure subroutine axangrot(v, th, rot)
Definition: axangrot.f90:10
real(8), dimension(3, 3) avec
Definition: modmain.f90:12
real(8) rndbfcmt
Definition: modmain.f90:277
real(8), dimension(2) wsfac
Definition: modmain.f90:1100
real(8) tempk
Definition: modmain.f90:682
integer, dimension(3, 27) optcomp
Definition: modmain.f90:1085
integer, dimension(2, maxkst) kstlist
Definition: modmain.f90:921
real(8) rgkmax
Definition: modmain.f90:493
integer, dimension(3) ngridk
Definition: modmain.f90:448
real(8) maxforce
Definition: modmain.f90:992
logical cmagz
Definition: modmain.f90:242
integer lradstp
Definition: modmain.f90:171
integer, parameter maxatoms
Definition: modmain.f90:32
real(8), dimension(3) dafieldc
Definition: modmain.f90:325
integer, parameter maxtasks
Definition: modmain.f90:1283
real(8), parameter sol
Definition: modmain.f90:1243
logical tpdos
Definition: modmain.f90:1077
real(8), dimension(3) dvqlss
Definition: modmain.f90:293
real(8), dimension(3) afieldc
Definition: modmain.f90:325
integer mefvs
Definition: modmain.f90:871
integer latvopt
Definition: modmain.f90:1029
real(8), dimension(3, 0:2) vclp2d
Definition: modmain.f90:1120
real(8), dimension(3, maxatoms, maxspecies) mommtfix
Definition: modmain.f90:259
logical mixsave
Definition: modmain.f90:698
integer ip01d
Definition: modmain.f90:1110
logical tavref
Definition: modmain.f90:1024
integer npfftq
Definition: modmain.f90:535
real(8) solsc
Definition: modmain.f90:1245
real(8) tau0atp
Definition: modmain.f90:1003
integer, dimension(3) ngridq
Definition: modmain.f90:515
real(8) momfixm
Definition: modmain.f90:255
real(8) drgkmax
Definition: modmain.f90:493
Definition: modgw.f90:6
integer, dimension(maxspecies) natoms
Definition: modmain.f90:36
integer nwrite
Definition: modmain.f90:1051
logical ssdph
Definition: modmain.f90:285
integer noptcomp
Definition: modmain.f90:1083
real(8) epslat
Definition: modmain.f90:24
integer maxitoep
Definition: modmain.f90:1134
real(8), dimension(maxatoms, maxspecies) mommtfixm
Definition: modmain.f90:261
integer ndspem
Definition: modmain.f90:483
integer stype
Definition: modmain.f90:885
logical hxbse
Definition: modmain.f90:1201
Definition: modmpi.f90:6
subroutine readspecies
Definition: readspecies.f90:7
logical dosssum
Definition: modmain.f90:1081
real(8) mstar
Definition: modmain.f90:893
real(8) epspot
Definition: modmain.f90:1053
real(8) dsxcscf
Definition: modmain.f90:666
logical ptnucl
Definition: modmain.f90:83
Definition: modpw.f90:6
real(8), dimension(:,:), allocatable vvlp1d
Definition: modmain.f90:1112
real(8) dlefe
Definition: modmain.f90:828
logical autolinengy
Definition: modmain.f90:826
integer nswplot
Definition: modmain.f90:1069
logical tbdip
Definition: modmain.f90:641
real(8) deapw
Definition: modmain.f90:762
logical spinorb
Definition: modmain.f90:230
integer maxatpstp
Definition: modmain.f90:1001
real(8), dimension(3) bfieldc0
Definition: modmain.f90:271
real(8) tau0oep
Definition: modmain.f90:1136
integer, dimension(3) dngridk
Definition: modmain.f90:448
subroutine readinput
Definition: readinput.f90:10
integer nspecies
Definition: modmain.f90:34
integer, dimension(2) np2d
Definition: modmain.f90:1122
logical bsefull
Definition: modmain.f90:1198
real(8), dimension(3) dmomfix
Definition: modmain.f90:253
logical autodlefe
Definition: modmain.f90:831
integer, dimension(3) np3d
Definition: modmain.f90:1126
integer reducek
Definition: modmain.f90:455
real(8) emaxelnes
Definition: modmain.f90:1098
real(8) tau0latv
Definition: modmain.f90:1033
real(8) bdipscf
Definition: modmain.f90:643
logical autoswidth
Definition: modmain.f90:891
integer nkspolar
Definition: modmain.f90:485
integer reduceq
Definition: modmain.f90:519
character(256), dimension(:), allocatable notes
Definition: modmain.f90:1299
pure subroutine r3mv(a, x, y)
Definition: r3mv.f90:10
Definition: modulr.f90:6
integer nwplot
Definition: modmain.f90:1065
real(8) deltabf
Definition: modmain.f90:281
integer mrmtav
Definition: modmain.f90:156
real(8) epsstress
Definition: modmain.f90:1059
real(8) chgexs
Definition: modmain.f90:722
real(8) demaxbnd
Definition: modmain.f90:821
integer nvp1d
Definition: modmain.f90:1106
real(8), dimension(3, maxatoms, maxspecies) atposc
Definition: modmain.f90:54
subroutine genspecies(fnum)
Definition: genspecies.f90:7
logical hdbse
Definition: modmain.f90:1201
integer, dimension(maxtasks) tasks
Definition: modmain.f90:1289
integer npp1d
Definition: modmain.f90:1108
logical molecule
Definition: modmain.f90:47
real(8) epsengy
Definition: modmain.f90:1055
real(8), dimension(3) efieldc
Definition: modmain.f90:312
integer dlmaxapw
Definition: modmain.f90:197
integer npfftgc
Definition: modmain.f90:408
integer nxoapwlo
Definition: modmain.f90:837
real(8), dimension(3) vklem
Definition: modmain.f90:479
real(8) hybridc
Definition: modmain.f90:1146
logical dosocc
Definition: modmain.f90:1075
logical trmt0
Definition: modmain.f90:165
integer maxscl
Definition: modmain.f90:1041
real(8) dgmaxvr
Definition: modmain.f90:384
logical spincore
Definition: modmain.f90:935
integer npmae0
Definition: modmain.f90:302
integer isgkmax
Definition: modmain.f90:491
real(8) atdfc
Definition: modmain.f90:66
logical primcell
Definition: modmain.f90:49
integer atpopt
Definition: modmain.f90:999
real(8) fracinr
Definition: modmain.f90:209
real(8) deltast
Definition: modmain.f90:1014
integer fsmtype
Definition: modmain.f90:251
logical trdatdv
Definition: modmain.f90:62
real(8) gmaxvr
Definition: modmain.f90:384
integer lmaxi
Definition: modmain.f90:205
logical bfdmag
Definition: modmain.f90:236