nmwsharp / geometry-central

Applied 3D geometry in C++, with a focus on surface meshes.
https://geometry-central.net
MIT License
1.1k stars 150 forks source link

How to get discrete connection Laplacian mentioned in The Vector Heat Method ? #100

Closed Lancelot899 closed 2 years ago

Lancelot899 commented 2 years ago

How to get discrete connection Laplacian mentioned in The Vector Heat Method ? Can we get a document ?

MarkGillespie commented 2 years ago

The connection Laplacian is called vertexConnectionLaplacian and is described briefly in the docs here.

Below I've included some sample code which uses the connection Laplacian to compute a smooth vector field.

Finally, I should point out that many of the uses of the connection Laplacian mentioned in the vector heat method paper are already included in geometry-central!

I hope this helps!

#include "geometrycentral/surface/manifold_surface_mesh.h"
#include "geometrycentral/surface/meshio.h"
#include "geometrycentral/surface/vertex_position_geometry.h"

#include "geometrycentral/numerical/linear_solvers.h"

#include "polyscope/surface_mesh.h"

using namespace geometrycentral;
using namespace geometrycentral::surface;

int main() {

    // ==============================================================
    //                   Load a mesh
    // ==============================================================

    std::unique_ptr<ManifoldSurfaceMesh> mesh;
    std::unique_ptr<VertexPositionGeometry> geometry;
    std::tie(mesh, geometry) = readManifoldSurfaceMesh("filename.obj");

    // ==============================================================
    // Use the connection Laplacian to compute a smooth vector field
    // ==============================================================

    // Get connection laplacian
    geometry->requireVertexConnectionLaplacian();
    SparseMatrix<std::complex<double>> connectionLaplacian =
        geometry->vertexConnectionLaplacian;

    // Get vector mass matrix
    geometry->requireVertexGalerkinMassMatrix();
    SparseMatrix<std::complex<double>> vectorMassMatrix =
        geometry->vertexGalerkinMassMatrix.cast<std::complex<double>>();

    // Find the smallest eigenvector
    Vector<std::complex<double>> connectionLaplacianEigenvector =
        smallestEigenvectorSquare(connectionLaplacian, vectorMassMatrix);

    // ==============================================================
    //      Display mesh and smooth vector field with polyscope
    // ==============================================================

    // Initialize polyscope
    polyscope::init();

    // Register mesh with polyscope
    polyscope::SurfaceMesh* psMesh = polyscope::registerSurfaceMesh(
        "mesh", geometry->vertexPositions, mesh->getFaceVertexList());

    // Set vertex tangent spaces
    geometry->requireVertexTangentBasis();
    VertexData<Vector3> vBasisX(*mesh);
    for (Vertex v : mesh->vertices()) {
        vBasisX[v] = geometry->vertexTangentBasis[v][0];
    }
    psMesh->setVertexTangentBasisX(vBasisX);

    // Display vector field on surface
    psMesh->addVertexIntrinsicVectorQuantity("vector field",
                                             connectionLaplacianEigenvector);

    // Give control to the polyscope gui;
    polyscope::show();

    return 0;
}
Lancelot899 commented 2 years ago

The connection Laplacian is called vertexConnectionLaplacian and is described briefly in the docs here.

Below I've included some sample code which uses the connection Laplacian to compute a smooth vector field.

Finally, I should point out that many of the uses of the connection Laplacian mentioned in the vector heat method paper are already included in geometry-central!

* [Smooth vector fields/direction fields](http://geometry-central.net/surface/algorithms/direction_fields/)

* [The vector heat method itself](http://geometry-central.net/surface/algorithms/vector_heat_method/) (includes vector interpolation and computing logarithmic maps)

* [Surface centers](http://geometry-central.net/surface/algorithms/surface_centers/)

* [Voronoi tessellations](http://geometry-central.net/surface/algorithms/geodesic_voronoi_tessellations/)

I hope this helps!

#include "geometrycentral/surface/manifold_surface_mesh.h"
#include "geometrycentral/surface/meshio.h"
#include "geometrycentral/surface/vertex_position_geometry.h"

#include "geometrycentral/numerical/linear_solvers.h"

#include "polyscope/surface_mesh.h"

using namespace geometrycentral;
using namespace geometrycentral::surface;

int main() {

    // ==============================================================
    //                   Load a mesh
    // ==============================================================

    std::unique_ptr<ManifoldSurfaceMesh> mesh;
    std::unique_ptr<VertexPositionGeometry> geometry;
    std::tie(mesh, geometry) = readManifoldSurfaceMesh("filename.obj");

    // ==============================================================
    // Use the connection Laplacian to compute a smooth vector field
    // ==============================================================

    // Get connection laplacian
    geometry->requireVertexConnectionLaplacian();
    SparseMatrix<std::complex<double>> connectionLaplacian =
        geometry->vertexConnectionLaplacian;

    // Get vector mass matrix
    geometry->requireVertexGalerkinMassMatrix();
    SparseMatrix<std::complex<double>> vectorMassMatrix =
        geometry->vertexGalerkinMassMatrix.cast<std::complex<double>>();

    // Find the smallest eigenvector
    Vector<std::complex<double>> connectionLaplacianEigenvector =
        smallestEigenvectorSquare(connectionLaplacian, vectorMassMatrix);

    // ==============================================================
    //      Display mesh and smooth vector field with polyscope
    // ==============================================================

    // Initialize polyscope
    polyscope::init();

    // Register mesh with polyscope
    polyscope::SurfaceMesh* psMesh = polyscope::registerSurfaceMesh(
        "mesh", geometry->vertexPositions, mesh->getFaceVertexList());

    // Set vertex tangent spaces
    geometry->requireVertexTangentBasis();
    VertexData<Vector3> vBasisX(*mesh);
    for (Vertex v : mesh->vertices()) {
        vBasisX[v] = geometry->vertexTangentBasis[v][0];
    }
    psMesh->setVertexTangentBasisX(vBasisX);

    // Display vector field on surface
    psMesh->addVertexIntrinsicVectorQuantity("vector field",
                                             connectionLaplacianEigenvector);

    // Give control to the polyscope gui;
    polyscope::show();

    return 0;
}

Thinks a lot, the code and doc are very useful, and may I get some derivation nodes of vertexConnectionLaplacian or some literature about it?

keenancrane commented 2 years ago

The vector heat method paper is probably the best reference for the discussion/derivation of that version of the discrete connection Laplacian. It was also used in this paper, albeit using a connection derived from an arbitrary vector field rather than the surface geometry: https://www.cs.cmu.edu/~kmcrane/Projects/StripePatterns/paper.pdf

Lancelot899 commented 2 years ago

The vector heat method paper is probably the best reference for the discussion/derivation of that version of the discrete connection Laplacian. It was also used in this paper, albeit using a connection derived from an arbitrary vector field rather than the surface geometry: https://www.cs.cmu.edu/~kmcrane/Projects/StripePatterns/paper.pdf

Oh, thinks, I seem to get the idea now. We need to calculate the energy of two vector, but this two vectors are not in the save vector space, so the need do Parallel Transport to determine the difference. Then the the vector Dirichlet energy need to formulate as \sum w_{ij}|vj - r{ij}v_i|^2, and finally we get the vertexConnectionLaplacian. May I get the right point?

keenancrane commented 2 years ago

Yep! That’s it. :-)