anchor.diff:      type(parameters)                 :: para
   
       c2 = c**2
  -    
  +
       call get_parameters(para)
       rs = 2._x_precision*G*para%M/c2
  -    
  +
       state_0%omega_0 = sqrt( G*para%M/(rs)**3 )
       state_0%temps_0 = 2._x_precision / state_0%omega_0
       state_0%x_0     = sqrt(rs)
  
positive.code:       if (i .eq. 1) then
           state_out%v(i)  = 0._x_precision
        else if (i .eq. n_cell) then
           state_out%v(i)  = - 1._x_precision / state_out%S(i) / state_out%x(i) 
        else
           state_out%v(i)  = - 1._x_precision / state_out%S(i) / state_out%x(i) * ( ((state_out%nu(i+1) * &
                             state_out%S(i+1)) - (state_out%nu(i-1) * state_out%S(i-1))) /  &
  
negative.code:   private
  
    public :: compute_variables, dim_adim, init_variable_0, state_0
  
  contains
    !--------------------------------------------------------------
    subroutine compute_variables(state_out)
  