Triplets of integers of (start vertex, end vertex, lagtime) representing the path of particles are written to a binary file and suitably sorted. By performing the sorting, all paths starting or ending at a given vertex form a contiguous block, irrespective of the lagtime, ∆. Then, a single call of the VTP algorithm on the source vertex i is sufficient to evaluate one complete block, as all the paths that start or end at vertex i are grouped together, independently of the lagtime. Finally, the calculated geodesic distances are assigned both to the source and end vertices of the path, while keeping track of the lagtime. Once properly processed, the mesh contains the discretized spatial distributions of the geodesic MSD (gMSD) values at various lagtimes. Importantly, our approach is flexible enough to handle surface meshes of arbitrary shapes subject to the condition that they do not change in time. It is also modular, so that the distance calculation algorithm can readily be swapped with other existing methods. An illustration of the major steps of the algorithm is presented in Fig. . Details on the creation of the surface meshes, the magnitude of the discretization error, the handling of periodic boundary conditions and information about the scaling of the algorithm can be found in the Supporting Information. Here, we used meshes of 40,000 points with a spacing of ≈ 0.4 nm to cover four periodic images, yet the algorithm scales reasonably well to meshes at least twice this size, corresponding to membranes with over 10,000 lipids. Independent calls of the VPT algorithm for different source vertices allow for a further increase in studied system size trough the implemented trivial parallelization. For the sake of simplicity, the computed gMSD curves are interpreted by assuming a free diffusion model at each point of the surface, D geo = gMSD/4∆. The free diffusion model neglects the appearance of driftterms due to the non-flat free energy landscape, that is, the non-uniform lipid density.