anchor.diff:    ! interval [Smin,Smax] with an epsilon precision
     !-------------------------------------------------------------------------
     real(kind=x_precision) function dichotomy(Smin, Smax, eps, T, omega, optical_depth)
  -    use mod_read_parameters
  -    use mod_constants
  -    use mod_variables
       implicit none
   
       real(kind=x_precision),intent(inout)                     :: Smin,Smax
  
positive.code:     real(kind = x_precision)                 ,intent(out):: sigma_c_thick
      real(kind = x_precision)                 ,intent(out):: temp_c_thick
      integer                                              :: i
      i = 1
      do while (sigma_real_thick(i) < sigma_real_thick(i+1) .and. i < nb_it - 1)
         index_fcp     = i
         sigma_c_thick = sigma_real_thick(i)
         i             = i + 1
      enddo
      temp_c_thick  = temp_real_thick(i)
  
      endsubroutine first_critical_point
  
  
  
  
  
      !-------------------------------------------------------------------------
      ! Subroutine in order to find the second critical point
      !-------------------------------------------------------------------------
      subroutine second_critical_point(sigma_real_thick, sigma_real_thin, temp_real_thick,&
           index_fcp, index_scp, sigma_c_thin, temp_c_thin, nb_it)
      implicit none
      integer,intent(in)                                   :: nb_it
  
      real(kind = x_precision),dimension(nb_it),intent(in) :: sigma_real_thick
      real(kind = x_precision),dimension(nb_it),intent(in) :: temp_real_thick
      real(kind = x_precision),dimension(nb_it),intent(in) :: sigma_real_thin
  
      integer                                  ,intent(in) :: index_fcp
      integer                                  ,intent(out):: index_scp
  
  
negative.code:     real(kind=x_precision)                                   :: f_center
      real(kind=x_precision)                                   :: S_center = 0._x_precision
      integer                                                  :: j = 0
  
      !-------------------------------------------------------------------------
      ! N-> Number of iterations for the dichotomy
  