TABLE OF CONTENTS
ABINIT/ujdet [ Programs ]
NAME
ujdet
FUNCTION
Main routine. Determines U from inputfile ujdet.in (reduced inputfile containing mainly the atomic occupancies, atomic positions, and potential shifts.
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (DJA) 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 .
NOTES
INPUTS
OUTPUT
PARENTS
CHILDREN
abi_io_redirect,abimem_init,abinit_doctor,destroy_mpi_enreg,herald initmpi_seq,intagm,isfile,parsefile,pawuj_det,pawuj_free,pawuj_ini wrtout,xmpi_end,xmpi_init
SOURCE
29 #if defined HAVE_CONFIG_H 30 #include "config.h" 31 #endif 32 33 #include "abi_common.h" 34 35 36 program ujdet 37 38 use defs_basis 39 use defs_abitypes 40 use defs_datatypes 41 use m_xmpi 42 use m_abicore 43 use m_build_info 44 use m_errors 45 46 use m_paw_uj, only : pawuj_ini,pawuj_free,pawuj_det 47 48 use m_specialmsg, only : specialmsg_getcount, herald 49 use m_io_tools, only : open_file 50 use m_parser, only : intagm, parsefile 51 use m_mpinfo, only : destroy_mpi_enreg, initmpi_seq 52 use m_dtfil, only : isfile 53 54 !This section has been created automatically by the script Abilint (TD). 55 !Do not modify the following lines by hand. 56 #undef ABI_FUNC 57 #define ABI_FUNC 'ujdet' 58 !End of the abilint section 59 60 implicit none 61 62 !Local variables------------------------------- 63 !scalars 64 integer,parameter :: ndtpawuj=4,nwfchr=6,master=0 65 integer :: ii,jdtset,lenstr,marr,ndtset,tread,comm 66 integer :: nat,nproc,my_rank,nspden,macro_uj,pawujat,pawprtvol 67 logical :: iam_master 68 character(len=24) :: codename 69 character(len=30) :: token 70 character(len=fnlen) :: filnam(5) 71 character(len=strlen) :: string 72 character(len=500) :: message 73 type(MPI_type) :: mpi_enreg 74 real(dp) :: ures 75 !arrays 76 integer,allocatable :: idumar1(:),intarr(:),jdtset_(:) 77 real(dp),allocatable :: dpdumar(:),dprarr(:) 78 type(macro_uj_type),allocatable :: dtpawuj(:) 79 80 ! ********************************************************************* 81 82 !Change communicator for I/O (mandatory!) 83 call abi_io_redirect(new_io_comm=xmpi_world) 84 85 !Initialize MPI (not used but necessary) 86 call xmpi_init() 87 comm = xmpi_world 88 nproc = xmpi_comm_size(comm) 89 my_rank = xmpi_comm_rank(comm) 90 iam_master = (my_rank == master) 91 92 !Initialize memory profiling if it is activated 93 !if a full abimem.mocc report is desired, set the argument of abimem_init to "2" instead of "0" 94 !note that abimem.mocc files can easily be multiple GB in size so don't use this option normally 95 #ifdef HAVE_MEM_PROFILING 96 call abimem_init(0) 97 #endif 98 99 !Default for sequential use 100 call initmpi_seq(mpi_enreg) 101 102 !Signal MPI I/O compilation has been activated 103 #if defined HAVE_MPI_IO 104 if(xmpi_paral==0)then 105 write(message,'(3a)')& 106 & 'In order to use MPI_IO, you must compile with the MPI flag ',ch10,& 107 & 'Action : recompile your code with different CPP flags.' 108 MSG_ERROR(message) 109 end if 110 #endif 111 112 !---------------------------------------------------------------------- 113 114 !Check that old input file exists (18seqpar/iofn1.F90) 115 filnam(1)='ujdet.in' 116 call isfile(filnam(1),'old') 117 118 codename='UJDET'//repeat(' ',18) 119 filnam(2)='ujdet.out' 120 121 !Check that new output file does NOT exist (18seqpar/iofn1.F90) 122 if (iam_master) then 123 call isfile(filnam(2),'new') 124 if (open_file(filnam(2),message,unit=ab_out,form="formatted",status="new") /= 0) then 125 MSG_ERROR(message) 126 end if 127 rewind(unit=ab_out) 128 ! Print header 129 call herald(codename,abinit_version,ab_out) 130 end if 131 132 !Print header 133 call herald(codename,abinit_version,std_out) 134 135 !Read the file, stringify it and return the number of datasets. (main/abinit.F90) 136 call parsefile(filnam(1), lenstr, ndtset, string, comm) 137 138 !DEBUG 139 !write(std_out,*)'ujdet: trim(string):',trim(adjustl(string)) 140 !write(std_out,*)'ujdet: lenstr: ',lenstr 141 !write(std_out,*)'ujdet: ndtset: ',ndtset 142 !ENDDEBUG 143 144 ABI_DATATYPE_ALLOCATE(dtpawuj,(0:ndtpawuj)) 145 146 call pawuj_ini(dtpawuj,ndtset) 147 148 dtpawuj(1:ndtset)%ndtset=ndtset 149 150 marr=ndtset 151 ABI_ALLOCATE(idumar1,(marr)) 152 ABI_ALLOCATE(dpdumar,(marr)) 153 ABI_ALLOCATE(jdtset_,(ndtset)) 154 jdtset_=(/ ( ii,ii=1,ndtset )/) 155 156 !Read integers (main dimensions) 157 do jdtset=1,ndtset 158 ! Integer 159 token = 'pawprtvol' 160 call intagm(dpdumar,idumar1,jdtset,marr,1,string(1:lenstr),token,tread,'INT') 161 if(tread==1) dtpawuj(1:ndtset)%pawprtvol=idumar1(1) 162 163 token = 'nat' 164 call intagm(dpdumar,idumar1,jdtset,marr,1,string(1:lenstr),token,tread,'INT') 165 if(tread==1) dtpawuj(1:ndtset)%nat=idumar1(1) 166 167 token = 'nspden' 168 call intagm(dpdumar,idumar1,jdtset,marr,1,string(1:lenstr),token,tread,'INT') 169 if(tread==1) dtpawuj(1:ndtset)%nspden=idumar1(1) 170 171 token = 'macro_uj' 172 call intagm(dpdumar,idumar1,jdtset,marr,1,string(1:lenstr),token,tread,'INT') 173 if(tread==1) dtpawuj(1:ndtset)%macro_uj=idumar1(1) 174 175 token = 'pawujat' 176 call intagm(dpdumar,idumar1,jdtset,marr,1,string(1:lenstr),token,tread,'INT') 177 if(tread==1) dtpawuj(1:ndtset)%pawujat=idumar1(1) 178 179 token = 'dmatpuopt' 180 call intagm(dpdumar,idumar1,jdtset,marr,1,string(1:lenstr),token,tread,'INT') 181 if(tread==1) dtpawuj(1:ndtset)%dmatpuopt=idumar1(1) 182 183 ! DEBUG 184 ! write(std_out,*)'ujdet: read pawujat ; jdtset, idumar1(1)',jdtset, idumar1(1) 185 ! write(std_out,*)'ujdet: read pawujat ; dtpawuj(1:ndtset)%pawujat ', dtpawuj(1:ndtset)%pawujat 186 ! END DEBUG 187 188 token = 'pawujopt' 189 call intagm(dpdumar,idumar1,jdtset,marr,1,string(1:lenstr),token,tread,'INT') 190 if(tread==1) dtpawuj(1:ndtset)%option=idumar1(1) 191 192 ! Real 193 194 token = 'diemix' 195 call intagm(dpdumar,idumar1,jdtset,marr,1,string(1:lenstr),token,tread,'DPR') 196 if(tread==1) dtpawuj(1:ndtset)%diemix=idumar1(1) 197 198 token = 'mdist' 199 call intagm(dpdumar,idumar1,jdtset,marr,1,string(1:lenstr),token,tread,'DPR') 200 if(tread==1) dtpawuj(1:ndtset)%mdist=idumar1(1) 201 202 token = 'pawujga' 203 call intagm(dpdumar,idumar1,jdtset,marr,1,string(1:lenstr),token,tread,'DPR') 204 if(tread==1) dtpawuj(1:ndtset)%pawujga=dpdumar(1) 205 206 token = 'ph0phiint' 207 call intagm(dpdumar,idumar1,jdtset,marr,1,string(1:lenstr),token,tread,'DPR') 208 if(tread==1) dtpawuj(1:ndtset)%ph0phiint=dpdumar(1) 209 210 token = 'pawrad' 211 call intagm(dpdumar,idumar1,jdtset,marr,1,string(1:lenstr),token,tread,'LEN') 212 if(tread==1) dtpawuj(1:ndtset)%pawrad=dpdumar(1) 213 214 token = 'pawujrad' 215 call intagm(dpdumar,idumar1,jdtset,marr,1,string(1:lenstr),token,tread,'LEN') 216 if(tread==1) dtpawuj(1:ndtset)%pawujrad=dpdumar(1) 217 218 end do !ndtset 219 220 nat=dtpawuj(1)%nat 221 nspden=dtpawuj(1)%nspden 222 macro_uj=dtpawuj(1)%macro_uj 223 pawujat=dtpawuj(1)%pawujat 224 pawprtvol=dtpawuj(1)%pawprtvol 225 226 !Minimal tests 227 if (.not.all(dtpawuj(1:ndtset)%nat==nat).and.all(dtpawuj(1:ndtset)%pawujat==pawujat)& 228 & .and.all(dtpawuj(1:ndtset)%macro_uj==macro_uj)) then 229 write(message,'(a)')'Problems with nat, pawujat or nspden. Values for datasets differ.' 230 call wrtout(std_out,message,'COLL') 231 write(message,*) ' ujdet: nat ', all(dtpawuj(1:ndtset)%nat==nat) 232 call wrtout(std_out,message,'COLL') 233 write(message,*) ' ujdet: pawujat ', all(dtpawuj(1:ndtset)%pawujat==pawujat) 234 call wrtout(std_out,message,'COLL') 235 write(message,*) ' ujdet: macro_uj', all(dtpawuj(1:ndtset)%macro_uj==macro_uj) 236 call wrtout(std_out,message,'COLL') 237 end if 238 239 if (pawprtvol>1) then 240 write(message,fmt='(a,150i3)')' ujdet: dtpawuj(0:ndtset)%macro_uj ',dtpawuj(0:ndtset)%macro_uj 241 call wrtout(std_out,message,'COLL') 242 write(message,fmt='(a,150i3)')' ujdet: dtpawuj(0:ndtset)%pawujat ',dtpawuj(0:ndtset)%pawujat 243 call wrtout(std_out,message,'COLL') 244 write(message,fmt='(a,150i3)')' ujdet: dtpawuj(0:ndtset)%nspden ',dtpawuj(0:ndtset)%nspden 245 call wrtout(std_out,message,'COLL') 246 write(message,fmt='(a,150i3)')' ujdet: nat*nspden ',nat*nspden 247 call wrtout(std_out,message,'COLL') 248 end if 249 250 !Read arrays 251 252 marr=maxval((/nspden*nat*nat, 3*3 ,nat*3/)) 253 ABI_ALLOCATE(intarr,(marr)) 254 ABI_ALLOCATE(dprarr,(marr)) 255 256 do jdtset=0,ndtset 257 ! DEBUG 258 ! write(std_out,*)'ujdet 2a, jdtset ',jdtset 259 ! END DEBUG 260 ABI_ALLOCATE(dtpawuj(jdtset)%vsh,(nspden,nat)) 261 ABI_ALLOCATE(dtpawuj(jdtset)%occ,(nspden,nat)) 262 ABI_ALLOCATE(dtpawuj(jdtset)%xred,(3,nat)) 263 264 dtpawuj(jdtset)%iuj=jdtset 265 dtpawuj(jdtset)%vsh=zero 266 dtpawuj(jdtset)%occ=zero 267 dtpawuj(jdtset)%xred=zero 268 end do 269 270 do jdtset=1,ndtset 271 dprarr=-1_dp 272 intarr=-1 273 274 ! Integer arrays 275 276 ! scdim, wfchr and rprimd allocated in pawuj_ini 277 token = 'scdim' 278 call intagm(dprarr,intarr,jdtset,marr,3,string(1:lenstr),token,tread,'INT') 279 if(tread==1) dtpawuj(jdtset)%scdim(1:3)=intarr(1:3) 280 281 token = 'wfchr' 282 call intagm(dprarr,intarr,jdtset,marr,nwfchr,string(1:lenstr),token,tread,'DPR') 283 if(tread==1) dtpawuj(jdtset)%wfchr(1:nwfchr)=dprarr(1:nwfchr) 284 285 write(std_out,*)'ujdet: wfchr ',dtpawuj(jdtset)%wfchr 286 287 token = 'wfchr' 288 call intagm(dprarr,intarr,jdtset,marr,nwfchr,string(1:lenstr),token,tread,'DPR') 289 if(tread==1) dtpawuj(jdtset)%wfchr(1:nwfchr)=dprarr(1:nwfchr) 290 291 write(std_out,*)'ujdet: wfchr ',dtpawuj(jdtset)%wfchr 292 293 token = 'rprimd' 294 call intagm(dprarr,intarr,jdtset,marr,3*3,string(1:lenstr),token,tread,'DPR') 295 if(tread==1) dtpawuj(jdtset)%rprimd(1:3,1:3)=reshape(dprarr(1:3*3),(/ 3,3/)) 296 297 token = 'occ' 298 call intagm(dprarr,intarr,jdtset,marr,nspden*nat,string(1:lenstr),token,tread,'DPR') 299 if(tread==1) dtpawuj(jdtset)%occ(1:nspden,1:nat)=reshape(dprarr(1:nspden*nat),(/ nspden,nat/)) 300 301 token = 'vsh' 302 call intagm(dprarr,intarr,jdtset,marr,nat*nspden,string(1:lenstr),token,tread,'ENE') 303 if(tread==1) dtpawuj(jdtset)%vsh(1:nspden,1:nat)=reshape(dprarr(1:nspden*nat),(/ nspden,nat/)) 304 305 token = 'xred' 306 call intagm(dprarr,intarr,jdtset,marr,nat*3,string(1:lenstr),token,tread,'DPR') 307 if(tread==1) dtpawuj(jdtset)%xred(1:3,1:nat)=reshape(dprarr(1:3*nat),(/ 3, nat/)) 308 309 ! DEBUG 310 ! write(std_out,*)'dtpawuj(',jdtset,')%occ ',dtpawuj(jdtset)%occ,ch10 311 ! write(std_out,*)'dtpawuj(',jdtset,')%rprimd ',dtpawuj(jdtset)%rprimd,ch10 312 ! write(std_out,*)'dtpawuj(',jdtset,')%vsh ',dtpawuj(jdtset)%vsh,ch10 313 ! write(std_out,*)'dtpawuj(',jdtset,')%xred ',dtpawuj(jdtset)%xred,ch10,ch10 314 ! END DEBUG 315 316 end do ! jdtset 317 318 !Not yet really parallelized ... 319 if (iam_master) then 320 call pawuj_det(dtpawuj,ndtset, filnam(2)//"_UJDET.nc",ures) 321 end if 322 323 !Close files 324 if (iam_master) then 325 close(ab_out) 326 end if 327 328 !Deallocations 329 do jdtset=0,ndtset 330 call pawuj_free(dtpawuj(jdtset)) 331 end do 332 ABI_DATATYPE_DEALLOCATE(dtpawuj) 333 334 ABI_DEALLOCATE(intarr) 335 ABI_DEALLOCATE(dprarr) 336 ABI_DEALLOCATE(idumar1) 337 ABI_DEALLOCATE(dpdumar) 338 ABI_DEALLOCATE(jdtset_) 339 340 call destroy_mpi_enreg(mpi_enreg) 341 342 !Write information on file about the memory before ending mpi module, if memory profiling is enabled 343 call abinit_doctor("__ujdet") 344 345 call xmpi_end() 346 347 end program ujdet