TABLE OF CONTENTS
ABINIT/indefo [ Functions ]
NAME
indefo
FUNCTION
Initialisation phase: default values for most input variables (some are initialized earlier, see indefo1 routine, or even at the definition of the input variables (m_dtset.F90))
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. NOTE that Scalars and static arrays can be initialized directly at the level of the datatype declaration provided the value does not depend on runtime conditions.
SOURCE
2072 subroutine indefo(dtsets, ndtset_alloc, nprocs) 2073 2074 use m_gwdefs 2075 #if defined DEV_YP_VDWXC 2076 use m_xc_vdw 2077 #endif 2078 2079 use m_fftcore, only : get_cache_kb, fftalg_for_npfft 2080 2081 !Arguments ------------------------------------ 2082 !scalars 2083 integer,intent(in) :: ndtset_alloc,nprocs 2084 !arrays 2085 type(dataset_type),intent(inout) :: dtsets(0:ndtset_alloc) !vz_i 2086 2087 !Local variables ------------------------------- 2088 !scalars 2089 integer :: idtset,ii,jdtset,paral_atom_default 2090 logical :: wvl_bigdft 2091 #if defined DEV_YP_VDWXC 2092 type(xc_vdw_type) :: vdw_defaults 2093 #endif 2094 2095 !****************************************************************** 2096 2097 DBG_ENTER("COLL") 2098 2099 !Set up default values. All variables to be output in outvars.f 2100 !should have a default, even if a nonsensible one can be chosen to garantee print in that routine. 2101 !Some default values are also set at the definition of the input variables (m_dtset.F90). 2102 2103 !These variables have already been initialized, for idtset/=0 2104 dtsets(0)%istatr=0 2105 dtsets(0)%istatshft=1 2106 dtsets(0)%kptrlatt(1:3,1:3)=0 2107 !dtsets(0)%kptrlatt_orig=0 2108 dtsets(0)%qptrlatt(1:3,1:3)=0 2109 dtsets(0)%ptgroupma=0 2110 dtsets(0)%spgroup=0 2111 dtsets(0)%shiftk(:,:)=half 2112 !XG20200801 Changed the default value. This default value is also defined in m_ingeo.F90 . Must be coherent ! 2113 !dtsets(0)%tolsym=tol8 2114 dtsets(0)%tolsym=tol5 2115 dtsets(0)%znucl(:)=zero 2116 dtsets(0)%ucrpa=0 2117 dtsets(0)%usedmft=0 2118 2119 paral_atom_default=0 2120 if (nprocs>1.and.maxval(dtsets(:)%usepaw)>0) paral_atom_default=1 2121 2122 !WARNING: set default in all datasets, including idtset=0 !!! 2123 !Use alphabetic order 2124 2125 do idtset=0,ndtset_alloc 2126 jdtset=dtsets(idtset)%jdtset 2127 2128 wvl_bigdft=.false. 2129 if(dtsets(idtset)%usewvl==1 .and. dtsets(idtset)%wvl_bigdft_comp==1) wvl_bigdft=.true. 2130 ! Special case of use_gpu_cuda (can be undertermined at this point) 2131 ! use_gpu_cuda=-1 means undetermined ; here impose its value due to some restrictions 2132 if (dtsets(idtset)%use_gpu_cuda==-1) then 2133 if (dtsets(idtset)%optdriver/=0.or. dtsets(idtset)%tfkinfunc/=0.or. dtsets(idtset)%nspinor/=1) then 2134 dtsets(idtset)%use_gpu_cuda=0 2135 else 2136 dtsets(idtset)%use_gpu_cuda=1 2137 end if 2138 end if 2139 2140 ! A 2141 ! Here we change the default value of iomode according to the configuration options. 2142 ! Ideally, all the sequential tests should pass independently of the default value. 2143 ! The parallel tests may require IO_MODE_MPI or, alternatively, IO_MODE_ETSF with HDF5 support. 2144 ! MG FIXME Sun Sep 6 2015: Many tests fail if IO_MODE_MPI is used as default. IO errors in v1, v2 ... 2145 ! with np=1 and wonderful deadlocks if np>1. 2146 2147 ! Note that this default value might be overriden for specific datasets later, in case of parallelism 2148 dtsets(idtset)%iomode=IO_MODE_FORTRAN 2149 #ifdef HAVE_NETCDF_DEFAULT 2150 dtsets(idtset)%iomode=IO_MODE_ETSF 2151 #endif 2152 #ifdef HAVE_MPI_IO_DEFAULT 2153 dtsets(idtset)%iomode=IO_MODE_MPI 2154 #endif 2155 2156 dtsets(idtset)%adpimd=0 2157 dtsets(idtset)%adpimd_gamma=one 2158 dtsets(idtset)%accuracy=0 2159 dtsets(idtset)%atvshift(:,:,:)=zero 2160 dtsets(idtset)%auxc_ixc=11 2161 dtsets(idtset)%auxc_scal=one 2162 ! B 2163 dtsets(idtset)%bdberry(1:4)=0 2164 dtsets(idtset)%bdeigrf=-1 2165 dtsets(idtset)%bdgw=0 2166 dtsets(idtset)%berrystep=1 2167 dtsets(idtset)%bmass=ten 2168 dtsets(idtset)%boxcenter(1:3)=half 2169 dtsets(idtset)%boxcutmin=two 2170 dtsets(idtset)%brvltt=0 2171 dtsets(idtset)%bs_nstates=0 2172 dtsets(idtset)%bs_hayd_term=1 2173 dtsets(idtset)%builtintest=0 2174 dtsets(idtset)%bxctmindg=two 2175 ! C 2176 dtsets(idtset)%cd_halfway_freq=3.674930883_dp !(100 eV) 2177 dtsets(idtset)%cd_imfrqs(:) = zero 2178 dtsets(idtset)%cd_max_freq=36.74930883_dp !(1000 eV) 2179 dtsets(idtset)%cd_subset_freq(1:2)=0 2180 dtsets(idtset)%cd_frqim_method=1 2181 dtsets(idtset)%cd_full_grid=0 2182 dtsets(idtset)%cellcharge(:)=zero 2183 dtsets(idtset)%chempot(:,:,:)=zero 2184 dtsets(idtset)%chkdilatmx=1 2185 dtsets(idtset)%chkexit=0 2186 dtsets(idtset)%chksymbreak=1 2187 dtsets(idtset)%chksymtnons=1 2188 dtsets(idtset)%cineb_start=7 2189 dtsets(idtset)%corecs(:) = zero 2190 dtsets(idtset)%cprj_update_lvl=0 2191 ! D 2192 dtsets(idtset)%ddamp=0.1_dp 2193 dtsets(idtset)%delayperm=0 2194 dtsets(idtset)%densfor_pred=2 2195 if (dtsets(idtset)%paral_kgb>0.and.idtset>0) dtsets(idtset)%densfor_pred=6 ! Recommended for band-FFT parallelism 2196 dtsets(idtset)%dfpt_sciss=zero 2197 dtsets(idtset)%diecut=2.2_dp 2198 dtsets(idtset)%dielng=1.0774841_dp 2199 dtsets(idtset)%diemac=1.0d6 2200 if (dtsets(idtset)%usepaw==0) then 2201 dtsets(idtset)%diemix=one 2202 else 2203 dtsets(idtset)%diemix=0.7_dp 2204 end if 2205 dtsets(idtset)%diemixmag=dtsets(idtset)%diemix 2206 dtsets(idtset)%diegap=0.1_dp 2207 dtsets(idtset)%dielam=half 2208 dtsets(idtset)%diismemory=8 2209 dtsets(idtset)%dilatmx=one 2210 dtsets(idtset)%dmatpuopt=2 2211 if (size(dtsets(idtset)%dmatpawu,4)>0) dtsets(idtset)%dmatpawu=-10._dp 2212 dtsets(idtset)%dmatudiag=0 2213 dtsets(idtset)%dmft_entropy=0 2214 dtsets(idtset)%dmft_dc =1 2215 dtsets(idtset)%dmft_iter=0 2216 dtsets(idtset)%dmft_kspectralfunc=0 2217 dtsets(idtset)%dmft_nlambda=6 2218 dtsets(idtset)%dmft_nwli=0 2219 dtsets(idtset)%dmft_nwlo=0 2220 dtsets(idtset)%dmft_mxsf=0.3_dp 2221 dtsets(idtset)%dmft_occnd_imag=1 2222 dtsets(idtset)%dmft_read_occnd=0 2223 dtsets(idtset)%dmft_rslf=0 2224 dtsets(idtset)%dmft_solv=5 2225 if(dtsets(idtset)%ucrpa>0.and.dtsets(idtset)%usedmft==1) dtsets(idtset)%dmft_solv=0 2226 dtsets(idtset)%dmft_t2g=0 2227 ! dtsets(idtset)%dmft_x2my2d=0 2228 dtsets(idtset)%dmft_tolfreq=tol4 2229 dtsets(idtset)%dmft_tollc=tol5 2230 dtsets(idtset)%dmft_charge_prec=tol6 2231 dtsets(idtset)%dmft_wanorthnorm=3 2232 dtsets(idtset)%dmftbandi=0 2233 dtsets(idtset)%dmftbandf=0 2234 dtsets(idtset)%dmftcheck=0 2235 dtsets(idtset)%dmftctqmc_basis =1 2236 dtsets(idtset)%dmftctqmc_check =0 2237 dtsets(idtset)%dmftctqmc_correl=0 2238 dtsets(idtset)%dmftctqmc_grnns =0 2239 dtsets(idtset)%dmftctqmc_config =0 2240 dtsets(idtset)%dmftctqmc_meas =1 2241 dtsets(idtset)%dmftctqmc_mrka =0 2242 dtsets(idtset)%dmftctqmc_mov =0 2243 dtsets(idtset)%dmftctqmc_order =0 2244 dtsets(idtset)%dmftctqmc_triqs_nleg=30 2245 dtsets(idtset)%dmftqmc_l=0 2246 dtsets(idtset)%dmftqmc_n=0.0_dp 2247 dtsets(idtset)%dmftqmc_seed=jdtset 2248 dtsets(idtset)%dmftqmc_therm=1000 2249 dtsets(idtset)%dmftctqmc_gmove = dtsets(idtset)%dmftqmc_therm / 10 2250 dtsets(idtset)%dosdeltae=0.0 2251 dtsets(idtset)%dtion=100.0_dp 2252 dtsets(idtset)%dtele=0.1_dp 2253 dtsets(idtset)%d3e_pert1_atpol(1:2)=-1 2254 dtsets(idtset)%d3e_pert1_dir(1:3)=1 2255 dtsets(idtset)%d3e_pert1_elfd=0 2256 dtsets(idtset)%d3e_pert1_phon=0 2257 dtsets(idtset)%d3e_pert2_atpol(1:2)=-1 2258 dtsets(idtset)%d3e_pert2_dir(1:3)=1 2259 dtsets(idtset)%d3e_pert2_elfd=0 2260 dtsets(idtset)%d3e_pert2_phon=0 2261 dtsets(idtset)%d3e_pert2_strs=0 2262 dtsets(idtset)%d3e_pert3_atpol(1:2)=-1 2263 dtsets(idtset)%d3e_pert3_dir(1:3)=1 2264 dtsets(idtset)%d3e_pert3_elfd=0 2265 dtsets(idtset)%d3e_pert3_phon=0 2266 ! E 2267 dtsets(idtset)%ecut=-one 2268 dtsets(idtset)%ecuteps=zero 2269 dtsets(idtset)%ecutsigx=zero ! If ecutsigx is not defined explicitly, npwsigx will be initialized from ecutwfn. 2270 dtsets(idtset)%ecutsm=zero 2271 dtsets(idtset)%ecutwfn=zero ! The true default value is ecut . This is defined in invars2.F90 2272 dtsets(idtset)%effmass_free=one 2273 dtsets(idtset)%efmas=0 2274 dtsets(idtset)%efmas_bands=0 ! The true default is nband. This is defined in invars2.F90 2275 dtsets(idtset)%efmas_deg=1 2276 dtsets(idtset)%efmas_deg_tol=tol5 2277 dtsets(idtset)%efmas_dim=3 2278 dtsets(idtset)%efmas_dirs=zero 2279 dtsets(idtset)%efmas_ntheta=1000 2280 dtsets(idtset)%elph2_imagden=zero 2281 dtsets(idtset)%enunit=0 2282 dtsets(idtset)%eshift=zero 2283 dtsets(idtset)%esmear=0.01_dp 2284 dtsets(idtset)%exchn2n3d=0 2285 dtsets(idtset)%extrapwf=0 2286 dtsets(idtset)%exchmix=quarter 2287 dtsets(idtset)%expert_user=0 2288 ! F 2289 dtsets(idtset)%focktoldfe=zero 2290 dtsets(idtset)%fockoptmix=0 2291 dtsets(idtset)%fockdownsampling(:)=1 2292 dtsets(idtset)%fock_icutcoul=3 2293 dtsets(idtset)%freqim_alpha=five 2294 dtsets(idtset)%friction=0.001_dp 2295 dtsets(idtset)%frzfermi=0 2296 dtsets(idtset)%fxcartfactor=one ! Should be adjusted to the H2 conversion factor 2297 ! G 2298 dtsets(idtset)%ga_algor =1 2299 dtsets(idtset)%ga_fitness =1 2300 dtsets(idtset)%ga_opt_percent =0.2_dp 2301 dtsets(idtset)%ga_rules(:) =1 2302 dtsets(idtset)%goprecon =0 2303 dtsets(idtset)%goprecprm(:)=0 2304 dtsets(idtset)%gpu_devices=(/-1,-1,-1,-1,-1/) 2305 dtsets(idtset)%gpu_linalg_limit=2000000 2306 if (dtsets(idtset)%gw_customnfreqsp/=0) dtsets(idtset)%gw_freqsp(:) = zero 2307 if ( dtsets(idtset)%gw_nqlwl > 0 ) then 2308 dtsets(idtset)%gw_qlwl(:,:)=zero 2309 dtsets(idtset)%gw_qlwl(1,1)=0.00001_dp 2310 dtsets(idtset)%gw_qlwl(2,1)=0.00002_dp 2311 dtsets(idtset)%gw_qlwl(3,1)=0.00003_dp 2312 end if 2313 dtsets(idtset)%gw_frqim_inzgrid=0 2314 dtsets(idtset)%gw_frqre_inzgrid=0 2315 dtsets(idtset)%gw_frqre_tangrid=0 2316 dtsets(idtset)%gw_invalid_freq=0 2317 dtsets(idtset)%gw_icutcoul=6 2318 dtsets(idtset)%gw_qprange=0 2319 dtsets(idtset)%gw_sigxcore=0 2320 dtsets(idtset)%gwls_stern_kmax=1 2321 dtsets(idtset)%gwls_model_parameter=1.0_dp 2322 dtsets(idtset)%gwls_npt_gauss_quad=10 2323 dtsets(idtset)%gwls_diel_model=2 2324 dtsets(idtset)%gwls_print_debug=0 2325 if (dtsets(idtset)%gwls_n_proj_freq/=0) dtsets(idtset)%gwls_list_proj_freq(:) = zero 2326 dtsets(idtset)%gwls_nseeds=1 2327 dtsets(idtset)%gwls_recycle=2 2328 dtsets(idtset)%gwls_kmax_complement=1 2329 dtsets(idtset)%gwls_kmax_poles=4 2330 dtsets(idtset)%gwls_kmax_analytic=8 2331 dtsets(idtset)%gwls_kmax_numeric=16 2332 dtsets(idtset)%gwls_band_index=1 2333 dtsets(idtset)%gwls_exchange=1 2334 dtsets(idtset)%gwls_correlation=3 2335 dtsets(idtset)%gwls_first_seed=0 2336 ! H 2337 dtsets(idtset)%hmcsst=3 2338 dtsets(idtset)%hmctt=4 2339 dtsets(idtset)%hyb_mixing=-999.0_dp 2340 dtsets(idtset)%hyb_mixing_sr=-999.0_dp 2341 dtsets(idtset)%hyb_range_dft=-999.0_dp 2342 dtsets(idtset)%hyb_range_fock=-999.0_dp 2343 ! I 2344 if(dtsets(idtset)%natsph/=0) then 2345 ! do not use iatsph(:) but explicit boundaries 2346 ! to avoid to read to far away in the built array (/ ... /) 2347 dtsets(idtset)%iatsph(1:dtsets(idtset)%natsph)=(/ (ii,ii=1,dtsets(idtset)%natsph) /) 2348 else 2349 dtsets(idtset)%iatsph(:)=0 2350 end if 2351 dtsets(idtset)%iboxcut=0 2352 dtsets(idtset)%icutcoul=3 2353 dtsets(idtset)%ieig2rf=0 2354 dtsets(idtset)%imgwfstor=0 2355 dtsets(idtset)%intxc=0 2356 ! if (dtsets(idtset)%paral_kgb>0.and.idtset>0) dtsets(idtset)%intxc=0 2357 dtsets(idtset)%ionmov=0 2358 dtsets(idtset)%densfor_pred=2 2359 if (dtsets(idtset)%paral_kgb>0.and.idtset>0) dtsets(idtset)%densfor_pred=6 ! Recommended for band-FFT parallelism 2360 dtsets(idtset)%iprcel=0 2361 dtsets(idtset)%iprcfc=0 2362 dtsets(idtset)%irandom=3 2363 !iscf 2364 if(wvl_bigdft) then 2365 dtsets(idtset)%iscf=0 2366 else 2367 if(dtsets(idtset)%usepaw==0) then 2368 dtsets(idtset)%iscf=7 2369 else 2370 dtsets(idtset)%iscf=17 2371 end if 2372 end if 2373 dtsets(idtset)%isecur=0 2374 dtsets(idtset)%istatimg = 1 2375 dtsets(idtset)%istwfk(:)=0 2376 dtsets(idtset)%ixc=1 2377 dtsets(idtset)%ixc_sigma=1 2378 dtsets(idtset)%ixcpositron=1 2379 dtsets(idtset)%ixcrot=1 2380 ! J 2381 dtsets(idtset)%f4of2_sla(:)=-one 2382 dtsets(idtset)%f6of2_sla(:)=-one 2383 dtsets(idtset)%jpawu(:,:)=zero 2384 ! K 2385 dtsets(idtset)%kberry(1:3,:)=0 2386 dtsets(idtset)%kpt(:,:)=zero 2387 dtsets(idtset)%kptgw(:,:)=zero 2388 dtsets(idtset)%kptnrm=one 2389 dtsets(idtset)%kptns_hf(:,:)=zero 2390 dtsets(idtset)%kptopt=1 2391 if(dtsets(idtset)%nspden==4)dtsets(idtset)%kptopt=4 2392 dtsets(idtset)%kptrlen=30.0_dp 2393 ! L 2394 2395 #if defined HAVE_LOTF 2396 dtsets(idtset)%lotf_classic=5 2397 dtsets(idtset)%lotf_nitex=10 2398 dtsets(idtset)%lotf_nneigx=40 2399 dtsets(idtset)%lotf_version=2 2400 #endif 2401 dtsets(idtset)%lambsig(:) = zero 2402 dtsets(idtset)%lw_qdrpl=0 2403 dtsets(idtset)%lw_flexo=0 2404 dtsets(idtset)%lw_natopt=0 2405 ! M 2406 dtsets(idtset)%magconon = 0 2407 dtsets(idtset)%magcon_lambda = 0.01_dp 2408 dtsets(idtset)%mband = -1 2409 dtsets(idtset)%mdtemp(:)=300.0_dp 2410 dtsets(idtset)%mdwall=10000_dp 2411 dtsets(idtset)%mep_mxstep=100._dp 2412 dtsets(idtset)%mep_solver=0 2413 dtsets(idtset)%mffmem=1 2414 dtsets(idtset)%mgfft = -1 2415 dtsets(idtset)%mgfftdg = -1 2416 dtsets(idtset)%mixesimgf(:)=zero 2417 dtsets(idtset)%mpw = -1 2418 dtsets(idtset)%mqgrid=0 2419 dtsets(idtset)%mqgriddg=0 2420 ! N 2421 dtsets(idtset)%natrd = -1 2422 dtsets(idtset)%nband(:)=0 2423 dtsets(idtset)%nbandhf=0 2424 dtsets(idtset)%nbdblock=1 2425 dtsets(idtset)%nbdbuf=0 2426 dtsets(idtset)%nberry=1 2427 if (dtsets(idtset)%usepaw == 0) then 2428 dtsets(idtset)%nc_xccc_gspace = 0 2429 else 2430 dtsets(idtset)%nc_xccc_gspace = 1 2431 end if 2432 dtsets(idtset)%nctime=0 2433 dtsets(idtset)%ndtset = -1 2434 dtsets(idtset)%neb_algo=1 2435 dtsets(idtset)%neb_spring(1:2)=(/0.05_dp,0.05_dp/) 2436 dtsets(idtset)%nfft = -1 2437 dtsets(idtset)%nfftdg = -1 2438 2439 dtsets(idtset)%npulayit=7 2440 2441 ! ngfft is a special case 2442 dtsets(idtset)%ngfft(1:8)=0 2443 dtsets(idtset)%ngfft(7) = fftalg_for_npfft(1) 2444 ! fftcache=ngfft(8) is machine-dependent. 2445 dtsets(idtset)%ngfft(8) = get_cache_kb() 2446 2447 dtsets(idtset)%ngfftdg(:)=dtsets(idtset)%ngfft(:) 2448 ! 2449 !nline 2450 dtsets(idtset)%nline=4 2451 2452 if(dtsets(idtset)%usewvl==1 .and. .not. wvl_bigdft) then 2453 if(dtsets(idtset)%usepaw==1) then 2454 dtsets(idtset)%nline=4 2455 else 2456 dtsets(idtset)%nline=2 2457 end if 2458 end if 2459 2460 ! nloalg is also a special case 2461 dtsets(idtset)%nloalg(1)=4 2462 dtsets(idtset)%nloalg(2)=1 2463 dtsets(idtset)%nloalg(3)=dtsets(idtset)%usepaw 2464 !if (dtsets(idtset)%optdriver == RUNL_EPH) dtsets(idtset)%nloalg(3) = 1 2465 dtsets(idtset)%ngkpt=0 2466 dtsets(idtset)%nnsclo=0 2467 dtsets(idtset)%nnsclohf=0 2468 dtsets(idtset)%nonlinear_info=0 2469 dtsets(idtset)%noseinert=1.0d5 2470 dtsets(idtset)%npvel=0 2471 dtsets(idtset)%nqpt=0 2472 dtsets(idtset)%nscforder=16 2473 dtsets(idtset)%nshiftk=1 2474 dtsets(idtset)%nshiftk_orig=1 2475 dtsets(idtset)%nstep=30 2476 dtsets(idtset)%ntime=1 2477 dtsets(idtset)%nwfshist=0 2478 if(dtsets(idtset)%usewvl==1 .and. .not. wvl_bigdft) then 2479 if(dtsets(idtset)%usepaw==1) then 2480 dtsets(idtset)%nwfshist=4 2481 else 2482 dtsets(idtset)%nwfshist=2 2483 end if 2484 end if 2485 ! O 2486 dtsets(idtset)%occopt=1 2487 dtsets(idtset)%occ_orig(:,:)=zero 2488 dtsets(idtset)%optcell=0 2489 dtsets(idtset)%optforces=2 2490 if(dtsets(idtset)%usedmft>0) dtsets(idtset)%optforces=0 2491 dtsets(idtset)%optstress=1 2492 dtsets(idtset)%optnlxccc=1 2493 dtsets(idtset)%orbmag=0 2494 if (dtsets(idtset)%usepaw==0) then 2495 dtsets(idtset)%ortalg=2 2496 else 2497 dtsets(idtset)%ortalg=-2 2498 end if 2499 ! P 2500 dtsets(idtset)%paral_atom=paral_atom_default 2501 dtsets(idtset)%pawcpxocc=1 2502 dtsets(idtset)%pawcross=0 2503 dtsets(idtset)%pawecutdg=-one 2504 dtsets(idtset)%pawfatbnd=0 2505 dtsets(idtset)%pawlcutd=10 2506 dtsets(idtset)%pawlmix=10 2507 dtsets(idtset)%pawmixdg=0 ! Will be set to 1 when npfft>1 2508 dtsets(idtset)%pawnhatxc=1 2509 dtsets(idtset)%pawntheta=12 2510 dtsets(idtset)%pawnphi=13 2511 dtsets(idtset)%pawnzlm=1 2512 dtsets(idtset)%pawoptmix=0 2513 dtsets(idtset)%pawoptosc=0 2514 dtsets(idtset)%pawovlp=5._dp 2515 dtsets(idtset)%pawprtdos=0 2516 dtsets(idtset)%pawprtvol=0 2517 dtsets(idtset)%pawprtwf=0 2518 dtsets(idtset)%pawprt_k=0 2519 dtsets(idtset)%pawprt_b=0 2520 dtsets(idtset)%pawstgylm=1 2521 dtsets(idtset)%pawsushat=0 2522 dtsets(idtset)%pawujat=1 2523 dtsets(idtset)%pawujrad=20.0_dp 2524 dtsets(idtset)%pawujv=0.1_dp/Ha_eV 2525 dtsets(idtset)%pawusecp=1 2526 dtsets(idtset)%pawxcdev=1 2527 dtsets(idtset)%pimd_constraint=0 2528 dtsets(idtset)%pitransform=0 2529 dtsets(idtset)%ptcharge(:) = zero 2530 !dtsets(idtset)%plowan_compute=0 2531 dtsets(idtset)%plowan_bandi=0 2532 dtsets(idtset)%plowan_bandf=0 2533 !if(dtsets(idtset)%plowan_compute>0) then 2534 dtsets(idtset)%plowan_it(:)=0 2535 dtsets(idtset)%plowan_iatom(:)=0 2536 dtsets(idtset)%plowan_lcalc(:)=-1 2537 dtsets(idtset)%plowan_projcalc(:)=0 2538 dtsets(idtset)%plowan_nbl(:)=0 2539 !end if 2540 dtsets(idtset)%plowan_natom=0 2541 dtsets(idtset)%plowan_nt=0 2542 dtsets(idtset)%plowan_realspace=0 2543 dtsets(idtset)%pol(:)=zero 2544 dtsets(idtset)%polcen(:)=zero 2545 dtsets(idtset)%posdoppler=0 2546 dtsets(idtset)%positron=0 2547 dtsets(idtset)%posnstep=50 2548 dtsets(idtset)%posocc=one 2549 dtsets(idtset)%postoldfe=0.000001_dp 2550 dtsets(idtset)%postoldff=zero 2551 dtsets(idtset)%prepalw=0 2552 dtsets(idtset)%prepanl=0 2553 dtsets(idtset)%prtden=1 ; if (dtsets(idtset)%nimage>1) dtsets(idtset)%prtden=0 2554 dtsets(idtset)%prtebands=1 ; if (dtsets(idtset)%nimage>1) dtsets(idtset)%prtebands=0 2555 dtsets(idtset)%prteig=1 ; if (dtsets(idtset)%nimage>1) dtsets(idtset)%prteig=0 2556 dtsets(idtset)%prtgsr=1 ; if (dtsets(idtset)%nimage>1) dtsets(idtset)%prtgsr=0 2557 dtsets(idtset)%prtkpt = -1 2558 dtsets(idtset)%prtocc=0 2559 dtsets(idtset)%prtwf=1 ; if (dtsets(idtset)%nimage>1) dtsets(idtset)%prtwf=0 2560 !if (dtsets%(idtset)%optdriver == RUNL_RESPFN and all(dtsets(:)%optdriver /= RUNL_NONLINEAR) dtsets(idtset)%prtwf = -1 2561 do ii=1,dtsets(idtset)%natom,1 2562 dtsets(idtset)%prtatlist(ii)=ii 2563 end do 2564 dtsets(idtset)%pvelmax(:)=one 2565 dtsets(idtset)%pw_unbal_thresh=40._dp 2566 ! Q 2567 dtsets(idtset)%qmass(:)=ten 2568 dtsets(idtset)%qprtrb(1:3)=0 2569 dtsets(idtset)%qptdm(:,:)=zero 2570 dtsets(idtset)%quadmom(:) = zero 2571 ! R 2572 dtsets(idtset)%random_atpos=0 2573 dtsets(idtset)%ratsm=zero 2574 if (any(dtsets(idtset)%constraint_kind(1:dtsets(idtset)%ntypat)>0)) dtsets(idtset)%ratsm=0.05_dp 2575 dtsets(idtset)%ratsph_extra=two 2576 dtsets(idtset)%recefermi=zero 2577 dtsets(idtset)%recgratio=1 2578 dtsets(idtset)%recnpath=500 2579 dtsets(idtset)%recnrec=10 2580 dtsets(idtset)%recrcut=zero 2581 dtsets(idtset)%recptrott=0 2582 dtsets(idtset)%rectesteg=0 2583 dtsets(idtset)%rectolden=zero 2584 dtsets(idtset)%rcut=zero 2585 dtsets(idtset)%restartxf=0 2586 dtsets(idtset)%rfasr=0 2587 dtsets(idtset)%rfatpol(1:2)=-1 2588 dtsets(idtset)%rfddk=0 2589 dtsets(idtset)%rfdir(1:3)=1 2590 dtsets(idtset)%rfelfd=0 2591 dtsets(idtset)%rfmagn=0 2592 dtsets(idtset)%rfmeth=1 2593 dtsets(idtset)%rfphon=0 2594 dtsets(idtset)%rfstrs=0 2595 dtsets(idtset)%rfstrs_ref=0 2596 dtsets(idtset)%rfuser=0 2597 dtsets(idtset)%rf2_dkdk=0 2598 dtsets(idtset)%rf2_dkde=0 2599 dtsets(idtset)%rf2_pert1_dir(1:3)=1 2600 dtsets(idtset)%rf2_pert2_dir(1:3)=1 2601 dtsets(idtset)%rhoqpmix=one 2602 ! S 2603 dtsets(idtset)%shiftk_orig(:,:)=one 2604 dtsets(idtset)%signperm=1 2605 dtsets(idtset)%slabwsrad=zero 2606 dtsets(idtset)%slk_rankpp=1000 2607 dtsets(idtset)%smdelta=0 2608 dtsets(idtset)%spbroad=0.1_dp 2609 dtsets(idtset)%spgaxor = -1 2610 dtsets(idtset)%spgorig = -1 2611 dtsets(idtset)%spinmagntarget=-99.99_dp 2612 dtsets(idtset)%spnorbscl=one 2613 dtsets(idtset)%stmbias=zero 2614 dtsets(idtset)%strfact=100.0_dp 2615 dtsets(idtset)%string_algo=1 2616 dtsets(idtset)%strprecon=one 2617 dtsets(idtset)%strtarget(1:6)=zero 2618 ! T 2619 dtsets(idtset)%td_exp_order=4 2620 dtsets(idtset)%td_maxene=zero 2621 dtsets(idtset)%td_mexcit=0 2622 dtsets(idtset)%td_scnmax=3 2623 dtsets(idtset)%td_prtstr=10 2624 dtsets(idtset)%td_restart=0 2625 dtsets(idtset)%td_propagator=1 2626 dtsets(idtset)%td_scthr=1e-7_dp 2627 dtsets(idtset)%tfw_toldfe=0.000001_dp 2628 dtsets(idtset)%tim1rev = 1 2629 dtsets(idtset)%tl_nprccg = 30 2630 dtsets(idtset)%tl_radius = zero 2631 dtsets(idtset)%tphysel=zero 2632 dtsets(idtset)%toldfe=zero 2633 dtsets(idtset)%tolmxde=zero 2634 dtsets(idtset)%toldff=zero 2635 dtsets(idtset)%tolimg=5.0d-5 2636 dtsets(idtset)%tolrde=0.005_dp 2637 dtsets(idtset)%tolrff=zero 2638 dtsets(idtset)%tolmxf=5.0d-5 2639 dtsets(idtset)%tolvrs=zero 2640 dtsets(idtset)%tolwfr=zero 2641 2642 dtsets(idtset)%tsmear=0.01_dp 2643 ! U 2644 dtsets(idtset)%ucrpa_bands(:)=-1 2645 dtsets(idtset)%ucrpa_window(:)=-1.0_dp 2646 dtsets(idtset)%upawu(:,:)=zero 2647 dtsets(idtset)%usepead=1 2648 dtsets(idtset)%usefock=0 2649 dtsets(idtset)%usekden=0 2650 dtsets(idtset)%use_gemm_nonlop=0 2651 dtsets(idtset)%use_nonscf_gkk=0 !1 ! deactivate by default, for now 6 Oct 2013 2652 dtsets(idtset)%userec=0 2653 dtsets(idtset)%usexcnhat_orig=-1 2654 dtsets(idtset)%useylm=0 2655 ! V 2656 dtsets(idtset)%vacnum = -1 2657 dtsets(idtset)%vcutgeo(3)=zero 2658 dtsets(idtset)%vdw_nfrag = 1 2659 #if defined DEV_YP_VDWXC 2660 dtsets(idtset)%vdw_df_acutmin = vdw_defaults%acutmin 2661 dtsets(idtset)%vdw_df_aratio = vdw_defaults%aratio 2662 dtsets(idtset)%vdw_df_damax = vdw_defaults%damax 2663 dtsets(idtset)%vdw_df_damin = vdw_defaults%damin 2664 dtsets(idtset)%vdw_df_dcut = vdw_defaults%dcut 2665 dtsets(idtset)%vdw_df_dratio = vdw_defaults%dratio 2666 dtsets(idtset)%vdw_df_dsoft = vdw_defaults%dsoft 2667 dtsets(idtset)%vdw_df_gcut = vdw_defaults%gcut 2668 dtsets(idtset)%vdw_df_ndpts = vdw_defaults%ndpts 2669 dtsets(idtset)%vdw_df_ngpts = vdw_defaults%ngpts 2670 dtsets(idtset)%vdw_df_nqpts = vdw_defaults%nqpts 2671 dtsets(idtset)%vdw_df_nrpts = vdw_defaults%nrpts 2672 dtsets(idtset)%vdw_df_nsmooth = vdw_defaults%nsmooth 2673 dtsets(idtset)%vdw_df_phisoft = vdw_defaults%phisoft 2674 dtsets(idtset)%vdw_df_qcut = vdw_defaults%qcut 2675 dtsets(idtset)%vdw_df_qratio = vdw_defaults%qratio 2676 dtsets(idtset)%vdw_df_rcut = vdw_defaults%rcut 2677 dtsets(idtset)%vdw_df_rsoft = vdw_defaults%rsoft 2678 dtsets(idtset)%vdw_df_tolerance = vdw_defaults%tolerance 2679 dtsets(idtset)%vdw_df_tweaks = vdw_defaults%tweaks 2680 dtsets(idtset)%vdw_df_zab = vdw_defaults%zab 2681 dtsets(idtset)%vdw_df_threshold = 1.0d-2 2682 #endif 2683 dtsets(idtset)%vdw_supercell(:) = 0 2684 dtsets(idtset)%vdw_tol = tol10 2685 dtsets(idtset)%vdw_tol_3bt = -1 2686 dtsets(idtset)%vdw_typfrag(:) = 1 2687 dtsets(idtset)%vdw_xc = 0 2688 dtsets(idtset)%vis=100.0_dp 2689 dtsets(idtset)%vprtrb(1:2)=zero 2690 ! W 2691 dtsets(idtset)%wtatcon(:,:,:)=zero 2692 dtsets(idtset)%wfmix=one 2693 dtsets(idtset)%wfk_task=0 2694 dtsets(idtset)%wtk=one 2695 dtsets(idtset)%wvl_crmult = 6._dp 2696 dtsets(idtset)%wvl_frmult = 10._dp 2697 dtsets(idtset)%wvl_hgrid = 0.5_dp 2698 dtsets(idtset)%wvl_ngauss =(/1,100/) 2699 dtsets(idtset)%wvl_nprccg = 10 2700 dtsets(idtset)%w90iniprj = 1 2701 dtsets(idtset)%w90prtunk = 0 2702 2703 ! X 2704 dtsets(idtset)%xclevel = 0 2705 dtsets(idtset)%xc_denpos = tol14 2706 dtsets(idtset)%xc_tb09_c = 99.99_dp 2707 dtsets(idtset)%xredsph_extra(:,:)=zero 2708 ! Y 2709 ! Z 2710 dtsets(idtset)%zcut=3.67493260d-03 ! = 0.1eV 2711 if(dtsets(idtset)%optdriver == RUNL_GWLS) dtsets(idtset)%zcut=zero 2712 !if(dtsets(idtset)%optdriver == RUNL_EPH) dtsets(idtset)%zcut = 0.01 * eV_Ha 2713 dtsets(idtset)%ziontypat(:)=zero 2714 2715 dtsets(idtset)%bs_loband=0 2716 2717 !if (dtsets(idtset)%optdriver == RUNL_EPH) then 2718 ! dtsets(idtset)%mixprec = 1 2719 ! dtsets(idtset)%boxcutmin = 1.1_dp 2720 !end if 2721 2722 ! JB:UNINITIALIZED VALUES (not found in this file neither indefo1) 2723 ! They might be initialized somewhereelse, I don't know. 2724 ! That might cause unitialized error with valgrind depending on the compiler 2725 ! chkprim 2726 ! maxnsym 2727 ! nsym 2728 ! macro_uj 2729 ! prtpmp 2730 ! timopt 2731 ! useria 2732 ! userib 2733 ! useric 2734 ! userid 2735 ! userie 2736 ! bravais 2737 ! symafm 2738 ! symrel 2739 ! fband 2740 ! nelect 2741 ! userra 2742 ! userrb 2743 ! userrc 2744 ! userrd 2745 ! userre 2746 ! vacwidth 2747 ! genafm 2748 ! kptns 2749 ! rprimd_orig 2750 ! tnons 2751 2752 end do 2753 2754 DBG_EXIT("COLL") 2755 2756 end subroutine indefo
ABINIT/indefo1 [ 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.
SOURCE
918 subroutine indefo1(dtset) 919 920 !Arguments ------------------------------------ 921 !scalars 922 type(dataset_type),intent(inout) :: dtset 923 924 !Local variables ------------------------------- 925 !scalars 926 !integer :: ii 927 928 !****************************************************************** 929 930 !Set up default values. All variables to be output in outvars.f 931 !should have a default, even if a nonsensible one can be chosen to garantee print in that routine. 932 933 DBG_ENTER("COLL") 934 935 !Use alphabetic order 936 937 !A 938 dtset%acell_orig(:,:)=zero 939 dtset%algalch(:)=1 940 dtset%amu_orig(:,:)=-one 941 dtset%autoparal=0 942 !B 943 dtset%bandpp=1 944 dtset%berryopt=0 945 dtset%berrysav=0 946 dtset%bfield(:)=zero 947 !C 948 dtset%cd_customnimfrqs=0 949 dtset%chkprim=1 950 dtset%chrgat(:)=zero 951 dtset%constraint_kind(:)=0 952 !D 953 dtset%densty(:,:)=zero 954 dtset%dfield(:)=zero !!HONG 955 dtset%dynimage(:)=1 956 !E 957 dtset%efield(:)=zero 958 dtset%efmas_calc_dirs=0 959 dtset%efmas_n_dirs=0 960 !F 961 !G 962 dtset%ga_n_rules=1 963 dtset%gw_customnfreqsp=0 964 dtset%gw_nqlwl=0 965 dtset%gwls_n_proj_freq=0 966 !I 967 dtset%iatfix(:,:)=0 968 dtset%icoulomb=0 969 dtset%imgmov=0 970 dtset%ivalence=0 971 !J 972 dtset%jellslab=0 973 dtset%jfielddir(:)=0 974 !K 975 dtset%kptopt=0 976 !L 977 dtset%lexexch(:)=-1 978 dtset%ldaminushalf(:)=0 979 dtset%lpawu(:)=-1 980 !M 981 dtset%maxestep=0.005d0 982 dtset%mixalch_orig(:,:,:)=zero 983 dtset%mkmem=-1 984 dtset%mkqmem=-1 985 dtset%mk1mem=-1 986 !N 987 dtset%natpawu=0 988 dtset%natsph=0 989 dtset%natsph_extra=0 990 dtset%natvshift=0 991 dtset%nconeq=0 992 dtset%ndynimage=1 993 dtset%ne_qFD=zero 994 dtset%nh_qFD=zero 995 dtset%nkpt=-1 996 dtset%nkptgw=0 997 dtset%nkpthf=0 998 dtset%nnos=0 999 dtset%npband=1 1000 dtset%npfft=1 1001 dtset%nphf=1 1002 dtset%npimage=1 1003 dtset%np_spkpt=1 1004 dtset%nppert=1 1005 dtset%npspalch=0 1006 dtset%npspinor=1 1007 dtset%np_slk=1000000 1008 dtset%nqptdm=0 1009 dtset%nspden=1 1010 dtset%nspinor=1 1011 dtset%nsppol=1 1012 dtset%nsym=0 ! Actually, this default value is not used : it is to be reimposed before each call to ingeo in invars1 1013 dtset%ntimimage=1 1014 dtset%ntypalch=0 1015 dtset%ntyppure=-1 1016 dtset%nucdipmom(:,:)=zero 1017 dtset%nzchempot=0 1018 !O 1019 dtset%optdriver=0 1020 !P 1021 dtset%paral_rf=0 1022 !dtset%paral_kgb ! Is even initialized earlier. 1023 dtset%pawspnorb=0 ! will be changed to 1 as soon as usepaw==1 and nspinor==2 1024 dtset%pimass(:)=-one 1025 !Q 1026 dtset%qptn=zero 1027 !R 1028 dtset%red_efield(:)=zero 1029 dtset%red_dfield(:)=zero 1030 dtset%red_efieldbar(:)=zero 1031 dtset%rprim_orig(:,:,:)=zero 1032 dtset%rprim_orig(1,1,:)=one 1033 dtset%rprim_orig(2,2,:)=one 1034 dtset%rprim_orig(3,3,:)=one 1035 !S 1036 dtset%slabzbeg=zero 1037 dtset%slabzend=zero 1038 dtset%so_psp(:)=1 1039 dtset%spinat(:,:)=zero 1040 !T 1041 dtset%tfkinfunc=0 1042 dtset%typat(:)=0 ! This init is important because dimension of typat is mx%natom (and not natom). 1043 !U 1044 dtset%usedmatpu=0 1045 dtset%usedmft=0 1046 dtset%useexexch=0 1047 dtset%usepawu=0 1048 dtset%usepotzero=0 1049 dtset%use_slk=0 1050 dtset%use_oldchi=1 1051 !V 1052 dtset%vel_orig(:,:,:)=zero 1053 dtset%vel_cell_orig(:,:,:)=zero 1054 !W 1055 dtset%wtq=zero 1056 if (dtset%usepaw==0) dtset%wfoptalg=0 1057 if (dtset%usepaw/=0) dtset%wfoptalg=10 1058 if (dtset%optdriver==RUNL_GSTATE.and.dtset%paral_kgb>0) dtset%wfoptalg=14 1059 dtset%wvl_bigdft_comp=1 1060 1061 !X 1062 dtset%xred_orig(:,:,:)=zero 1063 !Y 1064 !Z 1065 dtset%zeemanfield(:)=zero 1066 1067 DBG_EXIT("COLL") 1068 1069 end subroutine indefo1
ABINIT/invars0 [ 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 comm= MPI communicator
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 pseudo_paths(npsp): List of paths to pseudopotential files as read from input file. List of empty strings if we are legacy "files file" mode. Allocated here, caller should free memory.
SOURCE
95 subroutine invars0(dtsets, istatr, istatshft, lenstr, msym, mxnatom, mxnimage, mxntypat, ndtset, ndtset_alloc, & 96 npsp, pseudo_paths, papiopt, timopt, string, comm) 97 98 !Arguments ------------------------------------ 99 !scalars 100 integer,intent(in) :: lenstr,ndtset,ndtset_alloc, comm 101 integer,intent(out) :: istatr,istatshft,msym,mxnatom,mxnimage,mxntypat,npsp,papiopt 102 integer,intent(inout) :: timopt 103 character(len=*),intent(in) :: string 104 !arrays 105 type(dataset_type),intent(inout) :: dtsets(0:ndtset_alloc) !vz_i 106 character(len=fnlen),allocatable,intent(out) :: pseudo_paths(:) 107 108 !Local variables------------------------------- 109 !scalars 110 integer :: i1,i2,idtset,ii,jdtset,marr,multiplicity,tjdtset,tread,treadh,treadm,tread_pseudos,cnt, tread_geo 111 integer :: treads, use_gpu_cuda 112 real(dp) :: cpus 113 character(len=500) :: msg 114 character(len=fnlen) :: pp_dirpath 115 character(len=20*fnlen) :: pseudos_string ! DO NOT decrease len 116 character(len=len(string)) :: geo_string 117 type(geo_t) :: geo 118 !arrays 119 integer,allocatable :: intarr(:), sidx(:) 120 real(dp),allocatable :: dprarr(:) 121 122 !****************************************************************** 123 124 !write(std_out,"(3a)")"invars1 with string:", ch10, trim(string) 125 126 marr=max(9,ndtset_alloc,2) 127 ABI_MALLOC(dprarr,(marr)) 128 ABI_MALLOC(intarr,(marr)) 129 130 ! Set up jdtset 131 if (ndtset/=0) then 132 ! Default values 133 dtsets(0)%jdtset = -1 ! unused value 134 dtsets(1:ndtset_alloc)%jdtset=(/ (ii,ii=1,ndtset_alloc) /) 135 136 ! Read explicitly the jdtset array 137 call intagm(dprarr,intarr,0,marr,ndtset,string(1:lenstr),'jdtset',tjdtset,'INT') 138 if(tjdtset==1) dtsets(1:ndtset)%jdtset=intarr(1:ndtset) 139 140 ! Read the udtset array 141 call intagm(dprarr,intarr,0,marr,2,string(1:lenstr),'udtset',tread,'INT') 142 143 ! jdtset and udtset cannot be defined together 144 if(tjdtset==1 .and. tread==1)then 145 write(msg, '(3a)' )& 146 'jdtset and udtset cannot be defined both in the input file.',ch10,& 147 'Action: remove one of them from your input file.' 148 ABI_ERROR(msg) 149 end if 150 151 ! Check values of udtset 152 if(tread==1)then 153 if(intarr(1)<1 .or. intarr(1)>999)then 154 write(msg, '(a,i0,3a)' )& 155 'udtset(1) must be between 1 and 999, but it is ',intarr(1),'.',ch10,& 156 'Action: change the value of udtset(1) in your input file.' 157 ABI_ERROR(msg) 158 end if 159 if(intarr(2)<1 .or. intarr(2)>9)then 160 write(msg, '(a,i0,3a)' )& 161 'udtset(2) must be between 1 and 9, but it is ',intarr(2),'.',ch10,& 162 'Action: change the value of udtset(2) in your input file.' 163 ABI_ERROR(msg) 164 end if 165 if(intarr(1)*intarr(2) /= ndtset)then 166 write(msg, '(3a,i0,3a,i0,a,i0,3a,i0,3a)' )& 167 'udtset(1)*udtset(2) must be equal to ndtset,',ch10,& 168 'but it is observed that udtset(1) = ',intarr(1),',',ch10,& 169 'and udtset(2) = ',intarr(2),' so that their product is ',intarr(1)*intarr(2),',',ch10,& 170 'while ndtset is ',ndtset,'.',ch10,& 171 'Action: change udtset or ndtset in your input file.' 172 ABI_ERROR(msg) 173 end if 174 idtset=0 175 do i1=1,intarr(1) 176 do i2=1,intarr(2) 177 idtset=idtset+1 178 dtsets(idtset)%jdtset=i1*10+i2 179 end do 180 end do 181 end if 182 183 ! Final check on the jdtset values 184 do idtset=1,ndtset 185 if(dtsets(idtset)%jdtset<1 .or. dtsets(idtset)%jdtset>9999)then 186 write(msg, '(3a,i0,a,i0,a,a)' )& 187 'The components of jdtset must be between 1 and 9999.',ch10,& 188 'However, the input value of the component ',idtset,' of jdtset is ',dtsets(idtset)%jdtset,ch10,& 189 'Action: correct jdtset in your input file.' 190 ABI_ERROR(msg) 191 end if 192 end do 193 194 else 195 dtsets(1)%jdtset=0 196 end if 197 198 papiopt = 0 199 call intagm(dprarr,intarr,0,1,1,string(1:lenstr),'papiopt',tread,'INT') 200 if(tread==1) papiopt=intarr(1) 201 202 ! Read timopt and pass it to timab 203 call intagm(dprarr,intarr,0,1,1,string(1:lenstr),'timopt',tread,'INT') 204 if(tread==1) timopt=intarr(1) 205 206 istatr=0 207 dtsets(0)%istatr=istatr 208 call intagm(dprarr,intarr,0,marr,1,string(1:lenstr),'istatr',tread,'INT') 209 if(tread==1) istatr=intarr(1) 210 dtsets(1:)%istatr=istatr 211 212 istatshft=1 213 dtsets(0)%istatshft=istatshft 214 call intagm(dprarr,intarr,0,marr,1,string(1:lenstr),'istatshft',tread,'INT') 215 if(tread==1) istatshft=intarr(1) 216 dtsets(1:)%istatshft=istatshft 217 218 cpus=zero 219 call intagm(dprarr,intarr,0,marr,1,string(1:lenstr),'cpus ',treads,'DPR') 220 if(treads==1) cpus=dprarr(1) 221 call intagm(dprarr,intarr,0,marr,1,string(1:lenstr),'cpum ',treadm,'DPR') 222 if(treadm==1) cpus=dprarr(1)*60.0_dp 223 call intagm(dprarr,intarr,0,marr,1,string(1:lenstr),'cpuh ',treadh,'DPR') 224 225 if(treadh==1) cpus=dprarr(1)*3600.0_dp 226 if(treads+treadm+treadh>1)then 227 write(msg, '(5a)' )& 228 'More than one input variable is used to defined the CPU time limit.',ch10,& 229 'This is not allowed.',ch10,& 230 'Action: in the input file, suppress either cpus, cpum or cpuh.' 231 ABI_ERROR(msg) 232 end if 233 dtsets(:)%cpus=cpus 234 235 ! Default for natom, nimage, ntypat, useri and userr 236 dtsets(:)%natom=1 237 dtsets(:)%nimage=1 238 dtsets(:)%ntypat=1 ; dtsets(0)%ntypat=0 ! Will always echo ntypat 239 dtsets(:)%macro_uj=0 240 dtsets(:)%maxnsym=384 241 dtsets(:)%useria=0 242 dtsets(:)%userib=0 243 dtsets(:)%useric=0 244 dtsets(:)%userid=0 245 dtsets(:)%userie=0 246 dtsets(:)%userra=zero 247 dtsets(:)%userrb=zero 248 dtsets(:)%userrc=zero 249 dtsets(:)%userrd=zero 250 dtsets(:)%userre=zero 251 dtsets(:)%usewvl = 0 252 dtsets(:)%plowan_compute=0 253 254 ! Loop on datasets, to find natom and mxnatom, as well as useri and userr 255 do idtset=1,ndtset_alloc 256 jdtset=dtsets(idtset)%jdtset ; if(ndtset==0)jdtset=0 257 258 ! proposal: supercell generation in input string before it is read in 259 ! call expand_supercell_input(jdtset, lenstr, string) 260 ! find supercell, else exit 261 ! determinant = ncells 262 ! copy rprim, acell, xred, xcart, vel, typat, to 263 ! rprim_uc, acell_uc, xred_uc, xcart_uc, vel_uc, typat_uc 264 ! NB: also rprim and angdeg need to be updated in non diagonal case!!! 265 ! generate supercell info for each of these copying out with translation vectors etc... 266 ! set chkprim to 0 267 ! done! 268 269 ! Generate the supercell if supercell_latt is specified and update string 270 dtsets(idtset)%supercell_latt(:) = 0 271 do ii=1,3 272 dtsets(idtset)%supercell_latt(ii) = 1 273 end do 274 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),"supercell_latt",tread,'INT') 275 if (tread==1) dtsets(idtset)%supercell_latt(:)=intarr(1:3) 276 !This test should be update if in the future we allow non-diagonal supercell 277 if (any(dtsets(idtset)%supercell_latt(:) < tol10 )) then 278 write(msg, '(5a)' )& 279 'supercell_latt must have positive parameters and diagonal part',ch10,& 280 'This is not allowed. ',ch10,& 281 'Action: modify supercell_latt in the input file.' 282 ABI_ERROR(msg) 283 end if 284 ! Compute the multiplicity of the supercell 285 multiplicity=dtsets(idtset)%supercell_latt(1) & 286 & *dtsets(idtset)%supercell_latt(2) & 287 & *dtsets(idtset)%supercell_latt(3) 288 ! call mati3det(dtsets(idtset)%supercell_latt,multiplicity) 289 290 ! Read natom from string 291 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natom',tread,'INT') 292 293 ! or get it from the structure variable 294 call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'structure', tread_geo, & 295 'KEY', key_value=geo_string) 296 297 if (tread_geo /= 0) then 298 geo = geo_from_abivar_string(geo_string, comm) 299 if (tread /= 0) then 300 ABI_CHECK(intarr(1) == geo%natom, "natom from variable and from structure do not agree with each other") 301 end if 302 intarr(1) = geo%natom 303 tread = 1 304 end if 305 306 ! Might also initialize natom from XYZ file 307 if (tread==0) call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'_natom',tread,'INT') 308 309 if (tread==1) then 310 dtsets(idtset)%natom=intarr(1) 311 else 312 write(msg, '(a,i0,2a)' )& 313 'Input natom must be defined, but was absent for dataset ',jdtset,ch10,& 314 'Action: check the input file.' 315 ABI_ERROR(msg) 316 end if 317 318 ! Check that natom is greater than 0 319 if (dtsets(idtset)%natom<=0) then 320 write(msg, '(a,i0,2a,i0,3a)' )& 321 'Input natom must be > 0, but was ',dtsets(idtset)%natom,ch10,& 322 'for dataset ',jdtset,'. This is not allowed.',ch10,& 323 'Action: check the input file.' 324 ABI_ERROR(msg) 325 end if 326 327 if(multiplicity > 1)then 328 dtsets(idtset)%natom = dtsets(idtset)%natom * multiplicity 329 end if 330 331 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nimage',tread,'INT') 332 if(tread==1) dtsets(idtset)%nimage=intarr(1) 333 334 ! Check that nimage is greater than 0 335 if (dtsets(idtset)%nimage<=0) then 336 write(msg, '(a,i0,4a)' )& 337 'nimage must be > 0, but was ',dtsets(idtset)%nimage,ch10,& 338 'This is not allowed.',ch10,& 339 'Action: check the input file.' 340 ABI_ERROR(msg) 341 end if 342 343 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ntypat',tread,'INT') 344 if (tread==1) dtsets(idtset)%ntypat=intarr(1) 345 346 if (tread_geo /= 0) then 347 if (tread == 1) then 348 ABI_CHECK(geo%ntypat == dtsets(idtset)%ntypat, "ntypat and geo%ntypat do not agree with each other") 349 end if 350 dtsets(idtset)%ntypat = geo%ntypat 351 end if 352 353 ! Check that ntypat is greater than 0 354 if (dtsets(idtset)%ntypat<=0) then 355 write(msg, '(a,i0,2a,i0,3a)' )& 356 'Input ntypat must be > 0, but was ',dtsets(idtset)%ntypat,ch10,& 357 'for dataset ',jdtset,'. This is not allowed.',ch10,& 358 'Action: check the input file.' 359 ABI_ERROR(msg) 360 end if 361 362 ! Read msym from string 363 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'maxnsym',tread,'INT') 364 if(tread==1)dtsets(idtset)%maxnsym=intarr(1) 365 ! Check that maxnsym is greater than 1 366 if (dtsets(idtset)%maxnsym<1) then 367 write(msg, '(a,i0,2a,i0,3a)' )& 368 'Input maxnsym must be > 1, but was ',dtsets(idtset)%maxnsym,ch10,& 369 'for dataset ',jdtset,'. This is not allowed.',ch10,& 370 'Action: check the input file.' 371 ABI_ERROR(msg) 372 end if 373 374 ! Read plowan_compute 375 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'plowan_compute',tread,'INT') 376 if(tread==1) dtsets(idtset)%plowan_compute=intarr(1) 377 378 ! Read extfpmd calculations 379 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'useextfpmd',tread,'INT') 380 if(tread==1) dtsets(idtset)%useextfpmd=intarr(1) 381 382 ! Read user* variables 383 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'useria',tread,'INT') 384 if(tread==1) dtsets(idtset)%useria=intarr(1) 385 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userib',tread,'INT') 386 if(tread==1) dtsets(idtset)%userib=intarr(1) 387 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'useric',tread,'INT') 388 if(tread==1) dtsets(idtset)%useric=intarr(1) 389 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userid',tread,'INT') 390 if(tread==1) dtsets(idtset)%userid=intarr(1) 391 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userie',tread,'INT') 392 if(tread==1) dtsets(idtset)%userie=intarr(1) 393 394 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userra',tread,'DPR') 395 if(tread==1) dtsets(idtset)%userra=dprarr(1) 396 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userrb',tread,'DPR') 397 if(tread==1) dtsets(idtset)%userrb=dprarr(1) 398 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userrc',tread,'DPR') 399 if(tread==1) dtsets(idtset)%userrc=dprarr(1) 400 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userrd',tread,'DPR') 401 if(tread==1) dtsets(idtset)%userrd=dprarr(1) 402 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'userre',tread,'DPR') 403 if(tread==1) dtsets(idtset)%userre=dprarr(1) 404 405 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usewvl',tread,'INT') 406 if(tread==1) dtsets(idtset)%usewvl=intarr(1) 407 408 call geo%free() 409 end do ! idtset 410 411 !mxnatom =maxval(dtsets(1:ndtset_alloc)%natom) 412 !mxntypat =maxval(dtsets(1:ndtset_alloc)%ntypat) 413 !msym =maxval(dtsets(1:ndtset_alloc)%maxnsym) 414 !There is a bug in the HP compiler, the following should execute properly 415 mxnatom=dtsets(1)%natom ; mxnimage=dtsets(1)%nimage 416 mxntypat=dtsets(1)%ntypat ; msym=dtsets(1)%maxnsym 417 if(ndtset_alloc>1)then 418 do idtset=2,ndtset_alloc 419 mxnatom =max(dtsets(idtset)%natom,mxnatom) 420 mxnimage=max(dtsets(idtset)%nimage,mxnimage) 421 mxntypat=max(dtsets(idtset)%ntypat,mxntypat) 422 msym =max(dtsets(idtset)%maxnsym,msym) 423 end do 424 end if 425 426 if(mxnimage>1)then 427 do idtset=2,ndtset_alloc 428 if(mxnatom/=dtsets(idtset)%natom)then 429 write(msg,'(5a,i0,a,i0,3a,i0,a)')& 430 'When there exist one dataset with more than one image,',ch10,& 431 'the number of atoms in each dataset must be the same.',ch10,& 432 'However, it has been found that for dataset= ',idtset,ch10,& 433 'natom= ',dtsets(idtset)%natom,' differs from the maximum number',ch10,& 434 'of atoms, mxnatom= ',mxnatom,& 435 'Action: check the input variables natom for different datasets.' 436 ABI_ERROR(msg) 437 end if 438 end do 439 end if 440 441 ! Set up npsp 442 npsp=mxntypat ! Default value 443 call intagm(dprarr,intarr,0,marr,1,string(1:lenstr),'npsp',tread,'INT') 444 445 if(tread==1)then 446 npsp=intarr(1) 447 else 448 if(ndtset_alloc>1)then 449 do idtset=1,ndtset_alloc 450 if(dtsets(idtset)%ntypat/=mxntypat)then 451 write(msg, '(5a,i0,a,i0,2a,i0,2a)' )& 452 ' When npsp is not defined, the input variable ntypat must be',ch10,& 453 ' the same for all datasets. However, it has been found that for',ch10,& 454 ' jdtset: ',dtsets(idtset)%jdtset,', ntypat= ',dtsets(idtset)%ntypat,ch10,& 455 ' differs from the maximum value of ntypat= ',mxntypat,ch10,& 456 ' Action: check the input variables npsp and ntypat.' 457 ABI_ERROR(msg) 458 end if 459 if(dtsets(idtset)%ntypat>npsp)then 460 write(msg, '(5a,i0,a,i0,a,i0,2a)' )& 461 ' The number of pseudopotentials, npsp, must never be smaller than ntypat.',ch10,& 462 ' However, it has been found that for',ch10,& 463 ' jdtset: ',dtsets(idtset)%jdtset,', ntypat= ',dtsets(idtset)%ntypat,' and npsp=',npsp,ch10,& 464 ' Action: check the input variables npsp and ntypat.' 465 ABI_ERROR(msg) 466 endif 467 end do 468 end if 469 end if 470 dtsets(0)%npsp = mxntypat ! Default value 471 dtsets(1:ndtset_alloc)%npsp = npsp 472 473 ! Read pseudopotential directory and pseudo paths from input. 474 ! Remember that in "files file mode", this info is passed through the files file so these variables are optional 475 pp_dirpath = "" 476 call intagm(dprarr, intarr, 0, marr, 1, string(1:lenstr), 'pp_dirpath', tread, 'KEY', key_value=pp_dirpath) 477 if (tread == 1) then 478 !! XG2020_07_20 Now, the replacement of environment variables is done at the level of the parser 479 ! if (pp_dirpath(1:1) == "$") then 480 ! shell_var = pp_dirpath(2:) 481 ! call get_environment_variable(shell_var, pp_dirpath, status=ierr) 482 ! if (ierr == -1) ABI_ERROR(sjoin(shell_var, "is present but string too short for the environment variable")) 483 ! if (ierr == +1) ABI_ERROR(sjoin(shell_var, "variable is not defined!")) 484 ! if (ierr == +2) ABI_ERROR(sjoin(shell_var, "used in input file but processor does not support environment variables")) 485 ! call wrtout(std_out, sjoin(shell_var, "found in env. Assuming pseudos located in:", pp_dirpath)) 486 ! end if 487 if (.not. endswith(pp_dirpath, "/")) pp_dirpath = strcat(pp_dirpath, "/") 488 end if 489 490 ! String must be large enough to contain ntypat filepaths. 491 pseudos_string = "" 492 call intagm(dprarr, intarr, 0, marr, 1, string(1:lenstr), "pseudos", tread_pseudos, 'KEY', key_value=pseudos_string) 493 494 ABI_MALLOC(pseudo_paths, (npsp)) 495 pseudo_paths = "" 496 497 if (tread_pseudos == 1) then 498 ! Split pseudos_string using comma and transfer results to pseudos_paths 499 ! Make sure string length is large enough and input string is consistent with npsp 500 ! Lot of checks must be done here! 501 !print *, "pseudos_string: ", trim(pseudos_string) 502 ABI_ICALLOC(sidx, (npsp + 1)) 503 sidx(1) = 1; sidx(npsp + 1) = len(pseudos_string) 504 cnt = 1 505 do ii=1,len(pseudos_string) 506 if (pseudos_string(ii:ii) == ",") then 507 pseudos_string(ii:ii) = " " 508 cnt = cnt + 1 509 sidx(cnt) = ii 510 ABI_CHECK(cnt <= npsp, "Too many commas in pseudos string!") 511 end if 512 end do 513 if (cnt /= npsp) then 514 write(msg,'(4a)')& 515 & "Not enough pseudopotentials in input `pseudos` string, expecting npsp: ",itoa(npsp),ch10,& 516 & "Perhaps the separator (=a comma) is missing between pseudopotentials in input `pseudos` string." 517 ABI_ERROR(msg) 518 end if 519 520 do ii=1,npsp 521 i1 = sidx(ii) 522 i2 = sidx(ii + 1) 523 cnt = len(adjustl(trim(pseudos_string(i1:i2)))) 524 ABI_CHECK(cnt <= fnlen, "pseudo path too small, increase fnlen") 525 pseudo_paths(ii) = adjustl(trim(pseudos_string(i1:i2))) 526 if (len_trim(pp_dirpath) > 0) then 527 if (len_trim(pp_dirpath) + len_trim(pseudo_paths(ii)) > fnlen) then 528 ABI_ERROR(sjoin("String of len fnlen:", itoa(fnlen), " too small to contain full pseudo path")) 529 end if 530 pseudo_paths(ii) = strcat(pp_dirpath, pseudo_paths(ii)) 531 end if 532 end do 533 ABI_FREE(sidx) 534 !print *, "pp_dirpath: ", trim(pp_dirpath), "pseudos: ", trim(pseudos_string) 535 end if 536 537 ! KGB parallelism information (needed at this stage) 538 dtsets(:)%paral_kgb=0 539 do idtset=1,ndtset_alloc 540 jdtset=dtsets(idtset)%jdtset ; if(ndtset==0)jdtset=0 541 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'paral_kgb',tread,'INT') 542 if(tread==1)dtsets(idtset)%paral_kgb=intarr(1) 543 544 if (dtsets(idtset)%paral_kgb<0 .or. dtsets(idtset)%paral_kgb>1) then 545 write(msg,'(a,i0,2a,i0,3a)')& 546 'Input paral_kgb must be 0 or 1, but was ',dtsets(idtset)%paral_kgb,ch10,& 547 'for dataset ',jdtset,'. This is not allowed.',ch10,& 548 'Action: check the input file.' 549 ABI_ERROR(msg) 550 end if 551 552 553 end do 554 555 dtsets(:)%diago_apply_block_sliced=1 556 do idtset=1,ndtset_alloc 557 jdtset=dtsets(idtset)%jdtset ; if(ndtset==0)jdtset=0 558 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'diago_apply_block_sliced',tread,'INT') 559 if(tread==1)dtsets(idtset)%diago_apply_block_sliced=intarr(1) 560 end do 561 562 563 !GPU information 564 use_gpu_cuda=0 565 dtsets(:)%use_gpu_cuda=0 566 #if defined HAVE_GPU_CUDA && defined HAVE_GPU_CUDA_DP 567 call Get_ndevice(ii) 568 if (ii>0) then 569 do i1=1,ndtset_alloc 570 dtsets(i1)%use_gpu_cuda=-1 571 end do 572 end if 573 #endif 574 do idtset=1,ndtset_alloc 575 jdtset=dtsets(idtset)%jdtset ; if(ndtset==0)jdtset=0 576 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'use_gpu_cuda',tread,'INT') 577 if(tread==1)dtsets(idtset)%use_gpu_cuda=intarr(1) 578 if (dtsets(idtset)%use_gpu_cuda==1) use_gpu_cuda=1 579 end do 580 581 dtsets(:)%use_nvtx=0 582 #if defined HAVE_GPU_CUDA && defined HAVE_GPU_NVTX_V3 583 do idtset=1,ndtset_alloc 584 jdtset=dtsets(idtset)%jdtset ; if(ndtset==0)jdtset=0 585 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'use_nvtx',tread,'INT') 586 if(tread==1)dtsets(idtset)%use_nvtx=intarr(1) 587 end do 588 #endif 589 590 if (use_gpu_cuda==1) then 591 #if defined HAVE_GPU_CUDA && defined HAVE_GPU_CUDA_DP 592 if (ii<=0) then 593 write(msg,'(5a)')& 594 & 'Input variables use_gpu_cuda is on',ch10,& 595 & 'but no available GPU device has been detected !',ch10,& 596 & 'Action: change the input variable use_gpu_cuda.' 597 ABI_ERROR(msg) 598 end if 599 #else 600 write(msg,'(7a)')& 601 & 'Input variables use_gpu_cuda is on but abinit hasn''t been built',ch10,& 602 & 'with (double precision) gpu mode enabled !',ch10,& 603 & 'Action: change the input variable use_gpu_cuda',ch10,& 604 & ' or re-compile ABINIT with double-precision Cuda enabled.' 605 ABI_ERROR(msg) 606 #endif 607 end if 608 609 ABI_FREE(dprarr) 610 ABI_FREE(intarr) 611 612 ! We allocate the internal array, depending on the computed values. 613 ! WARNING: do not forget to deallocate these arrays in the routine dtset_free 614 ! (should make a separate subroutine for allocating/deallocating these records) 615 do idtset=0,ndtset_alloc 616 ABI_MALLOC(dtsets(idtset)%acell_orig,(3,mxnimage)) 617 ABI_MALLOC(dtsets(idtset)%algalch,(mxntypat)) 618 ABI_MALLOC(dtsets(idtset)%amu_orig,(mxntypat,mxnimage)) 619 ABI_MALLOC(dtsets(idtset)%cellcharge,(mxnimage)) 620 ABI_MALLOC(dtsets(idtset)%chrgat,(mxnatom)) 621 ABI_MALLOC(dtsets(idtset)%constraint_kind,(mxntypat)) 622 ABI_MALLOC(dtsets(idtset)%corecs,(mxntypat)) 623 ABI_MALLOC(dtsets(idtset)%densty,(mxntypat,4)) 624 ABI_MALLOC(dtsets(idtset)%dynimage,(mxnimage)) 625 ABI_MALLOC(dtsets(idtset)%iatfix,(3,mxnatom)) 626 ABI_MALLOC(dtsets(idtset)%f4of2_sla,(mxntypat)) 627 ABI_MALLOC(dtsets(idtset)%f6of2_sla,(mxntypat)) 628 ABI_MALLOC(dtsets(idtset)%jpawu,(mxntypat,mxnimage)) 629 ABI_MALLOC(dtsets(idtset)%kberry,(3,20)) 630 ABI_MALLOC(dtsets(idtset)%lambsig,(mxntypat)) 631 ABI_MALLOC(dtsets(idtset)%lexexch,(mxntypat)) 632 ABI_MALLOC(dtsets(idtset)%ldaminushalf,(mxntypat)) 633 ABI_MALLOC(dtsets(idtset)%lpawu,(mxntypat)) 634 ABI_MALLOC(dtsets(idtset)%mixalch_orig,(npsp,mxntypat,mxnimage)) 635 ABI_MALLOC(dtsets(idtset)%mixesimgf,(mxnimage)) 636 ABI_MALLOC(dtsets(idtset)%nucdipmom,(3,mxnatom)) 637 ABI_MALLOC(dtsets(idtset)%pimass,(mxntypat)) 638 ABI_MALLOC(dtsets(idtset)%ptcharge,(mxntypat)) 639 ABI_MALLOC(dtsets(idtset)%prtatlist,(mxnatom)) 640 ABI_MALLOC(dtsets(idtset)%quadmom,(mxntypat)) 641 ABI_MALLOC(dtsets(idtset)%ratsph,(mxntypat)) 642 ABI_MALLOC(dtsets(idtset)%rprim_orig,(3,3,mxnimage)) 643 ABI_MALLOC(dtsets(idtset)%rprimd_orig,(3,3,mxnimage)) 644 ABI_MALLOC(dtsets(idtset)%so_psp,(npsp)) 645 ABI_MALLOC(dtsets(idtset)%spinat,(3,mxnatom)) 646 ABI_MALLOC(dtsets(idtset)%shiftk,(3,MAX_NSHIFTK)) 647 ABI_MALLOC(dtsets(idtset)%typat,(mxnatom)) 648 ABI_MALLOC(dtsets(idtset)%upawu,(mxntypat,mxnimage)) 649 ABI_MALLOC(dtsets(idtset)%plowan_iatom,(mxnatom)) 650 ABI_MALLOC(dtsets(idtset)%plowan_it,(100*3)) 651 ABI_MALLOC(dtsets(idtset)%plowan_nbl,(mxnatom)) 652 ABI_MALLOC(dtsets(idtset)%plowan_lcalc,(12*mxnatom)) 653 ABI_MALLOC(dtsets(idtset)%plowan_projcalc,(12*mxnatom)) 654 ABI_MALLOC(dtsets(idtset)%vel_orig,(3,mxnatom,mxnimage)) 655 ABI_MALLOC(dtsets(idtset)%vel_cell_orig,(3,3,mxnimage)) 656 ABI_MALLOC(dtsets(idtset)%xred_orig,(3,mxnatom,mxnimage)) 657 ABI_MALLOC(dtsets(idtset)%ziontypat,(mxntypat)) 658 ABI_MALLOC(dtsets(idtset)%znucl,(npsp)) 659 end do 660 661 !DEBUG 662 !write(std_out,*)' invars0 : nimage, mxnimage = ',dtsets(:)%nimage, mxnimage 663 !write(std_out,*)' invars0 : natom = ',dtsets(:)%natom 664 !write(std_out,*)' invars0 : mxnatom = ',mxnatom 665 !write(std_out,*)' m_invars1%invars0 : exit ' 666 !call flush(std_out) 667 !ENDDEBUG 668 669 end subroutine invars0
ABINIT/invars1 [ Functions ]
NAME
invars1
FUNCTION
Initialize the dimensions needed to allocate the input arrays for one dataset characterized by jdtset, by taking from string the necessary data. Perform some preliminary checks and echo these dimensions.
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 comm= MPI communicator
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,chrgat,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, cellcharge 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, cellcharge They should be kept consistent with defaults of the same variables provided to the invars routines.
SOURCE
1125 subroutine invars1(bravais,dtset,iout,jdtset,lenstr,mband_upper,msym,npsp1,& 1126 & string,symafm,symrel,tnons,zionpsp, comm) 1127 1128 !Arguments ------------------------------------ 1129 !scalars 1130 integer,intent(in) :: iout,jdtset,lenstr,msym,npsp1, comm 1131 integer,intent(out) :: mband_upper 1132 character(len=*),intent(inout) :: string 1133 type(dataset_type),intent(inout) :: dtset 1134 !arrays 1135 integer,intent(inout) :: bravais(11),symafm(msym),symrel(3,3,msym) 1136 real(dp),intent(inout) :: tnons(3,msym) 1137 real(dp),intent(in) :: zionpsp(npsp1) 1138 1139 !Local variables------------------------------- 1140 !scalars 1141 integer,parameter :: master = 0 1142 integer :: chksymbreak,expert_user,found,ierr,iatom,ii,ikpt,iimage,index_blank,index_lower, tread_geo 1143 integer :: index_typsymb,index_upper,ipsp,iscf,intimage,itypat,leave,marr 1144 integer :: natom,nkpt,nkpthf,npsp,npspalch, ncid 1145 integer :: nqpt,nspinor,nsppol,ntypat,ntypalch,ntyppure,occopt,response 1146 integer :: rfddk,rfelfd,rfphon,rfstrs,rfuser,rf2_dkdk,rf2_dkde,rfmagn 1147 integer :: tfband,tnband,tread,tread_alt, my_rank, nprocs 1148 real(dp) :: cellcharge,cellcharge_min, fband,kptnrm,kptrlen,sum_spinat,zelect,zval 1149 character(len=1) :: blank=' ',string1 1150 character(len=2) :: string2,symbol 1151 character(len=500) :: msg 1152 type(atomdata_t) :: atom 1153 !arrays 1154 integer :: cond_values(4),vacuum(3) 1155 integer,allocatable :: iatfix(:,:),intarr(:),istwfk(:),nband(:),typat(:) 1156 real(dp) :: acell(3),rprim(3,3) 1157 !real(dp) :: field(3) 1158 real(dp),allocatable :: amu(:),chrgat(:),dprarr(:),kpt(:,:),kpthf(:,:),mixalch(:,:),nucdipmom(:,:) 1159 real(dp),allocatable :: ratsph(:),reaalloc(:),spinat(:,:) 1160 real(dp),allocatable :: vel(:,:),vel_cell(:,:),wtk(:),xred(:,:),znucl(:) 1161 character(len=32) :: cond_string(4) 1162 character(len=fnlen) :: key_value 1163 character(len=len(string)) :: geo_string 1164 type(geo_t) :: geo 1165 1166 !************************************************************************ 1167 1168 !DEBUG 1169 !write(std_out,'(a)')' m_invars1%invars1 : enter ' 1170 !call flush(std_out) 1171 !ENDDEBUG 1172 1173 my_rank = xmpi_comm_rank(comm); nprocs = xmpi_comm_size(comm) 1174 1175 ! This counter is incremented when we find a non-critical error. 1176 ! The code outputs a warning and stops at end. 1177 leave = 0 1178 1179 ! Some initialisations 1180 ierr=0 1181 cond_string(1:4)=' ' 1182 cond_values(1:4)=(/0,0,0,0/) 1183 1184 ! Read parameters 1185 marr=dtset%npsp;if (dtset%npsp<3) marr=3 1186 marr=max(marr,dtset%nimage) 1187 ABI_MALLOC(intarr,(marr)) 1188 ABI_MALLOC(dprarr,(marr)) 1189 1190 !--------------------------------------------------------------------------- 1191 1192 rfddk=0; rfelfd=0; rfphon=0; rfmagn=0; rfstrs=0; rfuser=0; rf2_dkdk=0; rf2_dkde=0 1193 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfddk',tread,'INT') 1194 if(tread==1) rfddk=intarr(1) 1195 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfelfd',tread,'INT') 1196 if(tread==1) rfelfd=intarr(1) 1197 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfmagn',tread,'INT') 1198 if(tread==1) rfmagn=intarr(1) 1199 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfphon',tread,'INT') 1200 if(tread==1) rfphon=intarr(1) 1201 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfstrs',tread,'INT') 1202 if(tread==1) rfstrs=intarr(1) 1203 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rfuser',tread,'INT') 1204 if(tread==1) rfuser=intarr(1) 1205 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rf2_dkdk',tread,'INT') 1206 if(tread==1) rf2_dkdk=intarr(1) 1207 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'rf2_dkde',tread,'INT') 1208 if(tread==1) rf2_dkde=intarr(1) 1209 1210 response=0 1211 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 1212 1213 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'optdriver',tread,'INT') 1214 if (tread==1) then 1215 dtset%optdriver=intarr(1) 1216 else 1217 ! If optdriver was not read, while response=1, set optdriver to 1 1218 if(response==1)dtset%optdriver=1 1219 end if 1220 1221 !--------------------------------------------------------------------------- 1222 !For now, waiting express parallelisation for recursion 1223 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'tfkinfunc',tread,'INT') 1224 if(tread==1) dtset%tfkinfunc=intarr(1) 1225 1226 !--------------------------------------------------------------------------- 1227 ! wvl_bigdft_comp, done here since default values of nline, nwfshist and iscf depend on its value (see indefo) 1228 if(dtset%usewvl==1) then 1229 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'wvl_bigdft_comp',tread,'INT') 1230 if(tread==1) dtset%wvl_bigdft_comp=intarr(1) 1231 end if 1232 1233 !--------------------------------------------------------------------------- 1234 1235 natom=dtset%natom 1236 npsp=dtset%npsp 1237 ntypat=dtset%ntypat 1238 1239 call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), 'structure', tread_geo, & 1240 'KEY', key_value=geo_string) 1241 1242 if (tread_geo == 0) then 1243 ! No default value for znucl 1244 call intagm(dprarr,intarr,jdtset,marr,dtset%npsp,string(1:lenstr),'znucl',tread,'DPR') 1245 if(tread==1) dtset%znucl(1:dtset%npsp)=dprarr(1:dtset%npsp) 1246 1247 if(tread/=1)then 1248 write(msg, '(3a)' )& 1249 'The array znucl MUST be initialized in the input file while this is not done.',ch10,& 1250 'Action: initialize znucl in your input file.' 1251 ABI_ERROR(msg) 1252 end if 1253 1254 else 1255 call wrtout(std_out, sjoin(" Initializing lattice and positions from:", geo_string)) 1256 geo = geo_from_abivar_string(geo_string, comm) 1257 dtset%znucl(1:dtset%ntypat) = geo%znucl 1258 call geo%free() 1259 end if 1260 1261 ! The default for ratsph has already been initialized 1262 call intagm(dprarr,intarr,jdtset,marr,dtset%ntypat,string(1:lenstr),'ratsph',tread,'LEN') 1263 if(tread==1)then 1264 do ii=1,dtset%ntypat 1265 dtset%ratsph(ii)=dprarr(ii) 1266 end do 1267 end if 1268 ABI_MALLOC(ratsph,(dtset%ntypat)) 1269 do ii=1,dtset%ntypat 1270 ratsph(ii)=dtset%ratsph(ii) 1271 end do 1272 1273 !Special treatment of _TYPAX (from a XYZ file), taking into account 1274 !the fact that znucl does NOT depend on the dataset 1275 !Examine all occurences of '_TYPAX' 1276 1277 do 1278 index_typsymb=index(string(1:lenstr),'_TYPAX') 1279 if(index_typsymb==0)exit 1280 ! Replace '_TYPAX' by '_TYPAT' 1281 string(index_typsymb:index_typsymb+5)='_TYPAT' 1282 index_upper=index_typsymb+5 1283 ! Must start from the first blank after the tag (including possible dtset_char) 1284 index_upper=index(string(index_upper:lenstr),blank)+index_upper-1 1285 index_lower=index_upper 1286 1287 ! Examine all atoms (the end of the symbol string is delimited by a XX ) 1288 do 1289 index_blank=index(string(index_upper:lenstr),blank)+index_upper-1 1290 string2=string(index_blank+1:index_blank+2) 1291 if(string2=="XX")exit 1292 found=0 1293 ! Find the matching symbol 1294 do ipsp=1,dtset%npsp 1295 call atomdata_from_znucl(atom,dtset%znucl(ipsp)) 1296 symbol = atom%symbol 1297 call inupper(symbol) 1298 call inupper(string2) 1299 ! write(std_out,'(a)')' invars1 : before test, trim(adjustl(symbol)),trim(adjustl(string2))' 1300 ! write(std_out,'(5a)' )'"',trim(adjustl(symbol)),'","',trim(adjustl(string2)),'"' 1301 if(trim(adjustl(symbol))==trim(adjustl(string2)))then 1302 found=1 1303 index_upper=index_blank+1 1304 ! Cannot deal properly with more that 9 psps 1305 if(ipsp>=10)then 1306 ABI_ERROR('Need to use a pseudopotential with number larger than 9. Not allowed yet.') 1307 end if 1308 1309 ! write(std_out,*)' invars1 : found ipsp=',ipsp 1310 write(string1,'(i1)')ipsp 1311 string(index_lower:index_lower+1)=blank//string1 1312 index_lower=index_lower+2 1313 end if 1314 end do ! ipsp 1315 ! if not found ... 1316 if(found==0)then 1317 write(msg,'(6a)' )& 1318 & 'Did not find matching pseudopotential for XYZ atomic symbol,',ch10,& 1319 & 'with value ',string2,ch10,& 1320 & 'Action: check that the atoms required by the XYZ file correspond to one psp file.' 1321 ABI_ERROR(msg) 1322 end if 1323 end do ! Loop on atoms 1324 ! One should find blanks after the last significant type value 1325 string(index_lower:index_blank+2)=blank 1326 end do ! loop to identify _TYPAX 1327 1328 !--------------------------------------------------------------------------- 1329 1330 ! Here, set up quantities that are related to geometrical description of the system (acell,rprim,xred), as well as 1331 ! initial velocity(vel), cellcharge (to compute mband_upper) and spin of atoms (chrgat,spinat), nuclear dipole moments of atoms (nucdipmom), 1332 ! the symmetries (symrel,symafm, and tnons) and the list of fixed atoms (iatfix,iatfixx,iatfixy,iatfixz). 1333 ! Arrays have already been dimensioned thanks to the knowledge of msym and mx%natom 1334 1335 !ji: We need to read the electric field before calling ingeo 1336 !****** Temporary ****** 1337 1338 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'berryopt',tread,'INT') 1339 if(tread==1) dtset%berryopt=intarr(1) 1340 1341 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'berrysav',tread,'INT') 1342 if(tread==1) dtset%berrysav=intarr(1) 1343 1344 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'bfield',tread,'DPR') 1345 if (tread==1) dtset%bfield(1:3) = dprarr(1:3) 1346 1347 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'dfield',tread,'DPR') 1348 if (tread==1) dtset%dfield(1:3) = dprarr(1:3) 1349 1350 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'efield',tread,'DPR') 1351 if (tread==1) dtset%efield(1:3) = dprarr(1:3) 1352 1353 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'red_dfield',tread,'DPR') 1354 if (tread==1) dtset%red_dfield(1:3) = dprarr(1:3) 1355 1356 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'red_efield',tread,'DPR') 1357 if (tread==1) dtset%red_efield(1:3) = dprarr(1:3) 1358 1359 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'red_efieldbar',tread,'DPR') 1360 if (tread==1) dtset%red_efieldbar(1:3) = dprarr(1:3) 1361 1362 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'jfielddir',tread,'INT') 1363 if(tread==1) dtset%jfielddir(1:3)=intarr(1:3) 1364 1365 ! We need to know nsppol/nspinor/nspden before calling ingeo 1366 nsppol=dtset%nsppol 1367 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nsppol',tread,'INT') 1368 if(tread==1) nsppol=intarr(1) 1369 1370 !Alternate SIESTA definition of nsppol 1371 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'SpinPolarized',tread_alt,'LOG') 1372 if(tread_alt==1)then 1373 if(tread==1)then 1374 msg = 'nsppol and SpinPolarized cannot be specified simultaneously for the same dataset.' 1375 ABI_ERROR_NOSTOP(msg, leave) 1376 else 1377 ! Note that SpinPolarized is a logical input variable 1378 nsppol=1 1379 if(intarr(1)==1)nsppol=2 1380 tread=1 1381 end if 1382 end if 1383 dtset%nsppol=nsppol 1384 1385 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nspinor',tread,'INT') 1386 if(tread==1) dtset%nspinor=intarr(1) 1387 1388 !Has to read pawspnorb now, in order to adjust nspinor 1389 !Also, if nspinor=1, turn on spin-orbit coupling by default, here for the PAW case. NC case is treated elsewhere. 1390 if (dtset%usepaw>0)then 1391 ! Change the default value 1392 if(dtset%nspinor==2)dtset%pawspnorb=1 1393 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'pawspnorb',tread,'INT') 1394 if(tread==1)then 1395 dtset%pawspnorb=intarr(1) 1396 if(dtset%pawspnorb>0) dtset%nspinor=2 1397 else 1398 if(dtset%nspinor==2)then 1399 write(msg, '(4a)' ) ch10,& 1400 & ' invars1: COMMENT -',ch10,& 1401 & ' With nspinor=2 and usepaw=1, pawspnorb=1 has been switched on by default.' 1402 call wrtout(iout, msg,'COLL') 1403 end if 1404 end if 1405 end if 1406 nspinor=dtset%nspinor 1407 1408 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nspden',tread,'INT') 1409 if(tread==1) then 1410 dtset%nspden=intarr(1) 1411 else 1412 dtset%nspden=dtset%nsppol 1413 end if 1414 1415 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ntypalch',tread,'INT') 1416 if(tread==1) dtset%ntypalch=intarr(1) 1417 1418 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nzchempot',tread,'INT') 1419 if(tread==1) dtset%nzchempot=intarr(1) 1420 1421 ntypalch=dtset%ntypalch 1422 if(ntypalch>ntypat)then 1423 write(msg, '(3a,i0,a,i0,a,a)' )& 1424 'The input variable ntypalch must be smaller than ntypat, while it is',ch10,& 1425 'ntypalch=',dtset%ntypalch,', and ntypat=',ntypat,ch10,& 1426 'Action: check ntypalch vs ntypat in your input file.' 1427 ABI_ERROR(msg) 1428 end if 1429 1430 ntyppure=ntypat-ntypalch 1431 dtset%ntyppure=ntyppure 1432 npspalch=npsp-ntyppure 1433 dtset%npspalch=npspalch 1434 if(npspalch<0)then 1435 write(msg, '(a,i0,2a,i0,a,a)' )& 1436 'The number of available pseudopotentials, npsp=',npsp,ch10,& 1437 'is smaller than the requested number of types of pure atoms, ntyppure=',ntyppure,ch10,& 1438 'Action: check ntypalch versus ntypat and npsp in your input file.' 1439 ABI_ERROR(msg) 1440 end if 1441 1442 if(ntypalch>0)then 1443 call intagm(dprarr,intarr,jdtset,marr,ntypalch,string(1:lenstr),'algalch',tread,'INT') 1444 if(tread==1) dtset%algalch(1:ntypalch)=intarr(1:ntypalch) 1445 if (tread_geo /= 0) then 1446 ABI_ERROR("Alchemical mixing cannot be used with geo variable, use typat, znucl etc.") 1447 end if 1448 end if 1449 1450 !Read the Zeeman field 1451 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),'zeemanfield',tread,'BFI') 1452 if(tread==1) then 1453 if(dtset%nspden == 2)then 1454 write(msg,'(7a)')& 1455 'A Zeeman field has been specified without noncollinear spins.',ch10,& 1456 'Only the z-component of the magnetic field will be used.' 1457 ABI_WARNING(msg) 1458 else if (dtset%nspden == 1)then 1459 write(msg, '(a,a,a)' )& 1460 'A Zeeman field has been specified for a non-spin-polarized calculation.',ch10,& 1461 'Action: check the input file.' 1462 ABI_ERROR(msg) 1463 end if 1464 1465 dtset%zeemanfield(1:3) = dprarr(1:3) 1466 end if 1467 1468 !Initialize geometry of the system, for different images. Also initialize cellcharge_min to be used later for estimating mband_upper.. 1469 ABI_MALLOC(amu,(ntypat)) 1470 ABI_MALLOC(mixalch,(npspalch,ntypalch)) 1471 ABI_MALLOC(vel,(3,natom)) 1472 ABI_MALLOC(vel_cell,(3,3)) 1473 ABI_MALLOC(xred,(3,natom)) 1474 !Only take into account negative cellcharge, to compute maximum number of bands, so initialize cellcharge_min to zero 1475 cellcharge_min=zero 1476 intimage=2 ; if(dtset%nimage==1)intimage=1 1477 do ii=1,dtset%nimage+1 1478 iimage=ii 1479 if(dtset%nimage==1 .and. ii==2)exit 1480 if(dtset%nimage==2 .and. ii==3)exit 1481 if(dtset%nimage> 2 .and. ii==intimage)cycle ! Will do the intermediate reference image at the last reading 1482 if(dtset%nimage>=2 .and. ii==dtset%nimage+1)iimage=intimage 1483 1484 if (dtset%nimage /= 1) call wrtout(std_out, sjoin(' invars1: treat image number: ',itoa(iimage))) 1485 1486 ! Need to reset nsym to default value for each image 1487 dtset%nsym=0 1488 1489 ! Call ingeo for each image in turn, with the possible default values 1490 acell=dtset%acell_orig(1:3,iimage) 1491 amu=dtset%amu_orig(1:ntypat,iimage) 1492 mixalch=dtset%mixalch_orig(1:npspalch,1:ntypalch,iimage) 1493 rprim=dtset%rprim_orig(1:3,1:3,iimage) 1494 vel=dtset%vel_orig(1:3,1:natom,iimage) 1495 vel_cell=dtset%vel_cell_orig(1:3,1:3,iimage) 1496 xred=dtset%xred_orig(1:3,1:natom,iimage) 1497 ABI_MALLOC(chrgat,(natom)) 1498 ABI_MALLOC(iatfix,(3,natom)) 1499 ABI_MALLOC(nucdipmom,(3,natom)) 1500 ABI_MALLOC(spinat,(3,natom)) 1501 ABI_MALLOC(typat,(natom)) 1502 ABI_MALLOC(znucl,(dtset%npsp)) 1503 chrgat(1:natom)=dtset%chrgat(1:natom) 1504 nucdipmom(1:3,1:natom)=dtset%nucdipmom(1:3,1:natom) 1505 spinat(1:3,1:natom)=dtset%spinat(1:3,1:natom) 1506 znucl(1:dtset%npsp)=dtset%znucl(1:dtset%npsp) 1507 1508 call ingeo(acell,amu,bravais,chrgat,dtset,dtset%genafm(1:3),iatfix,& 1509 dtset%icoulomb,iimage,iout,jdtset,dtset%jellslab,lenstr,mixalch,& 1510 msym,natom,dtset%nimage,dtset%npsp,npspalch,dtset%nspden,dtset%nsppol,& 1511 dtset%nsym,ntypalch,dtset%ntypat,nucdipmom,dtset%nzchempot,& 1512 dtset%pawspnorb,dtset%ptgroupma,ratsph,& 1513 rprim,dtset%slabzbeg,dtset%slabzend,dtset%spgroup,spinat,& 1514 string,dtset%supercell_latt,symafm,dtset%symmorphi,symrel,tnons,dtset%tolsym,& 1515 typat,vel,vel_cell,xred,znucl, comm) 1516 1517 dtset%chrgat(1:natom)=chrgat(1:natom) 1518 dtset%iatfix(1:3,1:natom)=iatfix(1:3,1:natom) 1519 dtset%nucdipmom(1:3,1:natom)=nucdipmom(1:3,1:natom) 1520 dtset%spinat(1:3,1:natom)=spinat(1:3,1:natom) 1521 dtset%typat(1:natom)=typat(1:natom) 1522 ABI_FREE(chrgat) 1523 ABI_FREE(iatfix) 1524 ABI_FREE(nucdipmom) 1525 ABI_FREE(spinat) 1526 ABI_FREE(typat) 1527 ABI_FREE(znucl) 1528 dtset%acell_orig(1:3,iimage)=acell 1529 dtset%amu_orig(1:ntypat,iimage)=amu 1530 dtset%mixalch_orig(1:npspalch,1:ntypalch,iimage)=mixalch 1531 dtset%rprim_orig(1:3,1:3,iimage)=rprim 1532 dtset%vel_orig(1:3,1:natom,iimage)=vel 1533 dtset%vel_cell_orig(1:3,1:3,iimage)=vel_cell 1534 dtset%xred_orig(1:3,1:natom,iimage)=xred 1535 call mkrdim(dtset%acell_orig(1:3,iimage),dtset%rprim_orig(1:3,1:3,iimage),dtset%rprimd_orig(1:3,1:3,iimage)) 1536 1537 ! Read cellcharge for each image, but use it only to initialize cellcharge_min 1538 ! The old name 'charge' is still tolerated. Will be removed in due time. 1539 cellcharge=zero 1540 ! Initialize cellcharge with the value for the first image 1541 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'cellcharge',tread,'DPR') 1542 if(tread==1)then 1543 cellcharge=dprarr(1) 1544 else 1545 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'charge',tread,'DPR') 1546 if(tread==1) cellcharge=dprarr(1) 1547 endif 1548 ! Possibly overwrite cellcharge from the first image with a specific value for the current image 1549 call intagm_img(dprarr,iimage,jdtset,lenstr,dtset%nimage,1,string,'cellcharge',tread_alt,'DPR') 1550 if(tread_alt==1)then 1551 cellcharge=dprarr(1) 1552 else 1553 call intagm_img(dprarr,iimage,jdtset,lenstr,dtset%nimage,1,string,'charge',tread_alt,'DPR') 1554 if(tread_alt==1) cellcharge=dprarr(1) 1555 endif 1556 1557 if(cellcharge < cellcharge_min)cellcharge_min=cellcharge 1558 1559 end do 1560 1561 ABI_FREE(amu) 1562 ABI_FREE(mixalch) 1563 ABI_FREE(vel) 1564 ABI_FREE(vel_cell) 1565 ABI_FREE(xred) 1566 1567 ! Examine whether there is some vacuum space in the unit cell 1568 call invacuum(jdtset,lenstr,natom,dtset%rprimd_orig(1:3,1:3,intimage),string,vacuum,& 1569 & dtset%xred_orig(1:3,1:natom,intimage)) 1570 1571 !write(std_out,*)' invars1: before inkpts, dtset%mixalch_orig(1:npspalch,1:ntypalch,:)=',& 1572 !dtset%mixalch_orig(1:npspalch,1:ntypalch,1:dtset%nimage) 1573 1574 !--------------------------------------------------------------------------- 1575 1576 !Set up k point grid number 1577 !First, get additional information 1578 dtset%kptopt=1 1579 if(dtset%nspden==4)dtset%kptopt=4 1580 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'kptopt',tread,'INT') 1581 if(tread==1) dtset%kptopt=intarr(1) 1582 1583 dtset%qptopt=1 1584 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'qptopt',tread,'INT') 1585 if(tread==1) dtset%qptopt=intarr(1) 1586 1587 iscf=5 1588 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'iscf',tread,'INT') 1589 if(tread==1) iscf=intarr(1) 1590 1591 dtset%natsph=dtset%natom 1592 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natsph',tread,'INT') 1593 if(tread==1) dtset%natsph=intarr(1) 1594 1595 dtset%natsph_extra=0 1596 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natsph_extra',tread,'INT') 1597 if(tread==1) dtset%natsph_extra=intarr(1) 1598 1599 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'natvshift',tread,'INT') 1600 if(tread==1) dtset%natvshift=intarr(1) 1601 1602 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nconeq',tread,'INT') 1603 if(tread==1) dtset%nconeq=intarr(1) 1604 1605 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nkptgw',tread,'INT') 1606 if(tread==1) dtset%nkptgw=intarr(1) 1607 if (dtset%nkptgw<0) then 1608 write(msg, '(a,i0,4a)' )& 1609 'Input nkptgw must be >= 0, but was ',dtset%nkptgw,ch10,& 1610 'This is not allowed.',ch10,'Action: check the input file.' 1611 ABI_ERROR(msg) 1612 end if 1613 1614 ! Number of points for long wavelength limit. Default is dtset%gw_nqlwl=0 1615 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_nqlwl',tread,'INT') 1616 if(tread==1) dtset%gw_nqlwl=intarr(1) 1617 if (dtset%gw_nqlwl<0) then 1618 write(msg, '(a,i0,4a)' )& 1619 'Input gw_nqlwl must be > 0, but was ',dtset%gw_nqlwl,ch10,& 1620 'This is not allowed.',ch10,'Action: check the input file.' 1621 ABI_ERROR(msg) 1622 end if 1623 1624 ! Read number of k-points from input file (if specified) 1625 nkpt=0 1626 if(dtset%kptopt==0)nkpt=1 1627 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nkpt',tread,'INT') 1628 if(tread==1) nkpt=intarr(1) 1629 1630 ! or from KERANGE file. 1631 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr), "getkerange_filepath", tread, 'KEY', key_value=key_value) 1632 if (tread==1) dtset%getkerange_filepath = key_value 1633 1634 #ifdef HAVE_NETCDF 1635 if (dtset%getkerange_filepath /= ABI_NOFILE) then 1636 ! Get number of k-points in sigma_erange energy windows. 1637 !dtset%kptopt = 0 1638 if (my_rank == master) then 1639 NCF_CHECK(nctk_open_read(ncid, dtset%getkerange_filepath, xmpi_comm_self)) 1640 NCF_CHECK(nctk_get_dim(ncid, "nkpt_inerange", nkpt, datamode=.True.)) 1641 NCF_CHECK(nf90_close(ncid)) 1642 end if 1643 call xmpi_bcast(nkpt, master, comm, ierr) 1644 end if 1645 #endif 1646 1647 dtset%nkpt = nkpt 1648 1649 call chkint_ge(0,0,cond_string,cond_values,ierr,'nkpt',nkpt,0,iout) 1650 if (dtset%kptopt==0) then 1651 cond_string(1)='kptopt'; cond_values(1)=0 1652 call chkint_ge(1,1,cond_string,cond_values,ierr,'nkpt',nkpt,1,iout) 1653 end if 1654 1655 nkpthf=nkpt 1656 dtset%nkpthf=nkpt 1657 1658 ! Will compute the actual value of nkpt, if needed. Otherwise, 1659 ! test that the value of nkpt is OK, if kptopt/=0 1660 ! Set up dummy arrays istwfk, kpt, wtk 1661 1662 if(nkpt/=0 .or. dtset%kptopt/=0)then 1663 ABI_MALLOC(istwfk,(nkpt)) 1664 ABI_MALLOC(kpt,(3,nkpt)) 1665 ABI_MALLOC(kpthf,(3,nkpthf)) 1666 ABI_MALLOC(wtk,(nkpt)) 1667 ! Here, occopt is also a dummy argument 1668 occopt=1; dtset%nshiftk=1; dtset%kptrlatt(:,:)=0 1669 1670 kptrlen=20.0_dp ; wtk(:)=1.0_dp 1671 dtset%shiftk(:,:)=half 1672 1673 nqpt=0 1674 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nqpt',tread,'INT') 1675 if(tread==1) nqpt=intarr(1) 1676 1677 expert_user=0 1678 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'expert_user',tread,'INT') 1679 if (tread==1) expert_user=intarr(1) 1680 1681 ! The default value of chksymbreak depends on expert_user but we still allow user to specify it. 1682 chksymbreak=1; if (expert_user > 0) chksymbreak = 0 1683 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'chksymbreak',tread,'INT') 1684 if(tread==1) chksymbreak=intarr(1) 1685 1686 ! Use the first image to predict k and/or q points, except if an intermediate image is available 1687 intimage=1; if(dtset%nimage>2)intimage=(1+dtset%nimage)/2 1688 1689 ! Find the q-point, if any. 1690 if(nqpt/=0)then 1691 call inqpt(chksymbreak,std_out,jdtset,lenstr,msym,natom,dtset%qptn,dtset%wtq,& 1692 dtset%rprimd_orig(1:3,1:3,intimage),dtset%spinat,string,dtset%typat,& 1693 vacuum,dtset%xred_orig(1:3,1:natom,intimage),dtset%qptrlatt) 1694 endif 1695 1696 ! Find the k point grid 1697 call inkpts(bravais,chksymbreak,dtset%fockdownsampling,iout,iscf,istwfk,jdtset,& 1698 kpt,kpthf,dtset%kptopt,kptnrm,dtset%kptrlatt_orig,dtset%kptrlatt,kptrlen,lenstr,msym, dtset%getkerange_filepath, & 1699 nkpt,nkpthf,nqpt,dtset%ngkpt,dtset%nshiftk,dtset%nshiftk_orig,dtset%shiftk_orig,dtset%nsym,& 1700 occopt,dtset%qptn,response,dtset%rprimd_orig(1:3,1:3,intimage),dtset%shiftk,& 1701 string,symafm,symrel,vacuum,wtk,comm) 1702 1703 ABI_FREE(istwfk) 1704 ABI_FREE(kpt) 1705 ABI_FREE(kpthf) 1706 ABI_FREE(wtk) 1707 1708 ! nkpt and nkpthf have been computed, as well as the k point grid, if needed 1709 dtset%nkpt=nkpt 1710 dtset%nkpthf=nkpthf 1711 end if 1712 1713 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nqptdm',tread,'INT') 1714 if(tread==1) dtset%nqptdm=intarr(1) 1715 1716 if (dtset%nqptdm<-1) then 1717 write(msg, '(a,i0,4a)' )& 1718 'Input nqptdm must be >= 0, but was ',dtset%nqptdm,ch10,& 1719 'This is not allowed.',ch10,'Action: check the input file.' 1720 ABI_ERROR(msg) 1721 end if 1722 1723 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nzchempot',tread,'INT') 1724 if(tread==1) dtset%nzchempot=intarr(1) 1725 1726 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'cd_customnimfrqs',tread,'INT') 1727 if(tread==1) dtset%cd_customnimfrqs=intarr(1) 1728 1729 if (dtset%cd_customnimfrqs<0) then 1730 write(msg, '(a,i0,4a)' )& 1731 'Input cd_customnimfrqs must be >= 0, but was ',dtset%cd_customnimfrqs,ch10,& 1732 'This is not allowed.',ch10,'Action: check the input file.' 1733 ABI_ERROR(msg) 1734 end if 1735 1736 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gw_customnfreqsp',tread,'INT') 1737 if(tread==1) dtset%gw_customnfreqsp=intarr(1) 1738 1739 if (dtset%gw_customnfreqsp<0) then 1740 write(msg, '(a,i0,4a)' )& 1741 'Input gw_customnfreqsp must be >= 0, but was ',dtset%gw_customnfreqsp,ch10,& 1742 'This is not allowed.',ch10,'Action: check the input file.' 1743 ABI_ERROR(msg) 1744 end if 1745 1746 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'gwls_n_proj_freq',tread,'INT') 1747 if(tread==1) dtset%gwls_n_proj_freq=intarr(1) 1748 1749 if (dtset%gwls_n_proj_freq<0) then 1750 write(msg, '(a,i0,4a)' )& 1751 'Input gwls_n_proj_freq must be >= 0, but was ',dtset%gwls_n_proj_freq,ch10,& 1752 'This is not allowed.',ch10,'Action: check the input file.' 1753 ABI_ERROR(msg) 1754 end if 1755 1756 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'efmas_calc_dirs',tread,'INT') 1757 if(tread==1) dtset%efmas_calc_dirs=intarr(1) 1758 1759 if (ABS(dtset%efmas_calc_dirs)>3) then 1760 write(msg, '(a,i0,4a)' )& 1761 'Input efmas_calc_dirs must be between -3 and 3, but was ',dtset%efmas_calc_dirs,ch10,& 1762 'This is not allowed.',ch10,'Action: check the input file.' 1763 ABI_ERROR(msg) 1764 end if 1765 1766 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'efmas_n_dirs',tread,'INT') 1767 if(tread==1) dtset%efmas_n_dirs=intarr(1) 1768 1769 if (dtset%efmas_n_dirs<0) then 1770 write(msg, '(a,i0,4a)' )& 1771 'Input efmas_n_dirs must be >= 0, but was ',dtset%efmas_n_dirs,ch10,& 1772 'This is not allowed.',ch10,'Action: check the input file.' 1773 ABI_ERROR(msg) 1774 end if 1775 1776 !--------------------------------------------------------------------------- 1777 1778 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nnos',tread,'INT') 1779 if(tread==1) dtset%nnos=intarr(1) 1780 1781 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ga_n_rules',tread,'INT') 1782 if(tread==1) dtset%ga_n_rules=intarr(1) 1783 1784 ! Perform the first checks 1785 ! Check that nkpt is greater than 0 1786 if (nkpt<=0) then 1787 write(msg, '(a,i0)' )'After inkpts, nkpt must be > 0, but was ',nkpt 1788 ABI_ERROR_NOSTOP(msg, leave) 1789 end if 1790 1791 ! Check that nsppol is 1 or 2 1792 if (nsppol/=1 .and. nsppol/=2) then 1793 write(msg, '(a,i0)' )'Input nsppol must be 1 or 2, but was ',nsppol 1794 ABI_ERROR_NOSTOP(msg, leave) 1795 end if 1796 1797 ! Check that nspinor is 1 or 2 1798 if (nspinor/=1 .and. nspinor/=2) then 1799 write(msg, '(a,i0)' )'Input nspinor must be 1 or 2, but was ',nspinor 1800 ABI_ERROR_NOSTOP(msg, leave) 1801 end if 1802 1803 ! Check that nspinor and nsppol are not 2 together 1804 if (nsppol==2 .and. nspinor==2) then 1805 ABI_ERROR_NOSTOP('nspinor and nsppol cannot be 2 together!', leave) 1806 end if 1807 1808 ! Here, leave if an error has been detected earlier 1809 if (leave /= 0) then 1810 ABI_ERROR('Errors are present in the input file. See ABOVE messages') 1811 end if 1812 1813 ! Now, take care of mband_upper 1814 mband_upper=1 1815 occopt=1 1816 fband=0.5_dp 1817 1818 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'occopt',tread,'INT') 1819 if(tread==1) occopt=intarr(1) 1820 1821 ! Also read fband, that is an alternative to nband. The default 1822 ! is different for occopt==1 and for metallic occupations. 1823 if(occopt==1)fband=0.125_dp 1824 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'fband',tfband,'DPR') 1825 if(tfband==1)fband=dprarr(1) 1826 1827 ! fband cannot be used when occopt==0 or occopt==2 1828 if(tfband==1 .and. (occopt==0 .or. occopt==2) )then 1829 write(msg, '(3a)' )& 1830 'fband cannot be used if occopt==0 or occopt==2 ',ch10,& 1831 'Action: correct your input file, suppress fband, or change occopt.' 1832 ABI_ERROR(msg) 1833 end if 1834 1835 ABI_MALLOC(nband,(nkpt*nsppol)) 1836 tnband=0 1837 1838 ! Compute ziontypat 1839 ! When the pseudo-atom is pure, simple copy 1840 if(ntyppure>0)then 1841 do itypat=1,ntyppure 1842 dtset%ziontypat(itypat)=zionpsp(itypat) 1843 end do 1844 end if 1845 1846 ! When the pseudo-atom is alchemical, must make mixing 1847 if(ntypalch>0)then 1848 do itypat=ntyppure+1,ntypat 1849 dtset%ziontypat(itypat)=zero 1850 do ipsp=ntyppure+1,npsp 1851 dtset%ziontypat(itypat)=dtset%ziontypat(itypat) & 1852 & +dtset%mixalch_orig(ipsp-ntyppure,itypat-ntyppure,1)*zionpsp(ipsp) 1853 end do 1854 end do 1855 end if 1856 1857 ! CP modified 1858 !if (occopt==0 .or. occopt==1 .or. (occopt>=3 .and. occopt<=8) ) then 1859 if (occopt==0 .or. occopt==1 .or. (occopt>=3 .and. occopt<=9) ) then 1860 ! End CP modified 1861 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'nband',tnband,'INT') 1862 ! Note: mband_upper is initialized, not nband 1863 if(tnband==1) mband_upper=intarr(1) 1864 1865 if(tfband==1 .and. tnband==1)then 1866 write(msg, '(3a)' )& 1867 'fband and nband cannot be used together. ',ch10,& 1868 'Action: correct your input file, suppress either fband or nband.' 1869 ABI_ERROR(msg) 1870 end if 1871 1872 ! In case nband was not read, use fband, either read, or the default, 1873 ! to provide an upper limit for mband_upper 1874 if(tnband==0)then 1875 1876 ! mband_upper=nspinor*((nint(zion_max)*natom+1)/2 - floor(cellcharge_min/2.0_dp)& 1877 !& + ceiling(fband*natom-1.0d-10)) 1878 zval=zero 1879 sum_spinat=zero 1880 do iatom=1,natom 1881 zval=zval+dtset%ziontypat(dtset%typat(iatom)) 1882 sum_spinat=sum_spinat+dtset%spinat(3,iatom) 1883 end do 1884 zelect=zval-cellcharge_min 1885 mband_upper=nspinor * ((ceiling(zelect-tol10)+1)/2 + ceiling( fband*natom - tol10 )) & 1886 & + (nsppol-1)*(ceiling(half*(sum_spinat -tol10))) 1887 nband(:)=mband_upper 1888 1889 ! write(std_out,*)' invars1 : zion_max,natom,fband,mband_upper ' 1890 ! write(std_out,*)zion_max,natom,fband,mband_upper 1891 end if 1892 1893 nband(:)=mband_upper 1894 1895 else if (occopt==2) then 1896 ABI_MALLOC(reaalloc,(nkpt*nsppol)) 1897 call intagm(reaalloc,nband,jdtset,nkpt*nsppol,nkpt*nsppol,string(1:lenstr),'nband',tnband,'INT') 1898 if(tnband==1)then 1899 do ikpt=1,nkpt*nsppol 1900 if (nband(ikpt)>mband_upper) mband_upper=nband(ikpt) 1901 end do 1902 end if 1903 ABI_FREE(reaalloc) 1904 else 1905 write(msg, '(a,i0,3a)' )'occopt=',occopt,' is not an allowed value.',ch10,'Action: correct your input file.' 1906 ABI_ERROR(msg) 1907 end if 1908 1909 ! Check that mband_upper is greater than 0 1910 if (mband_upper<=0) then 1911 write(msg, '(a,i0,4a)' )& 1912 'Maximal nband must be > 0, but was ',mband_upper,ch10,& 1913 'This is not allowed.',ch10,'Action: check the input file.' 1914 ABI_ERROR(msg) 1915 end if 1916 1917 ! The following 3 values are needed to dimension the parallelism over images 1918 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'imgmov',tread,'INT') 1919 if(tread==1) dtset%imgmov=intarr(1) 1920 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ntimimage',tread,'INT') 1921 if(tread==1) dtset%ntimimage=intarr(1) 1922 call intagm(dprarr,intarr,jdtset,marr,dtset%nimage,string(1:lenstr),'dynimage',tread,'INT') 1923 if(tread==1)then 1924 dtset%dynimage(1:dtset%nimage)=intarr(1:dtset%nimage) 1925 else if (dtset%imgmov==2.or.dtset%imgmov==5) then 1926 dtset%dynimage(1)=0;dtset%dynimage(dtset%nimage)=0 1927 end if 1928 dtset%ndynimage=count(dtset%dynimage(1:dtset%nimage)/=0) 1929 1930 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'wfoptalg',tread,'INT') 1931 if(tread==1) then 1932 dtset%wfoptalg=intarr(1) 1933 else 1934 if (dtset%usepaw==0) dtset%wfoptalg=0 1935 if (dtset%usepaw/=0) dtset%wfoptalg=10 1936 if (dtset%optdriver==RUNL_GSTATE) then 1937 if (dtset%paral_kgb/=0) dtset%wfoptalg=14 1938 end if 1939 end if 1940 1941 !--------------------------------------------------------------------------- 1942 !Some PAW+DMFT keywords 1943 dtset%usedmft=0 1944 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usedmft',tread,'INT') 1945 if(tread==1) dtset%usedmft=intarr(1) 1946 1947 !Some ucrpa keywords 1948 dtset%ucrpa=0 1949 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'ucrpa',tread,'INT') 1950 if(tread==1) dtset%ucrpa=intarr(1) 1951 1952 if (dtset%ucrpa > 0 .and. dtset%usedmft > 0) then 1953 write(msg, '(9a)' )& 1954 'usedmft and ucrpa are both activated in the input file ',ch10,& 1955 'In the following, abinit assume you are doing a ucrpa calculation and ',ch10,& 1956 'you define Wannier functions as in DFT+DMFT calculation',ch10,& 1957 'If instead, you want to do a full dft+dmft calculation and not only the Wannier construction, use ucrpa=0',ch10,& 1958 'This keywords are depreciated, please use the new keywords to perform cRPA calculation' 1959 ABI_WARNING(msg) 1960 end if 1961 1962 !Some PAW+U keywords 1963 dtset%usepawu=0 1964 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usepawu',tread,'INT') 1965 if(tread==1) dtset%usepawu=intarr(1) 1966 !if(dtset%usedmft>0.and.(dtset%usepawu==14.or.dtset%usepawu==4)) then 1967 ! dtset%usepawu=14 1968 !else if(dtset%usedmft>0.and.dtset%usepawu>=0) then 1969 ! dtset%usepawu=1 1970 !endif 1971 1972 1973 dtset%usedmatpu=0 1974 dtset%lpawu(1:dtset%ntypat)=-1 1975 dtset%optdcmagpawu=3 1976 if (dtset%usepawu/=0.or.dtset%usedmft>0) then 1977 call intagm(dprarr,intarr,jdtset,marr,dtset%ntypat,string(1:lenstr),'lpawu',tread,'INT') 1978 if(tread==1) dtset%lpawu(1:dtset%ntypat)=intarr(1:dtset%ntypat) 1979 1980 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usedmatpu',tread,'INT') 1981 if(tread==1) dtset%usedmatpu=intarr(1) 1982 if (dtset%nspden==4) then 1983 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'optdcmagpawu',tread,'INT') 1984 if(tread==1) dtset%optdcmagpawu=intarr(1) 1985 end if 1986 end if 1987 1988 !Some PAW+Exact exchange keywords 1989 dtset%useexexch=0 1990 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'useexexch',tread,'INT') 1991 if(tread==1) dtset%useexexch=intarr(1) 1992 1993 dtset%lexexch(1:dtset%ntypat)=-1 1994 1995 if (dtset%useexexch/=0) then 1996 call intagm(dprarr,intarr,jdtset,marr,dtset%ntypat,string(1:lenstr),'lexexch',tread,'INT') 1997 if(tread==1) dtset%lexexch(1:dtset%ntypat)=intarr(1:dtset%ntypat) 1998 end if 1999 2000 !LDA minus half keyword 2001 call intagm(dprarr,intarr,jdtset,marr,dtset%ntypat,string(1:lenstr),'ldaminushalf',tread,'INT') 2002 if(tread==1) dtset%ldaminushalf(1:dtset%ntypat)=intarr(1:dtset%ntypat) 2003 2004 !Some plowan data 2005 dtset%plowan_natom=0 2006 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'plowan_natom',tread,'INT') 2007 if(tread==1) dtset%plowan_natom=intarr(1) 2008 2009 dtset%plowan_nt=0 2010 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'plowan_nt',tread,'INT') 2011 if(tread==1) dtset%plowan_natom=intarr(1) 2012 2013 !if (dtset%ucrpa > 0 .and. dtset%plowan_compute==0) then 2014 !dtset%plowan_natom=1 2015 !dtset%plowan_nt=1 2016 !endif 2017 2018 !PAW potential zero keyword 2019 dtset%usepotzero=0 2020 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'usepotzero',tread,'INT') 2021 if(tread==1) dtset%usepotzero=intarr(1) 2022 2023 !Macro_uj (determination of U in PAW+U), governs also allocation of atvshift 2024 dtset%macro_uj = 0 2025 call intagm(dprarr,intarr,jdtset,marr,1,string(1:lenstr),'macro_uj',tread,'INT') 2026 if(tread==1) dtset%macro_uj=intarr(1) 2027 2028 !Constraint DFT keyword 2029 call intagm(dprarr,intarr,jdtset,marr,dtset%ntypat,string(1:lenstr),'constraint_kind',tread,'INT') 2030 if(tread==1) dtset%constraint_kind(1:dtset%ntypat)=intarr(1:dtset%ntypat) 2031 2032 ABI_FREE(nband) 2033 ABI_FREE(ratsph) 2034 ABI_FREE(intarr) 2035 ABI_FREE(dprarr) 2036 2037 !DEBUG 2038 !write(std_out,'(a)')' m_invars1%invars1 : exit ' 2039 !call flush(std_out) 2040 !ENDDEBUG 2041 2042 end subroutine invars1
ABINIT/invars1m [ 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 comm=MPI communicator
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
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) mx<ab_dimensions>=datatype storing the maximal dimensions. Partly initialized in input.
SOURCE
705 subroutine invars1m(dmatpuflag, dtsets, iout, lenstr, mband_upper_, mx,& 706 & msym, ndtset, ndtset_alloc, string, npsp, zionpsp, comm) 707 708 !Arguments ------------------------------------ 709 !scalars 710 integer,intent(in) :: iout,lenstr,msym,ndtset,ndtset_alloc,npsp, comm 711 integer,intent(out) :: dmatpuflag 712 character(len=*),intent(inout) :: string 713 type(ab_dimensions),intent(inout) :: mx 714 !arrays 715 integer,intent(out) :: mband_upper_(0:ndtset_alloc) 716 type(dataset_type),intent(inout) :: dtsets(0:ndtset_alloc) 717 real(dp),intent(in) :: zionpsp(npsp) 718 719 !Local variables------------------------------- 720 !scalars 721 integer :: idtset,ii,jdtset,lpawu,mband_upper,iatom,nat,nsp 722 !arrays 723 integer,allocatable :: symafm_(:,:),symrel_(:,:,:,:) 724 integer,allocatable :: symafm(:),symrel(:,:,:) 725 real(dp),allocatable :: tnons_(:,:,:),tnons(:,:) 726 727 !****************************************************************** 728 729 !DEBUG 730 !write(std_out,'(a)')' m_invars1%invars1m : enter ' 731 !call flush(std_out) 732 !ENDDEBUG 733 734 ! Here, allocation of the arrays that depend on msym. 735 ABI_MALLOC(symrel_,(3,3,msym,0:ndtset_alloc)) 736 ABI_MALLOC(symafm_,(msym,0:ndtset_alloc)) 737 ABI_MALLOC(tnons_,(3,msym,0:ndtset_alloc)) 738 ABI_MALLOC(symafm,(msym)) 739 ABI_MALLOC(symrel,(3,3,msym)) 740 ABI_MALLOC(tnons,(3,msym)) 741 742 ! Set up default values (note that the default acell, amu mkmem, mkmem1,mkqmem, and nkpt must be overcome 743 do idtset=0,ndtset_alloc 744 call indefo1(dtsets(idtset)) 745 end do 746 747 ! natom and nimage are already initialized in invars0 748 dtsets(0)%natom=-1 749 dtsets(0)%nimage=1 750 751 !Initialization for parallelization data has changed 752 !these lines aim to keep old original default values 753 dtsets(0)%npimage=1 754 dtsets(0)%np_spkpt=1 755 dtsets(0)%npspinor=1 756 dtsets(0)%npfft=1 757 dtsets(0)%npband=1 758 dtsets(0)%bandpp=1 759 760 symafm_(:,0)=1 761 symrel_(:,:,:,0)=0 762 symrel_(1,1,:,0)=1 ; symrel_(2,2,:,0)=1 ; symrel_(3,3,:,0)=1 763 tnons_(:,:,0)=0.0_dp 764 765 ! Loop on datasets 766 do idtset=1,ndtset_alloc 767 jdtset=dtsets(idtset)%jdtset ; if(ndtset==0)jdtset=0 768 write(std_out,'(2a)') ch10,'======================================================= ' 769 write(std_out,'(a,i0)') ' invars1m : enter jdtset= ',jdtset 770 771 ! Input default values 772 dtsets(idtset)%bravais(:)=0 773 symafm(:)=symafm_(:,0) 774 symrel(:,:,:)=symrel_(:,:,:,0) 775 tnons(:,:)=tnons_(:,:,0) 776 777 call invars1(dtsets(idtset)%bravais,dtsets(idtset),iout,jdtset,lenstr,& 778 & mband_upper,msym,npsp,string,symafm,symrel,tnons,zionpsp, comm) 779 780 mband_upper_ (idtset)=mband_upper 781 symafm_(:,idtset)=symafm(:) 782 symrel_(:,:,:,idtset)=symrel(:,:,:) 783 tnons_(:,:,idtset)=tnons(:,:) 784 end do 785 786 mx%mband_upper = maxval(mband_upper_ (1:ndtset_alloc)) 787 788 dmatpuflag = 0; mx%natpawu = 0; mx%lpawu = 0 789 mx%natsph = dtsets(1)%natsph 790 mx%natsph_extra = dtsets(1)%natsph_extra 791 mx%natvshift = dtsets(1)%natvshift 792 mx%nconeq = dtsets(1)%nconeq 793 mx%n_efmas_dirs=0 794 mx%ga_n_rules = dtsets(1)%ga_n_rules 795 mx%gw_nqlwl = dtsets(1)%gw_nqlwl 796 mx%nimfrqs = 0 797 mx%nfreqsp = 0 798 mx%n_projection_frequencies = 0 799 mx%nkpt = dtsets(1)%nkpt 800 mx%nkptgw = dtsets(1)%nkptgw 801 mx%nkpthf = dtsets(1)%nkpthf 802 mx%nnos = dtsets(1)%nnos 803 mx%nqptdm = dtsets(1)%nqptdm 804 mx%nspinor = dtsets(1)%nspinor 805 mx%nsppol = dtsets(1)%nsppol 806 mx%ntypat = dtsets(1)%ntypat 807 mx%nzchempot = dtsets(1)%nzchempot 808 mx%nberry = 20 ! This is presently a fixed value. Should be changed. 809 810 ! Get MAX dimension over datasets 811 do ii=1,ndtset_alloc 812 mx%natsph = max(dtsets(ii)%natsph, mx%natsph) 813 mx%natsph_extra=max(dtsets(ii)%natsph_extra, mx%natsph_extra) 814 mx%nconeq=max(dtsets(ii)%nconeq, mx%nconeq) 815 mx%n_efmas_dirs = max(dtsets(ii)%efmas_n_dirs, mx%n_efmas_dirs) 816 mx%ga_n_rules = max(dtsets(ii)%ga_n_rules,mx%ga_n_rules) 817 mx%gw_nqlwl = max(dtsets(ii)%gw_nqlwl,mx%gw_nqlwl) 818 mx%nimfrqs = max(dtsets(ii)%cd_customnimfrqs, mx%nimfrqs) 819 mx%nfreqsp = max(dtsets(ii)%gw_customnfreqsp, mx%nfreqsp) 820 mx%n_projection_frequencies = max(dtsets(ii)%gwls_n_proj_freq, mx%n_projection_frequencies) 821 mx%nkpt = max(dtsets(ii)%nkpt, mx%nkpt) 822 mx%nkptgw = max(dtsets(ii)%nkptgw, mx%nkptgw) 823 mx%nkpthf = max(dtsets(ii)%nkpthf, mx%nkpthf) 824 mx%nnos = max(dtsets(ii)%nnos, mx%nnos) 825 mx%nqptdm = max(dtsets(ii)%nqptdm, mx%nqptdm) 826 mx%nspinor = max(dtsets(ii)%nspinor, mx%nspinor) 827 mx%nsppol = max(dtsets(ii)%nsppol, mx%nsppol) 828 mx%ntypat = max(dtsets(ii)%ntypat, mx%ntypat) 829 mx%nzchempot = max(dtsets(ii)%nzchempot, mx%nzchempot) 830 if (dtsets(ii)%usepawu/=0) then 831 if (dtsets(ii)%usepawu>0.and.dtsets(ii)%usedmatpu/=0) dmatpuflag=1 832 lpawu=maxval(dtsets(ii)%lpawu(:)) 833 mx%lpawu=max(lpawu,mx%lpawu) 834 !dtsets(ii)%natpawu=count(dtsets(ii)%lpawu(dtsets(ii)%typat((/(i1,i1=1,dtsets(ii)%natom)/)))/=-1) 835 ! Old fashon way that should do fine 836 dtsets(ii)%natpawu = 0 837 do iatom=1, dtsets(ii)%natom 838 if (dtsets(ii)%lpawu(dtsets(ii)%typat(iatom)) /= -1 ) dtsets(ii)%natpawu = dtsets(ii)%natpawu + 1 839 end do 840 mx%natpawu = max(dtsets(ii)%natpawu, mx%natpawu) 841 if (dtsets(ii)%macro_uj/=0) dtsets(ii)%natvshift=lpawu*2+1 842 end if 843 mx%natvshift = max(dtsets(ii)%natvshift, mx%natvshift) 844 end do 845 846 !mx%nsym=maxval(dtsets(1:ndtset_alloc)%nsym) ! This might not work properly with HP compiler 847 mx%nsym=dtsets(1)%nsym 848 do idtset=1,ndtset_alloc 849 mx%nsym = max(dtsets(idtset)%nsym, mx%nsym) 850 end do 851 852 do idtset=0,ndtset_alloc 853 ABI_MALLOC(dtsets(idtset)%atvshift, (mx%natvshift, mx%nsppol, mx%natom)) 854 ABI_MALLOC(dtsets(idtset)%bs_loband,(mx%nsppol)) 855 ABI_MALLOC(dtsets(idtset)%bdgw,(2, mx%nkptgw, mx%nsppol)) 856 ABI_MALLOC(dtsets(idtset)%cd_imfrqs,(mx%nimfrqs)) 857 ABI_MALLOC(dtsets(idtset)%chempot,(3, mx%nzchempot, mx%ntypat)) 858 nsp = max(mx%nsppol, mx%nspinor); nat = mx%natpawu*dmatpuflag 859 ABI_MALLOC(dtsets(idtset)%dmatpawu,(2*mx%lpawu+1,2*mx%lpawu+1,nsp,nat, mx%nimage)) 860 ABI_MALLOC(dtsets(idtset)%efmas_bands,(2, mx%nkpt)) 861 ABI_MALLOC(dtsets(idtset)%efmas_dirs,(3, mx%n_efmas_dirs)) 862 ABI_MALLOC(dtsets(idtset)%gw_freqsp, (mx%nfreqsp)) 863 ABI_MALLOC(dtsets(idtset)%gwls_list_proj_freq, (mx%n_projection_frequencies)) 864 ABI_MALLOC(dtsets(idtset)%gw_qlwl,(3,mx%gw_nqlwl)) 865 ABI_MALLOC(dtsets(idtset)%kpt,(3, mx%nkpt)) 866 ABI_MALLOC(dtsets(idtset)%kptgw,(3, mx%nkptgw)) 867 ABI_MALLOC(dtsets(idtset)%kptns,(3, mx%nkpt)) 868 ABI_MALLOC(dtsets(idtset)%kptns_hf,(3, mx%nkpthf)) 869 ABI_MALLOC(dtsets(idtset)%iatsph,(mx%natsph)) 870 ABI_MALLOC(dtsets(idtset)%istwfk, (mx%nkpt)) 871 ABI_MALLOC(dtsets(idtset)%nband, (mx%nkpt*mx%nsppol)) 872 ABI_MALLOC(dtsets(idtset)%occ_orig,(mx%mband_upper*mx%nkpt*mx%nsppol, mx%nimage)) 873 ABI_MALLOC(dtsets(idtset)%qmass, (mx%nnos)) 874 ABI_MALLOC(dtsets(idtset)%qptdm,(3, mx%nqptdm)) 875 ABI_MALLOC(dtsets(idtset)%symafm, (mx%nsym)) 876 ABI_MALLOC(dtsets(idtset)%symrel,(3,3,mx%nsym)) 877 ABI_MALLOC(dtsets(idtset)%tnons,(3,mx%nsym)) 878 ABI_MALLOC(dtsets(idtset)%wtatcon,(3,mx%natom, mx%nconeq)) 879 ABI_MALLOC(dtsets(idtset)%wtk, (mx%nkpt)) 880 ABI_MALLOC(dtsets(idtset)%xredsph_extra,(3, mx%natsph_extra)) 881 dtsets(idtset)%symrel(:,:,:)=symrel_(:,:,1:mx%nsym,idtset) 882 dtsets(idtset)%symafm(:) =symafm_(1:mx%nsym,idtset) 883 dtsets(idtset)%tnons (:,:) =tnons_ (:,1:mx%nsym,idtset) 884 end do 885 886 ABI_FREE(symafm_) 887 ABI_FREE(symrel_) 888 ABI_FREE(tnons_) 889 ABI_FREE(symafm) 890 ABI_FREE(symrel) 891 ABI_FREE(tnons) 892 893 !DEBUG 894 !write(std_out,'(a)')' m_invars1%invars1m : exit ' 895 !call flush(std_out) 896 !ENDDEBUG 897 898 end subroutine invars1m
ABINIT/m_invars1 [ Modules ]
NAME
m_invars1
FUNCTION
COPYRIGHT
Copyright (C) 1998-2022 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 .
SOURCE
17 #if defined HAVE_CONFIG_H 18 #include "config.h" 19 #endif 20 21 #include "abi_common.h" 22 23 module m_invars1 24 25 use defs_basis 26 use m_abicore 27 use m_xmpi 28 use m_errors 29 use m_atomdata 30 use m_dtset 31 use m_nctk 32 #ifdef HAVE_NETCDF 33 use netcdf 34 #endif 35 36 use m_fstrings, only : inupper, itoa, endswith, strcat, sjoin, startswith 37 use m_geometry, only : mkrdim 38 use m_parser, only : intagm, intagm_img, chkint_ge, ab_dimensions, geo_t, geo_from_abivar_string 39 use m_inkpts, only : inkpts, inqpt 40 use m_ingeo, only : ingeo, invacuum 41 use m_symtk, only : mati3det 42 43 #if defined HAVE_GPU_CUDA 44 use m_gpu_toolbox 45 #endif 46 47 implicit none 48 49 private