TABLE OF CONTENTS
ABINIT/dtfil_init [ 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).
COPYRIGHT
Copyright (C) 2010-2018 ABINIT group (XG) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt . For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
INPUTS
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 image_index= -optional argument- 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 jdtset_(0:ndtset)=actual index of the datasets mpi_enreg=information about MPI parallelization ndtset=number of datasets
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.
PARENTS
driver,gstateimg
CHILDREN
appdig,int2char4,mkfilename
SOURCE
51 #if defined HAVE_CONFIG_H 52 #include "config.h" 53 #endif 54 55 #include "abi_common.h" 56 57 subroutine dtfil_init(dtfil,dtset,filnam,filstat,idtset,jdtset_,mpi_enreg,ndtset,& 58 & image_index) ! optional argument 59 60 use defs_basis 61 use defs_abitypes 62 use m_profiling_abi 63 use m_errors 64 use m_xmpi 65 66 use m_fstrings, only : int2char4 67 68 !This section has been created automatically by the script Abilint (TD). 69 !Do not modify the following lines by hand. 70 #undef ABI_FUNC 71 #define ABI_FUNC 'dtfil_init' 72 use interfaces_32_util 73 use interfaces_54_abiutil 74 !End of the abilint section 75 76 implicit none 77 78 !Arguments ------------------------------------ 79 !scalars 80 integer, intent(in) :: idtset,ndtset 81 integer, optional, intent(in) :: image_index 82 character(len=fnlen),intent(in) :: filstat 83 type(MPI_type),intent(in) :: mpi_enreg 84 type(datafiles_type),intent(inout) :: dtfil !vz_i 85 !arrays 86 integer :: jdtset_(0:ndtset) 87 character(len=fnlen),intent(in) :: filnam(5) 88 type(dataset_type),intent(in) :: dtset 89 90 !Local variables------------------------------- 91 !scalars 92 ! Define input and output unit numbers (do not forget, unit 5 and 6 are standard input and output) 93 ! Also, unit number 21, 22 and 23 are used in dfpt_nstdy, for the 3 dot wavefunctions. 94 ! Unit 50,51,52 and 53 are used in dfpt_looppert (for ipert=natom+2, ipert=natom+10 and ipert=natom+11). 95 ! Others unit numbers will be used in the case of the variational and 2n+1 expressions. 96 ! In defs_basis, one defines : 97 ! std_in=5, ab_in=5, std_out=6, ab_out=7, tmp_unit=9, tmp_unit2=10 98 integer,parameter :: unchi0=42,unddb=16,unddk=50,undkdk=54,undkde=55,unkg1=19,unkg=17,unkgq=18 99 integer,parameter :: unpaw=26,unpaw1=27,unpawq=28,unpos=30 100 integer,parameter :: unwff1=1,unwff2=2,unwffgs=3,unwfkq=4,unwft1=11 101 integer,parameter :: unwft2=12,unwftgs=13,unwftkq=14,unylm=24,unylm1=25 102 integer,parameter :: unkss=40,unscr=41,unqps=43 103 integer :: ii,iimage,ireadden,ireadkden,ireadwf,ixx,jdtset,will_read 104 character(len=10) :: appen,tag 105 character(len=9) :: stringvar 106 character(len=15) :: stringfile 107 character(len=500) :: message 108 character(len=fnlen) :: filsus,filddbsin,fildens1in,fildensin,filpawdensin,filkdensin,filqps,filscr 109 character(len=fnlen) :: fnamewff1,fnamewffddk,fnamewffdelfd,fnamewffdkdk,fnamewffdkde,fnamewffk,fnamewffq 110 character(len=fnlen) :: filbseig,filfft,filhaydock,fil_bsreso,fil_bscoup 111 character(len=fnlen) :: filwfkfine 112 character(len=fnlen) :: filnam_ds(5) 113 character(len=fnlen) :: tmpfil(14) 114 integer :: idtmpfil(14) 115 116 !****************************************************************** 117 118 DBG_ENTER("COLL") 119 120 iimage=0;if (present(image_index)) iimage=image_index 121 122 dtfil%unchi0 =unchi0 123 dtfil%unddb =unddb 124 dtfil%unddk =unddk 125 dtfil%undkde =undkde 126 dtfil%undkdk =undkdk 127 dtfil%unkg =unkg 128 dtfil%unkgq =unkgq 129 dtfil%unkg1 =unkg1 130 dtfil%unkss =unkss 131 dtfil%unqps =unqps 132 dtfil%unscr =unscr 133 dtfil%unwff1 =unwff1 134 dtfil%unwff2 =unwff2 135 dtfil%unwffgs=unwffgs 136 dtfil%unwffkq=unwfkq 137 dtfil%unwft1 =unwft1 138 dtfil%unwft2 =unwft2 139 dtfil%unwftgs=unwftgs 140 dtfil%unwftkq=unwftkq 141 dtfil%unylm =unylm 142 dtfil%unylm1 =unylm1 143 dtfil%unpaw =unpaw 144 dtfil%unpaw1 =unpaw1 145 dtfil%unpawq =unpawq 146 dtfil%unpos =unpos 147 filnam_ds(1:5)=filnam(1:5) 148 jdtset=dtset%jdtset 149 150 !If multi dataset mode, special treatment of filenames 3 and 4 (density and 151 !wavefunctions input and output, as well as other output files) 152 if(ndtset>0)then 153 call appdig(jdtset,'',appen) 154 filnam_ds(3)=trim(filnam(3))//'_DS'//trim(appen) 155 filnam_ds(4)=trim(filnam(4))//'_DS'//trim(appen) 156 end if 157 158 !If multi image mode (nimage>1), special treatment of filenames 4 and 5 159 if(iimage>0)then 160 call appdig(iimage,'',appen) 161 filnam_ds(4)=trim(filnam_ds(4))//'_IMG'//trim(appen) 162 filnam_ds(5)=trim(filnam_ds(5))//'_IMG'//trim(appen) 163 end if 164 165 !According to getwfk and irdwfk, build _WFK file name, referred as fnamewffk 166 if (iimage>0.and.dtfil%getwfk_from_image/=0) then 167 if (dtfil%getwfk_from_image==-1) then 168 call appdig(iimage,'',appen) 169 else 170 call appdig(dtfil%getwfk_from_image,'',appen) 171 end if 172 stringfile='_IMG'//trim(appen)//'_WFK' 173 else 174 stringfile='_WFK' 175 end if 176 stringvar='wfk' 177 call mkfilename(filnam,fnamewffk,dtset%getwfk,idtset,dtset%irdwfk,jdtset_,& 178 & ndtset,stringfile,stringvar,will_read) 179 180 if(dtset%optdriver/=RUNL_RESPFN)ireadwf=will_read 181 if(ndtset/=0 .and. dtset%optdriver==RUNL_RESPFN .and. will_read==0)then 182 write(message, '(5a,i3,3a,i3,a,i3,3a)' )& 183 & 'At least one of the input variables irdwfk and getwfk ',ch10,& 184 & 'must refer to a valid _WFK file, in the response function',ch10,& 185 & 'case, while for idtset=',idtset,',',ch10,& 186 & 'they are irdwfk=',dtset%irdwfk,', and getwfk=',dtset%getwfk,'.',ch10,& 187 & 'Action: correct irdwfk or getwfk in your input file.' 188 MSG_ERROR(message) 189 end if 190 191 !Treatment of the other get wavefunction variable, if response function case or nonlinear case 192 if ( ANY(dtset%optdriver == (/RUNL_RESPFN, RUNL_NONLINEAR, RUNL_EPH/)) ) then 193 194 ! According to getwfq and irdwfq, build _WFQ file name, referred as fnamewffq 195 stringfile='_WFQ' ; stringvar='wfq' 196 call mkfilename(filnam,fnamewffq,dtset%getwfq,idtset,dtset%irdwfq,jdtset_,& 197 & ndtset,stringfile,stringvar,will_read) 198 ! If fnamewffq is not initialized thanks to getwfq or irdwfq, use fnamewffk 199 if(will_read==0)fnamewffq=fnamewffk 200 201 ! According to get1wf and ird1wf, build _1WF file name, referred as fnamewff1 202 stringfile='_1WF' ; stringvar='1wf' 203 call mkfilename(filnam,fnamewff1,dtset%get1wf,idtset,dtset%ird1wf,jdtset_,& 204 & ndtset,stringfile,stringvar,will_read) 205 ireadwf=will_read 206 207 ! According to getddk and irdddk, build _1WF file name, referred as fnamewffddk 208 stringfile='_1WF' ; stringvar='ddk' 209 call mkfilename(filnam,fnamewffddk,dtset%getddk,idtset,dtset%irdddk,jdtset_,& 210 & ndtset,stringfile,stringvar,will_read) 211 212 ! According to getdelfd, build _1WF file name, referred as fnamewffdelfd 213 stringfile='_1WF' ; stringvar='delfd' 214 call mkfilename(filnam,fnamewffdelfd,dtset%getdelfd,idtset,0,jdtset_,& 215 & ndtset,stringfile,stringvar,will_read) 216 217 ! According to getdkdk, build _1WF file name, referred as fnamewffdkdk 218 stringfile='_1WF' ; stringvar='dkdk' 219 call mkfilename(filnam,fnamewffdkdk,dtset%getdkdk,idtset,0,jdtset_,& 220 & ndtset,stringfile,stringvar,will_read) 221 222 ! According to getdkde, build _1WF file name, referred as fnamewffdkde 223 stringfile='_1WF' ; stringvar='dkde' 224 call mkfilename(filnam,fnamewffdkde,dtset%getdkde,idtset,0,jdtset_,& 225 & ndtset,stringfile,stringvar,will_read) 226 227 end if 228 229 !------------------------------------------------------------------------------------------- 230 !Build name of files from dtfil%filnam_ds(3) 231 232 !SP :According to getddb, build _DDB file name, referred as filddbsin 233 stringfile='_DDB' 234 stringvar='ddb' 235 call mkfilename(filnam,filddbsin,dtset%getddb,idtset,dtset%irdddb,jdtset_,& 236 & ndtset,stringfile,stringvar,will_read) 237 238 !According to getden, build _DEN file name, referred as fildensin 239 !A default is available if getden is 0 240 if (iimage>0.and.dtfil%getden_from_image/=0) then 241 if (dtfil%getden_from_image==-1) then 242 call appdig(iimage,'',appen) 243 else 244 call appdig(dtfil%getden_from_image,'',appen) 245 end if 246 stringfile='_IMG'//trim(appen)//'_DEN' 247 else 248 stringfile='_DEN' 249 end if 250 stringvar='den' 251 call mkfilename(filnam,fildensin,dtset%getden,idtset,dtset%irdden,jdtset_,& 252 & ndtset,stringfile,stringvar,will_read) 253 254 if(will_read==0)fildensin=trim(filnam_ds(3))//'_DEN' 255 ireadden=will_read 256 257 if ((dtset%optdriver==RUNL_GWLS.or.dtset%optdriver==RUNL_GSTATE) & 258 & .and.dtset%iscf<0) ireadden=1 259 !if (optdriver==RUNL_GSTATE.and.ireadwf/=0) ireadden=0 260 261 !According to getpawden, build _PAWDEN file name, referred as filpawdensin 262 !A default is available if getden is 0 263 if (iimage>0.and.dtfil%getpawden_from_image/=0) then 264 if (dtfil%getpawden_from_image==-1) then 265 call appdig(iimage,'',appen) 266 else 267 call appdig(dtfil%getpawden_from_image,'',appen) 268 end if 269 stringfile='_IMG'//trim(appen)//'_PAWDEN' 270 else 271 stringfile='_PAWDEN' 272 end if 273 stringvar='pawden' 274 call mkfilename(filnam,filpawdensin,dtset%getpawden,idtset,dtset%irdden,jdtset_,& 275 & ndtset,stringfile,stringvar,will_read) 276 if(will_read==0)filpawdensin=trim(filnam_ds(3))//'_PAWDEN' 277 278 !According to getden and usekden, build _KDEN file name, referred as filkdensin 279 !A default is available if getden is 0 280 if(dtset%usekden==1)then 281 if (iimage>0.and.dtfil%getden_from_image/=0) then 282 if (dtfil%getden_from_image==-1) then 283 call appdig(iimage,'',appen) 284 else 285 call appdig(dtfil%getden_from_image,'',appen) 286 end if 287 stringfile='_IMG'//trim(appen)//'_KDEN' 288 else 289 stringfile='_KDEN' 290 end if 291 stringvar='kden' 292 call mkfilename(filnam,filkdensin,dtset%getden,idtset,dtset%irdden,jdtset_,& 293 & ndtset,stringfile,stringvar,will_read) 294 if(will_read==0)filkdensin=trim(filnam_ds(3))//'_KDEN' 295 ireadkden=will_read 296 if ((dtset%optdriver==RUNL_GSTATE.or.dtset%optdriver==RUNL_GWLS).and.dtset%iscf<0) ireadkden=1 297 else 298 ireadkden=0 299 end if 300 301 !According to get1den, build _DEN file name, referred as fildens1in 302 !A default is available if get1den is 0 303 stringfile='_DEN' ; stringvar='1den' 304 call mkfilename(filnam,fildens1in,dtset%get1den,idtset,dtset%ird1den,jdtset_,& 305 & ndtset,stringfile,stringvar,will_read) 306 if(will_read==0)fildens1in=trim(filnam_ds(3))//'_DEN' 307 308 !According to getscr and irdscr, build _SCR file name, referred as filscr 309 !A default is available if getscr is 0 310 stringfile='_SCR' ; stringvar='scr' 311 call mkfilename(filnam,filscr,dtset%getscr,idtset,dtset%irdscr,jdtset_,& 312 & ndtset,stringfile,stringvar,will_read) 313 if(will_read==0)filscr=trim(filnam_ds(3))//'_SCR' 314 315 !According to getsuscep and irdsuscep, build _SUS file name, referred as filsus 316 !A default is available if getsuscep is 0 317 stringfile='_SUS' ; stringvar='sus' 318 call mkfilename(filnam,filsus,dtset%getsuscep,idtset,dtset%irdsuscep,jdtset_,& 319 & ndtset,stringfile,stringvar,will_read) 320 if(will_read==0)filsus=TRIM(filnam_ds(3))//'_SUS' 321 322 !According to getqps and irdqps, build _QPS file name, referred as filqps 323 !A default is available if getqps is 0 324 stringfile='_QPS' ; stringvar='qps' 325 call mkfilename(filnam,filqps,dtset%getqps,idtset,dtset%irdqps,jdtset_,& 326 & ndtset,stringfile,stringvar,will_read) 327 if(will_read==0)filqps=trim(filnam_ds(3))//'_QPS' 328 329 !According to getbseig and irdbseig, build _BSEIG file name, referred as filbseig 330 !A default is available if getbseig is 0 331 stringfile='_BSEIG' ; stringvar='bseig' 332 call mkfilename(filnam,filbseig,dtset%getbseig,idtset,dtset%irdbseig,jdtset_,ndtset,stringfile,stringvar,will_read) 333 if(will_read==0)filbseig=trim(filnam_ds(3))//'_BSEIG' 334 335 !According to gethaydock and irdhaydock, build _HAYD file name, referred as filhaydock. 336 !A default is available if gethaydock is 0 337 stringfile='_HAYDR_SAVE' ; stringvar='haydock' 338 call mkfilename(filnam,filhaydock,dtset%gethaydock,idtset,dtset%irdhaydock,jdtset_,ndtset,stringfile,stringvar,will_read) 339 if(will_read==0)filhaydock=trim(filnam_ds(3))//'_HAYDR_SAVE' 340 341 !According to getbsr and irdbsr, build _BSR file name, referred as fil_bsreso 342 !A default is available if getbsr is 0 343 stringfile='_BSR' ; stringvar='bsreso' 344 call mkfilename(filnam,fil_bsreso,dtset%getbsreso,idtset,dtset%irdbsreso,jdtset_,ndtset,stringfile,stringvar,will_read) 345 if(will_read==0) fil_bsreso=trim(filnam_ds(3))//'_BSR' 346 347 !According to getbsc and irdbsc, build _BSC file name, referred as fil_bscoup 348 !A default is available if getbsc is 0 349 stringfile='_BSC' ; stringvar='bscoup' 350 call mkfilename(filnam,fil_bscoup,dtset%getbscoup,idtset,dtset%irdbscoup,jdtset_,ndtset,stringfile,stringvar,will_read) 351 if(will_read==0)fil_bscoup=trim(filnam_ds(3))//'_BSC' 352 353 !According to getwfkfine and irdwfkfine, build _WFK file name, referred as filwfkfine 354 !A default is avaible if getwfkfine is 0 355 stringfile='_WFK' ; stringvar='wfkfine' 356 call mkfilename(filnam,filwfkfine,dtset%getwfkfine,idtset,dtset%irdwfkfine,jdtset_,& 357 & ndtset,stringfile,stringvar,will_read) 358 if(will_read==0)filwfkfine=trim(filnam_ds(3))//'_WFK' 359 360 dtfil%ireadden =ireadden 361 dtfil%ireadkden =ireadkden 362 dtfil%ireadwf =ireadwf 363 dtfil%filnam_ds(1:5)=filnam_ds(1:5) 364 365 dtfil%fnameabi_bsham_reso=fil_bsreso 366 dtfil%fnameabi_bsham_coup=fil_bscoup 367 dtfil%fnameabi_bseig=filbseig 368 dtfil%fnameabi_haydock=filhaydock 369 dtfil%fnameabi_sus =filsus 370 dtfil%fnameabi_qps =filqps 371 dtfil%fnameabi_scr =filscr 372 dtfil%filddbsin =filddbsin 373 dtfil%fildensin =fildensin 374 dtfil%fildens1in =fildens1in 375 dtfil%filkdensin =filkdensin 376 dtfil%filpawdensin =filpawdensin 377 dtfil%fnameabi_wfkfine = filwfkfine 378 dtfil%filstat =filstat 379 dtfil%fnamewffk =fnamewffk 380 dtfil%fnamewffq =fnamewffq 381 dtfil%fnamewffddk =fnamewffddk 382 dtfil%fnamewffdelfd =fnamewffdelfd 383 dtfil%fnamewffdkdk =fnamewffdkdk 384 dtfil%fnamewffdkde =fnamewffdkde 385 dtfil%fnamewff1 =fnamewff1 386 387 dtfil%fnameabi_hes=trim(dtfil%filnam_ds(3))//'_HES' 388 dtfil%fnameabi_phfrq=trim(dtfil%filnam_ds(3))//'_PHFRQ' 389 dtfil%fnameabi_phvec=trim(dtfil%filnam_ds(3))//'_PHVEC' 390 391 !------------------------------------------------------------------------------------------- 392 !Build name of files from dtfil%filnam_ds(4) 393 394 dtfil%fnameabo_ddb=trim(dtfil%filnam_ds(4))//'_DDB' 395 dtfil%fnameabo_den=trim(dtfil%filnam_ds(4))//'_DEN' 396 dtfil%fnameabo_dos=trim(dtfil%filnam_ds(4))//'_DOS' 397 dtfil%fnameabo_eelf=trim(dtfil%filnam_ds(4))//'_EELF' 398 dtfil%fnameabo_eig=trim(dtfil%filnam_ds(4))//'_EIG' 399 dtfil%fnameabo_eigi2d=trim(dtfil%filnam_ds(4))//'_EIGI2D' 400 dtfil%fnameabo_eigr2d=trim(dtfil%filnam_ds(4))//'_EIGR2D' 401 dtfil%fnameabo_em1=trim(dtfil%filnam_ds(4))//'_EM1' 402 dtfil%fnameabo_em1_lf=trim(dtfil%filnam_ds(4))//'_EM1_LF' 403 dtfil%fnameabo_em1_nlf=trim(dtfil%filnam_ds(4))//'_EM1_NLF' 404 dtfil%fnameabo_fan=trim(dtfil%filnam_ds(4))//'_FAN' 405 dtfil%fnameabo_gkk=trim(dtfil%filnam_ds(4))//'_GKK' 406 dtfil%fnameabo_gw=trim(dtfil%filnam_ds(4))//'_GW' ! TODO change name 407 dtfil%fnameabo_gwdiag=trim(dtfil%filnam_ds(4))//'_GWDIAG' 408 dtfil%fnameabo_kss=trim(dtfil%filnam_ds(4))//'_KSS' 409 dtfil%fnameabo_moldyn=trim(dtfil%filnam_ds(4))//'_MOLDYN' 410 dtfil%fnameabo_pot=trim(dtfil%filnam_ds(4))//'_POT' 411 dtfil%fnameabo_qps=trim(dtfil%filnam_ds(4))//'_QPS' 412 dtfil%fnameabo_qp_den=trim(dtfil%filnam_ds(4))//'_QP_DEN' 413 dtfil%fnameabo_qp_pawden=trim(dtfil%filnam_ds(4))//'_QP_PAWDEN' 414 dtfil%fnameabo_qp_dos=trim(dtfil%filnam_ds(4))//'_QP_DOS' 415 dtfil%fnameabo_qp_eig=trim(dtfil%filnam_ds(4))//'_QP_DB.nc' ! TODO change name 416 dtfil%fnameabo_rpa=trim(dtfil%filnam_ds(4))//'_RPA' 417 dtfil%fnameabo_scr=trim(dtfil%filnam_ds(4))//'_SCR' 418 dtfil%fnameabo_sgm=trim(dtfil%filnam_ds(4))//'_SGM' 419 dtfil%fnameabo_sgr=trim(dtfil%filnam_ds(4))//'_SGR' 420 dtfil%fnameabo_sig=trim(dtfil%filnam_ds(4))//'_SIG' 421 dtfil%fnameabo_spcur=trim(dtfil%filnam_ds(4))//'_SPCUR' 422 dtfil%fnameabo_sus=trim(dtfil%filnam_ds(4))//'_SUS' 423 dtfil%fnameabo_vha=trim(dtfil%filnam_ds(4))//'_VHA' 424 dtfil%fnameabo_vpsp=trim(dtfil%filnam_ds(4))//'_VPSP' 425 dtfil%fnameabo_vso=trim(dtfil%filnam_ds(4))//'_VSO' 426 dtfil%fnameabo_vxc=trim(dtfil%filnam_ds(4))//'_VXC' 427 dtfil%fnameabo_wan=trim(dtfil%filnam_ds(4))//'_WAN' 428 dtfil%fnameabo_wfk=trim(dtfil%filnam_ds(4))//'_WFK' 429 dtfil%fnameabo_wfq=trim(dtfil%filnam_ds(4))//'_WFQ' 430 dtfil%fnameabo_w90=trim(dtfil%filnam_ds(4))//'_w90' 431 dtfil%fnameabo_1wf=trim(dtfil%filnam_ds(4))//'_1WF' 432 dtfil%fnameabo_nlcc_derivs=trim(dtfil%filnam_ds(4))//'_nlcc_derivs_' 433 dtfil%fnameabo_pspdata=trim(dtfil%filnam_ds(4))//'_pspdata_' 434 435 !------------------------------------------------------------------------------------------- 436 !Build name of files from dtfil%filnam_ds(5) 437 438 dtfil%fnametmp_eig=trim(dtfil%filnam_ds(5))//'_EIG' 439 dtfil%fnametmp_1wf1_eig=trim(dtfil%filnam_ds(5))//'_1WF1_EIG' ! This appendix should be changed ! 440 dtfil%fnametmp_kgs=trim(dtfil%filnam_ds(5))//'_KGS' 441 dtfil%fnametmp_sustr=trim(dtfil%filnam_ds(5))//'_SUSTR' 442 dtfil%fnametmp_tdexcit=trim(dtfil%filnam_ds(5))//'_TDEXCIT' 443 dtfil%fnametmp_tdwf=trim(dtfil%filnam_ds(5))//'_TDWF' 444 dtfil%fnametmp_cg=trim(dtfil%filnam_ds(5))//'_cg' 445 dtfil%fnametmp_cprj=trim(dtfil%filnam_ds(5))//'_cprj' 446 447 !'_WF1' -> dtfil%unwft1 448 !'_WF2' -> dtfil%unwft2 449 !'_KG' -> dtfil%unkg 450 !'_DUM' -> tmp_unit (real dummy name) 451 !'_YLM' -> dtfil%unylm 452 !'_PAW' -> dtfil%unpaw 453 454 tmpfil(1)=trim(dtfil%filnam_ds(5))//'_WF1' ! tmpfil(1) 455 tmpfil(2)=trim(dtfil%filnam_ds(5))//'_WF2' ! tmpfil(2) 456 457 tmpfil(3)=trim(dtfil%filnam_ds(5))//'_KG' ! tmpfil(3) 458 tmpfil(4)=trim(dtfil%filnam_ds(5))//'_DUM' ! tmpfil(4) 459 tmpfil(5)=' ' ! to avoid Valgrind complain. 460 tmpfil(6)=trim(dtfil%filnam_ds(5))//'_YLM' ! tmpfil(6) 461 tmpfil(7)=trim(dtfil%filnam_ds(5))//'_PAW' ! tmpfil(7) 462 463 if(xmpi_paral==1)then ! parallel case : the index of the processor must be appended 464 call int2char4(mpi_enreg%me,tag) 465 ABI_CHECK((tag(1:1)/='#'),'Bug: string length too short!') 466 ixx=1 467 if (xmpi_mpiio == 1 .and. dtset%iomode == IO_MODE_MPI ) ixx=3 468 do ii=ixx,7 469 tmpfil(ii)=trim(tmpfil(ii))//'_P-'//trim(tag) 470 end do 471 end if 472 473 dtfil%fnametmp_wf1=trim(tmpfil(1)) 474 dtfil%fnametmp_wf2=trim(tmpfil(2)) 475 476 dtfil%fnametmp_kg =trim(tmpfil(3)) 477 dtfil%fnametmp_dum=trim(tmpfil(4)) 478 dtfil%fnametmp_ylm=trim(tmpfil(6)) 479 dtfil%fnametmp_paw=trim(tmpfil(7)) 480 481 !Create names for the temporary files based on dtfil%filnam_ds(5) 482 !by appending adequate string. 483 !'_1WF1' -> dtfil%unwft1 484 !'_1WF2' -> dtfil%unwft2 485 !'_KG' -> dtfil%unkg 486 !'_KGQ' -> dtfil%unkgq (not used for the time being) 487 !'_KG1' -> dtfil%unkg1 488 !'_DUM' -> tmp_unit (real dummy name) 489 !'_WFGS' -> dtfil%unwftgs 490 !'_WFKQ' -> dtfil%unwftkq 491 !'_YLM' -> dtfil%unylm 492 !'_YLM1' -> dtfil%unylm1 493 !'_PAW' -> dtfil%unpaw 494 !'_PAW1' -> dtfil%unpaw1 495 !'_PAWQ' -> dtfil%unpawq 496 tmpfil(1) =trim(dtfil%filnam_ds(5))//'_1WF1' 497 tmpfil(2) =trim(dtfil%filnam_ds(5))//'_1WF2' 498 tmpfil(3) =trim(dtfil%filnam_ds(5))//'_KG' 499 tmpfil(4) =trim(dtfil%filnam_ds(5))//'_KGQ' 500 tmpfil(5) =trim(dtfil%filnam_ds(5))//'_KG1' 501 tmpfil(6) =trim(dtfil%filnam_ds(5))//'_DUM' 502 tmpfil(7) =trim(dtfil%filnam_ds(5))//'_WFGS' 503 tmpfil(8) =trim(dtfil%filnam_ds(5))//'_WFKQ' 504 tmpfil(9) =' ' ! for Valgrind, to avoid uninitialized 505 tmpfil(10)=trim(dtfil%filnam_ds(5))//'_YLM' 506 tmpfil(11)=trim(dtfil%filnam_ds(5))//'_YLM1' 507 tmpfil(12)=trim(dtfil%filnam_ds(5))//'_PAW' 508 tmpfil(13)=trim(dtfil%filnam_ds(5))//'_PAW1' 509 tmpfil(14)=trim(dtfil%filnam_ds(5))//'_PAWQ' 510 511 if(xmpi_paral==1) then 512 idtmpfil(:)=0 513 do ii=1,14 514 idtmpfil(ii)=ii 515 end do 516 if (xmpi_mpiio==1.and.dtset%iomode==IO_MODE_MPI)then 517 idtmpfil(1)=0 !_1wf1 518 idtmpfil(2)=0 ! s1wf2 519 idtmpfil(7)=0 ! WFGS 520 idtmpfil(8)=0 ! WFKQ 521 end if 522 call int2char4(mpi_enreg%me,tag) 523 ABI_CHECK((tag(1:1)/='#'),'Bug: string length too short!') 524 do ii=1,14 525 if(idtmpfil(ii) /= 0) tmpfil(ii)=trim(tmpfil(ii))//'_P-'//trim(tag) 526 end do 527 end if 528 529 dtfil%fnametmp_1wf1=trim(tmpfil(1)) 530 dtfil%fnametmp_1wf2=trim(tmpfil(2)) 531 dtfil%fnametmp_kg =trim(tmpfil(3)) 532 dtfil%fnametmp_kgq =trim(tmpfil(4)) 533 dtfil%fnametmp_kg1 =trim(tmpfil(5)) 534 dtfil%fnametmp_dum =trim(tmpfil(6)) 535 dtfil%fnametmp_wfgs=trim(tmpfil(7)) 536 dtfil%fnametmp_wfkq=trim(tmpfil(8)) 537 dtfil%fnametmp_ylm =trim(tmpfil(10)) 538 dtfil%fnametmp_ylm1=trim(tmpfil(11)) 539 dtfil%fnametmp_paw =trim(tmpfil(12)) 540 dtfil%fnametmp_paw1=trim(tmpfil(13)) 541 dtfil%fnametmp_pawq=trim(tmpfil(14)) 542 543 !Prepare the name of the _FFT file 544 filfft=trim(dtfil%filnam_ds(5))//'_FFT' 545 !There is a definite problem in the treatment of // by CPP ... 546 if(xmpi_paral==1 .or. mpi_enreg%paral_kgb==1)then 547 call int2char4(mpi_enreg%me,tag) 548 ABI_CHECK((tag(1:1)/='#'),'Bug: string length too short!') 549 filfft=trim(filfft)//'_P-'//trim(tag) 550 end if 551 dtfil%fnametmp_fft=filfft 552 553 !These keywords are only used in algorithms using images of the cell 554 if (iimage==0) then 555 dtfil%getwfk_from_image =0 556 dtfil%getden_from_image =0 557 dtfil%getpawden_from_image=0 558 end if 559 560 DBG_EXIT("COLL") 561 562 end subroutine dtfil_init