1!! NAME 2!! input_dataset_mod 3!! 4!! FUNCTION 5!! This module contains objects and routines used to parse the atompaw 6!! input file. The complete input dataset is stored in the `input_dataset` 7!! Fortran data-structure. 8 9#if defined HAVE_CONFIG_H 10#include "config.h" 11#endif 12 13MODULE input_dataset_mod 14 15 USE io_tools 16 USE tools 17 18 IMPLICIT NONE 19 20 PRIVATE 21 22!Public functions 23 PUBLIC :: input_dataset_read ! Initialize input_dataset by reading it from file 24 PUBLIC :: input_dataset_free ! Destroy input_dataset data-structure 25 PUBLIC :: input_dataset_copy ! Copy a input data-structure into another 26 PUBLIC :: input_dataset_read_occ ! Read and change modified occupations in input_dataset 27 PUBLIC :: input_dataset_read_abinit ! Read ABINIT options and store them in input_dataset 28 PUBLIC :: input_dataset_read_xml ! Read XML options and store them in input_dataset 29 PUBLIC :: input_dataset_read_upf ! Read UPF options and store them in input_dataset 30 31!Public structured datatype 32!All data from input file 33 TYPE,PUBLIC :: input_dataset_t 34 CHARACTER(2) :: atomic_symbol ! Atomic symbol 35 INTEGER :: atomic_charge ! Atomic charge 36 LOGICAL :: scalarrelativistic ! Flag activating the scalar relativistic scheme 37 LOGICAL :: diracrelativistic ! Flag activating the Dirac relativistic scheme 38 LOGICAL :: finitenucleus ! Flag activating finite nucleus model 39 INTEGER :: finitenucleusmodel ! Option for the finite nucleus model 40 LOGICAL :: HFpostprocess ! Option for the post-processing of Hartree-Fock 41 LOGICAL :: localizedcoreexchange ! Flag activating the localization of the core exchange (HF) 42 CHARACTER(10) :: gridkey ! Type of radial grid 43 REAL(8) :: gridrange ! Range of the radial grid 44 REAL(8) :: gridmatch ! A matching radius in the radial grid 45 INTEGER :: gridpoints ! Number of points of the radial grid 46 INTEGER :: nlogderiv ! Number of points for the logarithmic derivative calculation 47 REAL(8) :: minlogderiv ! Minimum energy value for the logarithmic derivative calculation 48 REAL(8) :: maxlogderiv ! Maximum energy value for the logarithmic derivative calculation 49 LOGICAL :: BDsolve ! Flag activating the Block-Davidson solver 50 LOGICAL :: fixed_zero ! Flag activating the "fixed zero" exact exchange potential calculation 51 INTEGER :: fixed_zero_index ! Option for "fixed zero" calculation in the exact exchange potential 52 CHARACTER(132) :: exctype ! Exchange-correlation type (string) 53 INTEGER :: np(5) ! Electronic configuration: number of s,p,d,f,g shells 54 INTEGER :: norbit ! Electronic configuration: number of orbitals 55 LOGICAL,ALLOCATABLE :: orbit_iscore(:) ! Electronic configuration: TRUE for the core orbitals 56 INTEGER :: norbit_mod ! Electronic configuration: number of orbitals with modified occupations 57 INTEGER,ALLOCATABLE :: orbit_mod_n(:) ! Electronic config.: n number of the modified orbital 58 INTEGER,ALLOCATABLE :: orbit_mod_l(:) ! Electronic config.: l number of the modified orbital 59 INTEGER,ALLOCATABLE :: orbit_mod_k(:) ! Electronic config.: kappa number of the modified orbital 60 REAL(8),ALLOCATABLE :: orbit_mod_occ(:) ! Electronic config.: occupation of the modified orbital 61 INTEGER :: norbit_val ! Electronic configuration: number of valence orbitals 62 INTEGER,ALLOCATABLE :: orbit_val_n(:) ! Electronic config.: n number of the valence orbital 63 INTEGER,ALLOCATABLE :: orbit_val_l(:) ! Electronic config.: l number of the valence orbital 64 INTEGER,ALLOCATABLE :: orbit_val_k(:) ! Electronic config.: kappa number of the valence orbital 65 INTEGER :: lmax=-1 ! PAW Basis: maximum l value 66 REAL(8) :: rc=0.d0 ! PAW basis: cut-off radius for the augmentation regions 67 REAL(8) :: rc_shap=0.d0 ! PAW basis: cut-off radius of the compensation charge shape function 68 REAL(8) :: rc_vloc=0.d0 ! PAW basis: matching radius for the local potential 69 REAL(8) :: rc_core=0.d0 ! PAW basis: matching radius for the pseudo-core density 70 INTEGER :: nbasis ! PAW basis : number of basis functions 71 INTEGER :: nbasis_add ! PAW basis: number of additional basis functions (unbound states) 72 INTEGER,ALLOCATABLE :: basis_add_l(:) ! PAW basis: l number for the additional basis func. 73 INTEGER,ALLOCATABLE :: basis_add_k(:) ! PAW basis: kappa number for the additional basis func. 74 REAL(8),ALLOCATABLE :: basis_add_energy(:) ! PAW basis: ref. energy for the additional basis func. 75 REAL(8),ALLOCATABLE :: basis_func_rc(:) ! PAW basis: matching radii for the pseudo partial waves 76 INTEGER :: projector_type ! Type of projectors (Bloechl, Vanderbilt, ...) 77 INTEGER :: pseudo_type ! Type of pseudization scheme (Bessel, polynom, ...) 78 INTEGER :: ortho_type ! Type of orthogonalization scheme (Gram-Schmidt, ...) 79 INTEGER :: pseudo_polynom2_pdeg ! Polynom2 projectors: degree of the polynom 80 REAL(8) :: pseudo_polynom2_qcut ! Polynom2 projectors: q-value for Fourier filtering 81 INTEGER :: shapefunc_type ! Compensation shape function type (sinc2, gaussian, ...) 82 REAL(8) :: shapefunc_gaussian_param ! Compensation shape function: parameter for gaussian type 83 REAL(8) :: hf_coretol ! Tolerance for core density (Hartree-Fock only) 84 INTEGER :: vloc_type ! Type of local potential pseudization 85 INTEGER :: vloc_l ! Local potential: l quantum number (MTrouillier, Ultrasoft) 86 REAL(8) :: vloc_ene ! Local potential: reference energy (MTrouillier, Ultrasoft) 87 REAL(8) :: vloc_setvloc_coef ! "SetVloc" local potential: coefficient 88 REAL(8) :: vloc_setvloc_rad ! "SetVloc" local potential: radius 89 INTEGER :: vloc_kerker_power(4) ! "Kerker" locazl potential: polynomial powers 90 LOGICAL :: abinit_usexcnhat ! ABINIT option: TRUE if nhat density (compensation) is included in XC 91 LOGICAL :: abinit_prtcorewf ! ABINIT option: TRUE is printing of core WF is required 92 CHARACTER(132) :: abinit_author ! ABINIT option: author to be printed in the PAW dataset file 93 LOGICAL :: abinit_uselog ! ABINIT option: TRUE if data are transfered on a log. grid 94 INTEGER :: abinit_log_meshsz ! ABINIT option: mesh size for the logarithmic grid 95 REAL(8) :: abinit_log_step ! ABINIT option: logarithmic step for the logarithmic grid 96 LOGICAL :: abinit_userso ! ABINIT option: TRUE if REAL Space Optimization is required 97 REAL(8) :: abinit_rso_ecut ! ABINIT option: Real Space Optimization parameter - plane wave cutoff 98 REAL(8) :: abinit_rso_gfact ! ABINIT option: Real Space Optimization parameter - Gamma/Gmax ratio 99 REAL(8) :: abinit_rso_werror ! ABINIT option: Real Space Optimization parameter: max. error 100 LOGICAL :: xml_usexcnhat ! XML option: TRUE if nhat density (compensation) is included in XC 101 LOGICAL :: xml_prtcorewf ! XML option: TRUE is printing of core WF is required 102 CHARACTER(132) :: xml_author ! XML option: author to be printed in the PAW dataset file 103 CHARACTER(132) :: xml_comment ! XML option: comment to be printed in the header of PAW dataset file 104 LOGICAL :: xml_usespl ! XML option: TRUE if data are interpolated on a log. grid 105 INTEGER :: xml_spl_meshsz ! XML option: mesh size for the reduced grid 106 LOGICAL :: xml_userso ! XML option: TRUE if REAL Space Optimization is required 107 REAL(8) :: xml_rso_ecut ! XML option: Real Space Optimization parameter - plane wave cutoff 108 REAL(8) :: xml_rso_gfact ! XML option: Real Space Optimization parameter - Gamma/Gmax ratio 109 REAL(8) :: xml_rso_werror ! XML option: Real Space Optimization parameter: max. error 110 LOGICAL :: xml_uselda12 ! XML option: TRUE if LDA-1/2 potential calculation is required 111 INTEGER :: xml_lda12_orb_l ! XML option: LDA-1/2 parameter - l quantum number of the ionized orbital 112 INTEGER :: xml_lda12_orb_n ! XML option: LDA-1/2 parameter - n quantum number of the ionized orbital 113 REAL(8) :: xml_lda12_ion ! XML option: LDA-1/2 parameter - amount of charge removed from the ionized orbital 114 REAL(8) :: xml_lda12_rcut ! XML option: LDA-1/2 parameter - cut-off radius (in bohr) 115 CHARACTER(15) :: xml_lda12_logfile ! XML option: LDA-1/2 parameter - name of the log file 116 REAL(8) :: upf_grid_xmin ! UPF option: minimum radius given by the grid 117 REAL(8) :: upf_grid_zmesh ! UPF option: inverse of the rdial step of the grid 118 REAL(8) :: upf_grid_dx ! UPF option: logarithmic step of the grid 119 REAL(8) :: upf_grid_range ! UPF option: range of the grid 120 END TYPE input_dataset_t 121 122!Public variable containing the complete input dataset 123 TYPE(input_dataset_t),PUBLIC,SAVE,TARGET :: input_dataset 124 125!Public parameters 126 INTEGER,PARAMETER,PUBLIC :: PROJECTOR_TYPE_BLOECHL = 1 127 INTEGER,PARAMETER,PUBLIC :: PROJECTOR_TYPE_VANDERBILT= 2 128 INTEGER,PARAMETER,PUBLIC :: PROJECTOR_TYPE_CUSTOM = 3 129 INTEGER,PARAMETER,PUBLIC :: PROJECTOR_TYPE_MODRRKJ = 4 130 INTEGER,PARAMETER,PUBLIC :: PROJECTOR_TYPE_HF = 5 131 INTEGER,PARAMETER,PUBLIC :: PSEUDO_TYPE_BLOECHL = 1 132 INTEGER,PARAMETER,PUBLIC :: PSEUDO_TYPE_POLYNOM = 2 133 INTEGER,PARAMETER,PUBLIC :: PSEUDO_TYPE_POLYNOM2 = 3 134 INTEGER,PARAMETER,PUBLIC :: PSEUDO_TYPE_RRKJ = 4 135 INTEGER,PARAMETER,PUBLIC :: PSEUDO_TYPE_BLOECHL_K = 5 136 INTEGER,PARAMETER,PUBLIC :: PSEUDO_TYPE_HF = 6 137 INTEGER,PARAMETER,PUBLIC :: ORTHO_TYPE_GRAMSCHMIDT = 1 138 INTEGER,PARAMETER,PUBLIC :: ORTHO_TYPE_VANDERBILT = 2 139 INTEGER,PARAMETER,PUBLIC :: ORTHO_TYPE_SVD = 3 140 INTEGER,PARAMETER,PUBLIC :: ORTHO_TYPE_HF = 4 141 INTEGER,PARAMETER,PUBLIC :: SHAPEFUNC_TYPE_GAUSSIAN = 1 142 INTEGER,PARAMETER,PUBLIC :: SHAPEFUNC_TYPE_SINC = 2 143 INTEGER,PARAMETER,PUBLIC :: SHAPEFUNC_TYPE_BESSEL = 3 144 INTEGER,PARAMETER,PUBLIC :: VLOC_TYPE_MTROULLIER = 1 145 INTEGER,PARAMETER,PUBLIC :: VLOC_TYPE_ULTRASOFT = 2 146 INTEGER,PARAMETER,PUBLIC :: VLOC_TYPE_BESSEL = 3 147 INTEGER,PARAMETER,PUBLIC :: VLOC_TYPE_SETVLOC = 4 148 INTEGER,PARAMETER,PUBLIC :: VLOC_TYPE_KERKER_EXPF = 5 149 INTEGER,PARAMETER,PUBLIC :: VLOC_TYPE_KERKER_POLY = 6 150 INTEGER,PARAMETER,PUBLIC :: UNKNOWN_TYPE =-1 151 152 !Default values 153 LOGICAL,PARAMETER,PRIVATE :: PRTCOREWF_DEF =.false. 154 LOGICAL,PARAMETER,PRIVATE :: USEXCNHAT_DEF =.false. 155 LOGICAL,PARAMETER,PRIVATE :: USELOG_DEF =.false. 156 LOGICAL,PARAMETER,PRIVATE :: USERSO_DEF =.false. 157 LOGICAL,PARAMETER,PRIVATE :: USELDA12_DEF =.false. 158 INTEGER,PARAMETER,PRIVATE :: LOGGRD_SIZE_DEF =350 159 REAL(8),PARAMETER,PRIVATE :: LOGGRD_STEP_DEF =0.035d0 160 REAL(8),PARAMETER,PRIVATE :: RSO_ECUT_DEF =10.0d0 161 REAL(8),PARAMETER,PRIVATE :: RSO_GFACT_DEF =2.d0 162 REAL(8),PARAMETER,PRIVATE :: RSO_WERROR_DEF =0.0001d0 163 INTEGER,PARAMETER,PRIVATE :: LDA12_ORB_L_DEF =-1 164 INTEGER,PARAMETER,PRIVATE :: LDA12_ORB_N_DEF =-1 165 REAL(8),PARAMETER,PRIVATE :: LDA12_ION_DEF =0.d0 166 REAL(8),PARAMETER,PRIVATE :: LDA12_RCUT_DEF =0.d0 167 CHARACTER(15),PARAMETER,PRIVATE :: lda12_logfile='lda-12.log' 168 REAL(8),PARAMETER,PRIVATE :: UPF_DX_DEF =0.005d0 169 REAL(8),PARAMETER,PRIVATE :: UPF_XMIN_DEF =-9.d0 170 REAL(8),PARAMETER,PRIVATE :: UPF_ZMESH_DEF =1.d0 171 REAL(8),PARAMETER,PRIVATE :: UPF_RANGE_DEF =15.d0 172 CHARACTER(132),PARAMETER,PRIVATE :: AUTHOR_DEF ='', COMMENT_DEF='' 173 174!Private parameters 175 INTEGER,PARAMETER,PRIVATE :: norbit_max=50,nbasis_add_max=25 176 REAL(8),PARAMETER,PRIVATE :: linrange=50.d0,mxgridlin=20001 177 REAL(8),PARAMETER,PRIVATE :: logrange=80.d0,v4logrange=100.d0,mxgridlog=2001 178 REAL(8),PARAMETER,PRIVATE :: logder_min=-5.d0,logder_max=4.95d0,logder_pts=200 179 REAL(8),PARAMETER,PRIVATE :: polynom2_pdeg_def=4 180 REAL(8),PARAMETER,PRIVATE :: polynom2_qcut_def=10.d0 181 REAL(8),PARAMETER,PRIVATE :: gausstol_def=1.d-4 182 REAL(8),PARAMETER,PRIVATE :: hf_coretol_def=1.d-4 183 184!Private variables 185 LOGICAL,PRIVATE,SAVE :: has_to_print=.true. 186 LOGICAL,PRIVATE,SAVE :: has_to_ask=.true. 187 188CONTAINS 189 190!!================================================================= 191!! NAME 192!! input_dataset_read 193!! 194!! FUNCTION 195!! Initialize an input_dataset datastructure by reading it from 196!! a file. If file is omitted, then read from standard input. 197!! Note: we only read here data used to generate the PAW dataset, 198!! not data used for the post-processing (output, explore, scfpaw, ...) 199!! 200!! INPUTS (all optionals) 201!! [inputfile]= name of input file to be read 202!! [echofile]= name of a file to echo input file content 203!! [read_global_data]= if TRUE, read global data (atom, XC, grid, ...) - Default TRUE 204!! [read_elec_data]= if TRUE, read electronic configuration (orbital & occupations) - Default TRUE 205!! [read_coreval_data]= if TRUE, read electronic config (core and valence) - Default TRUE 206!! [read_basis_data]= if TRUE, read basis data (radii, pseudo scheme, ...) - Default TRUE 207!! 208!! OUTPUT 209!! [input_dt]= datastructure containing the complete input file. 210!! If omitted, then the global public `input_dataset` 211!! is used. 212!! 213!!================================================================= 214 SUBROUTINE input_dataset_read(input_dt,inputfile,echofile,& 215& read_global_data,read_elec_data,read_coreval_data,read_basis_data) 216 217!---- Arguments 218 CHARACTER*(*),INTENT(IN),OPTIONAL :: inputfile,echofile 219 LOGICAL,INTENT(IN),OPTIONAL :: read_global_data,read_elec_data,& 220& read_coreval_data,read_basis_data 221 TYPE(input_dataset_t),INTENT(INOUT),OPTIONAL,TARGET :: input_dt 222 223!---- Local variables 224 INTEGER,PARAMETER :: ifunit=111,ecunit=222 225 INTEGER,PARAMETER :: nkappa(5)=(/1,2,2,2,2/) 226 INTEGER :: input_unit 227 INTEGER :: ii,io,nadd,norb,nval,nbl,nn,ik,kk 228 INTEGER :: ilin,ilog,inrl,iscl,ipnt,ifin,iend,ihfpp,ilcex 229 INTEGER :: igrid,irelat,ilogder,ilogv4,ibd,idirac,ifixz,ll,nstart 230 LOGICAL :: has_to_echo 231 LOGICAL :: read_global_data_,read_elec_data_,read_coreval_data_,read_basis_data_ 232 CHARACTER(200) :: inputline,inputword 233 CHARACTER(128) :: exchangecorrelationandgridline 234 CHARACTER(1) :: CHR 235 REAL(8) :: x1,x2 236 INTEGER :: basis_add_l(nbasis_add_max) 237 INTEGER :: basis_add_k(nbasis_add_max) 238 REAL(8) :: basis_add_energy(nbasis_add_max) 239 TYPE(input_dataset_t),POINTER :: dataset 240 241!------------------------------------------------------------------ 242 243!Select datastruture to be read 244 IF (PRESENT(input_dt)) THEN 245 dataset => input_dt 246 ELSE 247 dataset => input_dataset 248 ENDIF 249 250!Select input file logical unit 251 IF (PRESENT(inputfile)) THEN 252 input_unit=ifunit 253 OPEN(ifunit,file=trim(inputfile),form='formatted') 254 ELSE 255 input_unit=STD_IN 256 END IF 257 258!Check if we read a TTY or a FILE 259 has_to_ask=unit_isatty(input_unit) 260 261!Do we echo input file content? 262 has_to_echo=PRESENT(echofile) 263 IF (has_to_echo) THEN 264 OPEN(ecunit,file=trim(echofile),form='formatted') 265 END IF 266 267!Select which components have to be read 268 read_global_data_=.true.;if (PRESENT(read_global_data)) read_global_data_=read_global_data 269 read_elec_data_=.true.;if (PRESENT(read_elec_data)) read_elec_data_=read_elec_data 270 read_coreval_data_=.true.;if (PRESENT(read_coreval_data)) read_coreval_data_=read_coreval_data 271 read_basis_data_=.true.;if (PRESENT(read_basis_data)) read_basis_data_=read_basis_data 272 273!Print a title 274 IF (read_global_data_.OR.read_elec_data_.OR.read_coreval_data_.OR.read_basis_data_) THEN 275 IF (has_to_ask) THEN 276 WRITE(STD_OUT,'(/,1x,a)') "===== READING OF INPUT DATA =====" 277 ELSE 278 WRITE(STD_OUT,'(/,3x,a)') "===== READING OF INPUT FILE =====" 279 END IF 280 END IF 281 282 283!------------------------------------------------------------------ 284!Start reading of AE data 285 IF (read_global_data_) THEN 286 287!------------------------------------------------------------------ 288!=== 1st line: read atomic symbol, atomic number 289 290 IF (has_to_ask) THEN 291 WRITE(STD_OUT,*) 'Enter atomic symbol and atomic number' 292 END IF 293 294 READ(input_unit,'(a)') inputline 295 IF (has_to_echo) WRITE(ecunit,'(a)') TRIM(inputline) 296 CALL eliminate_comment(inputline) 297 298 READ(inputline,*) dataset%atomic_symbol,dataset%atomic_charge 299 300!Print read data 301 IF (has_to_print) THEN 302 WRITE(STD_OUT,'(3x,a,a2)') "Atomic symbol : ",dataset%atomic_symbol 303 WRITE(STD_OUT,'(3x,a,i0)') "Atomic charge : ",dataset%atomic_charge 304 END IF 305 306!------------------------------------------------------------------ 307!=== 2nd line: read XC type, grid data, relativistic,point-nucleus, 308! logderiv data, HF data, Block-Davidson keyword 309 310 IF (has_to_ask) THEN 311 WRITE(STD_OUT,*) 'Enter exchange-correlation type, among the following:' 312 WRITE(STD_OUT,*) ' * LDA-PW (default), GGA-PBE, GGA-PBESOL' 313 WRITE(STD_OUT,*) ' * LibXC keyword beginning with XC_' 314 WRITE(STD_OUT,*) ' * EXX, EXXOCC, EXXKLI, EXXCS' 315 WRITE(STD_OUT,*) ' * HF, HFV' 316 WRITE(STD_OUT,*) ' further optionally (space) "nonrelativistic/scalarrelativistic/diracrelativistic" keyword' 317 WRITE(STD_OUT,*) ' further optionally (space) "point-nucleus/finite-nucleus" keyword' 318 WRITE(STD_OUT,*) ' further optionally (space) "loggrid/lineargrid" keyword if appropriate' 319 WRITE(STD_OUT,*) ' Note: "loggridv4" puts more points near origin' 320 WRITE(STD_OUT,*) ' further optionally n (number of grid points)' 321 WRITE(STD_OUT,*) ' r_max (max. grid radius)' 322 WRITE(STD_OUT,*) ' r_match (exact value of r(n))' 323 WRITE(STD_OUT,*) ' further optionally (space) "logderivrange" keyword' 324 WRITE(STD_OUT,*) ' further optionally emin (minimum energy for log. deriv. plot)' 325 WRITE(STD_OUT,*) ' emax (maximum energy for log. deriv. plot)' 326 WRITE(STD_OUT,*) ' ne (# of energies for log. deriv. plot)' 327 WRITE(STD_OUT,*) ' and optionally "BDsolve" keyword for Block-Davidson solver' 328 END IF 329 330!Read full line 331 READ(input_unit,'(a)') exchangecorrelationandgridline 332 IF (has_to_echo) WRITE(ecunit,'(a)') TRIM(exchangecorrelationandgridline) 333 CALL eliminate_comment(inputline) 334 335 CALL Uppercase(exchangecorrelationandgridline) 336 exchangecorrelationandgridline=trim(exchangecorrelationandgridline) 337 338!Retrieve keyword indexes 339 ilin=0;ilin=0;ilog=0;ilogv4=0;inrl=0;iscl=0;ipnt=0;ifin=0 340 ihfpp=0;ilcex=0;igrid=0;irelat=0;ilogder=0;ibd=0;idirac=0 341 ilin=INDEX(exchangecorrelationandgridline,'LINEARGRID') 342 ilog=INDEX(exchangecorrelationandgridline,'LOGGRID') 343 ilogv4=INDEX(exchangecorrelationandgridline,'LOGGRIDV4') 344 ibd=INDEX(exchangecorrelationandgridline,'BDSOLVE') 345 inrl=INDEX(exchangecorrelationandgridline,'NONRELATIVISTIC') 346 iscl=INDEX(exchangecorrelationandgridline,'SCALARRELATIVISTIC') 347 idirac=INDEX(exchangecorrelationandgridline,'DIRACRELATIVISTIC') 348 ipnt=INDEX(exchangecorrelationandgridline,'POINT-NUCLEUS') 349 ifin=INDEX(exchangecorrelationandgridline,'FINITE-NUCLEUS') 350 ilogder=INDEX(exchangecorrelationandgridline,'LOGDERIVRANGE') 351 ihfpp=INDEX(exchangecorrelationandgridline,'HFPOSTPROCESS') 352 ilcex=INDEX(exchangecorrelationandgridline,'LOCALIZEDCOREEXCHANGE') 353 ifixz=INDEX(exchangecorrelationandgridline,'FIXED_ZERO') 354 igrid=max(ilin,ilog) !This line may need attention.... 355 irelat=max(inrl,iscl) !This line may need attention.... 356 357!Treat simple logical variables 358 dataset%scalarrelativistic=(iscl>0.and.inrl==0) 359 dataset%diracrelativistic=(idirac>0.and.inrl==0) 360 dataset%finitenucleus=(ifin>0.and.ipnt==0) 361 dataset%BDsolve=(ibd>0) 362 dataset%HFpostprocess=(ihfpp>0) 363 364!Treat finite nucleus option 365 dataset%finitenucleusmodel=-1 366 IF (dataset%finitenucleus) THEN 367 READ(exchangecorrelationandgridline(ifin+14:ifin+14),'(a)') CHR 368 IF (CHR=="2") dataset%finitenucleusmodel=2 369 IF (CHR=="3") dataset%finitenucleusmodel=3 370 IF (CHR=="4") dataset%finitenucleusmodel=4 371 IF (CHR=="5") dataset%finitenucleusmodel=5 372 END IF 373 374!Treat grid data 375 dataset%gridkey='LINEAR' 376 dataset%gridpoints=mxgridlin 377 dataset%gridrange=linrange 378 dataset%gridmatch=linrange 379 IF (ilog>0.and.ilin==0.and.ilogv4==0) THEN 380 dataset%gridkey='LOGGRID' 381 dataset%gridpoints=mxgridlog; 382 dataset%gridrange=logrange 383 dataset%gridmatch=logrange 384 END IF 385 IF (ilog>0.and.ilin==0.and.ilogv4>0) THEN 386 dataset%gridkey='LOGGRID4' 387 dataset%gridpoints=mxgridlog; 388 dataset%gridrange=v4logrange 389 dataset%gridmatch=v4logrange 390 END IF 391 IF (igrid>0) THEN 392 iend=128 393 IF (irelat >igrid.and.irelat-1 <iend) iend=irelat -1 394 IF (ilogder>igrid.and.ilogder-1<iend) iend=ilogder-1 395 IF (ibd>igrid.and.ibd-1<iend) iend=ibd-1 396 inputline="" 397 IF (ilog>0.and.ilogv4==0.and.iend>igrid+7) & 398& inputline=TRIM(exchangecorrelationandgridline(igrid+7:iend)) 399 IF (ilog>0.and.ilogv4>0.and.iend>igrid+9) & 400& inputline=TRIM(exchangecorrelationandgridline(igrid+9:iend)) 401 IF (ilin>0.and.iend>igrid+10) & 402& inputline=TRIM(exchangecorrelationandgridline(igrid+10:iend)) 403 IF (inputline/="") THEN 404 CALL extractword(1,inputline,inputword);inputword=trim(inputword) 405 IF (inputword/="") THEN 406 READ(inputword,*) dataset%gridpoints 407 CALL extractword(2,inputline,inputword);inputword=trim(inputword) 408 IF (inputword/="") THEN 409 READ(inputword,*) dataset%gridrange 410 dataset%gridmatch=dataset%gridrange 411 CALL extractword(3,inputline,inputword);inputword=trim(inputword) 412 IF (inputword/="") read(inputword,*) dataset%gridmatch 413 END IF 414 END IF 415 END IF 416 IF (dataset%gridpoints<=0) STOP "input_dataset: error -- number of grid points should be >0!" 417 END IF 418 419!Treat logderiv data 420 dataset%minlogderiv=logder_min 421 dataset%maxlogderiv=logder_max 422 dataset%nlogderiv=logder_pts 423 IF (ilogder>0) THEN 424 iend=128 425 IF (igrid >ilogder.and.igrid-1 <iend) iend=igrid -1 426 IF (irelat>ilogder.and.irelat-1<iend) iend=irelat-1 427 inputline="" 428 IF (iend>ilogder+13) inputline=trim(exchangecorrelationandgridline(ilogder+13:iend)) 429 IF (inputline/="") THEN 430 CALL extractword(1,inputline,inputword);inputword=trim(inputword) 431 IF (inputword/="") THEN 432 READ(inputword,*) dataset%minlogderiv 433 CALL extractword(2,inputline,inputword);inputword=trim(inputword) 434 IF (inputword/="") THEN 435 READ(inputword,*) dataset%maxlogderiv 436 CALL extractword(3,inputline,inputword);inputword=trim(inputword) 437 IF (inputword/="") READ(inputword,*) dataset%nlogderiv 438 END IF 439 END IF 440 END IF 441 END IF 442 443!Treat XC/HF 444 READ(unit=exchangecorrelationandgridline,fmt=*) dataset%exctype 445 dataset%localizedcoreexchange=(ilcex>0) 446 dataset%fixed_zero=(ifixz>0) ; dataset%fixed_zero_index=-1 447 IF (dataset%fixed_zero) & 448& READ(unit=exchangecorrelationandgridline(ifixz+10:),fmt=*) dataset%fixed_zero_index 449 450!Print read data 451 IF (has_to_print) THEN 452 WRITE(STD_OUT,'(3x,2a)') "Scalar-relativistic calculation : ",MERGE("YES"," NO",dataset%scalarrelativistic) 453 WRITE(STD_OUT,'(3x,2a)') "Dirac-relativistic calculation : ",MERGE("YES"," NO",dataset%diracrelativistic) 454 WRITE(STD_OUT,'(3x,2a)') "Finite-nucleus calculation : ",MERGE("YES"," NO",dataset%finitenucleus) 455 IF (dataset%finitenucleus) THEN 456 WRITE(STD_OUT,'(3x,a,i0)') " - Finite-nucleus model : ",dataset%finitenucleusmodel 457 END IF 458 WRITE(STD_OUT,'(3x,2a)') "Block-Davidson calculation : ",MERGE("YES"," NO",dataset%BDsolve) 459 WRITE(STD_OUT,'(3x,2a)') "Grid type : ",TRIM(dataset%gridkey) 460 WRITE(STD_OUT,'(3x,a,i0)') "Grid size : ",dataset%gridpoints 461 WRITE(STD_OUT,'(3x,a,f7.3)') "Grid maximum value : ",dataset%gridrange 462 WRITE(STD_OUT,'(3x,a,f7.3)') "Grid imposed value : ",dataset%gridmatch 463 WRITE(STD_OUT,'(3x,a,i0)') "Log. derivative, number of pts : ",dataset%nlogderiv 464 WRITE(STD_OUT,'(3x,a,f7.3)') "Log. derivative, min. energy : ",dataset%minlogderiv 465 WRITE(STD_OUT,'(3x,a,f7.3)') "Log. derivative, max. energy : ",dataset%maxlogderiv 466 WRITE(STD_OUT,'(3x,2a)') "Hartree-Fock, post-processing : ",MERGE("YES"," NO",dataset%HFpostprocess) 467 WRITE(STD_OUT,'(3x,2a)') "Hartree-Fock, localized core ex.: ",MERGE("YES"," NO",dataset%localizedcoreexchange) 468 WRITE(STD_OUT,'(3x,2a)') "Hartree-Fock, fixed zero : ",MERGE("YES"," NO",dataset%fixed_zero) 469 IF (dataset%fixed_zero) THEN 470 WRITE(STD_OUT,'(3x,a,i0)') " - HF fixed zero index : ",dataset%fixed_zero_index 471 END IF 472 IF (dataset%BDsolve.and.dataset%gridkey=='LINEAR') THEN 473 WRITE(STD_OUT,'(/,3x,a)') "WARNING: BlockDavidson solver works very slowly with linear grid!" 474 END IF 475END IF 476 477 478!------------------------------------------------------------------ 479!End reading of global data. Start reading of electronic configuration data 480 ENDIF 481 IF (read_elec_data_) THEN 482 483 484!------------------------------------------------------------------ 485!=== 3rd line and following: electronic configuration of atom 486 487 IF (has_to_ask) THEN 488 WRITE(STD_OUT,*) 'Enter maximum principal quantum numbers for s,p,d,f,g' 489 END IF 490 491 READ(input_unit,'(a)') inputline 492 IF (has_to_echo) WRITE(ecunit,'(a)') TRIM(inputline) 493 CALL eliminate_comment(inputline) 494 495 READ(inputline,*) dataset%np(1:5) 496 497 DO ll=1,5 498 IF(dataset%np(ll)<0) dataset%np(ll)=0 499 END DO 500 dataset%norbit=dataset%np(1)+max(dataset%np(2)-1,0)+max(dataset%np(3)-2,0) & 501& +max(dataset%np(4)-3,0)+max(dataset%np(5)-4,0) 502 IF (dataset%diracrelativistic) dataset%norbit=2*dataset%norbit-dataset%np(1) 503 504!Print read data 505 IF (has_to_print) THEN 506 WRITE(STD_OUT,'(3x,a,5(1x,i0))') "Max. quantum numbers (s,p,d,f,g): ",dataset%np(1:5) 507 WRITE(STD_OUT,'(3x,a,i0)') "Total number of orbitals: ",dataset%norbit 508 END IF 509 510 CALL input_dataset_read_occ(dataset%norbit_mod,dataset%orbit_mod_l,& 511& dataset%orbit_mod_n,dataset%orbit_mod_k,dataset%orbit_mod_occ,& 512& dataset%np,dataset%diracrelativistic,& 513& inputfile_unit=input_unit,echofile_unit=ecunit) 514 515 516!------------------------------------------------------------------ 517!End reading of electronic data. Start reading of core/valence data 518 ENDIF 519 IF (read_coreval_data_) THEN 520 521 522!------------------------------------------------------------------ 523!=== Core and valence states 524 525 IF (has_to_ask) THEN 526 WRITE(STD_OUT,*) 'For each state enter c for core or v for valence' 527 END IF 528 529!Read core and valence states 530 IF (ALLOCATED(dataset%orbit_iscore)) DEALLOCATE(dataset%orbit_iscore) 531 ALLOCATE(dataset%orbit_iscore(dataset%norbit)) 532 DO io=1,dataset%norbit 533 DO 534 READ(input_unit,'(a)') inputline 535 CALL eliminate_comment(inputline) 536 READ(inputline,*) CHR 537 IF (CHR=='c'.OR.CHR=='C'.OR.& 538& CHR=='v'.OR.CHR=='V') THEN 539 IF (has_to_echo) WRITE(ecunit,'(a)') TRIM(inputline) 540 EXIT 541 ELSE 542 WRITE(STD_OUT,*) ' >> Please input c or v!' 543 END IF 544 END DO 545 dataset%orbit_iscore(io)=(CHR=='c'.OR.CHR=='C') 546 END DO 547 548!Store valence states 549 dataset%norbit_val=dataset%norbit-COUNT(dataset%orbit_iscore(:)) 550 IF (ALLOCATED(dataset%orbit_val_n)) DEALLOCATE(dataset%orbit_val_n) 551 IF (ALLOCATED(dataset%orbit_val_l)) DEALLOCATE(dataset%orbit_val_l) 552 IF (ALLOCATED(dataset%orbit_val_k)) DEALLOCATE(dataset%orbit_val_k) 553 ALLOCATE(dataset%orbit_val_n(dataset%norbit_val)) 554 ALLOCATE(dataset%orbit_val_l(dataset%norbit_val)) 555 ALLOCATE(dataset%orbit_val_k(dataset%norbit_val)) 556 kk=0 557 io=0;nval=0 558 DO ll=0,4 559 nn=dataset%np(ll+1) 560 IF (nn>0) THEN 561 DO ik=1,MERGE(nkappa(ll+1),1,dataset%diracrelativistic) 562 kk=MERGE(ll,-(ll+1),ik==1);IF (ll==0) kk=-1 563 IF (.NOT.dataset%diracrelativistic) kk=0 564 DO ii=1+ll,nn 565 io=io+1 566 IF (.NOT.dataset%orbit_iscore(io)) THEN 567 nval=nval+1 568 dataset%orbit_val_n(nval)=ii 569 dataset%orbit_val_l(nval)=ll 570 dataset%orbit_val_k(nval)=kk 571 END IF 572 END DO 573 END DO 574 END IF 575 END DO 576 IF (dataset%norbit_val/=nval) STOP 'input_dataset: bug -- wrong nval!' 577 578!Print read data 579 IF (has_to_print) THEN 580 WRITE(STD_OUT,'(3x,a)') "Core and valence orbitals:" 581 IF (.NOT.dataset%diracrelativistic) WRITE(STD_OUT,'(7x,a)') "n l : type" 582 IF (dataset%diracrelativistic) WRITE(STD_OUT,'(7x,a)') "n l kappa : type" 583 io=0 584 DO ll=0,4 585 nn=dataset%np(ll+1) 586 IF (nn>0) THEN 587 IF (.NOT.dataset%diracrelativistic) THEN 588 DO ii=1+ll,nn 589 io=io+1 590 WRITE(STD_OUT,'(7x,i1,1x,i1,2a)') ii,ll," : ", & 591& MERGE("CORE ","VALENCE",dataset%orbit_iscore(io)) 592 END DO 593 ELSE 594 DO ik=1,nkappa(ll+1) 595 kk=MERGE(ll,-(ll+1),ik==1);IF (ll==0) kk=-1 596 DO ii=1+ll,nn 597 io=io+1 598 WRITE(STD_OUT,'(7x,i1,1x,i1,2x,i2,2x,2a)') ii,ll,kk," : ", & 599& MERGE("CORE ","VALENCE",dataset%orbit_iscore(io)) 600 END DO 601 END DO 602 END IF 603 END IF 604 END DO 605 END IF 606 607 608!------------------------------------------------------------------ 609!End reading of AE data. Start reading of basis data 610 ENDIF 611 IF (read_basis_data_) THEN 612 613 614!------------------------------------------------------------------ 615!=== Maximum L for basis functions 616 617 IF (has_to_ask) THEN 618 WRITE(STD_OUT,*) 'Enter maximum L for basis and projector functions' 619 END IF 620 621 READ(STD_IN,'(a)') inputline 622 IF (has_to_echo) WRITE(ecunit,'(a)') TRIM(inputline) 623 CALL eliminate_comment(inputline) 624 625 READ(inputline,*) dataset%lmax 626 627!Print read data 628 IF (has_to_print) THEN 629 WRITE(STD_OUT,'(3x,a,i0)') "Basis, maximum L : ",dataset%lmax 630 END IF 631 632 633!------------------------------------------------------------------ 634!=== Cut-off radii 635 636 IF (has_to_ask) THEN 637 WRITE(STD_OUT,*) 'Enter rc [and possibly: rc_shape, rc_vloc, rc_core]' 638 END IF 639 640 READ(STD_IN,'(a)') inputline 641 IF (has_to_echo) WRITE(ecunit,'(a)') TRIM(inputline) 642 CALL eliminate_comment(inputline) 643 644 CALL extractword(1,inputline,inputword);inputword=trim(inputword) 645 IF (inputword/="") READ(inputword,*) dataset%rc 646 IF (dataset%rc<=1.d-12) THEN 647 WRITE(STD_OUT,*) 'input_dataset: error -- rc too small ',dataset%rc,'!' 648 STOP 649 END IF 650 651 dataset%rc_shap=dataset%rc 652 dataset%rc_vloc=dataset%rc 653 dataset%rc_core=dataset%rc 654 655 CALL extractword(2,inputline,inputword);inputword=trim(inputword) 656 IF (inputword/="") THEN 657 READ(inputword,*) dataset%rc_shap 658 CALL extractword(3,inputline,inputword);inputword=trim(inputword) 659 IF (inputword/="") THEN 660 READ(inputword,*) dataset%rc_vloc 661 CALL extractword(4,inputline,inputword);inputword=trim(inputword) 662 IF (inputword/="") THEN 663 READ(inputword,*) dataset%rc_core 664 ELSE 665 WRITE(STD_OUT,*) 'input_dataset: error -- rc(core) is missing!' 666 STOP 667 END IF 668 ELSE 669 WRITE(STD_OUT,*) 'input_dataset: error -- rc(Vloc) is missing!' 670 STOP 671 END IF 672 IF (dataset%rc_shap<=1.d-12.OR.dataset%rc_vloc<=1.d-12.OR.& 673& dataset%rc_core<=1.d-12) THEN 674 WRITE(STD_OUT,*) 'input_dataset: error -- one rc is too small!' 675 STOP 676 END IF 677 IF (dataset%rc_shap>dataset%rc.OR.dataset%rc_vloc>dataset%rc.OR.& 678& dataset%rc_core>dataset%rc) THEN 679 WRITE(STD_OUT,*) 'input_dataset: error -- rc_shape, rc_vloc and rc_core must be <rc!' 680 STOP 681 END IF 682 ENDIF 683 684!Print read data 685 IF (has_to_print) THEN 686 WRITE(STD_OUT,'(3x,a,f7.4)') "Augmentation region radius : ",dataset%rc 687 WRITE(STD_OUT,'(3x,a,f7.4)') "Core dens. matching radius : ",dataset%rc_core 688 WRITE(STD_OUT,'(3x,a,f7.4)') "Local pot. matching radius : ",dataset%rc_vloc 689 WRITE(STD_OUT,'(3x,a,f7.4)') "Compens. shape func radius : ",dataset%rc_shap 690 END IF 691 692 693!------------------------------------------------------------------ 694!=== Additional basis functions 695 696 nstart=0 ; dataset%nbasis_add=0 ; basis_add_k(:)=0 697 DO ll=0,dataset%lmax 698 nbl=0 699 nadd = MERGE(nkappa(ll+1),1,dataset%diracrelativistic) 700 IF (dataset%np(ll+1)>0) THEN 701 nbl=COUNT(.NOT.dataset%orbit_iscore(nstart+1:nstart+dataset%np(ll+1)-ll)) 702 nstart=nstart+dataset%np(ll+1)-ll 703 END IF 704 DO 705 IF (has_to_ask) THEN 706 WRITE(STD_OUT,*) 'For l = ',ll,' there are currently ',nbl,' basis functions' 707 IF (.NOT.dataset%diracrelativistic) THEN 708 WRITE(STD_OUT,*) 'Enter y to add an additional function or n to go to next l' 709 ELSE 710 WRITE(STD_OUT,*) 'Enter y to add additional functions (one per kappa) or n to go to next l' 711 END IF 712 END IF 713 READ(STD_IN,'(a)') inputline 714 IF (has_to_echo) WRITE(ecunit,'(a)') TRIM(inputline) 715 CALL eliminate_comment(inputline) 716 READ(inputline,*) CHR 717 IF (CHR/='y'.AND.CHR/='Y') THEN 718 IF (CHR/='n'.AND.CHR/='N') STOP 'input_dataset: error -- Please enter Y or N!' 719 EXIT 720 END IF 721 dataset%nbasis_add=dataset%nbasis_add+nadd 722 IF (dataset%nbasis_add>nbasis_add_max) STOP 'Too many additional basis functions!' 723 basis_add_l(dataset%nbasis_add-nadd+1:dataset%nbasis_add)=ll 724 IF (dataset%diracrelativistic) THEN 725 basis_add_k(dataset%nbasis_add)=-1 726 IF (ll/=0) THEN 727 basis_add_k(dataset%nbasis_add-1)=ll 728 basis_add_k(dataset%nbasis_add)=-(ll+1) 729 END IF 730 END IF 731 IF (has_to_ask) THEN 732 IF (.NOT.dataset%diracrelativistic) THEN 733 WRITE(STD_OUT,*) 'Enter energy for generalized function' 734 ELSE 735 IF (ll==0) WRITE(STD_OUT,*) 'Enter 1 energy for generalized function (kappa=-1)' 736 IF (ll/=0) WRITE(STD_OUT,*) 'Enter 2 energies for generalized functions (one energy per kappa)' 737 END IF 738 END IF 739 READ(STD_IN,'(a)') inputline 740 IF (has_to_echo) WRITE(ecunit,'(a)') TRIM(inputline) 741 CALL eliminate_comment(inputline) 742 READ(inputline,*) basis_add_energy(dataset%nbasis_add-nadd+1:dataset%nbasis_add) 743 END DO 744 END DO 745 746 IF (ALLOCATED(dataset%basis_add_l)) DEALLOCATE(dataset%basis_add_l) 747 IF (ALLOCATED(dataset%basis_add_k)) DEALLOCATE(dataset%basis_add_k) 748 IF (ALLOCATED(dataset%basis_add_energy)) DEALLOCATE(dataset%basis_add_energy) 749 ALLOCATE(dataset%basis_add_l(dataset%nbasis_add)) 750 ALLOCATE(dataset%basis_add_k(dataset%nbasis_add)) 751 ALLOCATE(dataset%basis_add_energy(dataset%nbasis_add)) 752 IF (dataset%nbasis_add>0) THEN 753 dataset%basis_add_l(1:dataset%nbasis_add)=basis_add_l(1:dataset%nbasis_add) 754 dataset%basis_add_k(1:dataset%nbasis_add)=basis_add_k(1:dataset%nbasis_add) 755 dataset%basis_add_energy(1:dataset%nbasis_add)=basis_add_energy(1:dataset%nbasis_add) 756 END IF 757 758 dataset%nbasis=COUNT(.NOT.dataset%orbit_iscore(:))+dataset%nbasis_add 759 760!Print read data 761 IF (has_to_print) THEN 762 WRITE(STD_OUT,'(3x,a,i0)') "Initial number of basis functions : ",dataset%nbasis-dataset%nbasis_add 763 WRITE(STD_OUT,'(3x,a,i0)') "Number of additional basis functions : ",dataset%nbasis_add 764 WRITE(STD_OUT,'(3x,a,i0)') "Total number of basis functions : ",dataset%nbasis 765 WRITE(STD_OUT,'(3x,a)') "Additional basis functions:" 766 IF (.NOT.dataset%diracrelativistic) THEN 767 WRITE(STD_OUT,'(7x,a)') "l : energy" 768 DO io=1,dataset%nbasis_add 769 WRITE(STD_OUT,'(7x,i1,a,f7.4)') dataset%basis_add_l(io)," : " ,dataset%basis_add_energy(io) 770 END DO 771 ELSE 772 WRITE(STD_OUT,'(7x,a)') "l kappa : energy" 773 DO io=1,dataset%nbasis_add 774 WRITE(STD_OUT,'(7x,i1,2x,i2,2x,a,f7.4)') dataset%basis_add_l(io), & 775& dataset%basis_add_k(io)," : " ,dataset%basis_add_energy(io) 776 END DO 777 END IF 778 END IF 779 780 781!------------------------------------------------------------------ 782!=== Projectors, compensation charge shape function, core tolerance 783 784 IF (has_to_ask) THEN 785 WRITE(STD_OUT,*) 'Enter "Bloechl", "Vanderbilt", "modrrkj" or "custom" keywords',& 786& ' for projector generation method.' 787 WRITE(STD_OUT,*) ' In case of "custom" choice, enter additional (optional) keywords:' 788 WRITE(STD_OUT,*) ' - For partial waves pseudization scheme:' 789 WRITE(STD_OUT,*) ' "bloechlps", "polynom", "polynom2 p qcut" or "RRKJ"' 790 WRITE(STD_OUT,*) ' - For orthogonalization scheme:' 791 WRITE(STD_OUT,*) ' "gramschmidtortho" or "vanderbiltortho" or "svdortho"' 792 WRITE(STD_OUT,*) ' - Sinc^2 shape is set by default,' 793 WRITE(STD_OUT,*) ' Gaussian shape can be specified by adding "Gaussian" keyword and tol,' 794 WRITE(STD_OUT,*) ' Bessel shape can be specified by adding "Besselshape" keyword' 795 END IF 796 797 READ(STD_IN,'(a)') inputline 798 IF (has_to_echo) WRITE(ecunit,'(a)') TRIM(inputline) 799 CALL eliminate_comment(inputline) 800 801 CALL Uppercase(inputline) 802 inputline=TRIM(inputline) 803 804 dataset%pseudo_type=PSEUDO_TYPE_BLOECHL 805 dataset%ortho_type=ORTHO_TYPE_GRAMSCHMIDT 806 dataset%pseudo_polynom2_pdeg=polynom2_pdeg_def 807 dataset%pseudo_polynom2_qcut=polynom2_qcut_def 808 dataset%shapefunc_type=SHAPEFUNC_TYPE_SINC 809 dataset%shapefunc_gaussian_param=gausstol_def 810 dataset%hf_coretol=hf_coretol_def 811 812 READ(unit=inputline,fmt=*) inputword 813 814 IF (TRIM(inputword)=='BLOECHL'.OR.TRIM(inputword)=='VNCT') THEN 815 dataset%projector_type=PROJECTOR_TYPE_BLOECHL 816 dataset%pseudo_type=PSEUDO_TYPE_BLOECHL 817 dataset%ortho_type=ORTHO_TYPE_GRAMSCHMIDT 818 ELSE IF (TRIM(inputword)=='VNCK') THEN 819 dataset%projector_type=PROJECTOR_TYPE_BLOECHL 820 dataset%pseudo_type=PSEUDO_TYPE_BLOECHL_K 821 dataset%ortho_type=ORTHO_TYPE_GRAMSCHMIDT 822 ELSE IF (TRIM(inputword)=='VANDERBILT'.OR.TRIM(inputword)=='VNCTV') THEN 823 dataset%projector_type=PROJECTOR_TYPE_VANDERBILT 824 dataset%pseudo_type=PSEUDO_TYPE_POLYNOM 825 dataset%ortho_type=ORTHO_TYPE_VANDERBILT 826 ELSE IF(TRIM(inputword)=='MODRRKJ') THEN 827 dataset%projector_type=PROJECTOR_TYPE_MODRRKJ 828 dataset%pseudo_type=PSEUDO_TYPE_RRKJ 829 dataset%ortho_type=ORTHO_TYPE_VANDERBILT 830 IF (INDEX(inputline,'VANDERBILTORTHO')>0) dataset%ortho_type=ORTHO_TYPE_VANDERBILT 831 IF (INDEX(inputline,'GRAMSCHMIDTORTHO')>0) dataset%ortho_type=ORTHO_TYPE_GRAMSCHMIDT 832 IF (INDEX(inputline,'SVDORTHO')>0) dataset%ortho_type=ORTHO_TYPE_SVD 833 ELSE IF (TRIM(inputword)=='CUSTOM') THEN 834 dataset%projector_type=PROJECTOR_TYPE_CUSTOM 835 IF (INDEX(inputline,'BLOECHLPS')>0) THEN 836 dataset%pseudo_type=PSEUDO_TYPE_BLOECHL 837 dataset%ortho_type=ORTHO_TYPE_GRAMSCHMIDT 838 ELSE IF (INDEX(inputline,'POLYNOM2')>0) THEN 839 dataset%pseudo_type=PSEUDO_TYPE_POLYNOM2 840 nstart=INDEX(inputline,'POLYNOM2') 841 READ(unit=inputline(nstart+8:),fmt=*,err=111,end=111,iostat=io) & 842& dataset%pseudo_polynom2_pdeg,dataset%pseudo_polynom2_qcut 843111 CONTINUE 844 ELSE IF (INDEX(inputline,'POLYNOM')>0) THEN 845 dataset%pseudo_type=PSEUDO_TYPE_POLYNOM 846 ELSE IF (INDEX(inputline,'RRKJ')>0) THEN 847 dataset%pseudo_type=PSEUDO_TYPE_RRKJ 848 END IF 849 IF (INDEX(inputline,'VANDERBILTORTHO')>0) dataset%ortho_type=ORTHO_TYPE_VANDERBILT 850 IF (INDEX(inputline,'GRAMSCHMIDTORTHO')>0) dataset%ortho_type=ORTHO_TYPE_GRAMSCHMIDT 851 END IF 852 853 IF (TRIM(dataset%exctype)=='HF') THEN 854 dataset%projector_type=PROJECTOR_TYPE_HF 855 dataset%pseudo_type=PSEUDO_TYPE_HF 856 dataset%ortho_type=ORTHO_TYPE_HF 857 WRITE(STD_OUT,'(3x,a)') '>> You are using HF XC type: pseudo and orthogonalization line will be ignored!' 858 END IF 859 860 IF ((dataset%pseudo_type==PSEUDO_TYPE_BLOECHL.OR.dataset%pseudo_type==PSEUDO_TYPE_BLOECHL_K)& 861& .AND.dataset%ortho_type==ORTHO_TYPE_VANDERBILT) STOP & 862& 'input_dataset: error -- Vanderbilt orthogonalization not compatible with Bloechls projector scheme!' 863 864 IF (INDEX(inputline,'SINC2')>0) THEN 865 dataset%shapefunc_type=SHAPEFUNC_TYPE_SINC 866 ELSE IF (INDEX(inputline,'GAUSSIAN')>0) THEN 867 dataset%shapefunc_type=SHAPEFUNC_TYPE_GAUSSIAN 868 nstart=INDEX(inputline,'GAUSSIAN') 869 READ(unit=inputline(nstart+8:),fmt=*,err=222,end=222,iostat=io) & 870& dataset%shapefunc_gaussian_param 871222 CONTINUE 872 ELSE IF (INDEX(inputline,'BESSELSHAPE')>0) THEN 873 dataset%shapefunc_type=SHAPEFUNC_TYPE_BESSEL 874 END IF 875 876 nstart=INDEX(inputline,'CORETOL') 877 IF (nstart>0) THEN 878 READ(unit=inputline(nstart+7:),fmt=*) dataset%hf_coretol 879 END IF 880 881!Print read data 882 IF (has_to_print) THEN 883 WRITE(STD_OUT,'(3x,a)') "Projectors description:" 884 IF (dataset%projector_type==PROJECTOR_TYPE_BLOECHL) & 885& WRITE(STD_OUT,'(7x,a)') "Type : BLOECHL" 886 IF (dataset%projector_type==PROJECTOR_TYPE_VANDERBILT) & 887& WRITE(STD_OUT,'(7x,a)') "Type : VANDERBILT" 888 IF (dataset%projector_type==PROJECTOR_TYPE_MODRRKJ) & 889& WRITE(STD_OUT,'(7x,a)') "Type : MODRRKJ" 890 IF (dataset%projector_type==PROJECTOR_TYPE_CUSTOM) & 891& WRITE(STD_OUT,'(7x,a)') "Type : CUSTOM" 892 IF (dataset%projector_type==PROJECTOR_TYPE_HF) & 893& WRITE(STD_OUT,'(7x,a)') "Type : HARTREE-FOCK" 894 IF (dataset%projector_type/=PROJECTOR_TYPE_HF) THEN 895 IF (dataset%pseudo_type==PSEUDO_TYPE_BLOECHL) & 896& WRITE(STD_OUT,'(7x,a)') "Pseudization : BLOECHL" 897 IF (dataset%pseudo_type==PSEUDO_TYPE_POLYNOM) & 898& WRITE(STD_OUT,'(7x,a)') "Pseudization : POLYNOM" 899 IF (dataset%pseudo_type==PSEUDO_TYPE_RRKJ) & 900& WRITE(STD_OUT,'(7x,a)') "Pseudization : RRKJ" 901 IF (dataset%pseudo_type==PSEUDO_TYPE_BLOECHL_K) & 902& WRITE(STD_OUT,'(7x,a)') "Pseudization : BLOECHL KERKER" 903 IF (dataset%pseudo_type==PSEUDO_TYPE_POLYNOM2) & 904& WRITE(STD_OUT,'(7x,a,i0,a,es9.3)') "Pseudization : POLYNOM2, pdeg=", & 905& dataset%pseudo_polynom2_pdeg,", qcut=",dataset%pseudo_polynom2_qcut 906 IF (dataset%ortho_type==ORTHO_TYPE_GRAMSCHMIDT) & 907& WRITE(STD_OUT,'(7x,a)') "Orthogonalisation : GRAM-SCHMIDT" 908 IF (dataset%ortho_type==ORTHO_TYPE_VANDERBILT) & 909& WRITE(STD_OUT,'(7x,a)') "Orthogonalisation : VANDERBILT" 910 IF (dataset%ortho_type==ORTHO_TYPE_SVD) & 911& WRITE(STD_OUT,'(7x,a)') "Orthogonalisation : SVD" 912 END IF 913 IF (dataset%shapefunc_type==SHAPEFUNC_TYPE_GAUSSIAN) & 914& WRITE(STD_OUT,'(3x,a,es9.3)') "Compensation charge shape function : GAUSSIAN, tol=",& 915& dataset%shapefunc_gaussian_param 916 IF (dataset%shapefunc_type==SHAPEFUNC_TYPE_SINC) & 917& WRITE(STD_OUT,'(3x,a)') "Compensation charge shape function : SINC2" 918 IF (dataset%shapefunc_type==SHAPEFUNC_TYPE_BESSEL) & 919& WRITE(STD_OUT,'(3x,a)') "Compensation charge shape function : BESSEL" 920 IF (INDEX(inputline,'CORETOL')>0) & 921& WRITE(STD_OUT,'(3x,a,es9.3)') "Core tolerance for Hartree-Fock : ",dataset%hf_coretol 922 END IF 923 924 925!------------------------------------------------------------------ 926!=== Local pseudopotential 927 928 IF (has_to_ask) THEN 929 WRITE(STD_OUT,*) 'To generate the local pseudopotential, this code can use:' 930 WRITE(STD_OUT,*) ' 1- a Troullier-Martins scheme for specified l value and energy' 931 WRITE(STD_OUT,*) ' 2- a non norm-conserving pseudopotential scheme for specified l value and energy' 932 WRITE(STD_OUT,*) ' 3- a simple pseudization scheme using a single spherical Bessel function' 933 WRITE(STD_OUT,*) ' 4- Vloc == VlocCoef*Shapefunc' 934 WRITE(STD_OUT,*) 'For choice 1, enter (high) l value and energy e' 935 WRITE(STD_OUT,*) 'For choice 2, enter (high) l value, energy e and "ultrasoft"' 936 WRITE(STD_OUT,*) 'For choice 3, enter "bessel"' 937 WRITE(STD_OUT,*) 'For choice 4, enter "setvloc x y" - x is VlocCoef y is VlocRad' 938 END IF 939 940 READ(STD_IN,'(a)') inputline 941 IF (has_to_echo) WRITE(ecunit,'(a)') TRIM(inputline) 942 CALL eliminate_comment(inputline) 943 944 call Uppercase(inputline) 945 inputline=TRIM(inputline) 946 947 dataset%vloc_type=VLOC_TYPE_MTROULLIER 948 dataset%vloc_l=-1 949 dataset%vloc_ene=0.d0 950 dataset%vloc_setvloc_coef=0.d0 951 dataset%vloc_setvloc_rad=dataset%rc 952 dataset%vloc_kerker_power(:)=0 953 954 IF (INDEX(inputline,'MTROULLIER')>0) THEN 955 dataset%vloc_type=VLOC_TYPE_MTROULLIER 956 ELSE IF (INDEX(inputline,'ULTRASOFT')>0) THEN 957 dataset%vloc_type=VLOC_TYPE_ULTRASOFT 958 ELSE IF (INDEX(inputline,'BESSEL')>0) THEN 959 dataset%vloc_type=VLOC_TYPE_BESSEL 960 ELSE IF (INDEX(inputline,'SETVLOC')>0) THEN 961 dataset%vloc_type=VLOC_TYPE_SETVLOC 962 nstart=INDEX(inputline,'SETVLOC') 963 READ(unit=inputline(nstart+8:),fmt=*,err=333,end=333,iostat=io) x1,x2 964 IF (x1<1.d3.AND.x1>-1.d3) dataset%vloc_setvloc_coef=x1 965 IF (x2>1.d-8.AND.x2<dataset%rc) dataset%vloc_setvloc_rad=x2 966333 CONTINUE 967 ELSE IF (INDEX(inputline,'KERKER')>0.OR.dataset%pseudo_type==PSEUDO_TYPE_BLOECHL_K) THEN 968 IF (INDEX(inputline,'EXPF')>0) THEN 969 dataset%vloc_type=VLOC_TYPE_KERKER_EXPF 970 nstart=INDEX(inputline,'EXPF') 971 ELSE IF (INDEX(inputline,'POLY')>0) THEN 972 dataset%vloc_type=VLOC_TYPE_KERKER_POLY 973 nstart=INDEX(inputline,'POLY') 974 ELSE 975 STOP "EXPF or POLY keyword missing!" 976 END IF 977 READ(unit=inputline(nstart+5:),fmt=*,err=334,end=334,iostat=io) & 978& dataset%vloc_kerker_power(1:4) 979334 CONTINUE 980 END IF 981 982 IF (dataset%vloc_type==VLOC_TYPE_MTROULLIER.OR.dataset%vloc_type==VLOC_TYPE_ULTRASOFT) THEN 983 READ(unit=inputline,fmt=*,err=444,end=444,iostat=io) dataset%vloc_l,dataset%vloc_ene 984444 CONTINUE 985 IF (dataset%vloc_l<0.or.dataset%vloc_l>10) STOP 'input_dataset: error while reading Vloc parameters!' 986 END IF 987 988!Print read data 989 IF (has_to_print) THEN 990 IF (dataset%vloc_type==VLOC_TYPE_MTROULLIER) & 991& WRITE(STD_OUT,'(7x,a,i0,a,f7.4)') "Local pseudopotential type : MTROULLIER, l=",& 992& dataset%vloc_l,", energy=",dataset%vloc_ene 993 IF (dataset%vloc_type==VLOC_TYPE_ULTRASOFT) & 994& WRITE(STD_OUT,'(7x,a,i0,a,f7.4)') "Local pseudopotential type : ULTRASOFT, l=",& 995& dataset%vloc_l,", energy=",dataset%vloc_ene 996 IF (dataset%vloc_type==VLOC_TYPE_BESSEL) & 997& WRITE(STD_OUT,'(7x,a)') "Local pseudopotential type : BESSEL" 998 IF (dataset%vloc_type==VLOC_TYPE_SETVLOC) & 999& WRITE(STD_OUT,'(7x,a,es9.4,a,es9.4)') "Local pseudopotential type : SETVLOC, coef=",& 1000& dataset%vloc_setvloc_coef,", rad=",dataset%vloc_setvloc_rad 1001 IF (dataset%vloc_type==VLOC_TYPE_KERKER_EXPF) & 1002& WRITE(STD_OUT,'(7x,a,4(1x,i0))') "Local pseudopotential type : KERKER EXPF, powers=",& 1003& dataset%vloc_kerker_power(1:4) 1004 IF (dataset%vloc_type==VLOC_TYPE_KERKER_POLY) & 1005& WRITE(STD_OUT,'(7x,a,4(1x,i0))') "Local pseudopotential type : KERKER POLY, powers=",& 1006& dataset%vloc_kerker_power(1:4) 1007 END IF 1008 1009 1010!------------------------------------------------------------------ 1011!=== Matching radii for the basis functions 1012 1013!Not for all choice of projectors 1014 IF (dataset%projector_type==PROJECTOR_TYPE_CUSTOM.OR.& 1015& dataset%projector_type==PROJECTOR_TYPE_VANDERBILT.OR.& 1016& dataset%projector_type==PROJECTOR_TYPE_MODRRKJ.OR.& 1017& dataset%projector_type==PROJECTOR_TYPE_HF.AND.dataset%vloc_type==VLOC_TYPE_MTROULLIER) THEN 1018 1019 IF (ALLOCATED(dataset%basis_func_rc)) DEALLOCATE(dataset%basis_func_rc) 1020 ALLOCATE(dataset%basis_func_rc(dataset%nbasis)) 1021 1022 IF (has_to_ask) WRITE(STD_OUT,*) 'For each of the following basis functions enter rc' 1023 1024 norb=0 1025 DO ll=0,dataset%lmax 1026 DO ik=1,MERGE(nkappa(ll+1),1,dataset%diracrelativistic) 1027 kk=MERGE(ll,-(ll+1),ik==1);IF (ll==0) kk=-1 1028 IF (.NOT.dataset%diracrelativistic) kk=0 1029 DO io=1,dataset%norbit_val 1030 IF (dataset%orbit_val_l(io)==ll.AND. & 1031& ((.NOT.dataset%diracrelativistic).OR.dataset%orbit_val_k(io)==kk)) THEN 1032 norb=norb+1 1033 IF (has_to_ask.AND.(.NOT.dataset%diracrelativistic)) & 1034& WRITE(STD_OUT,'(a,i2,a,2i4)') ' rc for basis function ',norb,' - n,l= ',norb,ll 1035 IF (has_to_ask.AND.dataset%diracrelativistic) & 1036& WRITE(STD_OUT,'(a,i2,a,3i4)') ' rc for basis function ',norb,' - n,l,kappa= ',norb,ll,kk 1037 READ(STD_IN,'(a)') inputline 1038 IF (has_to_echo) WRITE(ecunit,'(a)') TRIM(inputline) 1039 CALL eliminate_comment(inputline) 1040 READ(inputline,*) dataset%basis_func_rc(norb) 1041 END IF 1042 END DO 1043 IF (dataset%nbasis_add>0) THEN 1044 DO io=1,dataset%nbasis_add 1045 IF (dataset%basis_add_l(io)==ll.AND. & 1046& ((.NOT.dataset%diracrelativistic).OR.dataset%basis_add_k(io)==kk)) THEN 1047 norb=norb+1 1048 IF (has_to_ask.AND.(.NOT.dataset%diracrelativistic)) & 1049& WRITE(STD_OUT,'(a,i2,a,2i4)') ' rc for basis function ',norb,' - n,l= ',999,ll 1050 IF (has_to_ask.AND.dataset%diracrelativistic) & 1051& WRITE(STD_OUT,'(a,i2,a,3i4)') ' rc for basis function ',norb,' - n,l,kappa= ',999,ll,kk 1052 READ(STD_IN,'(a)') inputline 1053 IF (has_to_echo) WRITE(ecunit,'(a)') TRIM(inputline) 1054 CALL eliminate_comment(inputline) 1055 READ(inputline,*) dataset%basis_func_rc(norb) 1056 END IF 1057 END DO 1058 END IF 1059 END DO 1060 END DO 1061 1062 IF (dataset%nbasis/=norb) STOP 'input_dataset: error -- inconsistency in the number of basis functions!' 1063 1064! Print read data 1065 IF (has_to_print) THEN 1066 WRITE(STD_OUT,'(3x,a)') "Matching radius for basis functions:" 1067 IF (.NOT.dataset%diracrelativistic) WRITE(STD_OUT,'(7x,a)') " # - n l : radius" 1068 IF (dataset%diracrelativistic) WRITE(STD_OUT,'(7x,a)') " # - n l kappa : radius" 1069 norb=0 1070 DO ll=0,dataset%lmax 1071 DO ik=1,MERGE(nkappa(ll+1),1,dataset%diracrelativistic) 1072 kk=MERGE(ll,-(ll+1),ik==1);IF (ll==0) kk=-1 1073 IF (.NOT.dataset%diracrelativistic) kk=0 1074 DO io=1,dataset%norbit_val 1075 IF (dataset%orbit_val_l(io)==ll.AND. & 1076& ((.NOT.dataset%diracrelativistic).OR.dataset%orbit_val_k(io)==kk)) THEN 1077 norb=norb+1 1078 IF (.NOT.dataset%diracrelativistic) & 1079& WRITE(STD_OUT,'(7x,i2,a,i1,1x,i1,a,f7.4)') & 1080& norb," - ",dataset%orbit_val_n(io),ll," : ",dataset%basis_func_rc(norb) 1081 IF (dataset%diracrelativistic) & 1082& WRITE(STD_OUT,'(7x,i2,a,i1,1x,i1,2x,i2,2x,a,f7.4)') & 1083& norb," - ",dataset%orbit_val_n(io),ll,kk," : ",dataset%basis_func_rc(norb) 1084 END IF 1085 END DO 1086 IF (dataset%nbasis_add>0) THEN 1087 DO io=1,dataset%nbasis_add 1088 IF (dataset%basis_add_l(io)==ll.AND. & 1089& ((.NOT.dataset%diracrelativistic).OR.dataset%basis_add_k(io)==kk)) THEN 1090 norb=norb+1 1091 IF (.NOT.dataset%diracrelativistic) & 1092& WRITE(STD_OUT,'(7x,i2,a,a1,1x,i1,a,f7.4)') & 1093& norb," - ",".",ll," : ",dataset%basis_func_rc(norb) 1094 IF (dataset%diracrelativistic) & 1095& WRITE(STD_OUT,'(7x,i2,a,a1,1x,i1,2x,i2,2x,a,f7.4)') & 1096& norb," - ",".",ll,kk," : ",dataset%basis_func_rc(norb) 1097 END IF 1098 END DO 1099 END IF 1100 END DO 1101 END DO 1102 END IF 1103 1104 ELSE ! Other projectors 1105 IF (ALLOCATED(dataset%basis_func_rc)) DEALLOCATE(dataset%basis_func_rc) 1106 ALLOCATE(dataset%basis_func_rc(0)) 1107 END IF 1108 1109 1110!------------------------------------------------------------------ 1111!End reading of basis data 1112 ENDIF 1113 1114!Final message 1115 IF (read_global_data_.OR.read_elec_data_.OR.read_coreval_data_.OR.read_basis_data_) THEN 1116 IF (has_to_ask) THEN 1117 WRITE(STD_OUT,'(1x,a)') "===== END READING OF INPUT DATA =====" 1118 ELSE 1119 WRITE(STD_OUT,'(3x,a)') "===== END READING OF INPUT FILE =====" 1120 END IF 1121 END IF 1122 WRITE(STD_OUT,'(2/)') 1123 1124 1125!------------------------------------------------------------------ 1126!Close files 1127 IF (PRESENT(inputfile)) THEN 1128 CLOSE(ifunit) 1129 END IF 1130 IF (has_to_echo) THEN 1131 CLOSE(ecunit) 1132 END IF 1133 1134 END SUBROUTINE input_dataset_read 1135 1136 1137!!================================================================= 1138!! NAME 1139!! input_dataset_free 1140!! 1141!! FUNCTION 1142!! Free and destroy a input_dataset datastructure 1143!! 1144!! SIDE EFFECT 1145!! [input_dt]= datastructure to be destroyed (optional) 1146!! If omitted, then the global public `input_dataset` 1147!! is used. 1148!! 1149!!================================================================= 1150 SUBROUTINE input_dataset_free(input_dt) 1151 1152!---- Arguments 1153 TYPE(input_dataset_t),INTENT(inout),OPTIONAL,TARGET :: input_dt 1154 1155!---- Local variables 1156 TYPE(input_dataset_t),POINTER :: dataset 1157 1158!------------------------------------------------------------------ 1159 1160!Select datastruture to be read 1161 IF (PRESENT(input_dt)) THEN 1162 dataset => input_dt 1163 ELSE 1164 dataset => input_dataset 1165 ENDIF 1166 1167!Integers 1168 dataset%atomic_charge =-1 1169 dataset%finitenucleusmodel =-1 1170 dataset%gridpoints = 0 1171 dataset%nlogderiv = 0 1172 dataset%Fixed_zero_index = 0 1173 dataset%np(1:5) = 0 1174 dataset%norbit = 0 1175 dataset%norbit_mod = 0 1176 dataset%norbit_val = 0 1177 dataset%lmax =-1 1178 dataset%nbasis = 0 1179 dataset%nbasis_add = 0 1180 dataset%projector_type = UNKNOWN_TYPE 1181 dataset%pseudo_type = UNKNOWN_TYPE 1182 dataset%ortho_type = UNKNOWN_TYPE 1183 dataset%vloc_type = UNKNOWN_TYPE 1184 dataset%vloc_l = 0 1185 dataset%pseudo_polynom2_pdeg = 0 1186 dataset%shapefunc_type = UNKNOWN_TYPE 1187 dataset%abinit_log_meshsz = 0 1188 dataset%xml_spl_meshsz = 0 1189 dataset%xml_lda12_orb_l = -1 1190 dataset%xml_lda12_orb_n = -1 1191 dataset%vloc_kerker_power(1:4) = 0 1192 1193!Logicals 1194 dataset%scalarrelativistic = .false. 1195 dataset%diracrelativistic = .false. 1196 dataset%finitenucleus = .false. 1197 dataset%HFpostprocess = .false. 1198 dataset%localizedcoreexchange = .false. 1199 dataset%BDsolve = .false. 1200 dataset%Fixed_zero = .false. 1201 dataset%abinit_usexcnhat = .false. 1202 dataset%abinit_prtcorewf = .false. 1203 dataset%abinit_uselog = .false. 1204 dataset%abinit_userso = .false. 1205 dataset%xml_usexcnhat = .false. 1206 dataset%xml_prtcorewf = .false. 1207 dataset%xml_usespl = .false. 1208 dataset%xml_userso = .false. 1209 dataset%xml_uselda12 = .false. 1210 1211!Reals 1212 dataset%gridrange = 0.d0 1213 dataset%gridmatch = 0.d0 1214 dataset%minlogderiv = 0.d0 1215 dataset%maxlogderiv = 0.d0 1216 dataset%rc = 0.d0 1217 dataset%rc_shap = 0.d0 1218 dataset%rc_vloc = 0.d0 1219 dataset%rc_core = 0.d0 1220 dataset%pseudo_polynom2_qcut = 0.d0 1221 dataset%shapefunc_gaussian_param = 0.d0 1222 dataset%hf_coretol = 0.d0 1223 dataset%vloc_ene = 0.d0 1224 dataset%vloc_setvloc_coef = 0.d0 1225 dataset%vloc_setvloc_rad = 0.d0 1226 dataset%abinit_log_step = 0.d0 1227 dataset%abinit_rso_ecut = 0.d0 1228 dataset%abinit_rso_gfact = 0.d0 1229 dataset%abinit_rso_werror = 0.d0 1230 dataset%xml_rso_ecut = 0.d0 1231 dataset%xml_rso_gfact = 0.d0 1232 dataset%xml_rso_werror = 0.d0 1233 dataset%xml_lda12_ion = 0.d0 1234 dataset%xml_lda12_rcut = 0.d0 1235 dataset%upf_grid_xmin = 0.d0 1236 dataset%upf_grid_zmesh = 0.d0 1237 dataset%upf_grid_dx = 0.d0 1238 dataset%upf_grid_range = 0.d0 1239 1240!Character strings 1241 dataset%atomic_symbol = '' 1242 dataset%gridkey = '' 1243 dataset%exctype = '' 1244 dataset%abinit_author = '' 1245 dataset%xml_author = '' 1246 dataset%xml_comment = '' 1247 dataset%xml_lda12_logfile = '' 1248 1249!Integer allocatable arrays 1250 IF (ALLOCATED(dataset%orbit_mod_l)) DEALLOCATE(dataset%orbit_mod_l) 1251 IF (ALLOCATED(dataset%orbit_mod_n)) DEALLOCATE(dataset%orbit_mod_n) 1252 IF (ALLOCATED(dataset%orbit_mod_k)) DEALLOCATE(dataset%orbit_mod_k) 1253 IF (ALLOCATED(dataset%orbit_val_l)) DEALLOCATE(dataset%orbit_val_l) 1254 IF (ALLOCATED(dataset%orbit_val_n)) DEALLOCATE(dataset%orbit_val_n) 1255 IF (ALLOCATED(dataset%orbit_val_k)) DEALLOCATE(dataset%orbit_val_k) 1256 IF (ALLOCATED(dataset%basis_add_l)) DEALLOCATE(dataset%basis_add_l) 1257 IF (ALLOCATED(dataset%basis_add_k)) DEALLOCATE(dataset%basis_add_k) 1258 1259!Logical allocatable arrays 1260 IF (ALLOCATED(dataset%orbit_iscore)) DEALLOCATE(dataset%orbit_iscore) 1261 1262!Real allocatable arrays 1263 IF (ALLOCATED(dataset%orbit_mod_occ)) DEALLOCATE(dataset%orbit_mod_occ) 1264 IF (ALLOCATED(dataset%basis_add_energy)) DEALLOCATE(dataset%basis_add_energy) 1265 IF (ALLOCATED(dataset%basis_func_rc)) DEALLOCATE(dataset%basis_func_rc) 1266 1267 END SUBROUTINE input_dataset_free 1268 1269 1270!!================================================================= 1271!! NAME 1272!! input_dataset_copy 1273!! 1274!! FUNCTION 1275!! Copy a input_dataset datastructure into another 1276!! 1277!! INPUTS 1278!! input_dt= datastructure to copy 1279!! 1280!! OUTPUT 1281!! output_dt= copied datastructure 1282!! 1283!!================================================================= 1284 SUBROUTINE input_dataset_copy(input_dt,output_dt) 1285 1286!---- Arguments 1287 TYPE(input_dataset_t),INTENT(in) :: input_dt 1288 TYPE(input_dataset_t),INTENT(inout) :: output_dt 1289 1290!------------------------------------------------------------------ 1291 1292!Integers 1293 output_dt%atomic_charge =input_dt%atomic_charge 1294 output_dt%finitenucleusmodel =input_dt%finitenucleusmodel 1295 output_dt%gridpoints =input_dt%gridpoints 1296 output_dt%nlogderiv =input_dt%nlogderiv 1297 output_dt%Fixed_zero_index =input_dt%Fixed_zero_index 1298 output_dt%np(1:5) =input_dt%np(1:5) 1299 output_dt%norbit =input_dt%norbit 1300 output_dt%norbit_mod =input_dt%norbit_mod 1301 output_dt%norbit_val =input_dt%norbit_val 1302 output_dt%lmax =input_dt%lmax 1303 output_dt%nbasis =input_dt%nbasis 1304 output_dt%nbasis_add =input_dt%nbasis_add 1305 output_dt%projector_type =input_dt%projector_type 1306 output_dt%pseudo_type =input_dt%pseudo_type 1307 output_dt%ortho_type =input_dt%ortho_type 1308 output_dt%pseudo_polynom2_pdeg =input_dt%pseudo_polynom2_pdeg 1309 output_dt%shapefunc_type =input_dt%shapefunc_type 1310 output_dt%vloc_type =input_dt%vloc_type 1311 output_dt%vloc_l =input_dt%vloc_l 1312 output_dt%abinit_log_meshsz =input_dt%abinit_log_meshsz 1313 output_dt%xml_spl_meshsz =input_dt%xml_spl_meshsz 1314 output_dt%xml_lda12_orb_l =input_dt%xml_lda12_orb_l 1315 output_dt%xml_lda12_orb_n =input_dt%xml_lda12_orb_n 1316 output_dt%vloc_kerker_power(1:4) =input_dt%vloc_kerker_power(1:4) 1317 1318!Logicals 1319 output_dt%scalarrelativistic =input_dt%scalarrelativistic 1320 output_dt%diracrelativistic =input_dt%diracrelativistic 1321 output_dt%finitenucleus =input_dt%finitenucleus 1322 output_dt%HFpostprocess =input_dt%HFpostprocess 1323 output_dt%localizedcoreexchange =input_dt%localizedcoreexchange 1324 output_dt%BDsolve =input_dt%BDsolve 1325 output_dt%Fixed_zero =input_dt%Fixed_zero 1326 output_dt%abinit_usexcnhat =input_dt%abinit_usexcnhat 1327 output_dt%abinit_prtcorewf =input_dt%abinit_prtcorewf 1328 output_dt%abinit_uselog =input_dt%abinit_uselog 1329 output_dt%abinit_userso =input_dt%abinit_userso 1330 output_dt%xml_usexcnhat =input_dt%xml_usexcnhat 1331 output_dt%xml_prtcorewf =input_dt%xml_prtcorewf 1332 output_dt%xml_usespl =input_dt%xml_usespl 1333 output_dt%xml_userso =input_dt%xml_userso 1334 output_dt%xml_uselda12 =input_dt%xml_uselda12 1335 1336!Reals 1337 output_dt%gridrange =input_dt%gridrange 1338 output_dt%gridmatch =input_dt%gridmatch 1339 output_dt%minlogderiv =input_dt%minlogderiv 1340 output_dt%maxlogderiv =input_dt%maxlogderiv 1341 output_dt%rc =input_dt%rc 1342 output_dt%rc_shap =input_dt%rc_shap 1343 output_dt%rc_vloc =input_dt%rc_vloc 1344 output_dt%rc_core =input_dt%rc_core 1345 output_dt%pseudo_polynom2_qcut =input_dt%pseudo_polynom2_qcut 1346 output_dt%shapefunc_gaussian_param=input_dt%shapefunc_gaussian_param 1347 output_dt%hf_coretol =input_dt%hf_coretol 1348 output_dt%vloc_ene =input_dt%vloc_ene 1349 output_dt%vloc_setvloc_coef =input_dt%vloc_setvloc_coef 1350 output_dt%vloc_setvloc_rad =input_dt%vloc_setvloc_rad 1351 output_dt%abinit_log_step =input_dt%abinit_log_step 1352 output_dt%abinit_rso_ecut =input_dt%abinit_rso_ecut 1353 output_dt%abinit_rso_gfact =input_dt%abinit_rso_gfact 1354 output_dt%abinit_rso_werror =input_dt%abinit_rso_werror 1355 output_dt%xml_rso_ecut =input_dt%xml_rso_ecut 1356 output_dt%xml_rso_gfact =input_dt%xml_rso_gfact 1357 output_dt%xml_rso_werror =input_dt%xml_rso_werror 1358 output_dt%xml_lda12_ion =input_dt%xml_lda12_ion 1359 output_dt%xml_lda12_rcut =input_dt%xml_lda12_rcut 1360 output_dt%upf_grid_xmin =input_dt%upf_grid_xmin 1361 output_dt%upf_grid_zmesh =input_dt%upf_grid_zmesh 1362 output_dt%upf_grid_dx =input_dt%upf_grid_dx 1363 output_dt%upf_grid_range =input_dt%upf_grid_range 1364 1365!Character strings 1366 output_dt%atomic_symbol =TRIM(input_dt%atomic_symbol) 1367 output_dt%gridkey =TRIM(input_dt%gridkey) 1368 output_dt%exctype =TRIM(input_dt%exctype) 1369 output_dt%abinit_author =TRIM(input_dt%abinit_author) 1370 output_dt%xml_author =TRIM(input_dt%xml_author) 1371 output_dt%xml_comment =TRIM(input_dt%xml_comment) 1372 output_dt%xml_lda12_logfile =TRIM(input_dt%xml_lda12_logfile) 1373 1374!Integer allocatable arrays 1375 IF (ALLOCATED(output_dt%orbit_mod_l)) DEALLOCATE(output_dt%orbit_mod_l) 1376 IF (ALLOCATED(output_dt%orbit_mod_n)) DEALLOCATE(output_dt%orbit_mod_n) 1377 IF (ALLOCATED(output_dt%orbit_mod_k)) DEALLOCATE(output_dt%orbit_mod_k) 1378 IF (ALLOCATED(output_dt%orbit_val_l)) DEALLOCATE(output_dt%orbit_val_l) 1379 IF (ALLOCATED(output_dt%orbit_val_n)) DEALLOCATE(output_dt%orbit_val_n) 1380 IF (ALLOCATED(output_dt%orbit_val_k)) DEALLOCATE(output_dt%orbit_val_k) 1381 IF (ALLOCATED(output_dt%basis_add_l)) DEALLOCATE(output_dt%basis_add_l) 1382 IF (ALLOCATED(output_dt%basis_add_k)) DEALLOCATE(output_dt%basis_add_k) 1383 ALLOCATE(output_dt%orbit_mod_l(input_dt%norbit_mod)) 1384 ALLOCATE(output_dt%orbit_mod_n(input_dt%norbit_mod)) 1385 ALLOCATE(output_dt%orbit_mod_k(input_dt%norbit_mod)) 1386 ALLOCATE(output_dt%orbit_val_l(input_dt%norbit_val)) 1387 ALLOCATE(output_dt%orbit_val_n(input_dt%norbit_val)) 1388 ALLOCATE(output_dt%orbit_val_k(input_dt%norbit_val)) 1389 ALLOCATE(output_dt%basis_add_l(input_dt%nbasis_add)) 1390 ALLOCATE(output_dt%basis_add_k(input_dt%nbasis_add)) 1391 output_dt%orbit_mod_l(1:input_dt%norbit_mod)=input_dt%orbit_mod_l(1:input_dt%norbit_mod) 1392 output_dt%orbit_mod_n(1:input_dt%norbit_mod)=input_dt%orbit_mod_n(1:input_dt%norbit_mod) 1393 output_dt%orbit_mod_k(1:input_dt%norbit_mod)=input_dt%orbit_mod_k(1:input_dt%norbit_mod) 1394 output_dt%orbit_val_l(1:input_dt%norbit_val)=input_dt%orbit_val_l(1:input_dt%norbit_val) 1395 output_dt%orbit_val_n(1:input_dt%norbit_val)=input_dt%orbit_val_n(1:input_dt%norbit_val) 1396 output_dt%orbit_val_k(1:input_dt%norbit_val)=input_dt%orbit_val_k(1:input_dt%norbit_val) 1397 output_dt%basis_add_l(1:input_dt%nbasis_add)=input_dt%basis_add_l(1:input_dt%nbasis_add) 1398 output_dt%basis_add_k(1:input_dt%nbasis_add)=input_dt%basis_add_k(1:input_dt%nbasis_add) 1399 1400!Logical allocatable arrays 1401 IF (ALLOCATED(output_dt%orbit_iscore)) DEALLOCATE(output_dt%orbit_iscore) 1402 ALLOCATE(output_dt%orbit_iscore(input_dt%norbit)) 1403 output_dt%orbit_iscore(1:input_dt%norbit)=input_dt%orbit_iscore(1:input_dt%norbit) 1404 1405!Real allocatable arrays 1406 IF (ALLOCATED(output_dt%orbit_mod_occ)) DEALLOCATE(output_dt%orbit_mod_occ) 1407 IF (ALLOCATED(output_dt%basis_add_energy)) DEALLOCATE(output_dt%basis_add_energy) 1408 IF (ALLOCATED(output_dt%basis_func_rc)) DEALLOCATE(output_dt%basis_func_rc) 1409 ALLOCATE(output_dt%orbit_mod_occ(input_dt%norbit_mod)) 1410 ALLOCATE(output_dt%basis_add_energy(input_dt%nbasis_add)) 1411 ALLOCATE(output_dt%basis_func_rc(input_dt%nbasis)) 1412 output_dt%orbit_mod_occ(1:input_dt%norbit_mod)=input_dt%orbit_mod_occ(1:input_dt%norbit_mod) 1413 output_dt%basis_add_energy(1:input_dt%nbasis_add)=input_dt%basis_add_energy(1:input_dt%nbasis_add) 1414 output_dt%basis_func_rc(1:input_dt%nbasis)=input_dt%basis_func_rc(1:input_dt%nbasis) 1415 1416 END SUBROUTINE input_dataset_copy 1417 1418 1419!!================================================================= 1420!! NAME 1421!! input_dataset_read_occ 1422!! 1423!! FUNCTION 1424!! Initialize only modified occupations in a data-structure by reading them 1425!! from a file. If file is omitted, then read from standard input. 1426!! Note: file is supposed to be opened. 1427!! 1428!! INPUTS 1429!! np(5)= electronic configuration: number of s,p,d,f,g shells 1430!! diracrel= TRUE if we perform a Dirac relativistic calculation 1431!! [inputfile_unit]= logical unit of input file to be read (optional) 1432!! [echofile_unit]= logical unit of a file to echo input file content (optional) 1433!! 1434!! OUTPUT 1435!! norbit_mod= number of orbital with modified occupations 1436!! orbit_mod_l(norbit_mod)= l quantum numbers for the modifed orbitals 1437!! orbit_mod_n(norbit_mod)= n quantum numbers for the modifed orbitals 1438!! orbit_mod_k(norbit_mod)= kappa quantum numbers for the modifed orbitals 1439!! orbit_mod_occ(norbit_mod)= occupations for the modifed orbitals 1440!! [input_dt]= data-structure containing the complete input file (optional). 1441!! If omitted, then the global public `input_dataset` 1442!! is used. 1443!! 1444!!================================================================= 1445 SUBROUTINE input_dataset_read_occ(norbit_mod,orbit_mod_l,orbit_mod_n,orbit_mod_k,& 1446 & orbit_mod_occ,np,diracrel,inputfile_unit,echofile_unit) 1447 1448!---- Arguments 1449 INTEGER,INTENT(IN) :: np(5) 1450 LOGICAL,INTENT(IN) :: diracrel 1451 INTEGER,INTENT(IN),OPTIONAL :: inputfile_unit,echofile_unit 1452 INTEGER,INTENT(OUT) :: norbit_mod 1453 INTEGER,ALLOCATABLE,INTENT(INOUT) :: orbit_mod_l(:),orbit_mod_n(:),orbit_mod_k(:) 1454 REAL(8),ALLOCATABLE,INTENT(INOUT) :: orbit_mod_occ(:) 1455 1456!---- Local variables 1457 INTEGER,PARAMETER :: nkappa(5)=(/1,2,2,2,2/) 1458 INTEGER :: input_unit,echo_unit,ll,nn,io,ii,kk,ik 1459 LOGICAL :: has_to_echo 1460 CHARACTER(200) :: inputline 1461 REAL(8) :: xocc 1462 INTEGER :: tmp_n(norbit_max),tmp_l(norbit_max),tmp_k(norbit_max) 1463 REAL(8) :: tmp_occ(norbit_max),basis_add_energy(nbasis_add_max) 1464 1465!------------------------------------------------------------------ 1466 1467!Select input file logical unit 1468 input_unit=STD_IN;IF (PRESENT(inputfile_unit)) input_unit=inputfile_unit 1469 1470!Do we echo input file content? 1471 has_to_echo=PRESENT(echofile_unit) 1472 echo_unit=-1;IF (has_to_echo) echo_unit=echofile_unit 1473 1474 IF (has_to_ask) THEN 1475 IF (.NOT.diracrel) THEN 1476 WRITE(STD_OUT,*)'Enter np l occ for all occupations for all revisions' 1477 WRITE(STD_OUT,*)'Enter 0 0 0 to end' 1478 ELSE 1479 WRITE(STD_OUT,*)'Enter np l kappa occ for all occupations for all revisions' 1480 WRITE(STD_OUT,*)'Enter 0 0 0 0 to end' 1481 END IF 1482 END IF 1483 1484 norbit_mod=0 ; kk=0 1485 DO 1486 READ(input_unit,'(a)') inputline 1487 IF (has_to_echo) WRITE(echo_unit,'(a)') TRIM(inputline) 1488 if (.not.diracrel) READ(inputline,*) nn,ll,xocc 1489 if ( diracrel) READ(inputline,*) nn,ll,kk,xocc 1490 IF (nn<=0) EXIT 1491 IF (xocc<0.d0.OR.& 1492& ((.NOT.diracrel).AND.(xocc>2.d0*(2*ll+1))).OR.& 1493& (( diracrel).AND.(xocc>2.d0*ABS(kk)))) THEN 1494 WRITE(STD_OUT,*) 'input_dataset: error in occupations -- ip,l,kap,xocc: ',nn,ll,kk,xocc,'!' 1495 STOP 1496 END IF 1497 norbit_mod=norbit_mod+1 1498 if (norbit_mod>norbit_max) stop 'input_dataset_occ: error -- to many occupation lines!' 1499 tmp_l(norbit_mod)=ll 1500 tmp_n(norbit_mod)=nn 1501 tmp_k(norbit_mod)=kk 1502 tmp_occ(norbit_mod)=xocc 1503 END DO 1504 1505 IF (ALLOCATED(orbit_mod_l)) DEALLOCATE(orbit_mod_l) 1506 IF (ALLOCATED(orbit_mod_n)) DEALLOCATE(orbit_mod_n) 1507 IF (ALLOCATED(orbit_mod_n)) DEALLOCATE(orbit_mod_k) 1508 IF (ALLOCATED(orbit_mod_occ)) DEALLOCATE(orbit_mod_occ) 1509 ALLOCATE(orbit_mod_l(norbit_mod)) 1510 ALLOCATE(orbit_mod_n(norbit_mod)) 1511 ALLOCATE(orbit_mod_k(norbit_mod)) 1512 ALLOCATE(orbit_mod_occ(norbit_mod)) 1513 orbit_mod_l(1:norbit_mod)=tmp_l(1:norbit_mod) 1514 orbit_mod_n(1:norbit_mod)=tmp_n(1:norbit_mod) 1515 orbit_mod_k(1:norbit_mod)=tmp_k(1:norbit_mod) 1516 orbit_mod_occ(1:norbit_mod)=tmp_occ(1:norbit_mod) 1517 1518!Print read data 1519 IF (has_to_print) THEN 1520 WRITE(STD_OUT,'(3x,a)') "Occupations of orbitals:" 1521 IF (.NOT.diracrel) THEN 1522 WRITE(STD_OUT,'(7x,a)') "n l : occ" 1523 DO ll=0,4 1524 IF (np(ll+1)>0) THEN 1525 DO ii=1+ll,np(ll+1) 1526 nn=-1 1527 DO io=1,norbit_mod 1528 IF (orbit_mod_l(io)==ll.AND.orbit_mod_n(io)==ii) THEN 1529 nn=io ; EXIT 1530 END IF 1531 END DO 1532 IF (nn<=0) WRITE(STD_OUT,'(7x,i1,1x,i1,a,f4.1)') ii,ll," : ",real(2*(2*ll+1)) 1533 IF (nn >0) WRITE(STD_OUT,'(7x,i1,1x,i1,a,f4.1)') ii,ll," : ",orbit_mod_occ(nn) 1534 END DO 1535 END IF 1536 END DO 1537 ELSE 1538 WRITE(STD_OUT,'(7x,a)') "n l kappa : occ" 1539 DO ll=0,4 1540 IF (np(ll+1)>0) THEN 1541 DO ik=1,nkappa(ll+1) 1542 kk=MERGE(ll,-(ll+1),ik==1);IF (ll==0) kk=-1 1543 DO ii=1+ll,np(ll+1) 1544 nn=-1 1545 DO io=1,norbit_mod 1546 IF (orbit_mod_l(io)==ll.and.orbit_mod_n(io)==ii.AND.orbit_mod_k(io)==kk) THEN 1547 nn=io ; EXIT 1548 END IF 1549 END DO 1550 IF (nn<=0) WRITE(STD_OUT,'(7x,i1,1x,i1,2x,i2,2x,a,f4.1)') ii,ll,kk," : ",real(2*abs(kk)) 1551 IF (nn >0) WRITE(STD_OUT,'(7x,i1,1x,i1,2x,i2,2x,a,f4.1)') ii,ll,kk," : ",orbit_mod_occ(nn) 1552 END DO 1553 END DO 1554 END IF 1555 END DO 1556 END IF 1557 END IF 1558 1559 END SUBROUTINE input_dataset_read_occ 1560 1561 1562!!================================================================= 1563!! NAME 1564!! input_dataset_read_abinit 1565!! 1566!! FUNCTION 1567!! Read the input file (or standard input) in order to get ABINIT options. 1568!! Put them in a input_dataset datastructure. 1569!! If file is omitted, then read from standard input. 1570!! File is supposed to be opened. 1571!! 1572!! INPUTS 1573!! [inputfile_unit]= logical unit of input file to be read (optional) 1574!! [echofile_unit]= logical unit of a file to echo input file content (optional) 1575!! 1576!! OUTPUT 1577!! [input_dt]= data-structure containing the complete input file (optional). 1578!! If omitted, then the global public `input_dataset` 1579!! %abinit_usexcnhat=TRUE if nhat density (compensation) is included in XC potential 1580!! %abinit_prtcorewf= TRUE is printing of core WF is required 1581!! %abinit_author= author to be printed in the PAW dataset file 1582!! %abinit_uselog=TRUE if data are transfered on a log. grid before being written 1583!! %abinit_log_meshsz=mesh size for the logarithmic grid 1584!! %abinit_log_step=logarithmic step for the logarithmic grid 1585!! %abinit_userso=TRUE if REAL Space Optimization is required 1586!! %abinit_rso_ecut=Real Space Optimization parameter: plane wave cutoff = 1/2 Gmax**2 1587!! %abinit_rso_gfact=Real Space Optimization parameter: Gamma/Gmax ratio 1588!! %abinit_rso_werror=Real Space Optimization parameter: max. error W_l allowed 1589!! [abinit_string]= character string containing the ABINIT options line from input file 1590!! 1591!!================================================================= 1592 SUBROUTINE input_dataset_read_abinit(input_dt,inputfile_unit,echofile_unit,abinit_string) 1593 1594!---- Arguments 1595 INTEGER,INTENT(IN),OPTIONAL :: inputfile_unit,echofile_unit 1596 CHARACTER(*),INTENT(OUT),OPTIONAL :: abinit_string 1597 TYPE(input_dataset_t),INTENT(INOUT),OPTIONAL,TARGET :: input_dt 1598 1599!---- Local variables 1600 INTEGER :: input_unit,echo_unit 1601 INTEGER :: i_author,i_usexcnhat,i_prtcorewf,i_logspline,i_rsoptim,iend 1602 LOGICAL :: has_to_echo 1603 CHARACTER(200) :: inputline,inputstring,inputword 1604 TYPE(input_dataset_t),POINTER :: dataset 1605 1606!------------------------------------------------------------------ 1607 1608!Select datastruture to be read 1609 IF (PRESENT(input_dt)) THEN 1610 dataset => input_dt 1611 ELSE 1612 dataset => input_dataset 1613 ENDIF 1614 1615!Select input file logical unit 1616 input_unit=STD_IN;IF (PRESENT(inputfile_unit)) input_unit=inputfile_unit 1617 1618!Do we echo input file content? 1619 has_to_echo=PRESENT(echofile_unit) 1620 echo_unit=-1;IF (has_to_echo) echo_unit=echofile_unit 1621 1622!Reading of keywords 1623 READ(input_unit,'(a)') inputline 1624 IF (has_to_echo) WRITE(echo_unit,'(a)') TRIM(inputline) 1625 IF (PRESENT(abinit_string)) abinit_string=TRIM(inputline) 1626 CALL Uppercase(inputline) 1627 i_usexcnhat=INDEX(inputline,'USEXCNHAT') 1628 i_prtcorewf=INDEX(inputline,'PRTCOREWF') 1629 i_logspline=INDEX(inputline,'LOGSPLINE') 1630 i_rsoptim =INDEX(inputline,'RSOPTIM') 1631 i_author =INDEX(inputline,'AUTHOR') 1632 1633!Option for core WF printing 1634 dataset%abinit_prtcorewf=MERGE(.true.,PRTCOREWF_DEF,i_prtcorewf>0) 1635 1636!Option for use of NHAT in XC 1637 dataset%abinit_usexcnhat=MERGE(.true.,USEXCNHAT_DEF,i_usexcnhat>0) 1638 1639!To be activated later: 1640!dataset%abinit_vbare=merge(.true.,vbare_def,i_vbare>0) 1641 1642!Options related to the use of REAL SPACE OPTIMIZATION 1643 dataset%abinit_userso=MERGE(.true.,USERSO_DEF,i_rsoptim>0) 1644 dataset%abinit_rso_ecut=RSO_ECUT_DEF 1645 dataset%abinit_rso_gfact=RSO_GFACT_DEF 1646 dataset%abinit_rso_werror=RSO_WERROR_DEF 1647 IF (dataset%abinit_userso) THEN 1648 iend=200 1649 IF (i_usexcnhat>i_rsoptim.AND.i_usexcnhat-1<iend) iend=i_usexcnhat-1 1650 IF (i_prtcorewf>i_rsoptim.AND.i_prtcorewf-1<iend) iend=i_prtcorewf-1 1651 IF (i_logspline>i_rsoptim.AND.i_logspline-1<iend) iend=i_logspline-1 1652 IF (i_author >i_rsoptim.AND.i_author -1<iend) iend=i_author -1 1653 inputstring="";IF (iend>i_rsoptim+7) inputstring=TRIM(inputline(i_rsoptim+7:iend)) 1654 IF (inputstring/="") THEN 1655 CALL extractword(1,inputstring,inputword);inputword=TRIM(inputword) 1656 IF (inputword/="") THEN 1657 READ(inputword,*) dataset%abinit_rso_ecut 1658 CALL extractword(2,inputstring,inputword);inputword=TRIM(inputword) 1659 IF (inputword/="") THEN 1660 READ(inputword,*) dataset%abinit_rso_gfact 1661 CALL extractword(3,inputstring,inputword);inputword=TRIM(inputword) 1662 IF (inputword/="") READ(inputword,*) dataset%abinit_rso_werror 1663 END IF 1664 END IF 1665 END IF 1666 END IF 1667 1668!Options for related to the transfer to a reduced logarihmic grid 1669 dataset%abinit_uselog=MERGE(.true.,USELOG_DEF,i_logspline>0) 1670 dataset%abinit_log_meshsz=LOGGRD_SIZE_DEF 1671 dataset%abinit_log_step=LOGGRD_STEP_DEF 1672 IF (dataset%abinit_uselog) THEN 1673 iend=200 1674 IF (i_usexcnhat>i_logspline.AND.i_usexcnhat-1<iend) iend=i_usexcnhat-1 1675 IF (i_prtcorewf>i_logspline.AND.i_prtcorewf-1<iend) iend=i_prtcorewf-1 1676 IF (i_rsoptim >i_logspline.AND.i_rsoptim -1<iend) iend=i_rsoptim -1 1677 IF (i_author >i_logspline.AND.i_author -1<iend) iend=i_author -1 1678 inputstring="";IF (iend>i_logspline+9) inputstring=TRIM(inputline(i_logspline+9:iend)) 1679 IF (inputstring/="") THEN 1680 CALL extractword(1,inputstring,inputword);inputword=TRIM(inputword) 1681 IF (inputword/="") THEN 1682 READ(inputword,*) dataset%abinit_log_meshsz 1683 CALL extractword(2,inputstring,inputword);inputword=TRIM(inputword) 1684 IF (inputword/="") READ(inputword,*) dataset%abinit_log_step 1685 END IF 1686 END IF 1687 END IF 1688 1689!Author to be mentioned in the ABINIT file 1690 dataset%abinit_author=TRIM(AUTHOR_DEF) 1691 IF (i_author>0) then 1692 iend=200 1693 IF (i_usexcnhat>i_author.AND.i_usexcnhat-1<iend) iend=i_usexcnhat-1 1694 IF (i_prtcorewf>i_author.AND.i_prtcorewf-1<iend) iend=i_prtcorewf-1 1695 IF (i_logspline>i_author.AND.i_logspline-1<iend) iend=i_logspline-1 1696 IF (i_rsoptim >i_author.AND.i_rsoptim -1<iend) iend=i_rsoptim -1 1697 inputstring=TRIM(inputline(i_author+6:iend)) 1698 READ(unit=inputstring,fmt=*) dataset%abinit_author 1699 dataset%abinit_author=TRIM(dataset%abinit_author) 1700 END IF 1701 1702!Print read data 1703 IF (has_to_print) THEN 1704 WRITE(STD_OUT,'(/,2x,a)')'Options for ABINIT file:' 1705 WRITE(STD_OUT,'(2x,2a)') 'ABINIT option : output of core wave funcs =',MERGE("YES"," NO",dataset%abinit_prtcorewf) 1706 WRITE(STD_OUT,'(2x,2a)') 'ABINIT option : transfer to a log grid =',MERGE("YES"," NO",dataset%abinit_uselog) 1707 IF (dataset%abinit_uselog) THEN 1708 WRITE(STD_OUT,'(5x,a,i5)') '- mesh size of the log grid =',dataset%abinit_log_meshsz 1709 WRITE(STD_OUT,'(5x,a,g8.2)') '- step of the log grid = ',dataset%abinit_log_step 1710 END IF 1711 WRITE(STD_OUT,'(2x,2a)') 'ABINIT option : Real Space Optimization =',MERGE("YES"," NO",dataset%abinit_userso) 1712 IF (dataset%abinit_userso) THEN 1713 WRITE(STD_OUT,'(5x,a,f7.2)') '- RSO plane wave cutoff =',dataset%abinit_rso_ecut 1714 WRITE(STD_OUT,'(5x,a,1x,f6.2)') '- RSO Gamma/Gmax ratio =',dataset%abinit_rso_gfact 1715 WRITE(STD_OUT,'(5x,a,g11.3)') '- RSO maximum error = ',dataset%abinit_rso_werror 1716 END IF 1717 WRITE(STD_OUT,'(2x,2a)') 'ABINIT option : use of compensation in XC =',MERGE("YES"," NO",dataset%abinit_usexcnhat) 1718 IF (dataset%abinit_usexcnhat) THEN 1719 WRITE(STD_OUT,'(5x,a)') '- Kresse local ionic potential output in XML file' 1720 ELSE 1721 WRITE(STD_OUT,'(5x,a)') '- Blochl local ionic potential output in XML file' 1722 END IF 1723 END IF 1724 1725 END SUBROUTINE input_dataset_read_abinit 1726 1727 1728!!================================================================= 1729!! NAME 1730!! input_dataset_read_xml 1731!! 1732!! FUNCTION 1733!! Read the input file (or standard input) in order to get UPF options. 1734!! Put them in a input_dataset datastructure. 1735!! If file is omitted, then read from standard input. 1736!! File is supposed to be opened. 1737!! 1738!! INPUTS 1739!! [inputfile_unit]= logical unit of input file to be read (optional) 1740!! [echofile_unit]= logical unit of a file to echo input file content (optional) 1741!! 1742!! OUTPUT 1743!! [input_dt]= data-structure containing the complete input file (optional). 1744!! If omitted, then the global public `input_dataset` 1745!! %_usexcnhat=TRUE if nhat density (compensation) is included in XC potential 1746!! %xml_prtcorewf= TRUE is printing of core WF is required 1747!! %xml_author= author to be printed in the PAW dataset file 1748!! %xml_comment= comment line to be printed in the heder of the PAW dataset file 1749!! %xml_uselog=TRUE if data are transfered on a log. grid before being written 1750!! %xml_log_meshsz=mesh size for the logarithmic grid 1751!! %xml_log_step=logarithmic step for the logarithmic grid 1752!! %xml_userso=TRUE if REAL Space Optimization is required 1753!! %xml_rso_ecut=Real Space Optimization parameter: plane wave cutoff = 1/2 Gmax**2 1754!! %xml_rso_gfact=Real Space Optimization parameter: Gamma/Gmax ratio 1755!! %xml_rso_werror=Real Space Optimization parameter: max. error W_l allowed 1756!! %xml_uselda12= TRUE if LDA-1/2 potential calculation is required 1757!! %xml_lda12_orb_l=LDA-1/2 parameter: l quantum number of the orbital to be ionized 1758!! %xml_lda12_orb_n=LDA-1/2 parameter: n quantum number of the orbital to be ionized 1759!! %xml_lda12_ion=LDA-1/2 parameter: amount of charge to be removed from the ionized orbital 1760!! %xml_lda12_rcut=LDA-1/2 parameter: cut-off radius (in bohr) 1761!! %xml_lda12_logfile=LDA-1/2 parameter: name of the log file 1762!! [xml_string]= character string containing the ABINIT options line from input file 1763!! 1764!!================================================================= 1765 SUBROUTINE input_dataset_read_xml(input_dt,inputfile_unit,echofile_unit,xml_string) 1766 1767!---- Arguments 1768 INTEGER,INTENT(IN),OPTIONAL :: inputfile_unit,echofile_unit 1769 CHARACTER(*),INTENT(OUT),OPTIONAL :: xml_string 1770 TYPE(input_dataset_t),INTENT(INOUT),OPTIONAL,TARGET :: input_dt 1771 1772!---- Local variables 1773 INTEGER :: input_unit,echo_unit 1774 INTEGER :: i_author,i_comment,i_usexcnhat,i_prtcorewf,i_logspline,i_rsoptim,i_lda12,iend 1775 LOGICAL :: has_to_echo 1776 CHARACTER(200) :: inputline,inputstring,inputword 1777 TYPE(input_dataset_t),POINTER :: dataset 1778 1779!------------------------------------------------------------------ 1780 1781!Select datastruture to be read 1782 IF (PRESENT(input_dt)) THEN 1783 dataset => input_dt 1784 ELSE 1785 dataset => input_dataset 1786 ENDIF 1787 1788!Select input file logical unit 1789 input_unit=STD_IN;IF (PRESENT(inputfile_unit)) input_unit=inputfile_unit 1790 1791!Do we echo input file content? 1792 has_to_echo=PRESENT(echofile_unit) 1793 echo_unit=-1;IF (has_to_echo) echo_unit=echofile_unit 1794 1795!Reading of keywords 1796 READ(input_unit,'(a)') inputline 1797 IF (has_to_echo) WRITE(echo_unit,'(a)') TRIM(inputline) 1798 IF (PRESENT(xml_string)) xml_string=TRIM(inputline) 1799 CALL Uppercase(inputline) 1800 i_usexcnhat=INDEX(inputline,'USEXCNHAT') 1801 i_prtcorewf=INDEX(inputline,'PRTCOREWF') 1802 i_logspline=INDEX(inputline,'WITHSPLGRID') 1803 i_rsoptim =INDEX(inputline,'RSOPTIM') 1804 i_lda12 =INDEX(inputline,'LDA12') 1805 i_author =INDEX(inputline,'AUTHOR') 1806 i_comment =INDEX(inputline,'COMMENT') 1807 1808!Option for core WF printing 1809 dataset%xml_prtcorewf=MERGE(.true.,PRTCOREWF_DEF,i_prtcorewf>0) 1810 1811!Option for use of NHAT in XC 1812 dataset%xml_usexcnhat=MERGE(.true.,USEXCNHAT_DEF,i_usexcnhat>0) 1813 1814!To be activated later: 1815!dataset%xml_vbare=merge(.true.,vbare_def,i_vbare>0) 1816 1817!Options related to the use of REAL SPACE OPTIMIZATION 1818 dataset%xml_userso=MERGE(.true.,USERSO_DEF,i_rsoptim>0) 1819 dataset%xml_rso_ecut=RSO_ECUT_DEF 1820 dataset%xml_rso_gfact=RSO_GFACT_DEF 1821 dataset%xml_rso_werror=RSO_WERROR_DEF 1822 IF (dataset%xml_userso) THEN 1823 iend=200 1824 IF (i_usexcnhat>i_rsoptim.AND.i_usexcnhat-1<iend) iend=i_usexcnhat-1 1825 IF (i_prtcorewf>i_rsoptim.AND.i_prtcorewf-1<iend) iend=i_prtcorewf-1 1826 IF (i_logspline>i_rsoptim.AND.i_logspline-1<iend) iend=i_logspline-1 1827 IF (i_lda12 >i_rsoptim.AND.i_lda12 -1<iend) iend=i_lda12 -1 1828 IF (i_author >i_rsoptim.AND.i_author -1<iend) iend=i_author -1 1829 IF (i_comment >i_rsoptim.AND.i_comment -1<iend) iend=i_comment -1 1830 inputstring="";IF (iend>i_rsoptim+7) inputstring=TRIM(inputline(i_rsoptim+7:iend)) 1831 IF (inputstring/="") THEN 1832 CALL extractword(1,inputstring,inputword);inputword=TRIM(inputword) 1833 IF (inputword/="") THEN 1834 READ(inputword,*) dataset%xml_rso_ecut 1835 CALL extractword(2,inputstring,inputword);inputword=TRIM(inputword) 1836 IF (inputword/="") THEN 1837 READ(inputword,*) dataset%xml_rso_gfact 1838 CALL extractword(3,inputstring,inputword);inputword=TRIM(inputword) 1839 IF (inputword/="") READ(inputword,*) dataset%xml_rso_werror 1840 END IF 1841 END IF 1842 END IF 1843 END IF 1844 1845!Options related to the use of LDA-1/2 technique 1846 dataset%xml_uselda12=MERGE(.true.,USELDA12_DEF,i_lda12>0) 1847 dataset%xml_lda12_orb_l=LDA12_ORB_L_DEF 1848 dataset%xml_lda12_orb_n=LDA12_ORB_N_DEF 1849 dataset%xml_lda12_ion=LDA12_ION_DEF 1850 dataset%xml_lda12_rcut=LDA12_RCUT_DEF 1851 dataset%xml_lda12_logfile=TRIM(LDA12_LOGFILE) 1852 IF (dataset%xml_uselda12) THEN 1853 iend=200 1854 IF (i_usexcnhat>i_lda12.AND.i_usexcnhat-1<iend) iend=i_usexcnhat-1 1855 IF (i_prtcorewf>i_lda12.AND.i_prtcorewf-1<iend) iend=i_prtcorewf-1 1856 IF (i_logspline>i_lda12.AND.i_logspline-1<iend) iend=i_logspline-1 1857 IF (i_rsoptim >i_lda12.AND.i_rsoptim -1<iend) iend=i_rsoptim -1 1858 IF (i_author >i_lda12.AND.i_author -1<iend) iend=i_author -1 1859 IF (i_comment >i_lda12.AND.i_comment -1<iend) iend=i_comment -1 1860 inputstring="";IF (iend>i_lda12+5) inputstring=TRIM(inputline(i_lda12+5:iend)) 1861 IF (inputstring/="") THEN 1862 CALL extractword(1,inputstring,inputword);inputword=TRIM(inputword) 1863 IF (inputword/="") THEN 1864 READ(inputword,*) dataset%xml_lda12_orb_n 1865 CALL extractword(2,inputstring,inputword);inputword=TRIM(inputword) 1866 IF (inputword/="") THEN 1867 READ(inputword,*) dataset%xml_lda12_orb_l 1868 CALL extractword(3,inputstring,inputword);inputword=TRIM(inputword) 1869 IF (inputword/="") THEN 1870 READ(inputword,*) dataset%xml_lda12_ion 1871 CALL extractword(4,inputstring,inputword);inputword=TRIM(inputword) 1872 IF (inputword/="") READ(inputword,*) dataset%xml_lda12_rcut 1873 END IF 1874 END IF 1875 END IF 1876 END IF 1877 END IF 1878 1879!Options for related to the transfer to a reduced logarihmic grid 1880 dataset%xml_usespl=MERGE(.true.,USELOG_DEF,i_logspline>0) 1881 dataset%xml_spl_meshsz=LOGGRD_SIZE_DEF 1882 IF (dataset%xml_usespl) THEN 1883 iend=200 1884 IF (i_usexcnhat>i_logspline.AND.i_usexcnhat-1<iend) iend=i_usexcnhat-1 1885 IF (i_prtcorewf>i_logspline.AND.i_prtcorewf-1<iend) iend=i_prtcorewf-1 1886 IF (i_rsoptim >i_logspline.AND.i_rsoptim -1<iend) iend=i_rsoptim -1 1887 IF (i_lda12 >i_logspline.AND.i_lda12 -1<iend) iend=i_lda12 -1 1888 IF (i_author >i_logspline.AND.i_author -1<iend) iend=i_author -1 1889 IF (i_comment >i_logspline.AND.i_comment -1<iend) iend=i_comment -1 1890 inputstring="";IF (iend>i_logspline+11) inputstring=TRIM(inputline(i_logspline+11:iend)) 1891 IF (inputstring/="") THEN 1892 CALL extractword(1,inputstring,inputword);inputword=TRIM(inputword) 1893 IF (inputword/="") READ(inputword,*) dataset%xml_spl_meshsz 1894 END IF 1895 END IF 1896 1897!Author to be mentioned in the XML file 1898 dataset%xml_author=TRIM(AUTHOR_DEF) 1899 IF (i_author>0) then 1900 iend=200 1901 IF (i_usexcnhat>i_author.AND.i_usexcnhat-1<iend) iend=i_usexcnhat-1 1902 IF (i_prtcorewf>i_author.AND.i_prtcorewf-1<iend) iend=i_prtcorewf-1 1903 IF (i_logspline>i_author.AND.i_logspline-1<iend) iend=i_logspline-1 1904 IF (i_rsoptim >i_author.AND.i_rsoptim -1<iend) iend=i_rsoptim -1 1905 IF (i_lda12 >i_author.AND.i_lda12 -1<iend) iend=i_lda12 -1 1906 IF (i_comment >i_author.AND.i_comment -1<iend) iend=i_comment -1 1907 inputstring=TRIM(inputline(i_author+6:)) 1908 READ(unit=inputstring,fmt=*) dataset%xml_author 1909 dataset%xml_author=TRIM(dataset%xml_author) 1910 END IF 1911 1912!Comment to be added in the XML file 1913 dataset%xml_comment=TRIM(COMMENT_DEF) 1914 IF (i_comment>0) then 1915 iend=200 1916 IF (i_usexcnhat>i_comment.AND.i_usexcnhat-1<iend) iend=i_usexcnhat-1 1917 IF (i_prtcorewf>i_comment.AND.i_prtcorewf-1<iend) iend=i_prtcorewf-1 1918 IF (i_logspline>i_comment.AND.i_logspline-1<iend) iend=i_logspline-1 1919 IF (i_rsoptim >i_comment.AND.i_rsoptim -1<iend) iend=i_rsoptim -1 1920 IF (i_lda12 >i_comment.AND.i_lda12 -1<iend) iend=i_lda12 -1 1921 IF (i_author >i_comment.AND.i_author -1<iend) iend=i_author -1 1922 inputstring=TRIM(inputline(i_comment+7:)) 1923 READ(unit=inputstring,fmt=*) dataset%xml_comment 1924 dataset%xml_comment=TRIM(dataset%xml_comment) 1925 END IF 1926 1927!Print read data 1928 IF (has_to_print) THEN 1929 WRITE(STD_OUT,'(/,2x,a)')'Options for XML file:' 1930 WRITE(STD_OUT,'(2x,2a)') 'XML option : output of core wave funcs =',MERGE("YES"," NO",dataset%xml_prtcorewf) 1931 WRITE(STD_OUT,'(2x,2a)') 'XML option : spline to a log grid =',MERGE("YES"," NO",dataset%xml_usespl) 1932 IF (dataset%xml_usespl) THEN 1933 WRITE(STD_OUT,'(5x,a,i5)') '- mesh size of the log grid =',dataset%xml_spl_meshsz 1934 END IF 1935 WRITE(STD_OUT,'(2x,2a)') 'XML option : Real Space Optimization =',MERGE("YES"," NO",dataset%xml_userso) 1936 IF (dataset%xml_userso) THEN 1937 WRITE(STD_OUT,'(5x,a,f7.2)') '- RSO plane wave cutoff =',dataset%xml_rso_ecut 1938 WRITE(STD_OUT,'(5x,a,1x,f6.2)') '- RSO Gamma/Gmax ratio =',dataset%xml_rso_gfact 1939 WRITE(STD_OUT,'(5x,a,g11.3)') '- RSO maximum error = ',dataset%xml_rso_werror 1940 END IF 1941 WRITE(STD_OUT,'(2x,2a)') 'XML option : output of LDA-1/2 pot. =',MERGE("YES"," NO",dataset%xml_uselda12) 1942 IF (dataset%xml_uselda12) THEN 1943 WRITE(STD_OUT,'(5x,a,i2)') '- LDA-1/2 ionized orbital n =',dataset%xml_lda12_orb_n 1944 WRITE(STD_OUT,'(5x,a,i2)') '- LDA-1/2 ionized orbital l =',dataset%xml_lda12_orb_l 1945 WRITE(STD_OUT,'(5x,a,f5.2)') '- LDA-1/2 ionization =',dataset%xml_lda12_ion 1946 WRITE(STD_OUT,'(5x,a,f7.4)') '- LDA-1/2 cutoff radius (au)=',dataset%xml_lda12_rcut 1947 WRITE(STD_OUT,'(5x,3a)') '- LDA-1/2: see ''',trim(dataset%xml_lda12_logfile),& 1948& ''' file to check convergence' 1949 END IF 1950 WRITE(STD_OUT,'(2x,2a)') 'XML option : use of compensation in XC =',MERGE("YES"," NO",dataset%xml_usexcnhat) 1951 IF (dataset%xml_usexcnhat) THEN 1952 WRITE(STD_OUT,'(5x,a)') '- Zero potential and Kresse local ionic potential output in XML file' 1953 ELSE 1954 WRITE(STD_OUT,'(5x,a)') '- Zero potential and Blochl local ionic potential output in XML file' 1955 END IF 1956 END IF 1957 1958 END SUBROUTINE input_dataset_read_xml 1959 1960 1961!!================================================================= 1962!! NAME 1963!! input_dataset_read_upf 1964!! 1965!! FUNCTION 1966!! Read the input file (or standard input) in order to get XML options. 1967!! Put them in a input_dataset datastructure. 1968!! If file is omitted, then read from standard input. 1969!! File is supposed to be opened. 1970!! 1971!! INPUTS 1972!! [inputfile_unit]= logical unit of input file to be read (optional) 1973!! [echofile_unit]= logical unit of a file to echo input file content (optional) 1974!! 1975!! OUTPUT 1976!! [input_dt]= data-structure containing the complete input file (optional). 1977!! If omitted, then the global public `input_dataset` 1978!! %upf_grid_xmin= minimum radius given by the grid 1979!! %upf_grid_zmesh= inverse of the radial step of the grid 1980!! %upf_grid_dx= logarithmic step of the grid 1981!! %upf_grid_range= range of the grid 1982!! [upf_string]= character string containing the ABINIT options line from input file 1983!! 1984!!================================================================= 1985 SUBROUTINE input_dataset_read_upf(input_dt,inputfile_unit,echofile_unit,upf_string) 1986 1987!---- Arguments 1988 INTEGER,INTENT(IN),OPTIONAL :: inputfile_unit,echofile_unit 1989 CHARACTER(*),INTENT(OUT),OPTIONAL :: upf_string 1990 TYPE(input_dataset_t),INTENT(INOUT),OPTIONAL,TARGET :: input_dt 1991 1992!---- Local variables 1993 INTEGER :: input_unit,echo_unit 1994 INTEGER :: i_all,i_dx,i_xmin,i_zmesh,i_range 1995 LOGICAL :: has_to_echo 1996 CHARACTER(200) :: inputline 1997 TYPE(input_dataset_t),POINTER :: dataset 1998 1999!------------------------------------------------------------------ 2000 2001!Select datastruture to be read 2002 IF (PRESENT(input_dt)) THEN 2003 dataset => input_dt 2004 ELSE 2005 dataset => input_dataset 2006 ENDIF 2007 2008!Select input file logical unit 2009 input_unit=STD_IN;IF (PRESENT(inputfile_unit)) input_unit=inputfile_unit 2010 2011!Do we echo input file content? 2012 has_to_echo=PRESENT(echofile_unit) 2013 echo_unit=-1;IF (has_to_echo) echo_unit=echofile_unit 2014 2015!Reading of keywords 2016 READ(input_unit,'(a)') inputline 2017 IF (has_to_echo) WRITE(echo_unit,'(a)') TRIM(inputline) 2018 IF (PRESENT(upf_string)) upf_string=TRIM(inputline) 2019 CALL Uppercase(inputline) 2020 2021!Default values 2022 dataset%upf_grid_dx=UPF_DX_DEF 2023 dataset%upf_grid_xmin=UPF_XMIN_DEF 2024 dataset%upf_grid_zmesh=UPF_ZMESH_DEF 2025 dataset%upf_grid_range=UPF_RANGE_DEF 2026 2027 i_all=INDEX(inputline,'UPFDX')+INDEX(inputline,'UPFXMIN')+& 2028& INDEX(inputline,'UPFZMESH')+INDEX(inputline,'UPFRANGE') 2029 IF (i_all>0) then 2030 2031! Logarithmic step 2032 i_dx=INDEX(inputline,'UPFDX') 2033 IF (i_dx>0) READ(inputline(i_dx+5:256),*) dataset%upf_grid_dx 2034 2035! Minimum radius 2036 i_xmin=INDEX(inputline,'UPFXMIN') 2037 IF (i_xmin>0) READ(inputline(i_xmin+7:256),*) dataset%upf_grid_xmin 2038 2039! Inverse of radial step 2040 i_zmesh=INDEX(inputline,'UPFZMESH') 2041 IF (i_zmesh>0) READ(inputline(i_zmesh+8:256),*) dataset%upf_grid_zmesh 2042 2043! Range 2044 i_range=INDEX(inputline,'UPFRANGE') 2045 IF (i_range>0) READ(inputline(i_range+8:256),*) dataset%upf_grid_range 2046 2047 END IF 2048 2049!Print read data 2050 IF (has_to_print) THEN 2051 WRITE(STD_OUT,'(/,2x,a)') 'Options for UPF file:' 2052 WRITE(STD_OUT,'(2x,a,es10.3)') 'UPF grid option : logarithmic step (upfdx) =',dataset%upf_grid_dx 2053 WRITE(STD_OUT,'(2x,a,es10.3)') 'UPF grid option : minimum grid x (upfxmin) =',dataset%upf_grid_xmin 2054 WRITE(STD_OUT,'(2x,a,es10.3)') 'UPF grid option : inv. of radial step (upfzmesh)=',dataset%upf_grid_zmesh 2055 WRITE(STD_OUT,'(2x,a,es10.3)') 'UPF grid option : grid range (upf_gridrange) =',dataset%upf_grid_range 2056 END IF 2057 2058 END SUBROUTINE input_dataset_read_upf 2059 2060!!================================================================= 2061 2062END MODULE input_dataset_mod 2063