songsenand
/
cgal
mirror of https://github.com/CGAL/cgal
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

basic cache functionality

This commit is contained in:
Christina Vaz 2018-06-19 15:11:50 -04:00
parent 46f11a0a29
commit 3ff55a3a29
2 changed files with 46 additions and 22 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

28
Heat_method_3/include/CGAL/Heat_method_3/Heat_method_3.h
View File

@ -124,6 +124,11 @@ namespace Heat_method_3 {
* add `vd` to the source set, returning `false` if `vd` is already in the set.
*/
const VertexDistanceMap& get_vertex_distance_map() const
{
return vdm;
}
const Matrix& mass_matrix() const
{
return m_mass_matrix;
@ -148,6 +153,7 @@ namespace Heat_method_3 {
bool add_source(vertex_descriptor vd)
{
source_change_flag = true;
return sources.insert(vd).second;
}
@ -156,6 +162,7 @@ namespace Heat_method_3 {
*/
bool remove_source(vertex_descriptor vd)
{
source_change_flag = true;
return (sources.erase(vd) == 1);
}
@ -172,6 +179,7 @@ namespace Heat_method_3 {
*/
void clear_sources()
{
source_change_flag = true;
sources.clear();
return;
}
@ -420,11 +428,20 @@ namespace Heat_method_3 {
void update()
{
double d=0;
build();
if(source_change_flag)
{
//don't need to recompute Mass matrix, cotan matrix or timestep reflect that in this function
sources = get_sources();
kronecker = kronecker_delta(sources);
solved_u = solve_cotan_laplace(m_mass_matrix, m_cotan_matrix, kronecker, m_time_step, dimension);
X = compute_unit_gradient(solved_u);
index_divergence = compute_divergence(X, dimension);
solved_phi = solve_phi(m_cotan_matrix, index_divergence, dimension);
source_change_flag = false;
}
BOOST_FOREACH(vertex_descriptor vd, vertices(tm)){
Index i_d = get(vertex_id_map, vd);
d = solved_phi(i_d,0);
std::cout<<d<<"\n";
put(vdm,vd,d);
}
}
@ -433,6 +450,7 @@ namespace Heat_method_3 {
void build()
{
source_change_flag = false;
CGAL_precondition(is_triangle_mesh(tm));
vertex_id_map = get(Vertex_property_tag(),const_cast<TriangleMesh&>(tm));
Index i = 0;
@ -445,7 +463,7 @@ namespace Heat_method_3 {
put(face_id_map, fd, face_i++);
}
int m = static_cast<int>(num_vertices(tm));
dimension = m;
//mass matrix entries
std::vector<triplet> A_matrix_entries;
//cotan matrix entries
@ -512,6 +530,7 @@ namespace Heat_method_3 {
index_divergence = compute_divergence(X, m);
solved_phi = solve_phi(m_cotan_matrix, index_divergence, m);
}
int dimension;
const TriangleMesh& tm;
VertexDistanceMap vdm;
VertexPointMap vpm;
@ -525,11 +544,10 @@ namespace Heat_method_3 {
Matrix index_divergence;
Eigen::VectorXd solved_phi;
std::set<Index> source_index;
bool source_change_flag;
};
} // namespace Heat_method_3
} // namespace CGAL
#include <CGAL/enable_warnings.h>
#endif CGAL_HEAT_METHOD_3_HEAT_METHOD_3_H

40
Heat_method_3/test/Heat_method_3/heat_method_surface_mesh_test.cpp
View File

@ -23,7 +23,7 @@ typedef boost::graph_traits<Mesh>::vertex_descriptor vertex_descriptor;
typedef CGAL::Heat_method_3::Heat_method_3<Mesh,Kernel,Vertex_distance_map> Heat_method;
typedef CGAL::Heat_method_3::Heat_method_Eigen_traits_3::SparseMatrix SparseMatrix;
void source_set_tests(Heat_method hm, Mesh sm)
void source_set_tests(Heat_method hm, const Mesh& sm)
{
vertex_descriptor source = *(vertices(sm).first);
hm.add_source(source);
@ -41,7 +41,7 @@ void source_set_tests(Heat_method hm, Mesh sm)
assert(hm.remove_source(*(std::next(vertices(sm).first,3))));
}
void cotan_matrix_test(SparseMatrix c)
void cotan_matrix_test(const SparseMatrix& c)
{
double sum = 0;
for(int k = 0; k<c.outerSize(); ++k)
@ -55,7 +55,7 @@ void cotan_matrix_test(SparseMatrix c)
assert(sum < 0.000000001);
}
void mass_matrix_test(SparseMatrix M)
void mass_matrix_test(const SparseMatrix& M)
{
double sum = 0;
for(int k = 0; k<M.outerSize(); ++k)
@ -66,12 +66,12 @@ void mass_matrix_test(SparseMatrix M)
}
}
//total Area matrix should be equal to the sum of all faces on the mesh
//have to allow for the error because of rounding issues: Andreas might be able to help with this?
//have to allow for the error because of rounding
//this will only work for the pyramid mesh
assert((sum-1.866025)<=0.000005);
}
void check_for_zero(Eigen::VectorXd u)
void check_for_zero(const Eigen::VectorXd& u)
{
for(int c_i = 0; c_i<4; c_i++)
{
@ -79,7 +79,7 @@ void check_for_zero(Eigen::VectorXd u)
}
}
void check_for_unit(Eigen::MatrixXd X, int dimension)
void check_for_unit(const Eigen::MatrixXd& X, int dimension)
{
for(int k = 0; k<dimension; k++)
{
@ -88,10 +88,13 @@ void check_for_unit(Eigen::MatrixXd X, int dimension)
}
}
void check_no_update(const Mesh& sm, const Vertex_distance_map& original, const Vertex_distance_map& updated)
{
BOOST_FOREACH(vertex_descriptor vd, vertices(sm))
{
assert(get(original, vd) == get(updated,vd));
}
}
int main()
{
@ -135,9 +138,6 @@ int main()
Eigen::VectorXd solved_dist = hm.solve_phi(c, XD,4);
std::cout<<"PHASE 1 DONE \n";
std::cout<<"PHASE 2 DONE \n";
std::cout<<"PHASE 3 DONE \n";
Mesh sm2;
Vertex_distance_map vertex_distance_map2 = get(Vertex_distance_tag(),sm2);
@ -169,7 +169,6 @@ int main()
//verified a few of the actual values against the estimated values, avg. error was 0.0001
//In future, want to check performance against other solver
Mesh sm3;
Vertex_distance_map vertex_distance_map3 = get(Vertex_distance_tag(),sm3);
@ -186,11 +185,18 @@ int main()
const SparseMatrix& c3 = hm3.cotan_matrix();
cotan_matrix_test(c3);
const SparseMatrix& K3= hm3.kronecker_delta();
// AF: I commented the assert as I commented in build()
assert(K3.nonZeros()==1);
hm3.add_source(*(++(++(vertices(sm3).first))));
hm3.add_source(*(vertices(sm3).first));
const Vertex_distance_map& old_vdm = hm3.get_vertex_distance_map();
hm3.update();
const Vertex_distance_map& original_vdm = hm3.get_vertex_distance_map();
hm.update();
const Vertex_distance_map& new_vdm = hm3.get_vertex_distance_map();
check_no_update(sm3, original_vdm, new_vdm);
const SparseMatrix& K4 = hm3.kronecker_delta();
assert(K4.nonZeros()==2);
return 0;
}