code: def _ndtri(p):
      """Implements ndtri core logic."""
  
      # Constants used in piece-wise rational approximations. Taken from the cephes
      # library:
      # https://root.cern.ch/doc/v608/SpecFuncCephesInv_8cxx_source.html
  
      p0 = [
          -1.23916583867381258016e0,
          1.39312609387279679503e1,
          -5.66762857469070293439e1,
          9.80010754185999661536e1,
          -5.99633501014107895267e1,
      ]
      q0 = [
          -1.18331621121330003142e0,
          1.59056225126211695515e1,
          -8.20372256168333339912e1,
          2.00260212380060660359e2,
          -2.25462687854119370527e2,
          8.63602421390890590575e1,
          4.67627912898881538453e0,
          1.95448858338141759834e0,
          1.0,
      ]
      p1 = [
          -8.57456785154685413611e-4,
          -3.50424626827848203418e-2,
          -1.40256079171354495875e-1,
          2.18663306850790267539e0,
          1.46849561928858024014e1,
          4.40805073893200834700e1,
          5.71628192246421288162e1,
          3.15251094599893866154e1,
          4.05544892305962419923e0,
      ]
      q1 = [
          -9.33259480895457427372e-4,
          -3.80806407691578277194e-2,
          -1.42182922854787788574e-1,
          2.50464946208309415979e0,
          1.50425385692907503408e1,
          4.13172038254672030440e1,
          4.53907635128879210584e1,
          1.57799883256466749731e1,
          1.0,
      ]
      p2 = [
          6.23974539184983293730e-9,
          2.65806974686737550832e-6,
          3.01581553508235416007e-4,
          1.23716634817820021358e-2,
          2.01485389549179081538e-1,
          1.33303460815807542389e0,
          3.93881025292474443415e0,
          6.91522889068984211695e0,
          3.23774891776946035970e0,
      ]
      q2 = [
          6.79019408009981274425e-9,
          2.89247864745380683936e-6,
          3.28014464682127739104e-4,
          1.34204006088543189037e-2,
          2.16236993594496635890e-1,
          1.37702099489081330271e0,
          3.67983563856160859403e0,
          6.02427039364742014255e0,
          1.0,
      ]
  
      def _create_polynomial(var, coeffs):
          """Compute n_th order polynomial via Horner's method."""
          coeffs = np.array(coeffs, var.dtype.as_numpy_dtype)
          if not coeffs.size:
              return array_ops.zeros_like(var)
          return coeffs[0] + _create_polynomial(var, coeffs[1:]) * var
  
      maybe_complement_p = array_ops.where_v2(
          p > -np.expm1(-2.0), 1.0 - p, p
      )
      # Write in an arbitrary value in place of 0 for p since 0 will cause NaNs
      # later on. The result from the computation when p == 0 is not used so any
      # number that doesn't result in NaNs is fine.
      sanitized_mcp = array_ops.where_v2(
          maybe_complement_p <= 0.0,
          array_ops.fill(
              array_ops.shape(p), np.array(0.5, p.dtype.as_numpy_dtype)
          ),
          maybe_complement_p,
      )
  
      # Compute x for p > exp(-2): x/sqrt(2pi) = w + w**3 P0(w**2)/Q0(w**2).
      w = sanitized_mcp - 0.5
      ww = w**2
      x_for_big_p = w + w * ww * (
          _create_polynomial(ww, p0) / _create_polynomial(ww, q0)
      )
      x_for_big_p *= -np.sqrt(2.0 * np.pi)
  
      # Compute x for p <= exp(-2): x = z - log(z)/z - (1/z) P(1/z) / Q(1/z),
      # where z = sqrt(-2. * log(p)), and P/Q are chosen between two different
      # arrays based on whether p < exp(-32).
      z = math_ops.sqrt(-2.0 * math_ops.log(sanitized_mcp))
      first_term = z - math_ops.log(z) / z
      second_term_small_p = (
          _create_polynomial(1.0 / z, p2)
          / _create_polynomial(1.0 / z, q2)
          / z
      )
      second_term_otherwise = (
          _create_polynomial(1.0 / z, p1)
          / _create_polynomial(1.0 / z, q1)
          / z
      )
      x_for_small_p = first_term - second_term_small_p
      x_otherwise = first_term - second_term_otherwise
  
      x = array_ops.where_v2(
          sanitized_mcp > np.exp(-2.0),
          x_for_big_p,
          array_ops.where_v2(z >= 8.0, x_for_small_p, x_otherwise),
      )
  
      x = array_ops.where_v2(p > 1.0 - np.exp(-2.0), x, -x)
      infinity_scalar = constant_op.constant(np.inf, dtype=p.dtype)
      infinity = array_ops.fill(array_ops.shape(p), infinity_scalar)
      x_nan_replaced = array_ops.where_v2(
          p <= 0.0, -infinity, array_ops.where_v2(p >= 1.0, infinity, x)
      )
      return x_nan_replaced
  
file_path: tensorflow/tensorflow/python/ops/distributions/special_math.py
imports: import numpy
  from tensorflow.python.framework import constant_op
  from tensorflow.python.framework import ops
  from tensorflow.python.ops import array_ops
  from tensorflow.python.ops import math_ops
label_desc: Use image node_id 4 for calling a global function with example usage: _ndtri(p) and returns: x_nan_replaced, unknown, array_ops