TABLE OF CONTENTS


ABINIT/m_GreenHyb [ Modules ]

[ Top ] [ Modules ]

NAME

  m_GreenHyb

FUNCTION

  Manage a green function for one orbital

COPYRIGHT

  Copyright (C) 2013-2018 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 .

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

28 #include "defs.h"
29 MODULE m_GreenHyb
30 USE m_Global
31 USE m_MatrixHyb
32 USE m_Vector
33 USE m_VectorInt
34 USE m_ListCdagC
35 USE m_MapHyb
36 #ifdef HAVE_MPI2
37 USE mpi
38 #endif
39 IMPLICIT NONE

ABINIT/m_GreenHyb/GreenHyb_backFourier [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyb_backFourier

FUNCTION

  perform back fourier transform

COPYRIGHT

  Copyright (C) 2013-2018 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

  this=Green
  dvgc=divergence parameter

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

817 SUBROUTINE GreenHyb_backFourier(this,dvgc)
818 
819 
820 !This section has been created automatically by the script Abilint (TD).
821 !Do not modify the following lines by hand.
822 #undef ABI_FUNC
823 #define ABI_FUNC 'GreenHyb_backFourier'
824 !End of the abilint section
825 
826 
827 #ifdef HAVE_MPI1
828 include 'mpif.h'
829 #endif
830 !Arguments ------------------------------------
831   TYPE(GreenHyb)            , INTENT(INOUT) :: this
832   DOUBLE PRECISION, OPTIONAL, INTENT(IN   ) :: dvgc
833 !Local variables ------------------------------
834   INTEGER :: itau
835   INTEGER :: iomega
836   INTEGER :: omegaSamples
837   INTEGER :: tauSamples
838   INTEGER :: tauBegin
839   INTEGER :: tauEnd
840   INTEGER :: delta
841   INTEGER :: residu
842   INTEGER, ALLOCATABLE, DIMENSION(:) :: counts
843   INTEGER, ALLOCATABLE, DIMENSION(:) :: displs
844   DOUBLE PRECISION :: A ! Correction factor
845   DOUBLE PRECISION :: inv_beta
846   DOUBLE PRECISION :: pi_invBeta
847   DOUBLE PRECISION :: two_invBeta
848   DOUBLE PRECISION :: minusDt
849   DOUBLE PRECISION :: minusOmegaTau
850   DOUBLE PRECISION :: omega
851   DOUBLE PRECISION :: minusTau
852   DOUBLE PRECISION :: sumTerm
853   DOUBLE PRECISION :: pi
854   DOUBLE PRECISION :: twoPi
855   DOUBLE PRECISION :: correction
856   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: Domega
857   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: A_omega
858 
859   IF ( this%set .EQV. .FALSE. ) &
860     CALL ERROR("GreenHyb_backFourier : Uninitialized GreenHyb structure")
861   IF ( this%setW .EQV. .FALSE. ) &
862     CALL ERROR("GreenHyb_backFourier : no G(iw)")
863   
864   inv_beta     = this%inv_beta
865   two_invBeta  = 2.d0 * inv_beta
866   minusDt      = - this%delta_t
867   omegaSamples = this%Wmax
868   tauSamples   = this%samples-1
869 
870   this%oper = 0.d0
871 
872   pi         = ACOS(-1.d0)
873   twoPi        = 2.d0 * pi
874   pi_invBeta = pi * inv_beta
875 
876   IF ( PRESENT(dvgc) ) THEN
877     A = dvgc
878   ELSE
879     A = AIMAG(this%oper_w(omegaSamples)) &  ! A = \lim_\infty G(w)
880              *(2.d0*DBLE(omegaSamples)-1.d0) * pi_invBeta
881   END IF
882 
883   correction = A*0.5d0
884 
885   MALLOC(Domega,(1:omegaSamples))
886   MALLOC(A_omega,(1:omegaSamples))
887   Domega = (/ ((2.d0 * DBLE(iomega) - 1.d0)*pi_invbeta, iomega=1, omegaSamples) /)
888   A_omega = A / Domega
889   IF (this%have_MPI .EQV. .TRUE.) THEN
890     delta = tauSamples / this%size
891     residu = tauSamples - this%size*delta
892     IF ( this%rank .LT. this%size - residu ) THEN
893       tauBegin = 1 + this%rank*delta
894       tauEnd   = (this%rank + 1)*delta
895     ELSE
896 !      tauBegin = (this%size-residu)*delta + 1 + (this%rank-this%size+residu)*(delta+1)
897       tauBegin = 1 + this%rank*(delta + 1) -this%size + residu
898       tauEnd = tauBegin + delta
899     END IF
900     MALLOC(counts,(1:this%size))
901     MALLOC(displs,(1:this%size))
902     counts = (/ (delta, iTau=1, this%size-residu), &
903                 (delta+1, iTau=this%size-residu+1, this%size) /)
904     displs(1)=0
905     DO iTau = 2, this%size
906       displs(iTau) = displs(iTau-1) + counts (iTau-1)
907     END DO
908   ELSE
909     tauBegin = 1
910     tauEnd   = tauSamples
911   END IF
912   DO iTau = tauBegin, tauEnd
913     minusTau = DBLE(itau -1) * minusDt
914     DO iomega = 1, omegaSamples
915       omega         = Domega(iomega)
916       minusOmegaTau = MOD(omega*minusTau, TwoPi)
917       sumTerm       = REAL(( this%oper_w(iomega) - CMPLX(0.d0, A_omega(iomega),8) ) &
918                       * EXP( CMPLX(0.d0, minusOmegaTau, 8)))
919       this%oper(itau)    = this%oper(itau) + sumTerm
920     END DO
921     this%oper(itau) = correction + two_invBeta*this%oper(itau)
922   END DO
923   IF ( this%have_MPI .EQV. .TRUE. ) THEN
924 ! rassembler les resultats
925 #ifdef HAVE_MPI
926     CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DATATYPE_NULL, &
927                       this%oper, counts, displs, &
928                       MPI_DOUBLE_PRECISION, this%MY_COMM, residu)
929 #endif
930     FREE(counts)
931     FREE(displs)
932   END IF
933   this%oper(tauSamples+1) = A - this%oper(1) !G(0+)-G(0-)=G(0+)+G(beta-)=A
934   this%setT = .TRUE.
935   FREE(Domega)
936   FREE(A_omega)
937 
938 END SUBROUTINE GreenHyb_backFourier

ABINIT/m_GreenHyb/GreenHyb_clear [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyb_clear

FUNCTION

  clear green function

COPYRIGHT

  Copyright (C) 2013-2018 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

  this=Green

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

262 SUBROUTINE GreenHyb_clear(this)
263 
264 !Arguments ------------------------------------
265 
266 !This section has been created automatically by the script Abilint (TD).
267 !Do not modify the following lines by hand.
268 #undef ABI_FUNC
269 #define ABI_FUNC 'GreenHyb_clear'
270 !End of the abilint section
271 
272   TYPE(GreenHyb)     , INTENT(INOUT) :: this
273 
274   !CALL Vector_clear(this%oper_old)
275   !CALL VectorInt_clear(this%index_old)
276   CALL MapHyb_clear(this%this)
277   this%measurements = 0
278   IF ( ALLOCATED(this%oper) ) &
279   this%oper         = 0.d0
280   IF ( this%iTech .EQ. GREENHYB_OMEGA ) THEN
281     IF ( ALLOCATED(this%oper_w) ) &
282     this%oper_w       = CMPLX(0.d0,0.d0,8)
283     IF ( ALLOCATED(this%oper_w_old) ) &
284     this%oper_w_old   = CMPLX(0.d0,0.d0,8)
285   END IF
286   this%factor       = 0
287 END SUBROUTINE GreenHyb_clear

ABINIT/m_GreenHyb/GreenHyb_destroy [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyb_destroy

FUNCTION

  destroy green function

COPYRIGHT

  Copyright (C) 2013-2018 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

  this=Green

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

1337 SUBROUTINE GreenHyb_destroy(this)
1338 
1339 !Arguments ------------------------------------
1340 
1341 !This section has been created automatically by the script Abilint (TD).
1342 !Do not modify the following lines by hand.
1343 #undef ABI_FUNC
1344 #define ABI_FUNC 'GreenHyb_destroy'
1345 !End of the abilint section
1346 
1347   TYPE(GreenHyb), INTENT(INOUT) :: this
1348 
1349   this%set          = .FALSE.
1350   this%setT         = .FALSE.
1351   this%setW         = .FALSE.
1352   this%samples      = 0
1353   this%measurements = 0
1354   this%beta         = 0.d0
1355   this%inv_beta     = 0.d0
1356   this%inv_dt       = 0.d0
1357   this%delta_t      = 0.d0
1358   !CALL VectorInt_destroy(this%index_old)
1359   !CALL Vector_destroy(this%oper_old)
1360   CALL MapHyb_destroy(this%this)
1361   FREEIF(this%oper)
1362   FREEIF(this%oper_w)
1363   FREEIF(this%oper_w_old)
1364   FREEIF(this%omega)
1365 END SUBROUTINE GreenHyb_destroy

ABINIT/m_GreenHyb/GreenHyb_forFourier [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyb_forFourier

FUNCTION

  perform forward fourier transform

COPYRIGHT

  Copyright (C) 2013-2018 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

  this=Green
  Wmax=linear maximum frequency

OUTPUT

  Gomega=Results for omega frequencies
  omega=ask frequencies

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

 974 SUBROUTINE GreenHyb_forFourier(this, Gomega, omega, Wmax)
 975 !Arguments ------------------------------------
 976 
 977 !This section has been created automatically by the script Abilint (TD).
 978 !Do not modify the following lines by hand.
 979 #undef ABI_FUNC
 980 #define ABI_FUNC 'GreenHyb_forFourier'
 981 !End of the abilint section
 982 
 983 
 984 #ifdef HAVE_MPI1
 985 include 'mpif.h'
 986 #endif
 987   TYPE(GreenHyb)             , INTENT(INOUT) :: this
 988   COMPLEX(KIND=8), DIMENSION(:), OPTIONAL, INTENT(INOUT) :: Gomega  ! INOUT for MPI
 989   COMPLEX(KIND=8), DIMENSION(:), OPTIONAL, INTENT(IN   ) :: omega  
 990   INTEGER                 , OPTIONAL, INTENT(IN   ) :: Wmax   
 991   INTEGER :: i
 992   INTEGER :: j
 993   INTEGER :: L
 994   INTEGER :: Lspline
 995   INTEGER :: Nom
 996   INTEGER :: omegaBegin
 997   INTEGER :: omegaEnd
 998   INTEGER :: deltaw
 999   INTEGER :: residu
1000   INTEGER, ALLOCATABLE, DIMENSION(:) :: counts
1001   INTEGER, ALLOCATABLE, DIMENSION(:) :: displs
1002   DOUBLE PRECISION :: beta
1003   DOUBLE PRECISION :: tau
1004   DOUBLE PRECISION :: delta
1005   DOUBLE PRECISION :: deltabis
1006   DOUBLE PRECISION :: inv_delta
1007   DOUBLE PRECISION :: inv_delta2
1008   DOUBLE PRECISION :: omdeltabis
1009   DOUBLE PRECISION :: tmp
1010   DOUBLE PRECISION :: xpi
1011   DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE ::  diag
1012   DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE ::  diagL
1013   DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE ::  lastR
1014   DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE ::  lastC
1015   DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE ::  XM
1016   DOUBLE PRECISION, DIMENSION(:), ALLOCATABLE ::  X2
1017   DOUBLE PRECISION :: iw
1018   COMPLEX(KIND=8) :: iwtau
1019   COMPLEX(KIND=8), ALLOCATABLE, DIMENSION(:) :: Gwtmp  
1020   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:) :: omegatmp
1021 
1022   IF ( this%set .EQV. .FALSE. ) &
1023     CALL ERROR("GreenHyb_forFourier : Uninitialized GreenHyb structure")
1024   IF ( this%setT .EQV. .FALSE. ) &
1025     CALL ERROR("GreenHyb_forFourier : no G(tau)")
1026   IF ( this%setMk .NE. 2 ) &
1027     CALL WARNALL("GreenHyb_forFourier : green does not have moments    ")
1028 
1029   L  = this%samples
1030 
1031   xpi=acos(-1.d0)                !!! XPI=PI
1032   beta = this%beta
1033   Nom  = this%Wmax
1034   IF ( PRESENT(Gomega) ) THEN
1035     Nom = SIZE(Gomega)    
1036     !IF ( this%rank .EQ. 0 ) &
1037       !write(6,*) "size Gomega", Nom
1038   END IF
1039   IF ( PRESENT(omega) ) THEN
1040     IF ( PRESENT(Gomega) .AND. SIZE(omega) .NE. Nom ) THEN
1041       CALL ERROR("GreenHyb_forFourier : sizes mismatch              ")               
1042     !ELSE 
1043       !Nom = SIZE(omega)
1044     END IF 
1045   END IF
1046   IF ( .NOT. PRESENT(Gomega) .AND. .NOT. PRESENT(omega) ) THEN
1047     IF ( PRESENT(Wmax) ) THEN
1048       Nom=Wmax
1049     ELSE
1050       CALL ERROR("GreenHyb_forFourier : Missing argument Wmax")
1051     END IF
1052   END IF
1053 
1054   IF ( ALLOCATED(this%oper_w) ) THEN
1055     IF ( SIZE(this%oper_w) .NE. Nom ) THEN
1056       FREE(this%oper_w)
1057       MALLOC(this%oper_w,(1:Nom))
1058     END IF
1059   ELSE
1060     MALLOC(this%oper_w,(1:Nom))
1061   END IF
1062 
1063   !write(6,*) "PRESENT(GOMEGA)", PRESENT(GOMEGA)
1064   !write(6,*) "PRESENT(OMEGA)", PRESENT(OMEGA)
1065   !call flush(6)
1066 
1067   delta=this%delta_t
1068   inv_delta = this%inv_dt
1069   inv_delta2 = inv_delta*inv_delta
1070  
1071   MALLOC(diagL,(L-1))
1072   MALLOC(lastR,(L-1))
1073   MALLOC(diag,(L))
1074   MALLOC(lastC,(L-1))
1075 
1076 !(cf Stoer) for the spline interpolation : 
1077 ! second derivatives XM solution of A*XM=B.
1078 !A=(2.4.2.11) of Stoer&Bulirsch + 2 limit conditions
1079 !The LU decomposition of A is known explicitly; 
1080 
1081   diag (1) = 4.d0 ! 1.d0 *4.d0 factor 4 added for conditionning
1082   diagL(1) = 0.25d0 !1.d0/4.d0
1083   lastR(1) = -0.5d0 ! -2.d0/4.d0
1084   lastC(1) = 4.d0 ! 1.d0*4.d0
1085   diag (2) = 4.d0
1086   diagL(2) = 0.25d0
1087   lastR(2) = -0.25d0
1088   lastC(2) = -1.d0
1089 
1090   DO i = 3, L-2
1091     tmp = 4.d0 - diagL(i-1)
1092     diagL(i) = 1.d0 / tmp
1093   END DO
1094   DO i = 3, L-2
1095     diag (i) = 1.d0 / diagL(i)
1096     lastR(i) = -(lastR(i-1)*diagL(i))
1097     lastC(i) = -(lastC(i-1)*diagL(i-1))
1098   END DO
1099   
1100   tmp = 1.d0/diag(L-2)
1101   diag (L-1) = 4.d0 - tmp
1102   lastR(L-1) = (1.d0 - lastR(L-2))/ diag(L-1)
1103   !diagL(L-1) = lastR(L-1)
1104   diagL(L-1) = 0.d0 ! for the Lq=B resolution
1105   !lastC(L-1) = 1.d0 - lastC(L-2)*diagL(L-1) ! equivalent to the next line
1106   lastC(L-1) = 1.d0 - (lastC(L-2)*lastR(L-1)) ! True value
1107   diag (L  ) = 2.d0! - DOT_PRODUCT( lastR , lastC )
1108   tmp = 0.d0
1109   DO i = 1, L-1
1110     tmp = tmp + lastR(i)*lastC(i)
1111   END DO
1112   diag (L  ) = diag (L  ) - tmp
1113   lastC(L-1) = lastC(L-1)-1.d0 ! 1 is removed for the u.XM=q resolution
1114   
1115 ! construct the B vector from A.Xm=B
1116   MALLOC(XM,(L))
1117   XM(1) = 4.d0*this%Mk(3)
1118   XM(L) = (6.d0 * inv_delta) * ( this%Mk(2) - ( &
1119           (this%oper(2)-this%oper(1)) + &
1120           (this%oper(L)-this%oper(L-1)) ) * inv_delta )
1121   DO i = 2, L-1
1122     XM(i) = (6.d0 * inv_delta2) * ( (this%oper(i+1) &
1123                           - 2.d0 * this%oper(i)) &
1124                           +        this%oper(i-1) )
1125   END DO
1126 
1127 ! Find second derivatives XM: Solve the system
1128 ! SOLVING Lq= XM 
1129 !  q = XM 
1130   do j=1,L-1
1131       XM(j+1)=XM(j+1)-(diagL(j)*XM(j))
1132       XM(L)  =XM(L)  -(lastR(j)*XM(j))
1133   end do
1134   FREE(diagL)
1135   FREE(lastR)
1136 
1137 ! SOLVING U.XM=q 
1138 !  XM = q
1139   do j=L-1,2,-1
1140    XM(j+1)  = XM(j+1) / diag(j+1)
1141    XM(j)= (XM(j)-(XM(L)*lastC(j)))-XM(j+1)
1142   end do
1143   XM(2)  = XM(2) / diag(2)
1144   XM(1) = (XM(1)-XM(L)*lastC(1)) / diag(1)
1145 
1146   FREE(diag)
1147   FREE(lastC)
1148 
1149   Lspline = L-1
1150   MALLOC(X2,(1:Lspline+1)) ! We impose L = Nom
1151   !Construct L2 second derivative from known derivatives XM
1152   deltabis = beta / DBLE(Lspline)
1153   DO i = 1, Lspline
1154     tau = deltabis * DBLE(i-1)
1155     j = ((L-1)*(i-1))/Lspline + 1!INT(tau * inv_delta) + 1
1156     X2(i) = inv_delta * ( XM(j)*(DBLE(j)*delta - tau ) + XM(j+1)*(tau - DBLE(j-1)*delta) )
1157   END DO
1158   X2(Lspline+1) = XM(L)
1159   FREE(XM)
1160 
1161   IF ( this%have_MPI .EQV. .TRUE. ) THEN  
1162     deltaw = Nom / this%size
1163     residu = Nom - this%size*deltaw
1164     IF ( this%rank .LT. this%size - residu ) THEN
1165       omegaBegin = 1 + this%rank*deltaw
1166       omegaEnd   = (this%rank + 1)*deltaw
1167     ELSE
1168   !    tauBegin = (this%size-residu)*deltaw + 1 + (this%rank-this%size+residu)*(deltaw+1)
1169       omegaBegin = 1 + this%rank*(deltaw + 1) -this%size + residu
1170       omegaEnd = omegaBegin + deltaw
1171     END IF
1172     MALLOC(counts,(1:this%size))
1173     MALLOC(displs,(1:this%size))
1174     counts = (/ (deltaw, i=1, this%size-residu), &
1175                 (deltaw+1, i=this%size-residu+1, this%size) /)
1176     displs(1)=0
1177     DO i = 2, this%size
1178       displs(i) = displs(i-1) + counts (i-1)
1179     END DO
1180   ELSE
1181     omegaBegin = 1
1182     omegaEnd   = Nom 
1183   END IF
1184 
1185   this%Mk(1) = -1.d0
1186   MALLOC(omegatmp,(omegaBegin:omegaEnd))
1187   IF ( PRESENT(omega) ) THEN
1188     omegatmp(omegaBegin:omegaEnd) = (/ (AIMAG(omega(i)),i=omegaBegin,omegaEnd) /)
1189   ELSE
1190     omegatmp(omegaBegin:omegaEnd) = (/ ((((2.d0*DBLE(i)-1.d0)*xpi)/Beta), i=omegaBegin,omegaEnd) /)
1191   END IF
1192   MALLOC(Gwtmp,(1:Nom))
1193   DO i = omegaBegin, omegaEnd
1194     iw = omegatmp(i)
1195     omdeltabis = iw*deltabis
1196     Gwtmp(i)=CMPLX(0.d0,0.d0,8)
1197     DO j=2, Lspline ! We impose  L+1 = Nom
1198       iwtau = CMPLX(0.d0,omdeltabis*DBLE(j-1),8)
1199       Gwtmp(i) = Gwtmp(i) + EXP(iwtau) * CMPLX((X2(j+1) + X2(j-1))-2.d0*X2(j),0.d0,8)
1200     END DO
1201     Gwtmp(i) = Gwtmp(i)/CMPLX(((iw*iw)*(iw*iw)*deltabis),0.d0,8) &
1202 
1203               + CMPLX( &
1204                 ( ((X2(2)-X2(1))+(X2(Lspline+1)-X2(Lspline)))/((iw*iw)*deltabis) -this%Mk(2) ) &
1205                 /(iw*iw) &
1206                , &
1207                 (this%Mk(1)-this%Mk(3)/(iw*iw))/iw &
1208                , 8) 
1209               !+ CMPLX( (X2(2)-X2(1))+(X2(Lspline+1)-X2(Lspline)), 0.d0, 8 ) ) &
1210               !   / (((iw*iw)*(iw*iw))*CMPLX(deltabis,0.d0,8)) &
1211               !- CMPLX(this%Mk(1),0.d0,8)/iw  &
1212               !+ CMPLX(this%Mk(2),0.d0,8)/(iw*iw) &
1213               !- CMPLX(this%Mk(3),0.d0,8)/((iw*iw)*iw) 
1214     !IF ( this%rank .EQ. 0 )  write(12819,*) iw,gwtmp(i)
1215   END DO
1216   FREE(omegatmp)
1217   !call flush(12819)
1218   FREE(X2)
1219   IF ( this%have_MPI .EQV. .TRUE. ) THEN
1220 #ifdef HAVE_MPI
1221     CALL MPI_ALLGATHERV(MPI_IN_PLACE, 0, MPI_DATATYPE_NULL, &
1222                       Gwtmp  , counts, displs, &
1223                       MPI_DOUBLE_COMPLEX, this%MY_COMM, residu)
1224 #endif
1225     FREE(counts)
1226     FREE(displs)
1227   END IF
1228   IF ( PRESENT(Gomega) ) THEN
1229     Gomega = Gwtmp
1230   END IF
1231   this%oper_w = Gwtmp
1232   this%setW = .TRUE.
1233   FREE(Gwtmp)
1234 END SUBROUTINE GreenHyb_forFourier

ABINIT/m_GreenHyb/GreenHyb_getHybrid [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyb_getHybrid

FUNCTION

  reduce green function

COPYRIGHT

  Copyright (C) 2013-2018 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

  this=Green

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

598 SUBROUTINE GreenHyb_getHybrid(this)
599 
600 !Arguments ------------------------------------
601 
602 !This section has been created automatically by the script Abilint (TD).
603 !Do not modify the following lines by hand.
604 #undef ABI_FUNC
605 #define ABI_FUNC 'GreenHyb_getHybrid'
606 !End of the abilint section
607 
608   TYPE(GreenHyb), INTENT(INOUT) :: this
609 
610   IF ( this%set .EQV. .FALSE. ) &
611     CALL ERROR("GreenHyb_getHybrid : green operator not set          ")
612 
613   SELECT CASE(this%iTech)
614   CASE (GREENHYB_TAU)
615     this%oper = -(this%oper * this%inv_beta) / (DBLE(this%measurements) * this%delta_t)
616     this%setT = .TRUE.
617   CASE (GREENHYB_OMEGA)
618     this%oper_w = -(this%oper_w * this%inv_beta) / (DBLE(this%measurements) * this%delta_t)
619     this%setW = .TRUE.
620     CALL GreenHyb_backFourier(this)
621   END SELECT
622 
623 END SUBROUTINE GreenHyb_getHybrid

ABINIT/m_GreenHyb/GreenHyb_init [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyb_init

FUNCTION

  Initialize and allocate

COPYRIGHT

  Copyright (C) 2013-2018 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

  this=Green
  samples=imaginary time slices
  beta=inverse temperature
  iTech=SHOULD NOT BE USED => BUGGY
  MY_COMM=mpi_communicator

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

138 SUBROUTINE GreenHyb_init(this, samples, beta,iTech,MY_COMM)
139 
140 
141 !This section has been created automatically by the script Abilint (TD).
142 !Do not modify the following lines by hand.
143 #undef ABI_FUNC
144 #define ABI_FUNC 'GreenHyb_init'
145 !End of the abilint section
146 
147 
148 #ifdef HAVE_MPI1
149 include 'mpif.h'
150 #endif
151 !Arguments ------------------------------------
152   TYPE(GreenHyb)     , INTENT(INOUT) :: this
153   INTEGER         , INTENT(IN   ) :: samples
154   DOUBLE PRECISION, INTENT(IN   ) :: beta
155   !INTEGER          , INTENT(IN   ) :: Wmax
156   INTEGER, OPTIONAL, INTENT(IN   ) :: iTech
157   INTEGER, OPTIONAL, INTENT(IN   ) :: MY_COMM
158 !Local variables ------------------------------
159   INTEGER                         :: sp1
160   DOUBLE PRECISION                :: dt
161 #ifdef HAVE_MPI
162   INTEGER          :: ierr
163 #endif
164 
165   IF ( PRESENT(MY_COMM)) THEN
166 #ifdef HAVE_MPI
167     this%have_MPI = .TRUE.
168     this%MY_COMM = MY_COMM
169     CALL MPI_Comm_rank(this%MY_COMM, this%rank, ierr)
170     CALL MPI_Comm_size(this%MY_COMM, this%size, ierr)
171 #else
172     CALL WARN("GreenHyb_init : MPI is not used                                    ")
173     this%have_MPI = .FALSE.
174     this%MY_COMM = -1
175     this%rank = 0
176     this%size = 1
177 #endif
178   ELSE
179     this%have_MPI = .FALSE.
180     this%MY_COMM = -1
181     this%rank = 0
182     this%size = 1
183   END IF
184 
185   sp1             = samples + 1
186   this%samples      = sp1
187   this%measurements = 0
188   this%beta         = beta
189   this%inv_beta     = 1.d0 / beta
190   this%inv_dt       = DBLE(samples) * this%inv_beta
191   dt              = 1.d0 / this%inv_dt
192   this%delta_t      = dt
193   !this%Wmax         = Wmax
194   this%Wmax         = -1
195   FREEIF(this%oper)
196   MALLOC(this%oper,(1:sp1))
197   ! If we want to measure in frequences
198   ! let assume we first have "samples" frequences
199   IF ( PRESENT(iTech) ) THEN
200     this%iTech = iTech
201     SELECT CASE (this%iTech)
202     CASE (GREENHYB_TAU)  ! omega
203       this%iTech = GREENHYB_TAU
204     CASE (GREENHYB_OMEGA)  ! omega
205       this%Wmax = samples
206       FREEIF(this%oper_w)
207       MALLOC(this%oper_w,(1:this%Wmax))
208       FREEIF(this%oper_w_old)
209       MALLOC(this%oper_w_old,(1:this%Wmax))
210       this%oper_w     = CMPLX(0.d0,0.d0,8)
211       this%oper_w_old = CMPLX(0.d0,0.d0,8)
212       FREEIF(this%omega)
213       MALLOC(this%omega,(1:this%Wmax))
214       this%omega = (/ ((2.d0 * DBLE(sp1) - 1.d0)*ACOS(-1.d0)*this%inv_beta, sp1=1, this%Wmax) /)
215     END SELECT
216   ELSE
217     this%iTech = GREENHYB_TAU
218   END IF
219   ! end if
220   !CALL Vector_init(this%oper_old,10000)
221   !CALL VectorInt_init(this%index_old,10000)
222   CALL MapHyb_init(this%this,10000)
223 
224   this%oper       = 0.d0
225   this%set        = .TRUE.
226   this%factor     = 1
227   this%setMk      = 0
228   this%Mk         = 0.d0
229 END SUBROUTINE GreenHyb_init

ABINIT/m_GreenHyb/GreenHyb_measHybrid [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyb_measHybrid

FUNCTION

  Measure Green function

COPYRIGHT

  Copyright (C) 2013-2018 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

  this=Green
  Mthis=M this for the current flavor
  ListCdagC_1=list of all creator and annhilator operators
  updated=should we accumulate or not

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

436 SUBROUTINE GreenHyb_measHybrid(this, Mthis, ListCdagC_1, updated)
437 
438 !Arguments ------------------------------------
439 
440 !This section has been created automatically by the script Abilint (TD).
441 !Do not modify the following lines by hand.
442 #undef ABI_FUNC
443 #define ABI_FUNC 'GreenHyb_measHybrid'
444 !End of the abilint section
445 
446   TYPE(GreenHyb)    , INTENT(INOUT) :: this
447   TYPE(MatrixHyb)   , INTENT(IN   ) :: Mthis
448   TYPE(ListcdagC), INTENT(IN   ) :: ListCdagC_1
449   LOGICAL        , INTENT(IN   ) :: updated
450 !Local variables ------------------------------
451   INTEGER                        :: iC
452   INTEGER                        :: iCdag
453   INTEGER                        :: tail
454   !INTEGER                        :: index
455   INTEGER                        :: idx_old
456   INTEGER                        :: old_size
457   INTEGER                        :: omegaSamples
458   INTEGER                        :: iomega
459   DOUBLE PRECISION               :: pi_invBeta
460   DOUBLE PRECISION               :: mbeta_two
461   DOUBLE PRECISION               :: beta 
462   DOUBLE PRECISION               :: beta_tc
463   DOUBLE PRECISION               :: tcbeta_tc
464   DOUBLE PRECISION               :: inv_dt 
465   DOUBLE PRECISION               :: tC
466   DOUBLE PRECISION               :: tCdag
467   DOUBLE PRECISION               :: time
468   DOUBLE PRECISION               :: signe
469   DOUBLE PRECISION               :: argument
470   !DOUBLE PRECISION               :: taupi_invbeta
471   !COMPLEX(KIND=8)                   :: cargument
472   !COMPLEX(2*8)                   :: base_exp
473   !COMPLEX(2*8)                   :: increm_exp
474 
475   IF ( this%set .EQV. .FALSE. ) &
476     CALL ERROR("GreenHyb_measHybrid : green operator not set         ")
477 
478   tail = ListCdagC_1%tail
479   IF ( tail .NE. Mthis%tail ) &
480     CALL ERROR("GreenHyb_measHybrid : ListCdagC & M unconsistent     ")
481 
482   IF ( updated .EQV. .TRUE. ) THEN ! NEW change in the configuration
483     ! FIXME SHOULD be much more faster
484     
485       old_size = this%this%tail
486     SELECT CASE(this%iTech)
487     CASE (GREENHYB_TAU)
488       argument = DBLE(this%factor)
489       DO iC = 1, old_size
490         this%oper(this%this%listINT(iC)) = this%oper(this%this%listINT(iC)) + this%this%listDBLE(iC) * argument
491       END DO
492       this%measurements = this%measurements + this%factor
493   
494       CALL MapHyb_setSize(this%this,tail*tail)
495       this%factor = 1
496       idx_old = 0
497       beta   =  this%beta
498       mbeta_two = -(beta*0.5d0)
499       inv_dt =  this%inv_dt
500       ! WARNING time is not the time but just a temporary variable.
501       ! Index Time has been calculated previously and is in mat_tau
502       DO iC  = 1, tail
503         tC   = ListCdagC_1%list(iC,C_)
504         beta_tc = beta - tC
505         tcbeta_tc = tC * beta_tc
506         DO iCdag = 1, tail
507           tCdag  = ListCdagC_1%list(iCdag,Cdag_)
508           time = tcbeta_tc - tCdag*beta_tc
509   
510           !signe = SIGN(1.d0,time)
511           !time = time + (signe-1.d0)*mbeta_two
512     
513           !signe = signe * SIGN(1.d0,beta-tC)
514 
515           !signe = SIGN(1.d0,time) * SIGN(1.d0,beta-tC)
516           signe = SIGN(1.d0,time)
517   
518           argument = signe*Mthis%mat(iCdag,iC)
519   
520           !index = INT( ( time * inv_dt ) + 1.5d0 )
521           !IF (index .NE. Mthis%mat_tau(iCdag,iC)) THEN
522           !  WRITE(*,*) index, Mthis%mat_tau(iCdag,iC)
523           !!  CALL ERROR("Plantage")
524           !END IF
525   
526           idx_old = idx_old + 1
527           this%this%listDBLE(idx_old) = argument
528           !this%this%listINT(idx_old)  = index
529           this%this%listINT(idx_old)  = Mthis%mat_tau(iCdag,iC)
530         END DO
531       END DO
532     CASE (GREENHYB_OMEGA)
533       omegaSamples = this%Wmax
534       argument = DBLE(this%factor)
535       DO iomega = 1, omegaSamples
536         this%oper_w(iomega) = this%oper_w(iomega) + this%oper_w_old(iomega) * argument
537       END DO
538       this%measurements = this%measurements + this%factor
539 
540       this%factor = 1
541       beta   =  this%beta
542       mbeta_two = -(beta*0.5d0)
543       pi_invBeta = ACOS(-1.d0)/beta
544       DO iC  = 1, tail
545         tC   = ListCdagC_1%list(iC,C_)
546         DO iCdag = 1, tail
547           tCdag  = ListCdagC_1%list(iCdag,Cdag_)
548           time = tC - tCdag
549 
550           signe = SIGN(1.d0,time)
551           time = time + (signe-1.d0)*mbeta_two
552           signe = signe * SIGN(1.d0,beta-tC)
553           argument = signe*Mthis%mat(iCdag,iC)
554 
555           DO iomega = 1, omegaSamples
556             !this%oper_w_old(iomega) = Mthis%mat_tau(iCdag,iC)*CMPLX(0.d0,argument)
557             this%oper_w_old(iomega) = EXP(CMPLX(0.d0,this%omega(iomega)*time,8))*CMPLX(0.d0,argument,8)
558           END DO
559         END DO
560       END DO
561     END SELECT
562   ELSE
563     this%factor = this%factor + 1
564   END IF
565 END SUBROUTINE GreenHyb_measHybrid

ABINIT/m_GreenHyb/GreenHyb_print [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyb_print

FUNCTION

  print Green function

COPYRIGHT

  Copyright (C) 2013-2018 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

  this=Green
  ostream=file stream

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

1268 SUBROUTINE GreenHyb_print(this, ostream)
1269 
1270 !Arguments ------------------------------------
1271 
1272 !This section has been created automatically by the script Abilint (TD).
1273 !Do not modify the following lines by hand.
1274 #undef ABI_FUNC
1275 #define ABI_FUNC 'GreenHyb_print'
1276 !End of the abilint section
1277 
1278   TYPE(GreenHyb), INTENT(IN) :: this
1279   INTEGER, OPTIONAL , INTENT(IN) :: ostream
1280 !Local variables ------------------------------
1281   INTEGER                        :: ostream_val
1282   INTEGER                        :: i
1283   INTEGER                        :: samples
1284   
1285 
1286   IF ( this%set .EQV. .FALSE. ) &
1287     CALL ERROR("GreenHyb_print : green this%operator not set              ")
1288 
1289   IF ( PRESENT(ostream) ) THEN
1290     ostream_val = ostream
1291   ELSE
1292     ostream_val = 66
1293     OPEN(UNIT=ostream_val,FILE="Green.dat")
1294   END IF
1295 
1296   samples =  this%samples 
1297 
1298   DO i = 1, samples
1299     WRITE(ostream_val,*) DBLE(i-1)*this%delta_t, this%oper(i)
1300   END DO
1301 
1302   IF ( .NOT. PRESENT(ostream) ) &
1303     CLOSE(ostream_val)
1304 END SUBROUTINE GreenHyb_print

ABINIT/m_GreenHyb/GreenHyb_reset [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyb_reset

FUNCTION

  reset green function

COPYRIGHT

  Copyright (C) 2013-2018 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

  this=Green

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

320 SUBROUTINE GreenHyb_reset(this)
321 
322 !Arguments ------------------------------------
323 
324 !This section has been created automatically by the script Abilint (TD).
325 !Do not modify the following lines by hand.
326 #undef ABI_FUNC
327 #define ABI_FUNC 'GreenHyb_reset'
328 !End of the abilint section
329 
330   TYPE(GreenHyb)     , INTENT(INOUT) :: this
331 
332   CALL GreenHyb_clear(this)
333   this%setMk        = 0
334   this%Mk           = 0.d0
335   this%setT         = .FALSE.
336   this%setW         = .FALSE.
337 END SUBROUTINE GreenHyb_reset

ABINIT/m_GreenHyb/GreenHyb_setMoments [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyb_setMoments

FUNCTION

  Compute full moments

COPYRIGHT

  Copyright (C) 2013-2018 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

  this=Green
  u1=interaction energi like
  u2=idem order2

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

762 SUBROUTINE GreenHyb_setMoments(this,u1,u2)
763 
764 !Arguments ------------------------------------
765 
766 !This section has been created automatically by the script Abilint (TD).
767 !Do not modify the following lines by hand.
768 #undef ABI_FUNC
769 #define ABI_FUNC 'GreenHyb_setMoments'
770 !End of the abilint section
771 
772   TYPE(GreenHyb)  , INTENT(INOUT) :: this
773   DOUBLE PRECISION, INTENT(IN   ) :: u1
774   DOUBLE PRECISION, INTENT(IN   ) :: u2
775   
776   this%Mk(1) = -1.d0
777   this%Mk(3) = this%Mk(3) - 2.d0*(this%Mk(2)*u1)
778   this%Mk(2) = this%Mk(2) + u1
779   this%Mk(3) = this%Mk(3) - u2
780 
781   this%setMk = this%setMk + 1
782 
783 END SUBROUTINE GreenHyb_setMoments

ABINIT/m_GreenHyb/GreenHyb_setMuD1 [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyb_setMuD1

FUNCTION

  Set first moments for G

COPYRIGHT

  Copyright (C) 2013-2018 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

  this=Green
  mu=energy level (irrespectige with fermi level)
  d1=firt moment of hybridization function

OUTPUT

  argout(sizeout)=description

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

710 SUBROUTINE GreenHyb_setMuD1(this,mu,d1)
711 
712 !Arguments ------------------------------------
713 
714 !This section has been created automatically by the script Abilint (TD).
715 !Do not modify the following lines by hand.
716 #undef ABI_FUNC
717 #define ABI_FUNC 'GreenHyb_setMuD1'
718 !End of the abilint section
719 
720   TYPE(GreenHyb)  , INTENT(INOUT) :: this
721   DOUBLE PRECISION, INTENT(IN   ) :: mu
722   DOUBLE PRECISION, INTENT(IN   ) :: d1
723 
724   this%Mk(3) = -d1-(mu*mu)
725   this%Mk(2) = -mu
726   this%setMk = this%setMk + 1
727 END SUBROUTINE GreenHyb_setMuD1

ABINIT/m_GreenHyb/GreenHyb_setN [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyb_setN

FUNCTION

  impose number of electrons for this flavor

COPYRIGHT

  Copyright (C) 2013-2018 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

  this=Green
  N=number of electrons

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

657 SUBROUTINE GreenHyb_setN(this,N)
658 
659 !Arguments ------------------------------------
660 
661 !This section has been created automatically by the script Abilint (TD).
662 !Do not modify the following lines by hand.
663 #undef ABI_FUNC
664 #define ABI_FUNC 'GreenHyb_setN'
665 !End of the abilint section
666 
667   TYPE(GreenHyb)     , INTENT(INOUT) :: this
668   DOUBLE PRECISION, INTENT(IN   ) :: N
669 
670   IF ( this%set .EQV. .FALSE. ) &
671     CALL ERROR("GreenHyb_setN: green this%operator not set                ")
672   this%oper(1) = N - 1.d0
673   this%oper(this%samples) = - N
674 END SUBROUTINE GreenHyb_setN

ABINIT/m_GreenHyb/GreenHyb_setOperW [ Functions ]

[ Top ] [ Functions ]

NAME

  GreenHyb_setOperW

FUNCTION

  set Green function in frequencies

COPYRIGHT

  Copyright (C) 2013-2018 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

  this=Green
  Gomega=Input values

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

371 SUBROUTINE GreenHyb_setOperW(this, Gomega)
372 
373 !Arguments ------------------------------------
374 
375 !This section has been created automatically by the script Abilint (TD).
376 !Do not modify the following lines by hand.
377 #undef ABI_FUNC
378 #define ABI_FUNC 'GreenHyb_setOperW'
379 !End of the abilint section
380 
381   TYPE(GreenHyb)          , INTENT(INOUT) :: this
382   COMPLEX(KIND=8), DIMENSION(:), INTENT(IN   ) :: Gomega
383 !Loval variables ------------------------------
384   INTEGER :: tail
385 
386   tail = SIZE(Gomega)
387   IF ( .NOT. this%set ) &
388     CALL ERROR("GreenHyb_setOperW : Uninitialized GreenHyb structure")
389   IF ( ALLOCATED(this%oper_w) ) THEN
390     IF ( SIZE(this%oper_w) .NE. tail ) THEN
391       FREE(this%oper_w)
392       MALLOC(this%oper_w,(1:tail))
393     END IF
394   ELSE
395     MALLOC(this%oper_w,(1:tail))
396   END IF
397   this%oper_w(:) = Gomega(:)
398   this%Wmax = tail
399   this%setW = .TRUE.
400 END SUBROUTINE GreenHyb_setOperW

m_GreenHyb/GreenHyb [ Types ]

[ Top ] [ m_GreenHyb ] [ Types ]

NAME

  GreenHyb

FUNCTION

  This structured datatype contains the necessary data

COPYRIGHT

  Copyright (C) 2013-2018 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

60 TYPE, PUBLIC :: GreenHyb
61   LOGICAL _PRIVATE :: set = .FALSE.
62   LOGICAL _PRIVATE :: setT = .FALSE.
63   LOGICAL _PRIVATE :: setW = .FALSE.
64   LOGICAL _PRIVATE :: have_MPI = .FALSE.
65   INTEGER _PRIVATE :: setMk = 0
66   INTEGER _PRIVATE :: samples
67   INTEGER _PRIVATE :: measurements    
68   INTEGER          :: factor
69   INTEGER _PRIVATE :: MY_COMM
70   INTEGER _PRIVATE :: size
71   INTEGER _PRIVATE :: rank
72   INTEGER _PRIVATE :: Wmax
73   INTEGER _PRIVATE :: iTech
74   DOUBLE PRECISION _PRIVATE :: beta
75   DOUBLE PRECISION _PRIVATE :: inv_beta
76   DOUBLE PRECISION _PRIVATE :: delta_t
77   DOUBLE PRECISION _PRIVATE :: inv_dt
78   DOUBLE PRECISION, ALLOCATABLE , DIMENSION(:)            :: oper
79   DOUBLE PRECISION, ALLOCATABLE , DIMENSION(:)   _PRIVATE :: omega
80   DOUBLE PRECISION              , DIMENSION(1:3) _PRIVATE :: Mk
81   COMPLEX(KIND=8)  , ALLOCATABLE, DIMENSION(:)   _PRIVATE :: oper_w 
82   COMPLEX(KIND=8)  , ALLOCATABLE, DIMENSION(:)   _PRIVATE :: oper_w_old
83   TYPE(MapHyb)          :: this
84 END TYPE GreenHyb