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-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

29 #include "defs.h"
30 MODULE m_MatrixHyb
31 USE m_Global
32 IMPLICIT NONE

ABINIT/m_MatrixHyb/MatrixHyb_assign [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_assign

FUNCTION

  assign this=matrix2

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=Matrix 1
  Matrix2=Matrix 2

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

397 SUBROUTINE MatrixHyb_assign(this, matrix)
398 
399 !Arguments ------------------------------------
400 
401 !This section has been created automatically by the script Abilint (TD).
402 !Do not modify the following lines by hand.
403 #undef ABI_FUNC
404 #define ABI_FUNC 'MatrixHyb_assign'
405 !End of the abilint section
406 
407   TYPE(MatrixHyb), INTENT(INOUT) :: this
408   TYPE(MatrixHyb), INTENT(IN   ) :: matrix
409 !Local variables ------------------------------
410   INTEGER                        :: tail
411 
412   tail = matrix%tail
413   CALL MatrixHyb_setSize(this, tail)
414   this%mat(1:tail,1:tail) = matrix%mat(1:tail,1:tail)
415   IF ( this%iTech .NE. matrix%iTech ) & 
416     CALL ERROR("MatrixHyb_assign : not compatible matrices")
417   SELECT CASE(this%iTech)
418   CASE (GREENHYB_TAU)
419     this%mat_tau(1:tail,1:tail) = matrix%mat_tau(1:tail,1:tail)
420   CASE (GREENHYB_OMEGA)
421     IF ( this%Wmax .NE. matrix%Wmax ) &
422       CALL ERROR("MatrixHyb_assig : different numbers of omega")
423     this%mat_omega(1:this%Wmax,1:tail,1:tail) = &
424     matrix%mat_omega(1:this%Wmax,1:tail,1:tail)
425   END SELECT
426 
427 END SUBROUTINE MatrixHyb_assign

ABINIT/m_MatrixHyb/MatrixHyb_clear [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_clear

FUNCTION

  Clear this

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=matrix

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

351 SUBROUTINE MatrixHyb_clear(this)
352 
353 !Arguments ------------------------------------
354 
355 !This section has been created automatically by the script Abilint (TD).
356 !Do not modify the following lines by hand.
357 #undef ABI_FUNC
358 #define ABI_FUNC 'MatrixHyb_clear'
359 !End of the abilint section
360 
361   TYPE(MatrixHyb), INTENT(INOUT) :: this
362   this%tail = 0 
363 END SUBROUTINE MatrixHyb_clear

ABINIT/m_MatrixHyb/MatrixHyb_destroy [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_destroy

FUNCTION

  Destroy

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=Matrix

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

831 SUBROUTINE MatrixHyb_destroy(this)
832 
833 !Arguments ------------------------------------
834 
835 !This section has been created automatically by the script Abilint (TD).
836 !Do not modify the following lines by hand.
837 #undef ABI_FUNC
838 #define ABI_FUNC 'MatrixHyb_destroy'
839 !End of the abilint section
840 
841   TYPE(MatrixHyb), INTENT(INOUT) :: this
842 
843   FREEIF(this%mat)
844   SELECT CASE(this%iTech)
845   CASE (GREENHYB_TAU)
846     FREEIF(this%mat_tau)
847   CASE (GREENHYB_OMEGA)
848     FREEIF(this%mat_omega)
849   END SELECT
850 
851   this%tail     = 0
852   this%size     = 0
853   this%iTech    = -1
854 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-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=matrix
  size=new memory size

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

259 SUBROUTINE MatrixHyb_enlarge(this, size)
260 
261 !Arguments ------------------------------------
262 
263 !This section has been created automatically by the script Abilint (TD).
264 !Do not modify the following lines by hand.
265 #undef ABI_FUNC
266 #define ABI_FUNC 'MatrixHyb_enlarge'
267 !End of the abilint section
268 
269   TYPE(MatrixHyb)     , INTENT(INOUT)          :: this
270   INTEGER, OPTIONAL, INTENT(IN   )          :: size
271 !Local variables ------------------------------
272   INTEGER                                   :: width
273   INTEGER                                   :: tail
274   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: this_temp 
275   INTEGER         , ALLOCATABLE, DIMENSION(:,:) :: this_temp_tau
276   COMPLEX(KIND=8) , ALLOCATABLE, DIMENSION(:,:,:) :: this_temp_omega
277   INTEGER                                   :: size_val
278 
279   IF ( ALLOCATED(this%mat) ) THEN
280     FREEIF(this_temp)
281     FREEIF(this_temp_tau)
282     width = this%size
283     tail  = this%tail 
284     size_val = width
285     IF ( PRESENT(size) ) size_val = size 
286 
287 !   change size of mat
288     MALLOC(this_temp,(1:tail,1:tail))
289     this_temp(1:tail,1:tail) = this%mat(1:tail,1:tail)
290     FREE(this%mat)
291     this%size = width + size_val
292     MALLOC(this%mat,(1:this%size,1:this%size))
293     this%mat(1:tail,1:tail) = this_temp(1:tail,1:tail)
294     FREE(this_temp)
295     SELECT CASE(this%iTech)
296     CASE (GREENHYB_TAU)
297 
298 !   change size of mat_tau
299       MALLOC(this_temp_tau,(1:tail,1:tail))
300       this_temp_tau(1:tail,1:tail) = this%mat_tau(1:tail,1:tail)
301       FREE(this%mat_tau)
302       MALLOC(this%mat_tau,(1:this%size,1:this%size))
303       this%mat_tau(1:tail,1:tail) = this_temp_tau(1:tail,1:tail)
304       FREE(this_temp_tau)
305     CASE (GREENHYB_OMEGA)
306 
307 !   change size of mat_omega
308       MALLOC(this_temp_omega,(1:this%Wmax,1:tail,1:tail))
309       this_temp_omega(1:this%Wmax,1:tail,1:tail) = this%mat_omega(1:this%Wmax,1:tail,1:tail)
310       FREE(this%mat_omega)
311       MALLOC(this%mat_omega,(1:this%Wmax,1:this%size,1:this%size))
312       this%mat_omega(1:this%Wmax,1:tail,1:tail) = this_temp_omega(1:this%Wmax,1:tail,1:tail)
313       FREE(this_temp_omega)
314     END SELECT
315   ELSE
316     CALL MatrixHyb_init(this, this%iTech, size=Global_SIZE, Wmax=this%Wmax)
317   END IF
318 END SUBROUTINE MatrixHyb_enlarge

ABINIT/m_MatrixHyb/MatrixHyb_getDet [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_getDet

FUNCTION

  Just get the determinant 

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=matrix

OUTPUT

  det=determinant

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

701 SUBROUTINE MatrixHyb_getDet(this,det)
702 
703 !Arguments ------------------------------------
704 
705 !This section has been created automatically by the script Abilint (TD).
706 !Do not modify the following lines by hand.
707 #undef ABI_FUNC
708 #define ABI_FUNC 'MatrixHyb_getDet'
709 !End of the abilint section
710 
711   TYPE(MatrixHyb) , INTENT(INOUT) :: this
712   DOUBLE PRECISION, INTENT(  OUT) :: det
713 
714   IF ( this%tail .EQ. 0 ) THEN
715     det = 1.d0
716     RETURN
717   END IF
718   CALL MatrixHyb_LU(this, determinant=det)
719 END SUBROUTINE MatrixHyb_getDet

ABINIT/m_MatrixHyb/MatrixHyb_init [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_init

FUNCTION

  initialize

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=matrix
  iTech=DO NOT USE => BUG
  size=memory size for initialization
  Wmax=maximum frequency number

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

129 SUBROUTINE MatrixHyb_init(this, iTech, size, Wmax)
130 
131 !Arguments ------------------------------------
132 
133 !This section has been created automatically by the script Abilint (TD).
134 !Do not modify the following lines by hand.
135 #undef ABI_FUNC
136 #define ABI_FUNC 'MatrixHyb_init'
137 !End of the abilint section
138 
139   TYPE(MatrixHyb)     , INTENT(INOUT) :: this
140   INTEGER             , INTENT(IN   ) :: iTech
141   INTEGER, OPTIONAL, INTENT(IN   ) :: size
142   INTEGER, OPTIONAL, INTENT(IN   ) :: Wmax
143 !Local variables ------------------------------
144   INTEGER                          :: size_val
145 
146   size_val = Global_SIZE
147   IF ( PRESENT(size) ) size_val = size
148   this%size = size_val
149   FREEIF(this%mat)
150   MALLOC(this%mat,(1:size_val,1:size_val))
151   this%tail  = 0 
152   this%mat   = 0.d0
153   this%iTech = iTech
154   SELECT CASE(this%iTech)
155   CASE (GREENHYB_TAU)
156     FREEIF(this%mat_tau)
157     MALLOC(this%mat_tau,(1:size_val,1:size_val))
158   CASE (GREENHYB_OMEGA)
159     IF ( PRESENT(Wmax) .AND. Wmax .GT. 0 ) THEN
160       FREEIF(this%mat_omega)
161       MALLOC(this%mat_omega,(1:Wmax,1:size_val,1:size_val))
162       this%Wmax=Wmax
163     ELSE
164       CALL ERROR("MatrixHyb_init : Missing argument Wmax for measurement in omega")
165     END IF
166   CASE DEFAULT
167       CALL WARNALL("MatrixHyb_init : Wrong input argument iTech")
168   END SELECT
169 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-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=this

OUTPUT

  determinant=determinant of the this before inversion

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

461 SUBROUTINE MatrixHyb_inverse(this,determinant)
462 
463 !Arguments ------------------------------------
464 
465 !This section has been created automatically by the script Abilint (TD).
466 !Do not modify the following lines by hand.
467 #undef ABI_FUNC
468 #define ABI_FUNC 'MatrixHyb_inverse'
469 !End of the abilint section
470 
471   TYPE(MatrixHyb), INTENT(INOUT) :: this
472   DOUBLE PRECISION, OPTIONAL, INTENT(OUT) :: determinant
473   DOUBLE PRECISION, ALLOCATABLE, DIMENSION(:,:) :: invMatrix
474 !Local variables ------------------------------
475   !DOUBLE PRECISION :: reste
476   INTEGER :: ligne
477   INTEGER :: ligne_virtuelle
478   INTEGER :: colonne
479   INTEGER :: sys
480   INTEGER :: sys_vir
481   INTEGER :: tail
482   INTEGER, DIMENSION(:), ALLOCATABLE :: pivot 
483   DOUBLE PRECISION :: det
484 
485   tail = this%tail
486   IF ( tail .EQ. 0 ) THEN
487     IF ( PRESENT(determinant) ) determinant = 1.d0
488     RETURN
489   END IF
490 
491   MALLOC(invMatrix,(1:tail,1:tail))
492 
493   FREEIF(pivot)
494   MALLOC(pivot,(1:tail))
495 
496   CALL MatrixHyb_LU(this, pivot=pivot, determinant=det)
497 
498   !det = 1.d0
499   DO sys = 1, tail
500     sys_vir = pivot(sys)
501 !!    DO ligne=1,tail
502 !!      ligne_virtuelle = pivot(ligne)
503 !!      reste = 0.d0
504 !!      DO colonne = 1, ligne-1
505 !!        reste = reste + this%mat(ligne_virtuelle,colonne)*invMatrix%mat(colonne,sys_vir)
506 !!      END DO
507 !!      IF ( ligne .EQ. sys ) THEN
508 !!        invMatrix%mat(ligne,sys_vir) = 1.d0 - reste
509 !!      ELSE
510 !!        invMatrix%mat(ligne,sys_vir) = 0.d0 - reste
511 !!      END IF
512 !!    END DO
513 !!    
514 !!    DO ligne = tail, 1, -1
515 !!      ligne_virtuelle=pivot(ligne)
516 !!      reste = 0.d0
517 !!        DO colonne = ligne+1, tail
518 !!          reste= reste + this%mat(ligne_virtuelle,colonne)*invMatrix%mat(colonne,sys_vir)
519 !!        END DO
520 !!        invMatrix%mat(ligne, sys_vir) = (invMatrix%mat(ligne, sys_vir)-reste)/this%mat(ligne_virtuelle,ligne)
521 !!    END DO
522     !det = det*this%mat(sys_vir,sys)
523     invMatrix(:,sys_vir) = 0.d0
524     invMatrix(sys,sys_vir) = 1.d0
525     DO colonne = 1, tail-1
526       DO ligne=colonne+1,tail
527         ligne_virtuelle = pivot(ligne)
528         invMatrix(ligne, sys_vir) = invMatrix(ligne,sys_vir) &
529                                                 - this%mat(ligne_virtuelle,colonne)  &
530                                                 * invMatrix(colonne,sys_vir)
531       END DO
532     END DO
533 
534     DO colonne = tail, 1, -1
535       invMatrix(colonne, sys_vir) = invMatrix(colonne, sys_vir) / this%mat(pivot(colonne),colonne)
536       DO ligne = 1, colonne-1
537         ligne_virtuelle = pivot(ligne)
538         invMatrix(ligne, sys_vir) = invMatrix(ligne,sys_vir) &
539                                                  - this%mat(ligne_virtuelle,colonne) &
540                                                  * invMatrix(colonne, sys_vir)
541       END DO
542     END DO
543 
544   END DO
545   FREE(pivot)
546   this%mat(1:tail,1:tail) = invMatrix(1:tail,1:tail)
547   IF ( PRESENT(determinant) ) THEN
548     determinant = det
549   END IF
550   FREE(invMatrix)
551 END SUBROUTINE MatrixHyb_inverse

ABINIT/m_MatrixHyb/MatrixHyb_LU [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_LU

FUNCTION

  LU decomposition

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=matrix

OUTPUT

  pivot=gauss pivot
  determinant=determinant of the this

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

586 SUBROUTINE MatrixHyb_LU(this,pivot,determinant)
587 
588 !Arguments ------------------------------------
589 
590 !This section has been created automatically by the script Abilint (TD).
591 !Do not modify the following lines by hand.
592 #undef ABI_FUNC
593 #define ABI_FUNC 'MatrixHyb_LU'
594 !End of the abilint section
595 
596   TYPE(MatrixHyb), INTENT(INOUT) :: this
597   INTEGER, DIMENSION(:), ALLOCATABLE, OPTIONAL, INTENT(INOUT) :: pivot
598   DOUBLE PRECISION, OPTIONAL, INTENT(OUT) :: determinant
599 !Local variables ------------------------------
600   INTEGER :: ligne
601   INTEGER :: colonne 
602   INTEGER :: colonne_de_1_ligne
603   INTEGER :: max_ind_lig
604   INTEGER :: ligne_virtuelle
605   INTEGER :: tail
606   INTEGER, DIMENSION(:), ALLOCATABLE :: pivot_tmp
607   DOUBLE PRECISION :: max_col
608   DOUBLE PRECISION :: inverse_pivot
609   DOUBLE PRECISION :: coef
610   DOUBLE PRECISION :: det
611   DOUBLE PRECISION, DIMENSION(:,:), ALLOCATABLE :: mat_tmp
612 
613   tail = this%tail
614   det = 1.d0
615   FREEIF(pivot_tmp)
616   MALLOC(pivot_tmp,(1:tail))
617 
618   DO ligne=1, tail
619     pivot_tmp(ligne)=ligne
620   END DO
621 
622   MALLOC(mat_tmp,(1:tail,1:tail))
623   mat_tmp(1:tail,1:tail) = this%mat(1:tail,1:tail)
624   DO colonne = 1, tail-1 
625     max_col = ABS(mat_tmp(pivot_tmp(colonne),colonne))
626     max_ind_lig = colonne
627 
628     DO ligne = colonne+1, tail 
629       ligne_virtuelle = pivot_tmp(ligne)
630       IF ( ABS( mat_tmp(ligne_virtuelle,colonne)).GT.max_col) then
631         max_col = ABS(mat_tmp(ligne_virtuelle,colonne))
632         max_ind_lig = ligne
633       ENDIF
634     END DO
635 
636     ligne              = pivot_tmp(colonne) 
637     pivot_tmp(colonne)     = pivot_tmp(max_ind_lig)
638     pivot_tmp(max_ind_lig) = ligne
639     IF ( pivot_tmp(colonne) .NE. pivot_tmp(max_ind_lig) ) det = det * (-1.d0)
640     inverse_pivot=1.d0 / mat_tmp(pivot_tmp(colonne),colonne)
641     det = det * mat_tmp(pivot_tmp(colonne),colonne)
642     DO ligne = colonne+1, tail 
643       ligne_virtuelle = pivot_tmp(ligne)
644       coef = mat_tmp(ligne_virtuelle,colonne)*inverse_pivot
645       mat_tmp(ligne_virtuelle,colonne) = coef
646       DO colonne_de_1_ligne = colonne+1, tail 
647         mat_tmp(ligne_virtuelle,colonne_de_1_ligne)= mat_tmp(ligne_virtuelle,colonne_de_1_ligne)&
648                                          -coef * mat_tmp(pivot_tmp(colonne) ,colonne_de_1_ligne)
649       END DO
650     END DO
651   END DO
652   det = det * mat_tmp(pivot_tmp(tail),tail)
653   IF ( PRESENT(determinant) ) &
654     determinant = det
655   IF ( PRESENT(pivot) ) THEN
656     this%mat(1:tail,1:tail) = mat_tmp(1:tail,1:tail)
657     IF ( ALLOCATED(pivot) .AND. SIZE(pivot) .NE. tail ) THEN
658       FREE(pivot)
659       MALLOC(pivot,(1:tail))
660     ELSE IF ( .NOT. ALLOCATED(pivot) ) THEN
661       MALLOC(pivot,(1:tail))
662     END IF
663     pivot(1:tail)=pivot_tmp(1:tail)
664   END IF
665   FREE(mat_tmp)
666   FREE(pivot_tmp)
667 END SUBROUTINE MatrixHyb_LU

ABINIT/m_MatrixHyb/MatrixHyb_print [ Functions ]

[ Top ] [ Functions ]

NAME

  MatrixHyb_print

FUNCTION

  print Matrix

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=Matrix
  ostream=file stream
  opt_print=0 mat
            1 mat_tau
            2 both

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

756 SUBROUTINE MatrixHyb_print(this,ostream,opt_print)
757 
758 !Arguments ------------------------------------
759 
760 !This section has been created automatically by the script Abilint (TD).
761 !Do not modify the following lines by hand.
762 #undef ABI_FUNC
763 #define ABI_FUNC 'MatrixHyb_print'
764 !End of the abilint section
765 
766   TYPE(MatrixHyb), INTENT(IN) :: this
767   INTEGER, OPTIONAL, INTENT(IN) :: ostream
768   INTEGER, OPTIONAL, INTENT(IN) :: opt_print
769 !Local variables ------------------------------
770   INTEGER                       :: ostream_val
771   INTEGER                       :: opt_val
772   INTEGER                       :: it1
773   INTEGER                       :: it2
774   CHARACTER(LEN=4 )             :: size
775   CHARACTER(LEN=22)             :: string
776 
777   ostream_val = 6
778   opt_val=0
779   IF ( PRESENT(ostream) ) ostream_val = ostream
780   IF ( PRESENT(opt_print) ) opt_val = opt_print
781   WRITE(size,'(I4)') this%tail
782   IF ( MOD(opt_val,2) .EQ. 0 ) THEN
783     string ='(1x,A,1x,'//TRIM(ADJUSTL(size))//'(ES10.2,1x),A)'
784     WRITE(ostream_val,'(A)') "["
785     DO it1 = 1, this%tail
786       WRITE(ostream_val,string) "[",(/ (this%mat(it1,it2),it2=1,this%tail)  /)," ]"
787     END DO
788     WRITE(ostream_val,'(A)') "]"
789   END IF  
790   IF ( opt_val .GE.1 ) THEN
791     string ='(1x,A,1x,'//TRIM(ADJUSTL(size))//'(I4,1x),A)'
792     WRITE(ostream_val,'(A)') "["
793     DO it1 = 1, this%tail
794       WRITE(ostream_val,string) "[",(/ (this%mat_tau(it1,it2),it2=1,this%tail)  /),"]"
795     END DO
796     WRITE(ostream_val,'(A)') "]"
797   END IF
798 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-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=matrix
  new_tail=new size

OUTPUT

SIDE EFFECTS

NOTES

PARENTS

  Will be filled automatically by the parent script

CHILDREN

  Will be filled automatically by the parent script

SOURCE

203 SUBROUTINE MatrixHyb_setSize(this,new_tail)
204 
205 !Arguments ------------------------------------
206 
207 !This section has been created automatically by the script Abilint (TD).
208 !Do not modify the following lines by hand.
209 #undef ABI_FUNC
210 #define ABI_FUNC 'MatrixHyb_setSize'
211 !End of the abilint section
212 
213   TYPE(MatrixHyb), INTENT(INOUT) :: this
214   INTEGER        , INTENT(IN   ) :: new_tail
215 !Local variables ------------------------------
216   INTEGER                        :: size
217 
218   IF ( .NOT. ALLOCATED(this%mat) ) &
219     CALL MatrixHyb_init(this,this%iTech,Wmax=this%Wmax)
220   size = this%size
221   IF( new_tail .GT. size ) THEN
222     CALL MatrixHyb_enlarge(this, MAX(Global_SIZE,new_tail-size))
223   END IF
224   this%tail = new_tail
225 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-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

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