26#include <DGtal/base/Common.h>
27#include <DGtal/helpers/StdDefs.h>
28#include <DGtal/helpers/Shortcuts.h>
29#include <DGtal/helpers/ShortcutsGeometry.h>
30#include <DGtal/shapes/SurfaceMesh.h>
31#include <DGtal/geometry/surfaces/DigitalSurfaceRegularization.h>
32#include <DGtal/dec/PolygonalCalculus.h>
33#include <DGtal/math/linalg/DirichletConditions.h>
35#include <polyscope/polyscope.h>
36#include <polyscope/surface_mesh.h>
37#include <polyscope/point_cloud.h>
41#include <Eigen/Sparse>
51typedef std::size_t
Index;
53polyscope::SurfaceMesh *psMesh;
62 PolyDEC calculus(surfmesh);
64 PolyDEC::Form g = calculus.
form0();
65 DC::IntegerVector b = DC::IntegerVector::Zero( g.rows() );
69 std::cout<<
"Number of boundary edges= "<<boundaryEdges.size()<<std::endl;
73 for(
auto &e: boundaryEdges)
76 g(adjVertices.first) = pihVertex(adjVertices.first);
77 g(adjVertices.second) = pihVertex(adjVertices.second);
78 b(adjVertices.first) = 1;
79 b(adjVertices.second) = 1;
83 PolyDEC::Solver solver;
84 PolyDEC::SparseMatrix L_dirichlet = DC::dirichletOperator( L, b );
85 solver.compute( L_dirichlet );
86 ASSERT(solver.info()==Eigen::Success);
87 PolyDEC::Form g_dirichlet = DC::dirichletVector( L, g, b, g );
88 PolyDEC::Form x_dirichlet = solver.solve( g_dirichlet );
89 PolyDEC::Form u = DC::dirichletSolution( x_dirichlet, b, g );
92 psMesh->addVertexScalarQuantity(
"g", g);
93 psMesh->addVertexScalarQuantity(
"u", u);
98 ImGui::SliderFloat(
"Phi scale", &scale, 0., 10.);
99 if(ImGui::Button(
"Compute Laplace problem"))
108 params(
"polynomial",
"0.1*y*y -0.1*x*x - 2.0*z" )(
"gridstep", h );
117 std::vector<std::vector<SH3::SurfaceMesh::Vertex>> faces;
118 std::vector<RealPoint> positions = primalSurface->positions();
119 for(
auto face= 0 ; face < primalSurface->nbFaces(); ++face)
120 faces.push_back(primalSurface->incidentVertices( face ));
122 surfmesh =
SurfMesh(positions.begin(),
126 psMesh = polyscope::registerSurfaceMesh(
"digital surface", positions, faces);
132 polyscope::state::userCallback = myCallback;
Aim: A helper class to solve a system with Dirichlet boundary conditions.
Implements differential operators on polygonal surfaces from .
SparseMatrix globalLaplaceBeltrami(const double lambda=1.0) const
Aim: This class is used to simplify shape and surface creation. With it, you can create new shapes an...
static Parameters parametersGeometryEstimation()
static Parameters defaultParameters()
Aim: This class is used to simplify shape and surface creation. With it, you can create new shapes an...
static KSpace getKSpace(const Point &low, const Point &up, Parameters params=parametersKSpace())
static CountedPtr< DigitizedImplicitShape3D > makeDigitizedImplicitShape3D(CountedPtr< ImplicitShape3D > shape, Parameters params=parametersDigitizedImplicitShape3D())
static CountedPtr< DigitalSurface > makeDigitalSurface(CountedPtr< TPointPredicate > bimage, const KSpace &K, const Parameters ¶ms=parametersDigitalSurface())
static Parameters defaultParameters()
static CountedPtr< SurfaceMesh > makePrimalSurfaceMesh(Cell2Index &c2i, CountedPtr< ::DGtal::DigitalSurface< TContainer > > aSurface)
static CountedPtr< BinaryImage > makeBinaryImage(Domain shapeDomain)
static CountedPtr< ImplicitShape3D > makeImplicitShape3D(const Parameters ¶ms=parametersImplicitShape3D())
DGtal is the top-level namespace which contains all DGtal functions and types.
Aim: Represents an embedded mesh as faces and a list of vertices. Vertices may be shared among faces ...
const VertexPair & edgeVertices(Edge e) const
Edges computeManifoldBoundaryEdges() const
RealPoint & position(Vertex v)
int main(int argc, char **argv)