TABLE OF CONTENTS
- ABINIT/m_dtfil
- m_dtfil/datafiles_type
- m_dtfil/dtfil_init
- m_dtfil/dtfil_init_img
- m_dtfil/dtfil_init_time
- m_dtfil/fappnd
- m_dtfil/iofn1
- m_dtfil/isfile
- m_dtfil/mkfilename
ABINIT/m_dtfil [ Modules ]
NAME
m_dtfil
FUNCTION
object and procedures dealing with input/output filenames
COPYRIGHT
Copyright (C) 2008-2024 ABINIT group (XG, MT) 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
16 #if defined HAVE_CONFIG_H 17 #include "config.h" 18 #endif 19 20 #include "abi_common.h" 21 22 module m_dtfil 23 24 use defs_basis 25 use m_abicore 26 use m_errors 27 use m_xmpi 28 use m_dtset 29 30 use m_build_info, only : abinit_version 31 use defs_abitypes, only : MPI_type 32 use m_clib, only : clib_rename 33 use m_fstrings, only : int2char4, rmquotes, sjoin, strcat, basename 34 use m_io_tools, only : open_file, file_exists 35 use m_libpaw_tools, only : libpaw_log_flag_set 36 use m_parser, only : parsefile, intagm 37 38 implicit none 39 40 private
m_dtfil/datafiles_type [ Types ]
NAME
datafiles_type
FUNCTION
The datafiles_type structures datatype gather all the variables related to files, such as filename, and file units. For one dataset, it is initialized in 95_drive/dtfil_init1.F90, and will not change at all during the treatment of the dataset.
SOURCE
57 type, public :: datafiles_type 58 59 ! WARNING : if you modify this datatype, please check whether there might be creation/destruction/copy routines, 60 ! declared in another part of ABINIT, that might need to take into account your modification. 61 62 ! These keywords are only used in algorithms using images of the cell 63 integer :: getwfk_from_image 64 ! index of image from which read WFK file (0 if standard WFK) 65 ! -1: the same image as current one 66 ! 0: no image 67 ! >0: index of an image 68 69 integer :: getden_from_image 70 ! index of image from which read DEN file (0 if standard DEN) 71 ! -1: the same image as current one 72 ! 0: no image 73 ! >0: index of an image 74 75 integer :: getkden_from_image 76 ! index of image from which read KDEN file (0 if standard KDEN) 77 ! -1: the same image as current one 78 ! 0: no image 79 ! >0: index of an image 80 81 integer :: getpawden_from_image 82 ! index of image from which read PAWDEN file (0 if standard PAWDEN) 83 ! -1: the same image as current one 84 ! 0: no image 85 ! >0: index of an image 86 87 integer :: ireadddb 88 ! ireadddb non-zero if the ddb file must be read 89 90 integer :: ireadden 91 ! ireadden non-zero if the den file must be read 92 93 integer :: ireadkden 94 ! ireadkden non-zero if the kden file must be read 95 96 integer :: ireadwf 97 ! if(optdriver/=1), that is, no response-function computation, 98 ! ireadwf non-zero if the wffk file must be read 99 ! (if irdwfk non-zero or getwfk non-zero) 100 ! if(optdriver==1), that is, response-function computation, 101 ! ireadwf non-zero if the wff1 file must be read 102 ! (if ird1wf non-zero or get1wf non-zero) 103 104 integer :: unchi0 ! unit number for chi0 files 105 integer :: unddb ! unit number for Derivative DataBase 106 integer :: unddk ! unit number for ddk 1WF file 107 integer :: undkdk ! unit number for 2WF file (dkdk) 108 integer :: undkde ! unit number for 2WF file (dkde) 109 integer :: unkg ! unit number for k+G data 110 integer :: unkgq ! unit number for k+G+q data 111 integer :: unkg1 ! unit number for first-order k+G+q data 112 integer :: unkss ! unit number for KSS file 113 integer :: unqps ! unit number for QPS file 114 integer :: unscr ! unit number for SCR file 115 integer :: unwff1 ! unit number for wavefunctions, number one 116 integer :: unwff2 ! unit number for wavefunctions, number two 117 integer :: unwff3 ! unit number for wavefunctions, number three 118 integer :: unwffgs ! unit number for ground-state wavefunctions 119 integer :: unwffkq ! unit number for k+q ground-state wavefunctions 120 integer :: unwft1 ! unit number for wavefunctions, temporary one 121 integer :: unwft2 ! unit number for wavefunctions, temporary two 122 integer :: unwft3 ! unit number for wavefunctions, temporary three 123 integer :: unwftgs ! unit number for ground-state wavefunctions, temporary 124 integer :: unwftkq ! unit number for k+q ground-state wavefunctions, temporary 125 integer :: unylm ! unit number for Ylm(k) data 126 integer :: unylm1 ! unit number for first-order Ylm(k+q) data 127 integer :: unpaw ! unit number for temporary PAW data (for ex. rhoij_nk) (Paw only) 128 integer :: unpaw1 ! unit number for temporary PAW first-order cprj1=<c1_k,q|p>(1) data 129 integer :: unpawq ! unit number for temporary PAW cprjq=<c+_k+q|p> at k+qdata 130 integer :: unpos ! unit number for restart molecular dynamics 131 132 ! TODO: All this strings should be initialized with ABI_NOFILE 133 ! so that we can easily test for path /= ABI_NOFILE instead of getwfk /= 0 or irdwfk /= 0 134 135 character(len=fnlen) :: filnam_ds(5) 136 ! if no dataset mode, the five names from the standard input: 137 ! ab_in, ab_out, abi, abo, tmp 138 ! if dataset mode, the same 5 filenames, appended with //'_DS'//trim(jdtset) 139 140 character(len=fnlen) :: filddbsin 141 ! if no dataset mode : abi//'DDB' 142 ! if dataset mode, and getddb==0 : abi//'_DS'//trim(jdtset)//'DDB' 143 ! if dataset mode, and getddb/=0 : abo//'_DS'//trim(jgetddb)//'DDB' 144 145 character(len=fnlen) :: fildensin 146 ! if no dataset mode : abi//'DEN' 147 ! if dataset mode, and getden==0 : abi//'_DS'//trim(jdtset)//'DEN' 148 ! if dataset mode, and getden/=0 : abo//'_DS'//trim(jgetden)//'DEN' 149 150 character(len=fnlen) :: fildvdbin 151 ! if no dataset mode : abi//'DVDB' 152 ! if dataset mode, and getdvdb==0 : abi//'_DS'//trim(jdtset)//'DVDB' 153 ! if dataset mode, and getdvdb/=0 : abo//'_DS'//trim(jgetdvdb)//'DVDB' 154 155 character(len=fnlen) :: filpotin 156 ! Filename used to read POT file. 157 ! Initialize via getpot_filepath 158 159 character(len=fnlen) :: filkdensin 160 ! if no dataset mode : abi//'KDEN' 161 ! if dataset mode, and getkden==0 : abi//'_DS'//trim(jdtset)//'KDEN' 162 ! if dataset mode, and getkden/=0 : abo//'_DS'//trim(jgetkden)//'KDEN' 163 164 character(len=fnlen) :: filpawdensin 165 ! if no dataset mode : abi//'PAWDEN' 166 ! if dataset mode, and getpawden==0 : abi//'_DS'//trim(jdtset)//'PAWDEN' 167 ! if dataset mode, and getpawden/=0 : abo//'_DS'//trim(jgetpawden)//'PAWDEN' 168 169 character(len=fnlen) :: filsigephin 170 ! Filename used to read SIGEPH.nc file. 171 ! Initialize via getsigeph_filepath 172 173 character(len=fnlen) :: filgstorein 174 ! Filename used to read GSTOR.ncE file. 175 ! Initialize via getgstore_filepath 176 177 character(len=fnlen) :: filstat 178 ! tmp//'_STATUS' 179 180 character(len=fnlen) :: fnamewffk 181 ! the name of the ground-state wavefunction file to be read (see driver.F90) 182 183 character(len=fnlen) :: fnamewffq 184 ! the name of the k+q ground-state wavefunction file to be read (see driver.F90) 185 ! only useful in the response-function case 186 187 character(len=fnlen) :: fnamewffddk 188 ! the generic name of the ddk response wavefunction file(s) to be read (see driver.F90) 189 ! (the final name is formed by appending the number of the perturbation) 190 ! only useful in the response-function case 191 192 character(len=fnlen) :: fnamewffdelfd 193 ! the generic name of the electric field response wavefunction file(s) to be read (see driver.F90) 194 ! (the final name is formed by appending the number of the perturbation) 195 ! only useful in the response-function case 196 197 character(len=fnlen) :: fnamewffdkdk 198 ! the generic name of the 2nd order dkdk response wavefunction file(s) to be read (see driver.F90) 199 ! (the final name is formed by appending the number of the perturbation) 200 ! only useful in the response-function case 201 202 character(len=fnlen) :: fnamewffdkde 203 ! the generic name of the 2nd order dkde response wavefunction file(s) to be read (see driver.F90) 204 ! (the final name is formed by appending the number of the perturbation) 205 ! only useful in the response-function case 206 207 character(len=fnlen) :: fnamewff1 208 ! the generic name of the first-order wavefunction file(s) to be read (see driver.F90) 209 ! (the final name is formed by appending the number of the perturbation) 210 ! only useful in the response-function case 211 212 character(len=fnlen) :: fildens1in ! to be described by MVeithen 213 character(len=fnlen) :: fname_tdwf 214 character(len=fnlen) :: fname_w90 215 216 character(len=fnlen) :: fnametmp_wf1 217 character(len=fnlen) :: fnametmp_wf2 218 character(len=fnlen) :: fnametmp_1wf1 219 character(len=fnlen) :: fnametmp_1wf2 220 character(len=fnlen) :: fnametmp_wfgs 221 character(len=fnlen) :: fnametmp_wfkq 222 ! Set of filenames formed from trim(dtfil%filnam_ds(5))//APPEN where APPEN is _WF1, _WF2 ... 223 ! See dtfil_init 224 225 character(len=fnlen) :: fnametmp_kg 226 character(len=fnlen) :: fnametmp_kgq 227 character(len=fnlen) :: fnametmp_kg1 228 character(len=fnlen) :: fnametmp_dum 229 character(len=fnlen) :: fnametmp_ylm 230 character(len=fnlen) :: fnametmp_ylm1 231 character(len=fnlen) :: fnametmp_paw 232 character(len=fnlen) :: fnametmp_paw1 233 character(len=fnlen) :: fnametmp_pawq 234 ! Set of filenames formed from trim(dtfil%filnam_ds(5))//APPEN where APPEN is _KG, _DUM, followed 235 ! by the index of the processor. 236 ! See dtfil_init 237 238 character(len=fnlen) :: fnametmp_cg 239 character(len=fnlen) :: fnametmp_cprj 240 character(len=fnlen) :: fnametmp_eig 241 character(len=fnlen) :: fnametmp_1wf1_eig 242 character(len=fnlen) :: fnametmp_fft 243 character(len=fnlen) :: fnametmp_fft_mgga 244 character(len=fnlen) :: fnametmp_kgs 245 character(len=fnlen) :: fnametmp_sustr 246 character(len=fnlen) :: fnametmp_tdexcit 247 character(len=fnlen) :: fnametmp_tdwf 248 249 !@Bethe-Salpeter 250 ! New files introduced for the Bethe-Salpeter part. 251 252 character(len=fnlen) :: fnameabi_bsham_reso 253 ! if no dataset mode : abi//'BSR' 254 ! if dataset mode, and getbsreso==0 : abi//'_DS'//trim(jdtset)//'BSR' 255 ! if dataset mode, and getbsreso/=0 : abo//'_DS'//trim(jget_reso_bsham)//'BSR' 256 257 character(len=fnlen) :: fnameabi_bsham_coup 258 ! if no dataset mode : abi//'BSC' 259 ! if dataset mode, and getbscoup==0 : abi//'_DS'//trim(jdtset)//'BSC' 260 ! if dataset mode, and getbscoup/=0 : abo//'_DS'//trim(jget_coup_bsham)//'BSC' 261 262 character(len=fnlen) :: fnameabi_bseig 263 ! The name of the file containing the eigenstates and eigenvalues of the Bethe-Salpeter Hamiltonian 264 ! if no dataset mode : abi//'BS_EIG' 265 ! if dataset mode, and getbseig==0 : abi//'_DS'//trim(jdtset)//'BS_EIG' 266 ! if dataset mode, and getbseig/=0 : abo//'_DS'//trim(jget_bseig)//'BS_EIG' 267 268 character(len=fnlen) :: fnameabi_haydock 269 ! The prefix used to construct the names of the files containing the coefficients of the 270 ! continued fractions produced by the Haydock iterative algorithm. 271 ! if no dataset mode : abi//'HAYDOCK' 272 ! if dataset mode, and gethaydock==0 : abi//'_DS'//trim(jdtset)//'HAYDOCK' 273 ! if dataset mode, and gethaydock/=0 : abo//'_DS'//trim(jget_bseig)//'HAYDOCK' 274 275 character(len=fnlen) :: fnameabi_wfkfine 276 ! The name of the file containing the wavefunctions on a fine grid 277 ! if no dataset mode : abi//'WFK' 278 ! if dataset mode, and gethaydock==0 : abi//'_DS'//trim(jdtset)//'WFK' 279 ! if dataset mode, and gethaydock/=0 : abo//'_DS'//trim(jget_bseig)//'WFK' 280 281 !END @BEthe-Salpeter 282 283 !The following filenames do not depend on itimimage, iimage and itime loops. 284 !Note the following convention: 285 ! fnameabo_* are filenames used for ouput results (writing) 286 ! fnameabi_* are filenames used for data that should be read by the code. 287 ! fnametmp_* are filenames used for temporary files that should be erased at the end of each dataset. 288 ! 289 !If a file does not have the corresponding "abi" or the corresponding "abo" name, that means that 290 !that particular file is only used for writing or for reading results, respectively. 291 292 character(len=fnlen) :: fnameabi_efmas 293 character(len=fnlen) :: fnameabi_hes 294 character(len=fnlen) :: fnameabi_phfrq 295 character(len=fnlen) :: fnameabi_phvec 296 character(len=fnlen) :: fnameabi_qps 297 character(len=fnlen) :: fnameabi_scr ! SCReening file (symmetrized inverse dielectric matrix) 298 character(len=fnlen) :: fnameabi_sus ! KS independent-particle polarizability file 299 character(len=fnlen) :: fnameabi_chkp_rdm ! Checkpoint for GW@DFA to read 300 character(len=fnlen) :: fnameabo_ddb 301 character(len=fnlen) :: fnameabo_den 302 character(len=fnlen) :: fnameabo_ks_den ! KS DEN file at Sigma level 303 character(len=fnlen) :: fnameabo_chkp_rdm ! Checkpoint for GW@DFA to write 304 character(len=fnlen) :: fnameabo_dos 305 character(len=fnlen) :: fnameabo_dvdb 306 character(len=fnlen) :: fnameabo_eelf 307 character(len=fnlen) :: fnameabo_eig 308 character(len=fnlen) :: fnameabo_eigi2d 309 character(len=fnlen) :: fnameabo_eigr2d 310 character(len=fnlen) :: fnameabo_em1 311 character(len=fnlen) :: fnameabo_em1_lf 312 character(len=fnlen) :: fnameabo_em1_nlf 313 character(len=fnlen) :: fnameabo_fan 314 character(len=fnlen) :: fnameabo_gkk 315 character(len=fnlen) :: fnameabo_gw 316 character(len=fnlen) :: fnameabo_gw_nlf_mdf 317 character(len=fnlen) :: fnameabo_kss 318 character(len=fnlen) :: fnameabo_moldyn 319 character(len=fnlen) :: fnameabo_pot 320 character(len=fnlen) :: fnameabo_qps ! Quasi-Particle band structure file. 321 character(len=fnlen) :: fnameabo_qp_den 322 character(len=fnlen) :: fnameabo_qp_pawden ! Full QP density 323 character(len=fnlen) :: fnameabo_qp_dos 324 character(len=fnlen) :: fnameabo_qp_eig 325 character(len=fnlen) :: fnameabo_rpa 326 character(len=fnlen) :: fnameabo_rpa_nlf_mdf 327 character(len=fnlen) :: fnameabo_scr 328 character(len=fnlen) :: fnameabo_sgm 329 character(len=fnlen) :: fnameabo_sgr 330 character(len=fnlen) :: fnameabo_sig 331 character(len=fnlen) :: fnameabo_spcur 332 character(len=fnlen) :: fnameabo_sus 333 character(len=fnlen) :: fnameabo_td_ener 334 character(len=fnlen) :: fnameabo_vha 335 character(len=fnlen) :: fnameabo_vpsp 336 character(len=fnlen) :: fnameabo_vso 337 character(len=fnlen) :: fnameabo_vxc 338 character(len=fnlen) :: fnameabo_wan 339 character(len=fnlen) :: fnameabo_wfk 340 character(len=fnlen) :: fnameabo_wfq 341 character(len=fnlen) :: fnameabo_w90 342 character(len=fnlen) :: fnameabo_1wf 343 character(len=fnlen) :: fnameabo_gwdiag 344 character(len=fnlen) :: fnameabo_nlcc_derivs 345 character(len=fnlen) :: fnameabo_pspdata 346 347 !The following filenames are initialized only iniside itimimage, iimage and itime loops, 348 !and are appended with the adequate specifier 'app'. 349 350 character(len=fnlen) :: fnameabo_app 351 character(len=fnlen) :: fnameabo_app_atmden_core 352 character(len=fnlen) :: fnameabo_app_atmden_full 353 character(len=fnlen) :: fnameabo_app_atmden_val 354 character(len=fnlen) :: fnameabo_app_n_tilde 355 character(len=fnlen) :: fnameabo_app_n_one 356 character(len=fnlen) :: fnameabo_app_nt_one 357 character(len=fnlen) :: fnameabo_app_bxsf 358 character(len=fnlen) :: fnameabo_app_cif 359 character(len=fnlen) :: fnameabo_app_den 360 character(len=fnlen) :: fnameabo_app_dos 361 character(len=fnlen) :: fnameabo_app_elf 362 character(len=fnlen) :: fnameabo_app_elf_down 363 character(len=fnlen) :: fnameabo_app_elf_up 364 character(len=fnlen) :: fnameabo_app_eig 365 character(len=fnlen) :: fnameabo_app_fatbands 366 character(len=fnlen) :: fnameabo_app_gden1 367 character(len=fnlen) :: fnameabo_app_gden2 368 character(len=fnlen) :: fnameabo_app_gden3 369 character(len=fnlen) :: fnameabo_app_geo 370 character(len=fnlen) :: fnameabo_app_kden 371 character(len=fnlen) :: fnameabo_app_lden 372 character(len=fnlen) :: fnameabo_app_nesting 373 character(len=fnlen) :: fnameabo_app_pawden 374 character(len=fnlen) :: fnameabo_app_pot 375 character(len=fnlen) :: fnameabo_app_opt 376 character(len=fnlen) :: fnameabo_app_opt2 377 character(len=fnlen) :: fnameabo_app_stm 378 character(len=fnlen) :: fnameabo_app_vclmb 379 character(len=fnlen) :: fnameabo_app_vha 380 character(len=fnlen) :: fnameabo_app_vhxc 381 character(len=fnlen) :: fnameabo_app_vhpsp 382 character(len=fnlen) :: fnameabo_app_vpsp 383 character(len=fnlen) :: fnameabo_app_vxc 384 character(len=fnlen) :: fnameabo_app_wfk 385 character(len=fnlen) :: fnameabo_app_1dm 386 character(len=fnlen) :: fnameabo_app_vha_1dm 387 character(len=fnlen) :: fnameabo_app_vclmb_1dm 388 character(len=fnlen) :: fnametmp_app_den 389 character(len=fnlen) :: fnametmp_app_kden 390 391 end type datafiles_type
m_dtfil/dtfil_init [ Functions ]
[ Top ] [ m_dtfil ] [ Functions ]
NAME
dtfil_init
FUNCTION
Initialize most of the dtfil structured variable (what is left should be initialized inside the itimimage, iimage and itime loops).
INPUTS
dtset=<type datasets_type>contain all input variables for the current dataset filnam(5)=character strings giving file names filstat=character strings giving name of status file idtset=number of the dataset jdtset_(0:ndtset)=actual index of the datasets mpi_enreg=information about MPI parallelization ndtset=number of datasets [image_index]= index of image to be used when appending "_IMGxxx" string to file names. To be used only when an algorithm using images of the cell is activated
OUTPUT
dtfil=<type datafiles_type>infos about file names, file unit numbers (part of which were initialized previously)
NOTES
The array filnam is used for the name of input and output files, and roots for generic input, output or temporary files. Pseudopotential file names are set in pspini and pspatm, using another name. The name filstat will be needed beyond gstate to check the appearance of the "exit" flag, to make a hasty exit, as well as in order to output the status of the computation.
SOURCE
441 subroutine dtfil_init(dtfil,dtset,filnam,filstat,idtset,jdtset_,mpi_enreg,ndtset,& 442 & image_index) ! optional argument 443 444 !Arguments ------------------------------------ 445 !scalars 446 integer, intent(in) :: idtset,ndtset 447 integer, optional, intent(in) :: image_index 448 character(len=fnlen),intent(in) :: filstat 449 type(MPI_type),intent(in) :: mpi_enreg 450 type(datafiles_type),intent(inout) :: dtfil !vz_i 451 !arrays 452 integer :: jdtset_(0:ndtset) 453 character(len=fnlen),intent(in) :: filnam(5) 454 type(dataset_type),intent(in) :: dtset 455 456 !Local variables------------------------------- 457 !scalars 458 ! Define input and output unit numbers (do not forget, unit 5 and 6 are standard input and output) 459 ! Also, unit number 21, 22 and 23 are used in dfpt_nstdy, for the 3 dot wavefunctions. 460 ! Unit 50,51,52 and 53 are used in dfpt_looppert (for ipert=natom+2, ipert=natom+10 and ipert=natom+11). 461 ! Others unit numbers will be used in the case of the variational and 2n+1 expressions. 462 ! In defs_basis, one defines : 463 ! std_in=5, ab_in=5, std_out=6, ab_out=7, tmp_unit=9, tmp_unit2=10 464 ! TODO: Remove all these units and use get_unit API 465 integer,parameter :: unchi0=42,unddb=16,unddk=50,undkdk=54,undkde=55,unkg1=19,unkg=17,unkgq=18 466 integer,parameter :: unpaw=26,unpaw1=27,unpawq=28,unpos=30 467 integer,parameter :: unwff1=1,unwff2=2,unwff3=8,unwffgs=3,unwfkq=4,unwft1=11 468 integer,parameter :: unwft2=12,unwft3=15,unwftgs=13,unwftkq=14,unylm=24,unylm1=25 469 integer,parameter :: unkss=40,unscr=41,unqps=43 470 integer :: ii,iimage,ireadden,ireadkden,ireadwf,ixx,jdtset,will_read 471 character(len=10) :: appen,tag 472 character(len=9) :: stringvar 473 character(len=15) :: stringfile 474 character(len=500) :: msg 475 character(len=fnlen) :: filsus,filddbsin,fildens1in,fildensin,filpawdensin,filkdensin,filqps,filscr,fil_efmas 476 character(len=fnlen) :: fnamewff1,fnamewffddk,fnamewffdelfd,fnamewffdkdk,fnamewffdkde,fnamewffk,fnamewffq 477 character(len=fnlen) :: filbseig,filfft,filhaydock,fil_bsreso,fil_bscoup 478 character(len=fnlen) :: filwfkfine 479 character(len=fnlen) :: filnam_ds(5) 480 character(len=fnlen) :: tmpfil(14) 481 integer :: idtmpfil(14) 482 483 !****************************************************************** 484 485 DBG_ENTER("COLL") 486 487 iimage=0;if (present(image_index)) iimage=image_index 488 489 dtfil%unchi0 =unchi0 490 dtfil%unddb =unddb 491 dtfil%unddk =unddk 492 dtfil%undkde =undkde 493 dtfil%undkdk =undkdk 494 dtfil%unkg =unkg 495 dtfil%unkgq =unkgq 496 dtfil%unkg1 =unkg1 497 dtfil%unkss =unkss 498 dtfil%unqps =unqps 499 dtfil%unscr =unscr 500 dtfil%unwff1 =unwff1 501 dtfil%unwff2 =unwff2 502 dtfil%unwff3 =unwff3 503 dtfil%unwffgs=unwffgs 504 dtfil%unwffkq=unwfkq 505 dtfil%unwft1 =unwft1 506 dtfil%unwft2 =unwft2 507 dtfil%unwft3 =unwft3 508 dtfil%unwftgs=unwftgs 509 dtfil%unwftkq=unwftkq 510 dtfil%unylm =unylm 511 dtfil%unylm1 =unylm1 512 dtfil%unpaw =unpaw 513 dtfil%unpaw1 =unpaw1 514 dtfil%unpawq =unpawq 515 dtfil%unpos =unpos 516 filnam_ds(1:5)=filnam(1:5) 517 jdtset=dtset%jdtset 518 519 ! If multi dataset mode, special treatment of filenames 3 and 4 (density and 520 ! wavefunctions input and output, as well as other output files) 521 if(ndtset>0)then 522 call appdig(jdtset,'',appen) 523 filnam_ds(3)=trim(filnam(3))//'_DS'//trim(appen) 524 filnam_ds(4)=trim(filnam(4))//'_DS'//trim(appen) 525 end if 526 527 ! If multi image mode (nimage>1), special treatment of filenames 4 and 5 528 if(iimage>0)then 529 call appdig(iimage,'',appen) 530 filnam_ds(4)=trim(filnam_ds(4))//'_IMG'//trim(appen) 531 filnam_ds(5)=trim(filnam_ds(5))//'_IMG'//trim(appen) 532 end if 533 534 ! According to getwfk and irdwfk, build _WFK file name, referred as fnamewffk 535 if (iimage >0 .and. dtfil%getwfk_from_image /= 0) then 536 if (dtfil%getwfk_from_image==-1) then 537 call appdig(iimage,'',appen) 538 else 539 call appdig(dtfil%getwfk_from_image,'',appen) 540 end if 541 stringfile='_IMG'//trim(appen)//'_WFK' 542 else 543 stringfile='_WFK' 544 end if 545 stringvar='wfk' 546 call mkfilename(filnam,fnamewffk,dtset%getwfk,idtset,dtset%irdwfk,jdtset_,ndtset,stringfile,stringvar,will_read, & 547 getpath=dtset%getwfk_filepath) 548 549 if (dtset%optdriver /= RUNL_RESPFN) ireadwf = will_read 550 if(ndtset/=0 .and. dtset%optdriver==RUNL_RESPFN .and. will_read==0)then 551 write(msg, '(5a,i0,3a,i0,a,i0,3a)' )& 552 'At least one of the input variables irdwfk and getwfk ',ch10,& 553 'must refer to a valid _WFK file, in the response function',ch10,& 554 'case, while for idtset = ',idtset,',',ch10,& 555 'they are irdwfk= ',dtset%irdwfk,', and getwfk= ',dtset%getwfk,'.',ch10,& 556 'Action: correct irdwfk or getwfk in your input file.' 557 ABI_ERROR(msg) 558 end if 559 560 !Treatment of the other get wavefunction variable, if response function case or nonlinear case 561 if (ANY(dtset%optdriver == [RUNL_RESPFN, RUNL_NONLINEAR, RUNL_EPH, RUNL_LONGWAVE])) then 562 563 ! According to getwfq and irdwfq, build _WFQ file name, referred as fnamewffq 564 stringfile='_WFQ' ; stringvar='wfq' 565 call mkfilename(filnam,fnamewffq,dtset%getwfq,idtset,dtset%irdwfq,jdtset_,ndtset,stringfile,stringvar,will_read, & 566 getpath=dtset%getwfq_filepath) 567 ! If fnamewffq is not initialized thanks to getwfq or irdwfq, use fnamewffk 568 if(will_read==0) fnamewffq = fnamewffk 569 570 ! According to get1wf and ird1wf, build _1WF file name, referred as fnamewff1 571 stringfile='_1WF' ; stringvar='1wf' 572 call mkfilename(filnam,fnamewff1,dtset%get1wf,idtset,dtset%ird1wf,jdtset_,ndtset,stringfile,stringvar,will_read) 573 ireadwf=will_read 574 575 ! According to getddk and irdddk, build _1WF file name, referred as fnamewffddk 576 stringfile='_1WF' ; stringvar='ddk' 577 call mkfilename(filnam,fnamewffddk,dtset%getddk,idtset,dtset%irdddk,jdtset_,ndtset,stringfile,stringvar,will_read) 578 579 ! According to getdelfd, build _1WF file name, referred as fnamewffdelfd 580 stringfile='_1WF' ; stringvar='delfd' 581 call mkfilename(filnam,fnamewffdelfd,dtset%getdelfd,idtset,0,jdtset_,ndtset,stringfile,stringvar,will_read) 582 583 ! According to getdkdk, build _1WF file name, referred as fnamewffdkdk 584 stringfile='_1WF' ; stringvar='dkdk' 585 call mkfilename(filnam,fnamewffdkdk,dtset%getdkdk,idtset,0,jdtset_,ndtset,stringfile,stringvar,will_read) 586 587 ! According to getdkde, build _1WF file name, referred as fnamewffdkde 588 stringfile='_1WF' ; stringvar='dkde' 589 call mkfilename(filnam,fnamewffdkde,dtset%getdkde,idtset,0,jdtset_,ndtset,stringfile,stringvar,will_read) 590 end if 591 592 !------------------------------------------------------------------------------------------- 593 ! Build name of files from dtfil%filnam_ds(3) 594 595 ! According to getddb, build _DDB file name, referred as filddbsin 596 stringfile='_DDB'; stringvar='ddb' 597 call mkfilename(filnam,filddbsin,dtset%getddb,idtset,dtset%irdddb,jdtset_,ndtset,stringfile,stringvar,will_read, & 598 getpath=dtset%getddb_filepath) 599 600 ! According to getpot, build _POT file name 601 stringfile='_POT'; stringvar='pot' 602 call mkfilename(filnam, dtfil%filpotin, 0, idtset, 0, jdtset_, ndtset, stringfile, stringvar, will_read, & 603 getpath=dtset%getpot_filepath) 604 605 ! According to getdvdb, build _DVDB file name 606 stringfile='_DVDB'; stringvar='dvdb' 607 call mkfilename(filnam,dtfil%fildvdbin,dtset%getdvdb,idtset,dtset%irddvdb,jdtset_,ndtset,stringfile,stringvar,will_read, & 608 getpath=dtset%getdvdb_filepath) 609 if (will_read == 0) dtfil%fildvdbin = ABI_NOFILE 610 611 ! According to getsigeph_filepath, build _SIGEPH file name 612 stringfile='_SIGEPH.nc'; stringvar='sigeph' 613 call mkfilename(filnam, dtfil%filsigephin, 0, idtset, 0, jdtset_, ndtset, stringfile, stringvar, will_read, & 614 getpath=dtset%getsigeph_filepath) 615 ! If getsigeph_filepath is not used, will read the output as assumed in the transport driver when called after sigeph 616 if (will_read == 0) dtfil%filsigephin = strcat(filnam_ds(4), "_SIGEPH.nc") 617 618 ! According to getgstore_filepath, build _GSTORE file name 619 stringfile='_GSTORE.nc'; stringvar='gstore' 620 call mkfilename(filnam, dtfil%filgstorein, 0, idtset, 0, jdtset_, ndtset, stringfile, stringvar, will_read, & 621 getpath=dtset%getgstore_filepath) 622 if (will_read == 0) dtfil%filgstorein = ABI_NOFILE 623 624 ! According to getden, build _DEN file name, referred as fildensin 625 ! A default is available if getden is 0 626 if (iimage>0.and.dtfil%getden_from_image/=0) then 627 if (dtfil%getden_from_image==-1) then 628 call appdig(iimage,'',appen) 629 else 630 call appdig(dtfil%getden_from_image,'',appen) 631 end if 632 stringfile='_IMG'//trim(appen)//'_DEN' 633 else 634 stringfile='_DEN' 635 end if 636 stringvar='den' 637 call mkfilename(filnam,fildensin,dtset%getden,idtset,dtset%irdden,jdtset_,ndtset,stringfile,stringvar, will_read, & 638 getpath=dtset%getden_filepath) 639 if(will_read==0)fildensin=trim(filnam_ds(3))//'_DEN' 640 ireadden=will_read 641 if ((dtset%optdriver==RUNL_GWLS.or.dtset%optdriver==RUNL_GSTATE) .and.dtset%iscf<0) ireadden=1 642 643 ! According to getkden and usekden, build _KDEN file name, referred as filkdensin 644 ! A default is available if getkden is 0 645 if(dtset%usekden==1)then 646 if (iimage>0.and.dtfil%getkden_from_image/=0) then 647 if (dtfil%getkden_from_image==-1) then 648 call appdig(iimage,'',appen) 649 else 650 call appdig(dtfil%getkden_from_image,'',appen) 651 end if 652 stringfile='_IMG'//trim(appen)//'_KDEN' 653 else 654 stringfile='_KDEN' 655 end if 656 stringvar='kden' 657 call mkfilename(filnam,filkdensin,dtset%getkden,idtset,dtset%irdkden,jdtset_,ndtset,stringfile,stringvar,will_read) 658 if(will_read==0)filkdensin=trim(filnam_ds(3))//'_KDEN' 659 ireadkden=will_read 660 if ((dtset%optdriver==RUNL_GSTATE.or.dtset%optdriver==RUNL_GWLS).and.dtset%iscf<0) ireadkden=1 661 else 662 ireadkden=0 663 end if 664 665 ! According to getpawden, build _PAWDEN file name, referred as filpawdensin 666 ! A default is available if getpawden is 0 667 if (iimage>0.and.dtfil%getpawden_from_image/=0) then 668 if (dtfil%getpawden_from_image==-1) then 669 call appdig(iimage,'',appen) 670 else 671 call appdig(dtfil%getpawden_from_image,'',appen) 672 end if 673 stringfile='_IMG'//trim(appen)//'_PAWDEN' 674 else 675 stringfile='_PAWDEN' 676 end if 677 stringvar='pawden' 678 call mkfilename(filnam,filpawdensin,dtset%getpawden,idtset,dtset%irdden,jdtset_,ndtset,stringfile,stringvar,will_read) 679 if(will_read==0)filpawdensin=trim(filnam_ds(3))//'_PAWDEN' 680 681 ! According to get1den, build _DEN file name, referred as fildens1in 682 ! A default is available if get1den is 0 683 stringfile='_DEN' ; stringvar='1den' 684 call mkfilename(filnam,fildens1in,dtset%get1den,idtset,dtset%ird1den,jdtset_,ndtset,stringfile,stringvar,will_read) 685 if(will_read==0)fildens1in=trim(filnam_ds(3))//'_DEN' 686 687 ! According to getefmas and irdefmas, build _EFMAS file name, referred as fil_efmas 688 ! A default is available if getefmas is 0 689 stringfile='_EFMAS.nc' ; stringvar='efmas' 690 call mkfilename(filnam,fil_efmas,dtset%getefmas,idtset,dtset%irdefmas,jdtset_,ndtset,stringfile,stringvar,will_read) 691 if(will_read==0)fil_efmas=trim(filnam_ds(3))//'_EFMAS.nc' 692 693 ! According to getscr and irdscr, build _SCR file name, referred as filscr 694 ! A default is available if getscr is 0 695 stringfile='_SCR' ; stringvar='scr' 696 call mkfilename(filnam,filscr,dtset%getscr,idtset,dtset%irdscr,jdtset_,ndtset,stringfile,stringvar,will_read, & 697 getpath=dtset%getscr_filepath) 698 if(will_read==0)filscr=trim(filnam_ds(3))//'_SCR' 699 700 ! According to getsuscep and irdsuscep, build _SUS file name, referred as filsus 701 ! A default is available if getsuscep is 0 702 stringfile='_SUS' ; stringvar='sus' 703 call mkfilename(filnam,filsus,dtset%getsuscep,idtset,dtset%irdsuscep,jdtset_,ndtset,stringfile,stringvar,will_read) 704 if(will_read==0)filsus=TRIM(filnam_ds(3))//'_SUS' 705 706 ! According to getqps and irdqps, build _QPS file name, referred as filqps 707 ! A default is available if getqps is 0 708 stringfile='_QPS' ; stringvar='qps' 709 call mkfilename(filnam,filqps,dtset%getqps,idtset,dtset%irdqps,jdtset_,ndtset,stringfile,stringvar,will_read) 710 if(will_read==0)filqps=trim(filnam_ds(3))//'_QPS' 711 712 ! According to getbseig and irdbseig, build _BSEIG file name, referred as filbseig 713 ! A default is available if getbseig is 0 714 stringfile='_BSEIG' ; stringvar='bseig' 715 call mkfilename(filnam,filbseig,dtset%getbseig,idtset,dtset%irdbseig,jdtset_,ndtset,stringfile,stringvar,will_read) 716 if(will_read==0)filbseig=trim(filnam_ds(3))//'_BSEIG' 717 718 ! According to gethaydock and irdhaydock, build _HAYD file name, referred as filhaydock. 719 ! A default is available if gethaydock is 0 720 stringfile='_HAYDR_SAVE' ; stringvar='haydock' 721 call mkfilename(filnam,filhaydock,dtset%gethaydock,idtset,dtset%irdhaydock,jdtset_,ndtset,stringfile,stringvar,will_read) 722 if(will_read==0)filhaydock=trim(filnam_ds(3))//'_HAYDR_SAVE' 723 724 ! According to getbsr and irdbsr, build _BSR file name, referred as fil_bsreso 725 ! A default is available if getbsr is 0 726 stringfile='_BSR' ; stringvar='bsreso' 727 call mkfilename(filnam,fil_bsreso,dtset%getbsreso,idtset,dtset%irdbsreso,jdtset_,ndtset,stringfile,stringvar,will_read) 728 if(will_read==0) fil_bsreso=trim(filnam_ds(3))//'_BSR' 729 730 ! According to getbsc and irdbsc, build _BSC file name, referred as fil_bscoup 731 ! A default is available if getbsc is 0 732 stringfile='_BSC' ; stringvar='bscoup' 733 call mkfilename(filnam,fil_bscoup,dtset%getbscoup,idtset,dtset%irdbscoup,jdtset_,ndtset,stringfile,stringvar,will_read) 734 if(will_read==0)fil_bscoup=trim(filnam_ds(3))//'_BSC' 735 736 ! According to getwfkfine and irdwfkfine, build _WFK file name, referred as filwfkfine 737 ! A default is avaible if getwfkfine is 0 738 stringfile='_WFK' ; stringvar='wfkfine' 739 call mkfilename(filnam,filwfkfine,dtset%getwfkfine,idtset,dtset%irdwfkfine,jdtset_,ndtset,stringfile,stringvar,will_read, & 740 getpath=dtset%getwfkfine_filepath) 741 if(will_read==0)filwfkfine=trim(filnam_ds(3))//'_WFK' 742 743 dtfil%ireadden =ireadden 744 dtfil%ireadkden =ireadkden 745 dtfil%ireadwf =ireadwf 746 dtfil%filnam_ds(1:5)=filnam_ds(1:5) 747 748 dtfil%fnameabi_bsham_reso=fil_bsreso 749 dtfil%fnameabi_bsham_coup=fil_bscoup 750 dtfil%fnameabi_bseig=filbseig 751 dtfil%fnameabi_haydock=filhaydock 752 dtfil%fnameabi_sus =filsus 753 dtfil%fnameabi_qps =filqps 754 dtfil%fnameabi_scr =filscr 755 dtfil%fnameabi_efmas=fil_efmas 756 dtfil%filddbsin =filddbsin 757 dtfil%fildensin =fildensin 758 dtfil%fildens1in =fildens1in 759 dtfil%filkdensin =filkdensin 760 dtfil%filpawdensin =filpawdensin 761 dtfil%fnameabi_wfkfine = filwfkfine 762 dtfil%filstat =filstat 763 dtfil%fnamewffk =fnamewffk 764 dtfil%fnamewffq =fnamewffq 765 dtfil%fnamewffddk =fnamewffddk 766 dtfil%fnamewffdelfd =fnamewffdelfd 767 dtfil%fnamewffdkdk =fnamewffdkdk 768 dtfil%fnamewffdkde =fnamewffdkde 769 dtfil%fnamewff1 =fnamewff1 770 771 dtfil%fnameabi_hes=trim(dtfil%filnam_ds(3))//'_HES' 772 dtfil%fnameabi_phfrq=trim(dtfil%filnam_ds(3))//'_PHFRQ' 773 dtfil%fnameabi_phvec=trim(dtfil%filnam_ds(3))//'_PHVEC' 774 dtfil%fnameabi_chkp_rdm =trim(dtfil%filnam_ds(3))//'_CHKP_RDM_' 775 776 !------------------------------------------------------------------------------------------- 777 ! Build name of files from dtfil%filnam_ds(4) 778 779 dtfil%fnameabo_ddb=trim(dtfil%filnam_ds(4))//'_DDB' 780 dtfil%fnameabo_den=trim(dtfil%filnam_ds(4))//'_DEN' 781 dtfil%fnameabo_dos=trim(dtfil%filnam_ds(4))//'_DOS' 782 dtfil%fnameabo_dvdb=trim(dtfil%filnam_ds(4))//'_DVDB' 783 dtfil%fnameabo_eelf=trim(dtfil%filnam_ds(4))//'_EELF' 784 dtfil%fnameabo_eig=trim(dtfil%filnam_ds(4))//'_EIG' 785 dtfil%fnameabo_eigi2d=trim(dtfil%filnam_ds(4))//'_EIGI2D' 786 dtfil%fnameabo_eigr2d=trim(dtfil%filnam_ds(4))//'_EIGR2D' 787 dtfil%fnameabo_em1=trim(dtfil%filnam_ds(4))//'_EM1' 788 dtfil%fnameabo_em1_lf=trim(dtfil%filnam_ds(4))//'_EM1_LF' 789 dtfil%fnameabo_em1_nlf=trim(dtfil%filnam_ds(4))//'_EM1_NLF' 790 dtfil%fnameabo_fan=trim(dtfil%filnam_ds(4))//'_FAN' 791 dtfil%fnameabo_gkk=trim(dtfil%filnam_ds(4))//'_GKK' 792 dtfil%fnameabo_gw=trim(dtfil%filnam_ds(4))//'_GW' ! TODO change name 793 dtfil%fnameabo_gwdiag=trim(dtfil%filnam_ds(4))//'_GWDIAG' 794 dtfil%fnameabo_kss=trim(dtfil%filnam_ds(4))//'_KSS' 795 dtfil%fnameabo_moldyn=trim(dtfil%filnam_ds(4))//'_MOLDYN' 796 dtfil%fnameabo_pot=trim(dtfil%filnam_ds(4))//'_POT' 797 dtfil%fnameabo_qps=trim(dtfil%filnam_ds(4))//'_QPS' 798 dtfil%fnameabo_qp_den=trim(dtfil%filnam_ds(4))//'_QP_DEN' 799 dtfil%fnameabo_qp_pawden=trim(dtfil%filnam_ds(4))//'_QP_PAWDEN' 800 dtfil%fnameabo_qp_dos=trim(dtfil%filnam_ds(4))//'_QP_DOS' 801 dtfil%fnameabo_qp_eig=trim(dtfil%filnam_ds(4))//'_QP_DB.nc' ! TODO change name 802 dtfil%fnameabo_rpa=trim(dtfil%filnam_ds(4))//'_RPA' 803 dtfil%fnameabo_scr=trim(dtfil%filnam_ds(4))//'_SCR' 804 dtfil%fnameabo_ks_den=trim(dtfil%filnam_ds(4))//'_KS_DEN' 805 dtfil%fnameabo_chkp_rdm=trim(dtfil%filnam_ds(4))//'_CHKP_RDM_' 806 dtfil%fnameabo_sgm=trim(dtfil%filnam_ds(4))//'_SGM' 807 dtfil%fnameabo_sgr=trim(dtfil%filnam_ds(4))//'_SGR' 808 dtfil%fnameabo_sig=trim(dtfil%filnam_ds(4))//'_SIG' 809 dtfil%fnameabo_spcur=trim(dtfil%filnam_ds(4))//'_SPCUR' 810 dtfil%fnameabo_sus=trim(dtfil%filnam_ds(4))//'_SUS' 811 dtfil%fnameabo_td_ener=trim(dtfil%filnam_ds(4))//'_TDENER' 812 dtfil%fnameabo_vha=trim(dtfil%filnam_ds(4))//'_VHA' 813 dtfil%fnameabo_vpsp=trim(dtfil%filnam_ds(4))//'_VPSP' 814 dtfil%fnameabo_vso=trim(dtfil%filnam_ds(4))//'_VSO' 815 dtfil%fnameabo_vxc=trim(dtfil%filnam_ds(4))//'_VXC' 816 dtfil%fnameabo_wan=trim(dtfil%filnam_ds(4))//'_WAN' 817 dtfil%fnameabo_wfk=trim(dtfil%filnam_ds(4))//'_WFK' 818 dtfil%fnameabo_wfq=trim(dtfil%filnam_ds(4))//'_WFQ' 819 dtfil%fnameabo_w90=trim(dtfil%filnam_ds(4))//'_w90' 820 dtfil%fnameabo_1wf=trim(dtfil%filnam_ds(4))//'_1WF' 821 dtfil%fnameabo_nlcc_derivs=trim(dtfil%filnam_ds(4))//'_nlcc_derivs_' 822 dtfil%fnameabo_pspdata=trim(dtfil%filnam_ds(4))//'_pspdata_' 823 824 !------------------------------------------------------------------------------------------- 825 ! Build name of files from dtfil%filnam_ds(5) 826 dtfil%fnametmp_eig=trim(dtfil%filnam_ds(5))//'_EIG' 827 dtfil%fnametmp_1wf1_eig=trim(dtfil%filnam_ds(5))//'_1WF1_EIG' ! This appendix should be changed ! 828 dtfil%fnametmp_kgs=trim(dtfil%filnam_ds(5))//'_KGS' 829 dtfil%fnametmp_sustr=trim(dtfil%filnam_ds(5))//'_SUSTR' 830 dtfil%fnametmp_tdexcit=trim(dtfil%filnam_ds(5))//'_TDEXCIT' 831 dtfil%fnametmp_tdwf=trim(dtfil%filnam_ds(5))//'_TDWF' 832 dtfil%fnametmp_cg=trim(dtfil%filnam_ds(5))//'_cg' 833 dtfil%fnametmp_cprj=trim(dtfil%filnam_ds(5))//'_cprj' 834 835 !'_WF1' -> dtfil%unwft1 836 !'_WF2' -> dtfil%unwft2 837 !'_KG' -> dtfil%unkg 838 !'_DUM' -> tmp_unit (real dummy name) 839 !'_YLM' -> dtfil%unylm 840 !'_PAW' -> dtfil%unpaw 841 842 tmpfil(1)=trim(dtfil%filnam_ds(5))//'_WF1' ! tmpfil(1) 843 tmpfil(2)=trim(dtfil%filnam_ds(5))//'_WF2' ! tmpfil(2) 844 845 tmpfil(3)=trim(dtfil%filnam_ds(5))//'_KG' ! tmpfil(3) 846 tmpfil(4)=trim(dtfil%filnam_ds(5))//'_DUM' ! tmpfil(4) 847 tmpfil(5)=' ' ! to avoid Valgrind complain. 848 tmpfil(6)=trim(dtfil%filnam_ds(5))//'_YLM' ! tmpfil(6) 849 tmpfil(7)=trim(dtfil%filnam_ds(5))//'_PAW' ! tmpfil(7) 850 851 if(xmpi_paral==1)then ! parallel case: the index of the processor must be appended 852 call int2char4(mpi_enreg%me,tag) 853 ABI_CHECK((tag(1:1)/='#'),'Bug: string length too short!') 854 ixx=1 855 if (xmpi_mpiio == 1 .and. dtset%iomode == IO_MODE_MPI ) ixx=3 856 do ii=ixx,7 857 tmpfil(ii)=trim(tmpfil(ii))//'_P-'//trim(tag) 858 end do 859 end if 860 861 dtfil%fnametmp_wf1=trim(tmpfil(1)) 862 dtfil%fnametmp_wf2=trim(tmpfil(2)) 863 864 dtfil%fnametmp_kg =trim(tmpfil(3)) 865 dtfil%fnametmp_dum=trim(tmpfil(4)) 866 dtfil%fnametmp_ylm=trim(tmpfil(6)) 867 dtfil%fnametmp_paw=trim(tmpfil(7)) 868 869 !Create names for the temporary files based on dtfil%filnam_ds(5) 870 !by appending adequate string. 871 !'_1WF1' -> dtfil%unwft1 872 !'_1WF2' -> dtfil%unwft2 873 !'_KG' -> dtfil%unkg 874 !'_KGQ' -> dtfil%unkgq (not used for the time being) 875 !'_KG1' -> dtfil%unkg1 876 !'_DUM' -> tmp_unit (real dummy name) 877 !'_WFGS' -> dtfil%unwftgs 878 !'_WFKQ' -> dtfil%unwftkq 879 !'_YLM' -> dtfil%unylm 880 !'_YLM1' -> dtfil%unylm1 881 !'_PAW' -> dtfil%unpaw 882 !'_PAW1' -> dtfil%unpaw1 883 !'_PAWQ' -> dtfil%unpawq 884 tmpfil(1) =trim(dtfil%filnam_ds(5))//'_1WF1' 885 tmpfil(2) =trim(dtfil%filnam_ds(5))//'_1WF2' 886 tmpfil(3) =trim(dtfil%filnam_ds(5))//'_KG' 887 tmpfil(4) =trim(dtfil%filnam_ds(5))//'_KGQ' 888 tmpfil(5) =trim(dtfil%filnam_ds(5))//'_KG1' 889 tmpfil(6) =trim(dtfil%filnam_ds(5))//'_DUM' 890 tmpfil(7) =trim(dtfil%filnam_ds(5))//'_WFGS' 891 tmpfil(8) =trim(dtfil%filnam_ds(5))//'_WFKQ' 892 tmpfil(9) =' ' ! for Valgrind, to avoid uninitialized 893 tmpfil(10)=trim(dtfil%filnam_ds(5))//'_YLM' 894 tmpfil(11)=trim(dtfil%filnam_ds(5))//'_YLM1' 895 tmpfil(12)=trim(dtfil%filnam_ds(5))//'_PAW' 896 tmpfil(13)=trim(dtfil%filnam_ds(5))//'_PAW1' 897 tmpfil(14)=trim(dtfil%filnam_ds(5))//'_PAWQ' 898 899 if(xmpi_paral==1) then 900 idtmpfil(:)=0 901 do ii=1,14 902 idtmpfil(ii)=ii 903 end do 904 if (xmpi_mpiio==1.and.dtset%iomode==IO_MODE_MPI)then 905 idtmpfil(1)=0 !_1wf1 906 idtmpfil(2)=0 ! s1wf2 907 idtmpfil(7)=0 ! WFGS 908 idtmpfil(8)=0 ! WFKQ 909 end if 910 call int2char4(mpi_enreg%me,tag) 911 ABI_CHECK((tag(1:1)/='#'),'Bug: string length too short!') 912 do ii=1,14 913 if(idtmpfil(ii) /= 0) tmpfil(ii)=trim(tmpfil(ii))//'_P-'//trim(tag) 914 end do 915 end if 916 917 dtfil%fnametmp_1wf1=trim(tmpfil(1)) 918 dtfil%fnametmp_1wf2=trim(tmpfil(2)) 919 dtfil%fnametmp_kg =trim(tmpfil(3)) 920 dtfil%fnametmp_kgq =trim(tmpfil(4)) 921 dtfil%fnametmp_kg1 =trim(tmpfil(5)) 922 dtfil%fnametmp_dum =trim(tmpfil(6)) 923 dtfil%fnametmp_wfgs=trim(tmpfil(7)) 924 dtfil%fnametmp_wfkq=trim(tmpfil(8)) 925 dtfil%fnametmp_ylm =trim(tmpfil(10)) 926 dtfil%fnametmp_ylm1=trim(tmpfil(11)) 927 dtfil%fnametmp_paw =trim(tmpfil(12)) 928 dtfil%fnametmp_paw1=trim(tmpfil(13)) 929 dtfil%fnametmp_pawq=trim(tmpfil(14)) 930 931 ! Prepare the name of the _FFT file 932 filfft=trim(dtfil%filnam_ds(5))//'_FFT' 933 if(xmpi_paral==1 .or. mpi_enreg%paral_kgb==1)then 934 call int2char4(mpi_enreg%me,tag) 935 ABI_CHECK((tag(1:1)/='#'),'Bug: string length too short!') 936 filfft=trim(filfft)//'_P-'//trim(tag) 937 end if 938 dtfil%fnametmp_fft=filfft 939 dtfil%fnametmp_fft_mgga=trim(filfft)//'_MGGA' 940 941 ! These keywords are only used in algorithms using images of the cell 942 if (iimage==0) then 943 dtfil%getwfk_from_image =0 944 dtfil%getden_from_image =0 945 dtfil%getkden_from_image =0 946 dtfil%getpawden_from_image=0 947 end if 948 949 DBG_EXIT("COLL") 950 951 end subroutine dtfil_init
m_dtfil/dtfil_init_img [ Functions ]
[ Top ] [ m_dtfil ] [ Functions ]
NAME
dtfil_init_img
FUNCTION
Initialize few scalars in the dtfil structured variable when an alogrithm using image of the cell is selected. (initialize index of images from which read files)
INPUTS
dtset=<type datasets_type>=input variables for the current dataset dtsets(0:ndtset_alloc)=<type datasets_type>=input variables for all datasets idtset=number of the dataset jdtset(0:ndtset)=actual index of the datasets ndtset=number of datasets ndtset_alloc=number of datasets, corrected for allocation of at least one data set
OUTPUT
SIDE EFFECTS
dtfil=<type datafiles_type>= only getxxx_from_image flags are modified
SOURCE
1156 subroutine dtfil_init_img(dtfil,dtset,dtsets,idtset,jdtset,ndtset,ndtset_alloc) 1157 1158 !Arguments ------------------------------------ 1159 !scalars 1160 integer, intent(in) :: idtset,ndtset,ndtset_alloc 1161 type(datafiles_type),intent(out) :: dtfil 1162 type(dataset_type),intent(in) :: dtset 1163 !arrays 1164 integer,intent(in) :: jdtset(0:ndtset) 1165 type(dataset_type),intent(in) :: dtsets(0:ndtset_alloc) 1166 1167 !Local variables ------------------------- 1168 !scalars 1169 integer :: iget 1170 1171 ! ********************************************************************* 1172 1173 DBG_ENTER("COLL") 1174 1175 !Default values 1176 dtfil%getwfk_from_image =0 ! Get standard WFK from previous dataset 1177 dtfil%getden_from_image =0 ! Get standard DEN from previous dataset 1178 dtfil%getkden_from_image =0 ! Get standard KDEN from previous dataset 1179 dtfil%getpawden_from_image=0 ! Get standard PAWDEN from previous dataset 1180 1181 if (dtset%optdriver==RUNL_GSTATE.and.dtset%nimage>1) then 1182 1183 ! Define getwfk_from_image 1184 if (dtset%getwfk/=0.or.dtset%irdwfk/=0) then 1185 iget=-1 1186 if(dtset%getwfk<0) iget=jdtset(idtset+dtset%getwfk) 1187 if(dtset%getwfk>0) iget=dtset%getwfk 1188 if(dtset%irdwfk>0) iget=0 1189 if (iget>=0) then 1190 if (iget==0.or.dtsets(iget)%nimage==dtset%nimage) then 1191 dtfil%getwfk_from_image=-1 ! Get WFK from the same image of previous dataset 1192 else if (dtsets(iget)%nimage>1) then 1193 dtfil%getwfk_from_image=1 ! Get WFK from the first image of previous dataset 1194 end if 1195 end if 1196 end if 1197 1198 ! Define getden_from_image 1199 if (dtset%getden/=0.or.dtset%irdden/=0) then 1200 iget=-1 1201 if(dtset%getden<0) iget=jdtset(idtset+dtset%getden) 1202 if(dtset%getden>0) iget=dtset%getden 1203 if(dtset%irdden>0) iget=0 1204 if (iget>=0) then 1205 if (iget==0.or.dtsets(iget)%nimage==dtset%nimage) then 1206 dtfil%getden_from_image=-1 ! Get DEN from the same image of previous dataset 1207 else if (dtsets(iget)%nimage>1) then 1208 dtfil%getden_from_image=1 ! Get DEN from the first image of previous dataset 1209 end if 1210 end if 1211 end if 1212 1213 ! Define getkden_from_image 1214 if (dtset%getkden/=0.or.dtset%irdkden/=0) then 1215 iget=-1 1216 if(dtset%getkden<0) iget=jdtset(idtset+dtset%getkden) 1217 if(dtset%getkden>0) iget=dtset%getkden 1218 if(dtset%irdkden>0) iget=0 1219 if (iget>=0) then 1220 if (iget==0.or.dtsets(iget)%nimage==dtset%nimage) then 1221 dtfil%getkden_from_image=-1 ! Get KDEN from the same image of previous dataset 1222 else if (dtsets(iget)%nimage>1) then 1223 dtfil%getkden_from_image=1 ! Get KDEN from the first image of previous dataset 1224 end if 1225 end if 1226 end if 1227 1228 ! Define getpawden_from_image 1229 if (dtset%getpawden/=0.or.dtset%irdpawden/=0) then 1230 iget=-1 1231 if(dtset%getpawden<0) iget=jdtset(idtset+dtset%getpawden) 1232 if(dtset%getpawden>0) iget=dtset%getpawden 1233 if(dtset%irdpawden>0) iget=0 1234 if (iget>=0) then 1235 if (iget==0.or.dtsets(iget)%nimage==dtset%nimage) then 1236 dtfil%getpawden_from_image=-1 ! Get PAWDEN from the same image of previous dataset 1237 else if (dtsets(iget)%nimage>1) then 1238 dtfil%getpawden_from_image=1 ! Get PAWDEN from the first image of previous dataset 1239 end if 1240 end if 1241 end if 1242 end if 1243 1244 DBG_EXIT("COLL") 1245 1246 end subroutine dtfil_init_img
m_dtfil/dtfil_init_time [ Functions ]
[ Top ] [ m_dtfil ] [ Functions ]
NAME
dtfil_init_time
FUNCTION
Inside the itimimage, iimage and itime loops (this is only needed for optdriver=0), initialize the remaining parts of dtfil.
INPUTS
iapp=indicates the eventual suffix to be appended to the generic output root if 0 : no suffix to be appended (called directly from gstate) if positive : append "_TIM//iapp" (called from move or brdmin) if -1 : append "_TIM0" (called from brdmin) if -2, -3, -4, -5: append "_TIMA", ... ,"_TIMD", (called from move)
OUTPUT
SIDE EFFECTS
dtfil=<type datafiles_type>infos about file names, file unit numbers (part of which were initialized previously)
SOURCE
978 subroutine dtfil_init_time(dtfil,iapp) 979 980 !Arguments ------------------------------------ 981 !scalars 982 integer, intent(in) :: iapp 983 type(datafiles_type),intent(inout) :: dtfil 984 985 !Local variables------------------------------- 986 !scalars 987 character(len=fnlen) :: filapp,filprot 988 989 !****************************************************************** 990 991 DBG_ENTER("COLL") 992 993 !-------------------------------------------------------- 994 !Names based on dtfil%filnam_ds(4)+iapp 995 996 !Prepare the name of the auxiliary files DOS, EIG... 997 call fappnd(filapp,dtfil%filnam_ds(4),iapp) 998 dtfil%fnameabo_app=trim(filapp) 999 dtfil%fnameabo_app_atmden_core=trim(filapp)//'_ATMDEN_CORE' 1000 dtfil%fnameabo_app_atmden_val=trim(filapp)//'_ATMDEN_VAL' 1001 dtfil%fnameabo_app_atmden_full=trim(filapp)//'_ATMDEN_FULL' 1002 dtfil%fnameabo_app_n_tilde=trim(filapp)//'_N_TILDE' 1003 dtfil%fnameabo_app_n_one=trim(filapp)//'_N_ONE' 1004 dtfil%fnameabo_app_nt_one=trim(filapp)//'_NT_ONE' 1005 dtfil%fnameabo_app_bxsf=trim(filapp)//'_BXSF' 1006 dtfil%fnameabo_app_cif=trim(filapp)//'.cif' 1007 dtfil%fnameabo_app_den=trim(filapp)//'_DEN' 1008 dtfil%fnameabo_app_dos=trim(filapp)//'_DOS' 1009 dtfil%fnameabo_app_eig=trim(filapp)//'_EIG' 1010 dtfil%fnameabo_app_elf=trim(filapp)//'_ELF' 1011 dtfil%fnameabo_app_elf_down=trim(filapp)//'_ELF_DOWN' 1012 dtfil%fnameabo_app_elf_up=trim(filapp)//'_ELF_UP' 1013 dtfil%fnameabo_app_fatbands=trim(filapp)//'_FATBANDS' 1014 dtfil%fnameabo_app_gden1=trim(filapp)//'_GDEN1' 1015 dtfil%fnameabo_app_gden2=trim(filapp)//'_GDEN2' 1016 dtfil%fnameabo_app_gden3=trim(filapp)//'_GDEN3' 1017 dtfil%fnameabo_app_geo=trim(filapp)//'_GEO' 1018 dtfil%fnameabo_app_kden=trim(filapp)//'_KDEN' 1019 dtfil%fnameabo_app_lden=trim(filapp)//'_LDEN' 1020 dtfil%fnameabo_app_nesting=trim(filapp)//'_NEST' 1021 dtfil%fnameabo_app_opt=trim(filapp)//'_OPT' 1022 dtfil%fnameabo_app_opt2=trim(filapp)//'_OPT2' 1023 dtfil%fnameabo_app_pawden=trim(filapp)//'_PAWDEN' 1024 dtfil%fnameabo_app_pot=trim(filapp)//'_POT' 1025 dtfil%fnameabo_app_stm=trim(filapp)//'_STM' 1026 dtfil%fnameabo_app_vclmb=trim(filapp)//'_VCLMB' 1027 dtfil%fnameabo_app_vha=trim(filapp)//'_VHA' 1028 dtfil%fnameabo_app_vhxc=trim(filapp)//'_VHXC' 1029 dtfil%fnameabo_app_vpsp=trim(filapp)//'_VPSP' 1030 dtfil%fnameabo_app_vxc=trim(filapp)//'_VXC' 1031 dtfil%fnameabo_app_wfk=trim(filapp)//'_WFK' 1032 dtfil%fnameabo_app_vha_1dm=trim(filapp)//'_VHA_1DM' 1033 dtfil%fnameabo_app_vclmb_1dm=trim(filapp)//'_VCLMB_1DM' 1034 dtfil%fnameabo_app_1dm=trim(filapp)//'_1DM' 1035 1036 !-------------------------------------------------------- 1037 !Names based on dtfil%filnam_ds(5)+iapp 1038 1039 !Prepare the name of the auxiliary files for protection 1040 call fappnd(filprot,dtfil%filnam_ds(5),iapp) 1041 dtfil%fnametmp_app_den=trim(filprot)//'_DEN' 1042 dtfil%fnametmp_app_kden=trim(filprot)//'_KDEN' 1043 1044 DBG_EXIT("COLL") 1045 1046 end subroutine dtfil_init_time
m_dtfil/fappnd [ Functions ]
[ Top ] [ m_dtfil ] [ Functions ]
NAME
fappnd
FUNCTION
Create the modified root name to be used for output of density, potential, and geometry files. See the description of the iapp input variable.
INPUTS
filnam= generic output root name iapp=indicates the eventual suffix to be appended to the generic output root (the suffixe depends on the presence of the suff (optional) argument. if 0 : no suffix to be appended (called directly from gstate) if positive : append "_SUF//iapp" (called from move or brdmin) if -1 : append "_SUF0" (called from brdmin) if -2, -3, -4, -5: append "_SUFA", ... ,"_SUFD", (called from move) SUF can be TIM (default) or IMG [suff]= --optional argument--indicates the suffixe to be appended: SUF=TIM (default) or SUF=IMG or ...
OUTPUT
filapp= filename with appended string
SOURCE
1075 subroutine fappnd(filapp,filnam,iapp,& 1076 & suff) ! optional argument 1077 1078 !Arguments ------------------------------------ 1079 !scalars 1080 integer,intent(in) :: iapp 1081 character(len=fnlen),intent(in) :: filnam 1082 character(len=fnlen),intent(out) :: filapp 1083 character(len=3),optional,intent(in) :: suff 1084 1085 !Local variables------------------------------- 1086 !scalars 1087 integer :: ndig 1088 character(len=3) :: suffixe 1089 character(len=8) :: nchar 1090 character(len=500) :: msg 1091 1092 ! ************************************************************************* 1093 1094 if(iapp==0)then 1095 filapp=trim(filnam) 1096 else 1097 suffixe="TIM" 1098 if (present(suff)) suffixe=trim(suff(1:3)) 1099 if(iapp>0)then 1100 ! Create character string for filename. Determine the number of digits in iapp. 1101 ndig=int(log10(dble(iapp)+0.5_dp))+1 1102 ! Make integer format field of exact size (internal write) 1103 ! for assumed nchar string of 8 characters 1104 write(nchar, '(i8)' ) iapp 1105 if (ndig>8) then 1106 write(msg,'(5a,i0,2a,i0,2a)')& 1107 'Requested file name extension has more than the allowed 8 digits.',ch10,& 1108 'Action: resubmit the job with smaller value for ntime.',ch10,& 1109 'Value computed here was ndig=',ndig,ch10,& 1110 'iapp= ',iapp,' filnam= ',trim(filnam) 1111 ABI_ERROR(msg) 1112 end if 1113 ! Concatenate into character string, picking off exact number of digits 1114 ! The potential or density label will be appended in ioarr 1115 filapp=trim(filnam)//'_'//suffixe(1:3)//nchar(9-ndig:8) 1116 else if(iapp==-1)then 1117 filapp=trim(filnam)//'_'//suffixe(1:3)//'0' 1118 else if(iapp==-2)then 1119 filapp=trim(filnam)//'_'//suffixe(1:3)//'A' 1120 else if(iapp==-3)then 1121 filapp=trim(filnam)//'_'//suffixe(1:3)//'B' 1122 else if(iapp==-4)then 1123 filapp=trim(filnam)//'_'//suffixe(1:3)//'C' 1124 else if(iapp==-5)then 1125 filapp=trim(filnam)//'_'//suffixe(1:3)//'D' 1126 end if 1127 end if 1128 1129 end subroutine fappnd
m_dtfil/iofn1 [ Functions ]
[ Top ] [ m_dtfil ] [ Functions ]
NAME
iofn1
FUNCTION
Begin by eventual redefinition of unit std_in and std_out Then, print greetings for interactive user. Next, read filenames from unit std_in, AND check that new output file does not already exist.
INPUTS
input_path: String with input file path. Empty string activates files file legacy mode. comm=MPI communicator.
OUTPUT
character(len=fnlen) :: filnam(5)=character strings giving file names character(len=fnlen) :: filstat=character strings giving name of status file
NOTES
If it does exist, isfile will create a new name to avoid overwriting the output file. Also create name of status file File names refer to following files, in order: (1) Formatted input file (std_in) (2) Formatted output file (std_out) (3) Root name for generic input files (wavefunctions, potential, density ...) (4) Root name for generic output files (wavefunctions, potential, density, DOS, hessian ...) (5) Root name for generic temporary files (wftmp1,wftmp2,kgunit,status ...)
SOURCE
1538 subroutine iofn1(input_path, filnam, filstat, comm) 1539 1540 !Arguments ------------------------------------ 1541 integer,intent(in) :: comm 1542 character(len=fnlen), intent(in) :: input_path 1543 character(len=fnlen), intent(out) :: filstat 1544 character(len=fnlen), intent(out) :: filnam(5) 1545 1546 !Local variables------------------------------- 1547 character(len=1) :: blank 1548 integer,parameter :: master = 0 1549 integer :: me, ios, nproc, ierr, ndtset, lenstr, marr, jdtset, tread, i1,i2 1550 logical :: ex 1551 character(len=fnlen) :: fillog, tmpfil, fname 1552 character(len=10) :: tag 1553 character(len=500) :: msg, errmsg 1554 character(len=strlen) :: string 1555 !arrays 1556 integer,allocatable :: intarr(:) 1557 real(dp),allocatable :: dprarr(:) 1558 !************************************************************************* 1559 1560 ! NOTE: In this routine it's very important to perform tests 1561 ! on possible IO errors (err=10, iomsg) because we are initializing the IO stuff 1562 ! It there's some problem with the hardware or some misconfiguration, 1563 ! it's very likely that the code will crash here and we should try to give useful error messages. 1564 1565 blank = ' '; tmpfil = '' 1566 1567 ! Determine who I am in comm 1568 me = xmpi_comm_rank(comm); nproc = xmpi_comm_size(comm) 1569 1570 ! Define values of do_write_log and do_write_status parameters 1571 ! if a _NOLOG file exists no LOG file and no STATUS file are created for each cpu core 1572 ! if a _LOG file exists, a LOG file and a STATUS file are created for each cpu core 1573 ! if the #_of_cpu_core>NPROC_NO_EXTRA_LOG OR presence of ABI_MAIN_LOG_FILE, LOG file is only created for master proc 1574 ! if the #_of_cpu_core>NPROC_NO_EXTRA_STATUS OR presence of ABI_MAIN_LOG_FILE, STATUS file is only created for master proc 1575 inquire(file=ABI_NO_LOG_FILE, iostat=ios, exist=ex) 1576 if (ios /= 0) ex=.false. 1577 if (ex) then 1578 do_write_log=.false. ; do_write_status=.false. 1579 call abi_log_status_state(new_do_write_log=.false.,new_do_write_status=.false.) 1580 call libpaw_log_flag_set(.false.) 1581 else 1582 inquire(file=ABI_ENFORCE_LOG_FILE, iostat=ios, exist=ex) 1583 if (ios/=0) ex=.false. 1584 if (ex) then 1585 do_write_log=.true. ; do_write_status=.true. 1586 call abi_log_status_state(new_do_write_log=.true.,new_do_write_status=.true.) 1587 call libpaw_log_flag_set(.true.) 1588 else 1589 inquire(file=ABI_MAIN_LOG_FILE, iostat=ios, exist=ex) 1590 if (ios /= 0) ex=.false. 1591 if (ex .and. me /= master) then 1592 do_write_log=.false. ; do_write_status=.false. 1593 call abi_log_status_state(new_do_write_log=.false.,new_do_write_status=.false.) 1594 call libpaw_log_flag_set(.false.) 1595 else 1596 if (me /= master) then 1597 do_write_log= (nproc<NPROC_NO_EXTRA_LOG) 1598 call abi_log_status_state(new_do_write_log=(nproc<NPROC_NO_EXTRA_LOG)) 1599 call libpaw_log_flag_set((nproc<NPROC_NO_EXTRA_LOG)) 1600 do_write_status= (nproc<NPROC_NO_EXTRA_STATUS) 1601 call abi_log_status_state(new_do_write_status=(nproc<NPROC_NO_EXTRA_STATUS)) 1602 end if 1603 end if 1604 end if 1605 end if 1606 !do_write_log = .True. 1607 1608 if (me == master) then 1609 ! Eventually redefine standard input and standard output 1610 if (do_write_log) then 1611 #if defined READ_FROM_FILE 1612 ! Take care of the output file 1613 tmpfil(1:fnlen)=blank 1614 tmpfil(1:3)='log' 1615 call isfile(tmpfil,'new') 1616 close(std_out, err=10, iomsg=errmsg) 1617 if (open_file(tmpfil,msg,unit=std_out,form='formatted',status='new',action="write") /= 0) then 1618 ABI_ERROR(msg) 1619 end if 1620 #endif 1621 else 1622 ! Redirect standard output to null 1623 close(std_out, err=10, iomsg=errmsg) 1624 if (open_file(NULL_FILE, msg, unit=std_out, action="write") /= 0) then 1625 ABI_ERROR(msg) 1626 end if 1627 end if 1628 1629 #if defined READ_FROM_FILE 1630 ! Now take care of the "files" file 1631 tmpfil(1:fnlen)=blank 1632 tmpfil(1:9)='ab.files' 1633 write(msg, '(4a)' )& 1634 'Because of CPP option READ_FROM_FILE,',ch10,& 1635 'read file "ab.files" instead of standard input ' ,ch10 1636 ABI_COMMENT(msg) 1637 call isfile(tmpfil,'old') 1638 close(std_in, err=10, iomsg=errmsg) 1639 if (open_file(tmpfil,msg,unit=std_in,form='formatted',status='old',action="read") /= 0) then 1640 ABI_ERROR(msg) 1641 end if 1642 #endif 1643 1644 ! Print greetings for interactive user 1645 write(std_out,*,err=10,iomsg=errmsg)' ABINIT ',trim(abinit_version) 1646 write(std_out,*,err=10,iomsg=errmsg)' ' 1647 1648 if (len_trim(input_path) == 0) then 1649 ! Legacy Files file mode. 1650 write(std_out, "(2a)")" DeprecationWarning: ",ch10 1651 write(std_out, "(a)") " The files file has been deprecated in Abinit9 and will be removed in Abinit10." 1652 write(std_out, "(2a)")" Use the syntax `abinit t01.abi` where t01.abi is an input with pseudopotentials e.g.",ch10 1653 write(std_out, "(3a)")' pseudos = "al.psp8, as.psp8"',ch10,ch10 1654 1655 write(std_out,*,err=10,iomsg=errmsg)' Give name for formatted input file: ' 1656 read(std_in, '(a)',err=10,iomsg=errmsg ) filnam(1) 1657 write(std_out, '(a)',err=10,iomsg=errmsg ) trim(filnam(1)) 1658 write(std_out,*)' Give name for formatted output file:' 1659 read (std_in, '(a)',err=10,iomsg=errmsg ) filnam(2) 1660 write (std_out, '(a)',err=10,iomsg=errmsg ) trim(filnam(2)) 1661 write(std_out,*)' Give root name for generic input files:' 1662 read (std_in, '(a)',err=10,iomsg=errmsg ) filnam(3) 1663 write (std_out, '(a)',err=10,iomsg=errmsg ) trim(filnam(3)) 1664 write(std_out,*, err=10, iomsg=errmsg )' Give root name for generic output files:' 1665 read (std_in, '(a)', err=10, iomsg=errmsg ) filnam(4) 1666 write (std_out, '(a)', err=10, iomsg=errmsg ) trim(filnam(4)) 1667 write(std_out,*, err=10, iomsg=errmsg)' Give root name for generic temporary files:' 1668 read (std_in, '(a)', err=10, iomsg=errmsg ) filnam(5) 1669 write (std_out, '(a)', err=10, iomsg=errmsg ) trim(filnam(5)) 1670 1671 else 1672 ! Get prefix from input file. Default values are provided 1673 filnam(1) = input_path 1674 filnam(2) = trim(input_path)//".abo" 1675 filnam(3) = "i" 1676 filnam(4) = "o" 1677 filnam(5) = "t" 1678 1679 fname = basename(input_path) 1680 i1 = index(fname, ".") 1681 !if (i1 /= 0) then 1682 if (i1 > 1) then 1683 ! file ext is present --> use prefix to initialize filnam 1684 i2 = index(input_path, ".", back=.True.) 1685 filnam(2) = input_path(:i2) // "abo" 1686 filnam(3) = fname(:i1-1) // "i" 1687 filnam(4) = fname(:i1-1) // "o" 1688 filnam(5) = fname(:i1-1) // "t" 1689 end if 1690 1691 ! Read the file, stringify it and return the number of datasets. 1692 call parsefile(input_path, lenstr, ndtset, string, xmpi_comm_self) 1693 1694 marr = max(1, ndtset) 1695 ABI_MALLOC(dprarr, (marr)) 1696 ABI_MALLOC(intarr, (marr)) 1697 jdtset = 0 1698 1699 ! Allow user to change default values 1700 call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), "output_file", tread, 'KEY', key_value=filnam(2)) 1701 call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), "indata_prefix", tread, 'KEY', key_value=filnam(3)) 1702 call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), "outdata_prefix", tread, 'KEY', key_value=filnam(4)) 1703 call intagm(dprarr, intarr, jdtset, marr, 1, string(1:lenstr), "tmpdata_prefix", tread, 'KEY', key_value=filnam(5)) 1704 1705 ABI_FREE(dprarr) 1706 ABI_FREE(intarr) 1707 end if 1708 1709 ! Check that old input file exists 1710 call isfile(filnam(1), 'old') 1711 ! Check that new output file does NOT exist 1712 call isfile(filnam(2), 'new') 1713 1714 ! Check that root name for generic input and output differ 1715 if ( trim(filnam(3)) == trim(filnam(4)) ) then 1716 write(msg, '(3a)' )& 1717 'Root name for generic input and output files must differ ',ch10,& 1718 'Action: correct your "file" file.' 1719 ABI_ERROR(msg) 1720 end if 1721 1722 ! Check that root names are at least 20 characters less than fnlen 1723 if ( len_trim(filnam(3)) >= (fnlen-20) ) then 1724 write(msg, '(a,a,a,a,a,i0,a,i0,a,a)' )& 1725 'Root name for generic input files is too long. ',ch10,& 1726 'It must be 20 characters less than the maximal allowed ',ch10,& 1727 'length of names, that is ',fnlen,', while it is: ',len_trim(filnam(3)),ch10,& 1728 'Action: correct your "file" file.' 1729 ABI_ERROR(msg) 1730 end if 1731 if ( len_trim(filnam(4)) >= (fnlen-20) ) then 1732 write(msg, '(a,a,a,a,a,i0,a,i0,a,a)' )& 1733 'Root name for generic output files is too long. ',ch10,& 1734 'It must be 20 characters less than the maximal allowed ',ch10,& 1735 'length of names, that is ',fnlen,', while it is: ',len_trim(filnam(4)),ch10,& 1736 'Action: correct your "file" file.' 1737 ABI_ERROR(msg) 1738 end if 1739 if ( len_trim(filnam(5)) >= (fnlen-20) ) then 1740 write(msg, '(a,a,a,a,a,i0,a,i0,a,a)' )& 1741 'Root name for generic temporary files is too long. ',ch10,& 1742 'It must be 20 characters less than the maximal allowed ',ch10,& 1743 'length of names, that is ',fnlen,', while it is: ',len_trim(filnam(5)),ch10,& 1744 'Action: correct your "file" file.' 1745 ABI_ERROR(msg) 1746 end if 1747 1748 ! TODO: Create directories if needed but I need C routines to be portable. 1749 !i = index(filnam(5), "/"); if (i > 0) call mkdir(filnam(5)(1:i-1), ierr) 1750 !i = index(filnam(6), "/"); if (i > 0) call mkdir(filnam(6)(1:i-1), ierr) 1751 1752 end if ! master only 1753 1754 ! Communicate filenames to all processors 1755 call xmpi_bcast(filnam,master,comm,ierr) 1756 1757 ! Create a name for the status file, based on filnam(5) 1758 filstat=trim(filnam(5))//'_STATUS' 1759 1760 ! Redefine the log unit if not the master 1761 if (me /= master) then 1762 call int2char4(me,tag) 1763 ABI_CHECK((tag(1:1)/='#'),'Bug: string length too short!') 1764 filstat=trim(filstat)//'_P-'//trim(tag) 1765 if (do_write_log) then 1766 fillog=trim(filnam(5))//'_LOG_'//trim(tag) 1767 close(std_out, err=10, iomsg=errmsg) 1768 if (open_file(fillog,msg,unit=std_out,status='unknown',action="write") /= 0) then 1769 ABI_ERROR(msg) 1770 end if 1771 else 1772 close(std_out, err=10, iomsg=errmsg) 1773 if (open_file(NULL_FILE,msg,unit=std_out,action="write") /= 0) then 1774 ABI_ERROR(msg) 1775 end if 1776 end if 1777 end if 1778 1779 call xmpi_barrier(comm) 1780 return 1781 1782 ! Handle possibe IO errors 1783 10 continue 1784 ABI_ERROR(errmsg) 1785 1786 end subroutine iofn1
m_dtfil/isfile [ Functions ]
[ Top ] [ m_dtfil ] [ Functions ]
NAME
isfile
FUNCTION
Inquire Status of FILE Checks that for status = 'old': file already exists 'new': file does not exist; if file exists, filnam is modified to filnam.A or filnam.B,....
INPUTS
filnam=character string to specify filename status='old' or 'new'
OUTPUT
stops processing if old file does not exist; changes name and returns new name in redefined filnam if new file already exists.
SOURCE
1410 subroutine isfile(filnam, status) 1411 1412 !Arguments ------------------------------------ 1413 !scalars 1414 character(len=3),intent(in) :: status 1415 character(len=fnlen),intent(inout) :: filnam 1416 1417 !Local variables------------------------------- 1418 !scalars 1419 logical :: ex,found 1420 integer :: ii,ios, ioserr 1421 character(len=500) :: msg 1422 character(len=fnlen) :: filnam_tmp, trialnam 1423 1424 ! ************************************************************************* 1425 1426 filnam_tmp=filnam 1427 1428 if (status=='old') then 1429 ! Check that old file exists 1430 inquire(file=filnam,iostat=ios,exist=ex) 1431 1432 if (ios/=0) then 1433 write(msg,'(4a,i0,2a)')& 1434 'Checks for existence of file: ',trim(filnam),ch10,& 1435 'but INQUIRE statement returns error code',ios,ch10,& 1436 'Action: identify which problem appears with this file.' 1437 ABI_ERROR(msg) 1438 else if (.not.ex) then 1439 write(msg, '(5a)' )& 1440 'Checks for existence of file: ',trim(filnam),ch10,& 1441 'but INQUIRE finds file does not exist.',& 1442 'Action: check file name and re-run.' 1443 ABI_ERROR(msg) 1444 end if 1445 1446 else if (status=='new') then 1447 1448 ! Check that new output file does NOT exist 1449 ioserr = 0 1450 trialnam = filnam 1451 ii = 0 1452 inquire(file=trim(trialnam),iostat=ios,exist=ex) 1453 if (ios /= 0) then 1454 write(msg,'(3a)') & 1455 'Something is wrong with permissions for reading/writing on this filesystem.',ch10,& 1456 'Action: Check permissions.' 1457 ABI_ERROR(msg) 1458 end if 1459 1460 if ( ex .eqv. .true. ) then 1461 write(msg,'(3a)')'Output file: ',trim(trialnam),' already exists.' 1462 ABI_COMMENT(msg) 1463 found=.false. 1464 1465 ii=1 1466 do while ( (found .eqv. .false.) .and. (ii < 10000) ) 1467 call int2char4(ii,msg) 1468 trialnam=trim(trim(filnam_tmp)//msg) 1469 inquire(file=trim(trialnam),iostat=ios,exist=ex) 1470 if ( (ex .eqv. .false.) .and. (ios == 0)) then 1471 found = .true. 1472 end if 1473 if ( ios /= 0 ) ioserr=ioserr+1 1474 if ( ioserr > 10 ) then 1475 ! There is a problem => stop 1476 write(msg, '(2a,i0,2a)' )& 1477 'Check for permissions of reading/writing files on the filesystem', & 1478 '10 INQUIRE statements returned an error code like ',ios,ch10,& 1479 'Action: Check permissions' 1480 ABI_ERROR(msg) 1481 end if 1482 ii=ii+1 1483 end do 1484 if ( found .eqv. .true. ) then 1485 write(msg,'(4a)') 'Renaming old: ',trim(filnam),' to: ',trim(trialnam) 1486 ABI_COMMENT(msg) 1487 ioserr = clib_rename(filnam, trialnam) 1488 if ( ioserr /= 0 ) then 1489 write(msg,'(4a)') 'Failed to rename file: ', trim(filnam),' to: ',trim(trialnam) 1490 ABI_ERROR(msg) 1491 end if 1492 else 1493 write(msg,'(3a)')& 1494 'Have used all names of the form filenameXXXX, X in [0-9]',ch10,& 1495 'Action: clean up your directory and start over.' 1496 ABI_ERROR(msg) 1497 end if 1498 end if 1499 ! if ii > 0 we iterated so rename abi_out to abi_outXXXX and just write to abi_out 1500 else 1501 ABI_BUG(sjoin('Input status:', status, ' not recognized.')) 1502 end if 1503 1504 end subroutine isfile
m_dtfil/mkfilename [ Functions ]
[ Top ] [ m_dtfil ] [ Functions ]
NAME
mkfilename
FUNCTION
From the root (input or output) file names, produce a real file name.
INPUTS
filnam(5)=the root file names (only filnam(3) and filnam(4) are really needed) get=input 'get variable', if 1, must get the file from another dataset idtset=number of the dataset ird=input 'iread variable', if 1, must get the file from the input root jdtset_(0:ndtset)=actual index of the dataset ndtset=number of datasets stringfil=the string of characters to be appended e.g. '_WFK' or '_DEN' stringvar=the string of characters to be appended that defines the 'get' or 'ird' variables, e.g. 'wfk' or 'ddk' [getpath]=String with filename to be used as input, exclude get and ird option.
OUTPUT
filnam_out=the new file name will_read=1 if the file must be read ; 0 otherwise (ird and get were zero)
SOURCE
1275 subroutine mkfilename(filnam,filnam_out,get,idtset,ird,jdtset_,ndtset,stringfil,stringvar,will_read, & 1276 getpath) ! Optional 1277 1278 !Arguments ------------------------------------ 1279 !scalars 1280 integer,intent(in) :: get,idtset,ird,ndtset 1281 integer,intent(out) :: will_read 1282 character(len=*),intent(in) :: stringfil 1283 character(len=*),intent(in) :: stringvar 1284 character(len=fnlen),intent(out) :: filnam_out 1285 character(len=fnlen),optional,intent(in) :: getpath 1286 !arrays 1287 integer,intent(in) :: jdtset_(0:ndtset) 1288 character(len=fnlen),intent(in) :: filnam(5) 1289 1290 !Local variables------------------------------- 1291 !scalars 1292 integer :: jdtset,jget 1293 character(len=4) :: appen 1294 character(len=500) :: msg 1295 character(len=fnlen) :: filnam_appen 1296 1297 ! ************************************************************************* 1298 1299 ! Here, defaults if no get variable 1300 will_read=ird 1301 1302 filnam_appen = trim(filnam(3)) 1303 if (ndtset > 0) then 1304 jdtset = jdtset_(idtset) 1305 call appdig(jdtset, '', appen) 1306 filnam_appen = trim(filnam_appen)//'_DS'//appen 1307 end if 1308 filnam_out = trim(filnam_appen)//trim(stringfil) 1309 1310 if (present(getpath)) then 1311 if (getpath /= ABI_NOFILE) then 1312 if (ird /= 0 .or. get /= 0)then 1313 write(msg, '(11a,i0,3a,i0,a,i0,7a)' ) & 1314 'When the input variable: ', trim(getpath), ' is used ',ch10, & 1315 'the input variables ird',trim(stringvar),' and get',trim(stringvar),' cannot be',ch10,& 1316 'simultaneously non-zero, while for idtset = ',idtset,',',ch10,& 1317 'they are ',ird,', and ',get,'.',ch10,& 1318 'Action: correct ird',trim(stringvar),' or get',trim(stringvar),' in your input file.' 1319 ABI_ERROR(msg) 1320 end if 1321 filnam_out = rmquotes(getpath) 1322 write(msg, '(5a)')' mkfilename: get', trim(stringvar), " from: ",trim(filnam_out), ch10 1323 call wrtout([std_out, ab_out], msg) 1324 ! Check whether file exists taking into account a possible NC file extension. 1325 if (xmpi_comm_rank(xmpi_world) == 0) then 1326 if (.not. file_exists(filnam_out) .and. .not. file_exists(strcat(filnam_out, ".nc"))) then 1327 ABI_ERROR(sjoin("Cannot find file:", filnam_out, "(with or without .nc extension)")) 1328 end if 1329 end if 1330 will_read = 1; return 1331 end if 1332 end if 1333 1334 ! Treatment of the multi-dataset case (get is not relevant otherwise) 1335 if (ndtset /= 0) then 1336 1337 if(ndtset==1 .and. get<0 .and. (jdtset_(1)+get>0)) then 1338 write(msg, '(7a,i0,a,i0,5a)' )& 1339 'You cannot use a negative value of get',trim(stringvar),' with only 1 dataset!',ch10, & 1340 'If you want to refer to a previously computed dataset,',ch10, & 1341 'you should give the absolute index of it (i.e. ', jdtset_(idtset)+get,' instead of ',get,').',ch10, & 1342 'Action: correct get',trim(stringvar),' in your input file.' 1343 ABI_ERROR(msg) 1344 end if 1345 1346 if (idtset + get < 0) then 1347 write(msg, '(5a,i0,3a,i0,4a)' )& 1348 'The sum of idtset and get',trim(stringvar),' cannot be negative,',ch10,& 1349 'while they are idtset = ',idtset,', and get',trim(stringvar),' = ',get,ch10,& 1350 'Action: correct get',trim(stringvar),' in your input file.' 1351 ABI_ERROR(msg) 1352 end if 1353 1354 if(get>0 .or. (get<0 .and. idtset+get>0) )then 1355 1356 if(ird/=0 .and. get/=0)then 1357 write(msg, '(7a,i0,3a,i0,a,i0,7a)' )& 1358 'The input variables ird',trim(stringvar),' and get',trim(stringvar),' cannot be',ch10,& 1359 'simultaneously non-zero, while for idtset = ',idtset,',',ch10,& 1360 'they are ',ird,', and ',get,'.',ch10,& 1361 'Action: correct ird',trim(stringvar),' or get',trim(stringvar),' in your input file.' 1362 ABI_ERROR(msg) 1363 end if 1364 1365 will_read=1 1366 1367 ! Compute the dataset from which to take the file, and the corresponding index 1368 if(get<0 .and. idtset+get>0) jget=jdtset_(idtset+get) 1369 if(get>0) jget=get 1370 call appdig(jget,'',appen) 1371 1372 ! Note use of output filename (filnam(4)) 1373 filnam_out=trim(filnam(4))//'_DS'//trim(appen)//trim(stringfil) 1374 1375 if(jdtset>=100)then 1376 write(msg, '(5a,i5,2a)' )& 1377 ' mkfilename : get',trim(stringvar) ,'/=0, take file ',trim(stringfil),' from output of DATASET ',jget,'.',ch10 1378 else 1379 write(msg, '(5a,i3,2a)' )& 1380 ' mkfilename : get',trim(stringvar) ,'/=0, take file ',trim(stringfil),' from output of DATASET ',jget,'.',ch10 1381 end if 1382 call wrtout([std_out, ab_out], msg) 1383 end if ! conditions on get and idtset 1384 end if ! ndtset/=0 1385 1386 end subroutine mkfilename