1 | include(dom.inc)
   2 | 
   3 |       SUBROUTINE COFG(F,ngbands,ph,cb,g,gmax,dk)
   4 | 
   5 |       IMPLICIT NONE
   6 | 
   7 | !     IN
   8 |       DOM_INT  :: ngbands
   9 |       DOM_REAL :: g, gmax
  10 |       DOM_REAL :: ph(ngbands), cb(ngbands)
  11 |       DOM_REAL :: F(ngbands)
  12 | 
  13 | !     OUT
  14 |       DOM_REAL :: dk
  15 | 
  16 | !     LOCAL
  17 |       DOM_REAL :: dksup, dkinf
  18 |       DOM_REAL :: gdek, gdek_i, err, eps
  19 |       DOM_REAL :: ph_bound, cb_bound
  20 |       DOM_INT  :: i, n_it, n_itmax
  21 | 
  22 | 
  23 | !     -----------------------------------!
  24 | !     Limits for dichotomy (dkinf,dksup) !
  25 | !     -----------------------------------!
  26 | 
  27 |       n_itmax   = 100
  28 |       dk        = 0. 
  29 |       ph_bound  = minval(ph)
  30 |       cb_bound  = maxval(cb)
  31 |       dksup     = cb_bound
  32 |       dkinf     = 0.
  33 |       
  34 |       CALL cdss(ph_bound,cb_bound,dksup,gdek)
  35 | !     print*, '       + Boundaries:'
  36 | !     print*, '         ph_bound:', ph_bound
  37 | !     print*, '         cb_bound:', cb_bound
  38 | !     print*, '         On a gdek = ', gdek
  39 | !     print*, '         On cherche g = ', g
  40 | 
  41 |       DO WHILE (gdek.lt.g)
  42 |         dkinf=dksup
  43 |         dksup=dksup*10
  44 |         CALL cdss(ph_bound,cb_bound,dksup,gdek)
  45 |       ENDDO
  46 | !     print*, "         dkinf, dksup: ", dkinf, dksup
  47 | !     print*, "         gdek(dksup)= ",gdek,"g=",g
  48 | 
  49 | !     -----------------------!
  50 | !     Dichotomy to obtain dk !
  51 | !     -----------------------!
  52 |         
  53 |       eps=gmax/1000
  54 |       err = eps + 1.
  55 |       n_it=0
  56 | 
  57 | 
  58 | !     -------------------!
  59 | !     TEST fonction gdek !
  60 | !     -------------------!
  61 | !      open(unit=696,file="gdek_pondere.dat")
  62 | !      DO dk=0.1,60,0.5
  63 | !        gdek = 0.
  64 | !        DO i=1, ngbands
  65 | !          CALL cdss(ph(i),cb(i),dk,gdek_i)
  66 | !          gdek = gdek + F(i)*gdek_i
  67 | !        ENDDO
  68 | !        write(696,*) dk, gdek
  69 | !      ENDDO
  70 | !      close(696)
  71 | 
  72 | 
  73 |       DO WHILE (err.gt.eps.and.n_it.lt.n_itmax)
  74 | 
  75 |         dk=(dksup+dkinf)/2.
  76 | 
  77 | !       ---------------!
  78 | !       Band groupment !
  79 | !       ---------------!
  80 |         gdek = 0.
  81 |         DO i=1, ngbands
  82 | 
  83 |           IF ( ph(i).gt.0. ) THEN
  84 |             CALL cdss(ph(i),cb(i),dk,gdek_i)
  85 |           ELSE
  86 |             gdek_i = 1.
  87 |           ENDIF
  88 | 
  89 |           gdek = gdek + F(i)*gdek_i
  90 | 
  91 |         ENDDO
  92 | 
  93 |         IF (gdek.lt.g) THEN
  94 |          dkinf=dk
  95 |         ELSE
  96 |          dksup=dk
  97 |         ENDIF
  98 | 
  99 |         err=abs(gdek-g)
 100 |   
 101 |         n_it=n_it+1
 102 | 
 103 | !       print*, " "
 104 | !       print*, "         err,eps : ", err,eps
 105 | !       print*, "         dk : ", dk
 106 | !       print*, "         gdek: ", gdek
 107 | !       print*, "         g   : ", g
 108 | !       print*, "         dksup,dkmin   : ", dksup,dkinf
 109 | 
 110 |         IF(dksup.lt.1e-4) THEN
 111 |           err = 0.d0
 112 |           dk  = 0.d0
 113 |         ENDIF
 114 | 
 115 |       ENDDO
 116 | 
 117 | 
 118 | !     IF (n_it.eq.n_itmax) THEN
 119 | !       PRINT*, 'Dichotomy assumed to be converged after MAX iterations'
 120 | !     ENDIF
 121 | 
 122 |       END SUBROUTINE COFG