From 470ce56df5acce8b39665aa7af811cde99d4126a Mon Sep 17 00:00:00 2001 From: Andreas Fabri Date: Tue, 30 Oct 2018 22:31:24 +0100 Subject: [PATCH] Make it work with Heat_method_traits_3 --- .../Concepts/HeatMethodTraits_3.h | 4 + .../doc/Heat_method_3/Heat_method_3.txt | 3 +- .../doc/Heat_method_3/PackageDescription.txt | 2 + .../examples/Heat_method_3/CMakeLists.txt | 1 + .../Surface_mesh_geodesic_distances_3.h | 121 ++++++++++------- .../include/CGAL/Heat_method_3/internal/V2V.h | 4 +- .../test/Heat_method_3/CMakeLists.txt | 1 + .../test/Heat_method_3/Heat_method_traits_3.h | 123 ++++++++++++++++++ .../Heat_method_3/heat_method_concept.cpp | 63 +++++++++ .../Display/Display_property_plugin.cpp | 10 +- 10 files changed, 280 insertions(+), 52 deletions(-) create mode 100644 Heat_method_3/test/Heat_method_3/Heat_method_traits_3.h create mode 100644 Heat_method_3/test/Heat_method_3/heat_method_concept.cpp diff --git a/Heat_method_3/doc/Heat_method_3/Concepts/HeatMethodTraits_3.h b/Heat_method_3/doc/Heat_method_3/Concepts/HeatMethodTraits_3.h index c2860c5c908..f8303a0c46b 100644 --- a/Heat_method_3/doc/Heat_method_3/Concepts/HeatMethodTraits_3.h +++ b/Heat_method_3/doc/Heat_method_3/Concepts/HeatMethodTraits_3.h @@ -32,7 +32,11 @@ public: typedef unspecified_type Compute_y_3; typedef unspecified_type Compute_z_3; typedef unspecified_type Construct_vector_3; + typedef unspecified_type Construct_sum_of_vectors_3; + typedef unspecified_type Construct_scaled_vector_3; typedef unspecified_type Construct_cross_product_vector_3; + typedef unspecified_type Compute_scalar_product_3; + typedef unspecified_type Compute_squared_distance_3; /// @} diff --git a/Heat_method_3/doc/Heat_method_3/Heat_method_3.txt b/Heat_method_3/doc/Heat_method_3/Heat_method_3.txt index 6c34aad4169..169fb1cccb2 100644 --- a/Heat_method_3/doc/Heat_method_3/Heat_method_3.txt +++ b/Heat_method_3/doc/Heat_method_3/Heat_method_3.txt @@ -33,7 +33,8 @@ step will roughly double the overall preprocessing cost. In the next section we give some examples. Section \ref sec_HM_definitions presents the mathematical theory of the Heat method. The last section is about the \ref sec_HM_history. -Note that this package requires the third party \ref thirdpartyEigen library. +Note that this package depends on the third party \ref thirdpartyEigen library, or another +model of the concept `SparseLinearAlgebraWithFactorTraits_d`. This implementation is based on \cgalCite{cgal:cww-ghnac-13} and \cgalCite{cgal:fsbs-acidt-06} diff --git a/Heat_method_3/doc/Heat_method_3/PackageDescription.txt b/Heat_method_3/doc/Heat_method_3/PackageDescription.txt index be0af7031f0..4bf66859754 100644 --- a/Heat_method_3/doc/Heat_method_3/PackageDescription.txt +++ b/Heat_method_3/doc/Heat_method_3/PackageDescription.txt @@ -37,6 +37,8 @@ source vertices. } ## Classes ## - `CGAL::Heat_method_3::Surface_mesh_geodesic_distances_3` +- `CGAL::Heat_method_3::Direct` +- `CGAL::Heat_method_3::Intrinsic_Delaunay` ## Functions ## diff --git a/Heat_method_3/examples/Heat_method_3/CMakeLists.txt b/Heat_method_3/examples/Heat_method_3/CMakeLists.txt index 2b6e1f2ebe6..ec735680137 100644 --- a/Heat_method_3/examples/Heat_method_3/CMakeLists.txt +++ b/Heat_method_3/examples/Heat_method_3/CMakeLists.txt @@ -54,6 +54,7 @@ include_directories( BEFORE include ) include( CGAL_CreateSingleSourceCGALProgram ) create_single_source_cgal_program( "heat_method.cpp" ) +create_single_source_cgal_program( "heat_method_concept.cpp" ) create_single_source_cgal_program( "heat_method_polyhedron.cpp" ) create_single_source_cgal_program( "heat_method_surface_mesh.cpp" ) create_single_source_cgal_program( "heat_method_surface_mesh_intrinsic.cpp" ) diff --git a/Heat_method_3/include/CGAL/Heat_method_3/Surface_mesh_geodesic_distances_3.h b/Heat_method_3/include/CGAL/Heat_method_3/Surface_mesh_geodesic_distances_3.h index 69815dbf231..64b20edf5dc 100644 --- a/Heat_method_3/include/CGAL/Heat_method_3/Surface_mesh_geodesic_distances_3.h +++ b/Heat_method_3/include/CGAL/Heat_method_3/Surface_mesh_geodesic_distances_3.h @@ -231,6 +231,7 @@ private: double summation_of_edges() const { + typename Traits::Compute_squared_distance_3 squared_distance = Traits().compute_squared_distance_3_object(); double edge_sum = 0; CGAL::Vertex_around_face_iterator vbegin, vend, vmiddle; BOOST_FOREACH(face_descriptor f, faces(tm)) { @@ -245,9 +246,9 @@ private: VertexPointMap_reference pi = get(vpm,current); VertexPointMap_reference pj = get(vpm, neighbor_one); VertexPointMap_reference pk = get(vpm, neighbor_two); - edge_sum += (CGAL::is_border(opposite(hd,tm),tm)?1.0:0.5) * std::sqrt(CGAL::squared_distance(pi,pj)); - edge_sum += (CGAL::is_border(opposite(hd2,tm),tm)?1.0:0.5) * std::sqrt(CGAL::squared_distance(pj,pk)) ; - edge_sum += (CGAL::is_border(opposite(hd3,tm),tm)?1.0:0.5) * std::sqrt(CGAL::squared_distance(pk,pi)) ; + edge_sum += (CGAL::is_border(opposite(hd,tm),tm)?1.0:0.5) * std::sqrt(squared_distance(pi,pj)); + edge_sum += (CGAL::is_border(opposite(hd2,tm),tm)?1.0:0.5) * std::sqrt(squared_distance(pj,pk)) ; + edge_sum += (CGAL::is_border(opposite(hd3,tm),tm)?1.0:0.5) * std::sqrt(squared_distance(pk,pi)) ; } return edge_sum; } @@ -312,6 +313,11 @@ private: void compute_unit_gradient() { + Traits::Construct_vector_3 construct_vector = Traits().construct_vector_3_object(); + Traits::Construct_sum_of_vectors_3 sum = Traits().construct_sum_of_vectors_3_object(); + Traits::Compute_scalar_product_3 scalar_product = Traits().compute_scalar_product_3_object(); + Traits::Construct_cross_product_vector_3 cross_product = Traits().construct_cross_product_vector_3_object(); + Traits::Construct_scaled_vector_3 scale = Traits().construct_scaled_vector_3_object(); if(X.empty()){ X.resize(num_faces(tm)); } @@ -335,9 +341,13 @@ private: //cross that with eij, ejk, eki //so (Ncross eij) *uk and so on //sum all of those then multiply by 1./(2a) - Vector_3 cross = CGAL::cross_product((p_j-p_i), (p_k-p_i)); - double N_cross = (CGAL::sqrt(cross*cross)); - Vector_3 unit_cross = cross/N_cross; + + Vector_3 v_ij = construct_vector(p_i,p_j); + Vector_3 v_ik = construct_vector(p_i,p_k); + + Vector_3 cross = cross_product(v_ij, v_ik); + double N_cross = (CGAL::sqrt(scalar_product(cross,cross))); + Vector_3 unit_cross = scale(cross, 1./N_cross); double area_face = N_cross * (1./2); double u_i = CGAL::abs(solved_u(i)); double u_j = CGAL::abs(solved_u(j)); @@ -349,25 +359,24 @@ private: u_j = u_j * r_Mag; u_k = u_k * r_Mag; } - Vector_3 edge_sums = u_k * CGAL::cross_product(unit_cross,(p_j-p_i)); - edge_sums = edge_sums + u_i * (CGAL::cross_product(unit_cross, (p_k-p_j))); - edge_sums = edge_sums + u_j * CGAL::cross_product(unit_cross, (p_i-p_k)); - edge_sums = edge_sums * (1./area_face); - double e_magnitude = CGAL::sqrt(edge_sums*edge_sums); - X[face_i] = edge_sums*(1./e_magnitude); + Vector_3 edge_sums = scale(cross_product(unit_cross,v_ij), u_k); + edge_sums = sum(edge_sums, scale(cross_product(unit_cross, construct_vector(p_j,p_k)), u_i)); + edge_sums = sum(edge_sums, scale(cross_product(unit_cross, construct_vector(p_k,p_i)), u_j)); + edge_sums = scale(edge_sums, (1./area_face)); + double e_magnitude = CGAL::sqrt(scalar_product(edge_sums,edge_sums)); + X[face_i] = scale(edge_sums,(1./e_magnitude)); } } - double - dot_eigen_vector(const Vector_3& a, const Vector_3& b) const - { - return (a.x()*b.x() + a.y()*b.y() + a.z()*b.z()); - } + void compute_divergence() { + Traits::Construct_cross_product_vector_3 cross_product = Traits().construct_cross_product_vector_3_object(); + Traits::Compute_scalar_product_3 scalar_product = Traits().compute_scalar_product_3_object(); + Traits::Construct_vector_3 construct_vector = Traits().construct_vector_3_object(); Matrix indexD(dimension,1); CGAL::Vertex_around_face_iterator vbegin, vend, vmiddle; BOOST_FOREACH(face_descriptor f, faces(tm)) { @@ -382,24 +391,32 @@ private: VertexPointMap_reference p_j = get(vpm, neighbor_one); VertexPointMap_reference p_k = get(vpm, neighbor_two); Index face_i = get(face_id_map, f); - - Vector_3 cross = CGAL::cross_product((p_j-p_i), (p_k-p_i)); - double norm_cross = (CGAL::sqrt(cross*cross)); - double dot = (p_j-p_i)*(p_k-p_i); + + Vector_3 v_ij = construct_vector(p_i,p_j); + Vector_3 v_ik = construct_vector(p_i,p_k); + Vector_3 cross = cross_product(v_ij, v_ik); + double norm_cross = CGAL::sqrt(scalar_product(cross,cross)); + double dot = scalar_product(v_ij, v_ik); double cotan_i = dot/norm_cross; - cross = CGAL::cross_product((p_i-p_j), (p_k-p_j)); - dot = to_double((p_i-p_j)*(p_k-p_j)); + Vector_3 v_ji = construct_vector(p_j, p_i); + Vector_3 v_jk = construct_vector(p_j, p_k); + + cross = cross_product(v_ji, v_jk); + dot = to_double(scalar_product(v_ji, v_jk)); double cotan_j = dot/norm_cross; - cross = CGAL::cross_product((p_i-p_k), (p_j-p_k)); - dot = to_double((p_i-p_k)*(p_j-p_k)); + Vector_3 v_ki = construct_vector(p_k,p_i); + Vector_3 v_kj = construct_vector(p_k,p_j); + + cross = cross_product(v_ki, v_kj); + dot = to_double(scalar_product(v_ki,v_kj)); double cotan_k = dot/norm_cross; const Vector_3& a = X[face_i]; - double i_entry = cotan_k*(dot_eigen_vector(a,(p_j-p_i))) + cotan_j*(dot_eigen_vector(a,(p_k-p_i))); - double j_entry = cotan_i*(dot_eigen_vector(a,(p_k-p_j))) + cotan_k*(dot_eigen_vector(a,(p_i-p_j))); - double k_entry = cotan_j*(dot_eigen_vector(a,(p_i-p_k))) + cotan_i*(dot_eigen_vector(a,(p_j-p_k))); + double i_entry = (scalar_product(a,v_ij) * cotan_k) + (scalar_product(a,v_ik) * cotan_j); + double j_entry = (scalar_product(a,v_jk) * cotan_i) + (scalar_product(a,v_ji) * cotan_k); + double k_entry = (scalar_product(a,v_ki) * cotan_j) + (scalar_product(a,v_kj) * cotan_i); indexD.add_coef(i, 0, (1./2)*i_entry); indexD.add_coef(j, 0, (1./2)*j_entry); @@ -499,6 +516,10 @@ private: void build() { + Traits::Construct_cross_product_vector_3 cross_product = Traits().construct_cross_product_vector_3_object(); + Traits::Compute_scalar_product_3 scalar_product = Traits().compute_scalar_product_3_object(); + Traits::Construct_vector_3 construct_vector = Traits().construct_vector_3_object(); + source_change_flag = false; CGAL_precondition(is_triangle_mesh(tm)); @@ -540,27 +561,35 @@ private: pj = p_j; pk = p_k; - Vector_3 cross = CGAL::cross_product((pj-pi), (pk-pi)); - double dot = (pj-pi)*(pk-pi); - - double norm_cross = (CGAL::sqrt(cross*cross)); - + Vector_3 v_ij = construct_vector(p_i,p_j); + Vector_3 v_ik = construct_vector(p_i,p_k); + + Vector_3 cross = cross_product(v_ij, v_ik); + double dot = scalar_product(v_ij,v_ik); + + double norm_cross = (CGAL::sqrt(scalar_product(cross,cross))); + double cotan_i = dot/norm_cross; m_cotan_matrix.add_coef(j,k ,-(1./2)*cotan_i); m_cotan_matrix.add_coef(k,j,-(1./2)* cotan_i); m_cotan_matrix.add_coef(j,j,(1./2)*cotan_i); m_cotan_matrix.add_coef(k,k,(1./2)* cotan_i); - cross = CGAL::cross_product((pi-pj), (pk-pj)); - dot = to_double((pi-pj)*(pk-pj)); + Vector_3 v_ji = construct_vector(p_j,p_i); + Vector_3 v_jk = construct_vector(p_j,p_k); + + cross = cross_product(v_ji, v_jk); + dot = to_double(scalar_product(v_ji, v_jk)); double cotan_j = dot/norm_cross; m_cotan_matrix.add_coef(i,k ,-(1./2)*cotan_j); m_cotan_matrix.add_coef(k,i,-(1./2)* cotan_j); m_cotan_matrix.add_coef(i,i,(1./2)* cotan_j); m_cotan_matrix.add_coef(k,k,(1./2)* cotan_j); - cross = CGAL::cross_product((pi-pk), (pj-pk)); - dot = to_double((pi-pk)*(pj-pk)); + Vector_3 v_ki = construct_vector(p_k,p_i); + Vector_3 v_kj = construct_vector(p_k,p_j); + cross = cross_product(v_ki, v_kj); + dot = to_double(scalar_product(v_ki,v_kj)); double cotan_k = dot/norm_cross; m_cotan_matrix.add_coef(i,j,-(1./2)*cotan_k); m_cotan_matrix.add_coef(j,i,-(1./2)* cotan_k); @@ -604,9 +633,7 @@ private: Vector solved_phi; std::set source_index; bool source_change_flag; - - std::vector > m_triplets; - + }; template ::vertex_descriptor` as key type and `double` as value type. /// \tparam Mode either the tag `Direct` or `Intrinsic_Delaunay`, which determines if the geodesic distance @@ -867,7 +895,7 @@ public: /// The default is `Direct`. /// /// \sa CGAL::Heat_method_3::Surface_mesh_geodesic_distances_3 -template +template void estimate_geodesic_distances(const TriangleMesh& tm, VertexDistanceMap vdm, @@ -879,7 +907,8 @@ estimate_geodesic_distances(const TriangleMesh& tm, hm.estimate_geodesic_distances(vdm); } - + +#ifndef DOXYGEN_RUNNING template void estimate_geodesic_distances(const TriangleMesh& tm, @@ -890,11 +919,13 @@ estimate_geodesic_distances(const TriangleMesh& tm, hm.add_source(source); hm.estimate_geodesic_distances(vdm); } - +#endif + - /// \ingroup PkgHeatMethod +/// \ingroup PkgHeatMethod /// computes for each vertex of the triangle mesh `tm` the estimated geodesic distance to a given set of source vertices. /// \tparam TriangleMesh a triangulated surface mesh, model of `FaceListGraph` and `HalfedgeListGraph` +/// It must have an internal vertex point property map with the value type being a 3D point from a cgal Kernel model /// \tparam VertexDistanceMap a property map model of `WritablePropertyMap` /// with `boost::graph_traits::vertex_descriptor` as key type and `double` as value type. /// \tparam VertexRange a range with value type `boost::graph_traits::vertex_descriptor` diff --git a/Heat_method_3/include/CGAL/Heat_method_3/internal/V2V.h b/Heat_method_3/include/CGAL/Heat_method_3/internal/V2V.h index 2380e4b8f7c..ddc9d22bc70 100644 --- a/Heat_method_3/include/CGAL/Heat_method_3/internal/V2V.h +++ b/Heat_method_3/include/CGAL/Heat_method_3/internal/V2V.h @@ -24,6 +24,8 @@ #include +#ifndef DOXYGEN_RUNNING + namespace CGAL { namespace Heat_method_3 { @@ -39,7 +41,7 @@ namespace Heat_method_3 { return t; } }; - +#endif } // namespace Heat_method_3 } // namespace CGAL #endif //CGAL_HEAT_METHOD_3_INTERNAL_V2V_H diff --git a/Heat_method_3/test/Heat_method_3/CMakeLists.txt b/Heat_method_3/test/Heat_method_3/CMakeLists.txt index 9251ad24227..56d28fc5269 100644 --- a/Heat_method_3/test/Heat_method_3/CMakeLists.txt +++ b/Heat_method_3/test/Heat_method_3/CMakeLists.txt @@ -53,6 +53,7 @@ include( CGAL_CreateSingleSourceCGALProgram ) create_single_source_cgal_program( "validate.cpp" ) create_single_source_cgal_program( "idt.cpp" ) +create_single_source_cgal_program( "heat_method_concept.cpp" ) create_single_source_cgal_program( "heat_method_surface_mesh_test.cpp" ) create_single_source_cgal_program( "heat_method_intrinsic_surface_mesh_test.cpp" ) create_single_source_cgal_program( "intrinsic_delaunay_triangulation_test.cpp" ) diff --git a/Heat_method_3/test/Heat_method_3/Heat_method_traits_3.h b/Heat_method_3/test/Heat_method_3/Heat_method_traits_3.h new file mode 100644 index 00000000000..1d67423e4f9 --- /dev/null +++ b/Heat_method_3/test/Heat_method_3/Heat_method_traits_3.h @@ -0,0 +1,123 @@ + +struct Heat_method_traits_3 +{ + typedef double FT; + + struct Point_3 { + Point_3() + {} + + Point_3(double, double,double) + {} + Point_3(const Point_3&) + {} + }; + + + struct Vector_3 { + Vector_3() + {} + + Vector_3(double, double, double) + {} + Vector_3(const Vector_3&) + {} + + + }; + + + struct Point_2 { + Point_2() + {} + + Point_2(double, double) + {} + Point_2(const Point_2&) + {} + }; + + struct Compute_x_2 { + double operator()(const Point_2&)const + {return 0;} + }; + + struct Compute_y_2 { + double operator()(const Point_2&)const + {return 0;} + }; + + + struct Construct_vector_3{ + Vector_3 operator()(const Point_3&, const Point_3&) const + { return Vector_3();} + }; + + struct Construct_scaled_vector_3{ + Vector_3 operator()(const Vector_3&, double) const + { return Vector_3();} + }; + + struct Construct_cross_product_vector_3 { + + Vector_3 operator()(const Vector_3&, const Vector_3&) const + { return Vector_3();} + }; + + struct Construct_sum_of_vectors_3 { + + Vector_3 operator()(const Vector_3&, const Vector_3&) const + { return Vector_3();} + }; + struct Compute_scalar_product_3 { + + double operator()(const Vector_3&, const Vector_3&) const + { return 0;} + }; + + struct Compute_squared_distance_3 { + + double operator()(const Point_3&, const Point_3&) const + { return 0;} + }; + + Compute_x_2 compute_x_2_object() const + { + return Compute_x_2(); + } + + Compute_y_2 compute_y_2_object() const + { + return Compute_y_2(); + } + + + Construct_vector_3 construct_vector_3_object() const + { + return Construct_vector_3(); + } + Construct_scaled_vector_3 construct_scaled_vector_3_object() const + { + return Construct_scaled_vector_3(); + } + + Construct_cross_product_vector_3 construct_cross_product_vector_3_object() const + { + return Construct_cross_product_vector_3(); + } + + Compute_squared_distance_3 compute_squared_distance_3_object() const + { + return Compute_squared_distance_3(); + } + + Compute_scalar_product_3 compute_scalar_product_3_object() const + { + return Compute_scalar_product_3(); + } + + Construct_sum_of_vectors_3 construct_sum_of_vectors_3_object() const + { + return Construct_sum_of_vectors_3(); + } +}; diff --git a/Heat_method_3/test/Heat_method_3/heat_method_concept.cpp b/Heat_method_3/test/Heat_method_3/heat_method_concept.cpp new file mode 100644 index 00000000000..1b79626c681 --- /dev/null +++ b/Heat_method_3/test/Heat_method_3/heat_method_concept.cpp @@ -0,0 +1,63 @@ +#include +#include +#include + +#include "Heat_method_3/Heat_method_traits_3.h" + +#include +#include + +#include + +typedef Heat_method_traits_3 Kernel; +typedef Kernel::Point_3 Point_3; +typedef CGAL::Surface_mesh Surface_mesh; + +typedef boost::graph_traits::vertex_descriptor vertex_descriptor; +typedef Surface_mesh::Property_map Vertex_distance_map; +typedef CGAL::Heat_method_3::Surface_mesh_geodesic_distances_3::EigenType > >, + boost::property_map< Surface_mesh, CGAL::vertex_point_t>::const_type, + Kernel> Heat_method; + + + +int main(int argc, char* argv[]) +{ + //read in mesh + Surface_mesh sm; + + //the heat intensity will hold the distance values from the source set + Vertex_distance_map heat_intensity = sm.add_property_map("v:heat_intensity", 0).first; + + Heat_method hm(sm); + + //add the first vertex as the source set + vertex_descriptor source = *(vertices(sm).first); + hm.add_source(source); + hm.estimate_geodesic_distances(heat_intensity); + + Point_3 sp = sm.point(source); + + vertex_descriptor far; + double sdistance = 0; + + BOOST_FOREACH(vertex_descriptor vd , vertices(sm)){ + std::cout << vd << " is at distance " << get(heat_intensity, vd) << " from " << source << std::endl; + /* + if(squared_distance(sp,sm.point(vd)) > sdistance){ + far = vd; + sdistance = squared_distance(sp,sm.point(vd)); + } + */ + } + hm.add_source(far); + hm.estimate_geodesic_distances(heat_intensity); + + BOOST_FOREACH(vertex_descriptor vd , vertices(sm)){ + std::cout << vd << " is at distance " << get(heat_intensity, vd) << " " << std::endl; + } + + return 0; +} diff --git a/Polyhedron/demo/Polyhedron/Plugins/Display/Display_property_plugin.cpp b/Polyhedron/demo/Polyhedron/Plugins/Display/Display_property_plugin.cpp index c99fbe8f3c9..6b4b774e992 100644 --- a/Polyhedron/demo/Polyhedron/Plugins/Display/Display_property_plugin.cpp +++ b/Polyhedron/demo/Polyhedron/Plugins/Display/Display_property_plugin.cpp @@ -10,7 +10,7 @@ #include #include #include -#include +#include #include "Scene_points_with_normal_item.h" #include "Messages_interface.h" @@ -309,8 +309,8 @@ class DisplayPropertyPlugin : Q_INTERFACES(CGAL::Three::Polyhedron_demo_plugin_interface) Q_PLUGIN_METADATA(IID "com.geometryfactory.PolyhedronDemo.PluginInterface/1.0") typedef SMesh::Property_map::vertex_descriptor, double> Vertex_distance_map; - typedef CGAL::Heat_method_3::Heat_method_3 Heat_method; - typedef CGAL::Heat_method_3::Heat_method_3 Heat_method_idt; + typedef CGAL::Heat_method_3::Surface_mesh_geodesic_distances_3 Heat_method; + typedef CGAL::Heat_method_3::Surface_mesh_geodesic_distances_3 Heat_method_idt; public: @@ -769,9 +769,9 @@ private Q_SLOTS: } if(iDT){ - hm_idt->fill_distance_map(heat_intensity); + hm_idt->estimate_geodesic_distances(heat_intensity); }else{ - hm->fill_distance_map(heat_intensity); + hm->estimate_geodesic_distances(heat_intensity); } double max = 0;