[Getfem-commits] (no subject)

[Liang Jin Lim]

branch: uniform_composite_fem
commit 6d4786734cc2fedc45f9b16c9db7e19952252088
Author: lj <address@hidden>
Date:   Mon Jun 4 17:02:07 2018 +0200

    Use map instead of vector to reduce the memory consumption on composite fem.
---
 src/bgeot_poly_composite.cc       | 35 ++++++++++++++++++++++++-----------
 src/getfem/bgeot_poly_composite.h | 13 ++++++++-----
 src/getfem_fem_composite.cc       |  8 ++++----
 3 files changed, 36 insertions(+), 20 deletions(-)

diff --git a/src/bgeot_poly_composite.cc b/src/bgeot_poly_composite.cc
index a5d8418..7ce7569 100644
--- a/src/bgeot_poly_composite.cc
+++ b/src/bgeot_poly_composite.cc
@@ -120,9 +120,12 @@ namespace bgeot {
             elt[ii] = false;
             p0 = pt; p0 -= mp->orgs[ii];
             gmm::mult(gmm::transposed(mp->gtrans[ii]), p0, p1);
-            if (mp->trans_of_convex(ii)->convex_ref()->is_in(p1) < 1E-10)
-              return local_coordinate ? to_scalar(polytab[ii].eval(p1.begin()))
-              : to_scalar(polytab[ii].eval(pt.begin()));
+            auto it = polytab.find(ii);
+            if (it == polytab.end()) return 0.0;
+            if (mp->trans_of_convex(ii)->convex_ref()->is_in(p1) < 1E-10) {
+              return local_coordinate ? to_scalar(it->second.eval(p1.begin()))
+                                      : to_scalar(it->second.eval(pt.begin()));
+            }
           }
         }
         ++it1; i1 = it1.index();
@@ -141,9 +144,12 @@ namespace bgeot {
             elt[ii] = false;
             p0 = pt; p0 -= mp->orgs[ii];
             gmm::mult(gmm::transposed(mp->gtrans[ii]), p0, p1);
-            if (mp->trans_of_convex(ii)->convex_ref()->is_in(p1) < 1E-10)
-              return  local_coordinate ? 
to_scalar(polytab[ii].eval(p1.begin()))
-              : to_scalar(polytab[ii].eval(pt.begin()));
+            auto it = polytab.find(ii);
+            if (it == polytab.end()) return 0.0;
+            if (mp->trans_of_convex(ii)->convex_ref()->is_in(p1) < 1E-10) {
+              return  local_coordinate ? to_scalar(it->second.eval(p1.begin()))
+                                       : 
to_scalar(it->second.eval(pt.begin()));
+            }
           }
         }
         --it2; i2 = it2.index();
@@ -153,11 +159,8 @@ namespace bgeot {
   }
 
 
-  polynomial_composite::polynomial_composite(const mesh_precomposite &m,
-    bool lc)
-    : mp(&m), polytab(m.nb_convex()), local_coordinate(lc) {
-      std::fill(polytab.begin(), polytab.end(), base_poly(m.dim(), 0));
-  }
+  polynomial_composite::polynomial_composite(
+    const mesh_precomposite &m, bool lc) : mp(&m), local_coordinate(lc) {}
 
   void polynomial_composite::derivative(short_type k) {
     if (local_coordinate) {
@@ -180,6 +183,16 @@ namespace bgeot {
         polytab[ic].derivative(k);
   }
 
+  base_poly &polynomial_composite::set_poly_of_subelt(size_type l) {
+    return polytab[l];
+  }
+
+  const base_poly *polynomial_composite::poly_of_subelt(size_type l) const {
+    auto it = polytab.find(l);
+    if (it == polytab.end()) return nullptr;
+    else return &it->second;
+  }
+
 
   /* build a regularly refined mesh for a simplex of dimension <= 3.
   All simplexes have the same orientation as the original simplex.
diff --git a/src/getfem/bgeot_poly_composite.h 
b/src/getfem/bgeot_poly_composite.h
index 5c43f65..2a83a70 100644
--- a/src/getfem/bgeot_poly_composite.h
+++ b/src/getfem/bgeot_poly_composite.h
@@ -95,7 +95,7 @@ namespace bgeot {
 
   protected :
     const mesh_precomposite *mp;
-    std::vector<base_poly> polytab;
+    std::map<size_type, base_poly> polytab;
     bool local_coordinate; // are the polynomials described on the local
     // coordinates of each sub-element or on global coordinates.
 
@@ -104,8 +104,8 @@ namespace bgeot {
     template <class ITER> scalar_type eval(const ITER &it) const;
     scalar_type eval(const base_node &pt) const;
     void derivative(short_type k);
-    base_poly &poly_of_subelt(size_type l) { return polytab[l]; }
-    const base_poly &poly_of_subelt(size_type l) const { return polytab[l]; }
+    base_poly &set_poly_of_subelt(size_type l);
+    const base_poly *poly_of_subelt(size_type l) const;
     size_type nb_subelt() const { return polytab.size(); }
 
     polynomial_composite(bool lc = true) : local_coordinate(lc) {}
@@ -116,8 +116,11 @@ namespace bgeot {
   inline std::ostream &operator <<
   (std::ostream &o, const polynomial_composite& P) {
     o << "poly_composite [";
-    for (size_type i = 0; i < P.nb_subelt(); ++i) 
-      { if (i != 0) o << ", ";  o << P.poly_of_subelt(i); }
+    for (size_type i = 0; i < P.nb_subelt(); ++i) {
+      if (i == 0) continue;
+      auto poly = P.poly_of_subelt(i);
+      o << ", ";  o << (poly) ? *poly : base_poly();
+    }
     o << "]";
     return o;
   }
diff --git a/src/getfem_fem_composite.cc b/src/getfem_fem_composite.cc
index 2f3c81e..90d9ed4 100644
--- a/src/getfem_fem_composite.cc
+++ b/src/getfem_fem_composite.cc
@@ -59,7 +59,7 @@ namespace getfem {
       for (size_type k = 0; k < pf->nb_dof(cv); ++k) {
        size_type igl = mf.ind_basic_dof_of_element(cv)[k];
        base_poly fu = pf->base()[k];
-       base[igl].poly_of_subelt(cv) = fu;
+       base[igl].set_poly_of_subelt(cv) = fu;
        dofd[igl] = pf->dof_types()[k];
       }
     }
@@ -189,7 +189,7 @@ namespace getfem {
       (4, bgeot::polynomial_composite(mp, false));
     for (size_type k = 0; k < 4; ++k)
       for (size_type ic = 0; ic < 3; ++ic)
-       base()[k].poly_of_subelt(ic) = bgeot::read_base_poly(2, s);
+        base()[k].set_poly_of_subelt(ic) = bgeot::read_base_poly(2, s);
 
     for (size_type i = 0; i < 3; ++i) {
       base_node pt(0.0, 0.0);
@@ -355,7 +355,7 @@ namespace getfem {
       (12, bgeot::polynomial_composite(mp, false));
     for (size_type k = 0; k < 12; ++k)
       for (size_type ic = 0; ic < 3; ++ic)
-       base()[k].poly_of_subelt(ic) = bgeot::read_base_poly(2, s);
+        base()[k].set_poly_of_subelt(ic) = bgeot::read_base_poly(2, s);
 
     for (size_type i = 0; i < 3; ++i) {
       base_node pt(0.0, 0.0);
@@ -621,7 +621,7 @@ namespace getfem {
       (16, bgeot::polynomial_composite(mp, false));
     for (size_type k = 0; k < 16; ++k)
       for (size_type ic = 0; ic < 4; ++ic)
-       base()[k].poly_of_subelt(ic) = bgeot::read_base_poly(2, s);
+        base()[k].set_poly_of_subelt(ic) = bgeot::read_base_poly(2, s);
 
     for (size_type i = 0; i < 4; ++i) {
       base_node pt(0.0, 0.0);