code: def _compute_C_matrix(self):
          
          self.C = np.repeat(np.diag(self.Lp)[:, np.newaxis],
                             self.Lp.shape[0], axis=1)
          self.C += np.repeat(np.diag(self.Lp)[np.newaxis, :],
                              self.Lp.shape[0], axis=0)
          self.C -= 2*self.Lp
          # the following is much slower
          # self.C = np.zeros(self.Lp.shape)
          # for i in range(self.Lp.shape[0]):
          #     for j in range(self.Lp.shape[1]):
          #         self.C[i, j] = self.Lp[i, i] + self.Lp[j, j] - 2*self.Lp[i, j]
          volG = np.sum(self.z)
          self.C *= volG
          settings.mt(0, 'computed commute distance matrix')
          self.distances_dpt = self.C
docs: See Fouss et al. (2006) and von Luxburg et al. (2007).
  
          This is the commute-time matrix. It's a squared-euclidian distance
          matrix in :math:`\\mathbb{R}^n`.