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 | !     Limits for dichotomy (dkinf,dksup) !
  24 | !     -----------------------------------!
  25 | 
  26 |       n_itmax   = 100
  27 |       dk        = 0. 
  28 |       ph_bound  = minval(ph)
  29 |       cb_bound  = maxval(cb)
  30 |       dksup     = cb_bound
  31 |       dkinf     = 0.
  32 |       
  33 |       CALL cdss(ph_bound,cb_bound,dksup,gdek)
  34 | 
  35 |       DO WHILE (gdek.lt.g)
  36 | 
  37 |         dkinf=dksup
  38 |         dksup=dksup*10
  39 |         CALL cdss(ph_bound,cb_bound,dksup,gdek)
  40 |       ENDDO
  41 | 
  42 | !     -----------------------!
  43 | !     Dichotomy to obtain dk !
  44 | !     -----------------------!
  45 | 
  46 |       eps=gmax/1000000
  47 |       err = eps + 1.
  48 |       n_it=0
  49 | 
  50 | !     -------------------!
  51 | !     TEST fonction gdek !
  52 | !     -------------------!
  53 | !     dk=0.00001
  54 | !     open(unit=696,file="gdek_pondere.dat")
  55 | !	DO WHILE (dk.lt.1000) 
  56 | !         gdek = 0.
  57 | !       DO i=1, ngbands 
  58 | !         IF (ph(i).gt.0.) THEN
  59 | !           CALL cdss(ph(i),cb(i),dk,gdek_i)
  60 | !	  ELSE
  61 | !           gdek_i = 1.
  62 | !         ENDIF
  63 | !         gdek = gdek + F(i)*gdek_i
  64 | !      ENDDO 
  65 | !      write(696,*)  gdek,dk
  66 | !      dk=dk*5
  67 | !      ENDDO
  68 | !      close(696)
  69 | 
  70 |       DO WHILE (err.gt.eps.and.n_it.lt.n_itmax)
  71 |         
  72 |         dk=(dksup+dkinf)/2.
  73 | 
  74 | !       ---------------!
  75 | !       Band groupment !
  76 | !       ---------------!
  77 |         gdek = 0.
  78 |         DO i=1, ngbands
  79 | 
  80 |           IF (ph(i).gt.0.) THEN
  81 |             CALL cdss(ph(i),cb(i),dk,gdek_i)
  82 |           ELSE
  83 |             gdek_i = 1.
  84 |           ENDIF
  85 |  
  86 |           gdek = gdek + F(i)*gdek_i
  87 | 
  88 |         ENDDO
  89 | 
  90 |         IF (gdek.lt.g) THEN
  91 |          dkinf=dk
  92 |         ELSE
  93 |          dksup=dk
  94 |         ENDIF
  95 | 
  96 |         err=abs(gdek-g)
  97 |         n_it=n_it+1
  98 | 
  99 |         IF(dksup.lt.1e-4) THEN
 100 |           err = 0.d0
 101 |           dk  = 0.d0
 102 |         ENDIF
 103 | 
 104 |       ENDDO
 105 | 
 106 |       IF (n_it.eq.n_itmax) THEN
 107 |        WRITE(*,*) 'Dichotomy assumed to be converged'
 108 |        WRITE(*,*) 'nit =',n_it
 109 |       ENDIF
 110 | 
 111 |       END SUBROUTINE COFG