code: def _psi_solve(y, eps=0.00001):
      """Solve psi(c)-log(c)=y by dichotomy"""
      if y > 0:
          print("y", y)
          raise ValueError("y>0, the problem cannot be solved")
      u = 1.0
      if y > sp.psi(u) - np.log(u):
          while sp.psi(u) - np.log(u) < y:
              u *= 2
          u /= 2
      else:
          while sp.psi(u) - np.log(u) > y:
              u /= 2
      return _dichopsi_log(u, 2 * u, y, eps)
  
imports: import numpy
  import scipy.special
  import scipy.stats
label_desc: Use image node_id 2 for calling a global function with example usage: _psi_solve(y, eps) and returns: _dichopsi_log