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