code: def fast_cov(x, y=None, destination=None):
      
      validate_inputs(x, y, destination)
  
      if y is None:
          y = x
  
      if destination is None:
          destination = numpy.zeros((x.shape[1], y.shape[1]))
  
      mean_x = numpy.mean(x, axis=0)
      mean_y = numpy.mean(y, axis=0)
  
      mean_centered_x = (x - mean_x).astype(destination.dtype)
      mean_centered_y = (y - mean_y).astype(destination.dtype)
      
      numpy.dot(mean_centered_x.T, mean_centered_y, out=destination)
      numpy.divide(destination, (x.shape[0] - 1), out=destination)
  
      return destination
docs: calculate the covariance matrix for the columns of x (MxN), or optionally, the covariance matrix between the
      columns of x and and the columns of y (MxP).  (In the language of statistics, the columns are variables, the rows
      are observations).
  
      Args:
          x (numpy array-like) MxN in shape
          y (numpy array-like) MxP in shape
          destination (numpy array-like) optional location where to store the results as they are calculated (e.g. a numpy
              memmap of a file)
  
          returns (numpy array-like) array of the covariance values
              for defaults (y=None), shape is NxN
              if y is provided, shape is NxP