anchor.diff:      rmin               = 3._x_precision * rs
       Mdot_0             = param%Mdot
       Omega_0            = sqrt(G * param%M / rs**3 )
  -    Sigma_0            = Mdot_0 /(Omega_0 * rs**2 * 2 * pi)
  +    Sigma_0            = Mdot_0 /(Omega_0 * rs**2 * 2_x_precision * pi)
       T_0                = (Mdot_0 * c**2 / (48._x_precision * pi * rs**2 * stefan * &
                            sqrt(27._x_precision) ) )**(1._x_precision/4._x_precision)
       rho_0              = Sigma_0 / (2._x_precision * rs)
  
positive.code:     !Select the case for the opticaly depth to compute Fz
      !-------------------------------------------------------------------------
  
      optical_depth        = 1
  
      !if (tau_eff .ge. 1.)  then
      !   optical_depth     = 1
      !else
      !   optical_depth     = 0
      !end if
      
       select case(optical_depth)
        case(1)
           Fz = 4._x_precision * c**2 * T**4/(27._x_precision * sqrt(3._x_precision) &
                * (K_ff + K_e) * Sigma * Sigma_0)
        case (0)
           !Fz = 4._x_precision * rs * E_ff * H / (Omega_0 * Sigma_0)
           Fz = 6.22d20 * 2._x_precision / Omega_0 * rho_0 * H * rho**2 * sqrt(T * T_0)
       end select
      Q_minus             = 2._x_precision * Fz /Sigma
      Q_plus              = 9._x_precision /4._x_precision * nu * Omega**2
      f                   = Q_plus - Q_minus
      
    end subroutine variables
    
    
    !-------------------------------------------------------------------------
    ! Function in order to determine the change of sign in a given interval
    ! [Smin,Smax] with an epsilon precision
    !-------------------------------------------------------------------------
    real(kind=x_precision) function dichotomy(Smin, Smax, eps, T, omega, sigma_0, Omega_0,rs, T_0, rho_0)
      use mod_read_parameters
  
negative.code:     write(*,"(' rmin        =',1p,E12.4)") rmin
      write(*,"(' rs          =',1p,E12.4)") rs
      write(*,*)'****************************************'
  
      read(*,*)
    end subroutine display_initial_variables
  
  