!{\src2tex{textfont=tt}}
!!****f* ABINIT/invars2
!!
!! NAME
!! invars2
!!
!! FUNCTION
!! Initialize variables for the ABINIT code, for one particular
!! dataset, characterized by jdtset.
!! Note : some parameters have already been read in invars0 and invars1,
!! and were used to dimension the arrays needed here.
!!
!! COPYRIGHT
!! Copyright (C) 1998-2014 ABINIT group (DCA, XG, GMR)
!! This file is distributed under the terms of the
!! GNU General Public License, see ~abinit/COPYING
!! or http://www.gnu.org/copyleft/gpl.txt .
!! For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
!!
!! INPUTS
!! bravais(11): bravais(1)=iholohedry
!!              bravais(2)=center
!!              bravais(3:11)=coordinates of rprim in the axes
!!               of the conventional bravais lattice (*2 if center/=0)
!! iout=unit number for echoed output
!! jdtset=number of the dataset looked for
!! lenstr=actual length of the string
!! mband=maximum number of bands for any k-point
!! msym=default maximal number of symmetries
!! npsp=number of pseudopotentials
!! string*(*)=character string containing all the input data.
!!  Initialized previously in instrng.
!! usepaw= 0 for non paw calculation; =1 for paw calculation
!! zionpsp(npsp)=valence charge of each type of atom (coming from the psp files)
!!
!! OUTPUT
!!  (see side effects)
!!
!! SIDE EFFECTS
!! Input/Output:
!!  dtset=<type datafiles_type>contains all input variables,
!!   some of which are initialized here, while other were already
!! All rest of arguments given alphabetically from acell (length scales)
!! to wtk (k-point weights), used to control running of the main routine.
!! See abinit_help for definitions of these variables.
!! These values become defined by being read from string,
!! that contains all information from the input file,
!! in a compressed, standardized, format
!! At the input, they already contain a default value.
!!
!! NOTES
!!
!! PARENTS
!!      invars2m,m_ab7_invars_f90
!!
!! CHILDREN
!!      chkneu,inkpts,intagm,invacuum,smpbz,sort_int,timab,wrtout
!!
!! SOURCE

#if defined HAVE_CONFIG_H
#include "config.h"
#endif

#include "abi_common.h"

subroutine invars2(bravais,dtset,iout,jdtset,lenstr,&
& mband,msym,npsp,string,usepaw,zionpsp)

 use defs_basis
 use defs_abitypes
 use m_profiling
 use m_errors
#if defined HAVE_DFT_LIBXC
 use xc_f90_lib_m
#endif

!This section has been created automatically by the script Abilint (TD).
!Do not modify the following lines by hand.
#undef ABI_FUNC
#define ABI_FUNC 'invars2'
 use interfaces_14_hidewrite
 use interfaces_18_timing
 use interfaces_28_numeric_noabirule
 use interfaces_42_parser
 use interfaces_56_recipspace
 use interfaces_57_iovars, except_this_one => invars2
!End of the abilint section

 implicit none

!Arguments ------------------------------------
!scalars
 integer,intent(in) :: iout,jdtset,lenstr,mband,msym,npsp,usepaw
 character(len=*),intent(in) :: string
 type(dataset_type),intent(inout) :: dtset
!arrays
 integer,intent(in) :: bravais(11)
 real(dp),intent(in) :: zionpsp(npsp)

!Local variables-------------------------------
!scalars
 integer :: bantot,berryopt,dmatsize,getocc,iatom,ii,ikpt,index_nat,ionmov,intimage
 integer :: iprcch,ipsp,iscf,itypat,jj,kptopt,lpawu,marr,natom,nband1,nberry,mkpt
 integer :: niatcon,nimage,nkpt,npspalch,nqpt,nspinor,nsppol,nsym,ntypalch,ntypat,ntyppure,nkpt_hf
 integer :: occopt,occopt_tmp,response,tfband,tnband,tread,tread_alt,tread_key
 integer :: itol, itol_gen, ds_input, ifreq
 real(dp) :: areaxy,charge,fband,kptrlen,nelectjell
 real(dp) :: rhoavg,zelect,zval
 real(dp) :: toldfe_, tolrff_, toldff_, tolwfr_, tolvrs_
 character(len=500) :: message
 character(len=fnlen) :: key_value
!arrays
 integer :: vacuum(3)
 integer,allocatable :: iatcon(:),natcon(:)
 integer,allocatable :: intarr(:) !Dummy arguments for subroutine 'intagm' to parse input file
 real(dp) :: tsec(2)
 real(dp),allocatable :: dprarr(:) !Dummy arguments for subroutine 'intagm' to parse input file
 real(dp),allocatable :: kptns_hf(:,:) !Dummy arguments for subroutine 'smpbz' to calculate the value of nkpthf

! *************************************************************************

 call timab(191,1,tsec)

!Compute the maximum size of arrays intarr and dprarr
 natom=dtset%natom
 nimage=dtset%nimage
 nkpt=dtset%nkpt
 npspalch=dtset%npspalch
 nspinor=dtset%nspinor
 nsppol=dtset%nsppol
 ntypat=dtset%ntypat
 ntypalch=dtset%ntypalch
 ntyppure=dtset%ntyppure

 dmatsize=0
 if (dtset%usepawu>0.and.dtset%usedmatpu/=0) then
   do iatom=1,natom
     lpawu=dtset%lpawu(dtset%typat(iatom))
     if (lpawu/=-1) dmatsize=dmatsize+nsppol*nspinor*(2*lpawu+1)**2
   end do
 end if

 marr=max(3*natom,&
& nkpt*nsppol*mband,&
& 2*dtset%nkptgw*nsppol,&
& dmatsize,&
& 3*nkpt,&
& npsp,&
& ntypat,&
& 9*msym,&
& 60,&
& 3*dtset%nconeq*natom,&
& nimage,&
& 3*dtset%nqptdm,&
& 3*dtset%natsph_extra,&
& dtset%natvshift*nsppol*natom)
 ABI_ALLOCATE(intarr,(marr))
 ABI_ALLOCATE(dprarr,(marr))

!----------------------------------------------------------------------------

!****   Read parameters which set remaining array dimensions ****

!Note : some parameters have already been read in invars0 and invars1
!Also, some checking is needed here.

!Read ngfft(1), ngfft(2), and ngfft(3),
!then ngfft(7)=fftalg and ngfft(8)=fftcache.
!Read ngfftdg(1:3)

 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'ngfft',tread,'INT')
 if(tread==1) dtset%ngfft(1:3)=intarr(1:3)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fftalg',tread,'INT')
 if(tread==1) then
   dtset%ngfft(7)=intarr(1)
   if (usepaw==1) dtset%ngfftdg(7)=intarr(1)
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fftcache',tread,'INT')

 if(tread==1) then
   dtset%ngfft(8)=intarr(1)
   if (usepaw==1) dtset%ngfftdg(8)=intarr(1)
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fft_opt_lob',tread,'INT')
 if(tread==1) then
   dtset%fft_opt_lob=intarr(1)
 else
   if(dtset%paral_kgb==1) dtset%fft_opt_lob=2
 end if

 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'ngfftdg',tread,'INT')
 if(tread==1) dtset%ngfftdg(1:3)=intarr(1:3)

 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'ng',tread,'INT')
 if(tread==1) then
   write(message, '(a,a,a)' )&
&   'The use of the "ng" input variable is forbidden since version 1.8.',ch10,&
&   'Action: take "ng" out of your input file.'
   MSG_ERROR(message)
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'mqgrid',tread,'INT')
 if(tread==1) dtset%mqgrid=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'mqgriddg',tread,'INT')
 if(tread==1) dtset%mqgriddg=intarr(1)

!Make special arrangements to check nband: may be a scalar
!(for occopt=0, 1 or 3, 4, 5, 6, 7, 8) or may be an array (for occopt=2)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'occopt',tread,'INT')
 if(tread==1) dtset%occopt=intarr(1)
 occopt=dtset%occopt

!check for variables related to genetic algorithm. ga_n_rules has been already read

 if (dtset%imgmov==4) then
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ga_algor',tread,'INT')
   if(tread==1) dtset%ga_algor=intarr(1)

   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ga_fitness',tread,'INT')
   if(tread==1) dtset%ga_fitness=intarr(1)

   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ga_opt_percent',tread,'ENE')
   if(tread==1) dtset%ga_opt_percent=dprarr(1)

   call intagm(dprarr,intarr,jdtset,marr,dtset%ga_n_rules,string(1:lenstr),'ga_rules',tread,'INT')
!  DEBUG AHR
!  write(std_out,*) ' write ga_n_rules ',dtset%ga_n_rules
!  call flush(std_out)
!  ENDDEBUG
   if(tread==1)then
     dtset%ga_rules(1:dtset%ga_n_rules)=intarr(1:dtset%ga_n_rules)
     do ii=1,dtset%ga_n_rules
       if(dtset%ga_rules(ii)<0)then
         write(message, '(a,i4,a)' )&
&         'All values for Genetic rules must be greater than 0.',ch10,&
&         'Action: check the values of ga_rules.'
         MSG_ERROR(message)
       end if
     end do
   end if
 end if


 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gwcalctyp',tread,'INT')
 if(tread==1) dtset%gwcalctyp=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'chkgwcomp',tread,'INT')
 if(tread==1) dtset%chkgwcomp=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gwcomp',tread,'INT')
 if(tread==1) dtset%gwcomp=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gwencomp',tread,'ENE')
 if(tread==1) dtset%gwencomp=dprarr(1)

 if (usepaw==1) then
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_sigxcore',tread,'INT')
   if(tread==1) dtset%gw_sigxcore=intarr(1)
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gwmem',tread,'INT')
 if(tread==1) dtset%gwmem=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_sctype',tread,'INT')
 if(tread==1) dtset%gw_sctype=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_nstep',tread,'INT')
 if(tread==1) dtset%gw_nstep=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_toldfeig',tread,'ENE')
 if(tread==1) dtset%gw_toldfeig=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rhoqpmix',tread,'DPR')
 if(tread==1) dtset%rhoqpmix=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nfreqim',tread,'INT')
 if(tread==1) dtset%nfreqim=intarr(1)
 if (dtset%cd_customnimfrqs/=0) then
   write(message, '(a,a,a)' )&
&   'cd_customnimfrqs not equal to zero and not equal to nfreqim',ch10,&
&   'setting nfreqim = cd_customnimfrqs'
   MSG_WARNING(message)
   dtset%nfreqim = dtset%cd_customnimfrqs
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'freqim_alpha',tread,'DPR')
 if(tread==1) dtset%freqim_alpha=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'freqremin',tread,'ENE')
 if(tread==1) dtset%freqremin=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'freqremax',tread,'ENE')
 if(tread==1) dtset%freqremax=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nfreqre',tread,'INT')
 if(tread==1) dtset%nfreqre=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nfreqsp',tread,'INT')
 if(tread==1) dtset%nfreqsp=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'freqspmax',tread,'ENE')
 if(tread==1) dtset%freqspmax=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'freqspmin',tread,'ENE')
 if(tread==1) then ! If found, set it
   dtset%freqspmin=dprarr(1)
 else ! Else give it the value -freqspmax
   dtset%freqspmin=-dtset%freqspmax
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_frqim_inzgrid',tread,'INT')
 if(tread==1) dtset%gw_frqim_inzgrid=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_frqre_inzgrid',tread,'INT')
 if(tread==1) dtset%gw_frqre_inzgrid=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_frqre_tangrid',tread,'INT')
 if(tread==1) dtset%gw_frqre_tangrid=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_use_pole_scr',tread,'INT')
 if(tread==1) dtset%gw_use_pole_scr=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_qprange',tread,'INT')
 if(tread==1) dtset%gw_qprange=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_npoles',tread,'INT')
 if(tread==1) dtset%gw_npoles=intarr(1)

 if(dtset%cd_customnimfrqs/=0) then
   call intagm(dprarr,intarr,jdtset,marr,dtset%cd_customnimfrqs,string(1:lenstr),'cd_imfrqs',tread,'ENE')
   if(tread==1) then
     dtset%cd_imfrqs(1:dtset%cd_customnimfrqs)=dprarr(1:dtset%cd_customnimfrqs)
     do ifreq=1,dtset%cd_customnimfrqs-1
       if (dtset%cd_imfrqs(ifreq+1)<dtset%cd_imfrqs(ifreq)) then
         write(message, '(a,a,a)' )&
&         'The frequencies specified in cd_imfrqs must be strictly increasing',ch10,&
&         'Action: Correct this in your input file'
         MSG_ERROR(message)
       end if
     end do
   end if
 end if

 if(dtset%gw_customnfreqsp/=0) then
   call intagm(dprarr,intarr,jdtset,marr,dtset%gw_customnfreqsp,string(1:lenstr),'gw_freqsp',tread,'ENE')
   if(tread==1) then
     dtset%gw_freqsp(1:dtset%gw_customnfreqsp)=dprarr(1:dtset%gw_customnfreqsp)
     do ifreq=1,dtset%gw_customnfreqsp-1
       if (dtset%gw_freqsp(ifreq+1)<dtset%gw_freqsp(ifreq)) then
         write(message, '(a,a,a)' )&
&         'The frequencies specified in gw_freqsp must be strictly increasing',ch10,&
&         'Action: Correct this in your input file'
         MSG_ERROR(message)
       end if
     end do
   end if
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'cd_full_grid',tread,'INT')
 if(tread==1) dtset%cd_full_grid=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'cd_frqim_method',tread,'INT')
 if(tread==1) dtset%cd_frqim_method=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'cd_halfway_freq',tread,'ENE')
 if(tread==1) dtset%cd_halfway_freq=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'cd_max_freq',tread,'ENE')
 if(tread==1) dtset%cd_max_freq=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_reconst_scr',tread,'INT')
 if(tread==1) dtset%gw_reconst_scr=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'cd_subset_freq',tread,'INT')
 if(tread==1) dtset%cd_subset_freq(1:2)=intarr(1:2)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gwrpacorr',tread,'INT')
 if(tread==1) dtset%gwrpacorr=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'supercell',tread,'INT')
 if(tread==1) dtset%supercell(1:3)=intarr(1:3)

!RESPFN integer input variables (needed here to get the value of response
!Presently, rfmeth is not used.
!Warning : rfddk,rfelfd,rfphon,rfstrs and rfuser are also read in invars1
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfasr',tread,'INT')
 if(tread==1) dtset%rfasr=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'rfatpol',tread,'INT')
 if(tread==1) dtset%rfatpol(1:2)=intarr(1:2)

 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'rfdir',tread,'INT')
 if(tread==1) dtset%rfdir(1:3)=intarr(1:3)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfddk',tread,'INT')
 if(tread==1) dtset%rfddk=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfelfd',tread,'INT')
 if(tread==1) dtset%rfelfd=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfmeth',tread,'INT')
 if(tread==1) dtset%rfmeth=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfphon',tread,'INT')
 if(tread==1) dtset%rfphon=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfstrs',tread,'INT')
 if(tread==1) dtset%rfstrs=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfuser',tread,'INT')
 if(tread==1) dtset%rfuser=intarr(1)

!set value of response to 1 and also set rfdir to 1 1 1 if we are doing
!response calculation but rfdir was left at default 0 0 0 value
 response=0
 if(dtset%rfddk/=0 .or. dtset%rfelfd/=0 .or. &
& dtset%rfphon/=0 .or. dtset%rfstrs/=0 .or. dtset%rfuser/=0 ) then
   response=1
   if( (dtset%rfdir(1) == 0) .AND. &
&   (dtset%rfdir(2) == 0) .AND. &
&   (dtset%rfdir(3) == 0) ) dtset%rfdir(1:3) = 1 
 end if

!NONLINEAR integer input variables (same definition as for rfarr)
!Presently, rf?asr, rf?meth,rf?strs and rf?thrd are not used
 call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'rf1atpol',tread,'INT')
 if(tread==1) dtset%rf1atpol(1:2)=intarr(1:2)

 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'rf1dir',tread,'INT')
 if(tread==1) dtset%rf1dir(1:3)=intarr(1:3)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rf1elfd',tread,'INT')
 if(tread==1) dtset%rf1elfd=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rf1phon',tread,'INT')
 if(tread==1) dtset%rf1phon=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'rf2atpol',tread,'INT')
 if(tread==1) dtset%rf2atpol(1:2)=intarr(1:2)

 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'rf2dir',tread,'INT')
 if(tread==1) dtset%rf2dir(1:3)=intarr(1:3)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rf2elfd',tread,'INT')
 if(tread==1) dtset%rf2elfd=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rf2phon',tread,'INT')
 if(tread==1) dtset%rf2phon=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'rf3atpol',tread,'INT')
 if(tread==1) dtset%rf3atpol(1:2)=intarr(1:2)

 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'rf3dir',tread,'INT')
 if(tread==1) dtset%rf3dir(1:3)=intarr(1:3)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rf3elfd',tread,'INT')
 if(tread==1) dtset%rf3elfd=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rf3phon',tread,'INT')
 if(tread==1) dtset%rf3phon=intarr(1)

 response=0
 if(dtset%rfddk/=0 .or. dtset%rfelfd/=0 .or. &
& dtset%rfphon/=0 .or. &
& dtset%rfstrs/=0 .or. dtset%rfuser/=0 .or. &
& dtset%rf1elfd/=0 .or. dtset%rf1phon/=0 .or. &
& dtset%rf2elfd/=0 .or. dtset%rf2phon/=0 .or. &
& dtset%rf3elfd/=0 .or. dtset%rf3phon/=0 ) response=1

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prepanl',tread,'INT')
 if(tread==1) dtset%prepanl=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prepgkk',tread,'INT')
 if(tread==1) dtset%prepgkk=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'use_nonscf_gkk',tread,'INT')
 if(tread==1) dtset%use_nonscf_gkk=intarr(1)

!real(dp) input variables
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'boxcutmin',tread,'DPR')
 if(tread==1) dtset%boxcutmin=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'charge',tread,'DPR')
 if(tread==1) dtset%charge=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dosdeltae',tread,'ENE')
 if(tread==1) dtset%dosdeltae=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dtion',tread,'DPR')
 if(tread==1) dtset%dtion=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ecut',tread,'ENE')
 if(tread==1) dtset%ecut=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'elph2_imagden',tread,'ENE')
 if(tread==1) dtset%elph2_imagden=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'esmear',tread,'ENE')
 if(tread==1) dtset%esmear=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fermie_nest',tread,'ENE')
 if(tread==1) dtset%fermie_nest=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'sciss',tread,'ENE')
 if(tread==1) dtset%sciss=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tsmear',tread,'ENE')
 if(tread==1) dtset%tsmear=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vis',tread,'DPR')
 if(tread==1) dtset%vis=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ecutsm',tread,'ENE')
 if(tread==1) dtset%ecutsm=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'exchmix',tread,'DPR')
 if(tread==1) dtset%exchmix=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dilatmx',tread,'DPR')
 if(tread==1) dtset%dilatmx=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fxcartfactor',tread,'DPR')
 if(tread==1) dtset%fxcartfactor=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'strfact',tread,'DPR')
 if(tread==1) dtset%strfact=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'freqsusin',tread,'DPR')
 if(tread==1) dtset%freqsusin=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'freqsuslo',tread,'DPR')
 if(tread==1) dtset%freqsuslo=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'optfreqsus',tread,'INT')
 if(tread==1) dtset%optfreqsus=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'effmass',tread,'DPR')
 if(tread==1) dtset%effmass=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'mdtemp',tread,'DPR')
 if(tread==1) dtset%mdtemp(1:2)=dprarr(1:2)

!Recursion input variables

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tfkinfunc',tread,'INT')
 if(tread==1) dtset%tfkinfunc=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'recgratio',tread,'INT')
 if(tread==1) dtset%recgratio=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'recefermi',tread,'ENE')
 if(tread==1) dtset%recefermi=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'recnpath',tread,'INT')
 if(tread==1) dtset%recnpath=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'recnrec',tread,'INT')
 if(tread==1) dtset%recnrec=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'recrcut',tread,'LEN')
 if(tread==1) dtset%recrcut=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'recptrott',tread,'INT')
 if(tread==1) dtset%recptrott=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rectesteg',tread,'INT')
 if(tread==1) dtset%rectesteg=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rectolden',tread,'DPR')
 if(tread==1) dtset%rectolden=dprarr(1)

!Constant NPT Molecular Dynamics Input variables
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'noseinert',tread,'DPR')
 if(tread==1) dtset%noseinert=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bmass',tread,'DPR')
 if(tread==1) dtset%bmass=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,dtset%nnos,string(1:lenstr),'qmass',tread,'DPR')
 if(tread==1) dtset%qmass(:)=dprarr(1:dtset%nnos)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tphysel',tread,'ENE')
 if(tread==1) dtset%tphysel=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,6,string(1:lenstr),'strtarget',tread,'DPR')
 if(tread==1) dtset%strtarget(1:6)=dprarr(1:6)

 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'vcutgeo',tread,'DPR')
 if(tread==1) dtset%vcutgeo(1:3)=dprarr(1:3)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'friction',tread,'DPR')
 if(tread==1) dtset%friction=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'mdwall',tread,'LEN')
 if(tread==1) dtset%mdwall=dprarr(1)


!Path-Integral Molecular Dynamics Input variables
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'adpimd',tread,'INT')
 if(tread==1) dtset%adpimd=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'adpimd_gamma',tread,'DPR')
 if(tread==1) dtset%adpimd_gamma=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pitransform',tread,'INT')
 if(tread==1) dtset%pitransform=intarr(1)

!Default for pimass is amu
 dtset%pimass(1:ntypat)=dtset%amu_orig(1:ntypat,1)   ! NOTE : initialisation with the first image only. TO BE MODIFIED ....
 call intagm(dprarr,intarr,jdtset,marr,ntypat,string(1:lenstr),'pimass',tread,'DPR')
 if(tread==1) dtset%pimass(1:ntypat)=dprarr(1:ntypat)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'spinmagntarget',tread,'DPR')
 if(tread==1) dtset%spinmagntarget=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'eshift',tread,'ENE')
 if(tread==1) dtset%eshift=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'boxcenter',tread,'DPR')
 if(tread==1) dtset%boxcenter(1:3)=dprarr(1:3)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ecuteps',tread,'ENE')
 if(tread==1) dtset%ecuteps=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ecutsigx',tread,'ENE')
 if(tread==1) then
   dtset%ecutsigx=dprarr(1)
!  else
!  if(dtset%optdriver==RUNL_SIGMA) dtset%ecutsigx=dtset%ecut
 end if
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ecutwfn',tread,'ENE')
 if(tread==1) then
   dtset%ecutwfn=dprarr(1)
 else
   if(dtset%optdriver==RUNL_SCREENING .or. dtset%optdriver==RUNL_SIGMA) dtset%ecutwfn=dtset%ecut
 end if
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'omegasimax',tread,'ENE')
 if(tread==1) dtset%omegasimax=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'omegasrdmax',tread,'ENE')
 if(tread==1) dtset%omegasrdmax=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'soenergy',tread,'ENE')
 if(tread==1) dtset%soenergy=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'spbroad',tread,'ENE')
 if(tread==1) dtset%spbroad=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'stmbias',tread,'ENE')
 if(tread==1) dtset%stmbias=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'awtr',tread,'INT')
 if(tread==1) dtset%awtr=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'inclvkb',tread,'INT')
 if(tread==1) dtset%inclvkb=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nomegasf',tread,'INT')
 if(tread==1) dtset%nomegasf=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'spmeth',tread,'INT')
 if(tread==1) dtset%spmeth=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'symchi',tread,'INT')
 if(tread==1) dtset%symchi=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getscr',tread,'INT')
 if(tread==1) dtset%getscr=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gwgamma',tread,'INT')
 if(tread==1) dtset%gwgamma=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdsuscep',tread,'INT')
 if(tread==1) dtset%irdsuscep=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdscr',tread,'INT')
 if(tread==1) dtset%irdscr=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nomegasi',tread,'INT')
 if(tread==1) dtset%nomegasi=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ppmodel',tread,'INT')
 if(tread==1) dtset%ppmodel=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'symsigma',tread,'INT')
 if(tread==1) dtset%symsigma=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fftgw',tread,'INT')
 if(tread==1) dtset%fftgw=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getsuscep',tread,'INT')
 if(tread==1) dtset%getsuscep=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getqps',tread,'INT')
 if(tread==1) dtset%getqps=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getkss',tread,'INT')
 if(tread==1) dtset%getkss=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gwpara',tread,'INT')
 if(tread==1) dtset%gwpara=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdqps',tread,'INT')
 if(tread==1) dtset%irdqps=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdkss',tread,'INT')
 if(tread==1) dtset%irdkss=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ppmfrq',tread,'ENE')
 if(tread==1) dtset%ppmfrq=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rcut',tread,'LEN')
 if(tread==1) dtset%rcut=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'zcut',tread,'ENE')
 if(tread==1) dtset%zcut=dprarr(1)

!q-points for long wave-length limit.
 if (dtset%gw_nqlwl>0) then
   call intagm(dprarr,intarr,jdtset,marr,3*dtset%gw_nqlwl,string(1:lenstr),'gw_qlwl',tread,'DPR')
   if(tread==1) dtset%gw_qlwl(1:3,1:dtset%gw_nqlwl) = reshape(dprarr(1:3*dtset%gw_nqlwl),(/3,dtset%gw_nqlwl/))
 end if

!@bethe_salpeter

!integer scalars
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bs_nstates',tread,'INT')
 if(tread==1) dtset%bs_nstates=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bs_algorithm',tread,'INT')
 if(tread==1) dtset%bs_algorithm=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bs_haydock_niter',tread,'INT')
 if(tread==1) dtset%bs_haydock_niter=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bs_hayd_term',tread,'INT')
 if(tread==1) dtset%bs_hayd_term=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bs_exchange_term',tread,'INT')
 if(tread==1) dtset%bs_exchange_term=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bs_calctype',tread,'INT')
 if(tread==1) dtset%bs_calctype=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bs_coulomb_term',tread,'INT')
 if(tread==1) dtset%bs_coulomb_term=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bs_coupling',tread,'INT')
 if(tread==1) dtset%bs_coupling=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdbseig',tread,'INT')
 if(tread==1) dtset%irdbseig=intarr(1)

!call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdbsr',tread,'INT')
!if(tread==1) dtset%irdbsr=intarr(1)

!call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdbsc',tread,'INT')
!if(tread==1) dtset%irdbsc=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getbseig',tread,'INT')
 if(tread==1) dtset%getbseig=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getbsreso',tread,'INT')
 if(tread==1) dtset%getbsreso=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getbscoup',tread,'INT')
 if(tread==1) dtset%getbscoup=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gethaydock',tread,'INT')
 if(tread==1) dtset%gethaydock=intarr(1)

!integer arrays
 call intagm(dprarr,intarr,jdtset,marr,dtset%nsppol,string(1:lenstr),'bs_loband',tread,'INT')
 if(tread==1) dtset%bs_loband=intarr(1:dtset%nsppol)

!real arrays
 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'bs_eh_cutoff',tread,'ENE')
 if(tread==1) dtset%bs_eh_cutoff(1:3)=dprarr(1:3)

 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'bs_freq_mesh',tread,'ENE')
 if(tread==1) dtset%bs_freq_mesh(1:3)=dprarr(1:3)

 call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'bs_haydock_tol',tread,'DPR')
 if(tread==1) dtset%bs_haydock_tol=dprarr(1:2)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_eet',tread,'INT')
 if(tread==1) dtset%gw_eet=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_eet_nband',tread,'INT')
 if(tread==1) dtset%gw_eet_nband=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_eet_inclvkb',tread,'INT')
 if(tread==1) dtset%gw_eet_inclvkb=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_eet_scale',tread,'ENE')
 if(tread==1) dtset%gw_eet_scale=dprarr(1)

 ntypalch=dtset%ntypalch
 npspalch=dtset%npspalch

!Compute ziontypat. When the pseudo-atom is pure, simple copy
 if(ntyppure>0)then
   do itypat=1,ntyppure
     dtset%ziontypat(itypat)=zionpsp(itypat)
   end do
 end if
!When the pseudo-atom is alchemical, must make mixing
 if(ntypalch>0)then
   do itypat=ntyppure+1,ntypat
     dtset%ziontypat(itypat)=zero
     do ipsp=ntyppure+1,npsp
       dtset%ziontypat(itypat)=dtset%ziontypat(itypat) &
&       +dtset%mixalch_orig(ipsp-ntyppure,itypat-ntyppure,1)*zionpsp(ipsp)
     end do
   end do
 end if

 charge=dtset%charge

 if (occopt==0 .or. occopt==1 .or. (occopt>=3 .and. occopt<=8) ) then

   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nband',tnband,'INT')
   if(tnband==1) then
     nband1=intarr(1)
   else
!    Default value in the metallic case, or in the insulating case
     fband=0.5_dp
     if(occopt==1)fband=0.125_dp
     if((occopt/=1).and.(dtset%accuracy==5.or.dtset%accuracy==6)) fband =0.75_dp
     if (dtset%usewvl == 1) fband = zero
     call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fband',tfband,'DPR')
     if(tfband==1)then
       fband=dprarr(1)
       write(message, '(a,es16.8,a)' )' invars2: read the value of fband=',fband,' from input file.'
     else
       write(message, '(a,es16.8)' )' invars2: take the default value of fband=',fband
     end if
     dtset%fband=fband
     call wrtout(std_out,message,'COLL')
!    First compute the total valence charge
     zval=0.0_dp
     do iatom=1,natom
       zval=zval+dtset%ziontypat(dtset%typat(iatom))
     end do
     zelect=zval-charge
!    Then select the minimum number of bands, and add the required number
!    Note that this number might be smaller than the one computed
!    by a slightly different formula in invars1
     nband1=dtset%nspinor * ((ceiling(zelect-1.0d-10)+1)/2 + ceiling( fband*natom - 1.0d-10 ))
   end if

!  Set nband to same input number for each k point and spin
!  where nband1 is the eventual input, computed value, or default
   do ikpt=1,nkpt*nsppol
     dtset%nband(ikpt)=nband1
   end do

 else if (occopt==2) then
!  Give nband explicitly for each k point and spin

   call intagm(dprarr,intarr,jdtset,nkpt*nsppol,nkpt*nsppol,string(1:lenstr),'nband',tnband,'INT')
   if(tnband==1) dtset%nband(1:nkpt*nsppol)=intarr(1:nkpt*nsppol)

 else
   write(message, '(a,i0,a,a,a)' )&
&   'occopt=',occopt,' not allowed.',ch10,&
&   'Action: correct your input file.'
   MSG_ERROR(message)
 end if

!----------------------------------------------------------------------------

!****   Read other parameters  ****

!All checking should be done in chkinp.f

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'builtintest',tread,'INT')
 if(tread==1) dtset%builtintest=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'chksymbreak',tread,'INT')
 if(tread==1) dtset%chksymbreak=intarr(1)

!Get array
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getocc',tread,'INT')
 if(tread==1) dtset%getocc=intarr(1)

 getocc=dtset%getocc
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getwfk',tread,'INT')
 if(tread==1) dtset%getwfk=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getwfkfine',tread,'INT')
 if(tread==1) dtset%getwfkfine=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getxcart',tread,'INT')
 if(tread==1) dtset%getxcart=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getxred',tread,'INT')
 if(tread==1) dtset%getxred=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getden',tread,'INT')
 if(tread==1) dtset%getden=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getpawden',tread,'INT')
 if(tread==1) dtset%getpawden=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getcell',tread,'INT')
 if(tread==1) dtset%getcell=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getgam_eig2nkq',tread,'INT')
 if(tread==1) dtset%getgam_eig2nkq=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getwfq',tread,'INT')
 if(tread==1) dtset%getwfq=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'get1wf',tread,'INT')
 if(tread==1) dtset%get1wf=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getddk',tread,'INT')
 if(tread==1) dtset%getddk=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getvel',tread,'INT')
 if(tread==1) dtset%getvel=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'get1den',tread,'INT')
 if(tread==1) dtset%get1den=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getwf',tread,'INT')
 if(tread==1) then
   write(message, '(a,a,a)' )&
&   'The use of the "getwf" variable is forbidden since version 2.0.',ch10,&
&   'Action: replace "getwf" by "getwfk" in your input file.'
   MSG_ERROR(message)
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'accesswff',tread,'INT')
 if(tread==1) dtset%accesswff=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'enunit',tread,'INT')
 if(tread==1) dtset%enunit=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'exchn2n3d',tread,'INT')
 if(tread==1) dtset%exchn2n3d=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'iboxcut',tread,'INT')
 if(tread==1) dtset%iboxcut=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'icutcoul',tread,'INT')
 if(tread==1) dtset%icutcoul=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'mdf_epsinf',tread,'DPR')
 if(tread==1) dtset%mdf_epsinf=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'mem_test',tread,'INT')
 if(tread==1) dtset%mem_test=intarr(1)

 if (dtset%imgmov==0.or.dtset%imgmov==2.or.dtset%imgmov==5) then
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'mep_mxstep',tread,'LEN')
   if(tread==1) then
     dtset%mep_mxstep=dprarr(1)
   else if(dtset%imgmov==5) then
     dtset%mep_mxstep=0.4_dp
   end if
 end if
 if (dtset%imgmov==2.or.dtset%imgmov==5) then
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'mep_solver',tread,'INT')
   if(tread==1) then
     dtset%mep_solver=intarr(1)
   else if(dtset%imgmov==2) then
     dtset%mep_solver=0
   else if(dtset%imgmov==5) then
     dtset%mep_solver=0
   end if
   if (dtset%imgmov==2) then
     call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'string_algo',tread,'INT')
     if(tread==1) dtset%string_algo=intarr(1)
   end if
   if (dtset%imgmov==5) then
     call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'neb_algo',tread,'INT')
     if(tread==1) dtset%neb_algo=intarr(1)
     call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'neb_spring',tread,'DPR')
     if(tread==1) then
       dtset%neb_spring(1:2)=dprarr(1:2)
     else if (dtset%neb_algo==2) then
       dtset%neb_spring(1:2)=(/0.02_dp,0.15_dp/)
     end if
     call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'cineb_start',tread,'INT')
     if(tread==1) dtset%cineb_start=intarr(1)
   end if
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'istatimg',tread,'INT')
 if(tread==1) dtset%istatimg=intarr(1)

!variables for random positions in unit cell
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'random_atpos',tread,'INT')
 if(tread==1) dtset%random_atpos=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ionmov',tread,'INT')
 if(tread==1) dtset%ionmov=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'intxc',tread,'INT')
 if(tread==1) dtset%intxc=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'iprcch',tread,'INT')
 if(tread==1) then
   dtset%iprcch=intarr(1)
 else
   if (dtset%paral_kgb==1) dtset%iprcch=6
 end if

!--For the moment LOTF does not use different from 2
 if(dtset%ionmov==23) then
#ifdef HAVE_LOTF
   if(dtset%iprcch/=2)then
     iprcch=2
     dtset%iprcch=iprcch
     write(message, '(4a)' )&
&     'When ionmov==23, iprcch must be 2.',ch10,&
&     'Set iprcch to 2.',ch10
     MSG_COMMENT(message)
   end if
#else
   dtset%ionmov=12
   write(message, '(4a)' )&
&   'LOTF is disabled, ionmov can not be 23.',ch10,&
&   'Set ionmov to 12.',ch10   
   MSG_COMMENT(message)
#endif
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'iprcel',tread,'INT')
 if(tread==1) dtset%iprcel=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'iprctfvw',tread,'INT')
 if(tread==1) dtset%iprctfvw=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'iprcfc',tread,'INT')
 if(tread==1) dtset%iprcfc=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irandom',tread,'INT')
 if(tread==1) dtset%irandom=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdden',tread,'INT')
 if(tread==1) dtset%irdden=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdpawden',tread,'INT')
 if(tread==1) dtset%irdpawden=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdhaydock',tread,'INT')
 if(tread==1) dtset%irdhaydock=intarr(1)
 
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdvdw',tread,'INT')
 if(tread==1) dtset%irdvdw=intarr(1)
 
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdwfk',tread,'INT')
 if(tread==1) dtset%irdwfk=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdwfkfine',tread,'INT')
 if(tread==1) dtset%irdwfkfine=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'iscf',tread,'INT')
 if(tread==1) then
   dtset%iscf=intarr(1)
 else if (dtset%optdriver==RUNL_RESPFN.and.dtset%iscf>=10) then
   dtset%iscf=dtset%iscf-10
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'isecur',tread,'INT')
 if(tread==1) dtset%isecur=intarr(1)

!Reading ixc must be immediately followed by reading xcname
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ixc',tread,'INT')
 if(tread==1) dtset%ixc=intarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),"xcname",tread_key,'KEY',key_value=key_value)
 if(tread_key==1)then
   if(tread==1)then
     write(message, '(5a)' )&
&     'ixc and xcname cannot be specified simultaneously',ch10,&
&     'for the same dataset.',ch10,&
&     'Action: check the input file.'
     MSG_ERROR(message)
   else
!    Note that xcname is a 'key' variable : its value is stored in keyw at output of intagm
     if(trim(key_value)=='PW92')dtset%ixc=7
     tread=1
   end if
 end if
 dtset%xclevel=0
 if( ( 1<=dtset%ixc .and. dtset%ixc<=10).or.(30<=dtset%ixc .and. dtset%ixc<=39) )dtset%xclevel=1 ! LDA
 if( (11<=dtset%ixc .and. dtset%ixc<=19).or.(23<=dtset%ixc .and. dtset%ixc<=29) )dtset%xclevel=2 ! GGA
 if( 20<=dtset%ixc .and. dtset%ixc<=22 )dtset%xclevel=3 ! ixc for TDDFT kernel tests
 if( dtset%ixc>=40 .and. dtset%ixc<=42 ) dtset%usefock=1 ! Hartree-Fock or hybrid functionals
 if( dtset%ixc>=41 .and. dtset%ixc<=42) dtset%xclevel=2 ! ixc for hybrids using GGA
#if defined HAVE_DFT_LIBXC
 if (dtset%ixc<0) then
   dtset%xclevel=1
   ii=-dtset%ixc/1000
   jj=xc_f90_family_from_id(ii)
   if (jj==XC_FAMILY_GGA.or.jj==XC_FAMILY_MGGA) dtset%xclevel=2
   ii=-dtset%ixc-ii*1000
   jj=xc_f90_family_from_id(ii)
   if (jj==XC_FAMILY_GGA.or.jj==XC_FAMILY_MGGA) dtset%xclevel=2
 end if
#endif

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_acutmin',tread,'DPR')
 if(tread==1) dtset%vdw_df_acutmin=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_aratio',tread,'DPR')
 if(tread==1) dtset%vdw_df_aratio=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_damax',tread,'DPR')
 if(tread==1) dtset%vdw_df_damax=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_damin',tread,'DPR')
 if(tread==1) dtset%vdw_df_damin=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_dcut',tread,'DPR')
 if(tread==1) dtset%vdw_df_dcut=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_dratio',tread,'DPR')
 if(tread==1) dtset%vdw_df_dratio=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_dsoft',tread,'DPR')
 if(tread==1) dtset%vdw_df_dsoft=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_gcut',tread,'DPR')
 if(tread==1) dtset%vdw_df_gcut=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_ndpts',tread,'INT')
 if(tread==1) dtset%vdw_df_ndpts=intarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_ngpts',tread,'INT')
 if(tread==1) dtset%vdw_df_ngpts=intarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_nqpts',tread,'INT')
 if(tread==1) dtset%vdw_df_nqpts=intarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_nrpts',tread,'INT')
 if(tread==1) dtset%vdw_df_nrpts=intarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_nsmooth',tread,'INT')
 if(tread==1) dtset%vdw_df_nsmooth=intarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_phisoft',tread,'DPR')
 if(tread==1) dtset%vdw_df_phisoft=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_qcut',tread,'DPR')
 if(tread==1) dtset%vdw_df_qcut=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_qratio',tread,'DPR')
 if(tread==1) dtset%vdw_df_qratio=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_rcut',tread,'DPR')
 if(tread==1) dtset%vdw_df_rcut=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_rsoft',tread,'DPR')
 if(tread==1) dtset%vdw_df_rsoft=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_tolerance',tread,'DPR')
 if(tread==1) dtset%vdw_df_tolerance=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_threshold',tread,'DPR',ds_input)
 if(tread==1) dtset%vdw_df_threshold=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_tweaks',tread,'INT')
 if(tread==1) dtset%vdw_df_tweaks=intarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_df_zab',tread,'DPR')
 if(tread==1) dtset%vdw_df_zab=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_xc',tread,'INT')
 if(tread==1) dtset%vdw_xc=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'positron',tread,'INT')
 if(tread==1) dtset%positron=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ixcpositron',tread,'INT')
 if(tread==1) dtset%ixcpositron=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'posnstep',tread,'INT')
 if(tread==1) dtset%posnstep=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'posocc',tread,'DPR')
 if(tread==1) dtset%posocc=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'frzfermi',tread,'INT')
 if(tread==1) dtset%frzfermi=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nqpt',tread,'INT')
 if(tread==1) dtset%nqpt=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ieig2rf',tread,'INT')
 if(tread==1) dtset%ieig2rf=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'smdelta',tread,'INT')
 if(tread==1) dtset%smdelta=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bdeigrf',tread,'INT')
 if(tread==1) dtset%bdeigrf=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'restartxf',tread,'INT')
 if(tread==1) dtset%restartxf=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'optcell',tread,'INT')
 if(tread==1) dtset%optcell=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdwfq',tread,'INT')
 if(tread==1) dtset%irdwfq=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ird1den',tread,'INT')
 if(tread==1) dtset%ird1den=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ird1wf',tread,'INT')
 if(tread==1) dtset%ird1wf=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdddk',tread,'INT')
 if(tread==1) dtset%irdddk=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'kptopt',tread,'INT')
 if(tread==1) dtset%kptopt=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'chkexit',tread,'INT')
 if(tread==1) dtset%chkexit=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ikhxc',tread,'INT')
 if(tread==1) dtset%ikhxc=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nbdblock',tread,'INT')
 if(tread==1) dtset%nbdblock=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nbdbuf',tread,'INT')
 if(tread==1)then
   dtset%nbdbuf=intarr(1)
 else
   if(response/=1 .and. dtset%iscf<0)dtset%nbdbuf=2*dtset%nspinor
   if(response==1 .and. 3<=occopt .and. occopt<=8 )dtset%nbdbuf=2*dtset%nspinor
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'localrdwf',tread,'INT')
 if(tread==1) dtset%localrdwf=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'magconon',tread,'INT')
 if(tread==1) dtset%magconon=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'magcon_lambda',tread,'DPR')
 if(tread==1) dtset%magcon_lambda=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'optforces',tread,'INT')
 if(tread==1) dtset%optforces=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'optstress',tread,'INT')
 if(tread==1) dtset%optstress=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'optnlxccc',tread,'INT')
 if(tread==1) dtset%optnlxccc=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nberry',tread,'INT')
 if(tread==1) dtset%nberry=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,2*nsppol,string(1:lenstr),'bdberry',tread,'INT')
 if(tread==1) then
   dtset%bdberry(1)=intarr(1); dtset%bdberry(2)=intarr(2)
   if(nsppol==2)then
     dtset%bdberry(3)=intarr(3); dtset%bdberry(4)=intarr(4)
   end if
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'berrystep',tread,'INT')
 if(tread==1) dtset%berrystep=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'delayperm',tread,'INT')
 if(tread==1) dtset%delayperm=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'signperm',tread,'INT')
 if(tread==1) dtset%signperm=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nbandkss',tread,'INT')
 if(tread==1) dtset%nbandkss=intarr(1)
 if ( dtset%usedmft > 0 .and. dtset%usepawu >= 0 .and. dtset%nbandkss==0) dtset%usepawu = 10

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'npwkss',tread,'INT')
 if(tread==1) dtset%npwkss=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'kssform',tread,'INT')
 if(tread==1) dtset%kssform=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'td_mexcit',tread,'INT')
 if(tread==1) dtset%td_mexcit=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'npweps',tread,'INT')
 if(tread==1) dtset%npweps=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'npvel',tread,'INT')
 if(tread==1) dtset%npvel=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'npulayit',tread,'INT')
 if(tread==1) dtset%npulayit=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'diismemory',tread,'INT')
 if(tread==1) dtset%diismemory=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'goprecon',tread,'INT')
 if(tread==1) dtset%goprecon=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'goprecprm',tread,'DPR')
 if(tread==1) dtset%goprecprm(1:3)=dprarr(1:3)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nwfshist',tread,'INT')
 if(tread==1) dtset%nwfshist=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'npwsigx',tread,'INT')
 if(tread==1) dtset%npwsigx=intarr(1)

 call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'nscforder', tread, 'INT')
 if (tread == 1) dtset%nscforder = intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nsheps',tread,'INT')
 if(tread==1) dtset%nsheps=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nshsigx',tread,'INT')
 if(tread==1) dtset%nshsigx=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nomegasrd',tread,'INT')
 if(tread==1) dtset%nomegasrd=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,ntypat,string(1:lenstr),'corecs',tread,'DPR')
 if(tread==1)then
   dtset%corecs(1:ntypat)=dprarr(1:ntypat)
 end if

 call intagm(dprarr,intarr,jdtset,marr,ntypat,string(1:lenstr),'ptcharge',tread,'DPR')
 if(tread==1)then
   dtset%ptcharge(1:ntypat)=dprarr(1:ntypat)
 end if

 call intagm(dprarr,intarr,jdtset,marr,ntypat,string(1:lenstr),'quadmom',tread,'DPR')
 if(tread==1)then
   dtset%quadmom(1:ntypat)=dprarr(1:ntypat)
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawcpxocc',tread,'INT')
 if(tread==1) then
   dtset%pawcpxocc=intarr(1)
 else if (dtset%nspinor==2.and.(dtset%usepawu>=1.or.dtset%usedmft>0)) then
   dtset%pawcpxocc=2
 else if (dtset%pawspnorb>0.and.(dtset%kptopt<=0.or.dtset%kptopt>=3)) then
   if (dtset%optdriver/=RUNL_GSTATE.or.dtset%ionmov<6.or.dtset%iscf<10) dtset%pawcpxocc=2
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawcross',tread,'INT')
 if(tread==1) dtset%pawcross=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawecutdg',tread,'ENE')
 if(tread==1) then
   dtset%pawecutdg=dprarr(1)
 else if (dtset%accuracy==1.or.dtset%accuracy==2) then                           
   dtset%pawecutdg=dtset%ecut
 else if (dtset%accuracy==3.) then                           
   dtset%pawecutdg=1.2_dp*dtset%ecut
 else if (dtset%accuracy==4) then                           
   dtset%pawecutdg=1.5_dp*dtset%ecut
 else if (dtset%accuracy==5.or.dtset%accuracy==6) then                           
   dtset%pawecutdg=two*dtset%ecut
! NOT YET ACTIVATED
! else if (dtset%accuracy==0.and.usepaw==1.and.dtset%usewvl==0) then                           
!  dtset%pawecutdg=1.5_dp*dtset%ecut
 end if
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawlcutd',tread,'INT')
 if(tread==1) dtset%pawlcutd=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawlmix',tread,'INT')
 if(tread==1) dtset%pawlmix=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawmixdg',tread,'INT')
 if(tread==1) then
   dtset%pawmixdg=intarr(1)
 else if (dtset%npfft>1.and.usepaw==1) then
   dtset%pawmixdg=1
 else if (dtset%usewvl==1.and.usepaw==1) then
   dtset%pawmixdg=1
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawnhatxc',tread,'INT')
 if(tread==1) dtset%pawnhatxc=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawntheta',tread,'INT')
 if(tread==1) dtset%pawntheta=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawnphi',tread,'INT')
 if(tread==1) dtset%pawnphi=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawnzlm',tread,'INT')
 if(tread==1) dtset%pawnzlm=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawoptmix',tread,'INT')
 if(tread==1) then
   dtset%pawoptmix=intarr(1)
 else
   if (usepaw==1.and.dtset%iscf<10.and.(dtset%positron==1.or.dtset%positron<0)) dtset%pawoptmix=1
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawoptosc',tread,'INT')
 if(tread==1) dtset%pawoptosc=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawovlp',tread,'DPR')
 if(tread==1) dtset%pawovlp=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawprtdos',tread,'INT')
 if(tread==1) dtset%pawprtdos=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawprtvol',tread,'INT')
 if(tread==1) dtset%pawprtvol=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawprtwf',tread,'INT')
 if(tread==1) dtset%pawprtwf=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawprt_k',tread,'INT')
 if(tread==1) dtset%pawprt_k=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawprt_b',tread,'INT')
 if(tread==1) dtset%pawprt_b=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawstgylm',tread,'INT')
 if(tread==1) dtset%pawstgylm=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawsushat',tread,'INT')
 if(tread==1) dtset%pawsushat=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawusecp',tread,'INT')
 if(tread==1) dtset%pawusecp=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawxcdev',tread,'INT')
 if(tread==1) dtset%pawxcdev=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'spnorbscl',tread,'DPR')
 if(tread==1) dtset%spnorbscl=dprarr(1)

 if (dtset%usedmft>0) then
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmft_dc',tread,'INT')
   if(tread==1) dtset%dmft_dc=intarr(1)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmft_iter',tread,'INT')
   if(tread==1) dtset%dmft_iter=intarr(1)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmft_mxsf',tread,'DPR')
   if(tread==1) dtset%dmft_mxsf=dprarr(1)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmft_nwli',tread,'INT')
   if(tread==1) dtset%dmft_nwli=intarr(1)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmft_nwlo',tread,'INT')
   if(tread==1) dtset%dmft_nwlo=intarr(1)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmft_read_occnd',tread,'INT')
   if(tread==1) dtset%dmft_read_occnd=intarr(1)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmft_rslf',tread,'INT')
   if(tread==1) dtset%dmft_rslf=intarr(1)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmft_solv',tread,'INT')
   if(tread==1) dtset%dmft_solv=intarr(1)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmft_t2g',tread,'INT')
   if(tread==1) dtset%dmft_t2g=intarr(1)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmft_tollc',tread,'DPR')
   if(tread==1) dtset%dmft_tollc=dprarr(1)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftbandi',tread,'INT')
   if(tread==1) dtset%dmftbandi=intarr(1)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftbandf',tread,'INT')
   if(tread==1) dtset%dmftbandf=intarr(1)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftcheck',tread,'INT')
   if(tread==1) dtset%dmftcheck=intarr(1)
   if(dtset%dmft_solv==4.or.dtset%dmft_solv==5) then
     call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftqmc_n',tread,'DPR')
     if(tread==1) then 
       dtset%dmftqmc_n=dprarr(1)
     else
       write(message, '(5a)' )&
&       'When DFT+DMFT is activated and one of QMC solvers is used,', ch10, &
&       'dmftqmc_n MUST be defined.',ch10,&
&       'Action: add dmftqmc_n keyword in input file.'
       MSG_ERROR(message)
     end if
     call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftqmc_l',tread,'INT')
     if(tread==1) then 
       dtset%dmftqmc_l=intarr(1)
     else
       write(message, '(5a)' )&
&       'When DFT+DMFT is activated and one of QMC solvers is used,', ch10, &
&       'dmftqmc_l MUST be defined.',ch10,&
&       'Action: add dmftqmc_l keyword in input file.'
       MSG_ERROR(message)
     end if
     call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftqmc_seed',tread,'INT')
     if(tread==1) dtset%dmftqmc_seed=intarr(1)
     call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftqmc_therm',tread,'INT')
     if(tread==1) then 
       dtset%dmftqmc_therm=intarr(1)
     else
       write(message, '(5a)' )&
&       'When DFT+DMFT is activated and one of QMC solvers is used,', ch10, &
&       'dmftqmc_therm MUST be defined.',ch10,&
&       'Action: add dmftqmc_therm keyword in input file.'
       MSG_ERROR(message)
     end if
     if(dtset%dmft_solv==5) then
       call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftctqmc_check',tread,'INT')
       if(tread==1) dtset%dmftctqmc_check  =intarr(1)
       call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftctqmc_correl',tread,'INT')
       if(tread==1) dtset%dmftctqmc_correl =intarr(1)
       call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftctqmc_gmove',tread,'INT')
       if(tread==1) dtset%dmftctqmc_gmove  =intarr(1)
       call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftctqmc_grnns',tread,'INT')
       if(tread==1) dtset%dmftctqmc_grnns  =intarr(1)
       call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftctqmc_meas',tread,'INT')
       if(tread==1) dtset%dmftctqmc_meas   =intarr(1)
       call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftctqmc_mrka',tread,'INT')
       if(tread==1) dtset%dmftctqmc_mrka   =intarr(1)
       call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftctqmc_mov',tread,'INT')
       if(tread==1) dtset%dmftctqmc_mov    =intarr(1)
       call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmftctqmc_order',tread,'INT')
       if(tread==1) dtset%dmftctqmc_order  =intarr(1)
     end if
   end if
 end if

 if (dtset%usepawu>0.or.dtset%usedmft>0) then
   call intagm(dprarr,intarr,jdtset,marr,ntypat,string(1:lenstr),'upawu',tread,'ENE')
   if(tread==1)then
     dtset%upawu(1:ntypat)=dprarr(1:ntypat)
   end if
   call intagm(dprarr,intarr,jdtset,marr,ntypat,string(1:lenstr),'f4of2_sla',tread,'ENE')
   if(tread==1)then
     dtset%f4of2_sla(1:ntypat)=dprarr(1:ntypat)
   end if
   call intagm(dprarr,intarr,jdtset,marr,ntypat,string(1:lenstr),'f6of2_sla',tread,'ENE')
   if(tread==1)then
     dtset%f6of2_sla(1:ntypat)=dprarr(1:ntypat)
   end if
   call intagm(dprarr,intarr,jdtset,marr,ntypat,string(1:lenstr),'jpawu',tread,'ENE')
   if(tread==1)then
     dtset%jpawu(1:ntypat)=dprarr(1:ntypat)
   end if
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmatpuopt',tread,'INT')
   if(tread==1) dtset%dmatpuopt=intarr(1)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dmatudiag',tread,'INT')
   if(tread==1) dtset%dmatudiag=intarr(1)

   if (dtset%usedmatpu/=0) then
     if (dmatsize>0) then
       call intagm(dprarr,intarr,jdtset,marr,dmatsize,string(1:lenstr),'dmatpawu',tread,'DPR')
       if(tread==1) then
         ii=1;jj=0
         do iatom=1,natom
           lpawu=dtset%lpawu(dtset%typat(iatom))
           if (lpawu/=-1) then
             dtset%dmatpawu(1:2*lpawu+1,1:2*lpawu+1,1:nsppol*nspinor,ii)= &
&             reshape(dprarr(jj+1:jj+nsppol*nspinor*(2*lpawu+1)**2),(/2*lpawu+1,2*lpawu+1,nsppol*nspinor/))
             ii=ii+1;jj=jj+nsppol*nspinor*(2*lpawu+1)**2
           end if
         end do
       else
         write(message, '(3a)' )&
&         'When LDA/GGA+U is activated and usedmatpu/=0, dmatpawu MUST be defined.',ch10,&
&         'Action: add dmatpawu keyword in input file.'
         MSG_ERROR(message)
       end if
     end if
   end if
 end if

 if (dtset%macro_uj>0) then
   dtset%pawujat=minval(pack((/ (ii,ii=1,dtset%natom) /) ,dtset%lpawu(dtset%typat(:))>0))
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawujat',tread,'INT')
   if(tread==1) dtset%pawujat=intarr(1)

   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawujrad',tread,'LEN')
   if(tread==1) dtset%pawujrad=dprarr(1)

   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawujv',tread,'ENE')
   if(tread==1) dtset%pawujv=dprarr(1)
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawfatbnd',tread,'INT')
 if(tread==1) dtset%pawfatbnd=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'bxctmindg',tread,'DPR')
 if(tread==1) then
   dtset%bxctmindg=dprarr(1)
 else if(usepaw==1.and.dtset%accuracy==0)then
   dtset%bxctmindg=dtset%boxcutmin
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usexcnhat',tread,'INT')
 if(tread==1) dtset%usexcnhat=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'useylm',tread,'INT')
 if(tread==1) then
   dtset%useylm=intarr(1)
   if ((usepaw==1).and.(dtset%useylm==0)) then
     write(message, '(a,a,a,a,a)' )&
&     'Pseudopotential file is PAW format (pspcod=7 or 17) while',ch10,&
&     'input variable "useylm" has the incompatible value 0 !',ch10,&
&     'Action: change psp format or "useylm" value in your input file.'
     MSG_ERROR(message)
   end if
   if ((dtset%tfkinfunc==2).and.(dtset%useylm==0)) then
     write(message, '(a,a,a,a,a)' )&
&     'You are using recursion method (tfkinfunc=2)  while',ch10,&
&     'input variable "useylm" has the incompatible value 0 !',ch10,&
&     'Action: change  "useylm" value in your input file.'
     MSG_ERROR(message)
   end if
 end if
 if (usepaw==1) dtset%useylm=1
 if (usepaw==1 .and. dtset%usewvl==1) then
   dtset%useylm=0
 end if
 if (dtset%use_gpu_cuda==1) dtset%useylm=1

 if(dtset%tfkinfunc==2 .and. dtset%usewvl==0 ) then
   dtset%useylm=1
   dtset%userec=1
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ireadc',tread,'INT')
 if(tread==1) then
   write(message, '(a,a,a)' )&
&   'The use of the "ireadc" variable is forbidden since version 2.0.',ch10,&
&   'Action: replace "ireadc" by "irdwfk" in your input file.'
   MSG_ERROR(message)
 end if

 ionmov=dtset%ionmov ; iprcch=dtset%iprcch ; iscf=dtset%iscf ; nqpt=dtset%nqpt
 kptopt=dtset%kptopt; nberry=dtset%nberry ; berryopt=dtset%berryopt

!Dielectric real(dp) input variables
!Reading of diemix/diemixmag must be inserted after iprcel
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'diecut',tread,'ENE')
 if(tread==1) dtset%diecut=dprarr(1)
!Special treatment if iscf==-1
 if(iscf==-1) dtset%diecut=four*dtset%ecut
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dielng',tread,'LEN')
 if(tread==1) dtset%dielng=dprarr(1)
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'diemac',tread,'DPR')
 if(tread==1) dtset%diemac=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'diemix',tread,'DPR')
 if(tread==1) then
   dtset%diemix=dprarr(1)
 else
   if (mod(dtset%iprcel,100)>19) dtset%diemix=one
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'diemixmag',tread,'DPR')
 if(tread==1) then
   dtset%diemixmag=dprarr(1)
 else
   if (dtset%iscf<10.or.dtset%nspden==1.or.dtset%iprcel==0.or.(dtset%iprcel>70.and.dtset%iprcel<80)) then
     dtset%diemixmag=dtset%diemix
   else
     dtset%diemixmag=-dtset%diemix
   end if
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'diegap',tread,'ENE')
 if(tread==1) dtset%diegap=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'dielam',tread,'DPR')
 if(tread==1) dtset%dielam=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'maxestep',tread,'DPR')
 if(tread==1) dtset%maxestep=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'td_maxene',tread,'DPR')
 if(tread==1) dtset%td_maxene=dprarr(1)

!ACFD input variables
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'idyson',tread,'INT')
 if(tread==1) dtset%idyson=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ndyson',tread,'INT')
 if(tread==1) dtset%ndyson=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'intexact',tread,'INT')
 if(tread==1) dtset%intexact=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nbandsus',tread,'INT')
 if(tread==1) dtset%nbandsus=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ldgapp',tread,'INT')
 if(tread==1) dtset%ldgapp=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'suskxcrs',tread,'INT')
 if(tread==1) dtset%suskxcrs=intarr(1)

 if((iscf==5.or.iscf==6) .and. ionmov==4 .and. iprcch/=3 )then
   iprcch=3
   dtset%iprcch=iprcch
   write(message, '(3a)' )&
&   '  When ionmov==4 and iscf==5 or 6, iprcch must be 3.',ch10,&
&   '  Set iprcch to 3.'
   MSG_COMMENT(message)
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'mffmem',tread,'INT')
 if(tread==1) dtset%mffmem=intarr(1)

!Set default values of 0 for occupation numbers and k pt wts
 bantot=0

!nkpt and nband must be defined to execute following loop
 if ( tnband == 1 ) then
   do ikpt=1,nkpt*nsppol
     do ii=1,dtset%nband(ikpt)
       bantot=bantot+1
       dtset%occ_orig(bantot)=0.0_dp
     end do
   end do
 end if

!NON DEVLOP CASE
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nloalg',tread,'INT')
 if(tread==1) then
   dtset%nloalg(1)=mod(intarr(1),10)
   dtset%nloalg(5)=intarr(1)/10
 end if
!DEVLOP CASE. Note that the fifth component is not used now.
!call intagm(dprarr,intarr,jdtset,marr,5,string(1:lenstr),'nloalg',tread,'INT')
!if(tread==1) then
!dtset%nloalg(1)=intarr(1) ; dtset%nloalg(2)=intarr(2)
!dtset%nloalg(3)=intarr(3) ; dtset%nloalg(4)=intarr(4) ; dtset%nloalg(4)=intarr(5)
!end if
!ENDDEVLOP

!LOOP variables
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nline',tread,'INT')
 if(tread==1) dtset%nline=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nnsclo',tread,'INT')
 if(tread==1) dtset%nnsclo=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nstep',tread,'INT')
 if(tread==1) dtset%nstep=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ntime',tread,'INT')
 if(tread==1) dtset%ntime=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nfreqsus',tread,'INT')
 if(tread==1) dtset%nfreqsus=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nctime',tread,'INT')
 if(tread==1) dtset%nctime=intarr(1)
 
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'orbmag',tread,'INT')
 if(tread==1) dtset%orbmag=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ortalg',tread,'INT')
 if(tread==1) then
   dtset%ortalg=intarr(1)
 else if (dtset%wfoptalg>=10 .and. dtset%ortalg>0) then
   dtset%ortalg=-dtset%ortalg
 end if

!Print variables
 call intagm(dprarr,intarr,jdtset,marr,natom,string(1:lenstr),'prtatlist',tread,'INT')
 if(tread==1) then
   dtset%prtatlist=intarr(1:natom)
   call sort_int(natom,dtset%prtatlist,intarr)
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtbbb',tread,'INT')
 if(tread==1) dtset%prtbbb=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtbltztrp',tread,'INT')
 if(tread==1) dtset%prtbltztrp=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtcif',tread,'INT')
 if(tread==1) dtset%prtcif=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtcml',tread,'INT')
 if(tread==1) dtset%prtcml=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtcs',tread,'INT')
 if(tread==1) dtset%prtcs=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtden',tread,'INT')
 if(tread==1) dtset%prtden=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtdipole',tread,'INT')
 if(tread==1) dtset%prtdipole=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtdos',tread,'INT')
 if(tread==1) dtset%prtdos=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtdosm',tread,'INT')
 if(tread==1) dtset%prtdosm=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtefg',tread,'INT')
 if(tread==1) dtset%prtefg=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prteig',tread,'INT')
 if(tread==1) dtset%prteig=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtelf',tread,'INT')
 if(tread==1) dtset%prtelf=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtfc',tread,'INT')
 if(tread==1) dtset%prtfc=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtfsurf',tread,'INT')
 if(tread==1) dtset%prtfsurf=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtgden',tread,'INT')
 if(tread==1) dtset%prtgden=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtgeo',tread,'INT')
 if(tread==1) dtset%prtgeo=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtgkk',tread,'INT')
 if(tread==1) dtset%prtgkk=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtkden',tread,'INT')
 if(tread==1) dtset%prtkden=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtlden',tread,'INT')
 if(tread==1) dtset%prtlden=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtnabla',tread,'INT')
 if(tread==1) dtset%prtnabla=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtnest',tread,'INT')
 if(tread==1) dtset%prtnest=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtposcar',tread,'INT')
 if(tread==1) dtset%prtposcar=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtpot',tread,'INT')
 if(tread==1) dtset%prtpot=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtspcur',tread,'INT')
 if(tread==1) dtset%prtspcur=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtsuscep',tread,'INT')
 if(tread==1) dtset%prtsuscep=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtstm',tread,'INT')
 if(tread==1) dtset%prtstm=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prt1dm',tread,'INT')
 if(tread==1) dtset%prt1dm=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtvdw',tread,'INT')
 if(tread==1) dtset%prtvdw=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtvha',tread,'INT')
 if(tread==1) dtset%prtvha=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtvhxc',tread,'INT')
 if(tread==1) dtset%prtvhxc=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtvol',tread,'INT')
 if(tread==1) dtset%prtvol=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtvolimg',tread,'INT')
 if(tread==1) dtset%prtvolimg=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtvxc',tread,'INT')
 if(tread==1) dtset%prtvxc=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtwant',tread,'INT')
 if(tread==1) dtset%prtwant=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtwf',tread,'INT')
 if(tread==1) then
   dtset%prtwf=intarr(1)
   if ((dtset%tfkinfunc==2).and.(dtset%prtwf==1)) then
     write(message, '(a,a,a,a,a)' )&
&     'You are using recursion method (tfkinfunc=2)  while',ch10,&
&     'input variable "prtwf" has the incompatible value 1 !',ch10,&
&     'Action: change  "prtwf=0" value in your input file.'
     MSG_ERROR(message)
   end if
 end if
 if (dtset%tfkinfunc==2) dtset%prtwf=0

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtxml',tread,'INT')
 if(tread==1) dtset%prtxml=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ratsph_extra',tread,'DPR')
 if(tread==1) dtset%ratsph_extra=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'qprtrb',tread,'INT')
 if(tread==1) dtset%qprtrb(:)=intarr(1:3)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'strprecon',tread,'DPR')
 if(tread==1) dtset%strprecon=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'xc_denpos',tread,'DPR')
 if(tread==1) dtset%xc_denpos=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'xc_tb09_c',tread,'DPR')
 if(tread==1) dtset%xc_tb09_c=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,3*dtset%natsph_extra,string(1:lenstr),'xredsph_extra',tread,'DPR')
 if(tread==1) dtset%xredsph_extra=reshape(dprarr(1:3*dtset%natsph_extra), (/3,dtset%natsph_extra/))

!WVL - Wavelets related values

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'wvl_hgrid',tread,'DPR')
 if(tread==1) dtset%wvl_hgrid=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'wvl_crmult',tread,'DPR')
 if(tread==1) dtset%wvl_crmult=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'wvl_frmult',tread,'DPR')
 if(tread==1) dtset%wvl_frmult=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'wvl_ngauss',tread,'INT')
 if(tread==1) dtset%wvl_ngauss(1:2)=intarr(1:2)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'wvl_nprccg',tread,'INT')
 if(tread==1) dtset%wvl_nprccg=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tl_nprccg',tread,'INT')
 if(tread==1) dtset%tl_nprccg=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tl_radius',tread,'DPR')
 if(tread==1) dtset%tl_radius=dprarr(1)

!Wannier90 interface related variables
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'w90iniprj',tread,'INT')
 if(tread==1) then
   dtset%w90iniprj=intarr(1)
   if ( usepaw == 0 .and. ( dtset%w90iniprj /= 2  .and. dtset%w90iniprj/=0 .and. dtset%w90iniprj /= 1 )) then
     write(message, '(a,i4,2a)' )&
&     'w90iniprj should be set to 0, 1 or 2, however, it was ',dtset%w90iniprj,ch10,&
&     'Action: check the values of w90iniprj.'
     MSG_ERROR(message)
   end if
   if ( usepaw == 1 .and. ( dtset%w90iniprj < 2 .or. dtset%w90iniprj>6 ) &
&   .and. ( dtset%w90iniprj /= 1 .and. dtset%w90iniprj/=0 )) then
     write(message, '(a,i4,2a)' )&
&     'In the PAW case, the only valid values for w90iniprj are 0, 1, 2, 5 and 6 however, it was ',dtset%w90iniprj,ch10,&
&     'Action: check the values of w90iniprj.'
     MSG_ERROR(message)
   end if
 end if
!w90prtunk
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'w90prtunk',tread,'INT')
 if(tread==1) then
   dtset%w90prtunk=intarr(1)
   if ( dtset%w90prtunk < 0 ) then
     write(message, '(4a)' )&
&     'w90prtunk should be greater or equal to zero, however, it was ',dtset%w90prtunk,ch10,&
&     'Action: check the values of w90prtunk.'
     MSG_ERROR(message)
   end if
 end if

!Wannier90 - GW quasiparticle interface
 if(dtset%prtwant==3) then
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'getqps',tread,'INT')
   if(tread==1) dtset%getqps=intarr(1)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'irdqps',tread,'INT')
   if(tread==1) dtset%irdqps=intarr(1)
 end if

!van der Waals with DFT-D2 approach
 if(dtset%vdw_xc==5) then
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_tol',tread,'DPR')
   if(tread==1) dtset%vdw_tol=dprarr(1)
 end if

!van der Waals with Wannier functions (Silvestrelli's approach versions 1,2 and 3)
 if(dtset%vdw_xc==10.or.dtset%vdw_xc==11.or.dtset%vdw_xc==12.or.dtset%vdw_xc==14) then
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'vdw_nfrag',tread,'INT')
   if(tread==1)then
     dtset%vdw_nfrag=intarr(1)
   end if

   call intagm(dprarr,intarr,jdtset,marr,dtset%natom,string(1:lenstr),'vdw_typfrag',tread,'INT')
   if(tread==1)then
     dtset%vdw_typfrag(1:dtset%natom)=intarr(1:dtset%natom)
     do ii=1,dtset%natom
       if(dtset%vdw_typfrag(ii)<0)then
         write(message, '(a,a,a,i0,a,i0,a,a)' )&
&         'All the components of vdw_typfrag must be greater than 0.',ch10,&
&         'The component',ii,' is equal to ',dtset%vdw_typfrag(ii),ch10,&
&         'Action: check the values of vdw_typfrag.'
         MSG_ERROR(message)
       end if
     end do
   end if
 end if
 if(dtset%vdw_xc==10.or.dtset%vdw_xc==11.or.dtset%vdw_xc==12.or.dtset%vdw_xc==14) then
   call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'vdw_supercell',tread,'INT')
   if(tread==1)then
     dtset%vdw_supercell(1:3)=intarr(1:3)
     if (dtset%supercell(1)<zero.and.dtset%supercell(2)<zero) then
       write(message, '(7a)' )&
&       ' only one component of vdw_supercell could be < 0, however, it was ',ch10,&
&       dtset%vdw_supercell(1),dtset%vdw_supercell(2),dtset%vdw_supercell(3),ch10,&
&       'Action: check the components of vdw_supercell.'
       MSG_ERROR(message)
     end if
     if (dtset%supercell(2)<zero.and.dtset%supercell(3)<zero) then
       write(message, '(7a)' )&
&       'only one component of vdw_supercell could be < 0, however, it was ',ch10,&
&       dtset%vdw_supercell(1),dtset%vdw_supercell(2),dtset%vdw_supercell(3),ch10,&
&       'Action: check the components of vdw_supercell.'
       MSG_ERROR(message)
     end if
     if (dtset%supercell(1)<zero.and.dtset%supercell(3)<zero) then
       write(message, '(7a)' )&
&       'only one component of vdw_supercell could be < 0, however, it was ',ch10,&
&       dtset%vdw_supercell(1),dtset%vdw_supercell(2),dtset%vdw_supercell(3),ch10,&
&       'Action: check the components of vdw_supercell.'
       MSG_ERROR(message)
     end if
   end if

 end if !vdw_xc==10 or vdw_xc==11 or vdw_xc==12 or vdw_xc==14

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ucrpa',tread,'INT')
 if(tread==1) dtset%ucrpa=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'ucrpa_bands',tread_key,'INT')
 if(tread_key==1) then
   if(dtset%nsym/=1) then
     write(message,"(a,i0,3a)")&
&     '  Value of nsym different from 1 when ucrpa_bands is used is under test ',dtset%nsym,&
&     '(because symmetry is not yet used)',ch10,&
&     '  Action : check your calculation  with nsym=1'
     MSG_WARNING(message)
   endif
   dtset%ucrpa_bands(1:2)=intarr(1:2)
 end if

 call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'ucrpa_window',tread,'ENE')
 if(tread==1) then
   if(dtset%nsym/=1) then
     write(message,*)&
&     '  Value of nsym different from 1 when ucrpa_windows is used is under test ',dtset%nsym,&
&     '  (because symmetry is not yet used)',ch10,&
&     '  Action : check your calculation  with nsym=1'
     MSG_WARNING(message)
   end if
   dtset%ucrpa_window(1:2)=dprarr(1:2)
 end if

 if(tread==1.and.tread_key==1) then
   write(message, '(a,a,a,a,a)' )&
&   ' ucrpa_bands and ucrpa_window cannot be specified simultaneously',ch10,&
&   ' for the same dataset.',ch10,&
&   ' Action : check the input file.'
   MSG_ERROR(message)
 endif

!Tolerance variables
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tolimg',tread,'ENE')
 if(tread==1) dtset%tolimg=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tolmxf',tread,'DPR')
 if(tread==1) dtset%tolmxf=dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tolrde',tread,'DPR')
 if(tread==1) dtset%tolrde=dprarr(1)

!find which tolXXX are defined generically and for this jdtset
 tolwfr_=zero
 toldfe_=zero
 toldff_=zero
 tolrff_=zero
 tolvrs_=zero
 itol=0
 itol_gen=0
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tolwfr',tread,'DPR',ds_input)
 if(tread==1) then
   if (ds_input == 0) then
     tolwfr_=dprarr(1)
     if(abs(dprarr(1))>tiny(0._dp))itol_gen=itol_gen+1
   else
     dtset%tolwfr=dprarr(1)
     if(abs(dprarr(1))>tiny(0._dp))itol=itol+1
   end if
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'toldff',tread,'DPR',ds_input)
 if(tread==1) then
   if (ds_input == 0) then
     toldff_=dprarr(1)
     if(abs(dprarr(1))>tiny(0._dp))itol_gen=itol_gen+1
   else
     dtset%toldff=dprarr(1)
     if(abs(dprarr(1))>tiny(0._dp))itol=itol+1
   end if
   dtset%optforces=1
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tolrff',tread,'DPR',ds_input)
 if(tread==1) then
   if (ds_input == 0) then
     tolrff_=dprarr(1)
     if(abs(dprarr(1))>tiny(0._dp))itol_gen=itol_gen+1
   else
     dtset%tolrff=dprarr(1)
     if(abs(dprarr(1))>tiny(0._dp))itol=itol+1
   end if
   dtset%optforces=1
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'toldfe',tread,'ENE',ds_input)
 if(tread==1) then
   if (ds_input == 0) then
     toldfe_=dprarr(1)
     if(abs(dprarr(1))>tiny(0._dp))itol_gen=itol_gen+1
   else
     dtset%toldfe=dprarr(1)
     if(abs(dprarr(1))>tiny(0._dp))itol=itol+1
   end if
 end if

 if((dtset%accuracy>0).and.(itol>0.or.itol_gen>0)) tolvrs_=zero

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tolvrs',tread,'DPR',ds_input)
 if(tread==1) then
   if (ds_input == 0) then
     tolvrs_=dprarr(1)
     if(abs(dprarr(1))>tiny(0._dp))itol_gen=itol_gen+1
   else
     dtset%tolvrs=dprarr(1)
     if(abs(dprarr(1))>tiny(0._dp))itol=itol+1
   end if
 end if

!check for multiple definitions of tolXXX for the present dataset
 if (itol > 1 .or. itol_gen > 1) then
   write(message, '(3a)' )&
&   'Only one of the tolXXX variables may be defined at once.',ch10,&
&   'Action: check values of tolvrs, toldfe, tolrff, tolwfr, and toldff.'
   MSG_ERROR(message)
 end if

!if no value is given for jdtset, use defaults
 if (itol == 0 .and. itol_gen == 1) then
   dtset%tolwfr=tolwfr_
   dtset%toldfe=toldfe_
   dtset%toldff=toldff_
   dtset%tolrff=tolrff_
   dtset%tolvrs=tolvrs_
 end if

!Tolerance variables for electrons-positron tc-dft
 toldfe_=zero;toldff_=zero
 itol=0;itol_gen=0
 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'postoldff',tread,'DPR',ds_input)
 if(tread==1) then
   if (ds_input == 0) then
     toldff_=dprarr(1)
     if(abs(dprarr(1))>tiny(0._dp))itol_gen=itol_gen+1
   else
     dtset%postoldff=dprarr(1)
     if(abs(dprarr(1))>tiny(0._dp))itol=itol+1
   end if
   dtset%postoldfe=zero
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'postoldfe',tread,'ENE',ds_input)
 if(tread==1) then
   if (ds_input == 0) then
     toldfe_=dprarr(1)
     if(abs(dprarr(1))>tiny(0._dp))itol_gen=itol_gen+1
   else
     dtset%postoldfe=dprarr(1)
     if(abs(dprarr(1))>tiny(0._dp))itol=itol+1
   end if
 end if
 if (itol > 1.or.itol_gen >1) then
   write(message, '(3a)' )&
&   'Only one of the postolXXX variables may be defined at once.',ch10,&
&   'Action: check values of postoldfe and postoldff.'
   MSG_ERROR(message)
 end if
 if (itol==0.and.itol_gen==1) then
   dtset%postoldfe=toldfe_
   dtset%postoldff=toldff_
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usekden',tread,'INT')
 if(tread==1) dtset%usekden=intarr(1)

 call intagm(dprarr,intarr,jdtset,marr,2,string(1:lenstr),'vprtrb',tread,'ENE')
 if(tread==1) dtset%vprtrb(:)=dprarr(1:2)

!TODO : should this test be just on natsph without caring about the others? natsph might be used for other purposes
 if((dtset%pawfatbnd>0.or.dtset%prtdos==3) .and. dtset%natsph>0)then
   call intagm(dprarr,intarr,jdtset,marr,dtset%natsph,string(1:lenstr),'iatsph',tread,'INT')
   if(tread==1) dtset%iatsph(1:dtset%natsph)=intarr(1:dtset%natsph)
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'prtdensph',tread,'INT')
 if(tread==1) then
   dtset%prtdensph=intarr(1)
 else if (dtset%nsppol==1.and.dtset%nspden==2) then
   dtset%prtdensph=1
 end if

!Initialize atvshift
 if(dtset%natvshift>0)then
   call intagm(dprarr,intarr,jdtset,marr,dtset%natvshift*nsppol*dtset%natom,string(1:lenstr),'atvshift',tread,'ENE')
   if(tread==1) dtset%atvshift(1:dtset%natvshift,1:nsppol,1:dtset%natom)=&
&   reshape(dprarr(1:dtset%natvshift*nsppol*dtset%natom),(/dtset%natvshift,nsppol,dtset%natom/))
 end if

!Initialize wtatcon
 if(dtset%nconeq>0)then

!  Read and check natcon
   ABI_ALLOCATE(natcon,(dtset%nconeq))
   call intagm(dprarr,intarr,jdtset,marr,dtset%nconeq,string(1:lenstr),'natcon',tread,'INT')
   if(tread==1)then
     natcon(:)=intarr(1:dtset%nconeq)
   else
     write(message, '(3a)' )&
&     'When nconeq is positive, natcon MUST be defined.',ch10,&
&     'Action: check the values of nconeq and natcon.'
     MSG_ERROR(message)
   end if
   do ii=1,dtset%nconeq
     if(natcon(ii)<0)then
       write(message, '(a,a,a,i4,a,i4,a,a)' )&
&       'All the components of natcon must be greater than 0.',ch10,&
&       'The component',ii,' is equal to ',natcon(ii),ch10,&
&       'Action: check the values of natcon.'
       MSG_ERROR(message)
     end if
   end do
   niatcon=sum(natcon(:))

!  Read and check iatcon
   ABI_ALLOCATE(iatcon,(niatcon))
   call intagm(dprarr,intarr,jdtset,marr,niatcon,string(1:lenstr),'iatcon',tread,'INT')
   if(tread==1)then
     iatcon(:)=intarr(1:niatcon)
   else
     write(message, '(3a)' )&
&     'When nconeq is positive, natcon MUST be defined.',ch10,&
&     'Action: check the values of nconeq and natcon.'
     MSG_ERROR(message)
   end if
   do ii=1,niatcon
     if(iatcon(ii)<0)then
       write(message, '(a,a,a,i4,a,i4,a,a)' )&
&       'All the components of iatcon must be greater than 0.',ch10,&
&       'The component',ii,' is equal to ',iatcon(ii),ch10,&
&       'Action: check the values of iatcon.'
       MSG_ERROR(message)
     end if
   end do

!  Read wtatcon, and unfold it.
   call intagm(dprarr,intarr,jdtset,marr,3*niatcon,string(1:lenstr),'wtatcon',tread,'DPR')
   if(tread/=1)then
     write(message, '(3a)' )&
&     'When nconeq is positive, wtatcon MUST be defined.',ch10,&
&     'Action: check the values of nconeq and wtatcon.'
     MSG_ERROR(message)
   end if
   index_nat=0
   do ii=1,dtset%nconeq
     do jj=1,natcon(ii)
       dtset%wtatcon(1:3,iatcon(jj+index_nat),ii)=dprarr(1+3*(jj+index_nat-1):3+3*(jj+index_nat-1))
     end do
     index_nat=index_nat+natcon(ii)
   end do

   ABI_DEALLOCATE(iatcon)
   ABI_DEALLOCATE(natcon)
 end if

!Initialize the list of k points, as well as wtk and istwfk
 intimage=1 ; if(nimage>2)intimage=(1+nimage)/2
 call invacuum(jdtset,lenstr,natom,dtset%rprimd_orig(1:3,1:3,intimage),string,vacuum,dtset%xred_orig(1:3,1:natom,intimage))

!reset rfdir to 1 1 1 in case it was left at the default of 0 0 0 for berryopt -1
 if( (dtset%berryopt == -1) .AND. &
& (dtset%rfdir(1) == 0) .AND. &
& (dtset%rfdir(2) == 0) .AND. &
& (dtset%rfdir(3) == 0) ) dtset%rfdir(1:3) = 1 

!In case of a Berryphase calculation or orbital magnetism, put response = 1 
!in order to set istwfk = 1 at all k-points
 if ( (dtset%berryopt/=0) .OR. (dtset%orbmag /= 0) ) then
   response = 1
 end if

 nsym=dtset%nsym 
 ii=0;if (mod(dtset%wfoptalg,10)==4) ii=2
 if (dtset%ngfft(7)==314) ii=1

 call inkpts(bravais,dtset%chksymbreak,iout,iscf,dtset%istwfk(1:nkpt),jdtset,&
& dtset%kpt(:,1:nkpt),kptopt,dtset%kptnrm,dtset%kptrlatt,kptrlen,lenstr,nsym,&
& nkpt,nqpt,dtset%ngkpt,dtset%nshiftk,dtset%nshiftk_orig,dtset%shiftk_orig,nsym,&
& occopt,dtset%qptn,response,dtset%rprimd_orig(1:3,1:3,intimage),dtset%shiftk,string,&
& dtset%symafm(1:nsym),dtset%symrel(:,:,1:nsym),vacuum,dtset%wtk(1:nkpt),&
& impose_istwf_1=ii)

 dtset%kptrlen=kptrlen

 dtset%kptns(:,1:nkpt)=dtset%kpt(:,1:nkpt)/dtset%kptnrm
 if(nqpt>=1 .and. dtset%optdriver/=RUNL_RESPFN)then
   dtset%kptns(1,1:nkpt)=dtset%kptns(1,1:nkpt)+dtset%qptn(1)
   dtset%kptns(2,1:nkpt)=dtset%kptns(2,1:nkpt)+dtset%qptn(2)
   dtset%kptns(3,1:nkpt)=dtset%kptns(3,1:nkpt)+dtset%qptn(3)
 end if

 if(dtset%nkptgw>0) then
   ! Read nkptgw and bdgw.
   call intagm(dprarr,intarr,jdtset,marr,2*dtset%nkptgw*dtset%nsppol,string(1:lenstr),'bdgw',tread,'INT')
   if(tread==1) then
     dtset%bdgw(1:2,1:dtset%nkptgw,1:dtset%nsppol)=reshape(intarr(1:2*dtset%nkptgw*dtset%nsppol),(/2,dtset%nkptgw,dtset%nsppol/))
   end if

   call intagm(dprarr,intarr,jdtset,marr,3*dtset%nkptgw,string(1:lenstr),'kptgw',tread,'DPR')
   if(tread==1) dtset%kptgw(1:3,1:dtset%nkptgw) = reshape(dprarr(1:3*dtset%nkptgw),(/3,dtset%nkptgw/))
 end if

 if(dtset%nqptdm>0) then
   call intagm(dprarr,intarr,jdtset,marr,3*dtset%nqptdm,string(1:lenstr),'qptdm',tread,'DPR')
   if(tread==1) dtset%qptdm(1:3,1:dtset%nqptdm)=reshape(dprarr(1:3*dtset%nqptdm),(/3,dtset%nqptdm/))
 end if

 if(dtset%npvel>0) then
   call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'pvelmax',tread,'DPR')
   if(tread==1) dtset%pvelmax(1:3)=dprarr(1:3)
 end if

!Only read occ if (iscf >0 or iscf=-1 or iscf=-3) and (occopt==0 or occopt==2)
 if  (iscf>=0.or.iscf==-1.or.iscf==-3)  then
   if (occopt==2 .and. getocc==0) then
!    Read occ(nband(kpt)*nkpt*nsppol) explicitly
     call wrtout(std_out,' invars2: reading occ(nband*nkpt*nsppol) explicitly','COLL')
     call intagm(dprarr,intarr,jdtset,marr,bantot,string(1:lenstr),'occ',tread,'DPR')
     if(tread==1) dtset%occ_orig(1:bantot)=dprarr(1:bantot)
   else if(occopt==0) then
!    Read usual occupancy--same for all k points and spins
     call intagm(dprarr,intarr,jdtset,marr,dtset%nband(1),string(1:lenstr),'occ',tread,'DPR')
     if(tread==1) dtset%occ_orig(1:dtset%nband(1))=dprarr(1:dtset%nband(1))
!    Fill in full occ array using input values for each k and spin
!    (make a separate copy for each k point and spin)
     do ikpt=1,nkpt*nsppol
       dtset%occ_orig(1+(ikpt-1)*dtset%nband(1):ikpt*dtset%nband(1))=&
&       dtset%occ_orig(1:dtset%nband(1))
     end do
   end if
 end if

 if(dtset%jellslab/=0)then
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'slabwsrad',tread,'LEN')
   if(tread==1) dtset%slabwsrad=dprarr(1)

!  Update number of electrons taking into account jellium
!  Suppose that the cell has z axis perpendicular to x and y axes.
!  This will be checked later
   areaxy=abs(dtset%rprimd_orig(1,1,1)*dtset%rprimd_orig(2,2,1)-dtset%rprimd_orig(1,2,1)*dtset%rprimd_orig(2,1,1))
   rhoavg=three/(four_pi*dtset%slabwsrad**3)
   nelectjell=areaxy*(dtset%slabzend-dtset%slabzbeg)*rhoavg
   charge=charge-nelectjell
 end if

!Initialize occ if occopt==1 or 3 ... 8,
!while if getocc/=0, make a fake initialization
!If iscf>0, check the charge of the system, and compute nelect.
 occopt_tmp=occopt
 if(getocc/=0)occopt_tmp=1
 call chkneu(charge,dtset,occopt_tmp)

!Initialize Berry phase vectors
!Should check that nberry is smaller than 20
 if(berryopt>0 .and. nberry>0)then
   call intagm(dprarr,intarr,jdtset,marr,3*nberry,string(1:lenstr),'kberry',tread,'INT')
   if(tread==1) dtset%kberry(1:3,1:nberry)=reshape(intarr(1:3*nberry), (/3,nberry/))
 end if

 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ddamp',tread,'DPR')
 if (tread==1) dtset%ddamp = dprarr(1)

 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'polcen',tread,'DPR')
 if (tread==1) dtset%polcen(1:3) = dprarr(1:3)

 call intagm(dprarr,intarr,jdtset,marr,ntypat,string(1:lenstr),'densty',tread,'DPR')
 if(tread==1) dtset%densty(1:ntypat,1)=dprarr(1:ntypat)

 call intagm(dprarr,intarr,jdtset,marr,npsp,string(1:lenstr),'so_psp',tread,'INT')
 if(tread==1)then
   dtset%so_psp(1:npsp)=intarr(1:npsp)
 end if

!@ HF or hybrid functionals
 if (dtset%usefock==1) then

! integer scalars (used to dimension specific arrays)
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'cgtyphf',tread,'INT')
   if(tread==1) dtset%cgtyphf=intarr(1)
   
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nbandhf',tread,'INT')
   if(tread==1) then 
     dtset%nbandhf=intarr(1)
   else
!  First compute the total valence charge
     zval=0.0_dp
     do iatom=1,natom
       zval=zval+dtset%ziontypat(dtset%typat(iatom))
     end do
!  Only take into account negative charge, to compute maximum number of bands
     if(charge > 0.0_dp)charge=0.0_dp
     zelect=zval-charge
     
!  Then select the minimum number of bands, and add the required number
!  Note that this number might be smaller than the one computed
!  by a slightly different formula in invars1
     dtset%nbandhf=dtset%nspinor*&
&     ((ceiling(zelect-1.0d-10)+1)/2 + ceiling( fband*natom - 1.0d-10 ))-1
!  More precisely, nbandhf = default value for nband (the number of occupied bands + 1) -1 = the number of occupied bands. 
   end if

   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nnsclohf',tread,'INT')
   if(tread==1) dtset%nnsclohf=intarr(1)

   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nkpthf',tread,'INT')
   if (dtset%kptopt==3) then 
     nkpt_hf=dtset%nkpt
   else 
!* Default value for nkpthf is the total number of k-points in BZ.

!* Compute the maximum number of k-points mkpt in the reciprocal space unit cell 
!* (same formula as used in getkgrid)
     mkpt=dtset%nshiftk*(dtset%kptrlatt(1,1)*dtset%kptrlatt(2,2)*dtset%kptrlatt(3,3) &
&     +dtset%kptrlatt(1,2)*dtset%kptrlatt(2,3)*dtset%kptrlatt(3,1) &
&     +dtset%kptrlatt(1,3)*dtset%kptrlatt(2,1)*dtset%kptrlatt(3,2) &
&     -dtset%kptrlatt(1,2)*dtset%kptrlatt(2,1)*dtset%kptrlatt(3,3) &
&     -dtset%kptrlatt(1,3)*dtset%kptrlatt(2,2)*dtset%kptrlatt(3,1) &
&     -dtset%kptrlatt(1,1)*dtset%kptrlatt(2,3)*dtset%kptrlatt(3,2))

     ABI_ALLOCATE(kptns_hf,(3,mkpt))
     kptns_hf=0.0_dp
!* Generate all the k-points in BZ (Monkhorst-Pack grid) 

!* brav=1 to treat all Bravais lattices ; iout=0 since we do not want any output ; option=0 since we consider k-points 
     call smpbz(1,0,dtset%kptrlatt,mkpt,nkpt_hf,dtset%nshiftk,0,dtset%shiftk,kptns_hf)
!* kptns_hf contains the special k points obtained by the Monkhorst & Pack method, in reduced coordinates. (output)
     ABI_DEALLOCATE(kptns_hf)

     if((tread==1).AND.(nkpt_hf/=intarr(1))) then 
       write(message, '(3a)' )&
&       '  This option is not yet fully implemented.',&
&       '  The default value should be equal to the total number of k-point in BZ.'
       MSG_ERROR(message)
     end if
     dtset%nkpthf=nkpt_hf
     
   end if 
 end if ! usefock

!LOTF variables
#if defined HAVE_LOTF
 if(dtset%ionmov==23) then
   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'lotf_classic',tread,'INT')
   if(tread==1) dtset%lotf_classic=intarr(1)

   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'lotf_nitex',tread,'INT')
   if(tread==1) dtset%lotf_nitex=intarr(1)

   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'lotf_nneigx',tread,'INT')
   if(tread==1) dtset%lotf_nneigx=intarr(1)

   call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'lotf_version',tread,'INT')
   if(tread==1) dtset%lotf_version=intarr(1)
 end if
#endif

!This is to be removed after a few version (from 5.4)
 call intagm(dprarr,intarr,jdtset,marr,ntypat,string(1:lenstr),'so_typat',tread_alt,'INT')
 if(tread_alt==1)then
   write(message, '(3a)' )&
&   'The use of the input variable so_typat has been discontinued in v5.4.',ch10,&
&   'Action: please switch to so_psp (dimensioned to npsp - with different meanings of 0 and 1.'
   MSG_ERROR(message)
 end if

 ABI_DEALLOCATE(intarr)
 ABI_DEALLOCATE(dprarr)

 call timab(191,2,tsec)

end subroutine invars2
!!***
