TABLE OF CONTENTS
ABINIT/prtrhomxmn [ Functions ]
NAME
prtrhomxmn
FUNCTION
If option==1, compute the maximum and minimum of the density (and spin-polarization if nspden==2), and print it. If option==2, also compute and print the second maximum or minimum
COPYRIGHT
Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
INPUTS
iout=unit for output file mpi_enreg=information about MPI parallelization nfft=(effective) number of FFT grid points (for this processor) ngfft(18)=contain all needed information about 3D FFT, see ~abinit/doc/variables/vargs.htm#ngfft nspden=number of spin-density components option, see above optrhor=option for rhor (If optrhor==0, rhor is expected to be the electron density) (If optrhor==1, rhor is expected to be the kinetic energy density (taur)) (If optrhor==2, rhor is expected to be the gradient of the electron density (grhor)) (If optrhor==3, rhor is expected to be the laplacian of the electron density (lrhor)) (If optrhor==4, rhor is expected to be the ELF (elfr)) rhor(nfft,nspden)=electron density (electrons/bohr^3)
OUTPUT
SIDE EFFECTS
NOTES
The tolerance tol12 aims at giving a machine-independent ordering. (this trick is used in bonds.f, listkk.f, prtrhomxmn.f and rsiaf9.f)
PARENTS
afterscfloop,bethe_salpeter,clnup1,mkrho,screening,sigma,vtorho
CHILDREN
wrtout,xmpi_sum
SOURCE
47 #if defined HAVE_CONFIG_H 48 #include "config.h" 49 #endif 50 51 #include "abi_common.h" 52 53 subroutine prtrhomxmn(iout,mpi_enreg,nfft,ngfft,nspden,option,rhor,optrhor,ucvol) 54 55 use defs_basis 56 use defs_abitypes 57 use m_profiling_abi 58 use m_xmpi 59 use m_errors 60 61 !This section has been created automatically by the script Abilint (TD). 62 !Do not modify the following lines by hand. 63 #undef ABI_FUNC 64 #define ABI_FUNC 'prtrhomxmn' 65 use interfaces_14_hidewrite 66 !End of the abilint section 67 68 implicit none 69 70 !Arguments ------------------------------------ 71 !scalars 72 integer,intent(in) :: iout,nfft,nspden,option 73 integer,intent(in),optional :: optrhor 74 real(dp),intent(in),optional :: ucvol 75 type(MPI_type),intent(in) :: mpi_enreg 76 !arrays 77 integer,intent(in) :: ngfft(18) 78 real(dp),intent(in) :: rhor(nfft,nspden) 79 80 !Local variables------------------------------- 81 !scalars 82 integer :: i1,i2,i3,ierr,ifft,ii,iisign,iitems,index1,ioptrhor 83 integer :: index2,indsign,iproc,istart,me,n1,n2,n3,nitems 84 integer :: nfft_,nfftot,nproc,spaceComm 85 real(dp) :: temp,value1,value2 86 character(len=500) :: message,txt1_in_mssg,txt2_in_mssg,txt3_in_mssg 87 logical :: reduce=.false. 88 !arrays 89 integer,allocatable :: iindex(:,:,:),index_fft(:,:,:,:) 90 real(dp) :: rhomn1(4),rhomn2(4),rhomx1(4),rhomx2(4),ri_rhomn1(3,4) 91 real(dp) :: ri_rhomn2(3,4),ri_rhomx1(3,4),ri_rhomx2(3,4),ri_zetmn1(3,2) 92 real(dp) :: ri_zetmn2(3,2),ri_zetmx1(3,2),ri_zetmx2(3,2),zetmn1(2) 93 real(dp) :: zetmn2(2),zetmx1(2),zetmx2(2) 94 real(dp),allocatable :: array(:),coord(:,:,:,:),value(:,:,:),integrated(:) 95 real(dp),allocatable :: value_fft(:,:,:) 96 97 ! ************************************************************************* 98 99 if(.not.(present(optrhor))) then 100 ioptrhor=0 101 else 102 ioptrhor=optrhor 103 end if 104 105 if(option/=1 .and. option/=2)then 106 write(message, '(a,i0)' )' Option must be 1 or 2, while it is ',option 107 MSG_BUG(message) 108 end if 109 110 if (mpi_enreg%nproc_wvl>1) then 111 ! nfft is always the potential size (in GGA, the density has buffers). 112 nfft_ = ngfft(1) * ngfft(2) * mpi_enreg%nscatterarr(mpi_enreg%me_wvl, 2) 113 n1 = ngfft(1) 114 n2 = ngfft(2) 115 n3 = sum(mpi_enreg%nscatterarr(:, 2)) 116 istart = mpi_enreg%nscatterarr(mpi_enreg%me_wvl, 4) 117 else 118 nfft_ = nfft 119 n1=ngfft(1) ; n2=ngfft(2) ; n3=ngfft(3) 120 istart = 0 121 end if 122 123 !-------------------------------------------------------------------------- 124 !One has to determine the maximum and minimum (etc...) values 125 !over all space, and then output it, as well as to identify 126 !the point at which it occurs ... 127 !This will require a bit of data exchange, and correct indirect indexing ... 128 129 !For the local processor, find different items : 130 !maximum and minimum total electron density and locations 131 !and also spin-polarisation and magnetization 132 !also keep the second maximal or minimal value 133 if(nspden==1)nitems=1 ! Simply the total density 134 if(nspden==2)nitems=5 ! Total density, spin up, spin down, magnetization, zeta 135 if(nspden==4)nitems=6 ! Total density, x, y, z, magnetization, zeta 136 137 ABI_ALLOCATE(value,(2,2,nitems)) 138 ABI_ALLOCATE(iindex,(2,2,nitems)) 139 ABI_ALLOCATE(array,(nfft)) 140 ABI_ALLOCATE(integrated,(nitems)) 141 142 do iitems=1,nitems 143 144 ! Copy the correct values into the array 145 ! First set of items : the density, for each spin component 146 if(iitems<=nspden)then 147 array(:)=rhor(:,iitems) 148 end if 149 ! Case nspden==2, some computation to be done 150 if(nspden==2)then 151 if(iitems==3)then ! Spin down 152 array(:)=rhor(:,1)-rhor(:,2) 153 else if(iitems==4)then ! Magnetization 154 array(:)=2*rhor(:,2)-rhor(:,1) 155 else if(iitems==5)then ! zeta = relative magnetization 156 ! Avoid 0/0: the limit of (x - y) / (x+ y) depends on the direction. 157 array(:)=zero 158 where (abs(rhor(:,1)) > tol12) array(:)=(2*rhor(:,2)-rhor(:,1))/rhor(:,1) 159 end if 160 ! Case nspden==4, some other computation to be done 161 else if(nspden==4)then 162 if(iitems==5)then ! Magnetization 163 array(:)=sqrt(rhor(:,2)**2+rhor(:,3)**2+rhor(:,4)**2) 164 else if(iitems==6)then ! zeta = relative magnetization 165 array(:)=(sqrt(rhor(:,2)**2+rhor(:,3)**2+rhor(:,4)**2))/rhor(:,1) 166 end if 167 end if 168 169 ! Zero all the absolute values that are lower than tol8, for portability reasons. 170 do ifft = 1, nfft_ 171 if(abs(array(ifft))<tol8)array(ifft)=zero 172 end do 173 174 ! DEBUG 175 ! write(std_out,*) ' iitems,array(1:2)=',iitems,array(1:2) 176 ! ENDDEBUG 177 178 do indsign=1,2 ! Find alternatively the maximum and the minimum 179 iisign=3-2*indsign 180 181 if (nfft_ > 1) then 182 ! Initialize the two first values 183 value1=array(istart + 1) ; value2=array(istart + 2) 184 index1=1 ; index2=2 185 186 ! Ordering, if needed 187 if( iisign*(value2+tol12) > iisign*(value1)) then 188 temp=value2 ; value2=value1 ; value1=temp 189 index1=2 ; index2=1 190 end if 191 192 ! Integration, if relevant 193 if(present(ucvol).and. indsign==1)then 194 integrated(iitems) = array(istart + 1)+array(istart + 2) 195 end if 196 else 197 value1 = zero; value2 = zero 198 index1 = 0; index2 = 0 199 end if 200 201 ! DEBUG 202 ! write(std_out,*) ' value1,value2,index1,index2=',value1,value2,index1,index2 203 ! ENDDEBUG 204 205 ! Loop over all points 206 do ifft = 3, nfft_ 207 208 temp=array(istart + ifft) 209 if(present(ucvol).and. indsign==1)integrated(iitems) = integrated(iitems)+temp 210 ! Compares it to the second value 211 if( iisign*(temp+tol12) > iisign*value2 ) then 212 ! Compare it to the first value 213 if( iisign*(temp+tol12) > iisign*value1 ) then 214 value2=value1 ; index2=index1 215 value1=temp ; index1=ifft 216 else 217 value2=temp ; index2=ifft 218 end if 219 end if 220 221 end do ! ifft 222 223 value(1,indsign,iitems)=value1 224 value(2,indsign,iitems)=value2 225 iindex(1,indsign,iitems)=index1 226 iindex(2,indsign,iitems)=index2 227 228 ! DEBUG 229 ! write(std_out,*) ' it,v1,i1=',iitems, value1,index1 230 ! write(std_out,*) ' it,v2,i2=',iitems, value2,index2 231 ! ENDDEBUG 232 233 end do ! indsign 234 235 if(present(ucvol))then 236 nfftot=ngfft(1) * ngfft(2) * ngfft(3) 237 integrated(iitems)=integrated(iitems)*ucvol/nfftot 238 end if 239 240 ! Integrate the array 241 ! integrated(iitems)=zero 242 ! do ifft=1,nfft_ 243 ! integrated(iitems) = integrated(iitems) + array(istart + ifft) 244 ! enddo 245 ! if(present(ucvol))integrated(iitems) = integrated(iitems)*ucvol/nfft_ 246 ! write(std_err,*)present(ucvol) 247 ! if(present(ucvol))then 248 ! write(std_err,*)ucvol 249 ! endif 250 251 end do ! iitems 252 253 ABI_DEALLOCATE(array) 254 255 !------------------------------------------------------------------- 256 !Enter section for FFT parallel case 257 !if(mpi_enreg%paral_kgb>1) spaceComm=mpi_enreg%comm_fft; reduce=.true. 258 spaceComm=mpi_enreg%comm_fft; reduce=.false. 259 if(mpi_enreg%nproc_fft>1) then 260 spaceComm=mpi_enreg%comm_fft; reduce=.true. 261 else if(mpi_enreg%nproc_wvl>1) then 262 spaceComm=mpi_enreg%comm_wvl; reduce=.true. 263 end if 264 nproc=xmpi_comm_size(spaceComm) 265 me=xmpi_comm_rank(spaceComm) 266 267 if (reduce) then 268 269 ! Communicate all data to all processors with only two global communications 270 ABI_ALLOCATE(value_fft,(5,nitems,nproc)) 271 ABI_ALLOCATE(index_fft,(2,2,nitems,nproc)) 272 value_fft(:,:,:)=zero 273 index_fft(:,:,:,:)=zero 274 value_fft(1,:,me + 1)=value(1,1,:) 275 value_fft(2,:,me + 1)=value(2,1,:) 276 value_fft(3,:,me + 1)=value(1,2,:) 277 value_fft(4,:,me + 1)=value(2,2,:) 278 if(present(ucvol))value_fft(5,:,me + 1)=integrated(:) 279 index_fft(:,:,:,me + 1)=iindex(:,:,:) 280 call xmpi_sum(value_fft,spaceComm,ierr) 281 call xmpi_sum(index_fft,spaceComm,ierr) 282 283 ! Determine the global optimum and second optimum for each item 284 ! Also, the integrated quantities, if relevant. 285 do iitems=1,nitems 286 287 if(present(ucvol))integrated(iitems)=sum(value_fft(5,iitems,1:nproc)) 288 289 do indsign=1,2 ! Find alternatively the maximum and the minimum 290 iisign=3-2*indsign 291 292 ! Initialisation 293 value1=value_fft(2*indsign-1,iitems,1) 294 value2=value_fft(2*indsign ,iitems,1) 295 index1=index_fft(1,indsign,iitems,1) 296 index2=index_fft(2,indsign,iitems,1) 297 298 ! Loop 299 do iproc=1, nproc, 1 300 do ii=1,2 301 if(iproc>1 .or. ii==2)then 302 303 temp=value_fft(ii+2*(indsign-1),iitems,iproc) 304 ! Compares it to the second value 305 if( iisign*(temp+tol12) > iisign*value2 ) then 306 ! Compare it to the first value 307 if( iisign*(temp+tol12) > iisign*value1 ) then 308 value2=value1 ; index2=index1 309 value1=temp ; index1=index_fft(ii,indsign,iitems,iproc) 310 else 311 value2=temp ; index2=index_fft(ii,indsign,iitems,iproc) 312 end if 313 end if 314 315 end if ! if(iproc>1 .or. ii==2) 316 end do ! ii 317 end do ! iproc 318 319 value(1,indsign,iitems)=value1 320 value(2,indsign,iitems)=value2 321 iindex(1,indsign,iitems)=index1 322 iindex(2,indsign,iitems)=index2 323 324 end do ! iisign 325 end do ! iitems 326 327 ABI_DEALLOCATE(value_fft) 328 ABI_DEALLOCATE(index_fft) 329 330 end if !if(reduce) 331 332 !------------------------------------------------------------------- 333 334 !Determines the reduced coordinates of the min and max for each item 335 ABI_ALLOCATE(coord,(3,2,2,nitems)) 336 do iitems=1,nitems 337 do indsign=1,2 338 do ii=1,2 339 index1=iindex(ii,indsign,iitems) 340 i3=(index1-1)/n1/n2 341 i2=(index1-1-i3*n1*n2)/n1 342 i1=index1-1-i3*n1*n2-i2*n1 343 coord(1,ii,indsign,iitems)=dble(i1)/dble(n1)+tol12 344 coord(2,ii,indsign,iitems)=dble(i2)/dble(n2)+tol12 345 coord(3,ii,indsign,iitems)=dble(i3)/dble(n3)+tol12 346 ! DEBUG 347 ! write(std_out,*)' ii,indsign,iitems,coord(1:3)=',ii,indsign,iitems,coord(:,ii,indsign,iitems) 348 ! write(std_out,*)' value ', value(ii, indsign, iitems) 349 ! ENDDEBUG 350 end do 351 end do 352 end do 353 354 !------------------------------------------------------------------------- 355 !Output 356 if (mpi_enreg%paral_kgb==0.or.mpi_enreg%me_fft==0) then 357 if(.true.)then 358 do iitems=1,nitems 359 360 if(ioptrhor==4 .and. iitems>2)exit 361 362 select case (ioptrhor) 363 case(0) 364 365 if(iitems==1) write(message,'(a)')' Total charge density [el/Bohr^3]' 366 if(nspden==2)then 367 if(iitems==2) write(message,'(a)')' Spin up density [el/Bohr^3]' 368 if(iitems==3) write(message,'(a)')' Spin down density [el/Bohr^3]' 369 if(iitems==4) write(message,'(a)')' Magnetization (spin up - spin down) [el/Bohr^3]' 370 if(iitems==5) write(message,'(a)')' Relative magnetization (=zeta, between -1 and 1) ' 371 else if(nspden==4)then 372 if(iitems==2) write(message,'(a)')' x component of magnetization [el/Bohr^3]' 373 if(iitems==3) write(message,'(a)')' y component of magnetization [el/Bohr^3]' 374 if(iitems==4) write(message,'(a)')' z component of magnetization [el/Bohr^3]' 375 if(iitems==5) write(message,'(a)')' Magnetization (absolute value) [el/Bohr^3]' 376 if(iitems==6) write(message,'(a)')' Relative magnetization (=zeta, between -1 and 1) ' 377 end if 378 379 case(1) 380 381 if(iitems==1) write(message,'(a)')' Total kinetic energy density [Ha/Bohr^3]' 382 if(nspden==2)then 383 if(iitems==2) write(message,'(a)')' Spin up density [Ha/Bohr^3]' 384 if(iitems==3) write(message,'(a)')' Spin down density [Ha/Bohr^3]' 385 if(iitems==4) write(message,'(a)')' Magnetization (spin up - spin down) [Ha/Bohr^3]' 386 if(iitems==5) write(message,'(a)')' Relative magnetization (=zeta, between -1 and 1) ' 387 else if(nspden==4)then 388 if(iitems==2) write(message,'(a)')' x component of magnetization [Ha/Bohr^3]' 389 if(iitems==3) write(message,'(a)')' y component of magnetization [Ha/Bohr^3]' 390 if(iitems==4) write(message,'(a)')' z component of magnetization [Ha/Bohr^3]' 391 if(iitems==5) write(message,'(a)')' Magnetization (absolute value) [Ha/Bohr^3]' 392 if(iitems==6) write(message,'(a)')' Relative magnetization (=zeta, between -1 and 1) ' 393 end if 394 395 case(2) 396 397 if(iitems==1) write(message,'(a)')' Gradient of the electronic density [el/Bohr^4]' 398 if(nspden==2)then 399 if(iitems==2) write(message,'(a)')' Spin up density [el/Bohr^4]' 400 if(iitems==3) write(message,'(a)')' Spin down density [el/Bohr^4]' 401 if(iitems==4) write(message,'(a)')' Magnetization (spin up - spin down) [el/Bohr^4]' 402 if(iitems==5) write(message,'(a)')' Relative magnetization (=zeta, between -1 and 1) ' 403 else if(nspden==4)then 404 if(iitems==2) write(message,'(a)')' x component of magnetization [el/Bohr^4]' 405 if(iitems==3) write(message,'(a)')' y component of magnetization [el/Bohr^4]' 406 if(iitems==4) write(message,'(a)')' z component of magnetization [el/Bohr^4]' 407 if(iitems==5) write(message,'(a)')' Magnetization (absolute value) [el/Bohr^4]' 408 if(iitems==6) write(message,'(a)')' Relative magnetization (=zeta, between -1 and 1) ' 409 end if 410 411 case(3) 412 413 if(iitems==1) write(message,'(a)')' Laplacian of the electronic density [el/Bohr^5]' 414 if(nspden==2)then 415 if(iitems==2) write(message,'(a)')' Spin up density [el/Bohr^5]' 416 if(iitems==3) write(message,'(a)')' Spin down density [el/Bohr^5]' 417 if(iitems==4) write(message,'(a)')' Magnetization (spin up - spin down) [el/Bohr^5]' 418 if(iitems==5) write(message,'(a)')' Relative magnetization (=zeta, between -1 and 1) ' 419 else if(nspden==4)then 420 if(iitems==2) write(message,'(a)')' x component of magnetization [el/Bohr^5]' 421 if(iitems==3) write(message,'(a)')' y component of magnetization [el/Bohr^5]' 422 if(iitems==4) write(message,'(a)')' z component of magnetization [el/Bohr^5]' 423 if(iitems==5) write(message,'(a)')' Magnetization (absolute value) [el/Bohr^5]' 424 if(iitems==6) write(message,'(a)')' Relative magnetization (=zeta, between -1 and 1) ' 425 end if 426 427 case(4) 428 429 if(iitems==1) write(message,'(a)')' Electron Localization Function (ELF) [min:0;max:1]' 430 if(nspden==2)then 431 if(iitems==2) write(message,'(a)')' Spin up ELF [min:0;max:1]' 432 ! if(iitems==3) write(message,'(a)')' Spin down ELF [min:0;max:1]' 433 ! if(iitems==4) write(message,'(a)')' Magnetization (spin up - spin down) [el/Bohr^4]' 434 ! if(iitems==5) write(message,'(a)')' Relative magnetization (=zeta, between -1 and 1) ' 435 else if(nspden==4)then 436 ! if(iitems==2) write(message,'(a)')' x component of magnetization [el/Bohr^4]' 437 ! if(iitems==3) write(message,'(a)')' y component of magnetization [el/Bohr^4]' 438 ! if(iitems==4) write(message,'(a)')' z component of magnetization [el/Bohr^4]' 439 ! if(iitems==5) write(message,'(a)')' Magnetization (spin up - spin down) [el/Bohr^4]' 440 ! if(iitems==6) write(message,'(a)')' Relative magnetization (=zeta, between -1 and 1) ' 441 end if 442 443 444 end select 445 446 call wrtout(iout,message,'COLL') 447 448 write(message,'(a,es13.4,a,3f10.4)') ' Maximum= ',& 449 & value(1,1,iitems),' at reduced coord.',coord(:,1,1,iitems) 450 call wrtout(iout,message,'COLL') 451 if(option==2)then 452 write(message,'(a,es13.4,a,3f10.4)')' Next maximum= ',& 453 & value(2,1,iitems),' at reduced coord.',coord(:,2,1,iitems) 454 call wrtout(iout,message,'COLL') 455 end if 456 write(message,'(a,es13.4,a,3f10.4)') ' Minimum= ',& 457 & value(1,2,iitems),' at reduced coord.',coord(:,1,2,iitems) 458 call wrtout(iout,message,'COLL') 459 if(option==2)then 460 write(message,'(a,es13.4,a,3f10.4)')' Next minimum= ',& 461 & value(2,2,iitems),' at reduced coord.',coord(:,2,2,iitems) 462 call wrtout(iout,message,'COLL') 463 end if 464 if(present(ucvol))then 465 if(.not.(nspden==2.and.iitems==5) .and. .not.(nspden==4.and.iitems==6))then 466 if(abs(integrated(iitems))<tol10)integrated(iitems)=zero 467 write(message,'(a,es13.4)')' Integrated= ',integrated(iitems) 468 call wrtout(iout,message,'COLL') 469 end if 470 end if 471 472 end do ! iitems 473 end if 474 475 if(.false.)then 476 477 select case(optrhor) 478 case(0) 479 write(txt1_in_mssg, '(a)')" Min el dens=" 480 write(txt2_in_mssg, '(a)')" el/bohr^3 at reduced coord." 481 write(txt3_in_mssg, '(a)')" Max el dens=" 482 case(1) 483 write(txt1_in_mssg, '(a)')" Min kin energy dens=" 484 write(txt2_in_mssg, '(a)')" bohr^(-5) at reduced coord." 485 write(txt3_in_mssg, '(a)')" Max kin energy dens=" 486 end select 487 488 write(message, '(a,a,1p,e12.4,a,0p,3f8.4)' ) ch10,& 489 & trim(txt1_in_mssg),value(1,2,1),& 490 & trim(txt2_in_mssg),coord(:,1,2,1) 491 call wrtout(iout,message,'COLL') 492 if(option==2)then 493 write(message, '(a,1p,e12.4,a,0p,3f8.4)' ) & 494 & ', next min=',value(2,2,1),& 495 & trim(txt2_in_mssg),coord(:,2,2,1) 496 call wrtout(iout,message,'COLL') 497 end if 498 write(message, '(a,1p,e12.4,a,0p,3f8.4)' )& 499 & trim(txt3_in_mssg),value(1,1,1),& 500 & trim(txt2_in_mssg),coord(:,1,1,1) 501 call wrtout(iout,message,'COLL') 502 if(option==2)then 503 write(message, '(a,1p,e12.4,a,0p,3f8.4)' )& 504 & ', next max=',value(2,1,1),& 505 & trim(txt2_in_mssg),coord(:,2,1,1) 506 call wrtout(iout,message,'COLL') 507 end if 508 509 if(nspden>=2)then 510 write(message, '(a,a,1p,e12.4,a,0p,3f8.4)' ) ch10,& 511 & ',Min spin pol zeta=',value(1,2,4+nspden/2),& 512 & ' at reduced coord.',coord(:,1,2,4+nspden/2) 513 call wrtout(iout,message,'COLL') 514 if(option==2)then 515 write(message, '(a,1p,e12.4,a,0p,3f8.4)' )& 516 & ', next min=',value(2,2,4+nspden/2),& 517 & ' at reduced coord.',coord(:,2,2,4+nspden/2) 518 call wrtout(iout,message,'COLL') 519 end if 520 write(message, '(a,1p,e12.4,a,0p,3f8.4)' )& 521 & ',Max spin pol zeta=',value(1,1,4+nspden/2),& 522 & ' at reduced coord.',coord(:,1,1,4+nspden/2) 523 call wrtout(iout,message,'COLL') 524 if(option==2)then 525 write(message, '(a,1p,e12.4,a,0p,3f8.4)' )& 526 & ', next max=',value(2,1,4+nspden/2),& 527 & ' at reduced coord.',coord(:,2,1,4+nspden/2) 528 call wrtout(iout,message,'COLL') 529 end if 530 end if ! nspden 531 532 end if ! second section always true 533 534 if(nspden==2 .and. .false.)then 535 write(message,'(a)')& 536 & ' Position in reduced coord. ( x y z )' 537 call wrtout(iout,message,'COLL') 538 write(message,'(a,es13.4,a,3f10.4)')' Minimum (Total el-den) : [el/Bohr^3]',& 539 & rhomn1(1),' at',ri_rhomn1(1,1),ri_rhomn1(2,1),ri_rhomn1(3,1) 540 call wrtout(iout,message,'COLL') 541 write(message,'(a,es13.4,a,3f10.4)')' Minimum (Spin-up den) : [el/Bohr^3]',& 542 & rhomn1(2),' at',ri_rhomn1(1,2),ri_rhomn1(2,2),ri_rhomn1(3,2) 543 call wrtout(iout,message,'COLL') 544 write(message,'(a,es13.4,a,3f10.4)')' Minimum (Spin-down den) : [el/Bohr^3]',& 545 & zetmn1(1),' at',ri_zetmn1(1,1),ri_zetmn1(2,1),ri_zetmn1(3,1) 546 call wrtout(iout,message,'COLL') 547 write(message,'(a,es13.4,a,3f10.4)')' Minimum (Spin pol zeta) : [m/|m|] ',& 548 & zetmn1(2),' at',ri_zetmn1(1,2),ri_zetmn1(2,2),ri_zetmn1(3,2) 549 call wrtout(iout,message,'COLL') 550 if(option==2)then 551 write(message,'(a,es13.4,a,3f10.4)')' Next minimum (Total el-den) : [el/Bohr^3]',& 552 & rhomn2(1),' at',ri_rhomn2(1,1),ri_rhomn2(2,1),ri_rhomn2(3,1) 553 call wrtout(iout,message,'COLL') 554 write(message,'(a,es13.4,a,3f10.4)')' Next minimum (Spin-up den) : [el/Bohr^3]',& 555 & rhomn2(2),' at',ri_rhomn2(1,2),ri_rhomn2(2,2),ri_rhomn2(3,2) 556 call wrtout(iout,message,'COLL') 557 write(message,'(a,es13.4,a,3f10.4)')' Next minimum (Spin-down den) : [el/Bohr^3]',& 558 & zetmn2(1),' at',ri_zetmn2(1,1),ri_zetmn2(2,1),ri_zetmn2(3,1) 559 call wrtout(iout,message,'COLL') 560 write(message,'(a,es13.4,a,3f10.4)')' Next minimum (Spin pol zeta) : [m/|m|] ',& 561 & zetmn2(2),' at',ri_zetmn2(1,2),ri_zetmn2(2,2),ri_zetmn2(3,2) 562 call wrtout(iout,message,'COLL') 563 end if 564 write(message,*)' ' 565 call wrtout(iout,message,'COLL') 566 write(message,'(a,es13.4,a,3f10.4)')' Maximum (Total el-den) : [el/Bohr^3]',& 567 & rhomx1(1),' at',ri_rhomx1(1,1),ri_rhomx1(2,1),ri_rhomx1(3,1) 568 call wrtout(iout,message,'COLL') 569 write(message,'(a,es13.4,a,3f10.4)')' Maximum (Spin-up den) : [el/Bohr^3]',& 570 & rhomx1(2),' at',ri_rhomx1(1,2),ri_rhomx1(2,2),ri_rhomx1(3,2) 571 call wrtout(iout,message,'COLL') 572 write(message,'(a,es13.4,a,3f10.4)')' Maximum (Spin-down den) : [el/Bohr^3]',& 573 & zetmx1(1),' at',ri_zetmx1(1,1),ri_zetmx1(2,1),ri_zetmx1(3,1) 574 call wrtout(iout,message,'COLL') 575 write(message,'(a,es13.4,a,3f10.4)')' Maximum (Spin pol zeta) : [m/|m|] ',& 576 & zetmx1(2),' at',ri_zetmx1(1,2),ri_zetmx1(2,2),ri_zetmx1(3,2) 577 call wrtout(iout,message,'COLL') 578 if(option==2)then 579 write(message,'(a,es13.4,a,3f10.4)')' Next maximum (Total el-den) : [el/Bohr^3]',& 580 & rhomx2(1),' at',ri_rhomx2(1,1),ri_rhomx2(2,1),ri_rhomx2(3,1) 581 call wrtout(iout,message,'COLL') 582 write(message,'(a,es13.4,a,3f10.4)')' Next maximum (Spin-up den) : [el/Bohr^3]',& 583 & rhomx2(2),' at',ri_rhomx2(1,2),ri_rhomx2(2,2),ri_rhomx2(3,2) 584 call wrtout(iout,message,'COLL') 585 write(message,'(a,es13.4,a,3f10.4)')' Next maximum (Spin-down den) : [el/Bohr^3]',& 586 & zetmx2(1),' at',ri_zetmx2(1,1),ri_zetmx2(2,1),ri_zetmx2(3,1) 587 call wrtout(iout,message,'COLL') 588 write(message,'(a,es13.4,a,3f10.4)')' Next maximum (Spin pol zeta) : [m/|m|] ',& 589 & zetmx2(2),' at',ri_zetmx2(1,2),ri_zetmx2(2,2),ri_zetmx2(3,2) 590 call wrtout(iout,message,'COLL') 591 end if 592 end if 593 594 if(nspden==4 .and. .false.)then 595 write(message,'(a)')& 596 & ' Position in reduced coord. ( x y z )' 597 call wrtout(iout,message,'COLL') 598 write(message,'(a,es13.4,a,3f10.4)')' Minimum (Total el-den) : [el/Bohr^3]',& 599 & rhomn1(1),' at',ri_rhomn1(1,1),ri_rhomn1(2,1),ri_rhomn1(3,1) 600 call wrtout(iout,message,'COLL') 601 write(message,'(a,es13.4,a,3f10.4)')' Minimum (Magnetizat.-x) : [m/|m|] ',& 602 & rhomn1(2),' at',ri_rhomn1(1,2),ri_rhomn1(2,2),ri_rhomn1(3,2) 603 call wrtout(iout,message,'COLL') 604 write(message,'(a,es13.4,a,3f10.4)')' Minimum (Magnetizat.-y) : [m/|m|] ',& 605 & rhomn1(3),' at',ri_rhomn1(1,3),ri_rhomn1(2,3),ri_rhomn1(3,3) 606 call wrtout(iout,message,'COLL') 607 write(message,'(a,es13.4,a,3f10.4)')' Minimum (Magnetizat.-z) : [m/|m|] ',& 608 & rhomn1(4),' at',ri_rhomn1(1,4),ri_rhomn1(2,4),ri_rhomn1(3,4) 609 call wrtout(iout,message,'COLL') 610 write(message,'(a,es13.4,a,3f10.4)')' Minimum (Spin pol zeta) : [m/|m|] ',& 611 & zetmn1(1),' at',ri_zetmn1(1,1),ri_zetmn1(2,1),ri_zetmn1(3,1) 612 call wrtout(iout,message,'COLL') 613 if(option==2)then 614 write(message,'(a,es13.4,a,3f10.4)')' Next-Minimum (Total el-den) : [el/Bohr^3]',& 615 & rhomn2(1),' at',ri_rhomn2(1,1),ri_rhomn2(2,1),ri_rhomn2(3,1) 616 call wrtout(iout,message,'COLL') 617 write(message,'(a,es13.4,a,3f10.4)')' Next-Minimum (Magnetizat.-x) : [m/|m|] ',& 618 & rhomn2(2),' at',ri_rhomn2(1,2),ri_rhomn2(2,2),ri_rhomn2(3,2) 619 call wrtout(iout,message,'COLL') 620 write(message,'(a,es13.4,a,3f10.4)')' Next-Minimum (Magnetizat.-y) : [m/|m|] ',& 621 & rhomn2(3),' at',ri_rhomn2(1,3),ri_rhomn2(2,3),ri_rhomn2(3,3) 622 call wrtout(iout,message,'COLL') 623 write(message,'(a,es13.4,a,3f10.4)')' Next-Minimum (Magnetizat.-z) : [m/|m|] ',& 624 & rhomn2(4),' at',ri_rhomn2(1,4),ri_rhomn2(2,4),ri_rhomn2(3,4) 625 call wrtout(iout,message,'COLL') 626 write(message,'(a,es13.4,a,3f10.4)')' Next-Minimum (Spin pol zeta) : [m/|m|] ',& 627 & zetmn2(1),' at',ri_zetmn2(1,1),ri_zetmn2(2,1),ri_zetmn2(3,1) 628 call wrtout(iout,message,'COLL') 629 end if 630 write(message,*)' ' 631 call wrtout(iout,message,'COLL') 632 write(message,'(a,es13.4,a,3f10.4)')' Maximum (Total el-den) : [el/Bohr^3]',& 633 & rhomx1(1),' at',ri_rhomx1(1,1),ri_rhomx1(2,1),ri_rhomx1(3,1) 634 call wrtout(iout,message,'COLL') 635 write(message,'(a,es13.4,a,3f10.4)')' Maximum (Magnetizat.-x) : [m/|m|] ',& 636 & rhomx1(2),' at',ri_rhomx1(1,2),ri_rhomx1(2,2),ri_rhomx1(3,2) 637 call wrtout(iout,message,'COLL') 638 write(message,'(a,es13.4,a,3f10.4)')' Maximum (Magnetizat.-y) : [m/|m|] ',& 639 & rhomx1(3),' at',ri_rhomx1(1,3),ri_rhomx1(2,3),ri_rhomx1(3,3) 640 call wrtout(iout,message,'COLL') 641 write(message,'(a,es13.4,a,3f10.4)')' Maximum (Magnetizat.-z) : [m/|m|] ',& 642 & rhomx1(4),' at',ri_rhomx1(1,4),ri_rhomx1(2,4),ri_rhomx1(3,4) 643 call wrtout(iout,message,'COLL') 644 write(message,'(a,es13.4,a,3f10.4)')' Maximum (Spin pol zeta) : [m/|m|] ',& 645 & zetmx1(1),' at',ri_zetmx1(1,1),ri_zetmx1(2,1),ri_zetmx1(3,1) 646 call wrtout(iout,message,'COLL') 647 if(option==2)then 648 write(message,'(a,es13.4,a,3f10.4)')' Next-Maximum (Total el-den) : [el/Bohr^3]',& 649 & rhomx2(1),' at',ri_rhomx2(1,1),ri_rhomx2(2,1),ri_rhomx2(3,1) 650 call wrtout(iout,message,'COLL') 651 write(message,'(a,es13.4,a,3f10.4)')' Next-Maximum (Magnetizat.-x) : [m/|m|] ',& 652 & rhomx2(2),' at',ri_rhomx2(1,2),ri_rhomx2(2,2),ri_rhomx2(3,2) 653 call wrtout(iout,message,'COLL') 654 write(message,'(a,es13.4,a,3f10.4)')' Next-Maximum (Magnetizat.-y) : [m/|m|] ',& 655 & rhomx2(3),' at',ri_rhomx2(1,3),ri_rhomx2(2,3),ri_rhomx2(3,3) 656 call wrtout(iout,message,'COLL') 657 write(message,'(a,es13.4,a,3f10.4)')' Next-Maximum (Magnetizat.-z) : [m/|m|] ',& 658 & rhomx2(4),' at',ri_rhomx2(1,4),ri_rhomx2(2,4),ri_rhomx2(3,4) 659 call wrtout(iout,message,'COLL') 660 write(message,'(a,es13.4,a,3f10.4)')' Next-Maximum (Spin pol zeta) : [m/|m|] ',& 661 & zetmx2(1),' at',ri_zetmx2(1,1),ri_zetmx2(2,1),ri_zetmx2(3,1) 662 call wrtout(iout,message,'COLL') 663 end if 664 end if 665 end if 666 667 ABI_DEALLOCATE(coord) 668 ABI_DEALLOCATE(value) 669 ABI_DEALLOCATE(iindex) 670 ABI_DEALLOCATE(integrated) 671 672 end subroutine prtrhomxmn