TABLE OF CONTENTS
ABINIT/m_paw_io [ Modules ]
NAME
m_paw_io
FUNCTION
PAW I/O related operations
COPYRIGHT
Copyright (C) 2012-2018 ABINIT group (MT, TR) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
NOTES
FOR DEVELOPPERS: in order to preserve the portability of libPAW library, please consult ~abinit/src/??_libpaw/libpaw-coding-rules.txt
SOURCE
21 #include "libpaw.h" 22 23 module m_paw_io 24 25 USE_DEFS 26 USE_MSG_HANDLING 27 USE_MEMORY_PROFILING 28 29 implicit none 30 31 private 32 33 public :: pawio_print_ij
m_paw_io/pawio_print_ij [ Functions ]
[ Top ] [ m_paw_io ] [ Functions ]
NAME
pawio_print_ij
FUNCTION
Print ij_ square matrixes in a "suitable" format. Data are "energy-like" in Hartree units. Devoted to the printing of rhoij, Dij -like PAW matrixes.
INPUTS
a_ij(cplex*adim)= input square matrix
asym_ij(cplex*adim)= -OPTIONAL ARGUMENT-
When present, A(j,i) is deduced from asym_ij
instead of a_ij
adim= dimension of array a_ij:
adim=ndim*(ndim+1)/2 if opt_pack= 0
adim=number of non-zero values of a_ij if opt_pack=+1
cplex=1 if a_ij is real, 2 if it is complex
[mode_paral]= --optional argument, default='COLL'--
'COLL' if all procs are calling the routine with the same message to be written once only.
'PERS' if the procs are calling the routine with different messages each to be written,
or if one proc is calling the routine
ndim= dimension of input square matrix
opt_l= if <0 all parts of a_ij are printed
if >=0 only parts of a_ij corresponding to li=lj=opt_l are printed
opt_l_index(ndim)= array giving l quantum number for each 1<=ilmn<=ndim
not used if opt_l<0
opt_pack= 0 if a_ij is given as A(j(j-1)/2+i), i<=j
+1 if a_ij is given as A(j(j-1)/2+i) and is in "packed storage"
(i.e. only non-zero values are stored)
opt_prtvol= >=0 if up to 12 components of _ij matrix have to be printed
<0 if all components of ij_ matrix have to be printed
opt_sym= -OPTIONAL ARGUMENT- (default if not present: opt_sym=2)
Define the symmetry of a_ij matrix:
opt_sym=1 : A(j,i)= A(i,j)
opt_sym=2 : A(j,i)= Conjg[A(i,j)]
opt_sym=3 : A(j,i)=-A(i,j)
opt_sym=4 : A(j,i)=-Conjg[A(i,j)]
When asym_ij argument is present, A[i,j] is taken from it.
pack2ij(adim)= gives the (i,j) index of of packed value of rhoij
used only if opt_packed=+1
test_value= (real number) if positive, print a warning when the
magnitude of a_ij is greater than opt_test
No test when test_value<0
unit=the unit number for output
Ha_or_eV= 1: output in hartrees, 2: output in eV
PARENTS
m_paw_ij,m_paw_slater,m_pawdij,m_pawrhoij,pawmkrhoij,pawprt setrhoijpbe0,wfd_pawrhoij
CHILDREN
wrtout
SOURCE
94 subroutine pawio_print_ij(unit,a_ij,adim,cplex,ndim,opt_l,opt_l_index, & 95 & opt_pack,opt_prtvol,pack2ij,test_value,Ha_or_eV, & 96 & mode_paral,opt_sym,asym_ij) ! Optional arguments 97 98 99 !This section has been created automatically by the script Abilint (TD). 100 !Do not modify the following lines by hand. 101 #undef ABI_FUNC 102 #define ABI_FUNC 'pawio_print_ij' 103 !End of the abilint section 104 105 implicit none 106 107 !Arguments --------------------------------------------- 108 !scalars 109 integer,intent(in) :: adim,cplex,ndim,opt_l,opt_pack,opt_prtvol,unit,Ha_or_eV 110 integer,intent(in),optional :: opt_sym 111 character(len=*),optional,intent(in) :: mode_paral 112 real(dp),intent(in) :: test_value 113 !arrays 114 integer,intent(in) :: opt_l_index(ndim*min(1+opt_l,1)),pack2ij(adim*opt_pack) 115 real(dp),intent(in) :: a_ij(cplex*adim) 116 real(dp),intent(in),optional :: asym_ij(cplex*adim) 117 118 !Local variables --------------------------------------- 119 ! Adjust format bellow according to maxprt 120 !scalars 121 integer,parameter :: maxprt_default=12 122 integer :: dplex,ilmn,ilmn1,j0lmn,jlmn,jlmn1,klmn,klmn1,klmn2,maxprt,nhigh 123 integer :: nmin,optsym 124 real(dp) :: testval 125 logical :: use_asym 126 character(len=4) :: mode_paral_ 127 character(len=500) :: msg='' 128 !arrays 129 real(dp),parameter :: fact_re(4)=(/one,one,-one,-one/),fact_im(4)=(/one,-one,-one,one/) 130 real(dp) :: tabmax(cplex),tabmin(cplex) 131 real(dp),allocatable :: b_ij(:),bsym_ij(:),prtab(:,:,:) 132 133 ! ************************************************************************* 134 135 10 format(100(1x,f9.5)) 136 11 format(12(1x,f9.5),a) !Change this format according to variable "maxprt" 137 138 139 !DEBUG 140 !write(std_out,*)' pawio_print_ij : enter ' 141 !ENDDEBUG 142 143 !Optional arguments 144 mode_paral_='COLL';if (present(mode_paral)) mode_paral_=mode_paral 145 use_asym=present(asym_ij) 146 if (present(opt_sym)) then 147 optsym=opt_sym 148 else 149 optsym=2 150 end if 151 152 !Define size of square matrix 153 if (opt_prtvol>=0) then 154 maxprt=maxprt_default 155 else 156 maxprt=ndim 157 end if 158 nmin=min(ndim,maxprt) 159 160 if (opt_l>=0) nmin=count(opt_l_index(:)==opt_l) 161 LIBPAW_ALLOCATE(prtab,(cplex,nmin,nmin)) 162 dplex=cplex-1 163 164 !Eventually unpack input matrix(es) 165 LIBPAW_ALLOCATE(b_ij,(cplex*ndim*(ndim+1)/2)) 166 if (opt_pack==0) then 167 b_ij=a_ij 168 else if (opt_pack==1) then 169 b_ij=zero 170 do klmn=1,adim 171 klmn1=cplex*klmn-dplex 172 klmn2=cplex*pack2ij(klmn)-dplex 173 b_ij(klmn2:klmn2+dplex)=a_ij(klmn1:klmn1+dplex) 174 end do 175 end if 176 if (opt_prtvol<0.and.opt_l<0) then 177 if (cplex==1) then 178 tabmax(1)=maxval(abs(b_ij)) 179 tabmin(1)=minval(abs(b_ij)) 180 else 181 tabmax(1:2)=zero;tabmin(1:2)=1.d20 182 do klmn=1,size(b_ij)/cplex 183 klmn2=2*klmn 184 tabmax(1)=max(tabmax(1),b_ij(klmn2-1)) 185 tabmin(1)=min(tabmin(1),b_ij(klmn2-1)) 186 tabmax(2)=max(tabmax(2),b_ij(klmn2 )) 187 tabmin(2)=min(tabmin(2),b_ij(klmn2 )) 188 end do 189 end if 190 end if 191 if (use_asym) then 192 LIBPAW_ALLOCATE(bsym_ij,(cplex*ndim*(ndim+1)/2)) 193 if (opt_pack==0) then 194 bsym_ij=asym_ij 195 else if (opt_pack==1) then 196 bsym_ij=zero 197 do klmn=1,adim 198 klmn1=cplex*klmn-dplex 199 klmn2=cplex*pack2ij(klmn)-dplex 200 bsym_ij(klmn2:klmn2+dplex)=asym_ij(klmn1:klmn1+dplex) 201 end do 202 end if 203 if (opt_prtvol<0.and.opt_l<0) then 204 if (cplex==1) then 205 tabmax(1)=max(tabmax(1),maxval(abs(bsym_ij))) 206 tabmin(1)=min(tabmin(1),minval(abs(bsym_ij))) 207 else 208 do klmn=1,ndim 209 klmn2=2*klmn 210 tabmax(1)=max(tabmax(1),bsym_ij(klmn2-1)) 211 tabmin(1)=min(tabmin(1),bsym_ij(klmn2-1)) 212 tabmax(2)=max(tabmax(2),bsym_ij(klmn2 )) 213 tabmin(2)=min(tabmin(2),bsym_ij(klmn2 )) 214 end do 215 end if 216 end if 217 end if 218 219 !Transfer triangular matrix to rectangular one 220 jlmn1=0 221 do jlmn=1,ndim 222 if (opt_l<0) then 223 jlmn1=jlmn;if (jlmn1>nmin) cycle 224 else if (opt_l_index(jlmn)==opt_l) then 225 jlmn1=jlmn1+1 226 else 227 cycle 228 end if 229 ilmn1=0;j0lmn=jlmn*(jlmn-1)/2 230 do ilmn=1,jlmn 231 if (opt_l<0) then 232 ilmn1=ilmn 233 else if (opt_l_index(ilmn)==opt_l) then 234 ilmn1=ilmn1+1 235 else 236 cycle 237 end if 238 klmn=j0lmn+ilmn 239 if (cplex==1) then 240 prtab(1,ilmn1,jlmn1)=b_ij(klmn) 241 if (use_asym) then 242 prtab(1,jlmn1,ilmn1)=fact_re(optsym)*bsym_ij(klmn) 243 else 244 prtab(1,jlmn1,ilmn1)=fact_re(optsym)*b_ij(klmn) 245 end if 246 else 247 klmn=2*klmn 248 prtab(1:2,ilmn1,jlmn1)=b_ij(klmn-1:klmn) 249 if (use_asym) then 250 prtab(1,jlmn1,ilmn1)=fact_re(optsym)*bsym_ij(klmn-1) 251 prtab(2,jlmn1,ilmn1)=fact_im(optsym)*bsym_ij(klmn ) 252 else 253 prtab(1,jlmn1,ilmn1)=fact_re(optsym)*b_ij(klmn-1) 254 prtab(2,jlmn1,ilmn1)=fact_im(optsym)*b_ij(klmn ) 255 end if 256 end if 257 end do 258 end do 259 LIBPAW_DEALLOCATE(b_ij) 260 261 if (use_asym) then 262 LIBPAW_DEALLOCATE(bsym_ij) 263 end if 264 265 if (Ha_or_eV==2) then 266 prtab=prtab*Ha_eV 267 if (opt_prtvol<0.and.opt_l<0) then 268 tabmax=tabmax*Ha_eV 269 tabmin=tabmin*Ha_eV 270 end if 271 end if 272 273 if (cplex==2) then 274 write(msg,'(3x,a)') '=== REAL PART:' 275 call wrtout(unit,msg,mode_paral_) 276 end if 277 278 if (ndim<=maxprt.or.opt_l>=0) then 279 do ilmn=1,nmin 280 write(msg,fmt=10) prtab(1,1:nmin,ilmn) 281 call wrtout(unit,msg,mode_paral_) 282 end do 283 else 284 do ilmn=1,nmin 285 write(msg,fmt=11) prtab(1,1:nmin,ilmn),' ...' 286 call wrtout(unit,msg,mode_paral_) 287 end do 288 write(msg,'(3x,a,i2,a)') '... only ',maxprt,' components have been written...' 289 call wrtout(unit,msg,mode_paral_) 290 end if 291 if (opt_prtvol<0.and.opt_l<0) then 292 write(msg,'(3x,2(a,es9.2))') 'max. value= ',tabmax(1),', min. value= ',tabmin(1) 293 call wrtout(unit,msg,mode_paral_) 294 end if 295 296 if (cplex==2) then 297 write(msg,'(3x,a)') '=== IMAGINARY PART:' 298 call wrtout(unit,msg,mode_paral_) 299 if (ndim<=maxprt.or.opt_l>=0) then 300 do ilmn=1,nmin 301 write(msg,fmt=10) prtab(2,1:nmin,ilmn) 302 call wrtout(unit,msg,mode_paral_) 303 end do 304 else 305 do ilmn=1,nmin 306 write(msg,fmt=11) prtab(2,1:nmin,ilmn),' ...' 307 call wrtout(unit,msg,mode_paral_) 308 end do 309 write(msg,'(3x,a,i2,a)') '... only ',maxprt,' components have been written...' 310 call wrtout(unit,msg,mode_paral_) 311 end if 312 if (opt_prtvol<0.and.opt_l<0) then 313 write(msg,'(3x,2(a,es9.2))') 'max. value= ',tabmax(2),', min. value= ',tabmin(2) 314 call wrtout(unit,msg,mode_paral_) 315 end if 316 end if 317 318 if (test_value>zero) then 319 testval=test_value;if (Ha_or_eV==2) testval=testval*Ha_eV 320 nhigh=0;nhigh=count(abs(prtab(:,:,:))>=testval) 321 if (nhigh>0) then 322 write(msg,'(5a,i3,a,f6.1,7a)')& 323 & ' pawio_print_ij: WARNING -',ch10,& 324 & ' The matrix seems to have high value(s) !',ch10,& 325 & ' (',nhigh,' components have a value greater than ',testval,').',ch10,& 326 & ' It can cause instabilities during SCF convergence.',ch10,& 327 & ' Action: you should check your atomic dataset (psp file)',ch10,& 328 & ' and look for "high" projector functions...' 329 call wrtout(unit,msg,mode_paral_) 330 end if 331 end if 332 333 LIBPAW_DEALLOCATE(prtab) 334 335 !DEBUG 336 !write(std_out,*)' pawio_print_ij : exit ' 337 !ENDDEBUG 338 339 end subroutine pawio_print_ij