TABLE OF CONTENTS
- ABINIT/m_self
- m_self/alloc_self
- m_self/dc_self
- m_self/destroy_self
- m_self/initialize_self
- m_self/kramerskronig_self
- m_self/make_qmcshift_self
- m_self/new_self
- m_self/print_self
- m_self/rw_self
- m_self/self_type
- m_self/selfreal2imag_self
ABINIT/m_self [ Modules ]
NAME
m_self
FUNCTION
COPYRIGHT
Copyright (C) 2006-2024 ABINIT group (BAmadon) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
SOURCE
15 #if defined HAVE_CONFIG_H 16 #include "config.h" 17 #endif 18 19 20 #include "abi_common.h" 21 22 MODULE m_self 23 24 use defs_basis 25 use m_xmpi 26 use m_errors 27 use m_abicore 28 29 use m_fstrings, only : int2char4 30 use m_crystal, only : crystal_t 31 use m_hu, only : hu_type 32 use m_io_tools, only : get_unit 33 use m_pawang, only : pawang_type 34 use m_paw_dmft, only : paw_dmft_type 35 use m_matlu, only : matlu_type, copy_matlu,shift_matlu,diag_matlu,rotate_matlu,init_matlu,destroy_matlu,print_matlu,zero_matlu 36 use m_oper, only : oper_type,init_oper,destroy_oper, loc_oper, print_oper 37 use m_datafordmft, only : compute_levels 38 39 implicit none 40 41 private 42 43 public :: alloc_self 44 public :: initialize_self 45 public :: destroy_self 46 public :: print_self 47 public :: rw_self 48 public :: dc_self 49 public :: new_self 50 public :: make_qmcshift_self 51 public :: selfreal2imag_self
m_self/alloc_self [ Functions ]
[ Top ] [ m_self ] [ Functions ]
NAME
alloc_self
FUNCTION
Allocate variables used in type self_type.
INPUTS
self <type(self_type)>= variables related to self-energy paw_dmft <type(paw_dmft_type)> = variables related to self-consistent DFT+DMFT calculations. opt_oper = 1 Allocate only quantities in the KS basis. 2 Allocate only quantities in the local basis. 3 Allocate quantities in both the KS and local basis. wtype = "real" Self energy will be computed for real frequencies = "imag" Self energy will be computed for imaginary frequencies OUTPUTS self <type(self_type)>= variables related to self-energy
SOURCE
126 subroutine alloc_self(self,paw_dmft,opt_oper,wtype) 127 128 !Arguments ------------------------------------ 129 !scalars 130 !type 131 type(self_type), intent(inout) :: self 132 type(paw_dmft_type), intent(in) :: paw_dmft 133 integer, optional, intent(in) :: opt_oper 134 character(len=4), optional :: wtype 135 !local variables ------------------------------------ 136 integer :: ifreq,optoper 137 138 !************************************************************************ 139 if(present(opt_oper)) then 140 optoper=opt_oper 141 else 142 optoper=2 143 endif 144 if(present(wtype)) then 145 self%w_type=wtype 146 else 147 self%w_type="imag" 148 endif 149 if(self%w_type=="imag") then 150 self%nw=paw_dmft%dmft_nwlo 151 self%omega=>paw_dmft%omega_lo 152 else if(self%w_type=="real") then 153 self%nw=size(paw_dmft%omega_r) 154 self%omega=>paw_dmft%omega_r 155 endif 156 157 self%dmft_nwlo=paw_dmft%dmft_nwlo 158 self%dmft_nwli=paw_dmft%dmft_nwli 159 self%iself_cv=0 160 161 call init_oper(paw_dmft,self%hdc,opt_ksloc=optoper) 162 ABI_MALLOC(self%oper,(self%nw)) 163 do ifreq=1,self%nw 164 call init_oper(paw_dmft,self%oper(ifreq),opt_ksloc=optoper) 165 enddo 166 167 if(paw_dmft%dmft_solv==4) then 168 ABI_MALLOC(self%qmc_shift,(paw_dmft%natom)) 169 ABI_MALLOC(self%qmc_xmu,(paw_dmft%natom)) 170 self%qmc_shift(:)=zero 171 self%qmc_xmu(:)=zero 172 endif 173 174 end subroutine alloc_self
m_self/dc_self [ Functions ]
[ Top ] [ m_self ] [ Functions ]
NAME
dc_self
FUNCTION
INPUTS
charge_loc(cryst_struc%natom,paw_dmft%nsppol+1)= total charge for correlated electrons on a given atom, and for spin cryst_struc <type(crystal_t)>=variables related to crystal structure hu <type(hu_type)>= variables related to the interaction between electrons self <type(self_type)>= variables related to self-energy dmft_dc = 1 Full localized Limit double counting. 2 Around Mean Field (without SO) 0 not double counting prtopt = integer which precises the amount of printing (not used here)
OUTPUT
self <type(self_type)>= variables related to self-energy hu <type(hu_type)>= variables related to the interaction between electrons
SOURCE
351 subroutine dc_self(charge_loc,cryst_struc,hu,self,dmft_dc,prtopt) 352 353 !Arguments ------------------------------------ 354 !type 355 type(crystal_t),intent(in) :: cryst_struc 356 type(self_type),intent(inout) :: self 357 real(dp), intent(in) :: charge_loc(cryst_struc%natom,self%hdc%nsppol+1) 358 type(hu_type),intent(inout) :: hu(cryst_struc%ntypat) 359 integer, intent(in) :: prtopt,dmft_dc 360 361 !Local variables------------------------------- 362 integer :: iatom,isppol,ispinor,lpawu,m1,nspinor,nsppol 363 real(dp) :: ntot,jpawu_dc 364 character(len=500) :: message 365 ! ********************************************************************* 366 nspinor = self%hdc%nspinor 367 nsppol = self%hdc%nsppol 368 369 do iatom=1,self%hdc%natom 370 lpawu=self%hdc%matlu(iatom)%lpawu 371 ntot=charge_loc(iatom,nsppol+1) 372 ! nup=charge_loc(iatom,1) 373 ! ndn=charge_loc(iatom,nsppol) 374 if(lpawu/=-1) then 375 self%hdc%matlu(iatom)%mat=czero 376 do isppol = 1 , nsppol 377 do ispinor = 1 , nspinor 378 do m1=1, 2*lpawu+1 379 if(dmft_dc==1.or.dmft_dc==4.or.dmft_dc==5) then ! FLL 380 if(dmft_dc==1.or.dmft_dc==5) then 381 jpawu_dc=hu(cryst_struc%typat(iatom))%jpawu 382 else if (dmft_dc==4) then ! FLL 383 jpawu_dc=zero 384 endif 385 if(nspinor==2.or.dmft_dc==5) then 386 self%hdc%matlu(iatom)%mat(m1,m1,isppol,ispinor,ispinor) = & 387 & hu(cryst_struc%typat(iatom))%upawu * ( ntot - one / two ) & 388 & - one/two * jpawu_dc * ( ntot - one ) 389 else 390 self%hdc%matlu(iatom)%mat(m1,m1,isppol,ispinor,ispinor)= & 391 & hu(cryst_struc%typat(iatom))%upawu * ( ntot - one / two ) & 392 !& - hu(cryst_struc%typat(iatom))%jpawu * ( charge_loc(iatom,2-nsppol+1) - one) 393 & - jpawu_dc * ( charge_loc(iatom,isppol) - one / two ) 394 endif 395 else if(dmft_dc==2.or.dmft_dc==6) then ! AMF 396 if(nspinor==2.or.dmft_dc==6) then 397 !write(message,'(a,i4,i4,2x,e20.10)') " AMF Double counting not implemented for SO" 398 !ABI_ERROR(message) 399 self%hdc%matlu(iatom)%mat(m1,m1,isppol,ispinor,ispinor)= & 400 hu(cryst_struc%typat(iatom))%upawu * ntot/two & 401 + ( hu(cryst_struc%typat(iatom))%upawu - hu(cryst_struc%typat(iatom))%jpawu )& 402 *ntot/two*(float(2*lpawu))/(float(2*lpawu+1)) 403 write(message,'(a,i4,i4,2x,e20.10)') " AMF Double counting is under test for SOC" 404 ABI_WARNING(message) 405 else 406 if(nsppol==2) then 407 self%hdc%matlu(iatom)%mat(m1,m1,isppol,ispinor,ispinor)= & 408 & hu(cryst_struc%typat(iatom))%upawu * charge_loc(iatom,2-isppol+1) & 409 & + (hu(cryst_struc%typat(iatom))%upawu - hu(cryst_struc%typat(iatom))%jpawu )& 410 & *charge_loc(iatom,isppol)*(float(2*lpawu))/(float(2*lpawu+1)) 411 else if(nsppol==1) then 412 self%hdc%matlu(iatom)%mat(m1,m1,isppol,ispinor,ispinor)= & 413 & hu(cryst_struc%typat(iatom))%upawu * charge_loc(iatom,isppol) & 414 & + (hu(cryst_struc%typat(iatom))%upawu - hu(cryst_struc%typat(iatom))%jpawu )& 415 & *charge_loc(iatom,isppol)*(float(2*lpawu))/(float(2*lpawu+1)) 416 endif 417 ! write(std_out,*) "AMF", charge_loc(iatom,2-isppol+1) 418 ! write(std_out,*) "AMF", charge_loc(iatom,isppol+1) 419 ! write(std_out,*) "AMF", lpawu 420 ! write(std_out,*) "AMF", hu(cryst_struc%typat(iatom))%upawu 421 ! write(std_out,*) "AMF", self%hdc%matlu(iatom)%mat(m1,m1,isppol,ispinor,ispinor) 422 endif 423 else 424 ABI_ERROR("not implemented") 425 endif 426 enddo ! m1 427 enddo ! ispinor 428 enddo ! isppol 429 endif 430 enddo ! iatom 431 432 if(prtopt>0) then 433 endif 434 435 436 end subroutine dc_self
m_self/destroy_self [ Functions ]
[ Top ] [ m_self ] [ Functions ]
NAME
destroy_self
FUNCTION
deallocate self
INPUTS
self <type(self_type)>= variables related to self-energy
OUTPUT
SOURCE
249 subroutine destroy_self(self) 250 251 !Arguments ------------------------------------ 252 !scalars 253 type(self_type),intent(inout) :: self 254 255 !local variables------------------------------- 256 integer :: ifreq 257 258 ! ********************************************************************* 259 if ( allocated(self%oper)) then 260 do ifreq=1,self%nw 261 call destroy_oper(self%oper(ifreq)) 262 enddo 263 ABI_FREE(self%oper) 264 end if 265 266 call destroy_oper(self%hdc) 267 if (allocated(self%qmc_shift)) then 268 ABI_FREE(self%qmc_shift) 269 end if 270 if (allocated(self%qmc_xmu)) then 271 ABI_FREE(self%qmc_xmu) 272 end if 273 self%omega => null() 274 275 end subroutine destroy_self
m_self/initialize_self [ Functions ]
[ Top ] [ m_self ] [ Functions ]
NAME
initialize_self
FUNCTION
Initialize self-energy.
INPUTS
cryst_struc <type(crystal_t)>=variables related to crystal structure self <type(self_type)>= variables related to self-energy paw_dmft <type(paw_dmft_type)> = variables related to self-consistent DFT+DMFT calculations. opt_read = not used for the moment wtype = "real" Self energy will be computed for real frequencies = "imag" Self energy will be computed for imaginary frequencies OUTPUTS self <type(self_type)>= variables related to self-energy
SOURCE
198 subroutine initialize_self(self,paw_dmft,wtype) 199 200 !Arguments ------------------------------------ 201 !scalars 202 !type 203 type(self_type), intent(inout) :: self 204 type(paw_dmft_type), intent(inout) :: paw_dmft 205 character(len=4), optional, intent(in) :: wtype 206 !local variables ------------------------------------ 207 ! character(len=500) :: message 208 integer :: iatom,ifreq 209 character(len=4) :: wtype2 210 !************************************************************************ 211 if(present(wtype)) then 212 wtype2=wtype 213 else 214 wtype2="imag" 215 endif 216 217 218 call alloc_self(self,paw_dmft,opt_oper=2,wtype=wtype2) ! opt_oper=1 is not useful and not implemented 219 do ifreq=1,self%nw 220 do iatom=1,paw_dmft%natom 221 self%oper(ifreq)%matlu(iatom)%mat=czero 222 enddo 223 enddo 224 ! if(paw_dmft%dmft_rslf==1.and.opt_read==1) then 225 ! call rw_self(cryst_struc,self,mpi_enreg,paw_dmft,pawtab,pawprtvol=2,opt_rw=1) 226 ! endif 227 ! write(message,'(a,a)') ch10," Self-energy for large frequency is" 228 ! call wrtout(std_out,message,'COLL') 229 ! call print_matlu(self%oper(paw_dmft%dmft_nwlo)%matlu, & 230 !& paw_dmft%natom,3) 231 232 end subroutine initialize_self
m_self/kramerskronig_self [ Functions ]
[ Top ] [ m_self ] [ Functions ]
NAME
kramerskronig_self
FUNCTION
INPUTS
hu <type(hu_type)> = U interaction paw_dmft <type(paw_dmft_type)> = paw+dmft related data
OUTPUT
self%qmc_shift in self <type(self_type)> = Self-energy
SOURCE
1475 subroutine kramerskronig_self(self,selflimit,selfhdc,filapp) 1476 1477 !Arguments ------------------------------------ 1478 !type 1479 type(self_type),intent(inout) :: self 1480 type(matlu_type),intent(in) :: selflimit(self%hdc%natom) 1481 type(matlu_type),intent(in) :: selfhdc(self%hdc%natom) 1482 character(len=fnlen), intent(in) :: filapp 1483 1484 !Local variables------------------------------- 1485 integer :: ifreq,jfreq,isppol,ispinor,ispinor1,im,im1,iatom 1486 real(dp), allocatable :: selftemp_re(:) 1487 real(dp), allocatable :: selftemp_imag(:) 1488 integer :: natom,ndim,nsppol,nspinor 1489 real(dp) :: delta 1490 character(len=500) :: message 1491 ! ********************************************************************* 1492 delta=0.0000000 1493 ABI_MALLOC(selftemp_re,(self%nw)) 1494 ABI_MALLOC(selftemp_imag,(self%nw)) 1495 natom=self%hdc%natom 1496 nsppol = self%hdc%nsppol 1497 nspinor=self%hdc%nspinor 1498 write(message,'(2a,i4)') ch10,' ------ Limit of real part of Self' 1499 call wrtout(std_out, message,'COLL') 1500 1501 call print_matlu(selflimit,natom,3) 1502 1503 !print norms 1504 write(message,'(2a,i4)') ch10,' ------ Double counting' 1505 call wrtout(std_out, message,'COLL') 1506 1507 call print_matlu(selfhdc,natom,3) 1508 open(unit=67,file=trim(filapp)//"_DFTDMFT_Self_realaxis_from_maxent_and_kramerskronig.dat", status='unknown',form='formatted') 1509 rewind(67) 1510 ! Compute limit of Real Part and put in double counting energy. 1511 ! call copy_matlu(selfhdc,self%hdc%matlu,natom) 1512 !!write(6,*) "selfhdc kramerskronig",selfhdc(1)%mat(1,1,1,1,1) 1513 !write(6,*) "selfr%hdc kramerskronig",self%hdc%matlu(1)%mat(1,1,1,1,1) 1514 do iatom=1,natom 1515 if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then 1516 ndim=2*self%oper(1)%matlu(iatom)%lpawu+1 1517 do isppol=1,nsppol 1518 do ispinor=1,nspinor 1519 do ispinor1=1,nspinor 1520 do im=1,ndim 1521 do im1=1,ndim 1522 !write(6,*) "realpart",real(selflimit(iatom)%mat(im,im1,isppol,ispinor,ispinor1)) 1523 do ifreq=1,self%nw 1524 selftemp_re(ifreq)=zero 1525 selftemp_imag(ifreq)=aimag(self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)) 1526 do jfreq=1,self%nw-1 1527 if(jfreq==ifreq) cycle 1528 ! selftemp_re(ifreq)=selftemp_re(ifreq) - & 1529 ! & aimag(self%oper(jfreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)) & 1530 ! & *(self%omega(ifreq)-self%omega(jfreq)) & 1531 ! & /((self%omega(ifreq)-self%omega(jfreq))**2+delta**2)& 1532 ! & *(self%omega(jfreq+1)-self%omega(jfreq)) 1533 selftemp_re(ifreq)=selftemp_re(ifreq) - & 1534 & aimag(self%oper(jfreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)) & 1535 & /(self%omega(ifreq)-self%omega(jfreq)) * (self%omega(jfreq+1)-self%omega(jfreq)) 1536 enddo 1537 selftemp_re(ifreq)=selftemp_re(ifreq)/pi 1538 !write(671,*) self%omega(ifreq),selftemp_re(ifreq),selftemp_imag(ifreq) 1539 enddo 1540 ! TEST************************* 1541 ! do ifreq=1,self%nw 1542 ! selftemp_imag(ifreq)=zero 1543 ! do jfreq=1,self%nw-1 1544 ! if(jfreq==ifreq) cycle 1545 !! selftemp_re(ifreq)=selftemp_re(ifreq) - & 1546 !! & aimag(self%oper(jfreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)) & 1547 !! & *(self%omega(ifreq)-self%omega(jfreq)) & 1548 !! & /((self%omega(ifreq)-self%omega(jfreq))**2+delta**2)& 1549 !! & *(self%omega(jfreq+1)-self%omega(jfreq)) 1550 ! selftemp_imag(ifreq)=selftemp_imag(ifreq) + & 1551 ! & selftemp_re(jfreq) & 1552 ! & /(self%omega(ifreq)-self%omega(jfreq)) * (self%omega(jfreq+1)-self%omega(jfreq)) 1553 ! enddo 1554 ! selftemp_imag(ifreq)=selftemp_imag(ifreq)/pi 1555 ! write(672,*) self%omega(ifreq),selftemp_re(ifreq),selftemp_imag(ifreq) 1556 ! enddo 1557 ! TEST************************* 1558 ! write(6,*) "TWO FACTOR IS PUT BECAUSE OF MAXENT CODE ??" 1559 do ifreq=1,self%nw 1560 ! write(68,*) self%omega(ifreq),selftemp_re(ifreq),selftemp_imag(ifreq) 1561 selftemp_re(ifreq)=selftemp_re(ifreq)+ & 1562 & real(selflimit(iatom)%mat(im,im1,isppol,ispinor,ispinor1)- & 1563 & selfhdc(iatom)%mat(im,im1,isppol,ispinor,ispinor1)) 1564 self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)& 1565 & =cmplx(selftemp_re(ifreq),selftemp_imag(ifreq),kind=dp)/two 1566 !& =cmplx(0.d0,selftemp_imag(ifreq),kind=dp)/two 1567 ! & =cmplx(selftemp_re(ifreq),0.d0,kind=dp)/two 1568 ! self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)& 1569 !& =cmplx(selftemp_re(ifreq),0.d0,kind=dp)/two 1570 !& =cmplx(0.d0,0.d0,kind=dp)/two 1571 ! The factor two is here to compensate for the factor two in OmegaMaxent.. 1572 ! & =cmplx(selftemp_re(ifreq),0.0,kind=dp) 1573 write(67,*) self%omega(ifreq),real(self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1))& 1574 ,aimag(self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)) 1575 enddo 1576 write(67,*) 1577 !write(68,*) 1578 !!!!!!!!!! Z renormalization 1579 ! i0=389 1580 ! slope=(selftemp_re(i0+1)-selftemp_re(i0))/& 1581 ! (self%omega(i0+1)-self%omega(i0)) 1582 ! y0= selftemp_re(i0) 1583 ! do ifreq=1,self%nw 1584 ! selftemp_re(ifreq)=slope * (self%omega(ifreq)-self%omega(i0)) + y0 1585 ! selftemp_imag(ifreq)=zero 1586 ! self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)& 1587 ! & =cmplx(selftemp_re(ifreq),selftemp_imag(ifreq),kind=dp)/two 1588 ! write(6777,*) self%omega(ifreq),real(self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)),aimag(self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)) 1589 ! enddo 1590 !!!!!!!!!! 1591 enddo 1592 enddo 1593 enddo 1594 enddo 1595 enddo ! isppol 1596 endif ! lpawu=/-1 1597 enddo ! iatom 1598 close(67) 1599 !write(6,*) "self1",aimag(self%oper(489)%matlu(1)%mat(1,1,1,1,1)) 1600 1601 ABI_FREE(selftemp_re) 1602 ABI_FREE(selftemp_imag) 1603 1604 1605 end subroutine kramerskronig_self
m_self/make_qmcshift_self [ Functions ]
[ Top ] [ m_self ] [ Functions ]
NAME
make_qmcshift_hu
FUNCTION
INPUTS
hu <type(hu_type)> = U interaction paw_dmft <type(paw_dmft_type)> = paw+dmft related data
OUTPUT
self%qmc_shift in self <type(self_type)> = Self-energy
SOURCE
1388 subroutine make_qmcshift_self(cryst_struc,hu,self,apply) 1389 1390 !Arguments ------------------------------------ 1391 !type 1392 type(crystal_t),intent(in) :: cryst_struc 1393 type(hu_type),intent(in) :: hu(cryst_struc%ntypat) 1394 type(self_type),intent(inout) :: self 1395 logical, optional :: apply 1396 1397 !Local variables------------------------------- 1398 integer :: im,iatom,ifreq,itypat,lpawu,tndim 1399 real(dp) :: hu_shift2 1400 character(len=500) :: message 1401 ! ********************************************************************* 1402 1403 do iatom = 1 , cryst_struc%natom 1404 lpawu=self%hdc%matlu(iatom)%lpawu 1405 tndim=2*lpawu+1 1406 itypat=cryst_struc%typat(iatom) 1407 if(lpawu/=-1) then 1408 self%qmc_shift(iatom) = zero 1409 do im =1, 2*tndim-1 1410 ! write(std_out,*)"make before",self%qmc_shift(iatom) 1411 self%qmc_shift(iatom) = self%qmc_shift(iatom) + hu(itypat)%uqmc(im) 1412 ! write(std_out,*)"make after",self%qmc_shift(iatom) 1413 enddo 1414 self%qmc_shift(iatom) = self%qmc_shift(iatom) / two 1415 hu_shift2 = hu(itypat)%uqmc(1) 1416 1417 do im = 2*tndim, 2*tndim + 2*tndim -3 1418 hu_shift2 = hu_shift2 + hu(itypat)%uqmc(im) 1419 enddo 1420 1421 hu_shift2 = hu_shift2 / two 1422 write(message,'(2a,i4)') ch10,' -------> For Correlated atom',iatom 1423 call wrtout(std_out, message,'COLL') 1424 1425 if(abs(self%qmc_shift(iatom)-hu_shift2)>tol6) then 1426 write(message,'(2a,2f16.7)') " Shift for QMC is not correctly"& 1427 & ," computed",self%qmc_shift(iatom),hu_shift2 1428 ABI_ERROR(message) 1429 endif ! shifts not equals 1430 1431 write(message,'(4x,a,f16.7)') & 1432 & " Shift for QMC (used to compute G(w)) is (in Ha) :",& 1433 & self%qmc_shift(iatom) 1434 call wrtout(std_out, message,'COLL') 1435 1436 self%qmc_xmu(iatom)=-self%qmc_shift(iatom) 1437 self%qmc_xmu(iatom)=zero 1438 write(message,'(4x,a,f16.7)') & 1439 & "Artificial Shift used in QMC AND to compute G is (in Ha) :",self%qmc_xmu(iatom) 1440 ABI_WARNING(message) 1441 1442 endif ! lpawu/=1 1443 enddo ! natom 1444 1445 if(present(apply)) then 1446 if (apply) then 1447 write(message,'(5x,a,f16.7,a)') " Shifts applied to self" 1448 call wrtout(std_out, message,'COLL') 1449 do ifreq=1,self%nw 1450 call shift_matlu(self%oper(ifreq)%matlu,cryst_struc%natom,cmplx(self%qmc_shift,0.d0,kind=dp),1) 1451 call shift_matlu(self%oper(ifreq)%matlu,cryst_struc%natom,cmplx(self%qmc_xmu,0.d0,kind=dp),1) 1452 enddo 1453 endif 1454 endif 1455 1456 1457 end subroutine make_qmcshift_self
m_self/new_self [ Functions ]
[ Top ] [ m_self ] [ Functions ]
NAME
new_self
FUNCTION
Mix Old and New self_energy with the mixing coefficient dmft_mxsf
INPUTS
self <type(self_type)>= variables related to self-energy self_new <type(self_type)>= variables related to the new self-energy paw_dmft <type(paw_dmft_type)>= paw+dmft related data opt_mix not used
OUTPUT
self <type(self_type)>= variables related to mixed self-energy
SOURCE
1295 subroutine new_self(self,self_new,paw_dmft,opt_mix) 1296 1297 !Arguments ------------------------------------ 1298 !type 1299 ! type(crystal_t),intent(in) :: cryst_struc 1300 type(self_type),intent(inout) :: self 1301 type(self_type),intent(in) :: self_new 1302 type(paw_dmft_type), intent(in) :: paw_dmft 1303 integer,intent(in) :: opt_mix 1304 1305 !local variables------------------------------- 1306 integer :: iatom,icount,ifreq,im,im1,ispinor,isppol,ispinor1 1307 integer :: natom,ndim,nspinor,nsppol 1308 real(dp) :: alpha,diff_self,sum_self 1309 complex(dpc) :: s1,s2 1310 character(len=500) :: message 1311 ! ********************************************************************* 1312 natom=self%hdc%natom 1313 nsppol=paw_dmft%nsppol 1314 nspinor=paw_dmft%nspinor 1315 alpha=paw_dmft%dmft_mxsf 1316 1317 if(opt_mix==1) then 1318 endif 1319 sum_self=zero 1320 diff_self=zero 1321 icount=0 1322 do iatom=1,natom 1323 if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then 1324 ndim=2*self%oper(1)%matlu(iatom)%lpawu+1 1325 do isppol=1,nsppol 1326 do ispinor=1,nspinor 1327 do ispinor1=1,nspinor 1328 do im=1,ndim 1329 do im1=1,ndim 1330 do ifreq=1,self%nw 1331 ! warning: self_new is the recent self-energy, which is mixed with self 1332 ! to give self= mixed self energy. self_new is deallocated just after. 1333 self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)= & 1334 & (one-alpha)*self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1) + & 1335 & (alpha)*self_new%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1) 1336 s1=self%hdc%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1) 1337 s2=self_new%hdc%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1) 1338 if((ispinor==ispinor1).and.(im==im1)) then 1339 diff_self=diff_self+dsqrt(real(s1-s2)**2+aimag(s1-s2)**2) 1340 sum_self=sum_self+dsqrt(real(s1)**2+aimag(s1)**2) 1341 icount=icount+1 1342 endif 1343 enddo 1344 self%hdc%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)= & 1345 & (one-alpha)*self%hdc%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1) + & 1346 & (alpha)*self_new%hdc%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1) 1347 enddo 1348 enddo 1349 enddo 1350 enddo 1351 enddo ! isppol 1352 endif ! lpawu=/-1 1353 enddo ! iatom 1354 diff_self=diff_self/icount 1355 1356 write(message,'(8x,a,e12.5)') "DMFT Loop: Precision on self-energy is",diff_self 1357 call wrtout(std_out,message,'COLL') 1358 if(diff_self<paw_dmft%dmft_fermi_prec.and.sum_self>tol6.and.paw_dmft%idmftloop>=2) then 1359 write(message,'(a,8x,a,e9.2,a,8x,a)') ch10, "Change of self =<", paw_dmft%dmft_fermi_prec,& 1360 & ch10,"DMFT Loop: Self Energy is converged" 1361 call wrtout(std_out,message,'COLL') 1362 self%iself_cv=1 1363 else 1364 write(message,'(a,8x,a)') ch10,"DMFT Loop: Self Energy is not converged" 1365 call wrtout(std_out,message,'COLL') 1366 self%iself_cv=0 1367 endif 1368 1369 1370 end subroutine new_self
m_self/print_self [ Functions ]
[ Top ] [ m_self ] [ Functions ]
NAME
print_self
FUNCTION
print occupations
INPUTS
self <type(self_type)>= variables related to self-energy option = 1 Do not print double counting. 2 Print double counting paw_dmft <type(paw_dmft_type)> = variables related to self-consistent DFT+DMFT calculations. prtopt = integer which precises the amount of printing in the subroutine called
OUTPUT
self <type(self_type)>= variables related to self-energy
SOURCE
297 subroutine print_self(self,prtdc,paw_dmft,prtopt) 298 299 !Arguments ------------------------------------ 300 !type 301 type(paw_dmft_type), intent(in) :: paw_dmft 302 type(self_type),intent(in) :: self 303 character(len=*), intent(in) :: prtdc 304 integer,intent(in) :: prtopt 305 306 307 !local variables------------------------------- 308 character(len=500) :: message 309 ! ********************************************************************* 310 311 write(message,'(2a)') ch10," == The self-energy for smallest frequency is == " 312 call wrtout(std_out,message,'COLL') 313 call print_oper(self%oper(1),1,paw_dmft,prtopt) 314 ! write(message,'(2a)') ch10," == The self-energy for small (3) frequency is == " 315 ! call wrtout(std_out,message,'COLL') 316 ! call print_oper(self%oper(3),1,paw_dmft,prtopt) 317 write(message,'(2a)') ch10," == The self-energy for large frequency is == " 318 call wrtout(std_out,message,'COLL') 319 call print_oper(self%oper(self%nw),1,paw_dmft,prtopt) 320 if(prtdc=="print_dc") then 321 write(message,'(2a)') ch10," == The double counting hamiltonian is (diagonal part) == " 322 call wrtout(std_out,message,'COLL') 323 call print_matlu(self%hdc%matlu,paw_dmft%natom,prtopt,opt_diag=1) 324 endif 325 326 end subroutine print_self
m_self/rw_self [ Functions ]
[ Top ] [ m_self ] [ Functions ]
NAME
rw_self
FUNCTION
INPUTS
self <type(self_type)>= variables related to self-energy paw_dmft <type(paw_dmft_type)>= paw+dmft related data prtopt = integer which precises the amount of printing opt_rw = 1 Read Self-Energy. 2 Write Self-Energy. 3 Impose Self-Energy.
OUTPUT
SOURCE
458 subroutine rw_self(self,paw_dmft,prtopt,opt_rw,istep_iter,opt_char,opt_imagonly,opt_selflimit,opt_hdc,opt_stop,pawang,cryst_struc) 459 460 !Arguments ------------------------------------ 461 !type 462 type(self_type),intent(inout) :: self 463 type(paw_dmft_type), intent(inout) :: paw_dmft 464 integer,intent(in) :: prtopt 465 integer,intent(in),optional :: opt_rw,istep_iter,opt_imagonly 466 character(len=4), optional :: opt_char 467 integer, intent(in), optional :: opt_stop 468 type(matlu_type), optional, intent(inout) :: opt_selflimit(paw_dmft%natom) 469 type(matlu_type), optional, intent(in) :: opt_hdc(paw_dmft%natom) 470 type(pawang_type), optional, intent(in) :: pawang 471 type(crystal_t), optional, intent(in) :: cryst_struc 472 473 !local variables------------------------------- 474 type(coeff2c_type), allocatable :: eigvectmatlu(:,:) 475 type(matlu_type), allocatable :: level_diag(:) 476 type(matlu_type), allocatable :: selfrotmatlu(:) 477 type(oper_type) :: energy_level 478 logical :: lexist 479 complex(dpc), allocatable :: buffer(:) 480 integer :: iall,iatom,iatu,ier,iexist2,ifreq,im,im1,ioerr,ispinor,ispinor1,isppol,istepiter,istep,istep_imp 481 integer :: icount,iexit,iter,iter_imp,master,mbandc,myproc,natom,ncount,ndim,nkpt,nproc,nrecl,nspinor,nsppol,spacecomm 482 integer :: natom_read,nsppol_read,nspinor_read,ndim_read,nw_read,optrw,readimagonly,tndim,iflavor,unitrot 483 logical :: nondiaglevels 484 character(len=30000) :: message ! Big buffer to avoid buffer overflow. 485 integer,allocatable :: unitselffunc_arr(:) 486 integer,allocatable :: unitselfrot(:,:,:,:) 487 character(len=fnlen) :: tmpfil,tmpfilrot,tmpmatrot 488 character(len=1) :: tag_is 489 character(len=10) :: tag_iflavor 490 character(len=10) :: tag_at 491 character(len=4) :: chtemp 492 real(dp):: xtemp,fermie_read,x_r,x_i 493 real(dp), allocatable:: s_r(:,:,:,:),s_i(:,:,:,:),fermie_read2(:) 494 ! ********************************************************************* 495 496 ! Initialise spaceComm, myproc, and nproc 497 istep=0 ; iter=0 ; istep_imp=0 ; iter_imp=0 498 if(present(opt_rw)) then 499 optrw=opt_rw 500 else 501 optrw=0 502 endif 503 readimagonly=0 504 if(present(opt_imagonly)) then 505 if(opt_imagonly==1.and.paw_dmft%dmft_solv>=5) then 506 readimagonly=opt_imagonly 507 write(message,*) 508 write(message,'(4x,2a)') "About to read imaginary part of Self energy" 509 call wrtout(std_out,message,'COLL') 510 else 511 readimagonly=0 512 write(message,'(4x,2a)') "About to read both real and imaginary part of Self energy" 513 call wrtout(std_out,message,'COLL') 514 endif 515 else 516 readimagonly=0 517 endif 518 if(present(istep_iter)) then 519 istepiter=istep_iter 520 else 521 istepiter=0 522 endif 523 if(paw_dmft%use_fixed_self>0) then 524 istep=istepiter/1000 525 iter=istepiter-(istepiter/1000)*1000 526 istep_imp=paw_dmft%use_fixed_self/1000 527 iter_imp=paw_dmft%use_fixed_self-(paw_dmft%use_fixed_self/1000)*1000 528 endif 529 530 if(paw_dmft%dmft_rslf<=0.and.optrw==1) optrw=0 531 iexit=0 532 ioerr=0 533 iexist2=1 534 lexist=.true. 535 spacecomm=paw_dmft%spacecomm 536 myproc=paw_dmft%myproc 537 nproc=paw_dmft%nproc 538 master=0 539 540 ! write(std_out,*) "myproc,master",myproc,master 541 if(prtopt>200) then 542 endif 543 natom=self%oper(1)%natom 544 nsppol=paw_dmft%nsppol 545 nspinor=paw_dmft%nspinor 546 mbandc=paw_dmft%mbandc 547 nkpt=paw_dmft%nkpt 548 549 ! - For the Tentative rotation of the self-energy file (begin init) 550 if(present(pawang)) then 551 ABI_MALLOC(unitselfrot,(natom,nsppol,nspinor,7)) ! 7 is the max ndim possible 552 if(optrw==2) then 553 write(message,'(a,2x,a,f13.5)') ch10,& 554 & " == About to print self-energy for MAXENT code " 555 else if (optrw==1) then 556 write(message,'(a,2x,a,f13.5)') ch10,& 557 & " == About to read self-energy from MAXENT code " 558 endif 559 call wrtout(std_out,message,'COLL') 560 561 ABI_MALLOC(eigvectmatlu,(natom,nsppol)) 562 do iatom=1,natom 563 if(paw_dmft%lpawu(iatom)/=-1) then 564 tndim=nspinor*(2*paw_dmft%lpawu(iatom)+1) 565 do isppol=1,nsppol 566 ABI_MALLOC(eigvectmatlu(iatom,isppol)%value,(tndim,tndim)) 567 end do 568 end if 569 end do 570 ABI_MALLOC(level_diag,(natom)) 571 ABI_MALLOC(selfrotmatlu,(natom)) 572 call init_matlu(natom,nspinor,nsppol,paw_dmft%lpawu,level_diag) 573 call init_matlu(natom,nspinor,nsppol,paw_dmft%lpawu,selfrotmatlu) 574 call init_oper(paw_dmft,energy_level,opt_ksloc=3) 575 endif 576 ! - For the Tentative rotation of the self-energy file (end init) 577 578 ! - For the Tentative rotation of the self-energy file (begin diag) 579 if(optrw==2.and.present(pawang)) then 580 call compute_levels(cryst_struc,energy_level,self%hdc,pawang,paw_dmft,nondiag=nondiaglevels) 581 write(message,'(a,2x,a,f13.5)') ch10,& 582 & " == Print not Diagonalized Self Energy for Fermi Level=",paw_dmft%fermie 583 call wrtout(std_out,message,'COLL') 584 call print_matlu(self%oper(2)%matlu,natom,1,compl=1,opt_exp=1) 585 call diag_matlu(energy_level%matlu,level_diag,natom,& 586 & prtopt=prtopt,eigvectmatlu=eigvectmatlu,& 587 & test=paw_dmft%dmft_solv) 588 write(message,'(a,2x,a,f13.5)') ch10,& 589 & " == Print Diagonalized levels for Fermi Level=",paw_dmft%fermie 590 call wrtout(std_out,message,'COLL') 591 call print_matlu(level_diag,natom,1,compl=1,opt_exp=1) 592 ! Create file for rotation 593 endif 594 if(present(pawang)) then 595 do iatom=1,natom 596 if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then 597 call int2char4(iatom,tag_at) 598 ABI_CHECK((tag_at(1:1)/='#'),'Bug: string length too short!') 599 if(optrw==2) then 600 tmpmatrot = trim(paw_dmft%filapp)//'.UnitaryMatrix_for_DiagLevel_iatom'//trim(tag_at) 601 else if (optrw==1) then 602 tmpmatrot = trim(paw_dmft%filnamei)//'.UnitaryMatrix_for_DiagLevel_iatom'//trim(tag_at) 603 endif 604 unitrot=3189+iatom 605 #ifdef FC_NAG 606 open (unit=unitrot,file=trim(tmpmatrot),status='unknown',form='formatted',recl=ABI_RECL) 607 #else 608 open (unit=unitrot,file=trim(tmpmatrot),status='unknown',form='formatted') 609 #endif 610 write(std_out,*) " Open file ",trim(tmpmatrot) 611 rewind(unitrot) 612 ndim=2*self%oper(1)%matlu(iatom)%lpawu+1 613 do isppol=1,nsppol 614 if(optrw==2) then 615 do im=1,ndim 616 do im1=1,ndim 617 write(message,*) real(eigvectmatlu(iatom,isppol)%value(im,im1)),aimag(eigvectmatlu(iatom,isppol)%value(im,im1)) 618 call wrtout(unitrot,message,'COLL') 619 enddo 620 enddo 621 else if (optrw==1) then 622 do im=1,ndim 623 do im1=1,ndim 624 read(unitrot,*) x_r,x_i 625 eigvectmatlu(iatom,isppol)%value(im,im1)=cmplx(x_r,x_i) 626 enddo 627 enddo 628 endif 629 enddo 630 close(unitrot) 631 endif 632 enddo 633 if(optrw==1) then 634 write(message,'(a,2x,a,i4)') ch10,& 635 & " == Print non rotated Self Limit read from Matsubara space=",ifreq 636 call wrtout(std_out,message,'COLL') 637 call print_matlu(opt_selflimit,natom,1,compl=1) 638 call rotate_matlu(opt_selflimit,eigvectmatlu,natom,3,1) 639 write(message,'(a,2x,a,i4)') ch10,& 640 & " == Print rotated Self Limit read from Matsubara space=",ifreq 641 call wrtout(std_out,message,'COLL') 642 call print_matlu(opt_selflimit,natom,1,compl=1) 643 endif 644 endif 645 ! - For the Tentative rotation of the self-energy file (end diag) 646 647 if((optrw==2.or.optrw==1).and.myproc==master)then 648 ABI_MALLOC(unitselffunc_arr,(natom*nsppol*nspinor)) 649 iall=0 650 do iatom=1,natom 651 if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then 652 ndim=2*self%oper(1)%matlu(iatom)%lpawu+1 653 ABI_MALLOC(s_r,(ndim,ndim,nspinor,nspinor)) 654 ABI_MALLOC(s_i,(ndim,ndim,nspinor,nspinor)) 655 ! write(std_out,*) "print_self",ndim 656 call int2char4(iatom,tag_at) 657 ABI_CHECK((tag_at(1:1)/='#'),'Bug: string length too short!') 658 do isppol=1,nsppol 659 write(tag_is,'(i1)')isppol 660 ! do ispinor=1,nspinor 661 iall=iall+1 662 ! write(tag_is2,'(i1)')ispinor 663 664 ! =========================== 665 ! == Create name for file 666 ! =========================== 667 if(self%w_type=="real") then 668 if (optrw==1) then 669 tmpfil = trim(paw_dmft%filnamei)//'Self_ra-omega_iatom'//trim(tag_at)//'_isppol'//tag_is 670 else 671 tmpfil = trim(paw_dmft%filapp)//'Self_ra-omega_iatom'//trim(tag_at)//'_isppol'//tag_is 672 endif 673 else 674 tmpfil = trim(paw_dmft%filapp)//'Self-omega_iatom'//trim(tag_at)//'_isppol'//tag_is 675 if(present(opt_char)) then 676 tmpfil = trim(paw_dmft%filapp)//'Self_ra-omega_iatom'//trim(tag_at)//'_isppol'//tag_is//opt_char 677 endif 678 endif 679 if(optrw==1) write(message,'(3a)') ch10," == Read self function and Fermi Level on file ",trim(tmpfil) 680 if(optrw==2) write(message,'(3a)') ch10," == Write self function and Fermi Level on file ",trim(tmpfil) 681 call wrtout(std_out,message,'COLL') 682 !unitselffunc_arr(iall)=300+iall-1 683 unitselffunc_arr(iall) = get_unit() 684 ABI_CHECK(unitselffunc_arr(iall) > 0, "Cannot find free IO unit!") 685 686 !- For the Tentative rotation of the self-energy file (create file) 687 if(optrw==2.and.present(pawang)) then 688 iflavor=0 689 do ispinor=1,nspinor 690 do im=1,ndim 691 iflavor=iflavor+1 692 call int2char4(iflavor,tag_iflavor) 693 unitselfrot(iatom,isppol,ispinor,im)=3000+iflavor 694 ABI_CHECK(unitselfrot(iatom,isppol,ispinor,im) > 0, "Cannot find free IO unit for unitselfrot!") 695 tmpfilrot = trim(paw_dmft%filapp)//'Selfrotformaxent'//& 696 & trim(tag_at)//'_isppol'//tag_is//'_iflavor'//trim(tag_iflavor) 697 write(std_out,*) "Create file ",trim(tmpfilrot)," unit ",unitselfrot(iatom,isppol,ispinor,im)," for flavor",iflavor 698 #ifdef FC_NAG 699 open (unit=unitselfrot(iatom,isppol,ispinor,im),file=trim(tmpfilrot),& 700 & status='unknown',form='formatted',recl=ABI_RECL) 701 #else 702 open (unit=unitselfrot(iatom,isppol,ispinor,im),file=trim(tmpfilrot),status='unknown',form='formatted') 703 #endif 704 rewind(unitselfrot(iatom,isppol,ispinor,im)) 705 enddo 706 enddo 707 endif 708 !- For the Tentative rotation of the self-energy file (create file) 709 710 ! write(std_out,*) "1" 711 712 ! =========================== 713 ! == Read: check that the file exists 714 ! =========================== 715 if(optrw==1) then 716 ! write(std_out,*) "3" 717 inquire(file=trim(tmpfil),exist=lexist,recl=nrecl) 718 if((.not.lexist)) then 719 ! write(std_out,*) "4" 720 iexist2=0 721 write(message,'(4x,a,i5,3a)') "File number",unitselffunc_arr(iall),& 722 & " called ",trim(tmpfil)," does not exist" 723 ! write(std_out,*) lexist,nrecl 724 call wrtout(std_out,message,'COLL') 725 endif 726 endif 727 !write(std_out,*) "2" 728 729 ! =========================== 730 ! == Open file 731 ! =========================== 732 if(optrw==2.or.(optrw==1.and.iexist2==1)) then 733 !write(std_out,*) "5" 734 #ifdef FC_NAG 735 open (unit=unitselffunc_arr(iall),file=trim(tmpfil),status='unknown',form='formatted',recl=ABI_RECL) 736 #else 737 open (unit=unitselffunc_arr(iall),file=trim(tmpfil),status='unknown',form='formatted') 738 #endif 739 rewind(unitselffunc_arr(iall)) 740 !write(std_out,*) "61",nrecl 741 if(prtopt>=3) then 742 write(message,'(a,a,a,i4)') ' opened file : ', trim(tmpfil), ' unit', unitselffunc_arr(iall) 743 call wrtout(std_out,message,'COLL') 744 endif 745 endif 746 !write(std_out,*) "6",nrecl 747 748 ! =========================== 749 ! == Check Header 750 ! =========================== 751 if(optrw==2) then 752 write(message,'(3a,5i5,2x,e25.17)') "# natom,nsppol,nspinor,ndim,nw,fermilevel",ch10& 753 & ,"####",natom,nsppol,nspinor,ndim,self%nw,paw_dmft%fermie 754 call wrtout(unitselffunc_arr(iall),message,'COLL') 755 else if(optrw==1.and.iexist2==1.and.readimagonly==0) then 756 read(unitselffunc_arr(iall),*) 757 read(unitselffunc_arr(iall),*,iostat=ioerr)& 758 & chtemp,natom_read,nsppol_read,nspinor_read,ndim_read,nw_read,fermie_read 759 if(ioerr<0) then 760 ! write(std_out,*)" HEADER IOERR" 761 ! write(std_out,'(a4,2x,31(e15.8,2x))') chtemp,natom_read,nsppol_read,nspinor_read,ndim_read,nw_read,fermie_read 762 endif 763 if(ioerr==0) then 764 write(message,'(a,3x,3a,i12,2a,i11,2a,i10,2a,i13,2a,i11,2a,e25.8)') ch10,"Data in Self file corresponds to",& 765 & ch10," natom",natom_read,& 766 & ch10," nsppol",nsppol_read,& 767 & ch10," nspinor",nspinor_read,& 768 & ch10," ndim",ndim_read, & 769 & ch10," nw",nw_read, & 770 & ch10," Fermi level",fermie_read 771 call wrtout(std_out,message,'COLL') 772 if((natom/=natom_read).or.(nsppol_read/=nsppol).or.& 773 & (nspinor/=nspinor_read).or.(nw_read/=self%nw)) then 774 write(message,'(a,3x,3a,i12,2a,i11,2a,i10,2a,i13,2a,i11,2a,e25.8)') ch10,"Data required is ",& 775 & ch10," natom",natom,& 776 & ch10," nsppol",nsppol,& 777 & ch10," nspinor",nspinor,& 778 & ch10," ndim",ndim, & 779 & ch10," nw",self%nw, & 780 & ch10," Fermi level",paw_dmft%fermie 781 call wrtout(std_out,message,'COLL') 782 message = "Dimensions in self are not correct" 783 if(readimagonly==1.or.present(opt_stop)) then 784 ABI_ERROR(message) 785 else 786 ABI_WARNING(message) 787 endif 788 iexist2=2 789 endif 790 else 791 ABI_WARNING("Self file is empty") 792 endif 793 endif 794 !write(std_out,*) "7" 795 796 ! =========================== 797 ! == Write/Read self in the file 798 ! =========================== 799 800 rewind(111) 801 do ifreq=1,self%nw 802 if(optrw==2) then 803 ! write(std_out,'(a,2x,31(e15.8,2x))') & 804 !& "SETEST",self%omega(ifreq),& 805 !& (self%oper(ifreq)%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)& 806 !& ,im=1,ndim) 807 ! write(std_out,*) self%omega(ifreq),& 808 !& ((self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor)& 809 !& ,im=1,ndim),im1=1,ndim) 810 ! write(message,'(2x,393(e25.17,2x))') self%omega(ifreq),& 811 !& ((((self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)& 812 !& ,im=1,ndim),im1=1,ndim),ispinor=1,nspinor),ispinor1=1,nspinor) 813 814 !MGNAG Runtime Error: wrtout_cpp.f90, line 896: Buffer overflow on output 815 !Is it possible to rewrite the code below to avoid such a long message 816 !What about Netcdf binary files ? 817 if(nspinor==1) then 818 write(message,'(2x,393(e25.17,2x))') self%omega(ifreq),& 819 & ((((real(self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)),& 820 & aimag(self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)),& 821 & im=1,ndim),im1=1,ndim),ispinor=1,nspinor),ispinor1=1,nspinor) 822 else 823 write(message,'(2x,393(e18.10,2x))') self%omega(ifreq),& 824 & ((((real(self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)),& 825 & aimag(self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)),& 826 & im=1,ndim),im1=1,ndim),ispinor=1,nspinor),ispinor1=1,nspinor) 827 endif 828 call wrtout(unitselffunc_arr(iall),message,'COLL') 829 830 !- For the Tentative rotation of the self-energy file (begin rot) 831 !---------------------------------------------------------------- 832 if(optrw==2.and.present(pawang)) then 833 call copy_matlu(self%oper(ifreq)%matlu,selfrotmatlu,natom) 834 if(ifreq<3) then 835 write(message,'(a,2x,a,i4)') ch10,& 836 & " == Print non Rotated Self Energy for freq=",ifreq 837 call wrtout(std_out,message,'COLL') 838 call print_matlu(selfrotmatlu,natom,1,compl=1) 839 endif 840 call rotate_matlu(selfrotmatlu,eigvectmatlu,natom,3,1) 841 if(ifreq<3) then 842 write(message,'(a,2x,a,i4)') ch10,& 843 & " == Print Rotated Self Energy for freq=",ifreq 844 call wrtout(std_out,message,'COLL') 845 call print_matlu(selfrotmatlu,natom,1,compl=1) 846 else if(ifreq==3) then 847 write(message,'(a,2x,a,i4)') ch10,& 848 & " (Other frequencies not printed)" 849 call wrtout(std_out,message,'COLL') 850 endif 851 !write(message,'(2x,393(e18.10,2x))') self%omega(ifreq),& 852 ! ((real(selfrotmatlu(iatom)%mat(im,im,isppol,ispinor,ispinor)),& 853 ! aimag(selfrotmatlu(iatom)%mat(im,im,isppol,ispinor,ispinor)),& 854 ! im=1,ndim),ispinor=1,nspinor) 855 iflavor=0 856 do ispinor=1,nspinor 857 do im=1,ndim 858 iflavor=iflavor+1 859 ! if(ifreq<5) then 860 ! write(std_out,*) "Write in file unit",unitselfrot(iatom,isppol,ispinor,im),"for flavor",iflavor 861 ! endif 862 write(message,'(2x,393(e18.10,2x))') self%omega(ifreq),& 863 & real(selfrotmatlu(iatom)%mat(im,im,isppol,ispinor,ispinor)),& 864 & aimag(selfrotmatlu(iatom)%mat(im,im,isppol,ispinor,ispinor)) 865 ! write(6,'(2x,393(e18.10,2x))') self%omega(ifreq),& 866 !& real(selfrotmatlu(iatom)%mat(im,im,isppol,ispinor,ispinor)),& 867 !& aimag(selfrotmatlu(iatom)%mat(im,im,isppol,ispinor,ispinor)) 868 ! if(iflavor==1) then 869 ! write(1024,*) iatom,isppol,ispinor,im,unitselfrot(iatom,isppol,ispinor,im) 870 ! write(1024,'(2x,393(e18.10,2x))') self%omega(ifreq),& 871 ! & real(selfrotmatlu(iatom)%mat(im,im,isppol,ispinor,ispinor)),& 872 ! & aimag(selfrotmatlu(iatom)%mat(im,im,isppol,ispinor,ispinor)) 873 ! endif 874 call wrtout(unitselfrot(iatom,isppol,ispinor,im),message,'COLL') 875 enddo 876 enddo 877 endif 878 !- For the Tentative rotation of the self-energy file (end rot) 879 880 ! write(std_out,*) unitselffunc_arr(iall) 881 else if(optrw==1.and.iexist2==1.and.ioerr==0) then 882 !write(std_out,*) "8" 883 ! read(unitselffunc_arr(iall),'(2x,31(e15.8,2x))',iostat=ioerr) & 884 !& xtemp,(s_r(im),s_i(im),im=1,ndim) 885 if(readimagonly==0) then 886 read(unitselffunc_arr(iall),*,iostat=ioerr) xtemp,& 887 & ((((s_r(im,im1,ispinor,ispinor1),s_i(im,im1,ispinor,ispinor1),im=1,ndim),im1=1,ndim),ispinor=1,nspinor),ispinor1=1,nspinor) 888 ! if(ioerr<0) then 889 ! write(std_out,*)" SELF IOERR<" 890 ! else if(ioerr>0) then 891 ! write(std_out,*)" SELF IOERR>" 892 ! write(std_out,'(a4,2x,31(e15.8,2x))') xtemp,(s_r(im),s_i(im),im=1,ndim) 893 ! endif 894 do im=1,ndim 895 do im1=1,ndim 896 do ispinor=1,nspinor 897 do ispinor1=1,nspinor 898 self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)& 899 & =cmplx(s_r(im,im1,ispinor,ispinor1),s_i(im,im1,ispinor,ispinor1),kind=dp) 900 !if((im1==im1).and.(ispinor==ispinor1)) then 901 ! write(6,*) "Self read",ifreq, self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1) 902 !endif 903 enddo 904 enddo 905 enddo 906 enddo 907 endif 908 endif 909 enddo ! ifreq 910 911 !- For the Tentative rotation of the self-energy file (begin close file) 912 if(optrw==2.and.present(pawang)) then 913 do ispinor=1,nspinor 914 do im=1,ndim 915 close(unitselfrot(iatom,isppol,ispinor,im)) 916 write(std_out,*) "Close file unit",unitselfrot(iatom,isppol,ispinor,im) 917 enddo 918 enddo 919 endif 920 !- For the Tentative rotation of the self-energy file (end close file) 921 922 923 if(optrw==1.and.iexist2==1.and.ioerr==0) then 924 if(readimagonly==1) then ! read from OmegaMaxent 925 s_r=zero 926 927 ! Read self energy from Maxent (imag part) on the real axis 928 !---------------------------------------------------------- 929 do ifreq=1,self%nw 930 call zero_matlu(self%oper(ifreq)%matlu,natom) 931 enddo 932 do ispinor=1,nspinor 933 do im=1,ndim 934 do ifreq=1,self%nw 935 read(unitselffunc_arr(iall),*,iostat=ioerr) xtemp,s_i(im,im,ispinor,ispinor) 936 ! minus sign because - Im Sigma is the output of OmegaMaxent 937 self%oper(ifreq)%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)& 938 & =cmplx(s_r(im,im,ispinor,ispinor),-s_i(im,im,ispinor,ispinor),kind=dp) 939 enddo ! ifreq 940 enddo 941 enddo 942 943 ! Kramers Kronig 944 !------------------- 945 write(message,'(4x,2a)') " Read only diagonal self energy from Maxent" 946 call wrtout(std_out,message,'COLL') 947 !write(6,*) "opt_hdc",opt_hdc(1)%mat(1,1,1,1,1) 948 call kramerskronig_self(self,opt_selflimit,opt_hdc,paw_dmft%filapp) 949 950 ! Rotate back rotate_matlu 951 !----------------------------- 952 if(present(pawang)) then 953 write(message,'(4x,2a)') " Rotate Back self in the original basis" 954 call wrtout(std_out,message,'COLL') 955 do ifreq=1,self%nw 956 if(ifreq<20) then 957 write(message,'(a,2x,a,i4)') ch10,& 958 & " == Print Rotated real axis Self Energy for freq=",ifreq 959 call wrtout(std_out,message,'COLL') 960 call print_matlu(self%oper(ifreq)%matlu,natom,1,compl=1) 961 endif 962 call rotate_matlu(self%oper(ifreq)%matlu,eigvectmatlu,natom,3,-1) 963 if(ifreq<20) then 964 write(message,'(a,2x,a,i4)') ch10,& 965 & " == Print Rotated back real axis Self Energy for freq=",ifreq 966 call wrtout(std_out,message,'COLL') 967 call print_matlu(self%oper(ifreq)%matlu,natom,1,compl=1) 968 endif 969 enddo 970 endif 971 972 endif 973 endif 974 975 ! =========================== 976 ! == Write/Read hdc in the file 977 ! =========================== 978 if(optrw==2) then 979 ! write(std_out,'(a,2x,31(e15.8,2x))') & 980 !& "SETEST #dc ",(self%hdc%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor),im=1,ndim) 981 write(message,'(a,2x,31(e25.17,2x))') & 982 & "#dc ",((self%hdc%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor),im=1,ndim),ispinor=1,nspinor) 983 call wrtout(unitselffunc_arr(iall),message,'COLL') 984 else if(optrw==1.and.iexist2==1.and.ioerr==0.and.readimagonly==0) then 985 !write(std_out,*) "8" 986 read(unitselffunc_arr(iall),*,iostat=ioerr) & 987 & chtemp,((s_r(im,1,ispinor,1),s_i(im,1,ispinor,1),im=1,ndim),ispinor=1,nspinor) 988 if(ioerr<0) then 989 ! write(std_out,*)" HDC IOERR<",ioerr 990 else if(ioerr>0) then 991 ! write(std_out,*)" HDC IOERR>",ioerr 992 endif 993 do ispinor=1,nspinor 994 do im=1,ndim 995 self%hdc%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)& 996 & =cmplx(s_r(im,1,ispinor,1),s_i(im,1,ispinor,1),kind=dp) 997 enddo 998 enddo 999 !write(6,*) "read selfhdc",self%hdc%matlu(1)%mat(1,1,1,1,1) 1000 else if(readimagonly==1.and..not.present(opt_hdc)) then 1001 do ispinor=1,nspinor 1002 do im=1,ndim 1003 self%hdc%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)=czero 1004 enddo 1005 enddo 1006 else 1007 write(std_out,*) " self%hdc fixed in kramerskronig_self" 1008 endif 1009 close(unitselffunc_arr(iall)) 1010 ! enddo ! ispinor 1011 enddo ! isppol 1012 ABI_FREE(s_r) 1013 ABI_FREE(s_i) 1014 endif ! lpawu=/-1 1015 enddo ! iatom 1016 ABI_FREE(unitselffunc_arr) 1017 endif ! optrw==2.or.myproc==master 1018 ! call xmpi_barrier(spacecomm) 1019 !write(std_out,*) "9" 1020 ! - For the Tentative rotation of the self-energy file (begin destroy) 1021 if(present(pawang)) then 1022 call destroy_oper(energy_level) 1023 call destroy_matlu(level_diag,natom) 1024 call destroy_matlu(selfrotmatlu,natom) 1025 ABI_FREE(level_diag) 1026 ABI_FREE(selfrotmatlu) 1027 do iatom=1,natom 1028 if(paw_dmft%lpawu(iatom)/=-1) then 1029 do isppol=1,nsppol 1030 ABI_FREE(eigvectmatlu(iatom,isppol)%value) 1031 end do 1032 end if 1033 end do 1034 ABI_FREE(eigvectmatlu) 1035 ABI_FREE(unitselfrot) ! 7 is the max ndim possible 1036 endif 1037 ! - For the Tentative rotation of the self-energy file (end destroy) 1038 1039 ! =========================== 1040 ! == Error messages 1041 ! =========================== 1042 if(optrw==1) then 1043 ! call xmpi_barrier(spacecomm) 1044 ncount=natom*nsppol*(nspinor**2)*(self%nw+1)*(maxval(self%oper(1)%matlu(:)%lpawu)*2+1)**2 1045 !write(std_out,*) ncount,maxval(pawtab(:)%lpawu)*2+1 1046 call xmpi_bcast(iexist2,master,spacecomm ,ier) 1047 call xmpi_bcast(ioerr,master,spacecomm ,ier) 1048 if(iexist2==0.or.ioerr<0.or.ioerr>0) then 1049 message = "Self file does not exist or is incomplete" 1050 if(readimagonly==1.or.present(opt_stop)) then 1051 if(readimagonly==1) then 1052 message = "Self file does not exist or is incomplete: check the number of self data in file" 1053 endif 1054 ABI_ERROR(message) 1055 else 1056 ABI_WARNING(message) 1057 endif 1058 if(iexist2==0) then 1059 write(message,'(4x,2a)') "File does not exist" 1060 call wrtout(std_out,message,'COLL') 1061 endif 1062 if(ioerr<0) then 1063 write(message,'(4x,2a)') "End of file reached" 1064 call wrtout(std_out,message,'COLL') 1065 endif 1066 if(ioerr>0) then 1067 write(message,'(4x,2a)') "Error during read statement" 1068 call wrtout(std_out,message,'COLL') 1069 endif 1070 if(paw_dmft%dmft_solv/=4) then 1071 write(message,'(4x,a,a,5i5,2x,e14.7)') "-> Put Self-Energy Equal to double counting term" 1072 else if(paw_dmft%dmft_solv==4) then 1073 write(message,'(4x,a,a,5i5,2x,e14.7)') "-> Put Self-Energy Equal to dc term - shift" 1074 call wrtout(std_out,message,'COLL') 1075 write(message,'(4x,a,a,5i5,2x,e14.7)') " No self energy is given, change dmft_rslf" 1076 ABI_ERROR(message) 1077 endif 1078 call wrtout(std_out,message,'COLL') 1079 do ifreq=1,self%nw 1080 ! write(std_out,*) "before",self%oper(1)%matlu(1)%mat(1,1,1,1,1) 1081 ! write(std_out,*) "before",self%hdc%matlu(1)%mat(1,1,1,1,1) 1082 call copy_matlu(self%hdc%matlu,self%oper(ifreq)%matlu,natom) 1083 ! write(std_out,*) "after",self%oper(1)%matlu(1)%mat(1,1,1,1,1) 1084 ! write(std_out,*) "before",self%hdc%matlu(1)%mat(1,1,1,1,1) 1085 if(paw_dmft%dmft_solv==4) then 1086 ! if(ifreq==1) write(std_out,*) "shift",self%qmc_shift(1) 1087 call shift_matlu(self%oper(ifreq)%matlu,natom,cmplx(self%qmc_shift,0.d0,kind=dp),1) 1088 ! if(ifreq==1) write(std_out,*) "self after dc and shift",self%oper(ifreq)%matlu(1)%mat(1,1,1,1,1) 1089 ! if(ifreq==1) write(std_out,*) "shift",self%qmc_shift(1) 1090 endif 1091 enddo 1092 else ! test read successfull 1093 ! call xmpi_barrier(spacecomm) 1094 !! lignes 924-928 semblent inutiles puisque la valeur de paw_dmft%fermie creee 1095 !! en ligne 927 est ecrasee en ligne 992. BA+jmb 1096 !! ABI_MALLOC(fermie_read2,(1)) 1097 !! fermie_read2(1)=fermie_read 1098 !! call xmpi_sum(fermie_read2,spacecomm ,ier) 1099 !! paw_dmft%fermie=fermie_read2(1) 1100 !! ABI_FREE(fermie_read2) 1101 1102 ! =========================== 1103 ! bcast to other proc 1104 ! =========================== 1105 ABI_MALLOC(buffer,(ncount)) 1106 ABI_MALLOC(fermie_read2,(1)) 1107 buffer(:)=czero 1108 !! BA+jmb 1109 fermie_read2=zero 1110 !write(std_out,*) self%nw 1111 if(myproc==master) then 1112 1113 ! == Send read data to all process 1114 icount=0 1115 fermie_read2(1)=fermie_read 1116 do iatom=1,natom 1117 if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then 1118 ndim=2*self%oper(1)%matlu(iatom)%lpawu+1 1119 do isppol=1,nsppol 1120 do ispinor=1,nspinor 1121 ! Self energy----------- 1122 do ifreq=1,self%nw 1123 do ispinor1=1,nspinor 1124 do im=1,ndim 1125 do im1=1,ndim 1126 icount=icount+1 1127 if(icount.gt.ncount) then 1128 write(message,'(2a,2i5)') ch10,"Error buffer",icount,ncount 1129 iexit=1 1130 ABI_ERROR(message) 1131 endif 1132 buffer(icount)=self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1) 1133 enddo 1134 enddo 1135 enddo ! ispinor1 1136 enddo 1137 ! Double counting------- 1138 do im=1,ndim 1139 icount=icount+1 1140 if(icount.gt.ncount) then 1141 write(message,'(2a,2i5)') ch10,"Error buffer",icount,ncount 1142 iexit=1 1143 ABI_ERROR(message) 1144 endif 1145 buffer(icount)=self%hdc%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor) 1146 enddo 1147 enddo ! ispinor 1148 enddo ! isppol 1149 endif ! lpawu=/-1 1150 enddo ! iatom 1151 endif 1152 ! call xmpi_bcast(buffer,master,spacecomm ,ier) 1153 ! call xmpi_sum(iexit,spacecomm ,ier) 1154 !!JB call xmpi_barrier(spacecomm) 1155 call xmpi_sum(buffer,spacecomm ,ier) 1156 !!JB call xmpi_barrier(spacecomm) 1157 1158 ! bcast fermi level 1159 call xmpi_sum(fermie_read2,spacecomm ,ier) 1160 1161 if(ier/=0) then 1162 message = "error in xmpi_sum in rw_self" 1163 ABI_ERROR(message) 1164 endif 1165 paw_dmft%fermie=fermie_read2(1) 1166 ! write(std_out,*) "Fermi level",paw_dmft%fermie 1167 icount=0 1168 do iatom=1,natom 1169 if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then 1170 ndim=2*self%oper(1)%matlu(iatom)%lpawu+1 1171 do isppol=1,nsppol 1172 do ispinor=1,nspinor 1173 ! self --------------- 1174 do ifreq=1,self%nw 1175 do ispinor1=1,nspinor 1176 do im=1,ndim 1177 do im1=1,ndim 1178 icount=icount+1 1179 self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)=buffer(icount) 1180 !write(6,*)'self procs', ifreq, self%oper(ifreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1) 1181 enddo 1182 enddo 1183 enddo 1184 enddo 1185 ! hdc --------------- 1186 do im=1,ndim 1187 icount=icount+1 1188 self%hdc%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)=buffer(icount) 1189 enddo 1190 enddo 1191 enddo ! isppol 1192 endif ! lpawu=/-1 1193 enddo ! iatom 1194 ABI_FREE(fermie_read2) 1195 ABI_FREE(buffer) 1196 endif ! test read successful 1197 endif ! optrw==1 1198 ! call flush_unit(std_out) 1199 ! ABI_ERROR("Aboring now") 1200 if(optrw==0) then 1201 if(paw_dmft%dmft_rslf==0) then 1202 if(paw_dmft%dmft_solv/=4) then 1203 write(message,'(4x,a,a,5i5,2x,e14.7)') "-> Put Self-Energy Equal to double counting term" 1204 else if(paw_dmft%dmft_solv==4) then 1205 write(message,'(4x,a,a,5i5,2x,e14.7)') "-> Put Self-Energy Equal to dc term - shift" 1206 endif 1207 else if (paw_dmft%dmft_rslf==-1) then 1208 write(message,'(4x,a,a,5i5,2x,e14.7)') "-> Put Self-Energy Equal to zero" 1209 endif 1210 call wrtout(std_out,message,'COLL') 1211 do ifreq=1,self%nw 1212 if(paw_dmft%dmft_rslf==0) then 1213 call copy_matlu(self%hdc%matlu,self%oper(ifreq)%matlu,natom) 1214 ! if(ifreq==1) write(std_out,*) "self after dc",self%oper(ifreq)%matlu(1)%mat(1,1,1,1,1) 1215 if(paw_dmft%dmft_solv==4) then 1216 ! if(ifreq==1) write(std_out,*) "shift",self%qmc_shift(1) 1217 call shift_matlu(self%oper(ifreq)%matlu,natom,cmplx(self%qmc_shift,0.d0,kind=dp),1) 1218 ! if(ifreq==1) write(std_out,*) "self after dc and shift",self%oper(ifreq)%matlu(1)%mat(1,1,1,1,1) 1219 ! if(ifreq==1) write(std_out,*) "shift",self%qmc_shift(1) 1220 endif 1221 else if (paw_dmft%dmft_rslf==-1) then 1222 do iatom=1,natom 1223 if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then 1224 ndim=2*self%oper(1)%matlu(iatom)%lpawu+1 1225 do isppol=1,nsppol 1226 do ispinor=1,nspinor 1227 do im=1,ndim 1228 self%oper(ifreq)%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)= czero 1229 enddo 1230 enddo 1231 enddo 1232 endif 1233 enddo 1234 endif 1235 enddo 1236 endif 1237 1238 ! write(std_out,*) "optrw,use_fixed_self,istep,iter,istep_imp,iter_imp" 1239 ! write(std_out,*) optrw,paw_dmft%use_fixed_self,istep,iter,istep_imp,iter_imp 1240 if((optrw==1.or.optrw==3).and.paw_dmft%use_fixed_self>0.and.istep<=istep_imp.and.iter<=iter_imp) then 1241 write(message,'(4x,a)') "-> Put Self-Energy Equal to imposed self-energy" 1242 call wrtout(std_out,message,'COLL') 1243 iatu=0 1244 do iatom=1,natom 1245 if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then 1246 iatu=iatu+1 1247 do ifreq=1,self%nw 1248 ndim=2*self%oper(1)%matlu(iatom)%lpawu+1 1249 do isppol=1,nsppol 1250 do ispinor=1,nspinor 1251 do im=1,ndim 1252 if(nspinor==1) then 1253 self%oper(ifreq)%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)= paw_dmft%fixed_self(im,im,isppol,iatu) 1254 ! write(std_out,*) paw_dmft%fixed_self(im,im,isppol,iatu) 1255 else 1256 self%oper(ifreq)%matlu(iatom)%mat(im,im,isppol,ispinor,ispinor)= paw_dmft%fixed_self(im,im,ispinor,iatu) 1257 ! write(message,'(a,i4,i4,2x,e20.10)') " Fixed self not implemented for nspinor==2" 1258 ! call wrtout(std_out, message,'COLL') 1259 ! ABI_ERROR("Aboring now") 1260 endif 1261 enddo 1262 enddo 1263 enddo 1264 enddo 1265 endif 1266 enddo 1267 1268 endif 1269 1270 1271 1272 1273 end subroutine rw_self
m_self/self_type [ Types ]
NAME
self_type
FUNCTION
This structured datatype contains the necessary data
SOURCE
64 type, public :: self_type ! for each atom 65 66 integer :: dmft_nwlo 67 ! dmft frequencies 68 69 character(len=4) :: w_type 70 ! type of frequencies used 71 72 integer :: nw 73 ! dmft frequencies 74 75 integer :: iself_cv 76 ! integer for self convergence 77 78 integer :: dmft_nwli 79 ! dmft frequencies 80 81 real(dp), pointer :: omega(:) => null() 82 ! value of frequencies 83 84 real(dp), allocatable :: qmc_shift(:) 85 ! value of frequencies 86 87 real(dp), allocatable :: qmc_xmu(:) 88 ! value of frequencies 89 90 type(oper_type), allocatable :: oper(:) 91 ! self function in different basis 92 93 type(oper_type):: hdc 94 ! self function in different basis 95 96 end type self_type
m_self/selfreal2imag_self [ Functions ]
[ Top ] [ m_self ] [ Functions ]
NAME
selfreal2imag_self
FUNCTION
INPUTS
hu <type(hu_type)> = U interaction paw_dmft <type(paw_dmft_type)> = paw+dmft related data
OUTPUT
self%qmc_shift in self <type(self_type)> = Self-energy
SOURCE
1623 subroutine selfreal2imag_self(selfr,self,filapp) 1624 1625 !Arguments ------------------------------------ 1626 !type 1627 type(self_type),intent(inout) :: selfr 1628 type(self_type),intent(inout) :: self 1629 character(len=fnlen), intent(in) :: filapp 1630 1631 !Local variables------------------------------- 1632 integer :: ifreq,jfreq,isppol,ispinor,ispinor1,im,im1,iatom 1633 complex(dpc), allocatable :: selftempmatsub(:) 1634 integer :: natom,ndim,nsppol,nspinor 1635 real(dp) :: delta 1636 ! ********************************************************************* 1637 delta=0.0000000 1638 ABI_MALLOC(selftempmatsub,(self%nw)) 1639 natom=self%hdc%natom 1640 nsppol = self%hdc%nsppol 1641 nspinor=self%hdc%nspinor 1642 ! Compute limit of Real Part and put in double counting energy. 1643 ! call copy_matlu(selfhdc,self%hdc%matlu,natom) 1644 !write(6,*) "self3",aimag(selfr%oper(489)%matlu(1)%mat(1,1,1,1,1)) 1645 1646 open(unit=672,file=trim(filapp)//"_DFTDMFT_Self_forcheck_imagaxis_from_realaxis.dat", status='unknown',form='formatted') 1647 do iatom=1,natom 1648 if(self%oper(1)%matlu(iatom)%lpawu.ne.-1) then 1649 ndim=2*self%oper(1)%matlu(iatom)%lpawu+1 1650 do isppol=1,nsppol 1651 do ispinor=1,nspinor 1652 do ispinor1=1,nspinor 1653 do im=1,ndim 1654 do im1=1,ndim 1655 do jfreq=1,selfr%nw-1 1656 ! write(6700,*) selfr%omega(jfreq),aimag(selfr%oper(jfreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)) 1657 enddo 1658 ! write(6700,*) 1659 do ifreq=1,self%nw 1660 selftempmatsub(ifreq)=czero 1661 do jfreq=1,selfr%nw-1 1662 selftempmatsub(ifreq)=selftempmatsub(ifreq) - & 1663 & aimag(selfr%oper(jfreq)%matlu(iatom)%mat(im,im1,isppol,ispinor,ispinor1)) & 1664 & /(cmplx(zero,self%omega(ifreq),kind=dp)-selfr%omega(jfreq)) & 1665 & * (selfr%omega(jfreq+1)-selfr%omega(jfreq)) 1666 enddo 1667 selftempmatsub(ifreq)=selftempmatsub(ifreq)/pi 1668 write(672,*) self%omega(ifreq),real(selftempmatsub(ifreq)),aimag(selftempmatsub(ifreq)) 1669 enddo 1670 write(672,*) 1671 enddo 1672 enddo 1673 enddo 1674 enddo 1675 enddo ! isppol 1676 endif ! lpawu=/-1 1677 enddo ! iatom 1678 close(672) 1679 1680 ABI_FREE(selftempmatsub) 1681 1682 1683 end subroutine selfreal2imag_self