def approx_geodesic_distance(self, verts, m=0.1):
        npt = len(self.pts)
        t = m * self.avg_edge_length ** 2 # time of heat evolution

        if m not in self._rlfac_solvers:
            B, D, W, V = self.laplace_operator
            nLC = W - V # negative laplace matrix
            spD = sparse.dia_matrix((D,[0]), (npt,npt)).tocsr() # lumped mass matrix
            
            lfac = spD - t * nLC # backward Euler matrix

            # Exclude rows with zero weight (these break the sparse LU, that finicky fuck)
            goodrows = np.nonzero(~np.array(lfac.sum(0) == 0).ravel())[0]
            self._goodrows = goodrows
            self._rlfac_solvers[m] = sparse.linalg.dsolve.factorized(lfac[goodrows][:,goodrows])

        # Solve system to get u, the heat values
        u0 = np.zeros((npt,)) # initial heat values
        u0[verts] = 1.0
        goodu = self._rlfac_solvers[m](u0[self._goodrows])
        u = np.zeros((npt,))
        u[self._goodrows] = goodu

        return -4 * t * np.log(u)
dataset/ETHPy150Open gallantlab/pycortex/cortex/polyutils.py/Surface.approx_geodesic_distance/original