code: def _kolmogn_DMTW(n, d, cdf=True):
      r"""Computes the Kolmogorov CDF:  Pr(D_n <= d) using the MTW approach to
      the Durbin matrix algorithm.
  
      Durbin (1968); Marsaglia, Tsang, Wang (2003). [1], [3].
      """
      # Write d = (k-h)/n, where k is positive integer and 0 <= h < 1
      # Generate initial matrix H of size m*m where m=(2k-1)
      # Compute k-th row of (n!/n^n) * H^n, scaling intermediate results.
      # Requires memory O(m^2) and computation O(m^2 log(n)).
      # Most suitable for small m.
  
      if d >= 1.0:
          return _select_and_clip_prob(1.0, 0.0, cdf)
      nd = n * d
      if nd <= 0.5:
          return _select_and_clip_prob(0.0, 1.0, cdf)
      k = int(np.ceil(nd))
      h = k - nd
      m = 2 * k - 1
  
      H = np.zeros([m, m])
  
      # Initialize: v is first column (and last row) of H
      #  v[j] = (1-h^(j+1)/(j+1)!  (except for v[-1])
      #  w[j] = 1/(j)!
      # q = k-th row of H (actually i!/n^i*H^i)
      intm = np.arange(1, m + 1)
      v = 1.0 - h**intm
      w = np.empty(m)
      fac = 1.0
      for j in intm:
          w[j - 1] = fac
          fac /= j  # This might underflow.  Isn't a problem.
          v[j - 1] *= fac
      tt = max(2 * h - 1.0, 0) ** m - 2 * h**m
      v[-1] = (1.0 + tt) * fac
  
      for i in range(1, m):
          H[i - 1 :, i] = w[: m - i + 1]
      H[:, 0] = v
      H[-1, :] = np.flip(v, axis=0)
  
      Hpwr = np.eye(np.shape(H)[0])  # Holds intermediate powers of H
      nn = n
      expnt = 0  # Scaling of Hpwr
      Hexpnt = 0  # Scaling of H
      while nn > 0:
          if nn % 2:
              Hpwr = np.matmul(Hpwr, H)
              expnt += Hexpnt
          H = np.matmul(H, H)
          Hexpnt *= 2
          # Scale as needed.
          if np.abs(H[k - 1, k - 1]) > _EP128:
              H /= _EP128
              Hexpnt += _E128
          nn = nn // 2
  
      p = Hpwr[k - 1, k - 1]
  
      # Multiply by n!/n^n
      for i in range(1, n + 1):
          p = i * p / n
          if np.abs(p) < _EM128:
              p *= _EP128
              expnt -= _E128
  
      # unscale
      if expnt != 0:
          p = np.ldexp(p, expnt)
  
      return _select_and_clip_prob(p, 1.0 - p, cdf)
  
file_path: scipy/scipy/stats/_ksstats.py
imports: import numpy
  import scipy.special
  import scipy.special._ufuncs
  from scipy._lib._finite_differences import _derivative
label_desc: Use image node_id 4 for calling a global function with example usage: _kolmogn_DMTW(n, d, cdf) and returns: _select_and_clip_prob, _select_and_clip_prob, _select_and_clip_prob