1 | c.......................................................................
   2 | c
   3 | c  inverse of the cumulative for the inverse transmission function
   4 | c  by dichotomy
   5 | c
   6 | c  in :  * $\phi$                ---> phcofg
   7 | c        * $\overline k$         ---> cbcofg
   8 | c        * $g$                   ---> gcofg
   9 | c
  10 | c  out : * $k$                   ---> dk
  11 | c
  12 | c.......................................................................
  13 |       subroutine cofg(phcofg,cbcofg,gcofg,gmax,dk)
  14 | c.......................................................................
  15 |       implicit double precision (a-h,o-z)
  16 | c.......................................................................
  17 |       parameter (icofg1=100)
  18 |       parameter (istop=100000)   !critere d'arret stupide si converge pas
  19 |                                  !amaury
  20 | c.......................................................................
  21 | c     Recherche des bornes pour la dichitomie (dkinf,dksup)
  22 | c.....................................................................
  23 | c            PRINT *,'cbcofg dans cofg en entree',cbcofg
  24 | 	dksup=cbcofg
  25 | 	dkinf=0.D+0
  26 | 10	CALL cdss(phcofg,cbcofg,0.D+0,dksup,gdek)
  27 | 	IF (gdek.lt.gcofg) THEN
  28 | 	    dkinf=dksup
  29 |             dksup=dksup*10
  30 | c           PRINT *,dkinf,dksup
  31 |             GO TO 10
  32 |         ENDIF
  33 | c.....................................................................
  34 | c     Dichotomie pour obtenir dk par evaluation de gdek
  35 | c.....................................................................
  36 |         index = 0   !amaury
  37 | 20	dk=(dksup+dkinf)/2
  38 |         index = index + 1
  39 |         CALL cdss(phcofg,cbcofg,0.D+0,dk,gdek)
  40 | 	IF (gdek.lt.gcofg) THEN
  41 | 	     dkinf=dk
  42 |         ELSE
  43 |              dksup=dk
  44 | 	ENDIF
  45 | c.....................................................................
  46 | 	err=dabs(gdek-gcofg)
  47 | 	eps=gmax/100000
  48 | !!!!!!!!!......!eps=1.d-4 !precision du fit sur g(a) A REVOIR PB AVEC CO,non H2O
  49 |                   ! ET PB NON Cest LONG
  50 |         if (index.gt.istop) then
  51 |             go to 30
  52 |         endif
  53 | !!!!!!!!...........   
  54 | 	IF (err.gt.eps) THEN
  55 | c           PRINT *,cbcofg,dk,eps,err
  56 |             GO TO 20
  57 | 	ENDIF
  58 | 30      RETURN
  59 | c.......................................................................
  60 |       end