TABLE OF CONTENTS


ABINIT/m_anaddb_dataset [ Modules ]

[ Top ] [ Modules ]

NAME

  m_anaddb_dataset

FUNCTION

COPYRIGHT

  Copyright (C) 2014-2018 ABINIT group (XG,JCC,CL,MVeithen,XW,MJV)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

CHILDREN

SOURCE

20 #if defined HAVE_CONFIG_H
21 #include "config.h"
22 #endif
23 
24 #include "abi_common.h"
25 
26 module m_anaddb_dataset
27 
28  use defs_basis
29  use m_abicore
30  use m_errors
31 
32  use m_fstrings,  only : next_token, rmquotes, sjoin, inupper
33  use m_symtk,     only : mati3det
34  use m_parser,    only : intagm, chkvars_in_string
35  use m_ddb,       only : DDB_QTOL
36 
37  implicit none
38 
39  private
40 
41  public :: anaddb_dataset_type
42  public :: anaddb_dtset_free
43  public :: anaddb_init
44  public :: outvars_anaddb
45  public :: invars9

m_anaddb_dataset/anaddb_chkvars [ Functions ]

[ Top ] [ m_anaddb_dataset ] [ Functions ]

NAME

 anaddb_chkvars

FUNCTION

  Examines the input string, to check whether all names are allowed.

INPUTS

  string*(*)=string of character
   the string (with upper case) from the input file, to which the XYZ data is (possibly) appended

OUTPUT

PARENTS

      m_anaddb_dataset

CHILDREN

      chkvars_in_string,inupper

SOURCE

2168 subroutine anaddb_chkvars(string)
2169 
2170 
2171 !This section has been created automatically by the script Abilint (TD).
2172 !Do not modify the following lines by hand.
2173 #undef ABI_FUNC
2174 #define ABI_FUNC 'anaddb_chkvars'
2175 !End of the abilint section
2176 
2177  implicit none
2178 
2179 !Arguments ------------------------------------
2180 !scalars
2181  character(len=*),intent(in) :: string
2182 
2183 !Local variables-------------------------------
2184 !scalars
2185  integer,parameter :: protocol0=0
2186  character(len=100) :: list_logicals,list_strings
2187  character(len=10000) :: list_vars
2188 
2189 !************************************************************************
2190 
2191 !Here, list all admitted variable names (max 10 per line, to fix the ideas)
2192 !Note: Do not use "double quotation mark" for the string since it triggers a bug in docchk.py (abirules script)
2193 !<ANADDB_VARS>
2194 !A
2195  list_vars=                 ' alphon asr a2fsmear atifc'
2196 !B
2197  list_vars=trim(list_vars)//' brav band_gap'
2198 !C
2199  list_vars=trim(list_vars)//' chneut'
2200 !D
2201  list_vars=trim(list_vars)//' dieflag dipdip dossum dosdeltae dossmear dostol'
2202 !E
2203  list_vars=trim(list_vars)//' ep_scalprod eivec elaflag elphflag enunit'
2204  list_vars=trim(list_vars)//' ep_b_min ep_b_max ep_int_gkk ep_keepbands ep_nqpt ep_nspline ep_prt_yambo'
2205  list_vars=trim(list_vars)//' elphsmear elph_fermie ep_extrael ep_qptlist'
2206 !F
2207  list_vars=trim(list_vars)//' freeze_displ frmax frmin'
2208 !G
2209  list_vars=trim(list_vars)//' gkk2write gkk_rptwrite gkqwrite gruns_nddbs'
2210 !H
2211 !I
2212  list_vars=trim(list_vars)//' ifcana ifcflag ifcout ifltransport instrflag istrfix iatfix iatprj_bs'
2213 !J
2214 !K
2215  list_vars=trim(list_vars)//' kptrlatt kptrlatt_fine'
2216 !L
2217 !M
2218  list_vars=trim(list_vars)//' mustar'
2219 !N
2220  list_vars=trim(list_vars)//' natfix natifc natom natprj_bs nchan ndivsm nfreq ngrids nlflag nph1l nph2l'
2221  list_vars=trim(list_vars)//' nqpath nqshft nsphere nstrfix ntemper nwchan ngqpt ng2qpt'
2222 !O
2223  list_vars=trim(list_vars)//' outboltztrap'
2224 !P
2225  list_vars=trim(list_vars)//' piezoflag polflag prtddb prtdos prt_ifc prtmbm prtfsurf'
2226  list_vars=trim(list_vars)//' prtnest prtphbands prtsrlr prtvol prtbltztrp'
2227 !Q
2228  list_vars=trim(list_vars)//' qrefine qgrid_type q1shft q2shft qnrml1 qnrml2 qpath qph1l qph2l'
2229 !R
2230  list_vars=trim(list_vars)//' ramansr relaxat relaxstr rfmeth rifcsph'
2231 !S
2232  list_vars=trim(list_vars)//' selectz symdynmat symgkq'
2233 !T
2234  list_vars=trim(list_vars)//' targetpol telphint thmflag temperinc tempermin thermal_supercell thmtol'
2235 !U
2236  list_vars=trim(list_vars)//' use_k_fine'
2237 !V
2238  list_vars=trim(list_vars)//' vs_qrad_tolkms'
2239 !W
2240 !X
2241 !Y
2242 !Z
2243 
2244 !Logical input variables
2245  list_logicals=' '
2246 
2247 !String input variables
2248  list_strings=' gruns_ddbs'
2249 !</ANADDB_VARS>
2250 
2251 !Extra token, also admitted:
2252 !<ANADDB_UNITS>
2253  list_vars=trim(list_vars)//' au Angstr Angstrom Angstroms Bohr Bohrs eV Ha'
2254  list_vars=trim(list_vars)//' Hartree Hartrees K Ry Rydberg Rydbergs T Tesla'
2255 !</ANADDB_UNITS>
2256 
2257 !<ANADDB_OPERATORS>
2258  list_vars=trim(list_vars)//' sqrt '
2259 !</ANADDB_OPERATORS>
2260 
2261 !Transform to upper case
2262  call inupper(list_vars)
2263  call inupper(list_logicals)
2264  call inupper(list_strings)
2265 
2266  call chkvars_in_string(protocol0, list_vars, list_logicals, list_strings, string)
2267 
2268 end subroutine anaddb_chkvars

m_anaddb_dataset/anaddb_dataset_type [ Types ]

[ Top ] [ m_anaddb_dataset ] [ Types ]

NAME

 anaddb_dataset_type

FUNCTION

 The anaddb_dataset_type structured datatype
 gather all the input variables for the anaddb code.

SOURCE

 60  type anaddb_dataset_type
 61 
 62 ! Variables should be declared on separated lines in order to reduce the occurence of git conflicts.
 63 ! Since all these input variables are described in the anaddb_help.html
 64 ! file, they are not described in length here ...
 65 ! Integer
 66   integer :: alphon
 67   integer :: asr
 68   integer :: brav
 69   integer :: chneut
 70   integer :: dieflag
 71   integer :: dipdip
 72   integer :: dossum
 73   integer :: ep_scalprod
 74   integer :: eivec
 75   integer :: elaflag
 76   integer :: elphflag
 77   integer :: enunit
 78   integer :: gkk2write
 79   integer :: gkk_rptwrite
 80   integer :: gkqwrite
 81   integer :: gruns_nddbs
 82   integer :: ifcana
 83   integer :: ifcflag
 84   integer :: ifcout
 85   integer :: ifltransport
 86   integer :: instrflag
 87   integer :: natfix
 88   integer :: natifc
 89   integer :: natom
 90   integer :: natprj_bs
 91   integer :: nchan
 92   integer :: ndivsm=20
 93   integer :: nfreq
 94   integer :: ngrids
 95   integer :: nlflag
 96   integer :: nph1l
 97   integer :: nph2l
 98   integer :: nqpath
 99   integer :: nqshft
100   integer :: nsphere
101   integer :: nstrfix
102   integer :: ntemper
103   integer :: nwchan
104   integer :: outboltztrap
105   integer :: piezoflag
106   integer :: polflag
107   integer :: prtdos
108   integer :: prt_ifc
109   integer :: prtddb
110   integer :: prtmbm
111   integer :: prtfsurf
112   integer :: prtnest
113   integer :: prtphbands
114   integer :: prtsrlr  ! print the short-range/long-range decomposition of phonon freq.
115   integer :: prtvol = 0
116   integer :: ramansr
117   integer :: relaxat
118   integer :: relaxstr
119   integer :: rfmeth
120   integer :: selectz
121   integer :: symdynmat
122   integer :: telphint
123   integer :: thmflag
124   integer :: qgrid_type
125   integer :: ep_b_min
126   integer :: ep_b_max
127   integer :: ep_int_gkk
128   integer :: ep_keepbands
129   integer :: ep_nqpt
130   integer :: ep_nspline
131   integer :: ep_prt_yambo
132   integer :: symgkq
133   integer :: use_k_fine
134   integer :: prtbltztrp
135 
136   integer :: ngqpt(9)             ! ngqpt(9) instead of ngqpt(3) is needed in wght9.f
137   integer :: istrfix(6)
138   integer :: ng2qpt(3)
139   integer :: qrefine(3)
140   integer :: kptrlatt(3,3)
141   integer :: kptrlatt_fine(3,3)
142   integer :: thermal_supercell(3,3)
143 
144 ! Real(dp)
145   real(dp) :: a2fsmear
146   real(dp) :: band_gap
147   real(dp) :: dosdeltae
148   real(dp) :: dossmear
149   real(dp) :: dostol
150   real(dp) :: elphsmear
151   real(dp) :: elph_fermie
152   real(dp) :: ep_extrael
153   real(dp) :: freeze_displ
154   real(dp) :: frmax
155   real(dp) :: frmin
156   real(dp) :: temperinc
157   real(dp) :: tempermin
158   real(dp) :: thmtol
159   real(dp) :: mustar
160   real(dp) :: rifcsph
161 
162   real(dp) :: q1shft(3,4)
163   real(dp) :: q2shft(3)
164   real(dp) :: targetpol(3)
165   real(dp) :: vs_qrad_tolkms(2) = 0
166 
167 ! Integer arrays
168   integer, allocatable :: atifc(:)
169    ! atifc(natom) WARNING : there is a transformation of this input variable, in chkin9
170    ! This should be changed ...
171 
172   integer, allocatable :: iatfix(:)
173   ! iatfix(natom)
174 
175   integer, allocatable :: iatprj_bs(:)
176 
177 ! Real arrays
178   real(dp), allocatable :: qnrml1(:)
179   ! qnrml1(nph1l)
180 
181   real(dp), allocatable :: qnrml2(:)
182   ! qnrml2(nph2l)
183 
184   real(dp), allocatable :: qpath(:,:)
185   ! qpath(3,nqpath)
186 
187   real(dp), allocatable :: qph1l(:,:)
188   ! qph1l(3,nph1l)
189 
190   real(dp), allocatable :: qph2l(:,:)
191   ! qph2l(3,nph2l)
192 
193   real(dp), allocatable :: ep_qptlist(:,:)
194   ! qph2l(3,ep_nqpt)
195 
196   character(len=fnlen), allocatable :: gruns_ddbs(:)
197   ! gruns_ddbs(gruns_nddbs)
198 
199  end type anaddb_dataset_type

m_anaddb_dataset/anaddb_dtset_free [ Functions ]

[ Top ] [ m_anaddb_dataset ] [ Functions ]

NAME

   anaddb_dtset_free

FUNCTION

   deallocate remaining arrays in the anaddb_dtset datastructure

INPUTS

  anaddb_dtset = anaddb datastructure

PARENTS

      anaddb

CHILDREN

      chkvars_in_string,inupper

NOTES

SOURCE

226 subroutine anaddb_dtset_free(anaddb_dtset)
227 
228 
229 !This section has been created automatically by the script Abilint (TD).
230 !Do not modify the following lines by hand.
231 #undef ABI_FUNC
232 #define ABI_FUNC 'anaddb_dtset_free'
233 !End of the abilint section
234 
235  implicit none
236 
237 !Arguments ------------------------------------
238 !scalars
239  type(anaddb_dataset_type), intent(inout) :: anaddb_dtset
240 
241 ! *************************************************************************
242 
243  if (allocated(anaddb_dtset%atifc))  then
244    ABI_DEALLOCATE(anaddb_dtset%atifc)
245  end if
246  if (allocated(anaddb_dtset%iatfix))  then
247    ABI_DEALLOCATE(anaddb_dtset%iatfix)
248  end if
249  if (allocated(anaddb_dtset%iatprj_bs))  then
250    ABI_DEALLOCATE(anaddb_dtset%iatprj_bs)
251  end if
252  if (allocated(anaddb_dtset%qnrml1))  then
253    ABI_DEALLOCATE(anaddb_dtset%qnrml1)
254  end if
255  if (allocated(anaddb_dtset%qnrml2))  then
256    ABI_DEALLOCATE(anaddb_dtset%qnrml2)
257  end if
258  if (allocated(anaddb_dtset%qpath))  then
259    ABI_DEALLOCATE(anaddb_dtset%qpath)
260  end if
261  if (allocated(anaddb_dtset%qph1l))  then
262    ABI_DEALLOCATE(anaddb_dtset%qph1l)
263  end if
264  if (allocated(anaddb_dtset%qph2l))  then
265    ABI_DEALLOCATE(anaddb_dtset%qph2l)
266  end if
267  if (allocated(anaddb_dtset%ep_qptlist))  then
268    ABI_DEALLOCATE(anaddb_dtset%ep_qptlist)
269  end if
270  if (allocated(anaddb_dtset%gruns_ddbs))  then
271    ABI_DEALLOCATE(anaddb_dtset%gruns_ddbs)
272  end if
273 
274 end subroutine anaddb_dtset_free

m_anaddb_dataset/anaddb_init [ Functions ]

[ Top ] [ m_anaddb_dataset ] [ Functions ]

NAME

 anaddb_init

FUNCTION

 Initialize the code ppddb9: write heading and make the first i/os

INPUTS

OUTPUT

 character(len=fnlen) filnam(7)=character strings giving file names

NOTES

 1. Should be executed by one processor only.
 2. File names refer to following files, in order:
     (1) Formatted input file
     (2) Formatted output file
     (3) Input Derivative Database
     (4) Output Molecular Dynamics
     (5) Input electron-phonon matrix elements
     (6) Root name for electron-phonon file names
     (7) Name of file containing the 3 ddk filenames and the GS wf file name

PARENTS

      anaddb

CHILDREN

      chkvars_in_string,inupper

SOURCE

2104 subroutine anaddb_init(filnam)
2105 
2106 
2107 !This section has been created automatically by the script Abilint (TD).
2108 !Do not modify the following lines by hand.
2109 #undef ABI_FUNC
2110 #define ABI_FUNC 'anaddb_init'
2111 !End of the abilint section
2112 
2113  implicit none
2114 
2115 !Arguments -------------------------------
2116 !arrays
2117  character(len=*),intent(out) :: filnam(7)
2118 
2119 ! *********************************************************************
2120 
2121 !Read the file names
2122  write(std_out,*)' Give name for formatted input file: '
2123  read(std_in, '(a)' ) filnam(1)
2124  write(std_out,'(a,a)' )'-   ',trim(filnam(1))
2125  write(std_out,*)' Give name for formatted output file: '
2126  read(std_in, '(a)' ) filnam(2)
2127  write(std_out,'(a,a)' )'-   ',trim(filnam(2))
2128  write(std_out,*)' Give name for input derivative database: '
2129  read(std_in, '(a)' ) filnam(3)
2130  write(std_out,'(a,a)' )'-   ',trim(filnam(3))
2131  write(std_out,*)' Give name for output molecular dynamics: '
2132  read(std_in, '(a)' ) filnam(4)
2133  write(std_out,'(a,a)' )'-   ',trim(filnam(4))
2134  write(std_out,*)' Give name for input elphon matrix elements (GKK file): '
2135  read(std_in, '(a)' ) filnam(5)
2136  write(std_out,'(a,a)' )'-   ',trim(filnam(5))
2137  write(std_out,*)' Give root name for elphon output files: '
2138  read(std_in, '(a)' ) filnam(6)
2139  write(std_out,'(a,a)' )'-   ',trim(filnam(6))
2140  write(std_out,*)' Give name for file containing ddk filenames for elphon/transport: '
2141  read(std_in, '(a)' ) filnam(7)
2142  write(std_out,'(a,a)' )'-   ',trim(filnam(7))
2143 
2144 end subroutine anaddb_init

m_anaddb_dataset/invars9 [ Functions ]

[ Top ] [ m_anaddb_dataset ] [ Functions ]

NAME

 invars9

FUNCTION

 Open input file for the anaddb code, then reads or echoes the input information.

INPUTS

 lenstr=actual length of string
 natom=number of atoms, needed for atifc
 string*(*)=string of characters containing all input variables and data

OUTPUT

 anaddb_dtset= (derived datatype) contains all the input variables

NOTES

 Should be executed by one processor only.

 27/01/2009: MJV: I have cleaned this routine extensively, putting all
  variables in alphabetical order, and in a second segment the dependent
  variables which need to be allocated depending on the dimensions read in.
  Could be divided into two routines as in abinit.
    FIXME: move checks to chkin9?

PARENTS

      anaddb

CHILDREN

      chkvars_in_string,inupper

SOURCE

 312 subroutine invars9 (anaddb_dtset,lenstr,natom,string)
 313 
 314 
 315 !This section has been created automatically by the script Abilint (TD).
 316 !Do not modify the following lines by hand.
 317 #undef ABI_FUNC
 318 #define ABI_FUNC 'invars9'
 319 !End of the abilint section
 320 
 321  implicit none
 322 
 323 !Arguments -------------------------------
 324 !scalars
 325  integer,intent(in) :: lenstr,natom
 326  character(len=*),intent(in) :: string
 327  type(anaddb_dataset_type),intent(inout) :: anaddb_dtset
 328 
 329 !Local variables -------------------------
 330 !scalars
 331  integer,parameter :: vrsddb=100401 !Set routine version number here:
 332  integer,parameter :: jdtset=1
 333  integer :: ii,iph1,iph2,marr,tread,start
 334  integer :: idet
 335  character(len=500) :: message
 336  character(len=fnlen) :: path
 337 !arrays
 338  integer,allocatable :: intarr(:)
 339  real(dp),allocatable :: dprarr(:)
 340 
 341 !*********************************************************************
 342  marr=3
 343  ABI_ALLOCATE(intarr,(marr))
 344  ABI_ALLOCATE(dprarr,(marr))
 345 
 346 !copy natom to anaddb_dtset
 347  anaddb_dtset%natom=natom
 348 
 349 !=====================================================================
 350 !start reading in dimensions and non-dependent variables
 351 !=====================================================================
 352 
 353 !A
 354 
 355 !typical value for gaussian smearing of a2F function
 356  anaddb_dtset%a2fsmear = 0.00002_dp
 357  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'a2fsmear',tread,'ENE')
 358  if(tread==1) anaddb_dtset%a2fsmear=dprarr(1)
 359  if (anaddb_dtset%a2fsmear < tol6) then
 360    write(message,'(a,f10.3,a,a,a,a,a)' )&
 361 &   'a2fsmear is',anaddb_dtset%a2fsmear,', but only values > 1.e-6 ',ch10,&
 362 &   'are allowed',ch10,'Action: correct a2fsmear in your input file.'
 363    MSG_ERROR(message)
 364  end if
 365 
 366  anaddb_dtset%alphon=0
 367  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'alphon',tread,'INT')
 368  if(tread==1) anaddb_dtset%alphon=intarr(1)
 369 !FIXME: need a test on input value
 370 
 371  anaddb_dtset%asr=1
 372  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'asr',tread,'INT')
 373  if(tread==1) anaddb_dtset%asr=intarr(1)
 374  if(anaddb_dtset%asr<-2.or.anaddb_dtset%asr>5)then
 375    write(message, '(a,i0,5a)' )&
 376 &   'asr is ',anaddb_dtset%asr,', but the only allowed values',ch10,&
 377 &   'are 0, 1, 2, 3, 4, 5, -1 or -2 .',ch10,'Action: correct asr in your input file.'
 378 !  Note : negative values are allowed when the acoustic sum rule
 379 !  is to be applied after the analysis of IFCs
 380 !  3,4 are for rotational invariance (under development)
 381 !  5 is for hermitian imposition of the ASR
 382    MSG_ERROR(message)
 383  end if
 384 
 385 !B
 386 
 387 !Target band gap in eV
 388 !The default value is just a very large number that will not be used in changing the band gap
 389  anaddb_dtset%band_gap = 999.0d0
 390  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'band_gap',tread,'DPR')
 391  if(tread==1) anaddb_dtset%band_gap=dprarr(1)
 392 
 393  anaddb_dtset%brav=1
 394  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'brav',tread,'INT')
 395  if(tread==1) anaddb_dtset%brav=intarr(1)
 396  if(anaddb_dtset%brav<=-2.or.anaddb_dtset%brav>=5 .or. anaddb_dtset%brav==0)then
 397    write(message, '(a,i0,a5)' )&
 398 &   'brav is ',anaddb_dtset%brav,', but the only allowed values',ch10,&
 399 &   'are -1, 1,2,3 or 4 .',ch10,'Action: correct brav in your input file.'
 400    MSG_ERROR(message)
 401  end if
 402 
 403 !C
 404 
 405  anaddb_dtset%chneut=0
 406  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'chneut',tread,'INT')
 407  if(tread==1) anaddb_dtset%chneut=intarr(1)
 408  if(anaddb_dtset%chneut<0.or.anaddb_dtset%chneut>2)then
 409    write(message, '(a,i0,5a)' )&
 410 &   'chneut is ',anaddb_dtset%chneut,', but the only allowed values',ch10,&
 411 &   'are 0, 1 or 2.',ch10,'Action: correct chneut in your input file.'
 412    MSG_ERROR(message)
 413  end if
 414 
 415 !D
 416 
 417  anaddb_dtset%dieflag=0
 418  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dieflag',tread,'INT')
 419  if(tread==1) anaddb_dtset%dieflag=intarr(1)
 420  if(anaddb_dtset%dieflag<0.or.anaddb_dtset%dieflag>4)then
 421    write(message, '(a,i0,5a)' )&
 422 &   'dieflag is ',anaddb_dtset%dieflag,', but the only allowed values',ch10,&
 423 &   'are 0, 1, 2, 3 or 4.',ch10,'Action: correct dieflag in your input file.'
 424    MSG_ERROR(message)
 425  end if
 426 
 427  anaddb_dtset%dipdip=1
 428  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dipdip',tread,'INT')
 429  if(tread==1) anaddb_dtset%dipdip=intarr(1)
 430  if(anaddb_dtset%dipdip<0.or.anaddb_dtset%dipdip>1)then
 431    write(message, '(a,i0,5a)' )&
 432 &   'dipdip is ',anaddb_dtset%dipdip,', but the only allowed values',ch10,&
 433 &   'are 0 or 1 .',ch10,'Action: correct dipdip in your input file.'
 434    MSG_ERROR(message)
 435  end if
 436 
 437  anaddb_dtset%ep_scalprod = 0
 438  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ep_scalprod',tread,'INT')
 439  if(tread==1) anaddb_dtset%ep_scalprod = intarr(1)
 440  if(anaddb_dtset%ep_scalprod < 0 .or. anaddb_dtset%ep_scalprod > 1) then
 441    write(message, '(a,i0,5a)' )&
 442 &   'ep_scalprod is ',anaddb_dtset%ep_scalprod,', but the only allowed values',ch10,&
 443 &   'are 0 or 1.',ch10,'Action: correct ep_scalprod in your input file.'
 444    MSG_ERROR(message)
 445  end if
 446 
 447  anaddb_dtset%dosdeltae=one/Ha_cmm1
 448  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dosdeltae',tread,'DPR')
 449  if(tread==1) anaddb_dtset%dosdeltae=dprarr(1)
 450  if(anaddb_dtset%dosdeltae<=zero)then
 451    write(message, '(a,es14.4,a,a,a)' )&
 452 &   'dosdeltae is ',anaddb_dtset%dosdeltae,', which is lower than 0 .',ch10,&
 453 &   'Action: correct dosdeltae in your input file.'
 454    MSG_ERROR(message)
 455  end if
 456 
 457 !FIXME : should probably be smaller
 458  anaddb_dtset%dossmear=5.0/Ha_cmm1
 459  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dossmear',tread,'DPR')
 460  if(tread==1) anaddb_dtset%dossmear=dprarr(1)
 461  if(anaddb_dtset%dossmear<=zero)then
 462    write(message, '(a,es14.4,3a)' )&
 463 &   'dossmear is ',anaddb_dtset%dossmear,', which is lower than 0 .',ch10,&
 464 &   'Action: correct dossmear in your input file.'
 465    MSG_ERROR(message)
 466  end if
 467 
 468  anaddb_dtset%dostol=0.25_dp
 469  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dostol',tread,'DPR')
 470  if(tread==1) anaddb_dtset%dostol=dprarr(1)
 471  if(anaddb_dtset%dostol<zero)then
 472    write(message, '(a,es14.4,3a)' )&
 473 &   'dostol is ',anaddb_dtset%dostol,', which is lower than 0 .',ch10,&
 474 &   'Action: correct dostol in your input file.'
 475    MSG_ERROR(message)
 476  end if
 477 
 478  anaddb_dtset%dossum=0
 479  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dossum',tread,'INT')
 480  if(tread==1) anaddb_dtset%dossum=intarr(1)
 481  if(anaddb_dtset%dossum < 0 .or. anaddb_dtset%dossum > one)then
 482    write(message, '(a,i0,5a)' )&
 483 &   'dossum is ',anaddb_dtset%dossum,', but the only allowed values',ch10,&
 484 &   'are 0, 1',ch10,'Action: correct dossum in your input file.'
 485    MSG_ERROR(message)
 486  end if
 487 
 488 !E
 489 
 490  anaddb_dtset%eivec=0
 491  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'eivec',tread,'INT')
 492  if(tread==1) anaddb_dtset%eivec=intarr(1)
 493  if(anaddb_dtset%eivec<0.or.anaddb_dtset%eivec>4)then
 494    write(message, '(a,i0,5a)' )&
 495 &   'eivec is ',anaddb_dtset%eivec,', but the only allowed values',ch10,&
 496 &   'are 0, 1, 2, 3 or 4.',ch10,'Action: correct eivec in your input file.'
 497    MSG_ERROR(message)
 498  end if
 499 
 500  anaddb_dtset%elaflag=0
 501  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'elaflag',tread,'INT')
 502  if(tread==1) anaddb_dtset%elaflag=intarr(1)
 503  if(anaddb_dtset%elaflag<0.or.anaddb_dtset%elaflag>5)then
 504    write(message,'(a,i0,5a)' )&
 505 &   'elaflag is ',anaddb_dtset%elaflag,', but the only allowed values',ch10,&
 506 &   'are 0,1,2,3,4 or 5 .',ch10,'Action: correct elaflag in your input file.'
 507    MSG_ERROR(message)
 508  end if
 509 
 510 !By default use the real fermie (tests for abs(elph_fermie) < tol10 in the code)
 511  anaddb_dtset%elph_fermie = zero
 512  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'elph_fermie',tread,'ENE')
 513  if(tread==1) anaddb_dtset%elph_fermie=dprarr(1)
 514 
 515 !extra charge in unit cell (number of electrons) wrt neutral cell
 516 !holes are negative values (reduce number of electrons)
 517  anaddb_dtset%ep_extrael = zero
 518  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ep_extrael',tread,'DPR')
 519  if(tread==1) anaddb_dtset%ep_extrael=dprarr(1)
 520 
 521 !number to control the spline interpolation in RTA
 522  anaddb_dtset%ep_nspline = 20
 523  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ep_nspline',tread,'INT')
 524  if(tread==1) anaddb_dtset%ep_nspline=intarr(1)
 525  if(anaddb_dtset%ep_nspline < 0 .or. anaddb_dtset%ep_nspline > 1000) then
 526    write(message, '(a,i0,5a)' )&
 527 &   'ep_nspline is ',anaddb_dtset%ep_nspline,', but this should not be ',ch10,&
 528 &   'negative or too large .',ch10,'Action: correct ep_nspline in your input file.'
 529    MSG_ERROR(message)
 530  end if
 531 
 532 !interpolate gkk or gamma. It should be better to interpolate gkk onto the
 533 !k_phon, since the integration weights will be treated the same way
 534  anaddb_dtset%ep_int_gkk = 0
 535  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ep_int_gkk',tread,'INT')
 536  if(tread==1) anaddb_dtset%ep_int_gkk = intarr(1)
 537  if(anaddb_dtset%ep_int_gkk < 0 .or. anaddb_dtset%ep_int_gkk > 1) then
 538    write(message, '(a,i0,5a)' )&
 539 &   'ep_int_gkk is ',anaddb_dtset%ep_int_gkk,', but the only allowed values',ch10,&
 540 &   'are 0 or 1.',ch10,'Action: correct ep_int_gkk in your input file.'
 541    MSG_ERROR(message)
 542  end if
 543 
 544  anaddb_dtset%elphflag=0
 545  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'elphflag',tread,'INT')
 546  if(tread==1) anaddb_dtset%elphflag=intarr(1)
 547  if(anaddb_dtset%elphflag<0.or.anaddb_dtset%elphflag>1)then
 548    write(message,'(a,i0,5a)' )&
 549 &   'elphflag = ',anaddb_dtset%elphflag,', but the allowed values',ch10,&
 550 &   'are 0, or 1.',ch10,'Action: correct elphflag in your input file.'
 551    MSG_ERROR(message)
 552  end if
 553 
 554 !typical value for gaussian smearing, but can vary sensibly with the metal
 555  anaddb_dtset%elphsmear = 0.01_dp
 556  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'elphsmear',tread,'ENE')
 557  if(tread==1) anaddb_dtset%elphsmear=dprarr(1)
 558  if (anaddb_dtset%elphsmear < tol6) then
 559    write(message,'(a,f10.3,5a)' )&
 560 &   'elphsmear is ',anaddb_dtset%elphsmear,'. Only values > 1.e-6 ',ch10,&
 561 &   'are allowed',ch10,'Action: correct elphsmear in your input file.'
 562    MSG_ERROR(message)
 563  end if
 564 
 565  anaddb_dtset%enunit=0
 566  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'enunit',tread,'INT')
 567  if(tread==1) anaddb_dtset%enunit=intarr(1)
 568  if(anaddb_dtset%enunit<0.or.anaddb_dtset%enunit>2)then
 569    write(message, '(a,i0,5a)' )&
 570 &   'enunit is ',anaddb_dtset%enunit,', but the only allowed values',ch10,&
 571 &   'are 0, 1 or 2.',ch10,'Action: correct enunit in your input file.'
 572    MSG_ERROR(message)
 573  end if
 574 
 575 !Default is 0 - not used unless telphint==2
 576  anaddb_dtset%ep_b_max = 0
 577  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ep_b_max',tread,'INT')
 578  if(tread==1) then
 579    anaddb_dtset%ep_b_max = intarr(1)
 580    if(anaddb_dtset%ep_b_max < 1) then
 581      write(message, '(a,i0,5a)' )&
 582 &     'ep_b_max is ',anaddb_dtset%ep_b_max,', but the only allowed values',ch10,&
 583 &     'are between 1 and nband.',ch10,'Action: correct ep_b_max in your input file.'
 584      MSG_ERROR(message)
 585    end if
 586  end if
 587 
 588 !Default is 0 - not used unless telphint==2
 589  anaddb_dtset%ep_b_min = 0
 590  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ep_b_min',tread,'INT')
 591  if(tread==1) then
 592    anaddb_dtset%ep_b_min = intarr(1)
 593    if(anaddb_dtset%ep_b_min < 1) then
 594      write(message, '(a,i0,5a)' )&
 595 &     'ep_b_min is ',anaddb_dtset%ep_b_min,', but the only allowed values',ch10,&
 596 &     'are between 1 and nband.',ch10,'Action: correct ep_b_min in your input file.'
 597      MSG_ERROR(message)
 598    end if
 599  end if
 600 
 601  anaddb_dtset%ep_keepbands = 0
 602  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ep_keepbands',tread,'INT')
 603  if(tread==1) anaddb_dtset%ep_keepbands = intarr(1)
 604  if(anaddb_dtset%ep_keepbands < 0 .or. anaddb_dtset%ep_keepbands > 1) then
 605    write(message, '(a,i0,5a)' )&
 606 &   'ep_keepbands is ',anaddb_dtset%ep_keepbands,', but the only allowed values',ch10,&
 607 &   'are 0 or 1 .',ch10,'Action: correct ep_keepbands in your input file.'
 608    MSG_ERROR(message)
 609  end if
 610 
 611  anaddb_dtset%ep_nqpt=0
 612  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ep_nqpt',tread,'INT')
 613  if(tread==1) anaddb_dtset%ep_nqpt = intarr(1)
 614  if(anaddb_dtset%ep_nqpt < 0) then
 615    write(message, '(a,i0,5a)' )&
 616 &   'ep_nqpt is ',anaddb_dtset%ep_nqpt,', but the only allowed values',ch10,&
 617 &   'are > 0.',ch10,'Action: correct ep_nqpt in your input file.'
 618    MSG_ERROR(message)
 619  end if
 620 
 621  if (anaddb_dtset%ep_nqpt > 0) then
 622    ABI_ALLOCATE(anaddb_dtset%ep_qptlist,(3,anaddb_dtset%ep_nqpt))
 623    if(3*anaddb_dtset%ep_nqpt>marr)then
 624      marr=3*anaddb_dtset%ep_nqpt
 625      ABI_DEALLOCATE(intarr)
 626      ABI_DEALLOCATE(dprarr)
 627      ABI_ALLOCATE(intarr,(marr))
 628      ABI_ALLOCATE(dprarr,(marr))
 629    end if
 630    anaddb_dtset%ep_qptlist(:,:)=zero
 631    call intagm(dprarr,intarr,jdtset,marr,3*anaddb_dtset%ep_nqpt,string(1:lenstr),'ep_qptlist',tread,'DPR')
 632    if(tread==1) then
 633      anaddb_dtset%ep_qptlist(1:3,1:anaddb_dtset%ep_nqpt)=&
 634 &     reshape(dprarr(1:3*anaddb_dtset%ep_nqpt),(/3,anaddb_dtset%ep_nqpt/))
 635    else
 636      write(message,'(3a)')&
 637 &     'ep_nqpt is non zero but ep_qptlist is absent ',ch10,&
 638 &     'Action: specify ep_qptlist in your input file.'
 639      MSG_ERROR(message)
 640    end if
 641  end if
 642 
 643 
 644 !F
 645  anaddb_dtset%freeze_displ = zero
 646  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'freeze_displ',tread,'DPR')
 647  if(tread==1) anaddb_dtset%freeze_displ=dprarr(1)
 648 
 649 
 650  anaddb_dtset%frmax=ten
 651  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'frmax',tread,'DPR')
 652  if(tread==1) anaddb_dtset%frmax=dprarr(1)
 653  if (anaddb_dtset%frmax < 0) then
 654    write(message,'(a,f10.3,5a)' )&
 655 &   'frmax is ',anaddb_dtset%frmax,'. Only values > 0 ',ch10,&
 656 &   'are allowed',ch10,'Action: correct frmax in your input file.'
 657    MSG_ERROR(message)
 658  end if
 659 
 660  anaddb_dtset%frmin=zero
 661  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'frmin',tread,'DPR')
 662  if(tread==1) anaddb_dtset%frmin=dprarr(1)
 663  if (anaddb_dtset%frmin < 0) then
 664    write(message,'(a,f10.3,5a)' )&
 665 &   'frmin is ',anaddb_dtset%frmin,'. Only values > 0 ',ch10,&
 666 &   'are allowed',ch10,'Action: correct frmin in your input file.'
 667    MSG_ERROR(message)
 668  end if
 669 
 670 !G
 671 
 672  anaddb_dtset%gkk2write = 0
 673 !call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gkk2write',&
 674 !& tread,'INT')
 675 !if(tread==1) anaddb_dtset%gkk2write = intarr(1)
 676 !if(anaddb_dtset%gkk2write < 0 .or. anaddb_dtset%gkk2write > 1) then
 677 !write(message, '(a,a,a,i8,a,a,a,a,a)' )&
 678 !&  ' invars9 : ERROR -',ch10,&
 679 !&  '  gkk2write is',anaddb_dtset%gkk2write,&
 680 !&  ', but the only allowed values',ch10,&
 681 !&  '  are 0 or 1 .',ch10,&
 682 !&  '  Action: correct gkk2write in your input file.'
 683 !MSG_ERROR(message)
 684 !end if
 685 
 686  anaddb_dtset%gkk_rptwrite = 0
 687 !call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gkk_rptwrite',&
 688 !& tread,'INT')
 689 !if(tread==1) anaddb_dtset%gkk_rptwrite = intarr(1)
 690 !if(anaddb_dtset%gkk_rptwrite < 0 .or. anaddb_dtset%gkk_rptwrite > 1) then
 691 !write(message, '(a,a,a,i8,a,a,a,a,a)' )&
 692 !&  ' invars9 : ERROR -',ch10,&
 693 !&  '  gkk_rptwrite is',anaddb_dtset%gkk_rptwrite,&
 694 !&  ', but the only allowed values',ch10,&
 695 !&  '  are 0 or 1 .',ch10,&
 696 !&  '  Action: correct gkk_rptwrite in your input file.'
 697 !MSG_ERROR(message)
 698 !end if
 699 
 700  anaddb_dtset%gkqwrite = 0
 701  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gkqwrite',tread,'INT')
 702  if(tread==1) anaddb_dtset%gkqwrite = intarr(1)
 703  if(anaddb_dtset%gkqwrite < 0 .or. anaddb_dtset%gkqwrite > 1) then
 704    write(message, '(a,i0,5a)' )&
 705 &   'gkqwrite is ',anaddb_dtset%gkqwrite,', but the only allowed values',ch10,&
 706 &   'are 0 or 1.',ch10,'Action: correct gkqwrite in your input file.'
 707    MSG_ERROR(message)
 708  end if
 709 
 710  anaddb_dtset%gruns_nddbs = 0
 711  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gruns_nddbs',tread,'INT')
 712  if (tread==1) anaddb_dtset%gruns_nddbs = intarr(1)
 713 
 714  if (anaddb_dtset%gruns_nddbs /= 0) then
 715    ! Read list of DDB paths.
 716    ABI_MALLOC(anaddb_dtset%gruns_ddbs, (anaddb_dtset%gruns_nddbs))
 717    start = index(string, "GRUNS_DDBS") + len("GRUNS_DDBS") + 1
 718    do ii=1,anaddb_dtset%gruns_nddbs
 719      if (next_token(string, start, path) /= 0) then
 720        MSG_ERROR(sjoin("Cannot find DDB path in input string:", ch10, string(start:)))
 721      end if
 722      anaddb_dtset%gruns_ddbs(ii) = rmquotes(path)
 723    end do
 724  end if
 725 
 726 !H
 727 
 728 !I
 729 
 730  anaddb_dtset%ifcana=0
 731  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ifcana',tread,'INT')
 732  if(tread==1) anaddb_dtset%ifcana=intarr(1)
 733  if(anaddb_dtset%ifcana<0.or.anaddb_dtset%ifcana>1)then
 734    write(message, '(a,i0,5a)' )&
 735 &   'ifcana is ',anaddb_dtset%ifcana,', but the only allowed values',ch10,&
 736 &   'are 0 or 1.',ch10,'Action: correct ifcana in your input file.'
 737    MSG_ERROR(message)
 738  end if
 739 
 740  anaddb_dtset%ifcflag=0
 741  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ifcflag',tread,'INT')
 742  if(tread==1) anaddb_dtset%ifcflag=intarr(1)
 743  if(anaddb_dtset%ifcflag<0.or.anaddb_dtset%ifcflag>1)then
 744    write(message, '(a,i0,5a)' )&
 745 &   'ifcflag is ',anaddb_dtset%ifcflag,', but the only allowed values',ch10,&
 746 &   'are 0 or 1.',ch10,'Action: correct ifcflag in your input file.'
 747    MSG_ERROR(message)
 748  end if
 749 
 750  anaddb_dtset%ifcout=0
 751  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ifcout',tread,'INT')
 752  if(tread==1) anaddb_dtset%ifcout=intarr(1)
 753  if(anaddb_dtset%ifcout<-1)then
 754    write(message, '(a,i0,3a)' )&
 755 &   'ifcout is ',anaddb_dtset%ifcout,', which is lower than -1.',ch10,&
 756 &   'Action: correct ifcout in your input file.'
 757    MSG_ERROR(message)
 758  end if
 759 
 760  anaddb_dtset%ifltransport = 0
 761  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ifltransport',tread,'INT')
 762  if(tread==1) anaddb_dtset%ifltransport = intarr(1)
 763  if(anaddb_dtset%ifltransport < 0 .or. anaddb_dtset%ifltransport > 3) then
 764    write(message, '(a,i0,5a)' )&
 765 &   'ifltransport is ',anaddb_dtset%ifltransport,', but the only allowed values',ch10,&
 766 &   'are 0 or 1 or 2 or 3.',ch10,'Action: correct ifltransport in your input file.'
 767    MSG_ERROR(message)
 768  end if
 769 
 770  anaddb_dtset%instrflag=0
 771  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'instrflag',tread,'INT')
 772  if(tread==1) anaddb_dtset%instrflag=intarr(1)
 773  if(anaddb_dtset%instrflag<0.or.anaddb_dtset%instrflag>1)then
 774    write(message,'(a,i0,5a)' )&
 775 &   'instrflag is ',anaddb_dtset%instrflag,', but the only allowed values',ch10,&
 776 &   'are 0, 1.',ch10,'Action: correct instrflag in your input file.'
 777    MSG_ERROR(message)
 778  end if
 779 
 780 !J
 781 
 782 !K
 783 
 784  anaddb_dtset%kptrlatt = 0
 785 !why this test on reading in kptrlatt?
 786  marr = 9
 787  ABI_DEALLOCATE(intarr)
 788  ABI_DEALLOCATE(dprarr)
 789  ABI_ALLOCATE(intarr,(marr))
 790  ABI_ALLOCATE(dprarr,(marr))
 791  call intagm(dprarr,intarr,jdtset,marr,9,string(1:lenstr),'kptrlatt',tread,'INT')
 792  if(tread==1)anaddb_dtset%kptrlatt(1:3,1:3)=reshape(intarr(1:9),(/3,3/))
 793 !NOTE: no a priori way to test the validity of the integers in kptrlatt
 794 
 795  anaddb_dtset%kptrlatt_fine(:,:)=0
 796  marr = 9
 797  ABI_DEALLOCATE(intarr)
 798  ABI_DEALLOCATE(dprarr)
 799  ABI_ALLOCATE(intarr,(marr))
 800  ABI_ALLOCATE(dprarr,(marr))
 801  call intagm(dprarr,intarr,jdtset,marr,9,string(1:lenstr),'kptrlatt_fine',tread,'INT')
 802  if(tread==1)anaddb_dtset%kptrlatt_fine(1:3,1:3)=reshape(intarr(1:9),(/3,3/))
 803 
 804 
 805 !L
 806 
 807 !M
 808 
 809 !typical value for mustar, but can vary sensibly with the metal
 810  anaddb_dtset%mustar = 0.1_dp
 811  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'mustar',tread,'DPR')
 812  if(tread==1) anaddb_dtset%mustar=dprarr(1)
 813  if (anaddb_dtset%mustar < zero) then
 814    write(message,'(a,f10.3,5a)' )&
 815 &   'mustar is ',anaddb_dtset%mustar,', but only positive values',ch10,&
 816 &   'are allowed',ch10,'Action: correct mustar in your input file.'
 817    MSG_ERROR(message)
 818  end if
 819 
 820 !N
 821 
 822  anaddb_dtset%natfix=0
 823  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natfix',tread,'INT')
 824  if(tread==1) anaddb_dtset%natfix=intarr(1)
 825  if(anaddb_dtset%natfix > natom)then
 826    write(message, '(a,i0,2a,i0,3a)' )&
 827 &   'natfix is ',anaddb_dtset%natfix,', which is larger than natom',' (=',natom,')',ch10,&
 828 &   'Action: correct natfix in your input file.'
 829    MSG_ERROR(message)
 830  end if
 831 
 832  if(anaddb_dtset%natfix < 0)then
 833    write(message, '(a,i0,3a)' )&
 834 &   'natfix is ',anaddb_dtset%natfix,', which is < 0',ch10,&
 835 &   'Action: correct natfix in your input file.'
 836    MSG_ERROR(message)
 837  end if
 838 
 839  anaddb_dtset%natifc=0
 840  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natifc',tread,'INT')
 841  if(tread==1) anaddb_dtset%natifc=intarr(1)
 842  if(anaddb_dtset%natifc<0)then
 843    write(message, '(a,i0,3a)' )&
 844 &   'natifc is ',anaddb_dtset%natifc,', which is lower than 0 .',ch10,&
 845 &   'Action: correct natifc in your input file.'
 846    MSG_ERROR(message)
 847  end if
 848 
 849  anaddb_dtset%natprj_bs=0
 850  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natprj_bs',tread,'INT')
 851  if(tread==1) anaddb_dtset%natprj_bs=intarr(1)
 852  if(anaddb_dtset%natprj_bs<0 .or. anaddb_dtset%natprj_bs > natom)then
 853    write(message, '(a,i0,a,i0,2a)' )&
 854 &   'natprj_bs is ',anaddb_dtset%natprj_bs,', but must be between 0 and natom = ',natom,ch10,&
 855 &   'Action: correct natprj_bs in your input file.'
 856    MSG_ERROR(message)
 857  end if
 858 
 859  anaddb_dtset%nchan=800
 860  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nchan',tread,'INT')
 861  if(tread==1) anaddb_dtset%nchan=intarr(1)
 862 !FIXME: check this - it should probably be .ge. 1, not 0
 863  if(anaddb_dtset%nchan <0)then
 864    write(message, '(a,i0,3a)' )&
 865 &   'nchan is ',anaddb_dtset%nchan,', which is lower than 0 .',ch10,&
 866 &   'Action: correct nchan in your input file.'
 867    MSG_ERROR(message)
 868  end if
 869 
 870  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ndivsm',tread,'INT')
 871  if(tread==1) anaddb_dtset%ndivsm=intarr(1)
 872  if(anaddb_dtset%ndivsm <=0)then
 873    write(message, '(a,i0,3a)' )&
 874 &   'ndivsm is ',anaddb_dtset%ndivsm,', which is <= 0 .',ch10,&
 875 &   'Action: correct ndivsm in your input file.'
 876    MSG_ERROR(message)
 877  end if
 878 
 879  anaddb_dtset%nfreq=1
 880  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nfreq',tread,'INT')
 881  if(tread==1) anaddb_dtset%nfreq=intarr(1)
 882  if(anaddb_dtset%nfreq<0)then
 883    write(message, '(a,i0,3a)' )&
 884 &   'nfreq is ',anaddb_dtset%nfreq,', which is lower than 0 .',ch10,&
 885 &   'Action: correct nfreq in your input file.'
 886    MSG_ERROR(message)
 887  end if
 888 
 889  anaddb_dtset%ng2qpt(:)=0
 890  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'ng2qpt',tread,'INT')
 891  if(tread==1) anaddb_dtset%ng2qpt(:)=intarr(1:3)
 892  do ii=1,3
 893    if(anaddb_dtset%ng2qpt(ii)<0)then
 894      write(message, '(a,i0,a,i0,3a,i0,a)' )&
 895 &     'ng2qpt(',ii,') is ',anaddb_dtset%ng2qpt(ii),', which is lower than 0 .',ch10,&
 896 &     'Action: correct ng2qpt(',ii,') in your input file.'
 897      MSG_ERROR(message)
 898    end if
 899  end do
 900 
 901  anaddb_dtset%ngqpt(:)=0
 902  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'ngqpt',tread,'INT')
 903  if(tread==1) anaddb_dtset%ngqpt(1:3)=intarr(1:3)
 904  do ii=1,3
 905    if(anaddb_dtset%ngqpt(ii)<0)then
 906      write(message, '(a,i0,a,i0,3a,i0,a)' )&
 907 &     'ngqpt(',ii,') is ',anaddb_dtset%ngqpt(ii),', which is lower than 0 .',ch10,&
 908 &     'Action: correct ngqpt(',ii,') in your input file.'
 909      MSG_ERROR(message)
 910    end if
 911  end do
 912 
 913  anaddb_dtset%ngrids=4
 914  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ngrids',tread,'INT')
 915  if(tread==1) anaddb_dtset%ngrids=intarr(1)
 916  if(anaddb_dtset%ngrids<0)then
 917    write(message, '(a,i0,3a)' )&
 918 &   'ngrids is ',anaddb_dtset%ngrids,', which is lower than 0 .',ch10,&
 919 &   'Action: correct ngrids in your input file.'
 920    MSG_ERROR(message)
 921  end if
 922 
 923  anaddb_dtset%nlflag=0
 924  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nlflag',tread,'INT')
 925  if(tread==1) anaddb_dtset%nlflag=intarr(1)
 926  if(anaddb_dtset%nlflag<0.or.anaddb_dtset%nlflag>3)then
 927    write(message, '(a,i0,5a)' )&
 928 &   'nlflag is ',anaddb_dtset%nlflag,', but the only allowed values',ch10,&
 929 &   'are 0, 1, 2 or 3.',ch10,'Action: correct nlflag in your input file.'
 930    MSG_ERROR(message)
 931  end if
 932 
 933  anaddb_dtset%nph1l=0
 934  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nph1l',tread,'INT')
 935  if(tread==1) anaddb_dtset%nph1l=intarr(1)
 936  if(anaddb_dtset%nph1l<0)then
 937    write(message, '(a,i0,3a)' )&
 938 &   'nph1l is ',anaddb_dtset%nph1l,', which is lower than 0 .',ch10,&
 939 &   'Action: correct nph1l in your input file.'
 940    MSG_ERROR(message)
 941  end if
 942 
 943  anaddb_dtset%nph2l=0
 944  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nph2l',tread,'INT')
 945  if(tread==1) anaddb_dtset%nph2l=intarr(1)
 946  if(anaddb_dtset%nph2l<0)then
 947    write(message, '(a,i0,3a)' )&
 948 &   'nph2l is ',anaddb_dtset%nph2l,', which is lower than 0 .',ch10,&
 949 &   'Action: correct nph2l in your input file.'
 950    MSG_ERROR(message)
 951  end if
 952 
 953  anaddb_dtset%nqpath=0
 954  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nqpath',tread,'INT')
 955  if(tread==1) anaddb_dtset%nqpath=intarr(1)
 956  if(anaddb_dtset%nqpath<0)then
 957    write(message,'(a,i0,3a)' )&
 958 &   'nqpath is ',anaddb_dtset%nqpath,', but must be positive',ch10,&
 959 &   'Action: correct elphflag in your input file.'
 960    MSG_ERROR(message)
 961  end if
 962 
 963  anaddb_dtset%nqshft=1
 964  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nqshft',tread,'INT')
 965  if(tread==1) anaddb_dtset%nqshft=intarr(1)
 966  if(anaddb_dtset%nqshft<0 .or. anaddb_dtset%nqshft==3 .or.&
 967 & anaddb_dtset%nqshft>=5 )then
 968    write(message, '(a,i0,5a)' )&
 969 &   'nqshft is ',anaddb_dtset%nqshft,', but the only allowed values',ch10,&
 970 &   'are 1, 2 or 4 .',ch10,'Action: correct nqshft in your input file.'
 971    MSG_ERROR(message)
 972  end if
 973 
 974  anaddb_dtset%nsphere=0
 975  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nsphere',tread,'INT')
 976  if(tread==1) anaddb_dtset%nsphere=intarr(1)
 977  if(anaddb_dtset%nsphere < -1)then
 978    write(message, '(a,i0,3a)' )&
 979 &   'nsphere is ',anaddb_dtset%nsphere,', while it must be >= 0 or equal to -1',ch10,&
 980 &   'Action: correct nsphere in your input file.'
 981    MSG_ERROR(message)
 982  end if
 983 
 984  anaddb_dtset%nstrfix=0
 985  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nstrfix',tread,'INT')
 986  if(tread==1) anaddb_dtset%nstrfix=intarr(1)
 987  if(anaddb_dtset%nstrfix > 6)then
 988    write(message, '(a,i0,3a)' )&
 989 &   'nstrfix is ',anaddb_dtset%nstrfix,', which is larger than 6',ch10,&
 990 &   'Action: correct nstrfix in your input file.'
 991    MSG_ERROR(message)
 992  end if
 993 
 994  if(anaddb_dtset%nstrfix < 0)then
 995    write(message, '(a,i0,3a)' )&
 996 &   'nstrfix is ',anaddb_dtset%nstrfix,', which is < 0',ch10,&
 997 &   'Action: correct nstrfix in your input file.'
 998    MSG_ERROR(message)
 999  end if
1000 
1001  anaddb_dtset%ntemper=10
1002  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ntemper',tread,'INT')
1003  if(tread==1) anaddb_dtset%ntemper=intarr(1)
1004  if(anaddb_dtset%ntemper <0)then
1005    write(message, '(a,i0,3a)' )&
1006 &   'ntemper is ',anaddb_dtset%ntemper,', which is lower than 0',ch10,&
1007 &   'Action: correct ntemper in your input file.'
1008    MSG_ERROR(message)
1009  end if
1010 
1011  anaddb_dtset%nwchan=10
1012  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nwchan',tread,'INT')
1013  if(tread==1) anaddb_dtset%nwchan=intarr(1)
1014 !FIXME: check this - it should probably be .ge. 1, not 0
1015  if(anaddb_dtset%nwchan<0)then
1016    write(message, '(a,i0,3a)' )&
1017 &   'nwchan is ',anaddb_dtset%nwchan,', which is lower than 0 .',ch10,&
1018 &   'Action: correct nwchan in your input file.'
1019    MSG_ERROR(message)
1020  end if
1021 
1022 !O
1023  anaddb_dtset%outboltztrap = 0
1024  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'outboltztrap',tread,'INT')
1025  if(tread==1) anaddb_dtset%outboltztrap=intarr(1)
1026  if(anaddb_dtset%outboltztrap<0.or.anaddb_dtset%outboltztrap>1)then
1027    write(message,'(a,i0,5a)' )&
1028 &   'outboltztrap is ',anaddb_dtset%outboltztrap,', but the only allowed values',ch10,&
1029 &   'are 0 or 1.',ch10,'Action: correct outboltztrap in your input file.'
1030    MSG_ERROR(message)
1031  end if
1032 
1033 
1034 !P
1035 
1036  anaddb_dtset%piezoflag=0
1037  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'piezoflag',tread,'INT')
1038  if(tread==1) anaddb_dtset%piezoflag=intarr(1)
1039  if(anaddb_dtset%piezoflag<0.or.anaddb_dtset%piezoflag>7)then
1040    write(message,'(3a,i0,5a)' )&
1041 &   ' piezoflag is ',anaddb_dtset%piezoflag,', but the only allowed values',ch10,&
1042 &   'are 0, 1,2,3,4,5,6,7  .',ch10,'Action: correct piezoflag in your input file.'
1043    MSG_ERROR(message)
1044  end if
1045 
1046  anaddb_dtset%polflag=0
1047  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'polflag',tread,'INT')
1048  if(tread==1) anaddb_dtset%polflag=intarr(1)
1049  if(anaddb_dtset%polflag<0.or.anaddb_dtset%polflag>1)then
1050    write(message, '(a,i0,5a)' )&
1051 &   'polflag is ',anaddb_dtset%polflag,', but the only allowed values',ch10,&
1052 &   'are 0 or 1.',ch10,'Action: correct polflag in your input file.'
1053    MSG_ERROR(message)
1054  end if
1055 
1056  anaddb_dtset%prtbltztrp=0
1057  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtbltztrp',tread,'INT')
1058  if(tread==1) anaddb_dtset%prtbltztrp=intarr(1)
1059  if(anaddb_dtset%prtbltztrp<0.or.anaddb_dtset%prtbltztrp>1)then
1060    write(message, '(a,i0,5a)' )&
1061 &   'prtbltztrp is ',anaddb_dtset%prtbltztrp,', but the only allowed values',ch10,&
1062 &   'are 0 or 1.',ch10,'Action: correct prtbltztrp in your input file.'
1063    MSG_ERROR(message)
1064  end if
1065 
1066 !Default is no output for PHDOS
1067  anaddb_dtset%prtdos=0
1068  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtdos',tread,'INT')
1069  if(tread==1) anaddb_dtset%prtdos = intarr(1)
1070  if(anaddb_dtset%prtdos < 0 .or. anaddb_dtset%prtdos > 2) then
1071    write(message, '(a,i0,5a)' )&
1072 &   'prtdos is ',anaddb_dtset%prtdos,', but the only allowed values',ch10,&
1073 &   'are 0 (no output) or 1 (gaussians) or 2 (tetrahedra) ',ch10,&
1074 &   'Action: correct prtdos in your input file.'
1075    MSG_ERROR(message)
1076  end if
1077 
1078 !Default is no output for the Fermi Surface
1079  anaddb_dtset%prtfsurf = 0
1080  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtfsurf',tread,'INT')
1081  if(tread==1) anaddb_dtset%prtfsurf = intarr(1)
1082  if(anaddb_dtset%prtfsurf < 0 .or. anaddb_dtset%prtfsurf > 2) then
1083    write(message, '(a,i0,5a)' )&
1084 &   'prtfsurf is ',anaddb_dtset%prtfsurf,'. The only allowed values',ch10,&
1085 &   'are 0 (no output) or 1 (Xcrysden bxsf format)',ch10,  &
1086 &   'Action: correct prtfsurf in your input file.'
1087    MSG_ERROR(message)
1088  end if
1089 
1090 !Default is no output of the real space IFC to file
1091  anaddb_dtset%prt_ifc = 0
1092  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prt_ifc',tread,'INT')
1093  if(tread==1) anaddb_dtset%prt_ifc = intarr(1)
1094  if(anaddb_dtset%prt_ifc < 0 .or. anaddb_dtset%prt_ifc > 1) then
1095    write(message, '(a,i0,5a)' )&
1096 &   'prtf_ifc is ',anaddb_dtset%prt_ifc,'. The only allowed values',ch10,&
1097 &   'are 0 (no output) or 1 (AI2PS format)',ch10,  &
1098 &   'Action: correct prt_ifc in your input file.'
1099    MSG_ERROR(message)
1100  end if
1101 ! check that ifcout is set
1102  if (anaddb_dtset%prt_ifc /= 0 .and. anaddb_dtset%ifcout == 0) then
1103    anaddb_dtset%ifcout = -1 ! this forces output of all IFC
1104  end if
1105 
1106 !Default is no output of the DDB to file
1107  anaddb_dtset%prtddb = 0
1108  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtddb',tread,'INT')
1109  if(tread==1) anaddb_dtset%prtddb = intarr(1)
1110  if(anaddb_dtset%prtddb < 0 .or. anaddb_dtset%prtddb > 1) then
1111    write(message, '(a,i0,5a)' )&
1112 &   'prtf_ddb is ',anaddb_dtset%prtddb,'. The only allowed values',ch10,&
1113 &   'are 0 (no output) or 1 (print DDB and DDB.nc files)',ch10,  &
1114 &   'Action: correct prtddb in your input file.'
1115    MSG_ERROR(message)
1116  end if
1117 ! check that ifcflag is set
1118  if (anaddb_dtset%prtddb /= 0 .and. anaddb_dtset%ifcflag == 0) then
1119    anaddb_dtset%ifcflag = 1 ! this forces the use of IFC
1120  end if
1121 
1122  anaddb_dtset%prtmbm=0
1123  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtmbm',tread,'INT')
1124  if(tread==1) anaddb_dtset%prtmbm=intarr(1)
1125 !FIXME: should check whether value of prtmbm is valid
1126 
1127 !Default is no output of the nesting factor
1128  anaddb_dtset%prtnest = 0
1129  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtnest',tread,'INT')
1130  if(tread==1) anaddb_dtset%prtnest = intarr(1)
1131  if(anaddb_dtset%prtnest < 0 .or. anaddb_dtset%prtnest > 2) then
1132    write(message, '(a,i0,5a)' )&
1133 &   'prtnest is ',anaddb_dtset%prtnest,' The only allowed values',ch10,&
1134 &   'are 0 (no nesting), 1 (XY format) or 2 (XY + Xcrysden format)',ch10,&
1135 &   'Action: correct prtnest in your input file.'
1136    MSG_ERROR(message)
1137  end if
1138 
1139  anaddb_dtset%prtphbands=1
1140  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtphbands',tread,'INT')
1141  if (tread==1) anaddb_dtset%prtphbands=intarr(1)
1142  if (all(anaddb_dtset%prtphbands /= [0,1,2])) then
1143    write(message, '(a,i0,a)' )&
1144 &   'prtphbands is ',anaddb_dtset%prtphbands,', but the only allowed values are [0, 1, 2].'
1145    MSG_ERROR(message)
1146  end if
1147 
1148  anaddb_dtset%prtsrlr=0
1149  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtsrlr',tread,'INT')
1150  if(tread==1) anaddb_dtset%prtsrlr=intarr(1)
1151  if(anaddb_dtset%prtsrlr<0.or.anaddb_dtset%prtsrlr>1)then
1152    write(message, '(a,i0,5a)' )&
1153 &   'prtsrlr is ',anaddb_dtset%prtsrlr,', but the only allowed values',ch10,&
1154 &   'are 0 or 1.',ch10,'Action: correct prtsrlr in your input file.'
1155    MSG_ERROR(message)
1156  end if
1157 
1158  anaddb_dtset%prtvol = 0
1159  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtvol',tread,'INT')
1160  if(tread==1) anaddb_dtset%prtvol = intarr(1)
1161 
1162 !Q
1163 
1164  anaddb_dtset%q2shft(:)=zero
1165  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'q2shft',tread,'DPR')
1166  if(tread==1) anaddb_dtset%q2shft(:)=dprarr(1:3)
1167 !FIXME: need a test on valid entries for q2shft
1168 
1169  anaddb_dtset%qgrid_type=1 ! default is uniform nqpt(:) grid
1170  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'qgrid_type',tread,'INT')
1171  if(tread==1) anaddb_dtset%qgrid_type = intarr(1)
1172  if(anaddb_dtset%qgrid_type < 1 .or. anaddb_dtset%qgrid_type > 2) then
1173    write(message, '(a,i0,5a)' )&
1174 &   'qgrid_type is ',anaddb_dtset%qgrid_type,' The only allowed values',ch10,&
1175 &   'are 1 (uniform grid from nqpt) or 2 (listed in ep_nqpt, ep_qptlist)',ch10,&
1176 &   'Action: correct qgrid_type in your input file.'
1177    MSG_ERROR(message)
1178  end if
1179 
1180  anaddb_dtset%qrefine=1 ! default is no refinement
1181  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'qrefine',tread,'INT')
1182  if(tread==1) anaddb_dtset%qrefine = intarr(1:3)
1183  do ii=1,3
1184    if(anaddb_dtset%qrefine(ii) < 1) then
1185      write(message, '(a,3i0,a,a,a,a,a)' )&
1186 &     'qrefine is',anaddb_dtset%qrefine,' The only allowed values',ch10,&
1187 &     'are integers >= 1 giving the refinement of the ngqpt grid',ch10,&
1188 &     'Action: correct qrefine in your input file.'
1189      MSG_ERROR(message)
1190    end if
1191  end do
1192 
1193 !R
1194 
1195  anaddb_dtset%ramansr=0
1196  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ramansr',tread,'INT')
1197  if(tread==1) anaddb_dtset%ramansr=intarr(1)
1198  if(anaddb_dtset%ramansr<0.or.anaddb_dtset%ramansr>2)then
1199    write(message, '(a,i0,5a)' )&
1200 &   'ramansr is ',anaddb_dtset%ramansr,', but the only allowed values',ch10,&
1201 &   'are 0, 1 or 2.',ch10,'Action: correct ramansr in your input file.'
1202    MSG_ERROR(message)
1203  end if
1204 
1205  anaddb_dtset%relaxat=0
1206  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'relaxat',tread,'INT')
1207  if(tread==1) anaddb_dtset%relaxat=intarr(1)
1208  if(anaddb_dtset%relaxat < 0.or.anaddb_dtset%relaxat > 1)then
1209    write(message, '(a,i0,5a)' )&
1210 &   'relaxat is ',anaddb_dtset%relaxat,', but the only allowed values',ch10,&
1211 &   'are 0 or 1.',ch10,'Action: correct relaxat in your input file.'
1212    MSG_ERROR(message)
1213  end if
1214 
1215  anaddb_dtset%relaxstr=0
1216  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'relaxstr',tread,'INT')
1217  if(tread==1) anaddb_dtset%relaxstr=intarr(1)
1218  if(anaddb_dtset%relaxstr<0.or.anaddb_dtset%relaxstr>1)then
1219    write(message, '(a,i0,5a)' )&
1220 &   'relaxstr is ',anaddb_dtset%relaxstr,'but the only allowed values',ch10,&
1221 &   'are 0 or 1.',ch10,'Action: correct relaxstr in your input file.'
1222    MSG_ERROR(message)
1223  end if
1224 
1225  anaddb_dtset%rfmeth=1
1226  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfmeth',tread,'INT')
1227  if(tread==1) anaddb_dtset%rfmeth=intarr(1)
1228  if(anaddb_dtset%rfmeth<1.or.anaddb_dtset%rfmeth>2)then
1229    write(message, '(a,i0,5a)' )&
1230 &   'rfmeth is ',anaddb_dtset%rfmeth,', but the only allowed values',ch10,&
1231 &   'are 1 or 2.',ch10,'Action: correct rfmeth in your input file.'
1232    MSG_ERROR(message)
1233  end if
1234 
1235  anaddb_dtset%rifcsph=zero
1236  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rifcsph',tread,'DPR')
1237  if(tread==1) anaddb_dtset%rifcsph=dprarr(1)
1238 ! if(anaddb_dtset%rifcsph<-tol12)then
1239 !   write(message, '(a,f10.3,3a)' )&
1240 !&   'rifcsph is ',anaddb_dtset%rifcsph,', which is lower than zero.',ch10,&
1241 !&   'Action: correct rifcsph in your input file.'
1242 !   MSG_ERROR(message)
1243 ! end if
1244 
1245 !S
1246 
1247  anaddb_dtset%selectz=0
1248  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'selectz',tread,'INT')
1249  if(tread==1) anaddb_dtset%selectz=intarr(1)
1250  if(anaddb_dtset%selectz<0.or.anaddb_dtset%selectz>2)then
1251    write(message, '(a,i0,5a)' )&
1252 &   'selectz is ',anaddb_dtset%selectz,', but the only allowed values',ch10,&
1253 &   'are 0, 1 or 2 .',ch10,'Action: correct selectz in your input file.'
1254    MSG_ERROR(message)
1255  end if
1256 
1257  anaddb_dtset%symdynmat=1
1258  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'symdynmat',tread,'INT')
1259  if(tread==1) anaddb_dtset%symdynmat=intarr(1)
1260  if(anaddb_dtset%symdynmat/=0.and.anaddb_dtset%symdynmat/=1)then
1261    write(message, '(a,i0,5a)' )&
1262 &   'symdynmat is ',anaddb_dtset%symdynmat,'. The only allowed values',ch10,&
1263 &   'are 0, or 1.',ch10,'Action: correct symdynmat in your input file.'
1264    MSG_ERROR(message)
1265  end if
1266 
1267 !T
1268 
1269  anaddb_dtset%targetpol(:) = 0._dp
1270  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'targetpol',tread,'DPR')
1271  if(tread==1) anaddb_dtset%targetpol(1:3) = dprarr(1:3)
1272 
1273 !Default is use gaussian integration
1274  anaddb_dtset%telphint = 1
1275  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'telphint',tread,'INT')
1276  if(tread==1) anaddb_dtset%telphint = intarr(1)
1277  if(anaddb_dtset%telphint < 0 .or. anaddb_dtset%telphint > 3) then
1278    write(message, '(a,i0,6a)' )&
1279 &   'telphint is ',anaddb_dtset%telphint,'. The only allowed values',ch10,&
1280 &   'are 0 (tetrahedron) or 1 (gaussian) or ','2 (set of bands occupied ep_b_min,ep_b_max) or 3 (Fermi Dirac).',ch10,&
1281 &   'Action: correct telphint in your input file.'
1282    MSG_ERROR(message)
1283  end if
1284 
1285  anaddb_dtset%temperinc=100.0_dp
1286  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'temperinc',tread,'DPR')
1287  if(tread==1) anaddb_dtset%temperinc=dprarr(1)
1288  if(anaddb_dtset%temperinc < zero)then
1289    write(message, '(a,f10.3,3a)' )&
1290 &   'temperinc is ',anaddb_dtset%temperinc,', which is lower than 0 .',ch10,&
1291 &   'Action: correct temperinc in your input file.'
1292    MSG_ERROR(message)
1293  end if
1294 
1295  anaddb_dtset%tempermin=100.0_dp
1296  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tempermin',tread,'DPR')
1297  if(tread==1) anaddb_dtset%tempermin=dprarr(1)
1298  if(anaddb_dtset%tempermin<-tol12)then
1299    write(message, '(a,f10.3,3a)' )&
1300 &   'tempermin is ',anaddb_dtset%tempermin,', which is lower than 0 .',ch10,&
1301 &   'Action: correct tempermin in your input file.'
1302    MSG_ERROR(message)
1303  end if
1304 
1305  anaddb_dtset%thermal_supercell(:,:)=0
1306  marr = 9
1307  ABI_DEALLOCATE(intarr)
1308  ABI_DEALLOCATE(dprarr)
1309  ABI_ALLOCATE(intarr,(marr))
1310  ABI_ALLOCATE(dprarr,(marr))
1311  call intagm(dprarr,intarr,jdtset,marr,9,string(1:lenstr),'thermal_supercell',tread,'INT')
1312  if(tread==1) anaddb_dtset%thermal_supercell(1:3,1:3)=reshape(intarr(1:9),(/3,3/))
1313  call mati3det(anaddb_dtset%thermal_supercell, idet)
1314  if(sum(abs(anaddb_dtset%thermal_supercell))>0 .and. idet == 0) then
1315    write(message, '(a,9I6,5a)' )&
1316 &   'thermal_supercell is ',anaddb_dtset%thermal_supercell,', but the matrix must be non singular',ch10,&
1317 &   'with a non zero determinant.',ch10,'Action: correct thermal_supercell in your input file.'
1318    MSG_ERROR(message)
1319  end if
1320 
1321  anaddb_dtset%thmflag=0
1322  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'thmflag',tread,'INT')
1323  if(tread==1) anaddb_dtset%thmflag=intarr(1)
1324  if(anaddb_dtset%thmflag<0.or.anaddb_dtset%thmflag>8)then
1325    write(message, '(a,i0,5a)' )&
1326 &   'thmflag is ',anaddb_dtset%thmflag,', but the only allowed values',ch10,&
1327 &   'are between 0 to 8 (included).',ch10,'Action: correct thmflag in your input file.'
1328    MSG_ERROR(message)
1329  end if
1330 
1331  anaddb_dtset%thmtol=0.25_dp
1332  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'thmtol',tread,'DPR')
1333  if(tread==1) anaddb_dtset%thmtol=dprarr(1)
1334  if(anaddb_dtset%thmtol<zero)then
1335    write(message, '(a,es14.4,3a)' )&
1336 &   'thmtol is ',anaddb_dtset%thmtol,', which is lower than 0 .',ch10,&
1337 &   'Action: correct thmtol in your input file.'
1338    MSG_ERROR(message)
1339  end if
1340 
1341  anaddb_dtset%ep_prt_yambo = 0
1342  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ep_prt_yambo',tread,'INT')
1343  if(tread==1) anaddb_dtset%ep_prt_yambo = intarr(1)
1344  if(anaddb_dtset%ep_prt_yambo< 0 .or. anaddb_dtset%ep_prt_yambo> 1) then
1345    write(message, '(a,i0,5a)' )&
1346 &   'ep_prt_yambo is ',anaddb_dtset%ep_prt_yambo,', but the only allowed values',ch10,&
1347 &   'are 0 or 1.',ch10,'Action: correct ep_prt_yambo in your input file.'
1348    MSG_ERROR(message)
1349  end if
1350 
1351 !default means _do_ symmetrize the ep coupling matrices over qpoints
1352  anaddb_dtset%symgkq = 1
1353  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'symgkq',tread,'INT')
1354  if(tread==1) anaddb_dtset%symgkq = intarr(1)
1355  if(anaddb_dtset%symgkq< 0 .or. anaddb_dtset%symgkq> 1) then
1356    write(message, '(a,i0,5a)' )&
1357 &   'symgkq is ',anaddb_dtset%symgkq,', but the only allowed values',ch10,&
1358 &   'are 0 or 1.',ch10,'Action: correct symgkq in your input file.'
1359    MSG_ERROR(message)
1360  else if (anaddb_dtset%symgkq == 0) then
1361    MSG_WARNING('You have turned off el-ph matrix symmetrization over q. Use at own risk')
1362  end if
1363 
1364 !U
1365 
1366  anaddb_dtset%use_k_fine = 0
1367  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'use_k_fine',tread,'INT')
1368  if(tread==1) anaddb_dtset%use_k_fine= intarr(1)
1369  if(anaddb_dtset%use_k_fine /= 1 .and. anaddb_dtset%use_k_fine /= 0) then
1370    write(message, '(a,i0,5a)' )&
1371 &   'use_k_fine is ',anaddb_dtset%use_k_fine,', but the only allowed values',ch10,&
1372 &   'are 1 or 0.',ch10,'Action: correct use_k_fine in your input file.'
1373    MSG_ERROR(message)
1374  end if
1375 
1376  if(anaddb_dtset%use_k_fine == 1) then
1377    if (sum(anaddb_dtset%kptrlatt) == 0 .or. sum(anaddb_dtset%kptrlatt_fine) == 0 ) then
1378      MSG_ERROR('If a finer k-grid is used, you must specify both kptrlatt and kptrlatt_fine')
1379    end if
1380  end if
1381 
1382 
1383 !V
1384  anaddb_dtset%vs_qrad_tolkms(:) = zero
1385  call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'vs_qrad_tolkms',tread,'DPR')
1386  if (tread==1) then
1387     anaddb_dtset%vs_qrad_tolkms(:) = dprarr(1:2)
1388     ABI_CHECK(dprarr(1) >= zero, "vs_qrad must be >= 0")
1389     ABI_CHECK(dprarr(2) > zero, "vs_tolkms must be > zero")
1390  end if
1391 !W
1392 
1393 !X
1394 
1395 !Y
1396 
1397 !Z
1398 
1399 !=====================================================================
1400 !end non-dependent variables
1401 !=====================================================================
1402 
1403 !=======================================================================
1404 !Read in dependent variables (dependent on dimensions above)
1405 !=======================================================================
1406 
1407 !A
1408 
1409  ABI_ALLOCATE(anaddb_dtset%atifc,(natom))
1410  anaddb_dtset%atifc(:) = 0
1411  if(anaddb_dtset%natifc>=1)then
1412    ! default to 1 for first natifc atoms
1413    anaddb_dtset%atifc(1:anaddb_dtset%natifc)=1
1414 
1415    if(anaddb_dtset%natifc>marr)then
1416      marr=anaddb_dtset%natifc
1417      ABI_DEALLOCATE(intarr)
1418      ABI_DEALLOCATE(dprarr)
1419      ABI_ALLOCATE(intarr,(marr))
1420      ABI_ALLOCATE(dprarr,(marr))
1421    end if
1422    call intagm(dprarr,intarr,jdtset,marr,anaddb_dtset%natifc,string(1:lenstr),'atifc',tread,'INT')
1423    if(tread==1) anaddb_dtset%atifc(1:anaddb_dtset%natifc)= intarr(1:anaddb_dtset%natifc)
1424 !  check of whether values of atifc are valid is done in chkin9
1425  end if
1426 
1427 !B
1428 
1429 !C
1430 
1431 !D
1432 
1433 !E
1434 
1435 !F
1436 
1437 !G
1438 
1439 !H
1440 
1441 !I
1442 
1443  ABI_ALLOCATE(anaddb_dtset%iatfix,(natom))
1444  anaddb_dtset%iatfix(:) = 0
1445  if ((anaddb_dtset%relaxat == 1).and.(anaddb_dtset%natfix > 0)) then
1446    if(natom > marr)then
1447      marr = natom
1448      ABI_DEALLOCATE(intarr)
1449      ABI_DEALLOCATE(dprarr)
1450      ABI_ALLOCATE(intarr,(marr))
1451      ABI_ALLOCATE(dprarr,(marr))
1452    end if
1453    call intagm(dprarr,intarr,jdtset,marr,anaddb_dtset%natfix,string(1:lenstr),'iatfix',tread,'INT')
1454    if(tread==1) anaddb_dtset%iatfix(1:anaddb_dtset%natfix) = intarr(1:anaddb_dtset%natfix)
1455  end if
1456 !FIXME: need a test on values of iatfix: are they just 1 or 0?
1457 
1458  if ((anaddb_dtset%relaxstr == 1).and.(anaddb_dtset%nstrfix > 0)) then
1459    anaddb_dtset%istrfix(:) = 0
1460    call intagm(dprarr,intarr,jdtset,marr,anaddb_dtset%nstrfix,string(1:lenstr),'istrfix',tread,'INT')
1461    if(tread==1) anaddb_dtset%istrfix(1:anaddb_dtset%nstrfix) = intarr(1:anaddb_dtset%nstrfix)
1462  end if
1463 !FIXME: need a test on values of istrfix
1464 
1465  if (anaddb_dtset%natprj_bs > 0) then
1466    ABI_ALLOCATE(anaddb_dtset%iatprj_bs,(anaddb_dtset%natprj_bs))
1467    if(anaddb_dtset%natprj_bs>marr)then
1468      marr=anaddb_dtset%natprj_bs
1469      ABI_DEALLOCATE(intarr)
1470      ABI_DEALLOCATE(dprarr)
1471      ABI_ALLOCATE(intarr,(marr))
1472      ABI_ALLOCATE(dprarr,(marr))
1473    end if
1474    anaddb_dtset%iatprj_bs(:)=0
1475    call intagm(dprarr,intarr,jdtset,marr,anaddb_dtset%natprj_bs,string(1:lenstr),'iatprj_bs',tread,'INT')
1476    if(tread==1) then
1477      anaddb_dtset%iatprj_bs(1:anaddb_dtset%natprj_bs)=intarr(1:anaddb_dtset%natprj_bs)
1478    else
1479      write(message,'(3a)')&
1480 &     'natprj_bs is non zero but iatprj_bs is absent ',ch10,&
1481 &     'Action: specify iatprj_bs in your input file.'
1482      MSG_ERROR(message)
1483    end if
1484  end if
1485 
1486 !J
1487 
1488 !K
1489 
1490 !L
1491 
1492 !M
1493 
1494 !N
1495 
1496 !O
1497 
1498 !P
1499 
1500 !Q
1501 
1502  if (anaddb_dtset%nqshft/=0)then
1503    if(3*anaddb_dtset%nqshft>marr)then
1504      marr=3*anaddb_dtset%nqshft
1505      ABI_DEALLOCATE(intarr)
1506      ABI_DEALLOCATE(dprarr)
1507      ABI_ALLOCATE(intarr,(marr))
1508      ABI_ALLOCATE(dprarr,(marr))
1509    end if
1510    anaddb_dtset%q1shft(:,:)=zero
1511    call intagm(dprarr,intarr,jdtset,marr,3*anaddb_dtset%nqshft, string(1:lenstr),'q1shft',tread,'DPR')
1512    if(tread==1) anaddb_dtset%q1shft(1:3,1:anaddb_dtset%nqshft)=&
1513 &   reshape(dprarr(1:3*anaddb_dtset%nqshft),(/3,anaddb_dtset%nqshft/))
1514  end if
1515 
1516  ABI_ALLOCATE(anaddb_dtset%qph1l,(3,anaddb_dtset%nph1l))
1517  ABI_ALLOCATE(anaddb_dtset%qnrml1,(anaddb_dtset%nph1l))
1518  if (anaddb_dtset%nph1l/=0)then
1519    if(4*anaddb_dtset%nph1l>marr)then
1520      marr=4*anaddb_dtset%nph1l
1521      ABI_DEALLOCATE(intarr)
1522      ABI_DEALLOCATE(dprarr)
1523      ABI_ALLOCATE(intarr,(marr))
1524      ABI_ALLOCATE(dprarr,(marr))
1525    end if
1526    anaddb_dtset%qph1l(:,:)=zero
1527    anaddb_dtset%qnrml1(:)=zero
1528    call intagm(dprarr,intarr,jdtset,marr,4*anaddb_dtset%nph1l,string(1:lenstr),'qph1l',tread,'DPR')
1529    if(tread==1)then
1530      do iph1=1,anaddb_dtset%nph1l
1531        do ii=1,3
1532          anaddb_dtset%qph1l(ii,iph1)=dprarr(ii+(iph1-1)*4)
1533        end do
1534        anaddb_dtset%qnrml1(iph1)=dprarr(4+(iph1-1)*4)
1535        if(abs(anaddb_dtset%qnrml1(iph1))<DDB_QTOL)then
1536          write(message, '(5a)' )&
1537 &         'The first list of wavevectors ','should not have non-analytical data.',ch10,&
1538 &         'Action: correct the first list',' of wavevectors in the input file.'
1539          MSG_ERROR(message)
1540        end if
1541      end do
1542    end if
1543  end if
1544 
1545  ABI_ALLOCATE(anaddb_dtset%qph2l,(3,anaddb_dtset%nph2l))
1546  ABI_ALLOCATE(anaddb_dtset%qnrml2,(anaddb_dtset%nph2l))
1547  if (anaddb_dtset%nph2l/=0)then
1548    if(4*anaddb_dtset%nph2l>marr)then
1549      marr=4*anaddb_dtset%nph2l
1550      ABI_DEALLOCATE(intarr)
1551      ABI_DEALLOCATE(dprarr)
1552      ABI_ALLOCATE(intarr,(marr))
1553      ABI_ALLOCATE(dprarr,(marr))
1554    end if
1555    anaddb_dtset%qph2l(:,:)=zero
1556    anaddb_dtset%qnrml2(:)=zero
1557    call intagm(dprarr,intarr,jdtset,marr,4*anaddb_dtset%nph2l,string(1:lenstr),'qph2l',tread,'DPR')
1558    if(tread==1)then
1559      do iph2=1,anaddb_dtset%nph2l
1560        do ii=1,3
1561          anaddb_dtset%qph2l(ii,iph2)=dprarr(ii+(iph2-1)*4)
1562        end do
1563        anaddb_dtset%qnrml2(iph2)=dprarr(4+(iph2-1)*4)
1564        if(abs(anaddb_dtset%qnrml2(iph2))>DDB_QTOL)then
1565          write(message, '(5a)' )&
1566 &         'The second list of wavevectors',' should have only non-analytical data.',ch10,&
1567 &         'Action: correct the second list','of wavevectors in the input file.'
1568          MSG_ERROR(message)
1569        end if
1570      end do
1571    end if
1572  end if
1573 
1574  if (anaddb_dtset%nqpath > 0) then
1575    ABI_ALLOCATE(anaddb_dtset%qpath,(3,anaddb_dtset%nqpath))
1576    if(3*anaddb_dtset%nqpath>marr)then
1577      marr=3*anaddb_dtset%nqpath
1578      ABI_DEALLOCATE(intarr)
1579      ABI_DEALLOCATE(dprarr)
1580      ABI_ALLOCATE(intarr,(marr))
1581      ABI_ALLOCATE(dprarr,(marr))
1582    end if
1583    anaddb_dtset%qpath(:,:)=zero
1584    call intagm(dprarr,intarr,jdtset,marr,3*anaddb_dtset%nqpath,string(1:lenstr),'qpath',tread,'DPR')
1585    if(tread==1) then
1586      anaddb_dtset%qpath(1:3,1:anaddb_dtset%nqpath)=&
1587 &     reshape(dprarr(1:3*anaddb_dtset%nqpath),(/3,anaddb_dtset%nqpath/))
1588    else
1589      write(message,'(3a)')&
1590 &     'nqpath is non zero but qpath is absent ',ch10,&
1591 &     'Action: specify qpath in your input file.'
1592      MSG_ERROR(message)
1593    end if
1594  end if
1595 
1596 !R
1597 
1598 !S
1599 
1600 !T
1601 
1602 !U
1603 
1604 !V
1605 
1606 !W
1607 
1608 !X
1609 
1610 !Y
1611 
1612 !Z
1613 
1614 !=======================================================================
1615 !Finished reading in variables - deallocate
1616 !=======================================================================
1617 
1618  ABI_DEALLOCATE(dprarr)
1619  ABI_DEALLOCATE(intarr)
1620 
1621 !=======================================================================
1622 !Check consistency of input variables:
1623 !=======================================================================
1624 
1625  if (anaddb_dtset%frmin > anaddb_dtset%frmax) then
1626    write(message,'(3a)' )&
1627 &   'frmax should be higher than frmin',ch10,&
1628 &   'Action: change frmax and/or frmin  in your input file.'
1629    MSG_ERROR(message)
1630  end if
1631 
1632  if (anaddb_dtset%nqpath==0 .and. anaddb_dtset%elphflag==1) then
1633    write(message,'(4a)' )&
1634 &   'elphflag is 1 but no nqpath has been specified','for phonon linewidths',ch10,&
1635 &   'Action: specify nqpath and qpath(3,nqpath) in your input file.'
1636    MSG_ERROR(message)
1637  end if
1638 
1639  if(anaddb_dtset%telphint /= 2 .and. (anaddb_dtset%ep_b_min /= 0 .or.&
1640 & anaddb_dtset%ep_b_max /= 0)) then
1641    write(message, '(a,i0,3a)' )&
1642 &   'telphint is ',anaddb_dtset%telphint,', but ep_b_min or ep_b_max',ch10,&
1643 &   'are set /= 1. They will not be used'
1644    call wrtout(std_out,message,'COLL')
1645    MSG_WARNING(message)
1646 
1647  else if(anaddb_dtset%telphint == 2 .and. &
1648 &   (anaddb_dtset%ep_b_min == 0 .or. anaddb_dtset%ep_b_max == 0)) then
1649    write(message, '(a,i0,4a)' )&
1650 &   'telphint is ',anaddb_dtset%telphint,', but ep_b_min or ep_b_max',ch10,&
1651 &   'are not both set. ',ch10,&
1652 &   'Action: set ep_b_min and ep_b_max in your input file.',ch10
1653    MSG_ERROR(message)
1654  end if
1655 
1656  if(anaddb_dtset%thmflag < 3) then
1657    if ((anaddb_dtset%telphint == 0 .or. anaddb_dtset%prtnest == 1 .or. &
1658 &   anaddb_dtset%prtnest == 2 .or. anaddb_dtset%prtfsurf== 1) .and.&
1659 &   sum(anaddb_dtset%kptrlatt) == 0 ) then
1660      write (message, '(3a)') &
1661 &     'if tetrahedron integration is used, ',&
1662 &     'or the output of the nesting function/Fermi surface is required, ',&
1663 &     'you must specify the kptrlatt'
1664      MSG_ERROR(message)
1665    end if
1666  end if
1667 
1668  if(anaddb_dtset%prtdos/=0 .and. anaddb_dtset%ifcflag/=1) then
1669    write(message, '(3a)' )&
1670 &   'ifcflag must be 1 when the calculation of the phonon DOS is required ',ch10,&
1671 &   'Action: correct ifcflag in your input file.'
1672    MSG_ERROR(message)
1673  end if
1674 
1675  if(anaddb_dtset%prtsrlr/=0 .and. anaddb_dtset%ifcflag/=1) then
1676    write(message, '(3a)' )&
1677 &   'ifcflag must be 1 for the SR/LR decomposition of the phonon frequencies',ch10,&
1678 &   'Action: correct ifcflag in your input file.'
1679    MSG_ERROR(message)
1680  end if
1681 
1682  if (anaddb_dtset%gruns_nddbs /=0 .and. anaddb_dtset%ifcflag /=1) then
1683    MSG_ERROR("ifcflag must be 1 for Grunesein calculation")
1684  end if
1685 
1686  if (anaddb_dtset%vs_qrad_tolkms(1) /= zero .and. anaddb_dtset%ifcflag /=1) then
1687    MSG_ERROR("ifcflag must be 1 to calculate speed of sound")
1688  end if
1689 
1690  if(anaddb_dtset%prtdos/=0 .and. sum(abs(anaddb_dtset%ng2qpt(:))) < 3 ) then
1691    write(message, '(3a)' )&
1692 &   'ng2qpt must be specified when the calculation of the phonon DOS is required ',ch10,&
1693 &   'Action: correct ng2qpt in your input file.'
1694    MSG_ERROR(message)
1695  end if
1696 
1697  if (anaddb_dtset%ifltransport /= 0 .and. anaddb_dtset%ep_keepbands /= 1) then
1698    write(message, '(3a)' )&
1699 &   'Band dependency of electron phonon matrix elements must be kept for transport ',ch10,&
1700 &   'Action: set ep_keepbands to 1 in your input file.'
1701    MSG_ERROR(message)
1702  end if
1703 
1704  if (anaddb_dtset%ifltransport > 1 .and. sum(abs(anaddb_dtset%kptrlatt)) == 0) then
1705    write(message, '(3a)' )&
1706 &   'For inelastic transport or electron lifetime calculations you must specify kprtlatt ',ch10,&
1707 &   'Action: copy kptrlatt from your abinit GS file to your anaddb input file.'
1708    MSG_ERROR(message)
1709  end if
1710 
1711 !FIXME: add check that if freeze_displ /= 0 then you need to be doing ifc and phonon interpolation
1712 
1713  if (anaddb_dtset%ifcflag > 0 .and. sum(abs(anaddb_dtset%ngqpt)) == 0) then
1714    write(message, '(3a)' )&
1715 &   'if you want interatomic force constant output, anaddb needs ngqpt input variable ',ch10,&
1716 &   'Action: set ngqpt in your input file.'
1717    MSG_ERROR(message)
1718  end if
1719 
1720 !check that q-grid refinement is a divisor of ngqpt in each direction
1721  if(any(anaddb_dtset%qrefine(1:3) > 1) .and. &
1722 &   any(abs(dmod(dble(anaddb_dtset%ngqpt(1:3))/dble(anaddb_dtset%qrefine(1:3)),one)) > tol10) ) then
1723    write(message, '(a,3i10,a,a,a,3i8,a,a)' )&
1724 &   'qrefine is',anaddb_dtset%qrefine(1:3),' The only allowed values',ch10,&
1725 &   'are integers which are divisors of the ngqpt grid', anaddb_dtset%ngqpt(1:3),ch10,&
1726 &   'Action: correct qrefine in your input file.'
1727    MSG_ERROR(message)
1728  end if
1729 
1730 !check that fermie and nelect are not both specified
1731  if(abs(anaddb_dtset%elph_fermie) > tol10 .and. &
1732 & abs(anaddb_dtset%ep_extrael) > tol10) then
1733    write(message, '(a,E10.2,a,E10.2,a,a,a)' )&
1734 &   'elph_fermie (',anaddb_dtset%elph_fermie,') and ep_extrael (',anaddb_dtset%ep_extrael, '), may not both be non 0',ch10,&
1735 &   'Action: remove one of the two in your input file.'
1736    MSG_ERROR(message)
1737  end if
1738 
1739  ! Check for possible typos.
1740  call anaddb_chkvars(string)
1741 
1742 end subroutine invars9

m_anaddb_dataset/outvars_anaddb [ Functions ]

[ Top ] [ m_anaddb_dataset ] [ Functions ]

NAME

 outvars_anaddb

FUNCTION

 Open input file for the anaddb code, then
 echoes the input information.

INPUTS

 anaddb_dtset= (derived datatype) contains all the input variables
 nunit=unit number for input or output

OUTPUT

  (only writing)

NOTES

 Should be executed by one processor only.

PARENTS

      anaddb

CHILDREN

      chkvars_in_string,inupper

SOURCE

1774 subroutine outvars_anaddb (anaddb_dtset,nunit)
1775 
1776 
1777 !This section has been created automatically by the script Abilint (TD).
1778 !Do not modify the following lines by hand.
1779 #undef ABI_FUNC
1780 #define ABI_FUNC 'outvars_anaddb'
1781 !End of the abilint section
1782 
1783  implicit none
1784 
1785 !Arguments -------------------------------
1786 !scalars
1787  integer,intent(in) :: nunit
1788  type(anaddb_dataset_type),intent(in) :: anaddb_dtset
1789 
1790 !Local variables -------------------------
1791 !scalars
1792  integer :: ii,iph1,iph2,iqpt,iqshft
1793 
1794 !*********************************************************************
1795 
1796 !Write the heading
1797  write(nunit,'(a,80a,a)') ch10,('=',ii=1,80),ch10
1798  write(nunit, '(2a)' )' -outvars_anaddb: echo values of input variables ----------------------',ch10
1799 
1800 !The flags
1801  if(anaddb_dtset%dieflag/=0 .or. anaddb_dtset%ifcflag/=0 .or. &
1802 & anaddb_dtset%nlflag/=0 .or. anaddb_dtset%thmflag/=0 .or. &
1803 & anaddb_dtset%elaflag/=0 .or. anaddb_dtset%elphflag/=0 .or. &
1804 & anaddb_dtset%polflag/=0 .or. anaddb_dtset%instrflag/=0 .or. &
1805 & anaddb_dtset%piezoflag/=0                                   &
1806 & )then
1807    write(nunit,'(a)')' Flags :'
1808    if(anaddb_dtset%dieflag/=0)write(nunit,'(3x,a9,3i10)')'  dieflag',anaddb_dtset%dieflag
1809    if(anaddb_dtset%ifcflag/=0)write(nunit,'(3x,a9,3i10)')'  ifcflag',anaddb_dtset%ifcflag
1810    if(anaddb_dtset%nlflag/=0)write(nunit,'(3x,a9,3i10)')'   nlflag',anaddb_dtset%nlflag
1811    if(anaddb_dtset%thmflag/=0)write(nunit,'(3x,a9,3i10)')'  thmflag',anaddb_dtset%thmflag
1812    if(anaddb_dtset%elaflag/=0)write(nunit,'(3x,a9,3i10)')'  elaflag',anaddb_dtset%elaflag
1813    if(anaddb_dtset%elphflag/=0)write(nunit,'(3x,a9,3i10)')' elphflag',anaddb_dtset%elphflag
1814    if(anaddb_dtset%polflag/=0)write(nunit,'(3x,a9,3i10)')'  polflag',anaddb_dtset%polflag
1815    if(anaddb_dtset%instrflag/=0)write(nunit,'(3x,a9,3i10)')'instrflag',anaddb_dtset%instrflag
1816    if(anaddb_dtset%piezoflag/=0)write(nunit,'(3x,a9,3i10)')'piezoflag',anaddb_dtset%piezoflag
1817  end if
1818 
1819 !Write the general information
1820  if( anaddb_dtset%rfmeth/=1 .or. &
1821 & anaddb_dtset%enunit/=0 .or. &
1822 & anaddb_dtset%eivec/=0 .or. &
1823 & anaddb_dtset%asr/=0 .or. &
1824 & anaddb_dtset%chneut/=0 .or. &
1825 & anaddb_dtset%selectz/=0  &
1826 & .or. anaddb_dtset%symdynmat/=1       &
1827 & )then
1828    write(nunit,'(a)')' Miscellaneous information :'
1829    if(anaddb_dtset%rfmeth/=1)write(nunit,'(3x,a9,3i10)')'   rfmeth',anaddb_dtset%rfmeth
1830    if(anaddb_dtset%enunit/=0)write(nunit,'(3x,a9,3i10)')'   enunit',anaddb_dtset%enunit
1831    if(anaddb_dtset%eivec/=0) write(nunit,'(3x,a9,3i10)')'    eivec',anaddb_dtset%eivec
1832    if(anaddb_dtset%asr/=0)   write(nunit,'(3x,a9,3i10)')'      asr',anaddb_dtset%asr
1833    if(anaddb_dtset%chneut/=0)write(nunit,'(3x,a9,3i10)')'   chneut',anaddb_dtset%chneut
1834    if(anaddb_dtset%selectz/=0)write(nunit,'(3x,a9,3i10)')'  selectz',anaddb_dtset%selectz
1835    if(anaddb_dtset%symdynmat/=1)write(nunit,'(3x,a9,3i10)')'symdynmat',anaddb_dtset%symdynmat
1836  end if
1837  if(anaddb_dtset%prtvol/=0) write(nunit,'(3x,a9,i10)')'   prtvol',anaddb_dtset%prtvol
1838 
1839 !Frequency information
1840  if(anaddb_dtset%dieflag==1)then
1841    write(nunit,'(a)')' Frequency information :'
1842    write(nunit,'(3x,a9,3i10)')'    nfreq',anaddb_dtset%nfreq
1843    write(nunit,'(3x,a9,7x,3es16.8)')'    frmin',anaddb_dtset%frmin
1844    write(nunit,'(3x,a9,7x,3es16.8)')'    frmax',anaddb_dtset%frmax
1845  end if
1846 
1847 !For interatomic force constant information
1848  if(anaddb_dtset%ifcflag/=0)then
1849    write(nunit,'(a)')' Interatomic Force Constants Inputs :'
1850    write(nunit,'(3x,a9,3i10)')'   dipdip',anaddb_dtset%dipdip
1851    if(anaddb_dtset%nsphere/=0)write(nunit,'(3x,a9,3i10)')'  nsphere',anaddb_dtset%nsphere
1852    if(abs(anaddb_dtset%rifcsph)>tol10)write(nunit,'(3x,a9,E16.6)')'  nsphere',anaddb_dtset%rifcsph
1853    write(nunit,'(3x,a9,3i10)')'   ifcana',anaddb_dtset%ifcana
1854    write(nunit,'(3x,a9,3i10)')'   ifcout',anaddb_dtset%ifcout
1855    if(anaddb_dtset%natifc>=1)then
1856      write(nunit,'(3x,a9,3i10)')'   natifc',anaddb_dtset%natifc
1857      write(nunit,'(3x,a9,8i10)')'    atifc',(anaddb_dtset%atifc(ii),ii=1,anaddb_dtset%natifc)
1858    end if
1859    write(nunit,'(a)')' Description of grid 1 :'
1860    write(nunit,'(3x,a9,3i10)')'     brav',anaddb_dtset%brav
1861    write(nunit,'(3x,a9,3i10)')'    ngqpt',anaddb_dtset%ngqpt(1:3)
1862    write(nunit,'(3x,a9,3i10)')'   nqshft',anaddb_dtset%nqshft
1863    if (anaddb_dtset%nqshft/=0)then
1864      write(nunit,'(3x,a9)')'   q1shft'
1865      do iqshft=1,anaddb_dtset%nqshft
1866        write(nunit,'(19x,4es16.8)') (anaddb_dtset%q1shft(ii,iqshft),ii=1,3)
1867      end do
1868    end if
1869    if (any(anaddb_dtset%qrefine(:) > 1)) then
1870      write(nunit,'(3x,a9,3i10)')'  qrefine', anaddb_dtset%qrefine
1871    end if
1872    ! Speed of sound
1873    if (anaddb_dtset%vs_qrad_tolkms(1) > zero) then
1874       write(nunit,'(a,2es16.8)')"vs_qrad_tolkms", (anaddb_dtset%vs_qrad_tolkms(:))
1875    end if
1876  end if
1877 
1878 !Phonon density of states with gaussian method
1879  if(anaddb_dtset%prtdos/=0)then
1880    write(nunit,'(a)')' Phonon DOS information :'
1881    write(nunit,'(3x,a9,es16.8)')'dosdeltae',anaddb_dtset%dosdeltae
1882    write(nunit,'(3x,a9,es16.8)')' dossmear',anaddb_dtset%dossmear
1883  end if
1884 
1885 !Thermal information
1886  if(anaddb_dtset%thmflag/=0)then
1887    write(nunit,'(a)')' Thermal information :'
1888    write(nunit,'(3x,a9,3i10)')'    nchan',anaddb_dtset%nchan
1889    write(nunit,'(3x,a9,3i10)')'   nwchan',anaddb_dtset%nwchan
1890    write(nunit,'(3x,a9,7x,3es16.8)')'   dostol',anaddb_dtset%dostol
1891    write(nunit,'(3x,a9,7x,3es16.8)')'   thmtol',anaddb_dtset%thmtol
1892    write(nunit,'(3x,a9,3i10)')'  ntemper',anaddb_dtset%ntemper
1893    write(nunit,'(3x,a9,7x,3es16.8)')'temperinc',anaddb_dtset%temperinc
1894    write(nunit,'(3x,a9,7x,3es16.8)')'tempermin',anaddb_dtset%tempermin
1895  endif
1896 
1897 !Grid 2 description
1898  if(anaddb_dtset%thmflag/=0 .or. anaddb_dtset%prtdos/=0)then
1899    write(nunit,'(a)')' Description of grid 2 (Fourier interp. or BZ sampling):'
1900    write(nunit,'(3x,a9,3i10)')'   ng2qpt',anaddb_dtset%ng2qpt(1:3)
1901    write(nunit,'(3x,a9,3i10)')'   ngrids',anaddb_dtset%ngrids
1902    write(nunit,'(3x,a9,7x,3es16.8)')'   q2shft',anaddb_dtset%q2shft(1:3)
1903  end if
1904 
1905 !Non-linear response information
1906  if (anaddb_dtset%nlflag /= 0) then
1907    write(nunit,'(a)')' Non-linear response information :'
1908    write(nunit,'(3x,a9,i10)') '   alphon',anaddb_dtset%alphon
1909    write(nunit,'(3x,a9,3i10)')'   prtmbm',anaddb_dtset%prtmbm
1910    write(nunit,'(3x,a9,3i10)')'  ramansr',anaddb_dtset%ramansr
1911  end if
1912 
1913 !Structural relaxation at fixed polarization
1914  if (anaddb_dtset%polflag /= 0) then
1915    write(nunit,'(a)')' Relaxation at fixed polarization :'
1916    if (anaddb_dtset%relaxat == 1) then
1917      write(nunit,'(3x,a9,i10)') '  relaxat',anaddb_dtset%relaxat
1918    end if
1919    if (anaddb_dtset%relaxstr == 1) then
1920      write(nunit,'(a12,i10)') ' relaxstr',anaddb_dtset%relaxstr
1921    end if
1922  end if
1923 
1924 !Elphon information
1925  if (anaddb_dtset%elphflag /= 0) then
1926    write(nunit,'(a)')' Elphon calculation will be carried out'
1927    write(nunit,'(a12,E16.6)') 'elphsmear', anaddb_dtset%elphsmear
1928    write(nunit,'(a12,E16.6)') 'a2fsmear', anaddb_dtset%a2fsmear
1929    write(nunit,'(a12,E16.6)') 'mustar', anaddb_dtset%mustar
1930    write(nunit,'(a12,i10)') 'nqpath', anaddb_dtset%nqpath
1931    write(nunit,'(a12)') 'qpath'
1932    do iqpt=1,anaddb_dtset%nqpath
1933      write(nunit,'(12x,3(E16.6,1x))') anaddb_dtset%qpath(:,iqpt)
1934    end do
1935    write(nunit,'(a12,i10)') 'telphint', anaddb_dtset%telphint
1936    if (anaddb_dtset%telphint == 0) then
1937      write(nunit,'(a)') ' Tetrahedron integration for elphon'
1938    else if (anaddb_dtset%telphint == 1) then
1939      write(nunit,'(a)') ' Smeared weight integration for elphon'
1940    else if (anaddb_dtset%telphint == 2) then
1941      write(nunit,'(a)') ' Band filtered integration for elphon'
1942    end if
1943    if (abs(anaddb_dtset%elph_fermie) > tol10) then
1944      write(nunit,'(a12,E16.6)')  'elph_fermie', anaddb_dtset%elph_fermie
1945    end if
1946    if (anaddb_dtset%ep_extrael /= 0) then
1947      if (abs(anaddb_dtset%ep_extrael) > 1.0d2) then
1948         write(nunit,'(a,E20.12)')' Doping set by the user is (negative for el doping) :',&
1949 &       anaddb_dtset%ep_extrael
1950      else
1951        write(nunit,'(a,E16.6)')  'Elphon: extra electrons per unit cell = ', anaddb_dtset%ep_extrael
1952      end if
1953    end if
1954    if (anaddb_dtset%ep_nspline /= 20) then
1955      write(nunit,'(a,I8)')  'Elphon: scale factor for spline interpolation in RTA = ', anaddb_dtset%ep_nspline
1956    end if
1957    if (anaddb_dtset%band_gap < 10.0d0) then
1958      write(nunit,'(a,E16.6)')  'Elphon: set band gap to (in eV) = ', anaddb_dtset%band_gap
1959    end if
1960 
1961    if (sum(abs(anaddb_dtset%kptrlatt)) > 0) then
1962      write(nunit,'(a12,3(3(i3,1x),2x))' ) 'kptrlatt',&
1963 &     reshape( anaddb_dtset%kptrlatt(:,:), (/9/) )
1964    end if
1965 
1966    if (sum(abs(anaddb_dtset%kptrlatt_fine)) > 0) then
1967      write(nunit,'(a12,3(3(i3,1x),2x))' ) 'kptrlatt_fine ',&
1968 &     reshape( anaddb_dtset%kptrlatt_fine(:,:), (/9/) )
1969    end if
1970 
1971    if (anaddb_dtset%ep_keepbands == 1) then
1972      write(nunit, '(a)') ' Will keep band dependency in gkk in memory.'
1973      write(nunit, '(a)') ' WARNING: the memory requirements will be multiplied by nbands**2 !!!'
1974    end if
1975 
1976    if (anaddb_dtset%ep_scalprod == 1) then
1977      write(nunit, '(a)') ' scalar product will be performed when assembling the gamma matrices.'
1978      write(nunit, '(a)') ' WARNING: with this option you can not distinguish which '
1979      write(nunit, '(a)') '    linewidth comes from which phonon mode !!!'
1980    end if
1981 
1982    if (anaddb_dtset%prtbltztrp== 1) then
1983      write(nunit, '(a)') ' Will output input files for BoltzTraP'
1984    end if
1985 
1986    if (anaddb_dtset%prtfsurf == 1) then
1987      write(nunit, '(a)') ' Will output fermi surface in XCrysDen format'
1988    end if
1989 
1990    if (anaddb_dtset%prt_ifc == 1) then
1991      write(nunit, '(a)') ' Will output real space IFC in AI2PS and TDEP format'
1992    end if
1993 
1994    if (anaddb_dtset%prtnest == 1) then
1995      write(nunit, '(a)') ' Will output nesting factor'
1996    end if
1997 
1998    if (anaddb_dtset%ifltransport == 1) then
1999      write(nunit, '(a)') ' Will perform transport calculation in elphon to get'
2000      write(nunit, '(a,a)') ' resistivity and thermal conductivity as a function of T',ch10
2001      write(nunit, '(a,es16.6,a)' ) ' Minimum temperature for transport outputs: ', &
2002 &     anaddb_dtset%tempermin, ' K'
2003      write(nunit, '(a,es16.6,a)' ) ' Maximum temperature for transport outputs: ', &
2004 &     anaddb_dtset%tempermin+anaddb_dtset%temperinc*anaddb_dtset%ntemper, ' K'
2005      write(nunit, '(a,i6)' ) ' Number of temperature points for transport outputs: ', anaddb_dtset%ntemper
2006      write(nunit, '(a)' )
2007    end if
2008 
2009    if (anaddb_dtset%gkqwrite == 1) then
2010      write(nunit,'(a,a)' ) 'Gkk matrix elements on input grid of ',&
2011      'qpoints will be written to disk. File gkqfile must be absent.'
2012    end if
2013    if (anaddb_dtset%gkk_rptwrite == 1) then
2014      write(nunit,'(a,a)' ) 'Gkk matrix elements in real space ',&
2015      'will be written to disk. File gkk_rpt_file must be absent.'
2016    end if
2017    if (anaddb_dtset%gkk2write == 1) then
2018      write(nunit,'(a,a)' ) 'Full grid gkk matrix elements ',&
2019      'will be written to disk. File gkk2file must be absent.'
2020    end if
2021  end if
2022 
2023  if (anaddb_dtset%gruns_nddbs /= 0) then
2024    write(nunit,'(a)' ) "Will compute Gruneisen parameters with finite difference method. DDB files:"
2025    do ii=1,anaddb_dtset%gruns_nddbs
2026      write(nunit, "(2a)")"    ",trim(anaddb_dtset%gruns_ddbs(ii))
2027    end do
2028  end if
2029 
2030 !List of vector 1  (reduced coordinates)
2031  if(anaddb_dtset%nph1l/=0)then
2032    write(nunit,'(a)')' First list of wavevector (reduced coord.) :'
2033    write(nunit,'(3x,a9,3i10)')'    nph1l',anaddb_dtset%nph1l
2034    write(nunit,'(3x,a9)')'    qph1l'
2035    do iph1=1,anaddb_dtset%nph1l
2036      write(nunit,'(19x,3es16.8,2x,es11.3)') &
2037 &     (anaddb_dtset%qph1l(ii,iph1),ii=1,3),anaddb_dtset%qnrml1(iph1)
2038    end do
2039  end if
2040 
2041 !List of vector 2  (cartesian coordinates)
2042  if(anaddb_dtset%nph2l/=0)then
2043    write(nunit,'(a)')' Second list of wavevector (cart. coord.) :'
2044    write(nunit,'(3x,a9,3i10)')'    nph2l',anaddb_dtset%nph2l
2045    write(nunit,'(3x,a9)')'    qph2l'
2046    do iph2=1,anaddb_dtset%nph2l
2047      write(nunit,'(19x,3es16.8,2x,es11.3)') &
2048 &     (anaddb_dtset%qph2l(ii,iph2),ii=1,3),anaddb_dtset%qnrml2(iph2)
2049    end do
2050  end if
2051 
2052 !phonon frozen in supercell
2053  if (abs(anaddb_dtset%freeze_displ) > tol10) then
2054    write(nunit,'(a)') 'Phonon displacements will be output, frozen into supercells'
2055    write(nunit,'(a,E20.10)') ' Chosen amplitude of frozen displacements = ', anaddb_dtset%freeze_displ
2056  end if
2057 
2058 !atom projected bs files
2059  if (abs(anaddb_dtset%natprj_bs) > 0) then
2060    write(nunit,'(a)') 'Phonon band structure files, with atomic projections, will be output '
2061    write(nunit,'(a)') ' Chosen atoms for projection = '
2062    write(nunit,'(10I6)') anaddb_dtset%iatprj_bs
2063  end if
2064 
2065  write(nunit,'(a,80a,a)') ch10,('=',ii=1,80),ch10
2066 
2067 end subroutine outvars_anaddb