TABLE OF CONTENTS


ABINIT/m_MatrixHyb [ Modules ]

[ Top ] [ Modules ]

NAME

  m_MatrixHyb

FUNCTION

  Module to deals with a matrix (mainly used for M matrix in m_BathOperatoroffdiag).
  Perform varius operation on matrices.

COPYRIGHT

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

SOURCE

23 #include "defs.h"
24 MODULE m_MatrixHyb
25 USE m_Global
26 IMPLICIT NONE

ABINIT/m_MatrixHyb/MatrixHyb_assign [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_assign

FUNCTION

  assign this=matrix2

COPYRIGHT

  Copyright (C) 2013-2022 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=Matrix 1
  Matrix2=Matrix 2

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

333 SUBROUTINE MatrixHyb_assign(this, matrix)
334 
335 !Arguments ------------------------------------
336   TYPE(MatrixHyb), INTENT(INOUT) :: this
337   TYPE(MatrixHyb), INTENT(IN   ) :: matrix
338 !Local variables ------------------------------
339   INTEGER                        :: tail
340 
341   tail = matrix%tail
342   CALL MatrixHyb_setSize(this, tail)
343   this%mat(1:tail,1:tail) = matrix%mat(1:tail,1:tail)
344   IF ( this%iTech .NE. matrix%iTech ) & 
345     CALL ERROR("MatrixHyb_assign : not compatible matrices")
346   SELECT CASE(this%iTech)
347   CASE (GREENHYB_TAU)
348     this%mat_tau(1:tail,1:tail) = matrix%mat_tau(1:tail,1:tail)
349   CASE (GREENHYB_OMEGA)
350     IF ( this%Wmax .NE. matrix%Wmax ) &
351       CALL ERROR("MatrixHyb_assig : different numbers of omega")
352     this%mat_omega(1:this%Wmax,1:tail,1:tail) = &
353     matrix%mat_omega(1:this%Wmax,1:tail,1:tail)
354   END SELECT
355 
356 END SUBROUTINE MatrixHyb_assign

ABINIT/m_MatrixHyb/MatrixHyb_clear [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_clear

FUNCTION

  Clear this

COPYRIGHT

  Copyright (C) 2013-2022 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=matrix

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

300 SUBROUTINE MatrixHyb_clear(this)
301 
302 !Arguments ------------------------------------
303   TYPE(MatrixHyb), INTENT(INOUT) :: this
304   this%tail = 0 
305 END SUBROUTINE MatrixHyb_clear

ABINIT/m_MatrixHyb/MatrixHyb_destroy [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_destroy

FUNCTION

  Destroy

COPYRIGHT

  Copyright (C) 2013-2022 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=Matrix

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

702 SUBROUTINE MatrixHyb_destroy(this)
703 
704 !Arguments ------------------------------------
705   TYPE(MatrixHyb), INTENT(INOUT) :: this
706 
707   FREEIF(this%mat)
708   SELECT CASE(this%iTech)
709   CASE (GREENHYB_TAU)
710     FREEIF(this%mat_tau)
711   CASE (GREENHYB_OMEGA)
712     FREEIF(this%mat_omega)
713   END SELECT
714 
715   this%tail     = 0
716   this%size     = 0
717   this%iTech    = -1
718 END SUBROUTINE MatrixHyb_destroy

ABINIT/m_MatrixHyb/MatrixHyb_enlarge [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_enlarge

FUNCTION

  This subroutine enlarges memory space

COPYRIGHT

  Copyright (C) 2013-2022 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=matrix
  size=new memory size

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

221 SUBROUTINE MatrixHyb_enlarge(this, size)
222 
223 !Arguments ------------------------------------
224   TYPE(MatrixHyb)     , INTENT(INOUT)          :: this
225   INTEGER, OPTIONAL, INTENT(IN   )          :: size
226 !Local variables ------------------------------
227   INTEGER                                   :: width
228   INTEGER                                   :: tail
229   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: this_temp 
230   INTEGER         , ALLOCATABLE, DIMENSION(:,:) :: this_temp_tau
231   COMPLEX(KIND=8) , ALLOCATABLE, DIMENSION(:,:,:) :: this_temp_omega
232   INTEGER                                   :: size_val
233 
234   IF ( ALLOCATED(this%mat) ) THEN
235     FREEIF(this_temp)
236     FREEIF(this_temp_tau)
237     width = this%size
238     tail  = this%tail 
239     size_val = width
240     IF ( PRESENT(size) ) size_val = size 
241 
242 !   change size of mat
243     MALLOC(this_temp,(1:tail,1:tail))
244     this_temp(1:tail,1:tail) = this%mat(1:tail,1:tail)
245     FREE(this%mat)
246     this%size = width + size_val
247     MALLOC(this%mat,(1:this%size,1:this%size))
248     this%mat(1:tail,1:tail) = this_temp(1:tail,1:tail)
249     FREE(this_temp)
250     SELECT CASE(this%iTech)
251     CASE (GREENHYB_TAU)
252 
253 !   change size of mat_tau
254       MALLOC(this_temp_tau,(1:tail,1:tail))
255       this_temp_tau(1:tail,1:tail) = this%mat_tau(1:tail,1:tail)
256       FREE(this%mat_tau)
257       MALLOC(this%mat_tau,(1:this%size,1:this%size))
258       this%mat_tau(1:tail,1:tail) = this_temp_tau(1:tail,1:tail)
259       FREE(this_temp_tau)
260     CASE (GREENHYB_OMEGA)
261 
262 !   change size of mat_omega
263       MALLOC(this_temp_omega,(1:this%Wmax,1:tail,1:tail))
264       this_temp_omega(1:this%Wmax,1:tail,1:tail) = this%mat_omega(1:this%Wmax,1:tail,1:tail)
265       FREE(this%mat_omega)
266       MALLOC(this%mat_omega,(1:this%Wmax,1:this%size,1:this%size))
267       this%mat_omega(1:this%Wmax,1:tail,1:tail) = this_temp_omega(1:this%Wmax,1:tail,1:tail)
268       FREE(this_temp_omega)
269     END SELECT
270   ELSE
271     CALL MatrixHyb_init(this, this%iTech, size=Global_SIZE, Wmax=this%Wmax)
272   END IF
273 END SUBROUTINE MatrixHyb_enlarge

ABINIT/m_MatrixHyb/MatrixHyb_getDet [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_getDet

FUNCTION

  Just get the determinant 

COPYRIGHT

  Copyright (C) 2013-2022 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=matrix

OUTPUT

  det=determinant

SIDE EFFECTS

NOTES

SOURCE

598 SUBROUTINE MatrixHyb_getDet(this,det)
599 
600 !Arguments ------------------------------------
601   TYPE(MatrixHyb) , INTENT(INOUT) :: this
602   DOUBLE PRECISION, INTENT(  OUT) :: det
603 
604   IF ( this%tail .EQ. 0 ) THEN
605     det = 1.d0
606     RETURN
607   END IF
608   CALL MatrixHyb_LU(this, determinant=det)
609 END SUBROUTINE MatrixHyb_getDet

ABINIT/m_MatrixHyb/MatrixHyb_init [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_init

FUNCTION

  initialize

COPYRIGHT

  Copyright (C) 2013-2022 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=matrix
  iTech=DO NOT USE => BUG
  size=memory size for initialization
  Wmax=maximum frequency number

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

117 SUBROUTINE MatrixHyb_init(this, iTech, size, Wmax)
118 
119 !Arguments ------------------------------------
120   TYPE(MatrixHyb)     , INTENT(INOUT) :: this
121   INTEGER             , INTENT(IN   ) :: iTech
122   INTEGER, OPTIONAL, INTENT(IN   ) :: size
123   INTEGER, OPTIONAL, INTENT(IN   ) :: Wmax
124 !Local variables ------------------------------
125   INTEGER                          :: size_val
126 
127   size_val = Global_SIZE
128   IF ( PRESENT(size) ) size_val = size
129   this%size = size_val
130   FREEIF(this%mat)
131   MALLOC(this%mat,(1:size_val,1:size_val))
132   this%tail  = 0 
133   this%mat   = 0.d0
134   this%iTech = iTech
135   SELECT CASE(this%iTech)
136   CASE (GREENHYB_TAU)
137     FREEIF(this%mat_tau)
138     MALLOC(this%mat_tau,(1:size_val,1:size_val))
139   CASE (GREENHYB_OMEGA)
140     IF ( PRESENT(Wmax) .AND. Wmax .GT. 0 ) THEN
141       FREEIF(this%mat_omega)
142       MALLOC(this%mat_omega,(1:Wmax,1:size_val,1:size_val))
143       this%Wmax=Wmax
144     ELSE
145       CALL ERROR("MatrixHyb_init : Missing argument Wmax for measurement in omega")
146     END IF
147   CASE DEFAULT
148       CALL WARNALL("MatrixHyb_init : Wrong input argument iTech")
149   END SELECT
150 END SUBROUTINE MatrixHyb_init

ABINIT/m_MatrixHyb/MatrixHyb_inverse [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_inverse

FUNCTION

  inverse the matrix and compute the determinant

COPYRIGHT

  Copyright (C) 2013-2022 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=this

OUTPUT

  determinant=determinant of the this before inversion

SIDE EFFECTS

NOTES

SOURCE

384 SUBROUTINE MatrixHyb_inverse(this,determinant)
385 
386 !Arguments ------------------------------------
387   TYPE(MatrixHyb), INTENT(INOUT) :: this
388   DOUBLE PRECISION, OPTIONAL, INTENT(OUT) :: determinant
389   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: invMatrix
390 !Local variables ------------------------------
391   !DOUBLE PRECISION :: reste
392   INTEGER :: ligne
393   INTEGER :: ligne_virtuelle
394   INTEGER :: colonne
395   INTEGER :: sys
396   INTEGER :: sys_vir
397   INTEGER :: tail
398   INTEGER, DIMENSION(:), ALLOCATABLE :: pivot 
399   DOUBLE PRECISION :: det
400 
401   tail = this%tail
402   IF ( tail .EQ. 0 ) THEN
403     IF ( PRESENT(determinant) ) determinant = 1.d0
404     RETURN
405   END IF
406 
407   MALLOC(invMatrix,(1:tail,1:tail))
408 
409   FREEIF(pivot)
410   MALLOC(pivot,(1:tail))
411 
412   CALL MatrixHyb_LU(this, pivot=pivot, determinant=det)
413 
414   !det = 1.d0
415   DO sys = 1, tail
416     sys_vir = pivot(sys)
417 !!    DO ligne=1,tail
418 !!      ligne_virtuelle = pivot(ligne)
419 !!      reste = 0.d0
420 !!      DO colonne = 1, ligne-1
421 !!        reste = reste + this%mat(ligne_virtuelle,colonne)*invMatrix%mat(colonne,sys_vir)
422 !!      END DO
423 !!      IF ( ligne .EQ. sys ) THEN
424 !!        invMatrix%mat(ligne,sys_vir) = 1.d0 - reste
425 !!      ELSE
426 !!        invMatrix%mat(ligne,sys_vir) = 0.d0 - reste
427 !!      END IF
428 !!    END DO
429 !!    
430 !!    DO ligne = tail, 1, -1
431 !!      ligne_virtuelle=pivot(ligne)
432 !!      reste = 0.d0
433 !!        DO colonne = ligne+1, tail
434 !!          reste= reste + this%mat(ligne_virtuelle,colonne)*invMatrix%mat(colonne,sys_vir)
435 !!        END DO
436 !!        invMatrix%mat(ligne, sys_vir) = (invMatrix%mat(ligne, sys_vir)-reste)/this%mat(ligne_virtuelle,ligne)
437 !!    END DO
438     !det = det*this%mat(sys_vir,sys)
439     invMatrix(:,sys_vir) = 0.d0
440     invMatrix(sys,sys_vir) = 1.d0
441     DO colonne = 1, tail-1
442       DO ligne=colonne+1,tail
443         ligne_virtuelle = pivot(ligne)
444         invMatrix(ligne, sys_vir) = invMatrix(ligne,sys_vir) &
445                                                 - this%mat(ligne_virtuelle,colonne)  &
446                                                 * invMatrix(colonne,sys_vir)
447       END DO
448     END DO
449 
450     DO colonne = tail, 1, -1
451       invMatrix(colonne, sys_vir) = invMatrix(colonne, sys_vir) / this%mat(pivot(colonne),colonne)
452       DO ligne = 1, colonne-1
453         ligne_virtuelle = pivot(ligne)
454         invMatrix(ligne, sys_vir) = invMatrix(ligne,sys_vir) &
455                                                  - this%mat(ligne_virtuelle,colonne) &
456                                                  * invMatrix(colonne, sys_vir)
457       END DO
458     END DO
459 
460   END DO
461   FREE(pivot)
462   this%mat(1:tail,1:tail) = invMatrix(1:tail,1:tail)
463   IF ( PRESENT(determinant) ) THEN
464     determinant = det
465   END IF
466   FREE(invMatrix)
467 END SUBROUTINE MatrixHyb_inverse

ABINIT/m_MatrixHyb/MatrixHyb_LU [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_LU

FUNCTION

  LU decomposition

COPYRIGHT

  Copyright (C) 2013-2022 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=matrix

OUTPUT

  pivot=gauss pivot
  determinant=determinant of the this

SIDE EFFECTS

NOTES

SOURCE

496 SUBROUTINE MatrixHyb_LU(this,pivot,determinant)
497 
498 !Arguments ------------------------------------
499   TYPE(MatrixHyb), INTENT(INOUT) :: this
500   INTEGER, DIMENSION(:), ALLOCATABLE, OPTIONAL, INTENT(INOUT) :: pivot
501   DOUBLE PRECISION, OPTIONAL, INTENT(OUT) :: determinant
502 !Local variables ------------------------------
503   INTEGER :: ligne
504   INTEGER :: colonne 
505   INTEGER :: colonne_de_1_ligne
506   INTEGER :: max_ind_lig
507   INTEGER :: ligne_virtuelle
508   INTEGER :: tail
509   INTEGER, DIMENSION(:), ALLOCATABLE :: pivot_tmp
510   DOUBLE PRECISION :: max_col
511   DOUBLE PRECISION :: inverse_pivot
512   DOUBLE PRECISION :: coef
513   DOUBLE PRECISION :: det
514   DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: mat_tmp
515 
516   tail = this%tail
517   det = 1.d0
518   FREEIF(pivot_tmp)
519   MALLOC(pivot_tmp,(1:tail))
520 
521   DO ligne=1, tail
522     pivot_tmp(ligne)=ligne
523   END DO
524 
525   MALLOC(mat_tmp,(1:tail,1:tail))
526   mat_tmp(1:tail,1:tail) = this%mat(1:tail,1:tail)
527   DO colonne = 1, tail-1 
528     max_col = ABS(mat_tmp(pivot_tmp(colonne),colonne))
529     max_ind_lig = colonne
530 
531     DO ligne = colonne+1, tail 
532       ligne_virtuelle = pivot_tmp(ligne)
533       IF ( ABS( mat_tmp(ligne_virtuelle,colonne)).GT.max_col) then
534         max_col = ABS(mat_tmp(ligne_virtuelle,colonne))
535         max_ind_lig = ligne
536       ENDIF
537     END DO
538 
539     ligne              = pivot_tmp(colonne) 
540     pivot_tmp(colonne)     = pivot_tmp(max_ind_lig)
541     pivot_tmp(max_ind_lig) = ligne
542     IF ( pivot_tmp(colonne) .NE. pivot_tmp(max_ind_lig) ) det = det * (-1.d0)
543     inverse_pivot=1.d0 / mat_tmp(pivot_tmp(colonne),colonne)
544     det = det * mat_tmp(pivot_tmp(colonne),colonne)
545     DO ligne = colonne+1, tail 
546       ligne_virtuelle = pivot_tmp(ligne)
547       coef = mat_tmp(ligne_virtuelle,colonne)*inverse_pivot
548       mat_tmp(ligne_virtuelle,colonne) = coef
549       DO colonne_de_1_ligne = colonne+1, tail 
550         mat_tmp(ligne_virtuelle,colonne_de_1_ligne)= mat_tmp(ligne_virtuelle,colonne_de_1_ligne)&
551                                          -coef * mat_tmp(pivot_tmp(colonne) ,colonne_de_1_ligne)
552       END DO
553     END DO
554   END DO
555   det = det * mat_tmp(pivot_tmp(tail),tail)
556   IF ( PRESENT(determinant) ) &
557     determinant = det
558   IF ( PRESENT(pivot) ) THEN
559     this%mat(1:tail,1:tail) = mat_tmp(1:tail,1:tail)
560     IF ( ALLOCATED(pivot) .AND. SIZE(pivot) .NE. tail ) THEN
561       FREE(pivot)
562       MALLOC(pivot,(1:tail))
563     ELSE IF ( .NOT. ALLOCATED(pivot) ) THEN
564       MALLOC(pivot,(1:tail))
565     END IF
566     pivot(1:tail)=pivot_tmp(1:tail)
567   END IF
568   FREE(mat_tmp)
569   FREE(pivot_tmp)
570 END SUBROUTINE MatrixHyb_LU

ABINIT/m_MatrixHyb/MatrixHyb_print [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_print

FUNCTION

  print Matrix

COPYRIGHT

  Copyright (C) 2013-2022 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=Matrix
  ostream=file stream
  opt_print=0 mat
            1 mat_tau
            2 both

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

640 SUBROUTINE MatrixHyb_print(this,ostream,opt_print)
641 
642 !Arguments ------------------------------------
643   TYPE(MatrixHyb), INTENT(IN) :: this
644   INTEGER, OPTIONAL, INTENT(IN) :: ostream
645   INTEGER, OPTIONAL, INTENT(IN) :: opt_print
646 !Local variables ------------------------------
647   INTEGER                       :: ostream_val
648   INTEGER                       :: opt_val
649   INTEGER                       :: it1
650   INTEGER                       :: it2
651   CHARACTER(LEN=4 )             :: size
652   CHARACTER(LEN=22)             :: string
653 
654   ostream_val = 6
655   opt_val=0
656   IF ( PRESENT(ostream) ) ostream_val = ostream
657   IF ( PRESENT(opt_print) ) opt_val = opt_print
658   WRITE(size,'(I4)') this%tail
659   IF ( MOD(opt_val,2) .EQ. 0 ) THEN
660     string ='(1x,A,1x,'//TRIM(ADJUSTL(size))//'(ES10.2,1x),A)'
661     WRITE(ostream_val,'(A)') "["
662     DO it1 = 1, this%tail
663       WRITE(ostream_val,string) "[",(/ (this%mat(it1,it2),it2=1,this%tail)  /)," ]"
664     END DO
665     WRITE(ostream_val,'(A)') "]"
666   END IF  
667   IF ( opt_val .GE.1 ) THEN
668     string ='(1x,A,1x,'//TRIM(ADJUSTL(size))//'(I4,1x),A)'
669     WRITE(ostream_val,'(A)') "["
670     DO it1 = 1, this%tail
671       WRITE(ostream_val,string) "[",(/ (this%mat_tau(it1,it2),it2=1,this%tail)  /),"]"
672     END DO
673     WRITE(ostream_val,'(A)') "]"
674   END IF
675 END SUBROUTINE MatrixHyb_print

ABINIT/m_MatrixHyb/MatrixHyb_setSize [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_setSize

FUNCTION

  impose size of the this

COPYRIGHT

  Copyright (C) 2013-2022 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=matrix
  new_tail=new size

OUTPUT

SIDE EFFECTS

NOTES

SOURCE

178 SUBROUTINE MatrixHyb_setSize(this,new_tail)
179 
180 !Arguments ------------------------------------
181   TYPE(MatrixHyb), INTENT(INOUT) :: this
182   INTEGER        , INTENT(IN   ) :: new_tail
183 !Local variables ------------------------------
184   INTEGER                        :: size
185 
186   IF ( .NOT. ALLOCATED(this%mat) ) &
187     CALL MatrixHyb_init(this,this%iTech,Wmax=this%Wmax)
188   size = this%size
189   IF( new_tail .GT. size ) THEN
190     CALL MatrixHyb_enlarge(this, MAX(Global_SIZE,new_tail-size))
191   END IF
192   this%tail = new_tail
193 END SUBROUTINE MatrixHyb_setSize  

m_MatrixHyb/MatrixHyb [ Types ]

[ Top ] [ m_MatrixHyb ] [ Types ]

NAME

  MatrixHyb

FUNCTION

  This structured datatype contains the necessary data

COPYRIGHT

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

47 TYPE, PUBLIC :: MatrixHyb
48   INTEGER _PRIVATE :: size
49   ! size of the matrix
50 
51   INTEGER          :: tail
52   ! the size of the matrix that is actually used.
53 
54   INTEGER _PRIVATE :: iTech = -1
55   ! precise if time or frequency values will be used.
56 
57   INTEGER _PRIVATE :: Wmax
58   ! size if the frequency grid for mat_omega
59 
60   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:)           :: mat
61   ! matrix of size "size"
62 
63   INTEGER         , ALLOCATABLE, DIMENSION(:,:)           :: mat_tau
64   ! matrix of size "size"
65 
66   COMPLEX(KIND=8) , ALLOCATABLE, DIMENSION(:,:,:)         :: mat_omega
67   ! array of size (Wmax, size, size)
68 
69 END TYPE MatrixHyb