anchor.diff:        write(fid_tot,'(3(A16))') 'Surface_density', 'Temperature', 'Optical_depth'
   
         do j = 1, nb_it
  -        call variables(k, temp_real(j), sigma_real(j), f, optical_depth, tau_eff)
  -
  -        temp_real(j)  =  temp_real(j) * state_0%T_0
  -        sigma_real(j) =  sigma_real(j) * state_0%S_0
  +        call variables(k, temperature(j), sigma(j), f, optical_depth, tau_eff)
  +
  +        temp_real(j)  = temperature(j) * state_0%T_0
  +        sigma_real(j) = sigma(j) * state_0%S_0
   
           write(fid_tot,fmt = '(3(e16.6e2))') sigma_real(j), temp_real(j), tau_eff
         enddo
  
positive.code:   real(kind = x_precision) function dichotomy(k, i, optical_depth)
      implicit none
  
      integer,                  intent(in) :: k ! Position in the disk, to get r/Omega
      integer,                  intent(in) :: i ! Iteration counter, to get the temperature
      integer,                  intent(in) :: optical_depth
      !-------------------------------------------------------------------------
      integer                  :: j
      real(kind = x_precision) :: T ! Temperature for which to solve the dichotomy
  
negative.code:   !-------------------------------------------------------------------------
    ! Subroutine in order to find the first critical point
    !-------------------------------------------------------------------------
    subroutine first_critical_point(sigma_real_thick, index_fcp, sigma_c_thick, temp_c_thick)
      implicit none
  
      real(kind = x_precision), dimension(nb_it), intent(in)  :: sigma_real_thick
  
      integer,                                    intent(out) :: index_fcp
      real(kind = x_precision),                   intent(out) :: sigma_c_thick
  