code: def get_uniformly_controlled_rotation_matrix(k):
      
      M = np.full((2 ** k, 2 ** k), 2 ** -k)
      for i in range(2 ** k):
          g_i = i ^ (i >> 1)  # Gray code for i
          for j in range(2 ** k):
              M[i, j] *= (-1) ** (bin(j & g_i).count("1"))
      return M
docs: Returns the matrix represented by :math:`M_{ij}` in arXiv:quant-ph/0407010.
  
      This matrix converts the angles of :math:`k`-fold uniformly
      controlled rotations to the angles of the efficient gate decomposition.
  
      :param int k: number of control qubits
      :return: the matrix :math:`M_{ij}`
      :rtype: 2darray