[Getfem-commits] r5324 - in /trunk/getfem: interface/tests/matlab/check_mesh_fem.m src/getfem/bgeot_kdtree.h src/getfem_mesh_fem.cc

[Yves . Renard]

Author: renard
Date: Sat May  7 14:33:14 2016
New Revision: 5324

URL: http://svn.gna.org/viewcvs/getfem?rev=5324&view=rev
Log:
optimization of enumerate_dof()

Modified:
    trunk/getfem/interface/tests/matlab/check_mesh_fem.m
    trunk/getfem/src/getfem/bgeot_kdtree.h
    trunk/getfem/src/getfem_mesh_fem.cc

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=5324&r1=5323&r2=5324&view=diff
==============================================================================
--- trunk/getfem/interface/tests/matlab/check_mesh_fem.m        (original)
+++ trunk/getfem/interface/tests/matlab/check_mesh_fem.m        Sat May  7 
14:33:14 2016
@@ -165,13 +165,13 @@
   gfassert('length(s2)>500');
   m2=gf_mesh('from string',s);
 
-  gf_mesh_fem_get(mf,'nbdof');
+  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 37 40 42 43 44 45 46 47 48 49 50 51 52 53 54 55 
56 57]']);
+  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]']);
   d=gf_mesh_fem_get(mf,'basic dof from cv',[1 5;1 2])
-  gfassert('d==[3 5 6 37 40 42 45 47 50 52 55 57]');
+  gfassert('d==[3 5 6 38 41 43 46 48 51 53 56 58]');
   d=gf_mesh_fem_get(mf,'basic dof from cvid',5)
-  gfassert('d==[37 43 44 45 46 47 40 48 49 50 51 52 42 53 54 55 56 57]');
+  gfassert('d==[38 44 45 46 47 48 41 49 50 51 52 53 43 54 55 56 57 58]');
   
   s2=gf_mesh_get(mf,'char');
   gfassert('length(s2)>500');

Modified: trunk/getfem/src/getfem/bgeot_kdtree.h
URL: 
http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/getfem/bgeot_kdtree.h?rev=5324&r1=5323&r2=5324&view=diff
==============================================================================
--- trunk/getfem/src/getfem/bgeot_kdtree.h      (original)
+++ trunk/getfem/src/getfem/bgeot_kdtree.h      Sat May  7 14:33:14 2016
@@ -107,7 +107,7 @@
   public:
     kdtree() : N(0) {}
     /// reset the tree, remove all points
-    void clear() { clear_tree(); pts.clear(); N = 0; }
+    void clear() { clear_tree(); pts = kdtree_tab_type(); N = 0; }
     void reserve(size_type n) { pts.reserve(n); }
     /// insert a new point
     size_type add_point(const base_node& n) { 

Modified: trunk/getfem/src/getfem_mesh_fem.cc
URL: 
http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/getfem_mesh_fem.cc?rev=5324&r1=5323&r2=5324&view=diff
==============================================================================
--- trunk/getfem/src/getfem_mesh_fem.cc (original)
+++ trunk/getfem/src/getfem_mesh_fem.cc Sat May  7 14:33:14 2016
@@ -305,11 +305,9 @@
 
     // Information stored per element
     size_type nb_max_cv = linked_mesh().convex_index().last_true()+1;
-    // std::vector<bgeot::node_tab> dof_nodes2(nb_max_cv, 
bgeot::node_tab(1.e5));
     std::vector<bgeot::kdtree> dof_nodes(nb_max_cv);
     std::vector<scalar_type> elt_car_sizes(nb_max_cv);
-    std::vector< std::map<fem_dof, size_type, dof_comp_> >
-      dof_sorts(nb_max_cv);
+    std::vector<std::map<fem_dof, size_type, dof_comp_>> dof_sorts(nb_max_cv);
 
     // Information for global dof
     dal::bit_vector encountered_global_dof, processed_elt;
@@ -327,18 +325,25 @@
     bgeot::pgeometric_trans pgt_old = 0;
     bgeot::pgeotrans_precomp pgp = 0;
 
+    for (size_type cv : cmk) {
+      if (fe_convex.is_in(cv)) {
+       gmm::copy(linked_mesh().points_of_convex(cv)[0], bmin);
+       gmm::copy(bmin, bmax);
+       for (size_type i = 0; i <  linked_mesh().nb_points_of_convex(cv); ++i) {
+         const base_node &pt = linked_mesh().points_of_convex(cv)[i];
+         for (size_type d = 1; d < bmin.size(); ++d) {
+           bmin[d] = std::min(bmin[d], pt[d]);
+           bmax[d] = std::max(bmax[d], pt[d]);
+         }
+       }       
+       elt_car_sizes[cv] = gmm::vect_dist2(bmin, bmax);
+      }
+    }
+
+    dal::bit_vector cv_done;
+
     for (size_type cv : cmk) { // Loop on elements
       if (!fe_convex.is_in(cv)) continue;
-      gmm::copy(linked_mesh().points_of_convex(cv)[0], bmin);
-      gmm::copy(bmin, bmax);
-      for (size_type i = 0; i <  linked_mesh().nb_points_of_convex(cv); ++i) {
-               const base_node &pt = linked_mesh().points_of_convex(cv)[i];
-               for (size_type d = 1; d < bmin.size(); ++d) {
-                 bmin[d] = std::min(bmin[d], pt[d]);
-                 bmax[d] = std::max(bmax[d], pt[d]);
-               }
-      }        
-      elt_car_sizes[cv] = gmm::vect_dist2(bmin, bmax);
       pfem pf = fem_of_element(cv);
       bgeot::pgeometric_trans pgt = linked_mesh().trans_of_convex(cv);
       bgeot::pstored_point_tab pspt = pf->node_tab(cv);
@@ -367,31 +372,40 @@
         } else {                            // For a standard linkable dof
           pgp->transform(linked_mesh().points_of_convex(cv), i, P);
           size_type idof = nbdof;
-          linked_mesh().neighbours_of_convex(cv, pf->faces_of_dof(cv, i), s);
-          for (size_type ncv : s) { // For each neighbour
-                                    // control if the dof already exists.
-           fd.ind_node = size_type(-1);
-            // size_type ii = fd.ind_node = dof_nodes2[ncv].search_node(P);
-           if (dof_nodes[ncv].nb_points() > 0) {
-              scalar_type dist = dof_nodes[ncv].nearest_neighbor(ipt, P);
-              if (gmm::abs(dist) <= 1e-6*elt_car_sizes[ncv]) fd.ind_node=ipt.i;
+
+         if (dof_nodes[cv].nb_points() > 0) {
+           scalar_type dist = dof_nodes[cv].nearest_neighbor(ipt, P);
+           if (gmm::abs(dist) <= 1e-6*elt_car_sizes[cv]) {
+             fd.ind_node=ipt.i;
+             auto it = dof_sorts[cv].find(fd);
+             if (it != dof_sorts[cv].end()) idof = it->second;
            }
-            if (fd.ind_node != size_type(-1)) {
-              auto it = dof_sorts[ncv].find(fd);
-              if (it != dof_sorts[ncv].end()) { idof = it->second; break; }
-            }
-          }
-          if (idof == nbdof) nbdof += Qdim / pf->target_dim();
+         }
+         
+          if (idof == nbdof) {
+           nbdof += Qdim / pf->target_dim();
+
+           linked_mesh().neighbours_of_convex(cv, pf->faces_of_dof(cv, i), s);
+           for (size_type ncv : s) { // For each unscanned neighbour
+             if (!cv_done[ncv] && fe_convex.is_in(ncv)) { // add the dof
+
+               fd.ind_node = size_type(-1);
+               if (dof_nodes[ncv].nb_points() > 0) {
+                 scalar_type dist = dof_nodes[ncv].nearest_neighbor(ipt, P);
+                 if (gmm::abs(dist) <= 1e-6*elt_car_sizes[ncv])
+                   fd.ind_node=ipt.i;
+               }
+               if (fd.ind_node == size_type(-1))
+                 fd.ind_node = dof_nodes[ncv].add_point(P);
+               dof_sorts[ncv][fd] = idof;
+             }
+           }
+         }
           itab[i] = idof;
-          // fd.ind_node = dof_nodes2[cv].add_node(P);
-         scalar_type dist = dof_nodes[cv].nearest_neighbor(ipt, P);
-         if (gmm::abs(dist) <= 1e-6*elt_car_sizes[cv])
-           fd.ind_node = ipt.i;
-         else
-           fd.ind_node = dof_nodes[cv].add_point(P);
-          dof_sorts[cv][fd] = idof;
         }
       }
+      cv_done.add(cv);
+      dof_sorts[cv].clear(); dof_nodes[cv].clear();
       dof_structure.add_convex_noverif(pf->structure(cv), itab.begin(), cv);
     }