code: def compute_permutation_for_rotation(positions_a,  # scaled positions
                                       positions_b,
                                       lattice,  # column vectors
                                       symprec):
      
  
      # Sort both sides by some measure which is likely to produce a small
      # maximum value of (sorted_rotated_index - sorted_original_index).
      # The C code is optimized for this case, reducing an O(n^2)
      # search down to ~O(n). (for O(n log n) work overall, including the sort)
      #
      # We choose distance from the nearest bravais lattice point as our measure.
      def sort_by_lattice_distance(fracs):
          carts = np.dot(fracs - np.rint(fracs), lattice.T)
          perm = np.argsort(np.sum(carts**2, axis=1))
          sorted_fracs = np.array(fracs[perm], dtype='double', order='C')
          return perm, sorted_fracs
  
      (perm_a, sorted_a) = sort_by_lattice_distance(positions_a)
      (perm_b, sorted_b) = sort_by_lattice_distance(positions_b)
  
      # Call the C code on our conditioned inputs.
      perm_between = _compute_permutation_c(sorted_a,
                                            sorted_b,
                                            lattice,
                                            symprec)
  
      # Compose all of the permutations for the full permutation.
      #
      # Note the following properties of permutation arrays:
      #
      # 1. Inverse:         if  x[perm] == y  then  x == y[argsort(perm)]
      # 2. Associativity:   x[p][q] == x[p[q]]
      return perm_a[perm_between][np.argsort(perm_b)]
docs: Get the overall permutation such that
  
          positions_a[perm[i]] == positions_b[i]   (modulo the lattice)
  
      or in numpy speak,
  
          positions_a[perm] == positions_b   (modulo the lattice)
  
      This version is optimized for the case where positions_a and positions_b
      are related by a rotation.