Author: renard Date: Fri Oct 14 09:37:49 2016 New Revision: 5417 URL: http://svn.gna.org/viewcvs/getfem?rev=5417&view=rev Log: transfer Cuthill Mc Kee ordering in mesh::optimize_structure. The dof numbering is modified
Modified: trunk/getfem/contrib/opt_assembly/opt_assembly.cc trunk/getfem/interface/src/gf_mesh.cc trunk/getfem/interface/src/gf_mesh_set.cc trunk/getfem/interface/tests/matlab/check_mesh_fem.m trunk/getfem/interface/tests/matlab/demo_laplacian.m trunk/getfem/src/getfem/bgeot_mesh_structure.h trunk/getfem/src/getfem/getfem_mesh.h trunk/getfem/src/getfem_mesh.cc trunk/getfem/src/getfem_mesh_fem.cc trunk/getfem/src/getfem_regular_meshes.cc trunk/getfem/src/gmm/gmm_vector.h trunk/getfem/tests/test_assembly.cc trunk/getfem/tests/test_mesh.cc Modified: trunk/getfem/contrib/opt_assembly/opt_assembly.cc URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/contrib/opt_assembly/opt_assembly.cc?rev=5417&r1=5416&r2=5417&view=diff ============================================================================== --- trunk/getfem/contrib/opt_assembly/opt_assembly.cc (original) +++ trunk/getfem/contrib/opt_assembly/opt_assembly.cc Fri Oct 14 09:37:49 2016 @@ -433,35 +433,35 @@ // - Instructions execution except for assembly ones // new | old | sto | asse | exec | J | Ins | test_new_assembly(2, 400, 1); // ndofu = 321602 ndofp = 160801 ndofchi = 1201 - // Mass (scalar) : 0.36 | 0.61 | 0.07 | 0.10 | 0.20 | 0.06 | 0.06 | - // Mass (vector) : 0.46 | 0.82 | 0.12 | 0.22 | 0.15 | 0.05 | 0.09 | - // Laplacian : 0.26 | 0.83 | 0.05 | 0.07 | 0.14 | 0.06 | 0.05 | - // Homogeneous elas : 0.46 | 1.88 | 0.14 | 0.22 | 0.13 | 0.05 | 0.11 | - // Non-homogeneous elast: 0.57 | 2.26 | 0.14 | 0.22 | 0.13 | 0.05 | 0.21 | + // Mass (scalar) : 0.34 | 0.61 | 0.05 | 0.08 | 0.20 | 0.06 | 0.06 | + // Mass (vector) : 0.43 | 0.82 | 0.10 | 0.20 | 0.14 | 0.05 | 0.09 | + // Laplacian : 0.23 | 0.83 | 0.02 | 0.04 | 0.14 | 0.06 | 0.05 | + // Homogeneous elas : 0.43 | 1.88 | 0.11 | 0.19 | 0.13 | 0.05 | 0.11 | + // Non-homogeneous elast: 0.53 | 2.26 | 0.10 | 0.18 | 0.13 | 0.05 | 0.21 | test_new_assembly(3, 36, 1); // ndofu = 151959 ndofp = 50653 ndofchi = 6553 - // Mass (scalar) : 0.44 | 0.77 | 0.08 | 0.12 | 0.24 | 0.11 | 0.08 | - // Mass (vector) : 0.94 | 1.57 | 0.18 | 0.39 | 0.20 | 0.11 | 0.35 | - // Laplacian : 0.41 | 1.38 | 0.06 | 0.09 | 0.19 | 0.11 | 0.14 | - // Homogeneous elas : 1.28 | 4.58 | 0.48 | 0.59 | 0.19 | 0.11 | 0.51 | - // Non-homogeneous elast: 1.42 | 6.81 | 0.49 | 0.63 | 0.18 | 0.11 | 0.61 | + // Mass (scalar) : 0.42 | 0.77 | 0.06 | 0.10 | 0.24 | 0.11 | 0.08 | + // Mass (vector) : 0.92 | 1.54 | 0.17 | 0.38 | 0.19 | 0.11 | 0.35 | + // Laplacian : 0.39 | 1.38 | 0.04 | 0.07 | 0.19 | 0.11 | 0.14 | + // Homogeneous elas : 1.26 | 4.58 | 0.46 | 0.57 | 0.19 | 0.11 | 0.51 | + // Non-homogeneous elast: 1.38 | 6.73 | 0.45 | 0.59 | 0.18 | 0.11 | 0.61 | test_new_assembly(2, 200, 2); // ndofu = 321602 ndofp = 160801 ndofchi = 1201 - // Mass (scalar) : 0.15 | 0.25 | 0.04 | 0.07 | 0.05 | 0.02 | 0.03 | - // Mass (vector) : 0.34 | 0.45 | 0.07 | 0.15 | 0.05 | 0.02 | 0.14 | + // Mass (scalar) : 0.14 | 0.25 | 0.04 | 0.07 | 0.04 | 0.02 | 0.03 | + // Mass (vector) : 0.34 | 0.44 | 0.07 | 0.15 | 0.05 | 0.02 | 0.14 | // Laplacian : 0.13 | 0.37 | 0.03 | 0.05 | 0.04 | 0.02 | 0.04 | // Homogeneous elas : 0.37 | 1.28 | 0.13 | 0.19 | 0.04 | 0.01 | 0.14 | - // Non-homogeneous elast: 0.45 | 2.40 | 0.13 | 0.18 | 0.04 | 0.01 | 0.23 | + // Non-homogeneous elast: 0.44 | 2.40 | 0.12 | 0.17 | 0.04 | 0.01 | 0.23 | test_new_assembly(3, 18, 2); // ndofu = 151959 ndofp = 50653 ndofchi = 6553 - // Mass (scalar) : 0.23 | 0.30 | 0.10 | 0.15 | 0.04 | 0.02 | 0.04 | + // Mass (scalar) : 0.23 | 0.29 | 0.10 | 0.15 | 0.04 | 0.02 | 0.04 | // Mass (vector) : 1.46 | 0.90 | 0.33 | 0.64 | 0.04 | 0.02 | 0.78 | // Laplacian : 0.18 | 0.55 | 0.08 | 0.10 | 0.03 | 0.02 | 0.05 | // Homogeneous elas : 2.37 | 3.47 | 1.20 | 1.37 | 0.03 | 0.02 | 0.97 | - // Non-homogeneous elast: 2.45 | 9.25 | 1.20 | 1.37 | 0.03 | 0.02 | 1.05 | + // Non-homogeneous elast: 2.44 | 9.25 | 1.20 | 1.37 | 0.03 | 0.02 | 1.04 | test_new_assembly(3, 9, 4); // ndofu = 151959 ndofp = 50653 ndofchi = 6553 - // Mass (scalar) : 0.75 | 0.38 | 0.25 | 0.39 | 0.01 | .005 | 0.35 | - // Mass (vector) : 5.10 | 1.33 | 0.42 | 1.86 | 0.01 | .005 | 3.23 | - // Laplacian : 0.61 | 0.79 | 0.25 | 0.34 | 0.01 | .005 | 0.26 | - // Homogeneous elas : 10.8 | 5.52 | 2.40 | 3.29 | 0.01 | .005 | 7.50 | - // Non-homogeneous elast: 11.0 | 48.0 | 2.35 | 3.19 | 0.01 | .005 | 7.80 | + // Mass (scalar) : 0.74 | 0.36 | 0.24 | 0.38 | 0.01 | .005 | 0.35 | + // Mass (vector) : 5.08 | 1.33 | 0.41 | 1.85 | 0.01 | .005 | 3.22 | + // Laplacian : 0.60 | 0.79 | 0.25 | 0.34 | 0.01 | .005 | 0.25 | + // Homogeneous elas : 10.7 | 5.47 | 2.40 | 3.29 | 0.01 | .005 | 7.40 | + // Non-homogeneous elast: 10.8 | 48.0 | 2.35 | 3.29 | 0.01 | .005 | 7.50 | // Conclusions : // - Desactivation of debug test has no sensible effect. Modified: trunk/getfem/interface/src/gf_mesh.cc URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/interface/src/gf_mesh.cc?rev=5417&r1=5416&r2=5417&view=diff ============================================================================== --- trunk/getfem/interface/src/gf_mesh.cc (original) +++ trunk/getfem/interface/src/gf_mesh.cc Fri Oct 14 09:37:49 2016 @@ -190,7 +190,7 @@ } pmesh->add_convex_by_points(pgt, pts.begin()); } - pmesh->optimize_structure(); + pmesh->optimize_structure(false); } getfem::base_small_vector diff(N); Modified: trunk/getfem/interface/src/gf_mesh_set.cc URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/interface/src/gf_mesh_set.cc?rev=5417&r1=5416&r2=5417&view=diff ============================================================================== --- trunk/getfem/interface/src/gf_mesh_set.cc (original) +++ trunk/getfem/interface/src/gf_mesh_set.cc Fri Oct 14 09:37:49 2016 @@ -372,7 +372,7 @@ ); - /*@SET ('optimize structure') + /*@SET ('optimize structure'[, @int with_renumbering]) Reset point and convex numbering. After optimisation, the points (resp. convexes) will @@ -381,8 +381,10 @@ (resp. MESH:GET('max cvid'))}@PYTHON{``0`` to ``MESH:GET('max pid')-1`` (resp. ``MESH:GET('max cvid')-1``)}.@*/ sub_command - ("optimize structure", 0, 0, 0, 0, - pmesh->optimize_structure(); + ("optimize structure", 0, 1, 0, 0, + bool with_renumbering = true; + if (in.remaining()) with_renumbering = (in.pop().to_integer(0,1) != 0); + pmesh->optimize_structure(with_renumbering); ); Modified: trunk/getfem/interface/tests/matlab/check_mesh_fem.m URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/interface/tests/matlab/check_mesh_fem.m?rev=5417&r1=5416&r2=5417&view=diff ============================================================================== --- trunk/getfem/interface/tests/matlab/check_mesh_fem.m (original) +++ trunk/getfem/interface/tests/matlab/check_mesh_fem.m Fri Oct 14 09:37:49 2016 @@ -167,11 +167,11 @@ gf_mesh_fem_get(mf,'nbdof') % should be 99 or 100 (element 0 and 1 are not really neigbhor but can be viewed as such) d=gf_mesh_fem_get(mf,'basic dof from cv',[1 5]) - gfassert(['d==[1 2 3 4 5 6 38 41 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58]']); + gfassert(['d==[1 2 3 4 5 6 29 32 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49]']); d=gf_mesh_fem_get(mf,'basic dof from cv',[1 5;1 2]) - gfassert('d==[3 5 6 38 41 43 46 48 51 53 56 58]'); + gfassert('d==[3 5 6 29 32 34 37 39 42 44 47 49]'); d=gf_mesh_fem_get(mf,'basic dof from cvid',5) - gfassert('d==[38 44 45 46 47 48 41 49 50 51 52 53 43 54 55 56 57 58]'); + gfassert('d==[29 35 36 37 38 39 32 40 41 42 43 44 34 45 46 47 48 49]'); s2=gf_mesh_get(mf,'char'); gfassert('length(s2)>500'); @@ -227,9 +227,9 @@ gf_mesh_set(m, 'boundary', 7, [3 4; 3 2]); cl=[1:5 7 8]; asserterr('gf_mesh_fem_get(mf_u, ''non conformal basic dof'')'); - d=gf_mesh_fem_get(mf_u, 'non conformal basic dof',cl); + d=gf_mesh_fem_get(mf_u, 'non conformal basic dof',cl) %gf_plot_mesh(mf_u, 'dof', 'on'); - gfassert('d==[11 12 13 14 15 16 21 22]'); + gfassert('d==[19 20 21 22 23 24 29 30]'); f=gf_mesh_fem_get(mf2, 'fem'); f5=gf_mesh_fem_get(mf2, 'fem',5); @@ -319,7 +319,7 @@ gfassert('ncv < maxcvid'); - gf_mesh_set(m,'optimize structure'); + gf_mesh_set(m,'optimize structure', false); maxpid=gf_mesh_get(m,'max pid'); maxcvid=gf_mesh_get(m,'max cvid'); Modified: trunk/getfem/interface/tests/matlab/demo_laplacian.m URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/interface/tests/matlab/demo_laplacian.m?rev=5417&r1=5416&r2=5417&view=diff ============================================================================== --- trunk/getfem/interface/tests/matlab/demo_laplacian.m (original) +++ trunk/getfem/interface/tests/matlab/demo_laplacian.m Fri Oct 14 09:37:49 2016 @@ -36,6 +36,7 @@ else m=gf_mesh('import','structured',sprintf('GT="GT_PK(2,1)";SIZES=[1,1];NOISED=0;NSUBDIV=[%d,%d];', NX, NX)); end +gf_mesh_set(m, 'optimize structure') % Create a mesh_fem of for a field of dimension 1 (i.e. a scalar field) mf = gf_mesh_fem(m,1); @@ -55,7 +56,8 @@ GAMMAD=1; gf_mesh_set(m, 'region', GAMMAD, border); if (draw) - gf_plot_mesh(m, 'regions', [GAMMAD]); % the boundary edges appear in red + figure(1); + gf_plot_mesh(m, 'regions', [GAMMAD], 'convexes', 'on'); % the boundary edges appear in red pause(1); end @@ -89,6 +91,7 @@ U = gf_model_get(md, 'variable', 'u'); if (draw) + figure(2); subplot(2,1,1); gf_plot(mf,U,'mesh','off'); colorbar; title('computed solution'); Modified: trunk/getfem/src/getfem/bgeot_mesh_structure.h URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/getfem/bgeot_mesh_structure.h?rev=5417&r1=5416&r2=5417&view=diff ============================================================================== --- trunk/getfem/src/getfem/bgeot_mesh_structure.h (original) +++ trunk/getfem/src/getfem/bgeot_mesh_structure.h Fri Oct 14 09:37:49 2016 @@ -253,7 +253,7 @@ /** Return the cuthill_mc_kee ordering on the convexes */ void APIDECL cuthill_mckee_on_convexes(const bgeot::mesh_structure &ms, - std::vector<size_type> &cmk); + std::vector<size_type> &cmk); template<class ITER> bool mesh_structure::is_convex_having_points(size_type ic, Modified: trunk/getfem/src/getfem/getfem_mesh.h URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/getfem/getfem_mesh.h?rev=5417&r1=5416&r2=5417&view=diff ============================================================================== --- trunk/getfem/src/getfem/getfem_mesh.h (original) +++ trunk/getfem/src/getfem/getfem_mesh.h Fri Oct 14 09:37:49 2016 @@ -443,7 +443,7 @@ void sup_convex_from_regions(size_type cv); /** Pack the mesh : renumber convexes and nodes such that there is no holes in their numbering. Do NOT do the Cuthill-McKee. */ - void optimize_structure(); + void optimize_structure(bool with_renumbering = true); /// Return the list of convex IDs for a Cuthill-McKee ordering const std::vector<size_type> &cuthill_mckee_ordering() const; /// Erase the mesh. Modified: trunk/getfem/src/getfem_mesh.cc URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/getfem_mesh.cc?rev=5417&r1=5416&r2=5417&view=diff ============================================================================== --- trunk/getfem/src/getfem_mesh.cc (original) +++ trunk/getfem/src/getfem_mesh.cc Fri Oct 14 09:37:49 2016 @@ -217,10 +217,9 @@ } #endif - void mesh::optimize_structure() { + void mesh::optimize_structure(bool with_renumbering) { pts.resort(); - size_type i, j; - j = nb_convex(); + size_type i, j = nb_convex(), nbc = j; for (i = 0; i < j; i++) if (!convex_tab.index_valid(i)) swap_convex(i, convex_tab.ind_last()); @@ -231,14 +230,20 @@ while (i < j && j != ST_NIL && !(pts.index()[j])) --j; if (i < j && j != ST_NIL ) swap_points(i, j); } - } - - const std::vector<size_type> &mesh::cuthill_mckee_ordering() const { - if (!cuthill_mckee_uptodate) { - bgeot::cuthill_mckee_on_convexes(*this, cmk_order); - cuthill_mckee_uptodate = true; - } - return cmk_order; + if (with_renumbering) { // Could be optimized no using only swap_convex + std::vector<size_type> cmk, iord(nbc), iordinv(nbc); + for (i = 0; i < nbc; ++i) iord[i] = iordinv[i] = i; + + bgeot::cuthill_mckee_on_convexes(*this, cmk); + for (i = 0; i < nbc; ++i) { + j = iordinv[cmk[i]]; + if (i != j) { + swap_convex(i, j); + std::swap(iord[i], iord[j]); + std::swap(iordinv[iord[i]], iordinv[iord[j]]); + } + } + } } void mesh::translation(const base_small_vector &V) Modified: trunk/getfem/src/getfem_mesh_fem.cc URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/getfem_mesh_fem.cc?rev=5417&r1=5416&r2=5417&view=diff ============================================================================== --- trunk/getfem/src/getfem_mesh_fem.cc (original) +++ trunk/getfem/src/getfem_mesh_fem.cc Fri Oct 14 09:37:49 2016 @@ -307,7 +307,15 @@ if (first_pf && first_pf->is_on_real_element()) is_uniform_ = false; // Gives the Cuthill McKee ordering to iterate on elements - const std::vector<size_type> &cmk = linked_mesh().cuthill_mckee_ordering(); + // const std::vector<size_type> &cmk = linked_mesh().cuthill_mckee_ordering(); + + + // std::vector<size_type> cmk; + // for (dal::bv_visitor cv(linked_mesh().convex_index()); !cv.finished(); ++cv) + // cmk.push_back(cv); + + + // Dof counter size_type nbdof = 0; @@ -334,7 +342,9 @@ bgeot::pgeometric_trans pgt_old = 0; bgeot::pgeotrans_precomp pgp = 0; - for (size_type cv : cmk) { + // for (size_type cv : cmk) { + for (dal::bv_visitor cv(linked_mesh().convex_index()); + !cv.finished(); ++cv) { if (fe_convex.is_in(cv)) { gmm::copy(linked_mesh().points_of_convex(cv)[0], bmin); gmm::copy(bmin, bmax); @@ -350,8 +360,10 @@ } dal::bit_vector cv_done; - - for (size_type cv : cmk) { // Loop on elements + + //for (size_type cv : cmk) { + for (dal::bv_visitor cv(linked_mesh().convex_index()); + !cv.finished(); ++cv) { // Loop on elements if (!fe_convex.is_in(cv)) continue; pfem pf = fem_of_element(cv); if (pf != first_pf) is_uniform_ = false; Modified: trunk/getfem/src/getfem_regular_meshes.cc URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/getfem_regular_meshes.cc?rev=5417&r1=5416&r2=5417&view=diff ============================================================================== --- trunk/getfem/src/getfem_regular_meshes.cc (original) +++ trunk/getfem/src/getfem_regular_meshes.cc Fri Oct 14 09:37:49 2016 @@ -234,7 +234,7 @@ /* apply a continuous deformation + some noise */ if (noised) noise_unit_mesh(m, nsubdiv, pgt); - m.optimize_structure(); + m.optimize_structure(false); } Modified: trunk/getfem/src/gmm/gmm_vector.h URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/gmm/gmm_vector.h?rev=5417&r1=5416&r2=5417&view=diff ============================================================================== --- trunk/getfem/src/gmm/gmm_vector.h (original) +++ trunk/getfem/src/gmm/gmm_vector.h Fri Oct 14 09:37:49 2016 @@ -1058,12 +1058,12 @@ iterator it = std::lower_bound(this->begin(), this->end(), ev); if (it != this->end() && it->c == c) it->e = e; else { - size_type ind = it - this->begin(); - if (this->nb_stored() - ind > 1100) + size_type ind = it - this->begin(), nb = this->nb_stored(); + if (nb - ind > 1100) GMM_WARNING2("Inefficient addition of element in rsvector with " << this->nb_stored() - ind << " non-zero entries"); - base_type_::resize(this->nb_stored()+1, ev); - if (ind != this->nb_stored() - 1) { + base_type_::resize(nb+1, ev); + if (ind != nb) { it = this->begin() + ind; iterator ite = this->end(); --ite; iterator itee = ite; for (; ite != it; --ite) { --itee; *ite = *itee; } @@ -1085,12 +1085,12 @@ iterator it = std::lower_bound(this->begin(), this->end(), ev); if (it != this->end() && it->c == c) it->e += e; else { - size_type ind = it - this->begin(); - if (this->nb_stored() - ind > 1100) + size_type ind = it - this->begin(), nb = this->nb_stored(); + if (nb - ind > 1100) GMM_WARNING2("Inefficient addition of element in rsvector with " << this->nb_stored() - ind << " non-zero entries"); - base_type_::resize(this->nb_stored()+1, ev); - if (ind != this->nb_stored() - 1) { + base_type_::resize(nb+1, ev); + if (ind != nb) { it = this->begin() + ind; iterator ite = this->end(); --ite; iterator itee = ite; for (; ite != it; --ite) { --itee; *ite = *itee; } Modified: trunk/getfem/tests/test_assembly.cc URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/tests/test_assembly.cc?rev=5417&r1=5416&r2=5417&view=diff ============================================================================== --- trunk/getfem/tests/test_assembly.cc (original) +++ trunk/getfem/tests/test_assembly.cc Fri Oct 14 09:37:49 2016 @@ -908,6 +908,7 @@ ((std::string("GT_PK(") + Ns + ",1)").c_str()); std::vector<size_type> nsubdiv(N, NX); getfem::regular_unit_mesh(m, nsubdiv, pgt); + m.optimize_structure(); const size_type NEUMANN_BOUNDARY_NUM = 1; const size_type DIRICHLET_BOUNDARY_NUM = 2; Modified: trunk/getfem/tests/test_mesh.cc URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/tests/test_mesh.cc?rev=5417&r1=5416&r2=5417&view=diff ============================================================================== --- trunk/getfem/tests/test_mesh.cc (original) +++ trunk/getfem/tests/test_mesh.cc Fri Oct 14 09:37:49 2016 @@ -198,7 +198,7 @@ m.sup_convex(0); m.sup_convex(1); assert(!m.region(3).is_in(0)); - m.optimize_structure(); + m.optimize_structure(false); m.write_to_file("test_mesh.mesh"); getfem::mesh m2; m2.read_from_file("test_mesh.mesh"); _______________________________________________ Getfem-commits mailing list Getfem-commits@gna.org https://mail.gna.org/listinfo/getfem-commits