TABLE OF CONTENTS


ABINIT/indefo [ Functions ]

[ Top ] [ Functions ]

NAME

 indefo

FUNCTION

 Initialisation phase : defaults values for most input variables
 (some are initialized earlier, see indefo1 routine)

INPUTS

  ndtset_alloc=number of datasets, corrected for allocation of at
               least one data set.
  nprocs=Number of MPI processors available.

OUTPUT

  dtsets(0:ndtset_alloc)=<type datafiles_type>contains all input variables,
   some of which are given a default value here.
   The dataset with number 0 should be the reference default value
   in the remaining of the code.

NOTES

 The outputs of this routine are the defaults values of input
 variables, stored at the index 0 of the last dimension of their
 multi-dataset representation.

TODO

  Scalars and static arrays can be initialized directly at the level of the datatype declaration
  provided the value does not depend on runtime conditions.

PARENTS

      m_ab7_invars_f90

CHILDREN

SOURCE

2001 subroutine indefo(dtsets,ndtset_alloc,nprocs)
2002 
2003  use defs_basis
2004  use m_abicore
2005  use m_errors
2006  use m_gwdefs
2007 #if defined DEV_YP_VDWXC
2008  use m_xc_vdw
2009 #endif
2010 
2011  use defs_abitypes,  only : dataset_type
2012  use m_fftcore,      only : get_cache_kb, fftalg_for_npfft
2013 
2014 !This section has been created automatically by the script Abilint (TD).
2015 !Do not modify the following lines by hand.
2016 #undef ABI_FUNC
2017 #define ABI_FUNC 'indefo'
2018 !End of the abilint section
2019 
2020  implicit none
2021 
2022 !Arguments ------------------------------------
2023 !scalars
2024  integer,intent(in) :: ndtset_alloc,nprocs
2025 !arrays
2026  type(dataset_type),intent(inout) :: dtsets(0:ndtset_alloc) !vz_i
2027 
2028 !Local variables -------------------------------
2029 !scalars
2030  integer :: idtset,ii,jdtset,paral_atom_default
2031  logical :: wvl_bigdft
2032 #if defined DEV_YP_VDWXC
2033  type(xc_vdw_type) :: vdw_defaults
2034 #endif
2035 
2036 !******************************************************************
2037 
2038  DBG_ENTER("COLL")
2039 
2040 !Set up default values. All variables to be output in outvars.f
2041 !should have a default, even if a nonsensible one can be
2042 !chosen to garantee print in that routine.
2043 
2044 !These variables have already been initialized, for idtset/=0
2045  dtsets(0)%istatr=0
2046  dtsets(0)%istatshft=1
2047  dtsets(0)%kptrlatt(1:3,1:3)=0
2048  !dtsets(0)%kptrlatt_orig=0
2049  dtsets(0)%qptrlatt(1:3,1:3)=0
2050  dtsets(0)%ptgroupma=0
2051  dtsets(0)%spgroup=0
2052  dtsets(0)%shiftk(:,:)=half
2053  dtsets(0)%tolsym=tol8
2054  dtsets(0)%znucl(:)=zero
2055  dtsets(0)%ucrpa=0
2056  dtsets(0)%usedmft=0
2057 
2058  paral_atom_default=0
2059  if (nprocs>1.and.maxval(dtsets(:)%usepaw)>0) paral_atom_default=1
2060 
2061 !WARNING : set default in all datasets, including idtset=0 !!!
2062 !Use alphabetic order
2063 
2064  do idtset=0,ndtset_alloc
2065    jdtset=dtsets(idtset)%jdtset
2066 
2067    wvl_bigdft=.false.
2068    if(dtsets(idtset)%usewvl==1 .and. dtsets(idtset)%wvl_bigdft_comp==1) then
2069      wvl_bigdft=.true.
2070    end if
2071 !  Special case of use_gpu_cuda (can be undertermined at this point)
2072 !  use_gpu_cuda=-1 means undetermined ; here impose its value due to some restrictions
2073    if (dtsets(idtset)%use_gpu_cuda==-1) then
2074      if (dtsets(idtset)%optdriver/=0.or.&
2075 &     dtsets(idtset)%tfkinfunc/=0.or.&
2076 &     dtsets(idtset)%nspinor/=1) then
2077        dtsets(idtset)%use_gpu_cuda=0
2078      else
2079        dtsets(idtset)%use_gpu_cuda=1
2080      end if
2081    end if
2082 
2083 !  A
2084 !  Here we change the default value of iomode according to the configuration options.
2085 !  Ideally, all the sequential tests should pass independently of the default value.
2086 !  The parallel tests may require IO_MODE_MPI or, alternatively, IO_MODE_ETSF with HDF5 support.
2087 !  MG FIXME Sun Sep 6 2015: Many tests fail if IO_MODE_MPI is used as default. IO errors in v1, v2 ...
2088 !  with np=1 and wonderful deadlocks if np>1.
2089 
2090 !  Note that this default value might be overriden for specific datasets later, in case of parallelism
2091    dtsets(idtset)%iomode=IO_MODE_FORTRAN
2092 #ifdef HAVE_NETCDF_DEFAULT
2093    dtsets(idtset)%iomode=IO_MODE_ETSF
2094 #endif
2095 #ifdef HAVE_MPI_IO_DEFAULT
2096    dtsets(idtset)%iomode=IO_MODE_MPI
2097 #endif
2098 
2099    dtsets(idtset)%adpimd=0
2100    dtsets(idtset)%adpimd_gamma=one
2101    dtsets(idtset)%accuracy=0
2102    dtsets(idtset)%atvshift(:,:,:)=zero
2103    dtsets(idtset)%auxc_ixc=11
2104    dtsets(idtset)%auxc_scal=one
2105    dtsets(idtset)%awtr=1
2106 !  B
2107    dtsets(idtset)%bdberry(1:4)=0
2108    dtsets(idtset)%bdeigrf=-1
2109    dtsets(idtset)%bdgw=0
2110    dtsets(idtset)%berrystep=1
2111    dtsets(idtset)%bmass=ten
2112    dtsets(idtset)%boxcenter(1:3)=half
2113    dtsets(idtset)%boxcutmin=two
2114    dtsets(idtset)%brvltt=0
2115    dtsets(idtset)%bs_nstates=0
2116 !  dtsets(idtset)%bs_hayd_term=0
2117    dtsets(idtset)%bs_hayd_term=1
2118    dtsets(idtset)%builtintest=0
2119    dtsets(idtset)%bxctmindg=two
2120 !  C
2121    dtsets(idtset)%cd_halfway_freq=3.674930883_dp !(100 eV)
2122    dtsets(idtset)%cd_imfrqs(:) = zero
2123    dtsets(idtset)%cd_max_freq=36.74930883_dp     !(1000 eV)
2124    dtsets(idtset)%cd_subset_freq(1:2)=0
2125    dtsets(idtset)%cd_frqim_method=1
2126    dtsets(idtset)%cd_full_grid=0
2127    dtsets(idtset)%charge=zero
2128    dtsets(idtset)%chempot(:,:,:)=zero
2129    dtsets(idtset)%chkdilatmx=1
2130    dtsets(idtset)%chkexit=0
2131    dtsets(idtset)%chksymbreak=1
2132    dtsets(idtset)%cineb_start=7
2133    dtsets(idtset)%corecs(:) = zero
2134 !  D
2135    dtsets(idtset)%ddamp=0.1_dp
2136    dtsets(idtset)%delayperm=0
2137    dtsets(idtset)%densfor_pred=2
2138    if (dtsets(idtset)%paral_kgb>0.and.idtset>0) dtsets(idtset)%densfor_pred=6 ! Recommended for band-FFT parallelism
2139 !XG170502 : This section is completely useless, as ionmov is NOT know at present !
2140 !#ifdef HAVE_LOTF
2141 !   if (dtsets(idtset)%ionmov==23) dtsets(idtset)%densfor_pred=2 ! Recommended for LOTF
2142 !#endif
2143    dtsets(idtset)%dfpt_sciss=zero
2144    dtsets(idtset)%diecut=2.2_dp
2145    dtsets(idtset)%dielng=1.0774841_dp
2146    dtsets(idtset)%diemac=1.0d6
2147    if (dtsets(idtset)%usepaw==0) then
2148      dtsets(idtset)%diemix=one
2149    else
2150      dtsets(idtset)%diemix=0.7_dp
2151    end if
2152    dtsets(idtset)%diemixmag=dtsets(idtset)%diemix
2153    dtsets(idtset)%diegap=0.1_dp
2154    dtsets(idtset)%dielam=half
2155    dtsets(idtset)%diismemory=8
2156    dtsets(idtset)%dilatmx=one
2157    dtsets(idtset)%dmatpuopt=2
2158    if (size(dtsets(idtset)%dmatpawu,4)>0) dtsets(idtset)%dmatpawu=-10._dp
2159    dtsets(idtset)%dmatudiag=0
2160    dtsets(idtset)%dmft_entropy=0
2161    dtsets(idtset)%dmft_dc  =1
2162    dtsets(idtset)%dmft_iter=0
2163    dtsets(idtset)%dmft_nlambda=6
2164    dtsets(idtset)%dmft_nwli=0
2165    dtsets(idtset)%dmft_nwlo=0
2166    dtsets(idtset)%dmft_mxsf=0.3_dp
2167    dtsets(idtset)%dmft_read_occnd=0
2168    dtsets(idtset)%dmft_rslf=0
2169    dtsets(idtset)%dmft_solv=5
2170    if(dtsets(idtset)%ucrpa>0.and.dtsets(idtset)%usedmft==1) dtsets(idtset)%dmft_solv=0
2171    dtsets(idtset)%dmft_t2g=0
2172    dtsets(idtset)%dmft_tolfreq=tol4
2173    dtsets(idtset)%dmft_tollc=tol5
2174    dtsets(idtset)%dmftbandi=0
2175    dtsets(idtset)%dmftbandf=0
2176    dtsets(idtset)%dmftcheck=0
2177    dtsets(idtset)%dmftctqmc_basis =1
2178    dtsets(idtset)%dmftctqmc_check =0
2179    dtsets(idtset)%dmftctqmc_correl=0
2180    dtsets(idtset)%dmftctqmc_grnns =0
2181    dtsets(idtset)%dmftctqmc_meas  =1
2182    dtsets(idtset)%dmftctqmc_mrka  =0
2183    dtsets(idtset)%dmftctqmc_mov   =0
2184    dtsets(idtset)%dmftctqmc_order =0
2185    dtsets(idtset)%dmftctqmc_triqs_nleg=30
2186    dtsets(idtset)%dmftqmc_l=0
2187    dtsets(idtset)%dmftqmc_n=0.0_dp
2188    dtsets(idtset)%dmftqmc_seed=jdtset
2189    dtsets(idtset)%dmftqmc_therm=1000
2190    dtsets(idtset)%dmftctqmc_gmove = dtsets(idtset)%dmftqmc_therm / 10
2191    dtsets(idtset)%dosdeltae=zero
2192    dtsets(idtset)%dtion=100.0_dp
2193    dtsets(idtset)%d3e_pert1_atpol(1:2)=1
2194    dtsets(idtset)%d3e_pert1_dir(1:3)=0
2195    dtsets(idtset)%d3e_pert1_elfd=0
2196    dtsets(idtset)%d3e_pert1_phon=0
2197    dtsets(idtset)%d3e_pert2_atpol(1:2)=1
2198    dtsets(idtset)%d3e_pert2_dir(1:3)=0
2199    dtsets(idtset)%d3e_pert2_elfd=0
2200    dtsets(idtset)%d3e_pert2_phon=0
2201    dtsets(idtset)%d3e_pert3_atpol(1:2)=1
2202    dtsets(idtset)%d3e_pert3_dir(1:3)=0
2203    dtsets(idtset)%d3e_pert3_elfd=0
2204    dtsets(idtset)%d3e_pert3_phon=0
2205 !  E
2206    dtsets(idtset)%ecut=-one
2207    dtsets(idtset)%ecuteps=zero
2208    dtsets(idtset)%ecutsigx=zero ! The true default value is ecut . This is defined in invars2.F90
2209    dtsets(idtset)%ecutsm=zero
2210    dtsets(idtset)%ecutwfn=zero ! The true default value is ecut . This is defined in invars2.F90
2211    dtsets(idtset)%effmass_free=one
2212    dtsets(idtset)%efmas=0
2213    dtsets(idtset)%efmas_bands=0 ! The true default is nband. This is defined in invars2.F90
2214    dtsets(idtset)%efmas_deg=1
2215    dtsets(idtset)%efmas_deg_tol=tol5
2216    dtsets(idtset)%efmas_dim=3
2217    dtsets(idtset)%efmas_dirs=zero
2218    dtsets(idtset)%efmas_ntheta=1000
2219    dtsets(idtset)%elph2_imagden=zero
2220    dtsets(idtset)%enunit=0
2221    dtsets(idtset)%eshift=zero
2222    dtsets(idtset)%esmear=0.01_dp
2223    dtsets(idtset)%exchn2n3d=0
2224    dtsets(idtset)%extrapwf=0
2225    dtsets(idtset)%exchmix=quarter
2226 !  F
2227    dtsets(idtset)%fermie_nest=zero
2228    dtsets(idtset)%fftgw=21
2229    dtsets(idtset)%focktoldfe=zero
2230    dtsets(idtset)%fockoptmix=0
2231    dtsets(idtset)%fockdownsampling(:)=1
2232    dtsets(idtset)%freqim_alpha=five
2233    dtsets(idtset)%freqremin=zero
2234    dtsets(idtset)%freqremax=zero
2235    dtsets(idtset)%freqspmin=zero
2236    dtsets(idtset)%freqspmax=zero
2237    dtsets(idtset)%friction=0.001_dp
2238    dtsets(idtset)%frzfermi=0
2239    dtsets(idtset)%fxcartfactor=one ! Should be adjusted to the H2 conversion factor
2240 !  G
2241    dtsets(idtset)%ga_algor =1
2242    dtsets(idtset)%ga_fitness =1
2243    dtsets(idtset)%ga_opt_percent =0.2_dp
2244    dtsets(idtset)%ga_rules(:) =1
2245    dtsets(idtset)%goprecon =0
2246    dtsets(idtset)%goprecprm(:)=0
2247    dtsets(idtset)%gpu_devices=(/-1,-1,-1,-1,-1/)
2248    dtsets(idtset)%gpu_linalg_limit=2000000
2249    if (dtsets(idtset)%gw_customnfreqsp/=0) dtsets(idtset)%gw_freqsp(:) = zero
2250    dtsets(idtset)%gw_nstep =30
2251    dtsets(idtset)%gwgamma =0
2252    if ( dtsets(idtset)%gw_nqlwl > 0 ) then
2253      dtsets(idtset)%gw_qlwl(:,:)=zero
2254      dtsets(idtset)%gw_qlwl(1,1)=0.00001_dp
2255      dtsets(idtset)%gw_qlwl(2,1)=0.00002_dp
2256      dtsets(idtset)%gw_qlwl(3,1)=0.00003_dp
2257    end if
2258    dtsets(idtset)%gw_frqim_inzgrid=0
2259    dtsets(idtset)%gw_frqre_inzgrid=0
2260    dtsets(idtset)%gw_frqre_tangrid=0
2261    dtsets(idtset)%gw_invalid_freq=0
2262    dtsets(idtset)%gw_qprange=0
2263    dtsets(idtset)%gw_sigxcore=0
2264    dtsets(idtset)%gw_sctype = GWSC_one_shot
2265    dtsets(idtset)%gw_toldfeig=0.1/Ha_eV
2266    dtsets(idtset)%getbseig=0
2267    dtsets(idtset)%getbsreso=0
2268    dtsets(idtset)%getbscoup=0
2269    dtsets(idtset)%getcell =0
2270    dtsets(idtset)%getddb  =0
2271    dtsets(idtset)%getddk  =0
2272    dtsets(idtset)%getdelfd=0
2273    dtsets(idtset)%getdkdk =0
2274    dtsets(idtset)%getdkde =0
2275    dtsets(idtset)%getden  =0
2276    dtsets(idtset)%getefmas=0
2277    dtsets(idtset)%getgam_eig2nkq  =0
2278    dtsets(idtset)%gethaydock=0
2279    dtsets(idtset)%getocc  =0
2280    dtsets(idtset)%getpawden=0
2281    dtsets(idtset)%getqps  =0
2282    dtsets(idtset)%getscr  =0
2283    dtsets(idtset)%getsuscep=0
2284    dtsets(idtset)%getvel  =0
2285    dtsets(idtset)%getwfk  =0
2286    dtsets(idtset)%getwfkfine = 0
2287    dtsets(idtset)%getwfq  =0
2288    dtsets(idtset)%getxcart=0
2289    dtsets(idtset)%getxred =0
2290    dtsets(idtset)%get1den =0
2291    dtsets(idtset)%get1wf  =0
2292    dtsets(idtset)%gwcalctyp=0
2293    dtsets(idtset)%gwcomp=0
2294    dtsets(idtset)%gwencomp=2.0_dp
2295    dtsets(idtset)%gwmem=11
2296    dtsets(idtset)%gwpara=2
2297    dtsets(idtset)%gwrpacorr=0
2298    dtsets(idtset)%gwls_stern_kmax=1
2299    dtsets(idtset)%gwls_model_parameter=1.0_dp
2300    dtsets(idtset)%gwls_npt_gauss_quad=10
2301    dtsets(idtset)%gwls_diel_model=2
2302    dtsets(idtset)%gwls_print_debug=0
2303    if (dtsets(idtset)%gwls_n_proj_freq/=0) dtsets(idtset)%gwls_list_proj_freq(:) = zero
2304    dtsets(idtset)%gwls_nseeds=1
2305    dtsets(idtset)%gwls_recycle=2
2306    dtsets(idtset)%gwls_kmax_complement=1
2307    dtsets(idtset)%gwls_kmax_poles=4
2308    dtsets(idtset)%gwls_kmax_analytic=8
2309    dtsets(idtset)%gwls_kmax_numeric=16
2310    dtsets(idtset)%gwls_band_index=1
2311    dtsets(idtset)%gwls_exchange=1
2312    dtsets(idtset)%gwls_correlation=3
2313    dtsets(idtset)%gwls_first_seed=0
2314 !  H
2315    dtsets(idtset)%hmcsst=3
2316    dtsets(idtset)%hmctt=4
2317    dtsets(idtset)%hyb_mixing=-999.0_dp
2318    dtsets(idtset)%hyb_mixing_sr=-999.0_dp
2319    dtsets(idtset)%hyb_range_dft=-999.0_dp
2320    dtsets(idtset)%hyb_range_fock=-999.0_dp
2321 !  I
2322    if(dtsets(idtset)%natsph/=0) then
2323 !    do not use iatsph(:) but explicit boundaries
2324 !    to avoid to read to far away in the built array (/ ... /)
2325      dtsets(idtset)%iatsph(1:dtsets(idtset)%natsph)=(/ (ii,ii=1,dtsets(idtset)%natsph) /)
2326    else
2327      dtsets(idtset)%iatsph(:)=0
2328    end if
2329    dtsets(idtset)%iboxcut=0
2330    dtsets(idtset)%icutcoul=6
2331    dtsets(idtset)%ieig2rf=0
2332    dtsets(idtset)%imgwfstor=0
2333    dtsets(idtset)%inclvkb=2
2334    dtsets(idtset)%intxc=0
2335 !  if (dtsets(idtset)%paral_kgb>0.and.idtset>0) dtsets(idtset)%intxc=0
2336    dtsets(idtset)%ionmov=0
2337    dtsets(idtset)%densfor_pred=2
2338    if (dtsets(idtset)%paral_kgb>0.and.idtset>0) dtsets(idtset)%densfor_pred=6 ! Recommended for band-FFT parallelism
2339 !This section is completely useless, as ionmov is NOT know at present !
2340 !#ifdef HAVE_LOTF
2341 !   if (dtsets(idtset)%ionmov==23) dtsets(idtset)%densfor_pred=2 ! Recommended for LOTF
2342 !#endif
2343    dtsets(idtset)%iprcel=0
2344    dtsets(idtset)%iprcfc=0
2345    dtsets(idtset)%irandom=3
2346    dtsets(idtset)%irdbseig=0
2347    dtsets(idtset)%irdbsreso=0
2348    dtsets(idtset)%irdbscoup=0
2349    dtsets(idtset)%irdddb=0
2350    dtsets(idtset)%irdddk=0
2351    dtsets(idtset)%irdden=0
2352    dtsets(idtset)%irdefmas=0
2353    dtsets(idtset)%irdhaydock=0
2354    dtsets(idtset)%irdpawden=0
2355    dtsets(idtset)%irdqps=0
2356    dtsets(idtset)%irdscr=0
2357    dtsets(idtset)%irdsuscep=0
2358    dtsets(idtset)%irdvdw=0
2359    dtsets(idtset)%irdwfk=0
2360    dtsets(idtset)%irdwfkfine=0
2361    dtsets(idtset)%irdwfq=0
2362    dtsets(idtset)%ird1den=0
2363    dtsets(idtset)%ird1wf=0
2364 !iscf
2365    if(wvl_bigdft) then
2366      dtsets(idtset)%iscf=0
2367    else
2368      if(dtsets(idtset)%usepaw==0) then
2369        dtsets(idtset)%iscf=7
2370      else
2371        dtsets(idtset)%iscf=17
2372      end if
2373    end if
2374    dtsets(idtset)%isecur=0
2375    dtsets(idtset)%istatimg = 1
2376    dtsets(idtset)%istwfk(:)=0
2377    dtsets(idtset)%ixc=1
2378    dtsets(idtset)%ixc_sigma=1
2379    dtsets(idtset)%ixcpositron=1
2380    dtsets(idtset)%ixcrot=3
2381 !  J
2382    dtsets(idtset)%f4of2_sla(:)=-one
2383    dtsets(idtset)%f6of2_sla(:)=-one
2384    dtsets(idtset)%jpawu(:,:)=zero
2385 !  K
2386    dtsets(idtset)%kberry(1:3,:)=0
2387    dtsets(idtset)%kpt(:,:)=zero
2388    dtsets(idtset)%kptgw(:,:)=zero
2389    dtsets(idtset)%kptnrm=one
2390    dtsets(idtset)%kptns_hf(:,:)=zero
2391    dtsets(idtset)%kptopt=1
2392    if(dtsets(idtset)%nspden==4)dtsets(idtset)%kptopt=4
2393    dtsets(idtset)%kptrlen=30.0_dp
2394    dtsets(idtset)%kssform=1
2395 !  L
2396    dtsets(idtset)%localrdwf=1
2397 
2398 #if defined HAVE_LOTF
2399    dtsets(idtset)%lotf_classic=5
2400    dtsets(idtset)%lotf_nitex=10
2401    dtsets(idtset)%lotf_nneigx=40
2402    dtsets(idtset)%lotf_version=2
2403 #endif
2404 !  M
2405    dtsets(idtset)%magconon = 0
2406    dtsets(idtset)%magcon_lambda = 0.01_dp
2407    dtsets(idtset)%max_ncpus = 0
2408    dtsets(idtset)%mbpt_sciss=zero
2409    dtsets(idtset)%mband = -1
2410    dtsets(idtset)%mdf_epsinf = zero
2411    dtsets(idtset)%mdtemp(:)=300.0_dp
2412    dtsets(idtset)%mdwall=10000_dp
2413    dtsets(idtset)%mep_mxstep=100._dp
2414    dtsets(idtset)%mep_solver=0
2415    dtsets(idtset)%mffmem=1
2416    dtsets(idtset)%mgfft = -1
2417    dtsets(idtset)%mgfftdg = -1
2418    dtsets(idtset)%mixesimgf(:)=zero
2419    dtsets(idtset)%mpw = -1
2420    dtsets(idtset)%mqgrid=0
2421    dtsets(idtset)%mqgriddg=0
2422 !  N
2423    dtsets(idtset)%natrd = -1
2424    dtsets(idtset)%nband(:)=0
2425    dtsets(idtset)%nbandhf=0
2426    dtsets(idtset)%nbdblock=1
2427    dtsets(idtset)%nbdbuf=0
2428    dtsets(idtset)%nberry=1
2429    if (dtsets(idtset)%usepaw==0) then
2430      dtsets(idtset)%nc_xccc_gspace=0
2431    else
2432      dtsets(idtset)%nc_xccc_gspace=1
2433    end if
2434    dtsets(idtset)%nbandkss=0
2435    dtsets(idtset)%nctime=0
2436    dtsets(idtset)%ndtset = -1
2437    dtsets(idtset)%neb_algo=1
2438    dtsets(idtset)%neb_spring(1:2)=(/0.05_dp,0.05_dp/)
2439    dtsets(idtset)%npwkss=0
2440    dtsets(idtset)%nfft = -1
2441    dtsets(idtset)%nfftdg = -1
2442 
2443    dtsets(idtset)%nfreqim=-1
2444    dtsets(idtset)%nfreqre=-1
2445    dtsets(idtset)%nfreqsp=0
2446 
2447    dtsets(idtset)%npulayit=7
2448 
2449 !  ngfft is a special case
2450    dtsets(idtset)%ngfft(1:8)=0
2451    dtsets(idtset)%ngfft(7) = fftalg_for_npfft(1)
2452 !  fftcache=ngfft(8) is machine-dependent.
2453    dtsets(idtset)%ngfft(8) = get_cache_kb()
2454 
2455    dtsets(idtset)%ngfftdg(:)=dtsets(idtset)%ngfft(:)
2456 !
2457    !nline
2458    dtsets(idtset)%nline=4
2459    if(dtsets(idtset)%usewvl==1 .and. .not. wvl_bigdft) then
2460      if(dtsets(idtset)%usepaw==1) then
2461        dtsets(idtset)%nline=4
2462      else
2463        dtsets(idtset)%nline=2
2464      end if
2465    end if
2466 
2467 !  nloalg is also a special case
2468    dtsets(idtset)%nloalg(1)=4
2469    dtsets(idtset)%nloalg(2)=1
2470    dtsets(idtset)%nloalg(3)=dtsets(idtset)%usepaw
2471    dtsets(idtset)%ngkpt=0
2472    dtsets(idtset)%nnsclo=0
2473    dtsets(idtset)%nnsclohf=0
2474    dtsets(idtset)%nomegasf=100
2475    dtsets(idtset)%nomegasrd=9
2476    dtsets(idtset)%nomegasi=12
2477    dtsets(idtset)%nonlinear_info=0
2478    dtsets(idtset)%noseinert=1.0d5
2479    dtsets(idtset)%npvel=0
2480    dtsets(idtset)%npweps=0
2481    dtsets(idtset)%npwsigx=0
2482    dtsets(idtset)%npwwfn=0
2483    dtsets(idtset)%nqpt=0
2484    dtsets(idtset)%nscforder=16
2485    dtsets(idtset)%nshiftk=1
2486    dtsets(idtset)%nshiftk_orig=1
2487    dtsets(idtset)%nstep=30
2488    dtsets(idtset)%ntime=1
2489    dtsets(idtset)%nwfshist=0
2490    if(dtsets(idtset)%usewvl==1 .and. .not. wvl_bigdft) then
2491      if(dtsets(idtset)%usepaw==1) then
2492        dtsets(idtset)%nwfshist=4
2493      else
2494        dtsets(idtset)%nwfshist=2
2495      end if
2496    end if
2497 !  O
2498    dtsets(idtset)%occopt=1
2499    dtsets(idtset)%occ_orig(:,:)=zero
2500    dtsets(idtset)%omegasrdmax=1.0_dp/Ha_eV  ! = 1eV
2501    dtsets(idtset)%omegasimax=50/Ha_eV
2502    dtsets(idtset)%optcell=0
2503    dtsets(idtset)%optforces=2
2504    if(dtsets(idtset)%usedmft>0) dtsets(idtset)%optforces=0
2505    dtsets(idtset)%optstress=1
2506    dtsets(idtset)%optnlxccc=1
2507    dtsets(idtset)%orbmag=0
2508    if (dtsets(idtset)%usepaw==0) then
2509      dtsets(idtset)%ortalg=2
2510 !    dtsets(idtset)%ortalg=999
2511    else
2512      dtsets(idtset)%ortalg=-2
2513 !    dtsets(idtset)%ortalg=999
2514    end if
2515 !  P
2516    dtsets(idtset)%paral_atom=paral_atom_default
2517    dtsets(idtset)%pawcpxocc=1
2518    dtsets(idtset)%pawcross=0
2519    dtsets(idtset)%pawecutdg=-one
2520    dtsets(idtset)%pawfatbnd=0
2521    dtsets(idtset)%pawlcutd=10
2522    dtsets(idtset)%pawlmix=10
2523    dtsets(idtset)%pawmixdg=0 ! Will be set to 1 when npfft>1
2524    dtsets(idtset)%pawnhatxc=1
2525    dtsets(idtset)%pawntheta=12
2526    dtsets(idtset)%pawnphi=13
2527    dtsets(idtset)%pawnzlm=1
2528    dtsets(idtset)%pawoptmix=0
2529    dtsets(idtset)%pawoptosc=0
2530    dtsets(idtset)%pawovlp=5._dp
2531    dtsets(idtset)%pawprtdos=0
2532    dtsets(idtset)%pawprtvol=0
2533    dtsets(idtset)%pawprtwf=0
2534    dtsets(idtset)%pawprt_k=0
2535    dtsets(idtset)%pawprt_b=0
2536    dtsets(idtset)%pawstgylm=1
2537    dtsets(idtset)%pawsushat=0
2538    dtsets(idtset)%pawujat=1
2539    dtsets(idtset)%pawujrad=20.0_dp
2540    dtsets(idtset)%pawujv=0.1_dp/Ha_eV
2541    dtsets(idtset)%pawusecp=1
2542    dtsets(idtset)%pawxcdev=1
2543    dtsets(idtset)%ph_nqshift = 0
2544    if(dtsets(idtset)%ph_nqshift > 0)then
2545      dtsets(idtset)%ph_qshift = zero
2546    end if
2547    dtsets(idtset)%pimd_constraint=0
2548    dtsets(idtset)%pitransform=0
2549    dtsets(idtset)%ptcharge(:) = zero
2550    dtsets(idtset)%plowan_bandi=0
2551    dtsets(idtset)%plowan_bandf=0
2552    if(dtsets(idtset)%plowan_compute>0) then
2553      dtsets(idtset)%plowan_it(:)=0
2554      dtsets(idtset)%plowan_iatom(:)=0
2555      dtsets(idtset)%plowan_lcalc(:)=-1
2556      dtsets(idtset)%plowan_projcalc(:)=0
2557      dtsets(idtset)%plowan_nbl(:)=0
2558    end if
2559    dtsets(idtset)%plowan_natom=0
2560    dtsets(idtset)%plowan_nt=0
2561    dtsets(idtset)%plowan_realspace=0
2562    dtsets(idtset)%pol(:)=zero
2563    dtsets(idtset)%polcen(:)=zero
2564    dtsets(idtset)%posdoppler=0
2565    dtsets(idtset)%positron=0
2566    dtsets(idtset)%posnstep=50
2567    dtsets(idtset)%posocc=one
2568    dtsets(idtset)%postoldfe=0.000001_dp
2569    dtsets(idtset)%postoldff=zero
2570    dtsets(idtset)%ppmodel=1
2571    dtsets(idtset)%ppmfrq=zero
2572    dtsets(idtset)%prepanl=0
2573    dtsets(idtset)%prepgkk=0
2574    dtsets(idtset)%prtbbb=0
2575    dtsets(idtset)%prtbltztrp=0
2576    dtsets(idtset)%prtcif=0
2577    dtsets(idtset)%prtden=1;if (dtsets(idtset)%nimage>1) dtsets(idtset)%prtden=0
2578    dtsets(idtset)%prtdensph=1
2579    dtsets(idtset)%prtdipole=0
2580    dtsets(idtset)%prtdos=0
2581    dtsets(idtset)%prtdosm=0
2582    dtsets(idtset)%prtebands=1;if (dtsets(idtset)%nimage>1) dtsets(idtset)%prtebands=0
2583    dtsets(idtset)%prtefg=0
2584    dtsets(idtset)%prtefmas=0
2585    dtsets(idtset)%prteig=1;if (dtsets(idtset)%nimage>1) dtsets(idtset)%prteig=0
2586    dtsets(idtset)%prtelf=0
2587    dtsets(idtset)%prtfc=0
2588    dtsets(idtset)%prtfull1wf=0
2589    dtsets(idtset)%prtfsurf=0
2590    dtsets(idtset)%prtgden=0
2591    dtsets(idtset)%prtgeo=0
2592    dtsets(idtset)%prtgkk=0
2593    dtsets(idtset)%prtkden=0
2594    dtsets(idtset)%prtkpt = -1
2595    dtsets(idtset)%prtlden=0
2596    dtsets(idtset)%prtnabla=0
2597    dtsets(idtset)%prtnest=0
2598    dtsets(idtset)%prtphdos=1
2599    dtsets(idtset)%prtphsurf=0
2600    dtsets(idtset)%prtposcar=0
2601    dtsets(idtset)%prtpot=0
2602    dtsets(idtset)%prtpsps=0
2603    dtsets(idtset)%prtspcur=0
2604    dtsets(idtset)%prtsuscep=0
2605    dtsets(idtset)%prtstm=0
2606    dtsets(idtset)%prtvclmb=0
2607    dtsets(idtset)%prtvdw=0
2608    dtsets(idtset)%prtvha=0
2609    dtsets(idtset)%prtvhxc=0
2610    dtsets(idtset)%prtvxc=0
2611    dtsets(idtset)%prtvol=0
2612    dtsets(idtset)%prtvolimg=0
2613    dtsets(idtset)%prtvpsp=0
2614    dtsets(idtset)%prtwant=0
2615    dtsets(idtset)%prtwf=1; if (dtsets(idtset)%nimage>1) dtsets(idtset)%prtwf=0
2616    dtsets(idtset)%prtwf_full=0
2617    dtsets(idtset)%prtxml = 0
2618    do ii=1,dtsets(idtset)%natom,1
2619      dtsets(idtset)%prtatlist(ii)=ii
2620    end do
2621    dtsets(idtset)%prt1dm=0
2622    dtsets(idtset)%pvelmax(:)=one
2623    dtsets(idtset)%pw_unbal_thresh=40.
2624 !  Q
2625    dtsets(idtset)%qmass(:)=ten
2626    dtsets(idtset)%qprtrb(1:3)=0
2627    dtsets(idtset)%qptdm(:,:)=zero
2628    dtsets(idtset)%quadmom(:) = zero
2629 !  R
2630    dtsets(idtset)%random_atpos=0
2631    dtsets(idtset)%ratsph_extra=two
2632    dtsets(idtset)%recefermi=zero
2633    dtsets(idtset)%recgratio=1
2634    dtsets(idtset)%recnpath=500
2635    dtsets(idtset)%recnrec=10
2636    dtsets(idtset)%recrcut=zero
2637    dtsets(idtset)%recptrott=0
2638    dtsets(idtset)%rectesteg=0
2639    dtsets(idtset)%rectolden=zero
2640    dtsets(idtset)%rcut=zero
2641    dtsets(idtset)%restartxf=0
2642    dtsets(idtset)%rfasr=0
2643    dtsets(idtset)%rfatpol(1:2)=1
2644    dtsets(idtset)%rfddk=0
2645    dtsets(idtset)%rfdir(1:3)=0
2646    dtsets(idtset)%rfelfd=0
2647    dtsets(idtset)%rfmagn=0
2648    dtsets(idtset)%rfmeth=1
2649    dtsets(idtset)%rfphon=0
2650    dtsets(idtset)%rfstrs=0
2651    dtsets(idtset)%rfuser=0
2652    dtsets(idtset)%rf2_dkdk=0
2653    dtsets(idtset)%rf2_dkde=0
2654    dtsets(idtset)%rf2_pert1_dir(1:3)=0
2655    dtsets(idtset)%rf2_pert2_dir(1:3)=0
2656    dtsets(idtset)%rhoqpmix=one
2657 !  S
2658    dtsets(idtset)%shiftk_orig(:,:)=one
2659    dtsets(idtset)%signperm=1
2660    dtsets(idtset)%slabwsrad=zero
2661    dtsets(idtset)%slk_rankpp=1000
2662    dtsets(idtset)%smdelta=0
2663    dtsets(idtset)%spbroad=0.1
2664    dtsets(idtset)%spgaxor = -1
2665    dtsets(idtset)%spgorig = -1
2666    dtsets(idtset)%spinmagntarget=-99.99_dp
2667    dtsets(idtset)%spmeth=0
2668    dtsets(idtset)%spnorbscl=one
2669    dtsets(idtset)%stmbias=zero
2670    dtsets(idtset)%strfact=100.0_dp
2671    dtsets(idtset)%string_algo=1
2672    dtsets(idtset)%strprecon=one
2673    dtsets(idtset)%strtarget(1:6)=zero
2674    dtsets(idtset)%symchi=1
2675    dtsets(idtset)%symsigma=0
2676 !  T
2677    dtsets(idtset)%td_maxene=zero
2678    dtsets(idtset)%td_mexcit=0
2679    dtsets(idtset)%tfw_toldfe=0.000001_dp
2680    dtsets(idtset)%tim1rev = 1
2681    dtsets(idtset)%tl_nprccg = 30
2682    dtsets(idtset)%tl_radius = zero
2683    dtsets(idtset)%tphysel=zero
2684    dtsets(idtset)%toldfe=zero
2685    dtsets(idtset)%tolmxde=zero
2686    dtsets(idtset)%toldff=zero
2687    dtsets(idtset)%tolimg=5.0d-5
2688    dtsets(idtset)%tolrde=0.005_dp
2689    dtsets(idtset)%tolrff=zero
2690    dtsets(idtset)%tolmxf=5.0d-5
2691    dtsets(idtset)%tolvrs=zero
2692    dtsets(idtset)%tolwfr=zero
2693    dtsets(idtset)%tmesh=[5._dp, 59._dp, 6._dp]
2694    dtsets(idtset)%tsmear=0.01_dp
2695 !  U
2696    dtsets(idtset)%ucrpa_bands(:)=-1
2697    dtsets(idtset)%ucrpa_window(:)=-1.0_dp
2698    dtsets(idtset)%upawu(:,:)=zero
2699    dtsets(idtset)%usepead=1
2700    dtsets(idtset)%usefock=0
2701    dtsets(idtset)%usekden=0
2702    dtsets(idtset)%use_gemm_nonlop=0
2703    dtsets(idtset)%use_nonscf_gkk=0 !1 ! deactivate by default, for now 6 Oct 2013
2704    dtsets(idtset)%userec=0
2705    dtsets(idtset)%usexcnhat_orig=-1
2706    dtsets(idtset)%useylm=0
2707 !  V
2708    dtsets(idtset)%vacnum = -1
2709    dtsets(idtset)%vcutgeo(:)=zero
2710    dtsets(idtset)%vdw_nfrag = 1
2711 #if defined DEV_YP_VDWXC
2712    dtsets(idtset)%vdw_df_acutmin = vdw_defaults%acutmin
2713    dtsets(idtset)%vdw_df_aratio = vdw_defaults%aratio
2714    dtsets(idtset)%vdw_df_damax = vdw_defaults%damax
2715    dtsets(idtset)%vdw_df_damin = vdw_defaults%damin
2716    dtsets(idtset)%vdw_df_dcut = vdw_defaults%dcut
2717    dtsets(idtset)%vdw_df_dratio = vdw_defaults%dratio
2718    dtsets(idtset)%vdw_df_dsoft = vdw_defaults%dsoft
2719    dtsets(idtset)%vdw_df_gcut = vdw_defaults%gcut
2720    dtsets(idtset)%vdw_df_ndpts = vdw_defaults%ndpts
2721    dtsets(idtset)%vdw_df_ngpts = vdw_defaults%ngpts
2722    dtsets(idtset)%vdw_df_nqpts = vdw_defaults%nqpts
2723    dtsets(idtset)%vdw_df_nrpts = vdw_defaults%nrpts
2724    dtsets(idtset)%vdw_df_nsmooth = vdw_defaults%nsmooth
2725    dtsets(idtset)%vdw_df_phisoft = vdw_defaults%phisoft
2726    dtsets(idtset)%vdw_df_qcut = vdw_defaults%qcut
2727    dtsets(idtset)%vdw_df_qratio = vdw_defaults%qratio
2728    dtsets(idtset)%vdw_df_rcut = vdw_defaults%rcut
2729    dtsets(idtset)%vdw_df_rsoft = vdw_defaults%rsoft
2730    dtsets(idtset)%vdw_df_tolerance = vdw_defaults%tolerance
2731    dtsets(idtset)%vdw_df_tweaks = vdw_defaults%tweaks
2732    dtsets(idtset)%vdw_df_zab = vdw_defaults%zab
2733    dtsets(idtset)%vdw_df_threshold = 1.0d-2
2734 #endif
2735    dtsets(idtset)%vdw_supercell(:) = 0
2736    dtsets(idtset)%vdw_tol = tol10
2737    dtsets(idtset)%vdw_tol_3bt = -1
2738    dtsets(idtset)%vdw_typfrag(:) = 1
2739    dtsets(idtset)%vdw_xc = 0
2740    dtsets(idtset)%vis=100.0_dp
2741    dtsets(idtset)%vprtrb(1:2)=zero
2742 !  W
2743    dtsets(idtset)%wtatcon(:,:,:)=zero
2744    dtsets(idtset)%wfmix=one
2745    dtsets(idtset)%wfk_task=0
2746    dtsets(idtset)%wtk=one
2747    dtsets(idtset)%wvl_crmult  = 6._dp
2748    dtsets(idtset)%wvl_frmult  = 10._dp
2749    dtsets(idtset)%wvl_hgrid   = 0.5_dp
2750    dtsets(idtset)%wvl_ngauss  =(1,100)
2751    dtsets(idtset)%wvl_nprccg  = 10
2752    dtsets(idtset)%w90iniprj   = 1
2753    dtsets(idtset)%w90prtunk   = 0
2754 
2755 !  X
2756    dtsets(idtset)%xclevel  = 0
2757    dtsets(idtset)%xc_denpos = tol14
2758    dtsets(idtset)%xc_tb09_c = 99.99_dp
2759    dtsets(idtset)%xredsph_extra(:,:)=zero
2760 !  Y
2761 !  Z
2762    dtsets(idtset)%zcut=3.67493260d-03  ! = 0.1eV
2763    if(dtsets(idtset)%optdriver == RUNL_GWLS) dtsets(idtset)%zcut=zero
2764    dtsets(idtset)%ziontypat(:)=zero
2765 
2766 !  BEGIN VARIABLES FOR @Bethe-Salpeter
2767    dtsets(idtset)%bs_algorithm    =2
2768    dtsets(idtset)%bs_haydock_niter=100
2769    dtsets(idtset)%bs_exchange_term=1
2770    dtsets(idtset)%bs_coulomb_term=11
2771    dtsets(idtset)%bs_calctype=1
2772    dtsets(idtset)%bs_coupling=0
2773 
2774    dtsets(idtset)%bs_haydock_tol=(0.02_dp,zero)
2775 
2776    dtsets(idtset)%bs_loband=0
2777 !  Take big absolute value numbers, but the the biggest ones, otherwise overflow can happen
2778    dtsets(idtset)%bs_eh_cutoff = [smallest_real*tol6,greatest_real*tol6]
2779    dtsets(idtset)%bs_freq_mesh = [zero,zero,0.01_dp/Ha_eV]
2780 
2781 !  Interpolation
2782    dtsets(idtset)%bs_interp_method=1 ! YG interpolation
2783    dtsets(idtset)%bs_interp_mode=0 ! No interpolation
2784    dtsets(idtset)%bs_interp_prep=0 ! Do not prepare interp
2785    dtsets(idtset)%bs_interp_kmult=(/zero,zero,zero/)
2786    dtsets(idtset)%bs_interp_m3_width=one
2787    dtsets(idtset)%bs_interp_rl_nb=1
2788 
2789 !  END VARIABLES FOR @Bethe-Salpeter.
2790 
2791 ! EPH variables
2792    dtsets(idtset)%asr = 1
2793    dtsets(idtset)%dipdip = 1
2794    dtsets(idtset)%chneut = 0
2795    dtsets(idtset)%symdynmat = 1
2796 
2797    dtsets(idtset)%ph_freez_disp_addStrain = 0
2798    dtsets(idtset)%ph_freez_disp_option = 0
2799    dtsets(idtset)%ph_freez_disp_nampl = 0
2800    if(dtsets(idtset)%ph_freez_disp_nampl>0)dtsets(idtset)%ph_freez_disp_ampl = zero
2801    dtsets(idtset)%ph_ndivsm = 20
2802    dtsets(idtset)%ph_nqpath = 0
2803    dtsets(idtset)%ph_ngqpt = [20, 20, 20]
2804 
2805    dtsets(idtset)%eph_mustar = 0.1_dp
2806    dtsets(idtset)%eph_intmeth = 2
2807    dtsets(idtset)%eph_extrael = zero
2808    dtsets(idtset)%eph_fermie = zero
2809    dtsets(idtset)%eph_frohlichm = 0
2810    dtsets(idtset)%eph_fsmear = 0.01
2811    dtsets(idtset)%eph_fsewin = 0.04
2812    dtsets(idtset)%eph_ngqpt_fine = [0, 0, 0]
2813    dtsets(idtset)%eph_task = 1
2814    dtsets(idtset)%eph_transport  = 0
2815 
2816    dtsets(idtset)%ph_wstep = 0.1/Ha_meV
2817    dtsets(idtset)%ph_intmeth = 2
2818    dtsets(idtset)%ph_nqshift = 1
2819    dtsets(idtset)%ph_smear = 0.00002_dp
2820    dtsets(idtset)%ddb_ngqpt = [0, 0, 0]
2821    dtsets(idtset)%ddb_shiftq(:) = zero
2822 
2823 ! JB:UNINITIALIZED VALUES (not found in this file neither indefo1)
2824 ! They might be initialized somewhereelse, I don't know.
2825 ! That might cause unitialized error with valgrind depending on the compilo
2826 ! chkprim
2827 ! maxnsym
2828 ! nsym
2829 ! macro_uj
2830 ! prtpmp
2831 ! timopt
2832 ! useria
2833 ! userib
2834 ! useric
2835 ! userid
2836 ! userie
2837 ! bravais
2838 ! symafm
2839 ! symrel
2840 ! fband
2841 ! nelect
2842 ! userra
2843 ! userrb
2844 ! userrc
2845 ! userrd
2846 ! userre
2847 ! vacwidth
2848 ! genafm
2849 ! kptns
2850 ! rprimd_orig
2851 ! tnons
2852 
2853  end do
2854 
2855  DBG_EXIT("COLL")
2856 
2857 end subroutine indefo

ABINIT/indefo1 [ Functions ]

[ Top ] [ Functions ]

NAME

 indefo1

FUNCTION

 Initialisation phase : defaults values for a first batch of input variables
 (especially dimensions, needed to allocate other parts of dtsets, as well
  as other input variables whose existence is needed for other initialisations to proceed).

INPUTS

OUTPUT

  dtset=<type datafiles_type>contains all input variables for one dataset,
   some of which are given a default value here.

PARENTS

      invars1m

CHILDREN

SOURCE

 841 subroutine indefo1(dtset)
 842 
 843 
 844 !This section has been created automatically by the script Abilint (TD).
 845 !Do not modify the following lines by hand.
 846 #undef ABI_FUNC
 847 #define ABI_FUNC 'indefo1'
 848 !End of the abilint section
 849 
 850  implicit none
 851 
 852 !Arguments ------------------------------------
 853 !scalars
 854 !arrays
 855  type(dataset_type),intent(inout) :: dtset
 856 
 857 !Local variables -------------------------------
 858 !scalars
 859  integer :: ii
 860 
 861 !******************************************************************
 862 !
 863 !Set up default values. All variables to be output in outvars.f
 864 !should have a default, even if a nonsensible one can be
 865 !chosen to garantee print in that routine.
 866 
 867  DBG_ENTER("COLL")
 868 
 869 !Use alphabetic order
 870 
 871 !A
 872  dtset%acell_orig(:,:)=zero
 873  dtset%algalch(:)=1
 874  dtset%amu_orig(:,:)=-one
 875  dtset%autoparal=0
 876 !B
 877  dtset%bandpp=1
 878  dtset%berryopt=0
 879  dtset%berrysav=0
 880  dtset%bfield(:)=zero
 881 !C
 882  dtset%cd_customnimfrqs=0
 883  dtset%chkprim=1
 884 !D
 885  dtset%densty(:,:)=zero
 886  dtset%dfield(:)=zero    !!HONG
 887  dtset%dynimage(:)=1
 888 !E
 889  dtset%efield(:)=zero
 890  dtset%efmas_calc_dirs=0
 891  dtset%efmas_n_dirs=0
 892 !F
 893 !G
 894  dtset%ga_n_rules=1
 895  dtset%gw_customnfreqsp=0
 896  dtset%gw_nqlwl=0
 897  dtset%gwls_n_proj_freq=0
 898 !I
 899  dtset%iatfix(:,:)=0
 900  dtset%icoulomb=0
 901  dtset%imgmov=0
 902 !J
 903  dtset%jellslab=0
 904  dtset%jfielddir(:)=0
 905 !K
 906  dtset%kptopt=0
 907 !L
 908  dtset%lexexch(:)=-1
 909  dtset%ldaminushalf(:)=0
 910  dtset%lpawu(:)=-1
 911 !M
 912  dtset%maxestep=0.005d0
 913  dtset%mixalch_orig(:,:,:)=zero
 914  dtset%mkmem=-1
 915  dtset%mkqmem=-1
 916  dtset%mk1mem=-1
 917 !N
 918  dtset%natpawu=0
 919  dtset%natsph=0
 920  dtset%natsph_extra=0
 921  dtset%natvshift=0
 922  dtset%nconeq=0
 923  dtset%ndynimage=1
 924  dtset%nkpt=-1
 925  dtset%nkptgw=0
 926  dtset%nkpthf=0
 927  dtset%nnos=0
 928  dtset%npband=1
 929  dtset%npfft=1
 930  dtset%nphf=1
 931  dtset%npimage=1
 932  dtset%npkpt=1
 933  dtset%nppert=1
 934  dtset%npspalch=0
 935  dtset%npspinor=1
 936  dtset%np_slk=1000000
 937  dtset%nqptdm=0
 938  dtset%nspden=1
 939  dtset%nspinor=1
 940  dtset%nsppol=1
 941  dtset%nsym=0     ! Actually, this default value is not used : it is to be reimposed before each call to ingeo in invars1
 942  dtset%ntimimage=1
 943  dtset%ntypalch=0
 944  dtset%ntyppure=-1
 945  dtset%nucdipmom(:,:)=zero
 946  dtset%nzchempot=0
 947 !O
 948  dtset%optdriver=0
 949 !P
 950  dtset%paral_rf=0
 951 !dtset%paral_kgb ! Is even initialized earlier.
 952  dtset%pawspnorb=0  ! will be changed to 1 as soon as usepaw==1 and nspinor==2
 953  dtset%pimass(:)=-one
 954 !Q
 955  dtset%qptn=zero
 956 !R
 957  dtset%red_efield(:)=zero
 958  dtset%red_dfield(:)=zero
 959  dtset%red_efieldbar(:)=zero
 960  dtset%rprim_orig(:,:,:)=zero
 961  dtset%rprim_orig(1,1,:)=one
 962  dtset%rprim_orig(2,2,:)=one
 963  dtset%rprim_orig(3,3,:)=one
 964 !S
 965  dtset%slabzbeg=zero
 966  dtset%slabzend=zero
 967  dtset%so_psp(:)=1
 968  dtset%spinat(:,:)=zero
 969  dtset%supercell_latt(:,:) = 0
 970  do ii=1,3
 971    dtset%supercell_latt(ii,ii) = 1
 972  end do
 973  dtset%symmorphi=1
 974 !T
 975  dtset%tfkinfunc=0
 976  dtset%typat(:)=0  ! This init is important because dimension of typat is mxnatom (and not natom).
 977 !U
 978  dtset%usedmatpu=0
 979  dtset%usedmft=0
 980  dtset%useexexch=0
 981  dtset%usepawu=0
 982  dtset%usepotzero=0
 983  dtset%use_slk=0
 984 !V
 985  dtset%vel_orig(:,:,:)=zero
 986  dtset%vel_cell_orig(:,:,:)=zero
 987 !W
 988  dtset%wtq=0
 989  if (dtset%usepaw==0) dtset%wfoptalg=0
 990  if (dtset%usepaw/=0) dtset%wfoptalg=10
 991  if (dtset%optdriver==RUNL_GSTATE.and.dtset%paral_kgb>0) dtset%wfoptalg=14
 992  dtset%wvl_bigdft_comp=1
 993 
 994 !X
 995  dtset%xred_orig(:,:,:)=zero
 996 !Y
 997 !Z
 998  dtset%zeemanfield(:)=zero
 999 
1000  DBG_EXIT("COLL")
1001 
1002 end subroutine indefo1

ABINIT/invars0 [ Functions ]

[ Top ] [ Functions ]

NAME

 invars0

FUNCTION

 Initialisation phase : prepare the main input subroutine call by
 reading most of the NO MULTI variables, as well as natom, nimage, and ntypat,
 needed for allocating some input arrays in abinit, and also useri
 and userr. The variable usewvl is also read here for later reading
 of input path for the atomic orbital file (if required).

INPUTS

  lenstr=actual length of string
  ndtset= number of datasets to be read; if 0, no multi-dataset mode
  ndtset_alloc=number of datasets, corrected for allocation of at least
               one data set.
  string*(*)=string of characters containing all input variables and data

OUTPUT

  dtsets(0:ndtset_alloc)=<type datafiles_type>contains all input variables,
   some of which are initialized here :
   cpus,jdtset,natom,nimage,npsp,ntypat,useri*,userr*
  istatr=repetition rate for status file
  istatshft=shift of the repetition rate for status file
  msym=maximal value of input msym for all the datasets
  mxnatom=maximal value of input natom for all the datasets
  mxnimage=maximal value of input nimage for all the datasets
  mxntypat=maximal value of input ntypat for all the datasets
  npsp=number of pseudopotentials

PARENTS

      m_ab7_invars_f90

CHILDREN

      get_ndevice,intagm

SOURCE

 99 subroutine invars0(dtsets,istatr,istatshft,lenstr,&
100 & msym,mxnatom,mxnimage,mxntypat,ndtset,ndtset_alloc,npsp,papiopt,timopt,string)
101 
102 
103 !This section has been created automatically by the script Abilint (TD).
104 !Do not modify the following lines by hand.
105 #undef ABI_FUNC
106 #define ABI_FUNC 'invars0'
107 !End of the abilint section
108 
109  implicit none
110 
111 !Arguments ------------------------------------
112 !scalars
113  integer,intent(in) :: lenstr,ndtset,ndtset_alloc
114  integer,intent(out) :: istatr,istatshft,msym,mxnatom,mxnimage,mxntypat,npsp,papiopt
115  integer,intent(inout) :: timopt
116  character(len=*),intent(in) :: string
117 !arrays
118  type(dataset_type),intent(inout) :: dtsets(0:ndtset_alloc) !vz_i
119 
120 !Local variables-------------------------------
121 !scalars
122  integer :: i1,i2,idtset,ii,jdtset,marr,multiplicity,tjdtset,tread,treadh,treadm
123  integer :: treads,use_gpu_cuda
124  real(dp) :: cpus
125  character(len=500) :: message
126 !arrays
127  integer,allocatable :: intarr(:)
128  real(dp),allocatable :: dprarr(:)
129 
130 !******************************************************************
131 
132 !Set ii to avoid warning of uninitialised variable
133  ii = 0
134 
135  marr=max(9,ndtset_alloc,2)
136  ABI_ALLOCATE(dprarr,(marr))
137  ABI_ALLOCATE(intarr,(marr))
138 
139 !Set up jdtset
140  if(ndtset/=0)then
141 
142 !  Default values
143    dtsets(0)%jdtset = -1 ! unused value
144    dtsets(1:ndtset_alloc)%jdtset=(/ (ii,ii=1,ndtset_alloc) /)
145 
146 !  Read explicitly the jdtset array
147    call intagm(dprarr,intarr,0,marr,ndtset,string(1:lenstr),'jdtset',tjdtset,'INT')
148    if(tjdtset==1) dtsets(1:ndtset)%jdtset=intarr(1:ndtset)
149 
150 !  Read the udtset array
151    call intagm(dprarr,intarr,0,marr,2,string(1:lenstr),'udtset',tread,'INT')
152 
153 !  jdtset and udtset cannot be defined together
154    if(tjdtset==1 .and. tread==1)then
155      write(message, '(3a)' )&
156 &     'jdtset and udtset cannot be defined both in the input file.',ch10,&
157 &     'Action: remove one of them from your input file.'
158      MSG_ERROR(message)
159    end if
160 
161 !  Check values of udtset
162    if(tread==1)then
163      if(intarr(1)<1 .or. intarr(1)>999)then
164        write(message, '(a,i0,3a)' )&
165 &       'udtset(1) must be between 1 and 999, but it is ',intarr(1),'.',ch10,&
166 &       'Action: change the value of udtset(1) in your input file.'
167        MSG_ERROR(message)
168      end if
169      if(intarr(2)<1 .or. intarr(2)>9)then
170        write(message, '(a,i0,3a)' )&
171 &       'udtset(2) must be between 1 and 9, but it is ',intarr(2),'.',ch10,&
172 &       'Action: change the value of udtset(2) in your input file.'
173        MSG_ERROR(message)
174      end if
175      if(intarr(1)*intarr(2) /= ndtset)then
176        write(message, '(3a,i0,3a,i0,a,i0,3a,i0,3a)' )&
177 &       'udtset(1)*udtset(2) must be equal to ndtset,',ch10,&
178 &       'but it is observed that udtset(1) = ',intarr(1),',',ch10,&
179 &       'and udtset(2) = ',intarr(2),' so that their product is ',intarr(1)*intarr(2),',',ch10,&
180 &       'while ndtset is ',ndtset,'.',ch10,&
181 &       'Action: change udtset or ndtset in your input file.'
182        MSG_ERROR(message)
183      end if
184      idtset=0
185      do i1=1,intarr(1)
186        do i2=1,intarr(2)
187          idtset=idtset+1
188          dtsets(idtset)%jdtset=i1*10+i2
189        end do
190      end do
191    end if
192 
193 !  Final check on the jdtset values
194    do idtset=1,ndtset
195      if(dtsets(idtset)%jdtset<1 .or. dtsets(idtset)%jdtset>9999)then
196        write(message, '(3a,i0,a,i0,a,a)' )&
197 &       'The components of jdtset must be between 1 and 9999.',ch10,&
198 &       'However, the input value of the component ',idtset,' of jdtset is ',dtsets(idtset)%jdtset,ch10,&
199 &       'Action : correct jdtset in your input file.'
200        MSG_ERROR(message)
201      end if
202    end do
203 
204  else
205    dtsets(1)%jdtset=0
206  end if
207 
208  papiopt = 0
209  call intagm(dprarr,intarr,0,1,1,string(1:lenstr),'papiopt',tread,'INT')
210  if(tread==1) papiopt=intarr(1)
211 
212 !Read timopt and pass it to timab
213  call intagm(dprarr,intarr,0,1,1,string(1:lenstr),'timopt',tread,'INT')
214  if(tread==1) timopt=intarr(1)
215 
216  istatr=0
217  dtsets(0)%istatr=istatr
218  call intagm(dprarr,intarr,0,marr,1,string(1:lenstr),'istatr',tread,'INT')
219  if(tread==1) istatr=intarr(1)
220  dtsets(1:)%istatr=istatr
221 
222  istatshft=1
223  dtsets(0)%istatshft=istatshft
224  call intagm(dprarr,intarr,0,marr,1,string(1:lenstr),'istatshft',tread,'INT')
225  if(tread==1) istatshft=intarr(1)
226  dtsets(1:)%istatshft=istatshft
227 
228  cpus=zero
229  call intagm(dprarr,intarr,0,marr,1,string(1:lenstr),'cpus ',treads,'DPR')
230  if(treads==1) cpus=dprarr(1)
231  call intagm(dprarr,intarr,0,marr,1,string(1:lenstr),'cpum ',treadm,'DPR')
232  if(treadm==1) cpus=dprarr(1)*60.0_dp
233  call intagm(dprarr,intarr,0,marr,1,string(1:lenstr),'cpuh ',treadh,'DPR')
234 
235  if(treadh==1) cpus=dprarr(1)*3600.0_dp
236  if(treads+treadm+treadh>1)then
237    write(message, '(5a)' )&
238 &   'More than one input variable is used to defined the CPU time limit.',ch10,&
239 &   'This is not allowed.',ch10,&
240 &   'Action: in the input file, suppress either cpus, cpum or cpuh.'
241    MSG_ERROR(message)
242  end if
243  dtsets(:)%cpus=cpus
244 
245 !Default for natom, nimage, ntypat, useri and userr
246  dtsets(:)%natom=1
247  dtsets(:)%nimage=1
248  dtsets(:)%ntypat=1 ; dtsets(0)%ntypat=0    ! Will always echo ntypat
249  dtsets(:)%macro_uj=0
250  dtsets(:)%maxnsym=384
251  dtsets(:)%useria=0
252  dtsets(:)%userib=0
253  dtsets(:)%useric=0
254  dtsets(:)%userid=0
255  dtsets(:)%userie=0
256  dtsets(:)%userra=zero
257  dtsets(:)%userrb=zero
258  dtsets(:)%userrc=zero
259  dtsets(:)%userrd=zero
260  dtsets(:)%userre=zero
261  dtsets(:)%usewvl = 0
262  dtsets(:)%plowan_compute=0
263 
264 !Loop on datasets, to find natom and mxnatom, as well as useri and userr
265  do idtset=1,ndtset_alloc
266    jdtset=dtsets(idtset)%jdtset ; if(ndtset==0)jdtset=0
267 
268 ! proposal: supercell generation in input string before it is read in
269 ! call expand_supercell_input(jdtset, lenstr, string)
270 !  find supercell, else exit
271 !  determinant = ncells
272 !  copy rprim,    acell,    xred,    xcart,    xangst,    vel,    typat,   to
273 !       rprim_uc, acell_uc, xred_uc, xcart_uc, xangst_uc, vel_uc, typat_uc
274 !     NB: also rprim and angdeg need to be updated in non diagonal case!!!
275 !  generate supercell info for each of these copying out with translation vectors etc...
276 !  set chkprim to 0
277 !  done!
278 
279 !  Generate the supercell if supercell_latt is specified and update string
280    dtsets(idtset)%supercell_latt(:,:) = 0
281    do ii=1,3
282      dtsets(idtset)%supercell_latt(ii,ii) = 1
283    end do
284    call intagm(dprarr,intarr,jdtset,marr,9,string(1:lenstr),"supercell_latt",tread,'INT')
285    if (tread==1) dtsets(idtset)%supercell_latt(:,:)=reshape(intarr(1:marr),(/3,3/))
286    !This test should be update if in the future we allow non-diagonal supercell
287    if (any(dtsets(idtset)%supercell_latt(:,:) < zero).or.&
288 &   (dtsets(idtset)%supercell_latt(1,1) < tol10 .or.&
289 &   dtsets(idtset)%supercell_latt(2,2) <tol10  .or.&
290 &   dtsets(idtset)%supercell_latt(3,3) < tol10 )) then
291      write(message, '(5a)' )&
292 &     'supercell_latt must have positive parameters and diagonal part',ch10,&
293 &     'This is not allowed.  ',ch10,&
294 &     'Action : modify supercell_latt in the input file.'
295      MSG_ERROR(message)
296    end if
297 !  Compute the multiplicity of the supercell
298    call mati3det(dtsets(idtset)%supercell_latt,multiplicity)
299 
300 !  Read natom from string
301    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natom',tread,'INT')
302 !  Might initialize natom from XYZ file
303    if(tread==0)then
304      call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'_natom',tread,'INT')
305    end if
306 
307    if(tread==1)then
308      dtsets(idtset)%natom=intarr(1)
309    else
310      write(message, '(a,i0,2a)' )&
311 &     'Input natom must be defined, but was absent for dataset ',jdtset,ch10,&
312 &     'Action: check the input file.'
313      MSG_ERROR(message)
314    end if
315 !  Check that natom is greater than 0
316    if (dtsets(idtset)%natom<=0) then
317      write(message, '(a,i0,2a,i0,3a)' )&
318 &     'Input natom must be > 0, but was ',dtsets(idtset)%natom,ch10,&
319 &     'for dataset ',jdtset,'. This is not allowed.',ch10,&
320 &     'Action: check the input file.'
321      MSG_ERROR(message)
322    end if
323 
324    if(multiplicity > 1)then
325      dtsets(idtset)%natom = dtsets(idtset)%natom * multiplicity
326    end if
327 
328    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nimage',tread,'INT')
329    if(tread==1) dtsets(idtset)%nimage=intarr(1)
330 
331 !  Check that nimage is greater than 0
332    if (dtsets(idtset)%nimage<=0) then
333      write(message, '(a,i0,4a)' )&
334 &     'nimage must be > 0, but was ',dtsets(idtset)%nimage,ch10,&
335 &     'This is not allowed.',ch10,&
336 &     'Action: check the input file.'
337      MSG_ERROR(message)
338    end if
339 
340    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ntypat',tread,'INT')
341    if(tread==1)dtsets(idtset)%ntypat=intarr(1)
342 !  Check that ntypat is greater than 0
343    if (dtsets(idtset)%ntypat<=0) then
344      write(message, '(a,i0,2a,i0,3a)' )&
345 &     'Input ntypat must be > 0, but was ',dtsets(idtset)%ntypat,ch10,&
346 &     'for dataset ',jdtset,'. This is not allowed.',ch10,&
347 &     'Action: check the input file.'
348      MSG_ERROR(message)
349    end if
350 
351 !  Read msym from string
352    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'maxnsym',tread,'INT')
353    if(tread==1)dtsets(idtset)%maxnsym=intarr(1)
354 !  Check that maxnsym is greater than 1
355    if (dtsets(idtset)%maxnsym<1) then
356      write(message, '(a,i0,2a,i0,3a)' )&
357 &     'Input maxnsym must be > 1, but was ',dtsets(idtset)%maxnsym,ch10,&
358 &     'for dataset ',jdtset,'. This is not allowed.',ch10,&
359 &     'Action: check the input file.'
360      MSG_ERROR(message)
361    end if
362 
363 ! Read plowan_compute
364    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'plowan_compute',tread,'INT')
365    if(tread==1) dtsets(idtset)%plowan_compute=intarr(1)
366 
367    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'useria',tread,'INT')
368    if(tread==1) dtsets(idtset)%useria=intarr(1)
369    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userib',tread,'INT')
370    if(tread==1) dtsets(idtset)%userib=intarr(1)
371    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'useric',tread,'INT')
372    if(tread==1) dtsets(idtset)%useric=intarr(1)
373    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userid',tread,'INT')
374    if(tread==1) dtsets(idtset)%userid=intarr(1)
375    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userie',tread,'INT')
376    if(tread==1) dtsets(idtset)%userie=intarr(1)
377 
378    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userra',tread,'DPR')
379    if(tread==1) dtsets(idtset)%userra=dprarr(1)
380    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userrb',tread,'DPR')
381    if(tread==1) dtsets(idtset)%userrb=dprarr(1)
382    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userrc',tread,'DPR')
383    if(tread==1) dtsets(idtset)%userrc=dprarr(1)
384    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userrd',tread,'DPR')
385    if(tread==1) dtsets(idtset)%userrd=dprarr(1)
386    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userre',tread,'DPR')
387    if(tread==1) dtsets(idtset)%userre=dprarr(1)
388 
389    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usewvl',tread,'INT')
390    if(tread==1) dtsets(idtset)%usewvl=intarr(1)
391 
392  end do
393 
394 !mxnatom =maxval(dtsets(1:ndtset_alloc)%natom)
395 !mxntypat =maxval(dtsets(1:ndtset_alloc)%ntypat)
396 !msym =maxval(dtsets(1:ndtset_alloc)%maxnsym)
397 !There is a bug in the HP compiler, the following should execute properly
398  mxnatom=dtsets(1)%natom ; mxnimage=dtsets(1)%nimage
399  mxntypat=dtsets(1)%ntypat ; msym=dtsets(1)%maxnsym
400  if(ndtset_alloc>1)then
401    do idtset=2,ndtset_alloc
402      mxnatom =max(dtsets(idtset)%natom,mxnatom)
403      mxnimage=max(dtsets(idtset)%nimage,mxnimage)
404      mxntypat=max(dtsets(idtset)%ntypat,mxntypat)
405      msym    =max(dtsets(idtset)%maxnsym,msym)
406    end do
407  end if
408 
409  if(mxnimage>1)then
410    do idtset=2,ndtset_alloc
411      if(mxnatom/=dtsets(idtset)%natom)then
412        write(message,'(5a,i0,a,i0,3a,i0,a)')&
413 &       'When there exist one dataset with more than one image,',ch10,&
414 &       'the number of atoms in each dataset must be the same.',ch10,&
415 &       'However, it has been found that for dataset= ',idtset,ch10,&
416 &       'natom= ',dtsets(idtset)%natom,' differs from the maximum number',ch10,&
417 &       'of atoms, mxnatom= ',mxnatom,&
418 &       'Action: check the input variables natom for different datasets.'
419        MSG_ERROR(message)
420      end if
421    end do
422  end if
423 
424 !Set up npsp
425  npsp=mxntypat   ! Default value
426  call intagm(dprarr,intarr,0,marr,1,string(1:lenstr),'npsp',tread,'INT')
427  if(tread==1)then
428    npsp=intarr(1)
429  else
430    if(ndtset_alloc>1)then
431      do idtset=1,ndtset_alloc
432        if(dtsets(idtset)%ntypat/=mxntypat)then
433          write(message, '(5a,i0,a,i0,2a,i0,2a)' )&
434 &         'When npsp is not defined, the input variable ntypat must be',ch10,&
435 &         'the same for all datasets. However, it has been found that for',ch10,&
436 &         'jdtset: ',dtsets(idtset)%jdtset,', ntypat= ',dtsets(idtset)%ntypat,ch10,&
437 &         'differs from the maximum value of ntypat= ',mxntypat,ch10,&
438 &         'Action: check the input variables npsp and ntypat.'
439          MSG_ERROR(message)
440        end if
441      end do
442    end if
443  end if
444  dtsets(0)%npsp = mxntypat   ! Default value
445  dtsets(1:ndtset_alloc)%npsp = npsp
446 
447 !KGB parallelism information (needed at this stage)
448  dtsets(:)%paral_kgb=0
449  do idtset=1,ndtset_alloc
450    jdtset=dtsets(idtset)%jdtset ; if(ndtset==0)jdtset=0
451    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'paral_kgb',tread,'INT')
452    if(tread==1)dtsets(idtset)%paral_kgb=intarr(1)
453 
454    if (dtsets(idtset)%paral_kgb<0 .or. dtsets(idtset)%paral_kgb>1) then
455      write(message,'(a,i0,2a,i0,3a)')&
456 &     'Input paral_kgb must be 0 or 1, but was ',dtsets(idtset)%paral_kgb,ch10,&
457 &     'for dataset',jdtset,'. This is not allowed.',ch10,&
458 &     'Action: check the input file.'
459      MSG_ERROR(message)
460    end if
461  end do
462 
463 !GPU information
464  use_gpu_cuda=0
465  dtsets(:)%use_gpu_cuda=0
466 #if defined HAVE_GPU_CUDA && defined HAVE_GPU_CUDA_DP
467  call Get_ndevice(ii)
468  if (ii>0) then
469    do i1=1,ndtset_alloc
470      dtsets(i1)%use_gpu_cuda=-1
471    end do
472  end if
473 #endif
474  do idtset=1,ndtset_alloc
475    jdtset=dtsets(idtset)%jdtset ; if(ndtset==0)jdtset=0
476    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'use_gpu_cuda',tread,'INT')
477    if(tread==1)dtsets(idtset)%use_gpu_cuda=intarr(1)
478    if (dtsets(idtset)%use_gpu_cuda==1) use_gpu_cuda=1
479  end do
480  if (use_gpu_cuda==1) then
481 #if defined HAVE_GPU_CUDA && defined HAVE_GPU_CUDA_DP
482    if (ii<=0) then
483      write(message,'(3a)')&
484 &     'Input variables use_gpu_cuda is on',ch10,&
485 &     'but no available GPU device has been detected !',ch10,&
486 &     'Action: change the input variable use_gpu_cuda.'
487      MSG_ERROR(message)
488    end if
489 #else
490    write(message,'(7a)')&
491 &   'Input variables use_gpu_cuda is on but abinit hasn''t been built',ch10,&
492 &   'with (double precision) gpu mode enabled !',ch10,&
493 &   'Action: change the input variable use_gpu_cuda',ch10,&
494 &   '        or re-compile ABINIT with double-precision Cuda enabled.'
495    MSG_ERROR(message)
496 #endif
497  end if
498 
499  ABI_DEALLOCATE(dprarr)
500  ABI_DEALLOCATE(intarr)
501 
502 !We allocate the internal array, depending on the computed values.
503 !WARNING : do not forget to deallocate these arrays in the routine dtset_free
504 !(should make a separate subroutine for allocating/deallocating these records)
505  do idtset=0,ndtset_alloc
506    ABI_ALLOCATE(dtsets(idtset)%acell_orig,(3,mxnimage))
507    ABI_ALLOCATE(dtsets(idtset)%algalch,(mxntypat))
508    ABI_ALLOCATE(dtsets(idtset)%amu_orig,(mxntypat,mxnimage))
509    ABI_ALLOCATE(dtsets(idtset)%corecs,(mxntypat))
510    ABI_ALLOCATE(dtsets(idtset)%densty,(mxntypat,4))
511    ABI_ALLOCATE(dtsets(idtset)%dynimage,(mxnimage))
512    ABI_ALLOCATE(dtsets(idtset)%iatfix,(3,mxnatom))
513    ABI_ALLOCATE(dtsets(idtset)%f4of2_sla,(mxntypat))
514    ABI_ALLOCATE(dtsets(idtset)%f6of2_sla,(mxntypat))
515    ABI_ALLOCATE(dtsets(idtset)%jpawu,(mxntypat,mxnimage))
516    ABI_ALLOCATE(dtsets(idtset)%kberry,(3,20))
517    ABI_ALLOCATE(dtsets(idtset)%lexexch,(mxntypat))
518    ABI_ALLOCATE(dtsets(idtset)%ldaminushalf,(mxntypat))
519    ABI_ALLOCATE(dtsets(idtset)%lpawu,(mxntypat))
520    ABI_ALLOCATE(dtsets(idtset)%mixalch_orig,(npsp,mxntypat,mxnimage))
521    ABI_ALLOCATE(dtsets(idtset)%mixesimgf,(mxnimage))
522    ABI_ALLOCATE(dtsets(idtset)%nucdipmom,(3,mxnatom))
523    ABI_ALLOCATE(dtsets(idtset)%pimass,(mxntypat))
524    ABI_ALLOCATE(dtsets(idtset)%ptcharge,(mxntypat))
525    ABI_ALLOCATE(dtsets(idtset)%prtatlist,(mxnatom))
526    ABI_ALLOCATE(dtsets(idtset)%quadmom,(mxntypat))
527    ABI_ALLOCATE(dtsets(idtset)%ratsph,(mxntypat))
528    ABI_ALLOCATE(dtsets(idtset)%rprim_orig,(3,3,mxnimage))
529    ABI_ALLOCATE(dtsets(idtset)%rprimd_orig,(3,3,mxnimage))
530    ABI_ALLOCATE(dtsets(idtset)%so_psp,(npsp))
531    ABI_ALLOCATE(dtsets(idtset)%spinat,(3,mxnatom))
532    ABI_ALLOCATE(dtsets(idtset)%shiftk,(3,210))
533    ABI_ALLOCATE(dtsets(idtset)%typat,(mxnatom))
534    ABI_ALLOCATE(dtsets(idtset)%upawu,(mxntypat,mxnimage))
535 !   if (dtsets(idtset)%plowan_compute>0) then
536    ABI_ALLOCATE(dtsets(idtset)%plowan_iatom,(mxnatom))
537    ABI_ALLOCATE(dtsets(idtset)%plowan_it,(100*3))
538    ABI_ALLOCATE(dtsets(idtset)%plowan_nbl,(mxnatom))
539    ABI_ALLOCATE(dtsets(idtset)%plowan_lcalc,(12*mxnatom))
540    ABI_ALLOCATE(dtsets(idtset)%plowan_projcalc,(12*mxnatom))
541 !   endif
542    ABI_ALLOCATE(dtsets(idtset)%vel_orig,(3,mxnatom,mxnimage))
543    ABI_ALLOCATE(dtsets(idtset)%vel_cell_orig,(3,3,mxnimage))
544    ABI_ALLOCATE(dtsets(idtset)%xred_orig,(3,mxnatom,mxnimage))
545    ABI_ALLOCATE(dtsets(idtset)%ziontypat,(mxntypat))
546    ABI_ALLOCATE(dtsets(idtset)%znucl,(npsp))
547  end do
548 
549 !DEBUG
550 !write(std_out,*)' invars0 : nimage, mxnimage = ',dtsets(:)%nimage, mxnimage
551 !write(std_out,*)' invars0 : natom = ',dtsets(:)%natom
552 !write(std_out,*)' invars0 : mxnatom = ',mxnatom
553 !ENDDEBUG
554 
555 end subroutine invars0

ABINIT/invars1 [ Functions ]

[ Top ] [ Functions ]

NAME

 invars1

FUNCTION

 Initialize the dimensions needed to allocate the input arrays
 for one dataset characterized by jdtset, by
 taking from string the necessary data.
 Perform some preliminary checks and echo these dimensions.

INPUTS

  iout=unit number of output file
  jdtset=number of the dataset looked for
  lenstr=actual length of string
  msym=default maximal number of symmetries
  npsp1= number of pseudopotential files
  zionpsp(npsp1)= valence charge over all psps

OUTPUT

  mband_upper=estimation of the maximum number of bands for any k-point

SIDE EFFECTS

 Input/Output (the default value is given in the calling routine)
  dtset=<type datafiles_type>contains all input variables,
   some of which are initialized here, while other were already
   initialized, while some others will still be initialized later.
   The list of records of dtset initialized in the present routine is:
   acell_orig,densty,iatfix,kptopt,kptrlatt,
   mkmem,mkqmem,mk1mem,natsph,natvshift,nconeq,nkpt,nkptgw,nkpthf,
   nqptdm,nshiftk,nucdipmom,nzchempot,optdriver,
   rprim_orig,rprimd_orig,shiftk,
   spgroup,spinat,typat,vel_orig,vel_cell_orig,xred_orig
  bravais(11)=characteristics of Bravais lattice (see symlatt.F90)
  symafm(1:msym)=(anti)ferromagnetic part of symmetry operations
  symrel(3,3,1:msym)=symmetry operations in real space in terms
   of primitive translations
  tnons(3,1:msym)=nonsymmorphic translations for symmetry operations
  string*(*)=string of characters containing all input variables and data

NOTES

 Must set up the geometry of the system, needed to compute
 k point grids in an automatic fashion.
 Treat separately mband_upper, since
 fband, charge and zionpsp must be known for being able to initialize it.

 Defaults are provided in the calling routine.
 Defaults are also provided here for the following variables :
 mband_upper, occopt, fband, charge
 They should be kept consistent with defaults of the same variables
 provided to the invars routines.

PARENTS

      invars1m

CHILDREN

      atomdata_from_znucl,chkint_ge,ingeo,inkpts,inqpt,intagm,inupper
      invacuum,mkrdim,wrtout

SOURCE

1065 subroutine invars1(bravais,dtset,iout,jdtset,lenstr,mband_upper,msym,npsp1,&
1066 & string,symafm,symrel,tnons,zionpsp)
1067 
1068 
1069 !This section has been created automatically by the script Abilint (TD).
1070 !Do not modify the following lines by hand.
1071 #undef ABI_FUNC
1072 #define ABI_FUNC 'invars1'
1073 !End of the abilint section
1074 
1075  implicit none
1076 
1077 !Arguments ------------------------------------
1078 !scalars
1079  integer,intent(in) :: iout,jdtset,lenstr,msym,npsp1
1080  integer,intent(out) :: mband_upper
1081  character(len=*),intent(inout) :: string
1082  type(dataset_type),intent(inout) :: dtset
1083 !arrays
1084  integer,intent(inout) :: bravais(11),symafm(msym),symrel(3,3,msym)
1085  real(dp),intent(inout) :: tnons(3,msym)
1086  real(dp),intent(in) :: zionpsp(npsp1)
1087 
1088 !Local variables-------------------------------
1089  character :: blank=' ',string1
1090 !scalars
1091  integer :: chksymbreak,found,ierr,iatom,ii,ikpt,iimage,index_blank,index_lower
1092  integer :: index_typsymb,index_upper,ipsp,iscf,intimage,itypat,leave,marr
1093  integer :: natom,nkpt,nkpthf,npsp,npspalch
1094  integer :: nqpt,nspinor,nsppol,ntypat,ntypalch,ntyppure,occopt,response
1095  integer :: rfddk,rfelfd,rfphon,rfstrs,rfuser,rf2_dkdk,rf2_dkde,rfmagn
1096  integer :: tfband,tnband,tread,tread_alt
1097  real(dp) :: charge,fband,kptnrm,kptrlen,zelect,zval
1098  character(len=2) :: string2,symbol
1099  character(len=500) :: message
1100  type(atomdata_t) :: atom
1101 !arrays
1102  integer :: cond_values(4),vacuum(3)
1103  integer,allocatable :: iatfix(:,:),intarr(:),istwfk(:),nband(:),typat(:)
1104  real(dp) :: acell(3),rprim(3,3)
1105 !real(dp) :: field(3)
1106  real(dp),allocatable :: amu(:),dprarr(:),kpt(:,:),kpthf(:,:),mixalch(:,:),nucdipmom(:,:)
1107  real(dp),allocatable :: ratsph(:),reaalloc(:),spinat(:,:)
1108  real(dp),allocatable :: vel(:,:),vel_cell(:,:),wtk(:),xred(:,:),znucl(:)
1109  character(len=32) :: cond_string(4)
1110 
1111 
1112 !************************************************************************
1113 
1114 !Some initialisations
1115  ierr=0
1116  cond_string(1:4)=' '
1117  cond_values(1:4)=(/0,0,0,0/)
1118 
1119 !Read parameters
1120  marr=dtset%npsp;if (dtset%npsp<3) marr=3
1121  marr=max(marr,dtset%nimage)
1122  ABI_ALLOCATE(intarr,(marr))
1123  ABI_ALLOCATE(dprarr,(marr))
1124 
1125 !---------------------------------------------------------------------------
1126 
1127  rfddk=0; rfelfd=0; rfphon=0; rfmagn=0; rfstrs=0; rfuser=0; rf2_dkdk=0; rf2_dkde=0
1128  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfddk',tread,'INT')
1129  if(tread==1) rfddk=intarr(1)
1130  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfelfd',tread,'INT')
1131  if(tread==1) rfelfd=intarr(1)
1132  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfmagn',tread,'INT')
1133  if(tread==1) rfmagn=intarr(1)
1134  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfphon',tread,'INT')
1135  if(tread==1) rfphon=intarr(1)
1136  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfstrs',tread,'INT')
1137  if(tread==1) rfstrs=intarr(1)
1138  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfuser',tread,'INT')
1139  if(tread==1) rfuser=intarr(1)
1140  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rf2_dkdk',tread,'INT')
1141  if(tread==1) rf2_dkdk=intarr(1)
1142  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rf2_dkde',tread,'INT')
1143  if(tread==1) rf2_dkde=intarr(1)
1144 
1145  response=0
1146  if(rfddk/=0.or.rf2_dkdk/=0.or.rf2_dkde/=0.or.rfelfd/=0.or.rfphon/=0.or.rfstrs/=0.or.rfuser/=0.or.rfmagn/=0)response=1
1147 
1148  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'optdriver',tread,'INT')
1149  if (tread==1) then
1150    dtset%optdriver=intarr(1)
1151  else
1152 !  If optdriver was not read, while response=1, set optdriver to 1
1153    if(response==1)dtset%optdriver=1
1154  end if
1155 
1156 !---------------------------------------------------------------------------
1157 !For now, waiting express parallelisation for recursion
1158  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tfkinfunc',tread,'INT')
1159  if(tread==1) dtset%tfkinfunc=intarr(1)
1160 
1161 !---------------------------------------------------------------------------
1162 ! wvl_bigdft_comp, done here since default values of nline, nwfshist and iscf
1163 ! depend on its value (see indefo)
1164  if(dtset%usewvl==1) then
1165    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'wvl_bigdft_comp',tread,'INT')
1166    if(tread==1) dtset%wvl_bigdft_comp=intarr(1)
1167  end if
1168 
1169 !---------------------------------------------------------------------------
1170 
1171  natom=dtset%natom
1172  npsp=dtset%npsp
1173  ntypat=dtset%ntypat
1174 
1175 !No default value for znucl
1176  call intagm(dprarr,intarr,jdtset,marr,dtset%npsp,string(1:lenstr),'znucl',tread,'DPR')
1177  if(tread==1)then
1178    dtset%znucl(1:dtset%npsp)=dprarr(1:dtset%npsp)
1179  end if
1180  if(tread/=1)then
1181    write(message, '(3a)' )&
1182 &   'The array znucl MUST be initialized in the input file while this is not done.',ch10,&
1183 &   'Action: initialize znucl in your input file.'
1184    MSG_ERROR(message)
1185  end if
1186 
1187 !The default for ratsph has already been initialized
1188  call intagm(dprarr,intarr,jdtset,marr,dtset%ntypat,string(1:lenstr),'ratsph',tread,'LEN')
1189  if(tread==1)then
1190    do ii=1,dtset%ntypat
1191      dtset%ratsph(ii)=dprarr(ii)
1192    end do
1193  end if
1194  ABI_ALLOCATE(ratsph,(dtset%ntypat))
1195  do ii=1,dtset%ntypat
1196    ratsph(ii)=dtset%ratsph(ii)
1197  end do
1198 
1199 !Special treatment of _TYPAX (from a XYZ file), taking into account
1200 !the fact that znucl does NOT depend on the dataset
1201 !Examine all occurences of '_TYPAX'
1202 
1203  do
1204    index_typsymb=index(string(1:lenstr),'_TYPAX')
1205    if(index_typsymb==0)exit
1206 !  Replace '_TYPAX' by '_TYPAT'
1207    string(index_typsymb:index_typsymb+5)='_TYPAT'
1208    index_upper=index_typsymb+5
1209 !  Must start from the first blank after the tag (including possible dtset_char)
1210    index_upper=index(string(index_upper:lenstr),blank)+index_upper-1
1211    index_lower=index_upper
1212 
1213 !  Examine all atoms (the end of the symbol string is delimited by a XX )
1214    do
1215      index_blank=index(string(index_upper:lenstr),blank)+index_upper-1
1216      string2=string(index_blank+1:index_blank+2)
1217      if(string2=="XX")exit
1218      found=0
1219 !    Find the matching symbol
1220      do ipsp=1,dtset%npsp
1221        call atomdata_from_znucl(atom,dtset%znucl(ipsp))
1222        symbol = atom%symbol
1223        call inupper(symbol)
1224        call inupper(string2)
1225 
1226 !      DEBUG
1227 !      write(std_out,'(a)')' invars1 : before test, trim(adjustl(symbol)),trim(adjustl(string2))'
1228 !      write(std_out,'(5a)' )'"',trim(adjustl(symbol)),'","',trim(adjustl(string2)),'"'
1229 !      ENDDEBUG
1230 
1231        if(trim(adjustl(symbol))==trim(adjustl(string2)))then
1232          found=1
1233          index_upper=index_blank+1
1234 !        Cannot deal properly with more that 9 psps
1235          if(ipsp>=10)then
1236            message = 'Need to use a pseudopotential with number larger than 9. Not allowed yet.'
1237            MSG_ERROR(message)
1238          end if
1239 
1240 !        DEBUG
1241 !        write(std_out,*)' invars1 : found ipsp=',ipsp
1242 !        ENDDEBUG
1243 
1244          write(string1,'(i1)')ipsp
1245          string(index_lower:index_lower+1)=blank//string1
1246          index_lower=index_lower+2
1247        end if
1248      end do ! ipsp
1249 !    if not found ...
1250      if(found==0)then
1251        write(message,'(6a)' )&
1252 &       'Did not find matching pseudopotential for XYZ atomic symbol,',ch10,&
1253 &       'with value ',string2,ch10,&
1254 &       'Action: check that the atoms required by the XYZ file correspond to one psp file.'
1255        MSG_ERROR(message)
1256      end if
1257    end do ! Loop on atoms
1258 !  One should find blanks after the last significant type value
1259    string(index_lower:index_blank+2)=blank
1260  end do ! loop to identify _TYPAX
1261 
1262 !---------------------------------------------------------------------------
1263 
1264 !Here, set up quantities that are related to geometrical description
1265 !of the system (acell,rprim,xred), as well as
1266 !initial velocity(vel), and spin of atoms (spinat), nuclear dipole moments
1267 ! of atoms (nucdipmom),
1268 !the symmetries (symrel,symafm, and tnons)
1269 !and the list of fixed atoms (iatfix,iatfixx,iatfixy,iatfixz).
1270 !Arrays have already been
1271 !dimensioned thanks to the knowledge of msym and mxnatom
1272 
1273 !ji: We need to read the electric field before calling ingeo
1274 !****** Temporary ******
1275 
1276  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'berryopt',tread,'INT')
1277  if(tread==1) dtset%berryopt=intarr(1)
1278 
1279  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'berrysav',tread,'INT')
1280  if(tread==1) dtset%berrysav=intarr(1)
1281 
1282  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'bfield',tread,'DPR')
1283  if (tread==1) dtset%bfield(1:3) = dprarr(1:3)
1284 
1285  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'dfield',tread,'DPR')
1286  if (tread==1) dtset%dfield(1:3) = dprarr(1:3)
1287 
1288  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'efield',tread,'DPR')
1289  if (tread==1) dtset%efield(1:3) = dprarr(1:3)
1290 
1291  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'red_dfield',tread,'DPR')
1292  if (tread==1) dtset%red_dfield(1:3) = dprarr(1:3)
1293 
1294  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'red_efield',tread,'DPR')
1295  if (tread==1) dtset%red_efield(1:3) = dprarr(1:3)
1296 
1297  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'red_efieldbar',tread,'DPR')
1298  if (tread==1) dtset%red_efieldbar(1:3) = dprarr(1:3)
1299 
1300  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'jfielddir',tread,'INT')
1301  if(tread==1) dtset%jfielddir(1:3)=intarr(1:3)
1302 
1303 !We need to know nsppol/nspinor/nspden before calling ingeo
1304  nsppol=dtset%nsppol
1305  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nsppol',tread,'INT')
1306  if(tread==1) nsppol=intarr(1)
1307 
1308 !Alternate SIESTA definition of nsppol
1309  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'SpinPolarized',tread_alt,'LOG')
1310  if(tread_alt==1)then
1311    if(tread==1)then
1312      write(message, '(5a)' )&
1313 &     'nsppol and SpinPolarized cannot be specified simultaneously',ch10,&
1314 &     'for the same dataset.',ch10,&
1315 &     'Action: check the input file.'
1316      call wrtout(std_out,  message,'COLL')
1317      leave=1
1318    else
1319 !    Note that SpinPolarized is a logical input variable
1320      nsppol=1
1321      if(intarr(1)==1)nsppol=2
1322      tread=1
1323    end if
1324  end if
1325  dtset%nsppol=nsppol
1326 
1327  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nspinor',tread,'INT')
1328  if(tread==1) dtset%nspinor=intarr(1)
1329 
1330 !Has to read pawspnorb now, in order to adjust nspinor
1331 !Also, if nspinor=1, turn on spin-orbit coupling by default, here for the PAW case. NC case is treated elsewhere.
1332  if (dtset%usepaw>0)then
1333 !  Change the default value
1334    if(dtset%nspinor==2)dtset%pawspnorb=1
1335    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawspnorb',tread,'INT')
1336    if(tread==1)then
1337      dtset%pawspnorb=intarr(1)
1338      if(dtset%pawspnorb>0) dtset%nspinor=2
1339    else
1340      if(dtset%nspinor==2)then
1341        write(message, '(4a)' ) ch10,&
1342 &       ' invars1: COMMENT -',ch10,&
1343 &       '  With nspinor=2 and usepaw=1, pawspnorb=1 has been switched on by default.'
1344        call wrtout(iout,  message,'COLL')
1345      end if
1346    end if
1347  end if
1348  nspinor=dtset%nspinor
1349 
1350  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nspden',tread,'INT')
1351  if(tread==1) then
1352    dtset%nspden=intarr(1)
1353  else
1354    dtset%nspden=dtset%nsppol
1355  end if
1356 
1357  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ntypalch',tread,'INT')
1358  if(tread==1) dtset%ntypalch=intarr(1)
1359 
1360  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nzchempot',tread,'INT')
1361  if(tread==1) dtset%nzchempot=intarr(1)
1362 
1363  ntypalch=dtset%ntypalch
1364  if(ntypalch>ntypat)then
1365    write(message, '(3a,i0,a,i0,a,a)' )&
1366 &   'The input variable ntypalch must be smaller than ntypat, while it is',ch10,&
1367 &   'ntypalch=',dtset%ntypalch,', and ntypat=',ntypat,ch10,&
1368 &   'Action: check ntypalch vs ntypat in your input file.'
1369    MSG_ERROR(message)
1370  end if
1371 
1372  ntyppure=ntypat-ntypalch
1373  dtset%ntyppure=ntyppure
1374  npspalch=npsp-ntyppure
1375  dtset%npspalch=npspalch
1376  if(npspalch<0)then
1377    write(message, '(a,i0,2a,i0,a,a)' )&
1378 &   'The number of available pseudopotentials, npsp=',npsp,ch10,&
1379 &   'is smaller than the requested number of types of pure atoms, ntyppure=',ntyppure,ch10,&
1380 &   'Action: check ntypalch versus ntypat and npsp in your input file.'
1381    MSG_ERROR(message)
1382  end if
1383 
1384  if(ntypalch>0)then
1385    call intagm(dprarr,intarr,jdtset,marr,ntypalch,string(1:lenstr),'algalch',tread,'INT')
1386    if(tread==1) dtset%algalch(1:ntypalch)=intarr(1:ntypalch)
1387  end if
1388 
1389 !Read the Zeeman field
1390  call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'zeemanfield',tread,'BFI')
1391  if(tread==1) then
1392    if(dtset%nspden == 2)then
1393      write(message,'(7a)')&
1394 &     'A Zeeman field has been specified without noncollinear spins.',ch10,&
1395 &     'Only the z-component of the magnetic field will be used.'
1396      MSG_WARNING(message)
1397    else if (dtset%nspden == 1)then
1398      write(message, '(a,a,a)' )&
1399 &     'A Zeeman field has been specified for a non-spin-polarized calculation.',ch10,&
1400 &     'Action: check the input file.'
1401      MSG_ERROR(message)
1402    end if
1403 
1404    dtset%zeemanfield(1:3) = dprarr(1:3)
1405  end if
1406 
1407 
1408  ABI_ALLOCATE(amu,(ntypat))
1409  ABI_ALLOCATE(mixalch,(npspalch,ntypalch))
1410  ABI_ALLOCATE(vel,(3,natom))
1411  ABI_ALLOCATE(vel_cell,(3,3))
1412  ABI_ALLOCATE(xred,(3,natom))
1413  intimage=2 ; if(dtset%nimage==1)intimage=1
1414  do ii=1,dtset%nimage+1
1415    iimage=ii
1416    if(dtset%nimage==1 .and. ii==2)exit
1417    if(dtset%nimage==2 .and. ii==3)exit
1418    if(dtset%nimage> 2 .and. ii==intimage)cycle ! Will do the intermediate reference image at the last reading
1419    if(dtset%nimage>=2 .and. ii==dtset%nimage+1)iimage=intimage
1420 
1421    write(message,'(a,i0)')' invars1 : treat image number: ',iimage
1422    call wrtout(std_out,message,'COLL')
1423 
1424 !  Need to reset nsym to default value for each image
1425    dtset%nsym=0
1426 
1427 !  Call ingeo for each image in turn, with the possible default values
1428    acell=dtset%acell_orig(1:3,iimage)
1429    amu=dtset%amu_orig(1:ntypat,iimage)
1430    mixalch=dtset%mixalch_orig(1:npspalch,1:ntypalch,iimage)
1431    rprim=dtset%rprim_orig(1:3,1:3,iimage)
1432    vel=dtset%vel_orig(1:3,1:natom,iimage)
1433    vel_cell=dtset%vel_cell_orig(1:3,1:3,iimage)
1434    xred=dtset%xred_orig(1:3,1:natom,iimage)
1435    ABI_ALLOCATE(iatfix,(3,natom))
1436    ABI_ALLOCATE(nucdipmom,(3,natom))
1437    ABI_ALLOCATE(spinat,(3,natom))
1438    ABI_ALLOCATE(typat,(natom))
1439    ABI_ALLOCATE(znucl,(dtset%npsp))
1440    nucdipmom(1:3,1:natom)=dtset%nucdipmom(1:3,1:natom)
1441    spinat(1:3,1:natom)=dtset%spinat(1:3,1:natom)
1442    znucl(1:dtset%npsp)=dtset%znucl(1:dtset%npsp)
1443 
1444    call ingeo(acell,amu,dtset,bravais,dtset%genafm(1:3),iatfix,&
1445 &   dtset%icoulomb,iimage,iout,jdtset,dtset%jellslab,lenstr,mixalch,&
1446 &   msym,natom,dtset%nimage,dtset%npsp,npspalch,dtset%nspden,dtset%nsppol,&
1447 &   dtset%nsym,ntypalch,dtset%ntypat,nucdipmom,dtset%nzchempot,&
1448 &   dtset%pawspnorb,dtset%ptgroupma,ratsph,&
1449 &   rprim,dtset%slabzbeg,dtset%slabzend,dtset%spgroup,spinat,&
1450 &   string,dtset%supercell_latt,symafm,dtset%symmorphi,symrel,tnons,dtset%tolsym,&
1451 &   typat,vel,vel_cell,xred,znucl)
1452    dtset%iatfix(1:3,1:natom)=iatfix(1:3,1:natom)
1453    dtset%nucdipmom(1:3,1:natom)=nucdipmom(1:3,1:natom)
1454    dtset%spinat(1:3,1:natom)=spinat(1:3,1:natom)
1455    dtset%typat(1:natom)=typat(1:natom)
1456    ABI_DEALLOCATE(iatfix)
1457    ABI_DEALLOCATE(nucdipmom)
1458    ABI_DEALLOCATE(spinat)
1459    ABI_DEALLOCATE(typat)
1460    ABI_DEALLOCATE(znucl)
1461    dtset%acell_orig(1:3,iimage)=acell
1462    dtset%amu_orig(1:ntypat,iimage)=amu
1463    dtset%mixalch_orig(1:npspalch,1:ntypalch,iimage)=mixalch
1464    dtset%rprim_orig(1:3,1:3,iimage)=rprim
1465    dtset%vel_orig(1:3,1:natom,iimage)=vel
1466    dtset%vel_cell_orig(1:3,1:3,iimage)=vel_cell
1467    dtset%xred_orig(1:3,1:natom,iimage)=xred
1468    call mkrdim(dtset%acell_orig(1:3,iimage),dtset%rprim_orig(1:3,1:3,iimage),dtset%rprimd_orig(1:3,1:3,iimage))
1469  end do
1470  ABI_DEALLOCATE(amu)
1471  ABI_DEALLOCATE(mixalch)
1472  ABI_DEALLOCATE(vel)
1473  ABI_DEALLOCATE(vel_cell)
1474  ABI_DEALLOCATE(xred)
1475 
1476 !Examine whether there is some vacuum space in the unit cell
1477  call invacuum(jdtset,lenstr,natom,dtset%rprimd_orig(1:3,1:3,intimage),string,vacuum,&
1478 & dtset%xred_orig(1:3,1:natom,intimage))
1479 
1480 !DEBUG
1481 !write(std_out,*)' invars1: before inkpts, dtset%mixalch_orig(1:npspalch,1:ntypalch,:)=',&
1482 !dtset%mixalch_orig(1:npspalch,1:ntypalch,1:dtset%nimage)
1483 !stop
1484 !ENDDEBUG
1485 
1486 !---------------------------------------------------------------------------
1487 
1488 !Set up k point grid number
1489 !First, get additional information
1490  dtset%kptopt=1
1491  if(dtset%nspden==4)dtset%kptopt=4
1492  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'kptopt',tread,'INT')
1493  if(tread==1) dtset%kptopt=intarr(1)
1494 
1495  dtset%qptopt=1
1496  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'qptopt',tread,'INT')
1497  if(tread==1) dtset%qptopt=intarr(1)
1498 
1499  iscf=5
1500  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'iscf',tread,'INT')
1501  if(tread==1) iscf=intarr(1)
1502 
1503 
1504  dtset%natsph=dtset%natom
1505  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natsph',tread,'INT')
1506  if(tread==1) dtset%natsph=intarr(1)
1507 
1508  dtset%natsph_extra=0
1509  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natsph_extra',tread,'INT')
1510  if(tread==1) dtset%natsph_extra=intarr(1)
1511 
1512  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natvshift',tread,'INT')
1513  if(tread==1) dtset%natvshift=intarr(1)
1514 
1515  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nconeq',tread,'INT')
1516  if(tread==1) dtset%nconeq=intarr(1)
1517 
1518  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nkptgw',tread,'INT')
1519  if(tread==1) dtset%nkptgw=intarr(1)
1520  if (dtset%nkptgw<0) then
1521    write(message, '(a,i0,a,a,a,a)' )&
1522 &   'Input nkptgw must be >= 0, but was ',dtset%nkptgw,ch10,&
1523 &   'This is not allowed.',ch10,&
1524 &   'Action: check the input file.'
1525    MSG_ERROR(message)
1526  end if
1527 
1528 !Number of points for long wavelength limit.
1529 !Default is dtset%gw_nqlwl=0
1530  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_nqlwl',tread,'INT')
1531  if(tread==1) dtset%gw_nqlwl=intarr(1)
1532  if (dtset%gw_nqlwl<0) then
1533    write(message, '(a,i12,a,a,a,a)' )&
1534 &   'Input gw_nqlwl must be > 0, but was ',dtset%gw_nqlwl,ch10,&
1535 &   'This is not allowed.',ch10,&
1536 &   'Action: check the input file.'
1537    MSG_ERROR(message)
1538  end if
1539 
1540 !Read number of k-points (if specified)
1541  nkpt=0
1542  if(dtset%kptopt==0)nkpt=1
1543  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nkpt',tread,'INT')
1544  if(tread==1) nkpt=intarr(1)
1545  dtset%nkpt=nkpt
1546 
1547  call chkint_ge(0,0,cond_string,cond_values,ierr,'nkpt',nkpt,0,iout)
1548  if(dtset%kptopt==0)then
1549    cond_string(1)='kptopt' ; cond_values(1)=0
1550    call chkint_ge(1,1,cond_string,cond_values,ierr,'nkpt',nkpt,1,iout)
1551  end if
1552 
1553  nkpthf=nkpt
1554  dtset%nkpthf=nkpt
1555 
1556 !Will compute the actual value of nkpt, if needed. Otherwise,
1557 !test that the value of nkpt is OK, if kptopt/=0
1558 !Set up dummy arrays istwfk, kpt, wtk
1559 
1560  if(nkpt/=0 .or. dtset%kptopt/=0)then
1561    ABI_ALLOCATE(istwfk,(nkpt))
1562    ABI_ALLOCATE(kpt,(3,nkpt))
1563    ABI_ALLOCATE(kpthf,(3,nkpthf))
1564    ABI_ALLOCATE(wtk,(nkpt))
1565 !  Here, occopt is also a dummy argument
1566    occopt=1 ; dtset%nshiftk=1 ; dtset%kptrlatt(:,:)=0
1567 
1568    kptrlen=20.0_dp ; wtk(:)=1.0_dp
1569    dtset%shiftk(:,:)=half
1570 
1571    nqpt=0
1572    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nqpt',tread,'INT')
1573    if(tread==1) nqpt=intarr(1)
1574 
1575    chksymbreak=1
1576    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'chksymbreak',tread,'INT')
1577    if(tread==1) chksymbreak=intarr(1)
1578 
1579 
1580 !  Use the first image to predict k and/or q points,
1581 !  except if an intermediate image is available
1582    intimage=1 ; if(dtset%nimage>2)intimage=(1+dtset%nimage)/2
1583 
1584 !  Find the q-point, if any.
1585    if(nqpt==1)then
1586      call inqpt(chksymbreak,std_out,jdtset,lenstr,msym,natom,dtset%qptn,dtset%wtq,&
1587 &     dtset%rprimd_orig(1:3,1:3,intimage),dtset%spinat,string,dtset%typat,&
1588 &     vacuum,dtset%xred_orig(1:3,1:natom,intimage),dtset%qptrlatt)
1589    end if
1590 
1591 !  Find the k point grid
1592    call inkpts(bravais,chksymbreak,dtset%fockdownsampling,iout,iscf,istwfk,jdtset,&
1593 &   kpt,kpthf,dtset%kptopt,kptnrm,dtset%kptrlatt_orig,dtset%kptrlatt,kptrlen,lenstr,msym,&
1594 &   nkpt,nkpthf,nqpt,dtset%ngkpt,dtset%nshiftk,dtset%nshiftk_orig,dtset%shiftk_orig,dtset%nsym,&
1595 &   occopt,dtset%qptn,response,dtset%rprimd_orig(1:3,1:3,intimage),dtset%shiftk,&
1596 &   string,symafm,symrel,vacuum,wtk)
1597 
1598    ABI_DEALLOCATE(istwfk)
1599    ABI_DEALLOCATE(kpt)
1600    ABI_DEALLOCATE(kpthf)
1601    ABI_DEALLOCATE(wtk)
1602 !  nkpt and nkpthf have been computed, as well as the k point grid, if needed
1603    dtset%nkpt=nkpt
1604    dtset%nkpthf=nkpthf
1605  end if
1606 
1607  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nqptdm',tread,'INT')
1608  if(tread==1) dtset%nqptdm=intarr(1)
1609 
1610  if (dtset%nqptdm<-1) then
1611    write(message, '(a,i12,a,a,a,a)' )&
1612 &   'Input nqptdm must be >= 0, but was ',dtset%nqptdm,ch10,&
1613 &   'This is not allowed.',ch10,&
1614 &   'Action: check the input file.'
1615    MSG_ERROR(message)
1616  end if
1617 
1618  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nzchempot',tread,'INT')
1619  if(tread==1) dtset%nzchempot=intarr(1)
1620 
1621  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'cd_customnimfrqs',tread,'INT')
1622  if(tread==1) dtset%cd_customnimfrqs=intarr(1)
1623 
1624  if (dtset%cd_customnimfrqs<0) then
1625    write(message, '(a,i0,a,a,a,a)' )&
1626 &   'Input cd_customnimfrqs must be >= 0, but was ',dtset%cd_customnimfrqs,ch10,&
1627 &   'This is not allowed.',ch10,&
1628 &   'Action: check the input file.'
1629    MSG_ERROR(message)
1630  end if
1631 
1632  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_customnfreqsp',tread,'INT')
1633  if(tread==1) dtset%gw_customnfreqsp=intarr(1)
1634 
1635  if (dtset%gw_customnfreqsp<0) then
1636    write(message, '(a,i0,a,a,a,a)' )&
1637 &   'Input gw_customnfreqsp must be >= 0, but was ',dtset%gw_customnfreqsp,ch10,&
1638 &   'This is not allowed.',ch10,&
1639 &   'Action: check the input file.'
1640    MSG_ERROR(message)
1641  end if
1642 
1643  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gwls_n_proj_freq',tread,'INT')
1644  if(tread==1) dtset%gwls_n_proj_freq=intarr(1)
1645 
1646  if (dtset%gwls_n_proj_freq<0) then
1647    write(message, '(a,i0,a,a,a,a)' )&
1648 &   'Input gwls_n_proj_freq must be >= 0, but was ',dtset%gwls_n_proj_freq,ch10,&
1649 &   'This is not allowed.',ch10,&
1650 &   'Action: check the input file.'
1651    MSG_ERROR(message)
1652  end if
1653 
1654  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'efmas_calc_dirs',tread,'INT')
1655  if(tread==1) dtset%efmas_calc_dirs=intarr(1)
1656 
1657  if (ABS(dtset%efmas_calc_dirs)>3) then
1658    write(message, '(a,i0,a,a,a,a)' )&
1659 &   'Input efmas_calc_dirs must be between -3 and 3, but was ',dtset%efmas_calc_dirs,ch10,&
1660 &   'This is not allowed.',ch10,&
1661 &   'Action: check the input file.'
1662    MSG_ERROR(message)
1663  end if
1664 
1665  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'efmas_n_dirs',tread,'INT')
1666  if(tread==1) dtset%efmas_n_dirs=intarr(1)
1667 
1668  if (dtset%efmas_n_dirs<0) then
1669    write(message, '(a,i0,a,a,a,a)' )&
1670 &   'Input efmas_n_dirs must be >= 0, but was ',dtset%efmas_n_dirs,ch10,&
1671 &   'This is not allowed.',ch10,&
1672 &   'Action: check the input file.'
1673    MSG_ERROR(message)
1674  end if
1675 !---------------------------------------------------------------------------
1676 
1677  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nnos',tread,'INT')
1678  if(tread==1) dtset%nnos=intarr(1)
1679 
1680  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ga_n_rules',tread,'INT')
1681  if(tread==1) dtset%ga_n_rules=intarr(1)
1682 
1683 !Perform the first checks
1684 
1685  leave=0
1686 
1687 !Check that nkpt is greater than 0
1688  if (nkpt<=0) then
1689    write(message, '(a,a,a,a,i12,a,a,a,a)' ) ch10,&
1690 &   ' invars1: ERROR -',ch10,&
1691 &   '  After inkpts, nkpt must be > 0, but was ',nkpt,ch10,&
1692 &   '  This is not allowed.',ch10,&
1693 &   '  Action: check the input file.'
1694    call wrtout(std_out,  message,'COLL')
1695    leave=1
1696  end if
1697 
1698 !Check that nsppol is 1 or 2
1699  if (nsppol/=1 .and. nsppol/=2) then
1700    write(message, '(a,a,a,a,i12,a,a,a,a)' ) ch10,&
1701 &   ' invars1: ERROR -',ch10,&
1702 &   '  Input nsppol must be 1 or 2, but was ',nsppol,ch10,&
1703 &   '  This is not allowed.',ch10,&
1704 &   '  Action: check the input file.'
1705    call wrtout(std_out,message,'COLL')
1706    leave=1
1707  end if
1708 
1709 !Check that nspinor is 1 or 2
1710  if (nspinor/=1 .and. nspinor/=2) then
1711    write(message, '(a,a,a,a,i12,a,a,a,a)' ) ch10,&
1712 &   ' invars1: ERROR -',ch10,&
1713 &   '  Input nspinor must be 1 or 2, but was ',nspinor,ch10,&
1714 &   '  This is not allowed.',ch10,&
1715 &   '  Action: check the input file.'
1716    call wrtout(std_out,message,'COLL')
1717    leave=1
1718  end if
1719 
1720 !Check that nspinor and nsppol are not 2 together
1721  if (nsppol==2 .and. nspinor==2) then
1722    write(message, '(8a)' ) ch10,&
1723 &   ' invars1: ERROR -',ch10,&
1724 &   '  nspinor and nsappol cannot be 2 together !',ch10,&
1725 &   '  This is not allowed.',ch10,&
1726 &   '  Action: check the input file.'
1727    call wrtout(std_out,message,'COLL')
1728    leave=1
1729  end if
1730 
1731 !Here, leave if an error has been detected earlier
1732  if(leave==1) then
1733    message = ' Other errors might be present in the input file. '
1734    MSG_ERROR(message)
1735  end if
1736 
1737 
1738 !Now, take care of mband_upper
1739 
1740  mband_upper=1
1741  occopt=1
1742  fband=0.5_dp
1743 
1744  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'occopt',tread,'INT')
1745  if(tread==1) occopt=intarr(1)
1746 
1747 !Also read fband, that is an alternative to nband. The default
1748 !is different for occopt==1 and for metallic occupations.
1749  if(occopt==1)fband=0.125_dp
1750  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fband',tfband,'DPR')
1751  if(tfband==1)fband=dprarr(1)
1752 
1753 !fband cannot be used when occopt==0 or occopt==2
1754  if(tfband==1 .and. (occopt==0 .or. occopt==2) )then
1755    write(message, '(3a)' )&
1756 &   'fband cannot be used if occopt==0 or occopt==2 ',ch10,&
1757 &   'Action: correct your input file, suppress fband, or change occopt.'
1758    MSG_ERROR(message)
1759  end if
1760 
1761  ABI_ALLOCATE(nband,(nkpt*nsppol))
1762  tnband=0
1763 
1764 !Compute ziontypat
1765 !When the pseudo-atom is pure, simple copy
1766  if(ntyppure>0)then
1767    do itypat=1,ntyppure
1768      dtset%ziontypat(itypat)=zionpsp(itypat)
1769    end do
1770  end if
1771 !When the pseudo-atom is alchemical, must make mixing
1772  if(ntypalch>0)then
1773    do itypat=ntyppure+1,ntypat
1774      dtset%ziontypat(itypat)=zero
1775      do ipsp=ntyppure+1,npsp
1776        dtset%ziontypat(itypat)=dtset%ziontypat(itypat) &
1777 &       +dtset%mixalch_orig(ipsp-ntyppure,itypat-ntyppure,1)*zionpsp(ipsp)
1778      end do
1779    end do
1780  end if
1781 
1782 
1783 
1784  if (occopt==0 .or. occopt==1 .or. (occopt>=3 .and. occopt<=8) ) then
1785    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nband',tnband,'INT')
1786 !  Note : mband_upper is initialized, not nband
1787    if(tnband==1) mband_upper=intarr(1)
1788 
1789    if(tfband==1 .and. tnband==1)then
1790      write(message, '(3a)' )&
1791 &     'fband and nband cannot be used together. ',ch10,&
1792 &     'Action: correct your input file, suppress either fband or nband.'
1793      MSG_ERROR(message)
1794    end if
1795 
1796 !  In case nband was not read, use fband, either read, or the default,
1797 !  to provide an upper limit for mband_upper
1798    if(tnband==0)then
1799 
1800      charge=0.0_dp
1801      call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'charge',tread,'DPR')
1802      if(tread==1) charge=dprarr(1)
1803 
1804 !    Only take into account negative charge, to compute maximum number of bands
1805      if(charge > 0.0_dp)charge=0.0_dp
1806 
1807 !     mband_upper=nspinor*((nint(zion_max)*natom+1)/2 - floor(charge/2.0_dp)&
1808 !&     + ceiling(fband*natom-1.0d-10))
1809      zval=0.0_dp
1810      do iatom=1,natom
1811        zval=zval+dtset%ziontypat(dtset%typat(iatom))
1812      end do
1813      zelect=zval-charge
1814      mband_upper=nspinor * ((ceiling(zelect-1.0d-10)+1)/2 + ceiling( fband*natom - 1.0d-10 ))
1815 
1816      nband(:)=mband_upper
1817 
1818 !    DEBUG
1819 !    write(std_out,*)' invars1 : zion_max,natom,fband,mband_upper '
1820 !    write(std_out,*)zion_max,natom,fband,mband_upper
1821 !    ENDDEBUG
1822 
1823    end if
1824 
1825    nband(:)=mband_upper
1826 
1827  else if (occopt==2) then
1828    ABI_ALLOCATE(reaalloc,(nkpt*nsppol))
1829    call intagm(reaalloc,nband,jdtset,nkpt*nsppol,nkpt*nsppol,string(1:lenstr),'nband',tnband,'INT')
1830    if(tnband==1)then
1831      do ikpt=1,nkpt*nsppol
1832        if (nband(ikpt)>mband_upper) mband_upper=nband(ikpt)
1833      end do
1834    end if
1835    ABI_DEALLOCATE(reaalloc)
1836  else
1837    write(message, '(a,i0,3a)' )&
1838 &   'occopt=',occopt,' is not an allowed value.',ch10,&
1839 &   'Action: correct your input file.'
1840    MSG_ERROR(message)
1841  end if
1842 
1843 !Check that mband_upper is greater than 0
1844  if (mband_upper<=0) then
1845    write(message, '(a,i0,4a)' )&
1846 &   'Maximal nband must be > 0, but was ',mband_upper,ch10,&
1847 &   'This is not allowed.',ch10,&
1848 &   'Action: check the input file.'
1849    MSG_ERROR(message)
1850  end if
1851 
1852 !The following 3 values are needed to dimension the parallelism over images
1853  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'imgmov',tread,'INT')
1854  if(tread==1) dtset%imgmov=intarr(1)
1855  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ntimimage',tread,'INT')
1856  if(tread==1) dtset%ntimimage=intarr(1)
1857  call intagm(dprarr,intarr,jdtset,marr,dtset%nimage,string(1:lenstr),'dynimage',tread,'INT')
1858  if(tread==1)then
1859    dtset%dynimage(1:dtset%nimage)=intarr(1:dtset%nimage)
1860  else if (dtset%imgmov==2.or.dtset%imgmov==5) then
1861    dtset%dynimage(1)=0;dtset%dynimage(dtset%nimage)=0
1862  end if
1863  dtset%ndynimage=count(dtset%dynimage(1:dtset%nimage)/=0)
1864 
1865  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'wfoptalg',tread,'INT')
1866  if(tread==1) then
1867    dtset%wfoptalg=intarr(1)
1868  else
1869    if (dtset%usepaw==0)    dtset%wfoptalg=0
1870    if (dtset%usepaw/=0)    dtset%wfoptalg=10
1871    if (dtset%optdriver==RUNL_GSTATE) then
1872      if (dtset%paral_kgb/=0) dtset%wfoptalg=14
1873    end if
1874  end if
1875 
1876 !---------------------------------------------------------------------------
1877 !Some PAW+DMFT keywords
1878  dtset%usedmft=0
1879  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usedmft',tread,'INT')
1880  if(tread==1) dtset%usedmft=intarr(1)
1881 
1882 !Some ucrpa keywords
1883  dtset%ucrpa=0
1884  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ucrpa',tread,'INT')
1885  if(tread==1) dtset%ucrpa=intarr(1)
1886 
1887  if (dtset%ucrpa > 0 .and. dtset%usedmft > 0) then
1888    write(message, '(7a)' )&
1889 &   'usedmft and ucrpa are both activated in the input file ',ch10,&
1890 &   'In the following, abinit assume you are doing a ucrpa calculation and ',ch10,&
1891 &   'you define Wannier functions as in DFT+DMFT calculation',ch10,&
1892 &   'If instead, you want to do a full dft+dmft calculation and not only the Wannier construction, use ucrpa=0'
1893    MSG_WARNING(message)
1894  end if
1895 
1896 !Some PAW+U keywords
1897  dtset%usepawu=0
1898  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usepawu',tread,'INT')
1899  if(tread==1) dtset%usepawu=intarr(1)
1900  if ( dtset%usedmft > 0 .and. dtset%usepawu >= 0 ) dtset%usepawu = 1
1901 
1902  if (dtset%usepawu > 0 ) then
1903    write(message, '(7a)' )&
1904 &   'usedmft and usepawu are both activated ',ch10,&
1905 &   'This is not an usual calculation:',ch10,&
1906 &   'usepawu will be put to a value >= 10:',ch10,&
1907 &   'LDA+U potential and energy will be put to zero'
1908    MSG_WARNING(message)
1909  end if
1910 
1911  dtset%usedmatpu=0
1912  dtset%lpawu(1:dtset%ntypat)=-1
1913  if (dtset%usepawu>0.or.dtset%usedmft>0) then
1914 
1915    call intagm(dprarr,intarr,jdtset,marr,dtset%ntypat,string(1:lenstr),'lpawu',tread,'INT')
1916    if(tread==1) dtset%lpawu(1:dtset%ntypat)=intarr(1:dtset%ntypat)
1917 
1918    call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usedmatpu',tread,'INT')
1919    if(tread==1) dtset%usedmatpu=intarr(1)
1920 
1921  end if
1922 
1923 !Some PAW+Exact exchange keywords
1924  dtset%useexexch=0
1925  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'useexexch',tread,'INT')
1926  if(tread==1) dtset%useexexch=intarr(1)
1927 
1928 
1929  dtset%lexexch(1:dtset%ntypat)=-1
1930 
1931  if (dtset%useexexch>0) then
1932    call intagm(dprarr,intarr,jdtset,marr,dtset%ntypat,string(1:lenstr),'lexexch',tread,'INT')
1933    if(tread==1) dtset%lexexch(1:dtset%ntypat)=intarr(1:dtset%ntypat)
1934  end if
1935 ! LDA minus half keyword
1936  call intagm(dprarr,intarr,jdtset,marr,dtset%ntypat,string(1:lenstr),'ldaminushalf',tread,'INT')
1937  if(tread==1) dtset%ldaminushalf(1:dtset%ntypat)=intarr(1:dtset%ntypat)
1938 !Some plowan data
1939  dtset%plowan_natom=0
1940  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'plowan_natom',tread,'INT')
1941  if(tread==1) dtset%plowan_natom=intarr(1)
1942 
1943  dtset%plowan_nt=0
1944  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'plowan_nt',tread,'INT')
1945  if(tread==1) dtset%plowan_natom=intarr(1)
1946 
1947 
1948 !PAW potential zero keyword
1949  dtset%usepotzero=0
1950  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usepotzero',tread,'INT')
1951  if(tread==1) dtset%usepotzero=intarr(1)
1952 
1953 !Macro_uj (determination of U in PAW+U), governs also allocation of atvshift
1954  dtset%macro_uj = 0
1955  call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'macro_uj',tread,'INT')
1956  if(tread==1) dtset%macro_uj=intarr(1)
1957 
1958  ABI_DEALLOCATE(nband)
1959  ABI_DEALLOCATE(ratsph)
1960  ABI_DEALLOCATE(intarr)
1961  ABI_DEALLOCATE(dprarr)
1962 
1963 end subroutine invars1

ABINIT/invars1m [ Functions ]

[ Top ] [ Functions ]

NAME

 invars1m

FUNCTION

 Initialisation phase : prepare the main input subroutine call by
 reading all the NO MULTI variables, as well as the dimensions
 needed for allocating the input arrays in abinit.

INPUTS

  iout=unit number of output file
  lenstr=actual length of string
  msym=default maximal number of symmetries
  mxnatom=maximal value of input natom for all the datasets
  mxnimage=maximal value of input nimage for all the datasets
  ndtset= number of datasets to be read; if 0, no multi-dataset mode
  ndtset_alloc=number of datasets, corrected for allocation of at least
               one data set.
  npsp= number of pseudopotential files
  string*(*)=string of characters containing all input variables and data
  zionpsp(npsp)= valence charge over all psps

OUTPUT

  dmatpuflag=flag controlling the use of an initial density matrix in PAW+U (max. value over datasets)
  mband_upper_(0:ndtset_alloc)=list of mband_upper values
  mxga_n_rules=maximal value of input ga_n_rules for all the datasets
  mxgw_nqlwl=maximal value of input gw_nqlwl for all the datasets
  mxlpawu=maximal value of input lpawu for all the datasets
  mxmband_upper=maximal value of input nband for all the datasets
  mxnatpawu=maximal value of number of atoms on which +U is applied for all the datasets
  mxnatsph=maximal value of input natsph for all the datasets
  mxnatsph_extra=maximal value of input natsph_extra for all the datasets
  mxnatvshift=maximal value of input natvshift for all the datasets
  mxnconeq=maximal value of input nconeq for all the datasets
  mxnkptgw=maximal value of input nkptgw for all the datasets
  mxnkpthf=maximal value of input nkpthf for all the datasets
  mxnkpt=maximal value of input nkpt for all the datasets
  mxnnos=maximal value of input nnos for all the datasets
  mxnqptdm=maximal value of input nqptdm for all the datasets
  mxnspinor=maximal value of input nspinor for all the datasets
  mxnsppol=maximal value of input nsppol for all the datasets
  mxnsym=maximum number of symmetries
  mxntypat=maximum number of types of atoms
  mxnzchempot=maximal value of input nzchempot for all the datasets

SIDE EFFECTS

  dtsets(0:ndtset_alloc)=<type datafiles_type>contains all input variables,
   some of which are initialized here (see invars1.f for more details on the
   initialized records)

PARENTS

      m_ab7_invars_f90

CHILDREN

      indefo1,invars1

SOURCE

617 subroutine invars1m(dmatpuflag,dtsets,iout,lenstr,mband_upper_,&
618 & msym,mxga_n_rules,mxgw_nqlwl,mxlpawu,mxmband_upper,mxnatom,&
619 & mxnatpawu,mxnatsph,mxnatsph_extra,mxnatvshift,mxnconeq,&
620 & mxnimage,mxn_efmas_dirs,mxnkpt,mxnkptgw,mxnkpthf,mxnnos,mxnqptdm,mxnspinor, &
621 & mxnsppol,mxnsym,mxntypat,mxnimfrqs,mxnfreqsp,mxnzchempot,&
622 & mxn_projection_frequencies,ndtset,ndtset_alloc,string,npsp,zionpsp)
623 
624 
625 !This section has been created automatically by the script Abilint (TD).
626 !Do not modify the following lines by hand.
627 #undef ABI_FUNC
628 #define ABI_FUNC 'invars1m'
629 !End of the abilint section
630 
631  implicit none
632 
633 !Arguments ------------------------------------
634 !scalars
635  integer,intent(in) :: iout,lenstr,msym,mxnatom,mxnimage,ndtset,ndtset_alloc,npsp
636  integer,intent(out) :: dmatpuflag,mxga_n_rules,mxgw_nqlwl,mxlpawu,mxmband_upper,mxnatpawu
637  integer,intent(out) :: mxnatsph, mxnatsph_extra
638  integer,intent(out) :: mxnatvshift,mxnconeq,mxn_efmas_dirs,mxnkpt,mxnkptgw,mxnkpthf,mxnnos
639  integer,intent(out) :: mxnqptdm,mxnspinor,mxnsppol,mxnsym,mxntypat
640  integer,intent(out) :: mxnimfrqs,mxnfreqsp,mxnzchempot,mxn_projection_frequencies
641  character(len=*),intent(inout) :: string
642 !arrays
643  integer,intent(out) :: mband_upper_(0:ndtset_alloc)
644  type(dataset_type),intent(inout) :: dtsets(0:ndtset_alloc)
645  real(dp),intent(in) :: zionpsp(npsp)
646 
647 !Local variables-------------------------------
648 !scalars
649  integer :: idtset,ii,jdtset,lpawu,mband_upper,iatom,nat,nsp
650 !arrays
651  integer,allocatable :: symafm_(:,:),symrel_(:,:,:,:)
652  integer,allocatable :: symafm(:),symrel(:,:,:)
653  real(dp),allocatable :: tnons_(:,:,:),tnons(:,:)
654 
655 !******************************************************************
656 
657 !Here, allocation of the arrays that depend on msym.
658  ABI_ALLOCATE(symrel_,(3,3,msym,0:ndtset_alloc))
659  ABI_ALLOCATE(symafm_,(msym,0:ndtset_alloc))
660  ABI_ALLOCATE(tnons_,(3,msym,0:ndtset_alloc))
661  ABI_ALLOCATE(symafm,(msym))
662  ABI_ALLOCATE(symrel,(3,3,msym))
663  ABI_ALLOCATE(tnons,(3,msym))
664 
665 !Set up default values (note that the default acell, amu
666 !mkmem, mkmem1,mkqmem, and nkpt must be overcome
667 
668  do idtset=0,ndtset_alloc
669    call indefo1(dtsets(idtset))
670  end do
671 
672 !natom and nimage are already initialized in invars0
673  dtsets(0)%natom=-1
674  dtsets(0)%nimage=1
675 
676 !Initialization for parallelization data has changed
677 !these lines aim to keep old original default values
678  dtsets(0)%npimage=1
679  dtsets(0)%npkpt=1
680  dtsets(0)%npspinor=1
681  dtsets(0)%npfft=1
682  dtsets(0)%npband=1
683  dtsets(0)%bandpp=1
684 
685  symafm_(:,0)=1
686  symrel_(:,:,:,0)=0
687  symrel_(1,1,:,0)=1 ; symrel_(2,2,:,0)=1 ; symrel_(3,3,:,0)=1
688  tnons_(:,:,0)=0.0_dp
689 
690 !Loop on datasets
691  do idtset=1,ndtset_alloc
692    !write(std_out,'(2a,i0)') ch10,' invars1m : enter jdtset= ',jdtset
693    jdtset=dtsets(idtset)%jdtset ; if(ndtset==0)jdtset=0
694 
695 !  Input default values
696    dtsets(idtset)%bravais(:)=0
697    symafm(:)=symafm_(:,0)
698    symrel(:,:,:)=symrel_(:,:,:,0)
699    tnons(:,:)=tnons_(:,:,0)
700 
701    call invars1(dtsets(idtset)%bravais,dtsets(idtset),iout,jdtset,lenstr,&
702 &   mband_upper,msym,npsp,string,symafm,symrel,tnons,zionpsp)
703 
704    mband_upper_ (idtset)=mband_upper
705    symafm_(:,idtset)=symafm(:)
706    symrel_(:,:,:,idtset)=symrel(:,:,:)
707    tnons_(:,:,idtset)=tnons(:,:)
708  end do
709 
710  mxmband_upper =maxval(mband_upper_ (1:ndtset_alloc))
711 
712  dmatpuflag=0;mxnatpawu=0;mxlpawu=0
713  mxnatsph=dtsets(1)%natsph
714  mxnatsph_extra=dtsets(1)%natsph_extra
715  mxnatvshift=dtsets(1)%natvshift
716  mxnconeq=dtsets(1)%nconeq
717  mxn_efmas_dirs=0
718  mxga_n_rules = dtsets(1)%ga_n_rules
719  mxgw_nqlwl = dtsets(1)%gw_nqlwl
720  mxnimfrqs = 0
721  mxnfreqsp = 0
722  mxn_projection_frequencies=0
723  mxnkpt  =dtsets(1)%nkpt
724  mxnkptgw=dtsets(1)%nkptgw
725  mxnkpthf=dtsets(1)%nkpthf
726  mxnnos  =dtsets(1)%nnos
727  mxnqptdm=dtsets(1)%nqptdm
728  mxnspinor=dtsets(1)%nspinor
729  mxnsppol=dtsets(1)%nsppol
730  mxntypat=dtsets(1)%ntypat
731  mxnzchempot=dtsets(1)%nzchempot
732 
733 !Get MAX dimension over datasets
734  do ii=1,ndtset_alloc
735    mxnatsph=max(dtsets(ii)%natsph,mxnatsph)
736    mxnatsph_extra=max(dtsets(ii)%natsph_extra,mxnatsph_extra)
737    mxnconeq=max(dtsets(ii)%nconeq,mxnconeq)
738    mxn_efmas_dirs=max(dtsets(ii)%efmas_n_dirs,mxn_efmas_dirs)
739    mxga_n_rules = max(dtsets(ii)%ga_n_rules,mxga_n_rules)
740    mxgw_nqlwl = max(dtsets(ii)%gw_nqlwl,mxgw_nqlwl)
741    mxnimfrqs = max(dtsets(ii)%cd_customnimfrqs,mxnimfrqs)
742    mxnfreqsp = max(dtsets(ii)%gw_customnfreqsp,mxnfreqsp)
743    mxn_projection_frequencies = max(dtsets(ii)%gwls_n_proj_freq,mxn_projection_frequencies)
744    mxnkpt  =max(dtsets(ii)%nkpt,mxnkpt)
745    mxnkptgw=max(dtsets(ii)%nkptgw,mxnkptgw)
746    mxnkpthf=max(dtsets(ii)%nkpthf,mxnkpthf)
747    mxnnos  =max(dtsets(ii)%nnos,mxnnos)
748    mxnqptdm=max(dtsets(ii)%nqptdm,mxnqptdm)
749    mxnspinor=max(dtsets(ii)%nspinor,mxnspinor)
750    mxnsppol=max(dtsets(ii)%nsppol,mxnsppol)
751    mxntypat=max(dtsets(ii)%ntypat,mxntypat)
752    mxnzchempot=max(dtsets(ii)%nzchempot,mxnzchempot)
753    if (dtsets(ii)%usepawu>0) then
754      if (dtsets(ii)%usedmatpu/=0) dmatpuflag=1
755      lpawu=maxval(dtsets(ii)%lpawu(:))
756      mxlpawu=max(lpawu,mxlpawu)
757      !dtsets(ii)%natpawu=count(dtsets(ii)%lpawu(dtsets(ii)%typat((/(i1,i1=1,dtsets(ii)%natom)/)))/=-1)
758      ! Old fashon way that should do fine
759      dtsets(ii)%natpawu = 0
760      do iatom=1, dtsets(ii)%natom
761        if (dtsets(ii)%lpawu(dtsets(ii)%typat(iatom)) /= -1 ) dtsets(ii)%natpawu = dtsets(ii)%natpawu + 1
762      end do
763      mxnatpawu=max(dtsets(ii)%natpawu,mxnatpawu)
764      if (dtsets(ii)%macro_uj/=0) dtsets(ii)%natvshift=lpawu*2+1
765    end if
766    mxnatvshift=max(dtsets(ii)%natvshift,mxnatvshift)
767  end do
768 
769 !mxnsym=maxval(dtsets(1:ndtset_alloc)%nsym) ! This might not work properly with HP compiler
770  mxnsym=dtsets(1)%nsym
771  do idtset=1,ndtset_alloc
772    mxnsym=max(dtsets(idtset)%nsym,mxnsym)
773  end do
774 
775  do idtset=0,ndtset_alloc
776    ABI_ALLOCATE(dtsets(idtset)%atvshift,(mxnatvshift,mxnsppol,mxnatom))
777    ABI_ALLOCATE(dtsets(idtset)%bs_loband,(mxnsppol))
778    ABI_ALLOCATE(dtsets(idtset)%bdgw,(2,mxnkptgw,mxnsppol))
779    ABI_ALLOCATE(dtsets(idtset)%cd_imfrqs,(mxnimfrqs))
780    ABI_ALLOCATE(dtsets(idtset)%chempot,(3,mxnzchempot,mxntypat))
781    nsp=max(mxnsppol,mxnspinor);nat=mxnatpawu*dmatpuflag
782    ABI_ALLOCATE(dtsets(idtset)%dmatpawu,(2*mxlpawu+1,2*mxlpawu+1,nsp,nat,mxnimage))
783    ABI_ALLOCATE(dtsets(idtset)%efmas_bands,(2,mxnkpt))
784    ABI_ALLOCATE(dtsets(idtset)%efmas_dirs,(3,mxn_efmas_dirs))
785    ABI_ALLOCATE(dtsets(idtset)%gw_freqsp,(mxnfreqsp))
786    ABI_ALLOCATE(dtsets(idtset)%gwls_list_proj_freq,(mxn_projection_frequencies))
787    ABI_ALLOCATE(dtsets(idtset)%gw_qlwl,(3,mxgw_nqlwl))
788    ABI_ALLOCATE(dtsets(idtset)%kpt,(3,mxnkpt))
789    ABI_ALLOCATE(dtsets(idtset)%kptgw,(3,mxnkptgw))
790    ABI_ALLOCATE(dtsets(idtset)%kptns,(3,mxnkpt))
791    ABI_ALLOCATE(dtsets(idtset)%kptns_hf,(3,mxnkpthf))
792    ABI_ALLOCATE(dtsets(idtset)%iatsph,(mxnatsph))
793    ABI_ALLOCATE(dtsets(idtset)%istwfk,(mxnkpt))
794    ABI_ALLOCATE(dtsets(idtset)%nband,(mxnkpt*mxnsppol))
795    ABI_ALLOCATE(dtsets(idtset)%occ_orig,(mxmband_upper*mxnkpt*mxnsppol,mxnimage))
796    ABI_ALLOCATE(dtsets(idtset)%qmass,(mxnnos))
797    ABI_ALLOCATE(dtsets(idtset)%qptdm,(3,mxnqptdm))
798    ABI_ALLOCATE(dtsets(idtset)%symafm,(mxnsym))
799    ABI_ALLOCATE(dtsets(idtset)%symrel,(3,3,mxnsym))
800    ABI_ALLOCATE(dtsets(idtset)%tnons,(3,mxnsym))
801    ABI_ALLOCATE(dtsets(idtset)%wtatcon,(3,mxnatom,mxnconeq))
802    ABI_ALLOCATE(dtsets(idtset)%wtk,(mxnkpt))
803    ABI_ALLOCATE(dtsets(idtset)%xredsph_extra,(3,mxnatsph_extra))
804    dtsets(idtset)%symrel(:,:,:)=symrel_(:,:,1:mxnsym,idtset)
805    dtsets(idtset)%symafm(:)    =symafm_(1:mxnsym,idtset)
806    dtsets(idtset)%tnons (:,:)  =tnons_ (:,1:mxnsym,idtset)
807  end do
808 
809  ABI_DEALLOCATE(symafm_)
810  ABI_DEALLOCATE(symrel_)
811  ABI_DEALLOCATE(tnons_)
812  ABI_DEALLOCATE(symafm)
813  ABI_DEALLOCATE(symrel)
814  ABI_DEALLOCATE(tnons)
815 
816 end subroutine invars1m

ABINIT/m_invars1 [ Modules ]

[ Top ] [ Modules ]

NAME

  m_invars1

FUNCTION

COPYRIGHT

 Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR, AR, MKV, FF, MM)
  This file is distributed under the terms of the
  GNU General Public License, see ~abinit/COPYING
  or http://www.gnu.org/copyleft/gpl.txt .

PARENTS

CHILDREN

SOURCE

21 #if defined HAVE_CONFIG_H
22 #include "config.h"
23 #endif
24 
25 #include "abi_common.h"
26 
27 module m_invars1
28 
29  use defs_basis
30  !use defs_abitypes
31  use m_abicore
32  use m_xmpi
33  use m_errors
34  use m_atomdata
35 
36  use defs_abitypes,  only : dataset_type
37  use m_fstrings, only : inupper
38  use m_geometry, only : mkrdim
39  use m_parser,   only : intagm, chkint_ge
40  use m_inkpts,   only : inkpts, inqpt
41  use m_ingeo,    only : ingeo, invacuum
42  use m_symtk,   only : mati3det
43 #if defined HAVE_GPU_CUDA
44  use m_initcuda, only : Get_ndevice
45 #endif
46 
47  implicit none
48 
49  private