TABLE OF CONTENTS
- ABINIT/m_BathOperatoroffdiag
- ABINIT/m_BathOperatoroffdiag/ BathOperatoroffdiag_destroy
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_activateParticle
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_checkM
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_doCheck
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_getDetAdd
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_getDetF
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_getDetRemove
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_getError
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_init
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_initF
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_printF
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_printM
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_printM_matrix
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_recomputeM
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_reset
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_setF
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_setMAdd
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_setMRemove
- ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_swap
- m_BathOperatoroffdiag/BathOperatoroffdiag
ABINIT/m_BathOperatoroffdiag [ Modules ]
NAME
m_BathOperatoroffdiag
FUNCTION
Manage all stuff related to the bath for the simgle Anderson Impurity Model
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder, J. Denier, B. Amadon) 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
SOURCE
24 #include "defs.h" 25 MODULE m_BathOperatoroffdiag 26 USE m_MatrixHyb 27 USE m_Vector 28 USE m_VectorInt 29 USE m_Global 30 USE m_ListCdagC 31 IMPLICIT NONE 32 33 ! subroutines 34 public :: BathOperatoroffdiag_init 35 public :: BathOperatoroffdiag_reset 36 public :: BathOperatoroffdiag_activateParticle 37 public :: BathOperatoroffdiag_setMAdd 38 public :: BathOperatoroffdiag_setMRemove 39 public :: BathOperatoroffdiag_swap 40 public :: BathOperatoroffdiag_initF 41 public :: BathOperatoroffdiag_setF 42 public :: BathOperatoroffdiag_printF 43 public :: BathOperatoroffdiag_printM 44 public :: BathOperatoroffdiag_destroy 45 public :: BathOperatoroffdiag_doCheck 46 public :: BathOperatoroffdiag_checkM 47 48 ! functions 49 ! public :: BathOperatoroffdiag_hybrid 50 public :: BathOperatoroffdiag_getDetAdd 51 public :: BathOperatoroffdiag_getDetRemove 52 public :: BathOperatoroffdiag_getDetF 53 public :: BathOperatoroffdiag_getError
ABINIT/m_BathOperatoroffdiag/ BathOperatoroffdiag_destroy [ Functions ]
NAME
BathOperatoroffdiag_destroy
FUNCTION
Deallocate and reset every thing
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1860 SUBROUTINE BathOperatoroffdiag_destroy(op) 1861 1862 TYPE(BathOperatoroffdiag), INTENT(INOUT) :: op 1863 1864 CALL MatrixHyb_destroy(op%M) 1865 CALL MatrixHyb_destroy(op%M_update) 1866 1867 CALL Vector_destroy(op%R) 1868 CALL Vector_destroy(op%Q) 1869 CALL Vector_destroy(op%Rtau) 1870 CALL Vector_destroy(op%Qtau) 1871 FREEIF(op%F) 1872 FREEIF(op%Fshift) 1873 FREEIF(op%tails) 1874 1875 op%MAddFlag = .FALSE. 1876 op%MRemoveFlag = .FALSE. 1877 op%flavors = 0 1878 op%beta = 0.d0 1879 op%dt = 0.d0 1880 op%inv_dt = 0.d0 1881 op%samples = 0 1882 op%sizeHybrid = 0 1883 op%activeFlavor = 0 1884 op%updatePosRow = 0 1885 op%updatePosCol = 0 1886 1887 END SUBROUTINE BathOperatoroffdiag_destroy
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_activateParticle [ Functions ]
NAME
BathOperatoroffdiag_activateParticle
FUNCTION
Just save on wicht flavor we are working It is better to use the macro defined in defs.h
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator flavor=the flavor to activate
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
371 SUBROUTINE BathOperatoroffdiag_activateParticle(op,flavor) 372 373 !Arguments ------------------------------------ 374 TYPE(BathOperatoroffdiag), INTENT(INOUT) :: op 375 !Local variables ------------------------------ 376 INTEGER , INTENT(IN ) :: flavor 377 378 IF ( flavor .GT. op%flavors ) & 379 CALL ERROR("BathOperatoroffdiag_activateParticle : out of range ") 380 IF ( op%set .EQV. .TRUE. ) THEN 381 op%activeFlavor = flavor 382 op%MAddFlag = .FALSE. 383 op%MRemoveFlag = .FALSE. 384 ELSE 385 CALL ERROR("BathOperatoroffdiag_activateParticle : not allocated ") 386 END IF 387 END SUBROUTINE BathOperatoroffdiag_activateParticle
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_checkM [ Functions ]
NAME
BathOperatoroffdiag_checkM
FUNCTION
compute from scratch the M matrix and compar it with the already computed M matrix
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator particle=list of all segments of the active flavor
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1952 SUBROUTINE BathOperatoroffdiag_checkM(op,particle) 1953 1954 !Arguments ------------------------------------ 1955 TYPE(BathOperatoroffdiag) , INTENT(INOUT) :: op 1956 TYPE(ListCdagC) , INTENT(IN ) :: particle(:) 1957 !Local variables ------------------------------ 1958 ! TYPE(MatrixHyb) :: checkMatrix 1959 LOGICAL :: checkTau 1960 INTEGER :: tail 1961 INTEGER :: iC 1962 INTEGER :: iCdag 1963 INTEGER :: aF 1964 INTEGER :: iflavora 1965 INTEGER :: iflavorb,it !,it1 1966 CHARACTER(LEN=6) :: a 1967 DOUBLE PRECISION :: time 1968 DOUBLE PRECISION :: beta 1969 DOUBLE PRECISION :: mbeta_two 1970 DOUBLE PRECISION :: errorabs 1971 DOUBLE PRECISION :: errormax 1972 DOUBLE PRECISION :: error1 1973 DOUBLE PRECISION :: errorrel 1974 DOUBLE PRECISION :: tc 1975 DOUBLE PRECISION :: tCdag 1976 DOUBLE PRECISION :: sumMmat 1977 DOUBLE PRECISION :: sumCheck 1978 #include "BathOperatoroffdiag_hybrid.h" 1979 1980 aF = op%activeFlavor 1981 !Construction de la matrix 1982 tail = op%sumtails 1983 ! CALL MatrixHyb_init(checkMatrix,op%iTech,size=tail,Wmax=op%samples) 1984 ! CALL MatrixHyb_setSize(checkMatrix,tail) 1985 1986 ! --- set size of the matrix 1987 CALL MatrixHyb_setSize(op%M_update,tail) 1988 1989 ! --- compute useful quantities 1990 beta = op%beta 1991 mbeta_two = -beta*0.5d0 1992 op%checkNumber = op%checkNumber + 1 1993 IF ( tail .NE. op%M%tail ) THEN 1994 CALL WARN("BathOperatoroffdiag_checkM : tails are different ") 1995 RETURN 1996 END IF 1997 1998 do it=1,op%sumtails 1999 !write(6,*) " checkM begin M_update%mat_tau",(op%M_update%mat_tau(it,it1),it1=1,op%sumtails) 2000 enddo 2001 ! --- build matrix 2002 !CALL ListCdagC_print(particle) 2003 DO iflavora = 1, op%flavors 2004 DO iCdag = 1, op%tails(iflavora) 2005 tCdag = particle(iflavora)%list(iCdag,Cdag_) 2006 !write(6,*) " checkM a",iflavora,tCdag 2007 DO iflavorb = 1, op%flavors 2008 DO iC = 1, op%tails(iflavorb) 2009 !tC = particle%list(C_,iC).MOD.beta 2010 MODCYCLE(particle(iflavorb)%list(iC,C_),beta,tC) ! tC is tC, or Tc-Beta if tc>beta 2011 !write(6,*) " checkM b",iflavorb,tC 2012 time = tC - tCdag ! time is positive or negative but lower than beta 2013 !write(6,*) " checkM time",time 2014 2015 #include "BathOperatoroffdiag_hybrid" 2016 2017 op%M_update%mat(op%Fshift(iflavorb)+iC,op%Fshift(iflavora)+iCdag) = hybrid 2018 2019 time = time + ( SIGN(1.d0,time) - 1.d0 )*mbeta_two 2020 op%M_update%mat_tau(op%Fshift(iflavora)+iCdag,op%Fshift(iflavorb)+iC) = INT ( (time*op%inv_dt) +1.5d0 ) 2021 !write(6,*) " checkM mat_tau",INT ( (time*op%inv_dt) +1.5d0 ) 2022 !write(6,*) " checkM shifts",op%Fshift(iflavorb),iCdag,op%Fshift(iflavora),iC 2023 END DO ! iC 2024 END DO ! iflavorb 2025 END DO ! iCdag 2026 END DO ! iflavora 2027 2028 ! CALL MatrixHyb_Print(checkMatrix) 2029 ! --- Inverse matrix 2030 CALL MatrixHyb_inverse(op%M_update) 2031 2032 ! CALL MatrixHyb_Print(checkMatrix) 2033 do it=1,op%sumtails 2034 !write(6,*) " checkM end M_update%mat_tau",(op%M_update%mat_tau(it,it1),it1=1,op%sumtails) 2035 enddo 2036 do it=1,op%sumtails 2037 !write(6,*) " checkM end M_update",(op%M%mat(it,it1),it1=1,op%sumtails) 2038 enddo 2039 2040 ! --- Compare M_update and M to check if calculation of M is correct 2041 sumMmat =0.d0 2042 sumCheck=0.d0 2043 error1 = 0.d0 2044 errormax = 0.d0 2045 checkTau = .FALSE. 2046 DO iCdag = 1, tail 2047 Do iC =1, tail 2048 errorrel= ABS((op%M_update%mat(iC, iCdag) - & 2049 op%M%mat(iC,iCdag))/op%M_update%mat(iC,iCdag)) 2050 errorabs= ABS(op%M_update%mat(iC, iCdag) - & 2051 op%M%mat(iC,iCdag)) 2052 IF ( errorrel .gt. errormax .and. errorabs .gt. 0.001d0 ) errormax = errorrel 2053 ! write(6,*) " checkM ", errorrel,errorabs 2054 IF ( op%M_update%mat_tau(iC,iCdag) .NE. op%M%mat_tau(iC,iCdag) ) then 2055 checkTau = .TRUE. 2056 !write(6,*) "op%M_update%mat_tau(iC,iCdag), op%M%mat_tau(iC,iCdag)",op%M_update%mat_tau(iC,iCdag), op%M%mat_tau(iC,iCdag) 2057 !call flush(6) 2058 CALL ERROR("BathOperatoroffdiag_checkM : "//a//"% ") 2059 ENDIF 2060 2061 END DO 2062 END DO 2063 2064 IF ( checkTau .EQV. .TRUE. ) THEN 2065 CALL WARN("BathOperatoroffdiag_checkM : mat_tau differs should be") 2066 CALL MatrixHyb_print(op%M_update,opt_print=1) 2067 CALL WARN("BathOperatoroffdiag_checkM : whereas it is") 2068 CALL MatrixHyb_print(op%M,opt_print=1) 2069 END IF 2070 op%meanError = op%meanError + errormax 2071 IF ( errormax .GT. 1.d0 ) THEN 2072 WRITE(a,'(I4)') INT(error1*100.d0) 2073 !write(6,'(I4)') INT(error1*100.d0) 2074 ! CALL MatrixHyb_Print(op%M) 2075 CALL WARN("BathOperatoroffdiag_checkM") 2076 END IF 2077 ! CALL MatrixHyb_destroy(checkMatrix) 2078 END SUBROUTINE BathOperatoroffdiag_checkM
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_doCheck [ Functions ]
NAME
BathOperatoroffdiag_doCheck
FUNCTION
Just store if we perfom check for updates of M
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator opt_check=second bit should be one
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1915 SUBROUTINE BathOperatoroffdiag_doCheck(op,opt_check) 1916 1917 !Arguments ------------------------------------ 1918 TYPE(BathOperatoroffdiag) , INTENT(INOUT) :: op 1919 INTEGER , INTENT(IN ) :: opt_check 1920 1921 IF ( opt_check .GE. 2 ) & 1922 op%doCheck = .TRUE. 1923 END SUBROUTINE BathOperatoroffdiag_doCheck
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_getDetAdd [ Functions ]
NAME
BathOperatoroffdiag_getDetAdd
FUNCTION
Compute the determinant ratio when a (anti)segment is trying to be added and store some array for setMAdd
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator CdagC_1=segment to be added position=ordered position of the Cdag time particle=full list of CdagC for activeFlavor
OUTPUT
BathOperatoroffdiag_getDetAdd=the det
SIDE EFFECTS
NOTES
SOURCE
418 DOUBLE PRECISION FUNCTION BathOperatoroffdiag_getDetAdd(op,CdagC_1, position, particle) 419 420 !Arguments ------------------------------------ 421 TYPE(BathOperatoroffdiag) , INTENT(INOUT) :: op 422 DOUBLE PRECISION, DIMENSION(1:2), INTENT(IN ) :: CdagC_1 423 INTEGER , INTENT(IN ) :: position 424 TYPE(ListCdagC), INTENT(IN ) :: particle(:) 425 !Local variables------------------------------- 426 INTEGER :: it1 427 INTEGER :: it2 428 INTEGER :: it3,iflavor,iflavora,iflavorb 429 INTEGER :: iflavorbegin,iflavorend 430 INTEGER :: tail,tailbegin,tailend 431 INTEGER :: tcheck 432 INTEGER :: new_tail 433 DOUBLE PRECISION :: C 434 DOUBLE PRECISION :: Cbeta 435 DOUBLE PRECISION :: Cibeta 436 DOUBLE PRECISION :: Cdag 437 DOUBLE PRECISION :: Cdagbeta 438 DOUBLE PRECISION :: beta 439 DOUBLE PRECISION :: ratio 440 DOUBLE PRECISION :: time 441 ! TYPE(CdagC) , POINTER, DIMENSION(:) :: list => NULL() 442 #include "BathOperatoroffdiag_hybrid.h" 443 444 op%antiShift = .FALSE. 445 beta = op%beta 446 C = CdagC_1(C_) 447 ! Cbeta = C.MOD.beta 448 MODCYCLE(C,beta,Cbeta) 449 Cdag = CdagC_1(Cdag_) 450 ! cdagbeta = Cdag.MOD.beta 451 MODCYCLE(Cdag,beta,Cdagbeta) 452 ! IF ( Cdag .GE. beta ) & 453 ! CALL ERROR("BathOperatoroffdiag_getDetAdd : bad case ... ") 454 IF ( op%activeFlavor .LE. 0 ) & 455 CALL ERROR("BathOperatoroffdiag_getDetAdd : no active hybrid function ") 456 457 IF ( size(particle)/=op%flavors ) & 458 CALL ERROR("BathOperatoroffdiag_getDetAdd : size of particle is erroneous ") 459 460 ! tail is now the complete size of the F matrix Fshift(nflavors+1) 461 tail = op%sumtails 462 new_tail = tail+1 463 ! list => particle%list 464 465 if(op%opt_nondiag==1) then 466 iflavorbegin = 1 467 iflavorend = op%flavors 468 tailbegin = 1 469 tailend = tail 470 else 471 !sui!write(6,*) "Bathoperator opt_nondiag=0" 472 iflavorbegin = op%activeflavor 473 iflavorend = op%activeflavor 474 tailbegin = op%Fshift(op%activeflavor)+1 475 tailend = op%Fshift(op%activeflavor)+op%tails(op%activeflavor) 476 endif 477 478 IF ( ((C .GT. Cdag) .AND. (position .EQ. -1)) & ! Segment added at the end of the segment 479 .OR. ((C .LT. Cdag) .AND. (tail .EQ. 0))) THEN ! empty orbital case: only adding a segment is possible 480 ! If ones add a segment to an empty orbital or a segment at the end 481 ! of a segment, then: 482 op%updatePosRow = op%tails(op%activeFlavor) + 1 483 op%updatePosCol = op%tails(op%activeFlavor) + 1 484 ELSE 485 ! For all the other cases, ABS(position) is the true position. 486 op%updatePosRow = ABS(position) 487 op%updatePosCol = ABS(position) 488 END IF 489 !write(6,*) " BathOperatoroffdiag_getDetAdd : op%updatePosRow",op%updatePosRow 490 !write(6,*) " BathOperatoroffdiag_getDetAdd : op%updatePosCol",op%updatePosCol 491 !write(6,*) " BathOperatoroffdiag_getDetAdd : C,Cdag",C,Cdag 492 493 IF ( C .LT. Cdag .AND. op%tails(op%activeFlavor) .GT. 0) THEN ! only if an antisegment is added 494 ! ratio = -ratio 495 op%updatePosRow = (op%updatePosRow + 1) !position in [1;tail] 496 ! If the antisegment created is such that a segment with tcdagger> tc 497 ! is suppressed 498 !write(6,*) " BathOperatoroffdiag_getDetAdd : op%updatePosRow",op%updatePosRow 499 !write(6,*) " BathOperatoroffdiag_getDetAdd : op%updatePosCol",op%updatePosCol 500 IF ( Cdagbeta .LT. particle(op%activeFlavor)%list(op%updatePosCol,Cdag_) ) op%antiShift = .TRUE. 501 END IF 502 503 ! CALL Vector_setSize(op%R,tail) 504 ! CALL Vector_setSize(op%Q,tail) 505 Vector_QuickResize(op%R,new_tail) 506 Vector_QuickResize(op%Q,new_tail) 507 Vector_QuickResize(op%Rtau,new_tail) 508 Vector_QuickResize(op%Qtau,new_tail) 509 510 ! This loop compute all Row and Col except op%updatePosRow 511 tcheck=0 512 DO iflavor = iflavorbegin,iflavorend 513 !write(6,*) " BathOperatoroffdiag_getDetAdd : tails(iflavor)",iflavor,op%tails(iflavor) 514 DO it1 = 1, op%tails(iflavor) 515 tcheck=tcheck+1 516 it2 = it1 517 it3 = it1 518 IF ( iflavor .GE. op%activeFlavor ) THEN 519 it2 = it1 + 1 520 it3 = it1 + 1 521 IF ( iflavor .EQ. op%activeFlavor .AND. it1 .LT. op%updatePosRow ) it2 = it1 522 IF ( iflavor .EQ. op%activeFlavor .AND. it1 .LT. op%updatePosCol ) it3 = it1 523 !it2 = it1 + ( 1+SIGN(1,it1-op%updatePosRow) )/2 524 !it3 = it1 + ( 1+SIGN(1,it1-op%updatePoscol) )/2 525 ! if it1>=op%updatePosRow and iflavor> activeflavor, then it2=it1+1 526 ! if it1< op%updatePosRow and iflavor> activeflavor, then it2=it1 527 END IF 528 529 !!write(6,*) size(op%Rtau%vec) 530 !!write(6,*) size(particle(iflavor)%list,1) 531 !!write(6,*) size(particle(iflavor)%list,2) 532 !!write(6,*) size(op%Fshift) 533 !!write(6,*) it1,Cdag_,op%Fshift(iflavor)+it2 534 op%Rtau%vec(op%Fshift(iflavor)+it2)= C - particle(iflavor)%list(it1,Cdag_) 535 ! the following line happend only for nondiag case 536 IF(op%Rtau%vec(op%Fshift(iflavor)+it2) .GT. beta) op%Rtau%vec(op%Fshift(iflavor)+it2)=op%Rtau%vec(op%Fshift(iflavor)+it2)-beta 537 !op%Rtau%vec(it1)= C - particle%list(it1,Cdag_) 538 time = Cbeta - particle(iflavor)%list(it1,Cdag_) 539 if(op%Rtau%vec(op%Fshift(iflavor)+it2)>beta) then 540 !write(6,*) "Rtau sup beta",op%Rtau%vec(op%Fshift(iflavor)+it2),C,particle(iflavor)%list(it1,Cdag_) 541 !write(6,*) time 542 stop 543 endif 544 545 ! "BathOperatoroffdiag_hybrid" interpolates between known values of F for the 546 ! selected time. 547 iflavora=iflavor 548 iflavorb=op%activeFlavor 549 #include "BathOperatoroffdiag_hybrid" 550 551 op%R%vec(op%Fshift(iflavor)+it1) = hybrid 552 ! op%R%vec(it) = BathOperatoroffdiag_hybrid(op, Cbeta - list(it)%Cdag) 553 ! Cibeta = list(it)%C.MOD.beta 554 MODCYCLE(particle(iflavor)%list(it1,C_),beta,Cibeta) 555 time = Cibeta - Cdagbeta 556 op%Qtau%vec(op%Fshift(iflavor)+it3)= time 557 !op%Qtau%vec(it1)= time 558 559 iflavora=op%activeFlavor 560 iflavorb=iflavor 561 #include "BathOperatoroffdiag_hybrid" 562 op%Q%vec(op%Fshift(iflavor)+it1) = hybrid 563 564 !op%Q%vec(it3) = hybrid 565 ! Q(it) = BathOperatoroffdiag_hybrid(op, Cibeta - Cdagbeta) 566 END DO 567 END DO 568 if(tcheck.ne.tail) then 569 !write(6,*) " PRB in the loop tail tcheck",tail,tcheck 570 stop 571 endif 572 573 ! Compute S 574 op%Stau = C - Cdagbeta 575 op%Rtau%vec(op%Fshift(op%activeFlavor)+op%updatePosRow) = op%Stau 576 if(op%Rtau%vec(op%Fshift(op%activeFlavor)+op%updatePosRow)>beta) then 577 !write(6,*) "Rtau sup beta", op%Stau,C,Cdagbeta 578 stop 579 endif 580 op%Qtau%vec(op%Fshift(op%activeFlavor)+op%updatePosCol) = op%Rtau%vec(op%Fshift(op%activeFlavor)+op%updatePosRow) 581 !write(6,*) " getdetAdd op%Stau",op%Stau 582 583 time = Cbeta-Cdagbeta 584 !write(6,*) " getdetAdd time",time 585 iflavora=op%activeFlavor 586 iflavorb=op%activeFlavor 587 !write(6,*) "time",time 588 #include "BathOperatoroffdiag_hybrid" 589 !write(6,*) "hybrid",hybrid 590 op%S = hybrid 591 !write(6,*) " getdetAdd hybrid=op%S",hybrid 592 593 !ratio = op%S - DOT_PRODUCT(MATMUL(op%R%vec(1:tail),op%M(op%activeFlavor)%mat(1:tail,1:tail)),op%Q%vec(1:tail)) 594 595 ! product of matrix R and M(k) is computed now: 596 ratio = 0.d0 597 DO it1 = tailbegin, tailend 598 time = 0.d0 599 DO it2 = tailbegin, tailend 600 time = time + op%R%vec(it2) * op%M%mat(it2,it1) 601 END DO 602 ratio = ratio + op%Q%vec(it1) * time 603 END DO 604 !sui!write(6,*) " = R Matrix",tail 605 !sui!write(6,*) " R ",(op%R%vec(it1),it1=1,tail) 606 !sui!write(6,*) " = Q Matrix",tail 607 !sui!do it1=1,tail 608 !sui!write(6,*) " Q ",op%Q%vec(it1) 609 !sui!enddo 610 !sui!write(6,*) " = M Matrix",tail 611 !sui!do it2=1,tail 612 !sui!write(6,*) " M ",(op%M%mat(it2,it1),it1=1,tail) 613 !sui!enddo 614 !sui!write(6,*) " RMQ =", ratio 615 !sui!write(6,*) " S =", op%S 616 ratio = op%S - ratio 617 !sui!write(6,*) " S-RMQ =", ratio 618 !sui!write(6,*) " getdetAdd ratio",ratio 619 620 op%Stilde = 1.d0 / ratio 621 ! If antisegment, the det ratio has to be multiplied by -1 ( sign of the signature of one 622 ! permutation line in the matrix) 623 IF ( C .LT. Cdag .AND. op%tails(op%activeFlavor) .GT. 0) THEN ! only if an antisegment is added 624 ratio=-ratio 625 ENDIF 626 627 ! This IF is the LAST "NON CORRECTION" in my opinion this should not appears. 628 ! IF ( MAX(C,Cdag) .GT. op%beta ) THEN 629 ! WRITE(*,*) op%Stilde 630 ! op%Stilde = - ABS(op%Stilde) 631 ! END IF 632 BathOperatoroffdiag_getDetAdd = ratio 633 !write(6,*) " getdetAdd",ratio,BathOperatoroffdiag_getDetAdd 634 op%MAddFlag = .TRUE. 635 !#ifdef CTQMC_CHECK 636 ! op%ListCdagC = particle 637 !!write(*,*) op%Stilde 638 !!write(*,*) op%antishift 639 !!write(*,*) op%updatePosRow 640 !!write(*,*) op%updatePosCol 641 !#endif 642 643 END FUNCTION BathOperatoroffdiag_getDetAdd
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_getDetF [ Functions ]
NAME
BathOperatoroffdiag_getDetF
FUNCTION
Compute the determinant of the F matrix using the hybridization of flavor and the segments of particle used for Gloval moves only
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator flavor=hybridization function to take particles=segments to use
OUTPUT
BathOperatoroffdiag_getDetF=the det
SIDE EFFECTS
NOTES
SOURCE
757 DOUBLE PRECISION FUNCTION BathOperatoroffdiag_getDetF(op,particle,option) 758 759 !Arguments ------------------------------------ 760 TYPE(BathOperatoroffdiag) , INTENT(INOUT) :: op 761 TYPE(ListCdagC), OPTIONAL, INTENT(IN ) :: particle(:) 762 INTEGER , optional :: option 763 !Local arguments------------------------------- 764 INTEGER :: iCdag 765 INTEGER :: iC 766 INTEGER :: tail 767 DOUBLE PRECISION :: time 768 DOUBLE PRECISION :: tC 769 DOUBLE PRECISION :: tCdag 770 DOUBLE PRECISION :: beta 771 DOUBLE PRECISION :: mbeta_two 772 DOUBLE PRECISION :: signe 773 DOUBLE PRECISION :: inv_dt 774 INTEGER :: iflavor,iflavora 775 INTEGER :: iflavordag,iflavorb 776 #include "BathOperatoroffdiag_hybrid.h" 777 778 BathOperatoroffdiag_getDetF = 1.d0 ! pour eviter des divisions par 0 779 IF ( PRESENT( particle ) ) THEN 780 tail = op%sumtails 781 beta = op%beta 782 mbeta_two = -beta*0.5d0 783 inv_dt = op%inv_dt 784 CALL MatrixHyb_setSize(op%M_update,tail) 785 DO iflavordag=1,op%flavors 786 DO iCdag = 1, op%tails(iflavordag) 787 tCdag = particle(iflavordag)%list(iCdag,Cdag_) 788 DO iflavor=1,op%flavors 789 DO iC = 1, op%tails(iflavor) 790 !tC = particle%list(C_,iC).MOD.beta 791 MODCYCLE(particle(iflavor)%list(iC,C_),beta,tC) 792 time = tC - tCdag 793 iflavora=iflavordag 794 iflavorb=iflavor 795 #include "BathOperatoroffdiag_hybrid" 796 op%M_update%mat(op%Fshift(iflavor)+iC,op%Fshift(iflavordag)+iCdag) = hybrid 797 END DO 798 END DO 799 END DO 800 END DO 801 ! mat_tau needs to be transpose of ordered time mat (way of measuring 802 ! G(tau)) 803 DO iflavor=1,op%flavors 804 DO iC = 1, tail 805 tC = particle(iflavor)%list(iC,C_) 806 DO iflavordag=1,op%flavors 807 DO iCdag = 1, tail 808 !sui!write(6,*) iCdag,Cdag_,size(particle(iflavordag)%list,1) 809 !stop 810 tCdag = particle(iflavordag)%list(iCdag,Cdag_) 811 time = tC - tCdag 812 signe = SIGN(1.d0,time) 813 time = time + (signe-1.d0)*mbeta_two 814 op%M_update%mat_tau(op%Fshift(iflavordag)+iCdag,op%Fshift(iflavor)+iC) = INT( ( time * inv_dt ) + 1.5d0 ) 815 END DO 816 END DO 817 END DO 818 END DO 819 CALL MatrixHyb_inverse(op%M_update,BathOperatoroffdiag_getDetF) ! calcul le det de la matrice et l'inverse 820 ELSE 821 if(present(option)) then 822 CALL MatrixHyb_getDet(op%M_update,BathOperatoroffdiag_getDetF) ! det M = 1/detF ! 823 else 824 CALL MatrixHyb_getDet(op%M,BathOperatoroffdiag_getDetF) ! det M = 1/detF ! 825 endif 826 BathOperatoroffdiag_getDetF = 1.d0 / BathOperatoroffdiag_getDetF 827 ENDIF 828 END FUNCTION BathOperatoroffdiag_getDetF
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_getDetRemove [ Functions ]
NAME
BathOperatoroffdiag_getDetRemove
FUNCTION
Compute the determinant ratio when a (anti)segment is trying to be removed
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator position=position of segment to be removed
OUTPUT
BathOperatoroffdiag_getDetRemove=the det
SIDE EFFECTS
NOTES
SOURCE
673 DOUBLE PRECISION FUNCTION BathOperatoroffdiag_getDetRemove(op,position) 674 675 !Arguments ------------------------------------ 676 TYPE(BathOperatoroffdiag), INTENT(INOUT) :: op 677 !Local arguments------------------------------- 678 INTEGER , INTENT(IN ) :: position 679 INTEGER :: ABSposition 680 INTEGER :: tail !,it,it1 681 682 IF ( op%activeFlavor .LE. 0 ) & 683 CALL ERROR("BathOperatoroffdiag_getDetRemove : no active hybrid fun ") 684 685 op%antiShift = .FALSE. 686 tail = op%sumtails 687 ABSposition = ABS(position) 688 IF ( ABSposition .GT. op%tails(op%activeFlavor) ) & 689 CALL ERROR("BathOperatoroffdiag_getDetRemove : position > M size ") 690 op%updatePosCol = ABSposition 691 op%antiShift = .FALSE. 692 IF ( position .GT. 0 ) THEN 693 op%updatePosRow = ABSposition 694 ELSE 695 op%updatePosRow = ABSposition+1 696 IF ( ABSposition .EQ. op%tails(op%activeFlavor) ) THEN 697 op%antiShift = .TRUE. 698 op%updatePosRow = 1 !ABSposition - 1 699 ! op%updatePosRow = ABSposition 700 ! IF ( op%updatePosCol .EQ. 0) op%updatePosCol = tail 701 END IF 702 ENDIF 703 op%Stilde = op%M%mat(op%Fshift(op%activeFlavor)+& 704 & op%updatePosRow,op%Fshift(op%activeFlavor)+op%updatePosCol) 705 !sui!write(6,*) "Fshift",op%Fshift(op%activeFlavor) 706 !sui!write(6,*) "updatepos",op%updatePosRow,op%updatePosCol 707 708 709 op%MRemoveFlag = .TRUE. 710 !write(6,*) " getdetRemove",op%Stilde 711 BathOperatoroffdiag_getDetRemove = op%Stilde 712 if(position<0.and.op%tails(op%activeFlavor)>1) then 713 BathOperatoroffdiag_getDetRemove = -op%Stilde 714 endif 715 !do it=1,op%sumtails 716 !!sui!write(6,*) " getdetRemove M",(op%M%mat(it,it1),it1=1,op%sumtails) 717 !enddo 718 !#ifdef CTQMC_CHECK 719 ! op%ListCdagC = particle 720 !!write(*,*) op%updatePosRow, op%updatePosCol, position 721 !!CALL ListCdagC_print(particle) 722 !#endif 723 724 END FUNCTION BathOperatoroffdiag_getDetRemove
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_getError [ Functions ]
NAME
BathOperatoroffdiag_getError
FUNCTION
compute a percentage error / checkM
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator
OUTPUT
BathOperatoroffdiag_getError=Error in percent
SIDE EFFECTS
NOTES
SOURCE
2227 DOUBLE PRECISION FUNCTION BathOperatoroffdiag_getError(op) 2228 2229 TYPE(BathOperatoroffdiag), INTENT(IN) :: op 2230 2231 IF ( op%doCheck .EQV. .TRUE. ) THEN 2232 BathOperatoroffdiag_getError = op%meanError / DBLE(op%checkNumber) 2233 ELSE 2234 BathOperatoroffdiag_getError = 0.d0 2235 END IF 2236 END FUNCTION BathOperatoroffdiag_getError
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_init [ Functions ]
NAME
BathOperatoroffdiag_init
FUNCTION
Initialize and allocate data
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath object flavors=numbers of flavors we have (including spin) samples=Time slices in the input file beta=inverse temperature iTech=imaginary time or frequencies It is imposes to imaginary time
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
226 SUBROUTINE BathOperatoroffdiag_init(op, flavors, samples, beta, iTech,opt_nondiag) 227 228 !Arguments ------------------------------------ 229 TYPE(BathOperatoroffdiag), INTENT(INOUT) :: op 230 INTEGER , INTENT(IN ) :: flavors 231 INTEGER , INTENT(IN ) :: samples 232 INTEGER , INTENT(IN ) :: opt_nondiag 233 DOUBLE PRECISION , INTENT(IN ) :: beta 234 !Local variables ------------------------------ 235 INTEGER , INTENT(IN ) :: iTech 236 !INTEGER :: it 237 238 op%MAddFlag = .FALSE. 239 op%MRemoveFlag = .FALSE. 240 op%flavors = flavors 241 op%opt_nondiag = opt_nondiag 242 op%beta = beta 243 op%samples = samples 244 op%sizeHybrid = samples + 1 245 op%dt = beta / DBLE(samples) 246 op%inv_dt = DBLE(samples) / beta 247 op%activeFlavor= 0 248 op%updatePosRow = 0 249 op%updatePosCol = 0 250 op%iTech = iTech 251 !#ifdef CTQMC_CHECK 252 op%checkNumber = 0 253 op%meanError = 0.d0 254 op%doCheck = .FALSE. 255 !#endif 256 257 FREEIF(op%F) 258 MALLOC(op%F,(1:op%sizeHybrid+1,1:flavors,1:flavors)) 259 DT_FREEIF(op%tails) 260 #ifdef FC_LLVM 261 ! LLVM 16 doesn't recognize this macro here 262 DT_MALLOC(op%tails, (1:op%flavors)) 263 #else 264 DT_MALLOC(op%tails,(1:op%flavors)) 265 #endif 266 op%tails=0 267 DT_FREEIF(op%Fshift) 268 #ifdef FC_LLVM 269 ! LLVM 16 doesn't recognize this macro here 270 DT_MALLOC(op%Fshift, (1:op%flavors+1)) 271 #else 272 DT_MALLOC(op%Fshift,(1:op%flavors+1)) 273 #endif 274 op%Fshift=0 275 276 CALL Vector_init(op%R,100*op%flavors) 277 CALL Vector_init(op%Q,100*op%flavors) 278 CALL Vector_init(op%Rtau,100*op%flavors) 279 CALL Vector_init(op%Qtau,100*op%flavors) 280 281 CALL MatrixHyb_init(op%M,op%iTech,size=Global_SIZE*op%flavors,Wmax=samples) !FIXME Should be consistent with ListCagC 282 CALL MatrixHyb_init(op%M_update,op%iTech,size=Global_SIZE*op%flavors,Wmax=samples) !FIXME Should be consistent with ListCagC 283 op%F = 0.d0 284 op%set = .TRUE. 285 286 END SUBROUTINE BathOperatoroffdiag_init
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_initF [ Functions ]
NAME
BathOperatoroffdiag_initF
FUNCTION
Copy input hybridization functions from a file
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator ifstream=file stream to read F
OUTPUT
argout(sizeout)=description
SIDE EFFECTS
NOTES
SOURCE
1601 SUBROUTINE BathOperatoroffdiag_initF(op,ifstream) 1602 1603 !Arguments ---------------------- 1604 TYPE(BathOperatoroffdiag), INTENT(INOUT) :: op 1605 INTEGER , INTENT(IN ) :: ifstream 1606 !Local variables ---------------- 1607 INTEGER :: iflavor1 1608 INTEGER :: iflavor2 1609 INTEGER :: sample 1610 1611 IF ( op%set .EQV. .FALSE. ) & 1612 CALL ERROR("BathOperatoroffdiag_initF : BathOperatoroffdiag not set ") 1613 1614 DO iflavor1=1,op%flavors 1615 DO iflavor2=2,op%flavors 1616 DO sample = 1, op%sizeHybrid 1617 READ(ifstream,*) op%F(sample,iflavor1,iflavor2) 1618 END DO 1619 END DO 1620 END DO 1621 END SUBROUTINE BathOperatoroffdiag_initF
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_printF [ Functions ]
NAME
BathOperatoroffdiag_printF
FUNCTION
print F function
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator ostream=file stream to write in
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1702 SUBROUTINE BathOperatoroffdiag_printF(op,ostream) 1703 1704 !Arguments ------------------------------------ 1705 TYPE(BathOperatoroffdiag), INTENT(INOUT) :: op 1706 INTEGER,OPTIONAL , INTENT(IN ) :: ostream 1707 !Local variables ------------------------------ 1708 CHARACTER(LEN=4) :: aflavor 1709 CHARACTER(LEN=50) :: string 1710 INTEGER :: iflavor1 1711 INTEGER :: iflavor2 1712 INTEGER :: sample 1713 INTEGER :: ostream_val 1714 1715 IF ( PRESENT(ostream) ) THEN 1716 ostream_val = ostream 1717 ELSE 1718 ostream_val = 65 1719 OPEN(UNIT=ostream_val, FILE="F.dat") 1720 END IF 1721 1722 WRITE(aflavor,'(I4)') (op%flavors*op%flavors+1) 1723 string = '(1x,'//TRIM(ADJUSTL(aflavor))//'E22.14)' 1724 DO sample = 1, op%sizeHybrid 1725 WRITE(ostream_val,string) (sample-1)*op%dt, ((op%F(sample,iflavor1,iflavor2),& 1726 iflavor1=1,op%flavors),iflavor2=1,op%flavors) 1727 END DO 1728 !CALL FLUSH(ostream_val) 1729 1730 IF ( .NOT. PRESENT(ostream) ) & 1731 CLOSE(ostream_val) 1732 1733 END SUBROUTINE BathOperatoroffdiag_printF
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_printM [ Functions ]
NAME
BathOperatoroffdiag_printM
FUNCTION
print M =F^{-1} matrix
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator ostream=file stream to write in
OUTPUT
argout(sizeout)=description
SIDE EFFECTS
NOTES
SOURCE
1762 SUBROUTINE BathOperatoroffdiag_printM(op,ostream) 1763 1764 !Arguments ------------------------------------ 1765 TYPE(BathOperatoroffdiag), INTENT(IN) :: op 1766 INTEGER, OPTIONAL , INTENT(IN) :: ostream 1767 !Local variables ------------------------------ 1768 INTEGER :: ostream_val 1769 1770 IF ( op%activeFlavor .LE. 0 ) & 1771 CALL ERROR("BathOperatoroffdiag_printM : no active hybrid function ") 1772 ostream_val = 6 1773 IF ( PRESENT(ostream) ) ostream_val = ostream 1774 CALL MatrixHyb_print(op%M,ostream_val) 1775 END SUBROUTINE BathOperatoroffdiag_printM
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_printM_matrix [ Functions ]
NAME
BathOperatoroffdiag_printM_matrix
FUNCTION
print M =F^{-1} matrix
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator ostream=file stream to write in
OUTPUT
argout(sizeout)=description
SIDE EFFECTS
NOTES
SOURCE
1804 SUBROUTINE BathOperatoroffdiag_printM_matrix(op,ostream) 1805 1806 !Arguments ------------------------------------ 1807 TYPE(BathOperatoroffdiag), INTENT(IN) :: op 1808 INTEGER, OPTIONAL , INTENT(IN) :: ostream 1809 !Local variables ------------------------------ 1810 INTEGER :: iflavor1 1811 INTEGER :: i1,it1,it2 1812 CHARACTER(LEN=22) :: string 1813 CHARACTER(LEN=22) :: string2 1814 CHARACTER(LEN=4 ) :: size 1815 1816 ABI_UNUSED(ostream) 1817 1818 WRITE(size,'(I4)') op%sumtails 1819 string ='(i2,x,i3,a,'//TRIM(ADJUSTL(size))//'(E5.2,1x))' 1820 string2 ='(6x,'//TRIM(ADJUSTL(size))//'(i6))' 1821 open(unit=222, file="M_matrix.dat") 1822 open(unit=223, file="M_matrix_tau.dat") 1823 it1=0 1824 write(222,string2) ((i1,i1=1,op%tails(iflavor1)),iflavor1=1,op%flavors) 1825 do iflavor1=1, op%flavors 1826 do i1=1, op%tails(iflavor1) 1827 it1=it1+1 1828 write(222,string) iflavor1,i1,'|',(op%M%mat(it1,it2),it2=1,op%sumtails) 1829 enddo 1830 enddo 1831 1832 1833 END SUBROUTINE BathOperatoroffdiag_printM_matrix
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_recomputeM [ Functions ]
NAME
BathOperatoroffdiag_recomputeM
FUNCTION
compute from scratch the M matrix
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (B. Amadon) 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
op=bath operator particle=list of all segments of the active flavor
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
2106 SUBROUTINE BathOperatoroffdiag_recomputeM(op,particle,flav_i,flav_j) 2107 2108 !Arguments ------------------------------------ 2109 TYPE(BathOperatoroffdiag) , INTENT(INOUT) :: op 2110 TYPE(ListCdagC) , INTENT(IN ) :: particle(:) 2111 INTEGER :: flav_i,flav_j 2112 !Local variables ------------------------------ 2113 ! TYPE(MatrixHyb) :: checkMatrix 2114 INTEGER :: tail 2115 INTEGER :: iC 2116 INTEGER :: iCdag 2117 INTEGER :: aF 2118 INTEGER :: iflavora 2119 INTEGER :: iflavorb,it !,it1 2120 INTEGER :: iflavora_imp 2121 INTEGER :: iflavorb_imp 2122 !CHARACTER(LEN=6) :: a 2123 DOUBLE PRECISION :: time 2124 DOUBLE PRECISION :: beta 2125 DOUBLE PRECISION :: mbeta_two 2126 DOUBLE PRECISION :: tc 2127 DOUBLE PRECISION :: tCdag 2128 !DOUBLE PRECISION :: sumMmat 2129 !DOUBLE PRECISION :: sumCheck 2130 #include "BathOperatoroffdiag_hybrid.h" 2131 2132 aF = op%activeFlavor 2133 !Construction de la matrix 2134 tail = op%sumtails 2135 ! CALL MatrixHyb_init(checkMatrix,op%iTech,size=tail,Wmax=op%samples) 2136 ! CALL MatrixHyb_setSize(checkMatrix,tail) 2137 2138 ! --- set size of the matrix 2139 CALL MatrixHyb_setSize(op%M_update,tail) 2140 2141 ! --- compute useful quantities 2142 beta = op%beta 2143 mbeta_two = -beta*0.5d0 2144 op%checkNumber = op%checkNumber + 1 2145 IF ( tail .NE. op%M%tail ) THEN 2146 CALL WARN("BathOperatoroffdiag_checkM : tails are different ") 2147 RETURN 2148 END IF 2149 2150 do it=1,op%sumtails 2151 !write(6,*) " checkM begin M_update%mat_tau",(op%M_update%mat_tau(it,it1),it1=1,op%sumtails) 2152 enddo 2153 ! --- build matrix 2154 !CALL ListCdagC_print(particle) 2155 DO iflavora = 1, op%flavors 2156 iflavora_imp=iflavora 2157 if(iflavora==flav_i) iflavora_imp=flav_j 2158 if(iflavora==flav_j) iflavora_imp=flav_i 2159 DO iCdag = 1, op%tails(iflavora_imp) 2160 tCdag = particle(iflavora_imp)%list(iCdag,Cdag_) 2161 !write(6,*) " checkM a",iflavora,tCdag 2162 DO iflavorb = 1, op%flavors 2163 iflavorb_imp=iflavorb 2164 if(iflavorb==flav_j) iflavorb_imp=flav_i 2165 if(iflavorb==flav_i) iflavorb_imp=flav_j 2166 DO iC = 1, op%tails(iflavorb_imp) 2167 !tC = particle%list(C_,iC).MOD.beta 2168 MODCYCLE(particle(iflavorb_imp)%list(iC,C_),beta,tC) ! tC is tC, or Tc-Beta if tc>beta 2169 !write(6,*) " checkM b",iflavorb,tC 2170 time = tC - tCdag ! time is positive or negative but lower than beta 2171 !write(6,*) " checkM time",time 2172 2173 #include "BathOperatoroffdiag_hybrid" 2174 2175 op%M_update%mat(op%Fshift(iflavorb_imp)+iC,op%Fshift(iflavora_imp)+iCdag) = hybrid 2176 2177 time = time + ( SIGN(1.d0,time) - 1.d0 )*mbeta_two 2178 op%M_update%mat_tau(op%Fshift(iflavora_imp)+iCdag,op%Fshift(iflavorb_imp)+iC) = INT ( (time*op%inv_dt) +1.5d0 ) 2179 !write(6,*) " checkM mat_tau",INT ( (time*op%inv_dt) +1.5d0 ) 2180 !write(6,*) " checkM shifts",op%Fshift(iflavorb),iCdag,op%Fshift(iflavora),iC 2181 END DO ! iC 2182 END DO ! iflavorb 2183 END DO ! iCdag 2184 END DO ! iflavora 2185 2186 ! CALL MatrixHyb_Print(checkMatrix) 2187 ! --- Inverse matrix 2188 CALL MatrixHyb_inverse(op%M_update) 2189 2190 ! CALL MatrixHyb_Print(checkMatrix) 2191 do it=1,op%sumtails 2192 !write(6,*) " checkM end M_update%mat_tau",(op%M_update%mat_tau(it,it1),it1=1,op%sumtails) 2193 enddo 2194 do it=1,op%sumtails 2195 !write(6,*) " checkM end M_update",(op%M%mat(it,it1),it1=1,op%sumtails) 2196 enddo 2197 2198 ! --- Compare M_update and M to check if calculation of M is correct 2199 END SUBROUTINE BathOperatoroffdiag_recomputeM
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_reset [ Functions ]
NAME
BathOperatoroffdiag_reset
FUNCTION
Reset all internal variables
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator to reset
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
313 SUBROUTINE BathOperatoroffdiag_reset(op) 314 315 !Arguments ------------------------------------ 316 TYPE(BathOperatoroffdiag), INTENT(INOUT) :: op 317 !Local variables ------------------------------ 318 INTEGER :: iflavor 319 op%MAddFlag = .FALSE. 320 op%MRemoveFlag = .FALSE. 321 op%activeFlavor = 0 322 op%updatePosRow = 0 323 op%updatePosCol = 0 324 !#ifdef CTQMC_CHECK 325 op%checkNumber = 0 326 op%meanError = 0.d0 327 op%sumtails = 0 328 !#endif 329 op%doCheck = .FALSE. 330 CALL Vector_clear(op%R) 331 CALL Vector_clear(op%Q) 332 CALL Vector_clear(op%Rtau) 333 CALL Vector_clear(op%Qtau) 334 335 CALL MatrixHyb_clear(op%M) !FIXME Should be consistent with ListCagC 336 op%F = 0.d0 337 do iflavor=1,op%flavors 338 op%tails(iflavor)=0 339 op%Fshift(iflavor)=0 340 enddo 341 342 END SUBROUTINE BathOperatoroffdiag_reset
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_setF [ Functions ]
NAME
BathOperatoroffdiag_setF
FUNCTION
Copy F from input array
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator F=array of the hybridization function
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1649 SUBROUTINE BathOperatoroffdiag_setF(op,F) 1650 1651 !Arguments ------------------------------------ 1652 TYPE(BathOperatoroffdiag) , INTENT(INOUT) :: op 1653 DOUBLE PRECISION, DIMENSION(:,:,:) , INTENT(IN ) :: F 1654 !Arguments ------------------------------------ 1655 INTEGER :: iflavor1 1656 INTEGER :: iflavor2 1657 INTEGER :: sample 1658 INTEGER :: length 1659 1660 IF ( op%set .EQV. .FALSE. ) & 1661 CALL ERROR("BathOperatoroffdiag_setF : BathOperatoroffdiag not set ") 1662 1663 length = SIZE(F) 1664 IF ( length .NE. (op%flavors * op%flavors * op%sizeHybrid) ) & 1665 CALL ERROR("BathOperatoroffdiag_setF : wrong input F ") 1666 1667 DO iflavor1=1,op%flavors 1668 DO iflavor2=1,op%flavors 1669 DO sample = 1, op%sizeHybrid 1670 op%F(sample,iflavor1,iflavor2) = F(sample,iflavor1,iflavor2) 1671 END DO 1672 END DO 1673 END DO 1674 END SUBROUTINE BathOperatoroffdiag_setF
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_setMAdd [ Functions ]
NAME
BathOperatoroffdiag_setMAdd
FUNCTION
Update de M matrix inserting a row and a column
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator particle=segments of active flavor
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
856 SUBROUTINE BathOperatoroffdiag_setMAdd(op,particle) 857 858 !Arguments ------------------------------------ 859 TYPE(BathOperatoroffdiag), INTENT(INOUT) :: op 860 TYPE(ListCdagC) , INTENT(IN ) :: particle(:) 861 !Local variables ------------------------------ 862 INTEGER :: tail 863 INTEGER :: new_tail 864 INTEGER :: col 865 INTEGER :: col_move 866 INTEGER :: row_move 867 INTEGER :: row 868 INTEGER :: positionRow 869 INTEGER :: positionCol 870 INTEGER :: aF,indice 871 INTEGER :: tailb,taile 872 DOUBLE PRECISION :: Stilde 873 DOUBLE PRECISION :: time 874 DOUBLE PRECISION :: mbeta_two 875 DOUBLE PRECISION :: inv_dt 876 TYPE(Vector) :: vec_tmp 877 TYPE(VectorInt) :: vecI_tmp 878 INTEGER :: m 879 INTEGER :: count 880 INTEGER :: i 881 INTEGER :: j 882 INTEGER :: p,it !,it1 883 884 885 886 ! --- op%MAddFlag is put to .TRUE. in BathOperatoroffdiag_getDetAdd. 887 IF ( op%MAddFlag .EQV. .FALSE. ) & 888 CALL ERROR("BathOperatoroffdiag_setMAdd : MAddFlag turn off ") 889 890 ! --- op%activeFlavor is put in ctqmc_loop 891 aF = op%activeFlavor 892 IF ( aF .LE. 0 ) & 893 CALL ERROR("BathOperatoroffdiag_setMAdd : no active hybrid function ") 894 895 !! do it=1,op%sumtails 896 !write(6,*) " setMAdd begin M",(op%M%mat(it,it1),it1=1,op%sumtails) 897 !! enddo 898 !! do it=1,op%sumtails 899 !write(6,*) " setMAdd begin M%mat_tau",(op%M%mat_tau(it,it1),it1=1,op%sumtails) 900 !! enddo 901 ! old tail 902 !write(6,*) " BathOperatoroffdiag_setMAdd op%sumtails",op%sumtails 903 tail = op%sumtails 904 new_tail = tail + 1 905 op%tails(aF)= op%tails(aF) + 1 906 DO indice = aF +1, op%flavors+1 907 op%Fshift(indice) = op%Fshift(indice) + 1 908 END DO 909 op%sumtails = op%Fshift(op%flavors) + op%tails(op%flavors) !last slot of Fshift is the tail of full matrix 910 !write(6,*) " BathOperatoroffdiag_setMAdd op%sumtails",op%sumtails 911 !write(6,*) " setMAdd actualized Fshift",(op%Fshift(it),it=1,op%flavors+1) 912 !write(6,*) " setMAdd actualized tails",(op%tails(it),it=1,op%flavors) 913 !CALL matrix_print(M) 914 915 if(op%opt_nondiag==1) then 916 tailb = 1 917 taile = tail 918 else 919 !sui!write(6,*) "Bathoperator a opt_nondiag=0" 920 tailb = op%Fshift(aF)+1 921 taile = op%Fshift(aF)+op%tails(aF) 922 endif 923 924 ! --- data obtained from BathOperatoroffdiag_getDetAdd 925 PositionRow = op%updatePosRow + op%Fshift(aF) ! position in the full matrix 926 PositionCol = op%updatePosCol + op%Fshift(aF) ! position in the full matrix 927 Stilde = op%Stilde 928 929 ! !write(6,*) "before", positionRow, positionCol 930 !CALL MatrixHyb_print(op%M(aF),opt_print=1) 931 ! --- MatrixHyb_setSize 932 !write(6,*) " BathOperatoroffdiag_setMAdd before setsize",size(op%M%mat,1) 933 CALL MatrixHyb_setSize(op%M,new_tail) 934 !write(6,*) " BathOperatoroffdiag_setMAdd after setsize",size(op%M%mat,1) 935 936 ! Compute Qtilde with Q 937 !op%Q%vec(1:tail) = (-1.d0) * MATMUL(op%M(aF)%mat(1:tail,1:tail),op%Q%vec(1:tail)) * Stilde 938 939 ! --- M*Q => Q 940 op%Q%vec(tailb:taile) = MATMUL(op%M%mat(tailb:taile,tailb:taile),op%Q%vec(tailb:taile)) 941 942 !op%Q%vec(PositionRow:new_tail) = EOSHIFT(op%Q%vec(PositionRow:new_tail), SHIFT=-1, BOUNDARY=-1.d0, DIM=1) 943 ! op%Qtau%vec(PositionCol:new_tail) = EOSHIFT(op%Qtau%vec(PositionCol:new_tail), SHIFT=-1, BOUNDARY=1.d0, DIM=1) 944 ! op%Qtau%vec(PositionCol) = op%Stau 945 946 !Compute Rtilde with R and without multiplying by Stilde 947 !op%R%vec(1:tail) = (-1.d0) * MATMUL(op%R%vec(1:tail),op%M(aF)%mat(1:tail,1:tail)) 948 949 ! --- R*M => R 950 op%R%vec(tailb:taile) = MATMUL(op%R%vec(tailb:taile),op%M%mat(tailb:taile,tailb:taile)) 951 952 !op%R%vec(PositionCol:new_tail) = EOSHIFT(op%R%vec(PositionCol:new_tail), SHIFT=-1, BOUNDARY=-1.d0, DIM=1) 953 ! op%Rtau%vec(PositionRow:new_tail) = EOSHIFT(op%Rtau%vec(PositionRow:new_tail), SHIFT=-1, BOUNDARY=1.d0, DIM=1) 954 ! op%Rtau%vec(PositionRow) = op%Stau 955 956 !Compute the new M matrix 957 !op%M(aF)%mat(PositionRow:new_tail,1:new_tail) = & 958 ! EOSHIFT(op%M(aF)%mat(PositionRow:new_tail,1:new_tail),SHIFT=-1, BOUNDARY=0.d0, DIM=1) 959 !op%M(aF)%mat(1:n12 characters (ABI_MALLOC) instead of 4 (FREE)ew_tail,PositionCol:new_tail) = & 960 ! EOSHIFT(op%M(aF)%mat(1:new_tail,PositionCol:new_tail),SHIFT=-1, BOUNDARY=0.d0, DIM=2) 961 ! ! op%M(aF)%mat(1:new_tail,1:new_tail) = op%M(aF)%mat(1:new_tail,1:new_tail) + & 962 ! ! Stilde * MATMUL(RESHAPE(op%Q%vec(1:new_tail),(/ new_tail,1 /)),RESHAPE(op%R%vec(1:new_tail),(/ 1,new_tail /))) 963 964 !op%M(aF)%mat_tau(PositionRow:new_tail,1:new_tail) = & 965 ! EOSHIFT(op%M(aF)%mat_tau(PositionRow:new_tail,1:new_tail),SHIFT=-1, BOUNDARY=0, DIM=1) 966 !op%M(aF)%mat_tau(1:new_tail,PositionCol:new_tail) = & 967 ! EOSHIFT(op%M(aF)%mat_tau(1:new_tail,PositionCol:new_tail),SHIFT=-1, BOUNDARY=0, DIM=2) 968 969 mbeta_two = -op%beta*0.5d0 970 inv_dt = op%inv_dt 971 972 ! ------ Shift mat_tau and update old M=Ptilde 973 DO col=tail,1,-1 ! decreasing order to avoid overwrite of data 974 col_move = col + ( 1+SIGN(1,col-PositionCol) )/2 975 ! if col>= PositionCol col_move=col+1 976 ! if col< PositionCol col_move=col 977 DO row=tail,1,-1 978 row_move = row + ( 1+SIGN(1,row-PositionRow) )/2 979 ! --- times for Ptilde are kept unchanged. But we have to copy it at the right place 980 op%M%mat_tau(row_move,col_move) = & 981 op%M%mat_tau(row,col) 982 ! --- Update Ptilde with the same indices as mat_tau 983 ! --- M + M*Q Stilde R*M => Ptilde => M 984 !if(row>=tailb.and.row<=taile.and.col>=tailb.and.col<=taile) then 985 op%M%mat(row_move,col_move) = & 986 op%M%mat(row,col) + op%Q%vec(row)*op%R%vec(col) * Stilde 987 !else 988 ! op%M%mat(row_move,col_move) = op%M%mat(row,col) 989 !endif 990 END DO 991 END DO 992 993 ! ------ Add new stuff for new row 994 DO row = 1, tail 995 row_move = row + ( 1+SIGN(1,row-PositionRow) )/2 996 ! --- M*Q Stilde => Qtilde => M with the good indices 997 !if(row>=tailb.and.row<=taile) then 998 op%M%mat(row_move,PositionCol) = -op%Q%vec(row)*Stilde 999 !else 1000 ! op%M%mat(row_move,PositionCol) = op%M%mat(row,PositionCol) 1001 !endif 1002 1003 time = op%Rtau%vec(row) ! pourquoi Rtau et pas Qtau ici ? 1004 time = time + ( SIGN(1.d0,time) - 1.d0 )*mbeta_two 1005 ! if time>=0 time=time 1006 ! if time< 0 time=time + beta 1007 ! --- mat_tau=int(time*L/beta+1.5) 1008 op%M%mat_tau(row,PositionCol) = INT ( (time*inv_dt) +1.5d0 ) 1009 !write(6,*) " setMadd new row", op%Rtau%vec(row),op%M%mat_tau(row,PositionCol) 1010 ! if(op%M%mat_tau(row,PositionCol)>301) then 1011 ! !write(6,*) ">301 a", time,inv_dt, op%M%mat_tau(row,PositionCol) 1012 ! time = op%Rtau%vec(row) ! pourquoi Rtau et pas Qtau ici ? 1013 ! !write(6,*) time,mbeta_two 1014 ! time = time + ( SIGN(1.d0,time) - 1.d0 )*mbeta_two 1015 ! !write(6,*) time 1016 ! !write(6,*) INT ( (time*inv_dt) +1.5d0 ) 1017 ! stop 1018 ! endif 1019 END DO 1020 ! Add last time missing in the loops 1021 time = op%Rtau%vec(new_tail) 1022 time = time + ( SIGN(1.d0,time) - 1.d0 )*mbeta_two 1023 op%M%mat_tau(new_tail,PositionCol) = INT ( (time*inv_dt) +1.5d0 ) 1024 !write(6,*) " setMadd last time", op%Rtau%vec(new_tail),op%M%mat_tau(new_tail,PositionCol) 1025 ! if(op%M%mat_tau(new_tail,PositionCol)>301) then 1026 ! !write(6,*) ">301 b", time,inv_dt, op%M%mat_tau(new_tail,PositionCol) 1027 ! time = op%Rtau%vec(new_tail) ! pourquoi Rtau et pas Qtau ici ? 1028 ! !write(6,*) time,mbeta_two 1029 ! time = time + ( SIGN(1.d0,time) - 1.d0 )*mbeta_two 1030 ! !write(6,*) time 1031 ! !write(6,*) INT ( (time*inv_dt) +1.5d0 ) 1032 ! stop 1033 ! endif 1034 1035 ! Add new stuff for new col 1036 DO col = 1, tail 1037 col_move = col + ( 1+SIGN(1,col-PositionCol) )/2 1038 ! --- Stilde RN => Rtilde => M 1039 !if(col>=tailb.and.col<=taile) then 1040 op%M%mat(PositionRow,col_move) = -op%R%vec(col)*Stilde 1041 !else 1042 ! op%M%mat(PositionRow,col_move) = op%M%mat(PositionRow,col) 1043 !endif 1044 time = op%Qtau%vec(col) 1045 time = time + ( SIGN(1.d0,time) - 1.d0 )*mbeta_two 1046 op%M%mat_tau(PositionRow,col) = INT ( (time*inv_dt) +1.5d0 ) 1047 !write(6,*) " setMadd new col", op%Qtau%vec(col),op%M%mat_tau(PositionRow,col) 1048 ! if(op%M%mat_tau(PositionRow,col)>301) then 1049 ! !write(6,*) ">301 c", time,inv_dt, op%M%mat_tau(PositionRow,col) 1050 ! time = op%Qtau%vec(col) ! pourquoi Rtau et pas Qtau ici ? 1051 ! !write(6,*) time,mbeta_two 1052 ! time = time + ( SIGN(1.d0,time) - 1.d0 )*mbeta_two 1053 ! !write(6,*) time 1054 ! !write(6,*) INT ( (time*inv_dt) +1.5d0 ) 1055 ! stop 1056 ! endif 1057 END DO 1058 ! Add last time missing in the loops 1059 time = op%Qtau%vec(new_tail) 1060 time = time + ( SIGN(1.d0,time) - 1.d0 )*mbeta_two 1061 op%M%mat_tau(PositionRow,new_tail) = INT ( (time*inv_dt) +1.5d0 ) 1062 !write(6,*) " setMadd last time", op%Qtau%vec(new_tail),op%M%mat_tau(PositionRow,new_tail) 1063 ! if(op%M%mat_tau(PositionRow,new_tail)>301) then 1064 ! !write(6,*) ">301 d", time,inv_dt, op%M%mat_tau(PositionRow,new_tail) 1065 ! time = op%Qtau%vec(new_tail) ! pourquoi Rtau et pas Qtau ici ? 1066 ! !write(6,*) time,mbeta_two 1067 ! time = time + ( SIGN(1.d0,time) - 1.d0 )*mbeta_two 1068 ! !write(6,*) time 1069 ! !write(6,*) INT ( (time*inv_dt) +1.5d0 ) 1070 ! stop 1071 ! endif 1072 1073 op%M%mat(PositionRow,PositionCol) = Stilde 1074 1075 !CALL MatrixHyb_print(op%M,opt_print=1) 1076 1077 ! DO col = 1, new_tail 1078 ! time = op%Rtau%vec(col) 1079 ! time = time + ( SIGN(1.d0,time) - 1.d0 )*mbeta_two 1080 ! op%M(aF)%mat_tau(col,PositionCol) = INT ( (time*inv_dt) +1.5d0 ) 1081 ! time = op%Qtau%vec(col) 1082 ! time = time + ( SIGN(1.d0,time) - 1.d0 )*mbeta_two 1083 ! op%M(aF)%mat_tau(PositionRow,Col) = INT ( (time*inv_dt) +1.5d0 ) 1084 ! time = op%R%vec(col)*Stilde 1085 ! DO row = 1, new_tail 1086 ! op%M(aF)%mat(row,col) = op%M(aF)%mat(row,col) + op%Q%vec(row)*time 1087 ! END DO 1088 ! END DO 1089 1090 !col_move = new_tail 1091 !col = tail 1092 !DO col_move = new_tail, 1, -1 1093 ! IF ( col_move .EQ. positionCol ) THEN 1094 ! ! on calcule rajoute Q tilde 1095 ! !row_move = new_tail 1096 ! row = tail 1097 ! DO row_move = new_tail, 1, -1 1098 ! ! calcul itau 1099 ! IF ( row_move .EQ. positionRow ) THEN 1100 ! op%M(aF)%mat(row_move,col_move) = Stilde 1101 ! !time = op%Stau 1102 ! ELSE 1103 ! op%M(aF)%mat(row_move,col_move) = -op%Q%vec(row)*Stilde 1104 ! !time = op%Rtau%vec(row_move) 1105 ! row = row - 1 1106 ! END IF 1107 ! !time = time + ( SIGN(1.d0,time) - 1.d0 )*mbeta_two 1108 ! !op%M(aF)%mat_tau(row_move,col_move) = INT ( (time*inv_dt) +1.5d0 ) 1109 ! END DO 1110 ! ! realignement des indices 1111 ! ELSE 1112 ! ! on calcule Ptilde 1113 ! !row_move = new_tail 1114 ! row = tail 1115 ! DO row_move = new_tail, 1, -1 1116 ! IF ( row_move .EQ. positionRow ) THEN 1117 ! op%M(aF)%mat(row_move,col_move) = -op%R%vec(col) * Stilde 1118 ! ! calcul itau 1119 ! !time = op%Qtau%vec(col_move) 1120 ! !time = time + ( SIGN(1.d0,time) - 1.d0 )*mbeta_two 1121 ! !op%M(aF)%mat_tau(row_move,col_move) = INT ( (time*inv_dt) +1.5d0 ) 1122 ! ELSE 1123 ! op%M(aF)%mat(row_move,col_move) = op%M(aF)%mat(row,col) + op%Q%vec(row)*op%R%vec(col)*Stilde 1124 ! ! copy itau 1125 ! !op%M(aF)%mat_tau(row_move,col_move) = op%M(aF)%mat_tau(row,col) 1126 ! row = row - 1 1127 ! END IF 1128 ! END DO 1129 ! col = col - 1 1130 ! END IF 1131 !END DO 1132 ! !write(6,*) "after" 1133 ! CALL MatrixHyb_print(op%M(aF),opt_print=1) 1134 !CALL matrix_inverse(M) 1135 !CALL MatrixHyb_print(M) 1136 !CALL matrix_inverse(M) 1137 1138 IF ( op%antiShift .EQV. .TRUE. ) THEN ! antisegment 1139 if(3==4) then 1140 CALL Vector_init(vec_tmp,new_tail) 1141 CALL VectorInt_init(vecI_tmp,new_tail) 1142 ! Shift if necessary according to op%antishift 1143 ! shift DIM=2 (col) 1144 1145 ! For new_tail=4, the following lines transform 1146 ! M=(a,b,c,d) vith a,b,c,d column vectors into 1147 ! M=(d,a,b,c) 1148 p = new_tail - 1 ! = tail 1149 m = 1 1150 ! count increases in the loop from 0 to new_tail-1 1151 count = 0 1152 DO WHILE ( count .NE. new_tail ) 1153 ! put column b in vec_tmp 1154 vec_tmp%vec(1:new_tail) = op%M%mat(1:new_tail,m) 1155 vecI_tmp%vec(1:new_tail) = op%M%mat_tau(1:new_tail,m) 1156 i = m 1157 !j = m+p 1158 MODCYCLE(m+p, new_tail, j) ! j=m+p modulo new_tail 1159 DO WHILE (j .NE. m) 1160 op%M%mat(1:new_tail,i) = op%M%mat(1:new_tail,j) 1161 op%M%mat_tau(1:new_tail,i) = op%M%mat_tau(1:new_tail,j) 1162 i = j 1163 MODCYCLE(j+p, new_tail, j) 1164 count = count+1 1165 END DO 1166 op%M%mat(1:new_tail,i) = vec_tmp%vec(1:new_tail) 1167 op%M%mat_tau(1:new_tail,i) = vecI_tmp%vec(1:new_tail) 1168 count = count+1 1169 m = m+1 1170 END DO 1171 ! shift DIM=1 (row) 1172 1173 ! below is similar to above but for rows instead of columns. 1174 p = new_tail - 1 1175 m = 1 1176 count = 0 1177 DO WHILE ( count .NE. new_tail) 1178 vec_tmp%vec(1:new_tail) = op%M%mat(m,1:new_tail) 1179 vecI_tmp%vec(1:new_tail) = op%M%mat_tau(m,1:new_tail) 1180 i = m 1181 !j = m+p 1182 MODCYCLE(m+p, new_tail, j) 1183 DO WHILE ( j .NE. m ) 1184 op%M%mat(i,1:new_tail) = op%M%mat(j,1:new_tail) 1185 op%M%mat_tau(i,1:new_tail) = op%M%mat_tau(j,1:new_tail) 1186 i = j 1187 MODCYCLE(j+p, new_tail, j) 1188 count = count+1 1189 END DO 1190 op%M%mat(i,1:new_tail) = vec_tmp%vec(1:new_tail) 1191 op%M%mat_tau(i,1:new_tail) = vecI_tmp%vec(1:new_tail) 1192 count = count+1 1193 m = m+1 1194 END DO 1195 CALL Vector_destroy(vec_tmp) 1196 CALL VectorInt_destroy(vecI_tmp) 1197 endif 1198 !op%M(aF)%mat(1:new_tail,1:new_tail) = CSHIFT(op%M(aF)%mat(1:new_tail,1:new_tail), SHIFT=-1, DIM=1) ! Shift to the bottom 1199 !op%M(aF)%mat(1:new_tail,1:new_tail) = CSHIFT(op%M(aF)%mat(1:new_tail,1:new_tail), SHIFT=-1, DIM=2) ! Shift to the right 1200 !op%M(aF)%mat_tau(1:new_tail,1:new_tail) = CSHIFT(op%M(aF)%mat_tau(1:new_tail,1:new_tail), SHIFT=-1, DIM=1) ! Shift to the bottom 1201 !op%M(aF)%mat_tau(1:new_tail,1:new_tail) = CSHIFT(op%M(aF)%mat_tau(1:new_tail,1:new_tail), SHIFT=-1, DIM=2) ! Shift to the right 1202 !write(6,*) " setMAdd size M%mat",size(op%M%mat,1),size(op%M%mat,2),new_tail 1203 !write(6,*) " setMAdd arguement M%mat",aF,op%Fshift(aF) 1204 !write(6,*) " setMAdd arguement M%mat",aF,op%Fshift(aF),op%Fshift(aF+1) 1205 do it=1,op%sumtails 1206 !write(6,*) " setMAdd before antishift M%mat_tau",(op%M%mat_tau(it,it1),it1=1,op%sumtails) 1207 enddo 1208 if (new_tail>0.and.op%Fshift(aF+1)>op%Fshift(aF)) then 1209 op%M%mat(op%Fshift(aF)+1:op%Fshift(aF+1) , 1:new_tail) = & 1210 CSHIFT( op%M%mat(op%Fshift(aF)+1:op%Fshift(aF+1) , 1:new_tail) , SHIFT=-1 , DIM=1) ! Shift to the bottom 1211 1212 op%M%mat(1:new_tail , op%Fshift(aF)+1:op%Fshift(aF+1)) = & 1213 CSHIFT( op%M%mat(1:new_tail , op%Fshift(aF)+1:op%Fshift(aF+1)) , SHIFT=-1 , DIM=2) ! Shift to the right 1214 1215 op%M%mat_tau(op%Fshift(aF)+1:op%Fshift(aF+1) , 1:new_tail) = & 1216 CSHIFT( op%M%mat_tau(op%Fshift(aF)+1:op%Fshift(aF+1) , 1:new_tail) , SHIFT=-1 , DIM=1) ! Shift to the bottom 1217 1218 op%M%mat_tau(1:new_tail , op%Fshift(aF)+1:op%Fshift(aF+1)) = & 1219 CSHIFT( op%M%mat_tau(1:new_tail , op%Fshift(aF)+1:op%Fshift(aF+1)) , SHIFT=-1 , DIM=2) ! Shift to the right 1220 end if 1221 !CALL matrix_print(M) 1222 END IF 1223 1224 !! do it=1,op%sumtails 1225 !write(6,*) " setMAdd end M",(op%M%mat(it,it1),it1=1,op%sumtails) 1226 !! enddo 1227 !! do it=1,op%sumtails 1228 !! !write(6,*) " setMAdd end M%mat_tau",(op%M%mat_tau(it,it1),it1=1,op%sumtails) 1229 !! enddo 1230 IF ( op%doCheck .EQV. .TRUE.) THEN 1231 !#ifdef CTQMC_CHECK 1232 CALL BathOperatoroffdiag_checkM(op,particle) 1233 !#endif 1234 END IF 1235 1236 op%MAddFlag = .FALSE. 1237 1238 END SUBROUTINE BathOperatoroffdiag_setMAdd
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_setMRemove [ Functions ]
NAME
BathOperatoroffdiag_setMRemove
FUNCTION
delete one row and one column of the M matrix
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator particle=segments of the active flavor
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1266 SUBROUTINE BathOperatoroffdiag_setMRemove(op,particle) 1267 1268 !Arguments ------------------------------------ 1269 TYPE(BathOperatoroffdiag), INTENT(INOUT) :: op 1270 TYPE(ListCdagC) , INTENT(IN ) :: particle(:) 1271 !Local variables ------------------------------ 1272 INTEGER :: tail,tailb,taile 1273 INTEGER :: new_tail 1274 INTEGER :: col 1275 INTEGER :: col_move 1276 INTEGER :: row_move 1277 INTEGER :: row 1278 INTEGER :: positionCol 1279 INTEGER :: positionRow 1280 INTEGER :: aF,iaf 1281 INTEGER :: m 1282 INTEGER :: count 1283 INTEGER :: i 1284 INTEGER :: j,it !,it1 1285 INTEGER :: p 1286 DOUBLE PRECISION :: invStilde 1287 DOUBLE PRECISION :: invStilde2 1288 TYPE(VectorInt) :: vecI_tmp 1289 TYPE(Vector) :: vec_tmp 1290 1291 IF ( op%MRemoveFlag .EQV. .FALSE. ) & 1292 CALL ERROR("BathOperatoroffdiag_setMRemove : MRemoveFlag turn off ") 1293 aF = op%activeFlavor 1294 IF ( aF .LE. 0 ) & 1295 CALL ERROR("BathOperatoroffdiag_setMRemove : no active hybrid func ") 1296 do it=1,op%sumtails 1297 !write(6,*) " setMRemove begin M",(op%M%mat(it,it1),it1=1,op%sumtails) 1298 enddo 1299 tail = op%sumtails 1300 new_tail = tail - 1 1301 op%tails(af)= op%tails(af) - 1 1302 DO iaf=af+1 , op%flavors+1 1303 op%Fshift(iaf) = op%Fshift(iaf) - 1 1304 END DO 1305 op%sumtails = op%Fshift(op%flavors) + op%tails(op%flavors) 1306 positionCol = op%updatePosCol + op%Fshift(af) 1307 positionRow = op%updatePosRow + op%Fshift(af) 1308 invStilde = 1.d0 / op%Stilde 1309 if(op%opt_nondiag==1) then 1310 tailb = 1 1311 taile = new_tail 1312 else 1313 !sui!write(6,*) "Bathoperator c opt_nondiag=0" 1314 tailb = op%Fshift(aF)+1 1315 taile = op%Fshift(aF)+op%tails(aF) 1316 endif 1317 1318 ! !write(6,*) "before", positionRow, positionCol 1319 ! CALL MatrixHyb_print(op%M(aF),opt_print=1) 1320 1321 ! IF ( new_tail .EQ. 0 ) THEN 1322 !! IF ( op%antiShift .EQV. .TRUE. ) THEN 1323 !! op%M(aF)%mat(1,1) = 1.d0/BathOperatoroffdiag_Hybrid(op, op%beta) 1324 !! op%MRemoveFlag = .FALSE. 1325 !! RETURN 1326 !! END IF 1327 ! CALL MatrixHyb_clear(op%M(aF)) 1328 ! op%MRemoveFlag = .FALSE. 1329 ! RETURN 1330 ! END IF 1331 1332 ! CALL Vector_setSize(op%Q,new_tail) 1333 ! CALL Vector_setSize(op%R,new_tail) 1334 Vector_QuickResize(op%Q,new_tail) 1335 Vector_QuickResize(op%R,new_tail) 1336 1337 ! We use R and Q as op%R%vec and op%Q%vec 1338 ! op%R%vec => op%R 1339 ! op%Q%vec => op%Q 1340 1341 row = 1 1342 !row_move = 1 1343 col = 1 1344 !col_move = 1 1345 DO row_move = 1, new_tail 1346 IF ( row .EQ. positionRow ) row = row + 1 1347 IF ( col .EQ. positionCol ) col = col + 1 1348 !col = row_move + (1+SIGN(1,row_move-positionCol))/2 1349 !row = row_move + (1+SIGN(1,row_move-positionRow))/2 1350 op%R%vec(row_move) = op%M%mat(positionRow,col) 1351 op%Q%vec(row_move) = op%M%mat(row,positionCol) 1352 row = row + 1 1353 col = col + 1 1354 END DO 1355 !! op%R%vec(1:positionCol-1) = op%M(aF)%mat(positionRow,1:positionCol-1) 1356 !! op%R%vec(positionCol:new_tail) = op%M(aF)%mat(positionRow,positionCol+1:tail) 1357 !! op%Q%vec(1:positionRow-1) = op%M(aF)%mat(1:positionRow-1,positionCol) 1358 !! op%Q%vec(positionRow:new_tail) = op%M(aF)%mat(positionRow+1:tail,positionCol) 1359 !write(*,*) positionRow, positionCol 1360 !CALL MatrixHyb_print(M) 1361 !CALL Vector_print(op%R) 1362 !CALL Vector_print(op%Q) 1363 !CALL ListCdagC_print(op%ListCdagC) 1364 1365 col = 1 1366 DO col_move = 1, new_tail 1367 IF ( col_move .EQ. positionCol ) col = col + 1 1368 !col = col_move + (1+SIGN(1,col_move-positionCol))/2 1369 row = 1 1370 invStilde2 = invStilde * op%R%vec(col_move) 1371 DO row_move = 1, new_tail 1372 IF ( row_move .EQ. positionRow ) row = row + 1 1373 !row = row_move + (1+SIGN(1,row_move-positionRow))/2 1374 ! Compute for all rows and cols M <= M - Q 1/S R 1375 !if(row_move>=tailb.and.row_move<=taile.and.col_move>=tailb.and.col_move<=taile) then 1376 op%M%mat(row_move,col_move) = op%M%mat(row,col) & 1377 - op%Q%vec(row_move)*invStilde2 1378 !else 1379 ! op%M%mat(row_move,col_move) = op%M%mat(row,col) 1380 !endif 1381 op%M%mat_tau(row_move,col_move) = op%M%mat_tau(row,col) 1382 row = row + 1 1383 END DO 1384 col = col + 1 1385 END DO 1386 CALL MatrixHyb_setSize(op%M,new_tail) 1387 1388 IF ( op%antiShift .EQV. .TRUE. ) THEN ! antisegment 1389 if(3==4) then 1390 ! Shift if necessary according to op%antishift 1391 ! shift DIM=2 (col) 1392 CALL Vector_init(vec_tmp,new_tail) 1393 CALL VectorInt_init(vecI_tmp,new_tail) 1394 p = 1 1395 m = 1 1396 count = 0 1397 DO WHILE ( count .NE. new_tail ) 1398 vec_tmp%vec(1:new_tail) = op%M%mat(1:new_tail,m) 1399 vecI_tmp%vec(1:new_tail) = op%M%mat_tau(1:new_tail,m) 1400 i = m 1401 !j = m+p 1402 MODCYCLE(m+p, new_tail, j) 1403 DO WHILE (j .NE. m) 1404 op%M%mat(1:new_tail,i) = op%M%mat(1:new_tail,j) 1405 op%M%mat_tau(1:new_tail,i) = op%M%mat_tau(1:new_tail,j) 1406 i = j 1407 MODCYCLE(j+p, new_tail, j) 1408 count = count+1 1409 END DO 1410 op%M%mat(1:new_tail,i) = vec_tmp%vec(1:new_tail) 1411 op%M%mat_tau(1:new_tail,i) = vecI_tmp%vec(1:new_tail) 1412 count = count+1 1413 m = m+1 1414 END DO 1415 CALL Vector_destroy(vec_tmp) 1416 CALL VectorInt_destroy(vecI_tmp) 1417 !op%M(aF)%mat(1:new_tail,1:new_tail) = & 1418 ! CSHIFT(op%M(aF)%mat(1:new_tail,1:new_tail), SHIFT=1, DIM=2) ! Shift to the top 1419 !op%M(aF)%mat_tau(1:new_tail,1:new_tail) = & 1420 ! CSHIFT(op%M(aF)%mat_tau(1:new_tail,1:new_tail), SHIFT=1, DIM=2) ! Shift to the top 1421 endif 1422 if (new_tail>0.and.op%Fshift(af+1)>op%Fshift(af)) then 1423 op%M%mat(1:new_tail,op%Fshift(af)+1:op%Fshift(af+1)) = & 1424 CSHIFT(op%M%mat(1:new_tail,op%Fshift(af)+1:op%Fshift(af+1)), SHIFT=1, DIM=2) ! Shift to the top 1425 op%M%mat_tau(1:new_tail,op%Fshift(af)+1:op%Fshift(af+1)) = & 1426 CSHIFT(op%M%mat_tau(1:new_tail,op%Fshift(af)+1:op%Fshift(af+1)), SHIFT=1, DIM=2) ! Shift to the top 1427 end if 1428 END IF 1429 ! !write(6,*) "after " 1430 ! CALL MatrixHyb_print(op%M(aF),opt_print=1) 1431 1432 IF ( op%doCheck .EQV. .TRUE. ) THEN 1433 !#ifdef CTQMC_CHECK 1434 CALL BathOperatoroffdiag_checkM(op,particle) 1435 !#endif 1436 END IF 1437 do it=1,op%sumtails 1438 !write(6,*) " setMRemove end M",(op%M%mat(it,it1),it1=1,op%sumtails) 1439 enddo 1440 1441 op%MRemoveFlag = .FALSE. 1442 1443 END SUBROUTINE BathOperatoroffdiag_setMRemove
ABINIT/m_BathOperatoroffdiag/BathOperatoroffdiag_swap [ Functions ]
NAME
BathOperatoroffdiag_swap
FUNCTION
Recompute 2 M matrix swaping the segments (used for Global moves)
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
op=bath operator iflavor1=flavor to swap with the next one iflavor2=favor to swap with the previous one
OUTPUT
SIDE EFFECTS
NOTES
SOURCE
1472 SUBROUTINE BathOperatoroffdiag_swap(op, flavor1, flavor2) 1473 1474 !Arguments ------------------------------------ 1475 TYPE(BathOperatoroffdiag), INTENT(INOUT) :: op 1476 INTEGER , INTENT(IN ) :: flavor1 1477 INTEGER , INTENT(IN ) :: flavor2 1478 INTEGER :: ii,iflavort,itmptail,flavora,flavorb 1479 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: mat_temp 1480 INTEGER , ALLOCATABLE, DIMENSION(:,:) :: mat_tau_temp 1481 1482 if(flavor1>flavor2) then 1483 flavora=flavor2 1484 flavorb=flavor1 1485 else 1486 flavora=flavor1 1487 flavorb=flavor2 1488 endif 1489 MALLOC(mat_temp,(1:op%sumtails,1:op%sumtails)) 1490 MALLOC(mat_tau_temp,(1:op%sumtails,1:op%sumtails)) 1491 !mat_temp= op%M%mat 1492 !mat_tau_temp= op%M%mat_tau 1493 !it1=0 1494 !do iflav1=1,op%flavors 1495 ! do ii1=1,op%tails(iflav1) 1496 ! it1=it1+1 1497 ! it2=0 1498 ! do iflav2=1,op%flavors 1499 ! do ii2=1,op%tails(iflav1) 1500 ! it2=it2+1 1501 ! op%M%mat(it1,it2)= 1502 ! enddo 1503 ! enddo 1504 ! enddo 1505 !enddo 1506 if(3==3) then 1507 op%M=op%M_update 1508 if (op%sumtails>0) then 1509 ! shift block flavorb at the place of flavora (column) 1510 if (op%Fshift(flavorb+1)>op%Fshift(flavora)) then 1511 do ii=1, op%tails(flavorb) 1512 op%M%mat(op%Fshift(flavora)+1:op%Fshift(flavorb+1) , 1:op%sumtails) = & 1513 CSHIFT( op%M%mat(op%Fshift(flavora)+1:op%Fshift(flavorb+1) , 1:op%sumtails) , SHIFT=-1 , DIM=1) 1514 op%M%mat_tau(op%Fshift(flavora)+1:op%Fshift(flavorb+1) , 1:op%sumtails) = & 1515 CSHIFT( op%M%mat_tau(op%Fshift(flavora)+1:op%Fshift(flavorb+1) , 1:op%sumtails) , SHIFT=-1 , DIM=1) 1516 enddo 1517 end if 1518 1519 ! shift block flavora at the place of flavorb (column) 1520 if (op%Fshift(flavorb)>op%Fshift(flavora)) then 1521 do ii=1, op%tails(flavora) 1522 op%M%mat(op%Fshift(flavora)+op%tails(flavorb)+& 1523 & 1:op%Fshift(flavorb)+op%tails(flavorb) , 1:op%sumtails) = & 1524 CSHIFT( op%M%mat( op%Fshift(flavora)+op%tails(flavorb)& 1525 & +1:op%Fshift(flavorb)+op%tails(flavorb) , 1:op%sumtails) , SHIFT=1 , DIM=1) 1526 op%M%mat_tau(op%Fshift(flavora)+op%tails(flavorb)+1:op%Fshift(flavorb)+& 1527 & op%tails(flavorb) , 1:op%sumtails) = & 1528 CSHIFT( op%M%mat_tau( op%Fshift(flavora)+op%tails(flavorb)+& 1529 & 1:op%Fshift(flavorb)+op%tails(flavorb) , 1:op%sumtails) , SHIFT=1 , DIM=1) 1530 enddo 1531 end if 1532 1533 ! shift block flavorb at the place of flavora (row) 1534 if (op%Fshift(flavorb+1)>op%Fshift(flavora)) then 1535 do ii=1, op%tails(flavorb) 1536 op%M%mat(1:op%sumtails , op%Fshift(flavora)+1:op%Fshift(flavorb+1)) = & 1537 CSHIFT( op%M%mat(1:op%sumtails , op%Fshift(flavora)+1:op%Fshift(flavorb+1)) , SHIFT=-1 , DIM=2) 1538 op%M%mat_tau(1:op%sumtails , op%Fshift(flavora)+1:op%Fshift(flavorb+1)) = & 1539 CSHIFT( op%M%mat_tau(1:op%sumtails , op%Fshift(flavora)+1:op%Fshift(flavorb+1)) , SHIFT=-1 , DIM=2) 1540 enddo 1541 end if 1542 1543 ! shift block flavora at the place of flavorb (row) 1544 if (op%Fshift(flavorb)>op%Fshift(flavora)) then 1545 do ii=1, op%tails(flavora) 1546 op%M%mat(1:op%sumtails , op%Fshift(flavora)+op%tails(flavorb)+1:op%Fshift(flavorb)+op%tails(flavorb)) = & 1547 CSHIFT( op%M%mat(1:op%sumtails ,op%Fshift(flavora)+op%tails(flavorb)& 1548 & +1:op%Fshift(flavorb)+op%tails(flavorb)) , SHIFT=1 , DIM=2) 1549 op%M%mat_tau(1:op%sumtails ,op%Fshift(flavora)+op%tails(flavorb)+& 1550 & 1:op%Fshift(flavorb)+op%tails(flavorb) ) = & 1551 CSHIFT( op%M%mat_tau(1:op%sumtails ,op%Fshift(flavora)+& 1552 & op%tails(flavorb)+1:op%Fshift(flavorb)+op%tails(flavorb) ) , SHIFT=1 , DIM=2) 1553 enddo 1554 end if 1555 end if 1556 endif 1557 if(3==4) then 1558 op%M=op%M_update 1559 endif 1560 1561 1562 do iflavort=flavora+1,flavorb 1563 op%Fshift(iflavort)=op%Fshift(iflavort)+op%tails(flavorb)-op%tails(flavora) 1564 enddo 1565 1566 itmptail=op%tails(flavora) 1567 op%tails(flavora)=op%tails(flavorb) 1568 op%tails(flavorb)=itmptail 1569 FREE(mat_temp) 1570 FREE(mat_tau_temp) 1571 1572 END SUBROUTINE BathOperatoroffdiag_swap
m_BathOperatoroffdiag/BathOperatoroffdiag [ Types ]
[ Top ] [ m_BathOperatoroffdiag ] [ Types ]
NAME
BathOperatoroffdiag
FUNCTION
This structured datatype contains the necessary data
COPYRIGHT
Copyright (C) 2013-2024 ABINIT group (J. Bieder) 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
71 TYPE BathOperatoroffdiag 72 LOGICAL :: set = .FALSE. 73 ! True if the BathOperatoroffdiag is initialized in BathOperatoroffdiag_init 74 75 LOGICAL :: MAddFlag = .FALSE. 76 ! Set to true if we can compute a new M (see updateDetXX) (ie in 77 ! BathOperatoroffdiag_getDetAdd) 78 79 LOGICAL :: MRemoveFlag = .FALSE. 80 ! Set to true if we can compute a new M (see updateDetXX) (ie in 81 ! BathOperatoroffdiag_getDetRemove) 82 83 LOGICAL :: antiShift = .FALSE. 84 ! shift when M is updated with antiseg 85 86 LOGICAL :: doCheck = .FALSE. 87 ! TRUE is checks are activated 88 89 INTEGER :: opt_nondiag = 0 90 ! if opt_nondiag = 1 F is non diagonal. 91 92 INTEGER :: flavors 93 ! number of flavors 94 ! if opt_nondiag = 0 , flavors= number of flavor 95 ! if opt_nondiag = 1 , flavors= 1 96 97 INTEGER :: activeFlavor 98 ! Active flavor on which a segment is added/suppressed... 99 100 INTEGER :: samples 101 ! Number of time slices (given in the input file) 102 103 INTEGER :: sizeHybrid 104 ! Number of time slices (given in the input file) + 1 (=qmc_l+1) 105 106 INTEGER :: updatePosRow 107 ! Gives the position of new Row to add 108 ! Modified in BathOperatoroffdiag_getDetAdd and BathOperatoroffdiag_getDetRemove 109 ! could be the Row in the Full matrix for the non diag implementation 110 111 INTEGER :: updatePosCol 112 ! Gives the position of new Col to add 113 ! Modified in BathOperatoroffdiag_getDetAdd and BathOperatoroffdiag_getDetRemove 114 115 INTEGER :: iTech 116 ! iTech is an integer which precise the technics used to compute the 117 ! Green's function (in time or frequency) 118 119 INTEGER :: sumtails 120 ! size of the full F matrix (sums of tails(iflavor) over iflavor) 121 122 INTEGER, ALLOCATABLE, DIMENSION(:) :: tails 123 ! tails(iflavor) is the current number of segments for the flavor iflavor 124 125 INTEGER, ALLOCATABLE, DIMENSION(:) :: Fshift 126 ! Fshift(iflavor) is the sum of number of segments for all flavors iflavor2 127 ! such that iflavor< iflavor 128 ! It is thus the shift in the F matrix to have the first segment of the flavor 129 ! iflavor 130 ! Fshift(nflavor+1) is the total nb of tails (=sumtails) 131 132 DOUBLE PRECISION :: beta 133 ! Inverse of Temperature 134 ! 135 136 DOUBLE PRECISION :: dt 137 ! dt=beta/samples 138 139 DOUBLE PRECISION :: inv_dt 140 ! inv_dt=1/dt 141 142 DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:,:) :: F ! qmc_l+2,Flavors 143 ! Hybridization function F(1:op%sizeHybrid+1,1:flavors,1:flavors) 144 145 DOUBLE PRECISION :: S 146 ! Sherman Morrison notations 147 148 DOUBLE PRECISION :: Stau 149 ! Sherman Morrison notations 150 151 DOUBLE PRECISION :: Stilde 152 ! Sherman Morrison notations 153 154 TYPE(Vector) :: R 155 ! Sherman Morrison notations R%vec(size). 156 ! computed for each flavor (As matrices are made of Blocks for each 157 ! flavor because the code is restricted to diagonal F matrices) 158 159 TYPE(Vector) :: Q 160 ! Sherman Morrison notations 161 ! computed for each flavor (As matrices are made of Blocks for each 162 ! flavor because the code is restricted to diagonal F matrices) 163 164 TYPE(Vector) :: Rtau 165 ! Sherman Morrison notations 166 ! Rtau gives the time length for each elements of R 167 ! computed for each flavor (As matrices are made of Blocks for each 168 ! flavor because the code is restricted to diagonal F matrices) 169 170 TYPE(Vector) :: Qtau 171 ! Sherman Morrison notations 172 ! Qtau gives the time length for each elements of Q 173 ! computed for each flavor (As matrices are made of Blocks for each 174 ! flavor because the code is restricted to diagonal F matrices) 175 176 TYPE(MatrixHyb) :: M ! Flavors 177 ! inverse of Hybridization matrix M%mat(global_size,global_size) 178 ! contains the value of the hybridization for all flavor and segments times, 179 ! the times (mat_tau), and possibly the 180 ! frequency 181 182 TYPE(MatrixHyb) :: M_update ! Flavors 183 ! used in BathOperatoroffdiag_getdetF and in BathOperatoroffdiag_checkM 184 ! for checks 185 186 !#ifdef CTQMC_CHECK 187 INTEGER :: checkNumber 188 DOUBLE PRECISION :: meanError 189 ! TYPE(ListCdagC) :: ListCdagC 190 !#endif 191 END TYPE BathOperatoroffdiag