getfem-commits
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

[Getfem-commits] (no subject)


From: Konstantinos Poulios
Subject: [Getfem-commits] (no subject)
Date: Thu, 27 Dec 2018 15:58:30 -0500 (EST)

branch: master
commit efdbaac130511c83be3d1015f3428e5bcc88fe59
Author: Konstantinos Poulios <address@hidden>
Date:   Thu Dec 27 21:58:20 2018 +0100

    Fix typos (also in class name) and code cleanup
---
 src/getfem/getfem_generic_assembly.h            |  12 +-
 src/getfem/getfem_generic_assembly_tree.h       |   4 +-
 src/getfem/getfem_models.h                      |  42 +-
 src/getfem_generic_assembly_compile_and_exec.cc |  84 +--
 src/getfem_generic_assembly_semantic.cc         |  13 +-
 src/getfem_generic_assembly_tree.cc             | 675 ++++++++++++------------
 src/getfem_generic_assembly_workspace.cc        |  21 +-
 src/getfem_models.cc                            | 131 ++---
 8 files changed, 482 insertions(+), 500 deletions(-)

diff --git a/src/getfem/getfem_generic_assembly.h 
b/src/getfem/getfem_generic_assembly.h
index 265ba11..5dd045c 100644
--- a/src/getfem/getfem_generic_assembly.h
+++ b/src/getfem/getfem_generic_assembly.h
@@ -183,10 +183,10 @@ namespace getfem {
   };
 
 
-  class ga_macro_dictionnary {
+  class ga_macro_dictionary {
 
   protected:
-    const ga_macro_dictionnary *parent;
+    const ga_macro_dictionary *parent;
     std::map<std::string, ga_macro> macros;
 
   public:
@@ -197,8 +197,8 @@ namespace getfem {
     void add_macro(const std::string &name, const std::string &expr);
     void del_macro(const std::string &name);
     
-    ga_macro_dictionnary() : parent(0) {}
-    ga_macro_dictionnary(bool, const ga_macro_dictionnary& gamd)
+    ga_macro_dictionary() : parent(0) {}
+    ga_macro_dictionary(bool, const ga_macro_dictionary& gamd)
       : parent(&gamd) {}
     
   };
@@ -351,7 +351,7 @@ namespace getfem {
 
     std::map<std::string, std::vector<std::string> > variable_groups;
 
-    ga_macro_dictionnary macro_dict;
+    ga_macro_dictionary macro_dict;
 
     struct m_tree {
       ga_tree *ptree;
@@ -505,7 +505,7 @@ namespace getfem {
 
     const std::string& get_macro(const std::string &name) const;
 
-    const ga_macro_dictionnary &macro_dictionnary() const { return macro_dict; 
}
+    const ga_macro_dictionary &macro_dictionary() const { return macro_dict; }
 
 
     // interpolate and elementary transformations
diff --git a/src/getfem/getfem_generic_assembly_tree.h 
b/src/getfem/getfem_generic_assembly_tree.h
index a2c9413..84b0dcb 100644
--- a/src/getfem/getfem_generic_assembly_tree.h
+++ b/src/getfem/getfem_generic_assembly_tree.h
@@ -480,9 +480,9 @@ namespace getfem {
   // Syntax analysis of an assembly string. Conversion to a tree.
   // No semantic analysis is done. The tree can be inconsistent.
   void ga_read_string(const std::string &expr, ga_tree &tree,
-                     const ga_macro_dictionnary &macro_dict);
+                      const ga_macro_dictionary &macro_dict);
   void ga_read_string_reg(const std::string &expr, ga_tree &tree,
-                         ga_macro_dictionnary &macro_dict);
+                          ga_macro_dictionary &macro_dict);
 
 
 } /* end of namespace */
diff --git a/src/getfem/getfem_models.h b/src/getfem/getfem_models.h
index 1b71bcc..9c3f01c 100644
--- a/src/getfem/getfem_models.h
+++ b/src/getfem/getfem_models.h
@@ -173,8 +173,8 @@ namespace getfem {
       std::vector<gmm::uint64_type> v_num_data;
 
       gmm::sub_interval I; // For a variable : indices on the whole system.
-      // For an affine dependent variable, should be the same than the
-      // orgininal variable
+      // For an affine dependent variable, should be the same as the
+      // original variable
       std::vector<model_real_plain_vector> real_value;
       std::vector<model_complex_plain_vector> complex_value;
       std::vector<gmm::uint64_type> v_num_var_iter;
@@ -235,7 +235,7 @@ namespace getfem {
       }
 
       size_type size() const // Should control that the variable is
-      // indeed initialized by actualize_sizes() ...
+                             // indeed initialized by actualize_sizes() ...
       { return is_complex ? complex_value[0].size() : real_value[0].size(); }
 
       void set_size();
@@ -278,8 +278,8 @@ namespace getfem {
 
   protected:
 
-    // rmatlist and cmatlist could be csc_matrix vectors to reduced the
-    // amount of memory (but this should add a supplementary copy).
+    // rmatlist and cmatlist could had been csc_matrix vectors to reduce the
+    // amount of memory (but this would require an additional copy).
     struct brick_description {
       mutable bool terms_to_be_computed;
       mutable gmm::uint64_type v_num;
@@ -292,23 +292,23 @@ namespace getfem {
       mimlist mims;              // List of integration methods.
       size_type region;          // Optional region size_type(-1) for all.
       bool is_update_brick;      // Flag for declaring special type of brick
-      // with no contributions.
+                                 // with no contributions.
       mutable scalar_type external_load; // External load computed in assembly
 
       mutable model_real_plain_vector coeffs;
       mutable scalar_type matrix_coeff = scalar_type(0);
-      mutable real_matlist rmatlist;    // Matrices the brick have to fill in
-      // (real version).
-      mutable std::vector<real_veclist> rveclist; // Rhs the brick have to
-      // fill in (real version).
-      mutable std::vector<real_veclist> rveclist_sym; // additional rhs for
-      //  symmetric terms (real version).
-      mutable complex_matlist cmatlist; // Matrices the brick have to fill in
-      // (complex version).
-      mutable std::vector<complex_veclist> cveclist; // Rhs the brick have to
-      // fill in (complex version).
-      mutable std::vector<complex_veclist> cveclist_sym;  // additional rhs
-      // for symmetric terms (real version).
+      // Matrices the brick has to fill in (real version)
+      mutable real_matlist rmatlist;
+      // Rhs the brick has to fill in (real version)
+      mutable std::vector<real_veclist> rveclist;
+      // additional rhs for symmetric terms (real version)
+      mutable std::vector<real_veclist> rveclist_sym;
+      // Matrices the brick has to fill in (complex version)
+      mutable complex_matlist cmatlist;
+      // Rhs the brick has to fill in (complex version)
+      mutable std::vector<complex_veclist> cveclist;
+      // additional rhs for symmetric terms (complex version)
+      mutable std::vector<complex_veclist> cveclist_sym;
 
       brick_description() : v_num(0) {}
 
@@ -368,17 +368,15 @@ namespace getfem {
 
     std::list<assignement_desc> assignments;
 
-    
     mutable std::list<gen_expr> generic_expressions;
 
     // Groups of variables for interpolation on different meshes
     // generic assembly
     std::map<std::string, std::vector<std::string> > variable_groups;
 
-    ga_macro_dictionnary macro_dict;
+    ga_macro_dictionary macro_dict;
     
 
-
     virtual void actualize_sizes() const;
     bool check_name_validity(const std::string &name, bool assert=true) const;
     void brick_init(size_type ib, build_version version,
@@ -842,7 +840,7 @@ namespace getfem {
                         size_type region, size_type niter = 1);
 
     /** Dictonnary of user defined macros. */
-    const ga_macro_dictionnary &macro_dictionnary() const { return macro_dict; 
}
+    const ga_macro_dictionary &macro_dictionary() const { return macro_dict; }
 
     /** Add a macro definition for the high generic assembly language.
         This macro can be used for the definition of generic assembly bricks.
diff --git a/src/getfem_generic_assembly_compile_and_exec.cc 
b/src/getfem_generic_assembly_compile_and_exec.cc
index baaa936..12b4ecc 100644
--- a/src/getfem_generic_assembly_compile_and_exec.cc
+++ b/src/getfem_generic_assembly_compile_and_exec.cc
@@ -220,6 +220,7 @@ namespace getfem {
       size_type s1 = pf1->nb_dof(cv_1) * Qmult1;
       size_type Qmult2 = qdim2 / pf2->target_dim();
       size_type s2 = pf2->nb_dof(cv_2) * Qmult2;
+      GMM_ASSERT1(s1 > 0 && s2 >0, "Element without degrees of freedom");
       if (t.sizes()[0] != s1 || t.sizes()[1] != s2) {
         bgeot::multi_index mi = t.sizes();
         mi[0] = s1; mi[1] = s2;
@@ -3378,15 +3379,16 @@ namespace getfem {
 
   // Performs Aij Bkl -> Cijkl, partially unrolled version
   template<int S1> struct ga_instruction_simple_tmult_unrolled
-    : public ga_instruction {
+  : public ga_instruction {
     base_tensor &t, &tc1, &tc2;
     virtual int exec() {
       size_type s2 = tc2.size();
+      GA_DEBUG_ASSERT(tc1.size() == S1,
+                      "Wrong sizes " << tc1.size() << " != " << S1);
       GA_DEBUG_INFO("Instruction: simple tensor product, unrolled with "
-                    << tc1.size() << " operations");
-      GA_DEBUG_ASSERT(t.size() == tc1.size() * s2, "Wrong sizes");
-      GA_DEBUG_ASSERT(tc1.size() == S1, "Wrong sizes");
-
+                    << S1 << " operations");
+      GA_DEBUG_ASSERT(t.size() == S1 * s2,
+                      "Wrong sizes " << t.size() << " != " << S1 << "*" << s2);
       base_tensor::iterator it = t.begin(), it2 = tc2.begin();
       for (size_type ii = 0; ii < s2; ++ii, ++it2) {
         base_tensor::iterator it1 = tc1.begin();
@@ -3987,7 +3989,7 @@ namespace getfem {
 
 
   struct ga_instruction_scalar_assembly : public ga_instruction {
-    base_tensor &t;
+    const base_tensor &t;
     scalar_type &E, &coeff;
      virtual int exec() {
       GA_DEBUG_INFO("Instruction: scalar term assembly");
@@ -4000,7 +4002,7 @@ namespace getfem {
   };
 
   struct ga_instruction_fem_vector_assembly : public ga_instruction {
-    base_tensor &t;
+    const base_tensor &t;
     base_vector &Vr, &Vn;
     const fem_interpolation_context &ctx;
     const gmm::sub_interval &Ir, &In;
@@ -4016,7 +4018,9 @@ namespace getfem {
         if (empty_weight) elem.resize(0);
         elem.resize(t.size());
         if (!empty_weight) {
-          auto itt = t.begin(); auto it = elem.begin(), ite = elem.end();
+          auto itt = t.begin();
+          auto it = elem.begin();
+          const auto ite = elem.end();
           size_type nd = ((t.size()) >> 2);
           for (size_type i = 0; i < nd; ++i) {
             *it++ = (*itt++) * coeff; *it++ = (*itt++) * coeff;
@@ -4025,7 +4029,9 @@ namespace getfem {
           for (; it != ite;) *it++ = (*itt++) * coeff;
         }
       } else if (!empty_weight) {
-        auto itt = t.begin(); auto it = elem.begin(), ite = elem.end();
+        auto itt = t.begin();
+        auto it = elem.begin();
+        const auto ite = elem.end();
         size_type nd = ((t.size()) >> 2);
         for (size_type i = 0; i < nd; ++i) {
           *it++ += (*itt++) * coeff; *it++ += (*itt++) * coeff;
@@ -4058,7 +4064,7 @@ namespace getfem {
       return 0;
     }
     ga_instruction_fem_vector_assembly
-    (base_tensor &t_, base_vector &Vr_, base_vector &Vn_,
+    (const base_tensor &t_, base_vector &Vr_, base_vector &Vn_,
      const fem_interpolation_context &ctx_,
      const gmm::sub_interval &Ir_, const gmm::sub_interval &In_,
      const mesh_fem *mfn_, const mesh_fem **mfg_,
@@ -4074,7 +4080,7 @@ namespace getfem {
     base_vector &V;
     const gmm::sub_interval &I;
     scalar_type &coeff;
-     virtual int exec() {
+    virtual int exec() {
       GA_DEBUG_INFO("Instruction: vector term assembly for "
                     "fixed size variable");
       gmm::add(gmm::scaled(t.as_vector(), coeff), gmm::sub_vector(V, I));
@@ -4268,7 +4274,7 @@ namespace getfem {
           } else {
             for (auto itt = ct1.begin(); itt != ct1.end(); ++itt)
               for (size_type q = 0; q < qmult1; ++q)
-                  *itd++ += *itt + q;
+                *itd++ += *itt + q;
           }
         } else
           for (size_type i=0; i < s1; ++i) dofs1[i] += i;
@@ -4712,9 +4718,9 @@ namespace getfem {
     } else {
       if (gis.var_intervals.find(varname) == gis.var_intervals.end()) {
         const mesh_fem *mf = workspace.associated_mf(varname);
-        size_type nd = mf ? mf->nb_basic_dof() :
-          gmm::vect_size(workspace.value(varname));
-        gis.var_intervals[varname]=gmm::sub_interval(gis.nb_dof, nd);
+        size_type nd = mf ? mf->nb_basic_dof()
+                          : gmm::vect_size(workspace.value(varname));
+        gis.var_intervals[varname] = gmm::sub_interval(gis.nb_dof, nd);
         gis.nb_dof += nd;
       }
       gis.max_dof = std::max(gis.max_dof,
@@ -6755,15 +6761,16 @@ namespace getfem {
                   GMM_ASSERT1(root->tensor_proper_size() == 1,
                               "Invalid vector or tensor quantity. An order 1 "
                               "weak form has to be a scalar quantity");
-                  const mesh_fem *mf=workspace.associated_mf(root->name_test1);
-                  const mesh_fem **mfg = 0;
+                  const mesh_fem *mf
+                    = workspace.associated_mf(root->name_test1);
                   add_interval_to_gis(workspace, root->name_test1, gis);
 
                   if (mf) {
-                    const std::string &intn1 = root->interpolate_name_test1;
+                    const mesh_fem **mfg = 0;
                     const gmm::sub_interval *Ir = 0, *In = 0;
+                    const std::string &intn1 = root->interpolate_name_test1;
                     bool secondary = intn1.size() &&
-                      workspace.secondary_domain_exists(intn1);
+                                     workspace.secondary_domain_exists(intn1);
                     if (intn1.size() && !secondary &&
                         workspace.variable_group_exists(root->name_test1)) {
                       ga_instruction_set::variable_group_info &vgi =
@@ -6778,12 +6785,13 @@ namespace getfem {
                       In = &(workspace.interval_of_variable(root->name_test1));
                     }
                     fem_interpolation_context &ctx
-                      = intn1.size() ?
-                      (secondary ? rmi.secondary_domain_infos.ctx
-                       : rmi.interpolate_infos[intn1].ctx) : gis.ctx;
-                    bool interpolate
-                      = (!intn1.empty() && intn1.compare("neighbour_elt")!=0 &&
-                         !secondary);
+                      = intn1.size()
+                      ? (secondary ? rmi.secondary_domain_infos.ctx
+                                   : rmi.interpolate_infos[intn1].ctx)
+                      : gis.ctx;
+                    bool interpolate = !(intn1.empty() ||
+                                         intn1 == "neighbour_elt" ||
+                                         secondary);
                     pgai = std::make_shared<ga_instruction_fem_vector_assembly>
                       (root->tensor(), workspace.unreduced_vector(),
                        workspace.assembled_vector(), ctx, *Ir, *In, mf, mfg,
@@ -6810,19 +6818,19 @@ namespace getfem {
                       workspace.secondary_domain_exists(intn1);
                   bool secondary2 = intn2.size() &&
                       workspace.secondary_domain_exists(intn2);
-                  fem_interpolation_context &ctx1
-                      = intn1.size() ?
-                      (secondary1 ? rmi.secondary_domain_infos.ctx
-                       : rmi.interpolate_infos[intn1].ctx) : gis.ctx;
-                  fem_interpolation_context &ctx2
-                      = intn2.size() ?
-                      (secondary2 ? rmi.secondary_domain_infos.ctx
-                       : rmi.interpolate_infos[intn2].ctx) : gis.ctx;
-                  bool interpolate
-                    = (!intn1.empty() && intn1.compare("neighbour_elt")!=0 &&
-                       !secondary1) ||
-                    (!intn2.empty() && intn2.compare("neighbour_elt")!=0 &&
-                     !secondary2);
+                  fem_interpolation_context
+                    &ctx1 = intn1.size()
+                          ? (secondary1 ? rmi.secondary_domain_infos.ctx
+                                        : rmi.interpolate_infos[intn1].ctx)
+                          : gis.ctx,
+                    &ctx2 = intn2.size()
+                          ? (secondary2 ? rmi.secondary_domain_infos.ctx
+                                        : rmi.interpolate_infos[intn2].ctx)
+                          : gis.ctx;
+                  bool interpolate = !(intn1.empty() || intn1 == 
"neighbour_elt"
+                                       || secondary1) ||
+                                     !(intn2.empty() || intn2 == 
"neighbour_elt"
+                                       || secondary2);
 
                   add_interval_to_gis(workspace, root->name_test1, gis);
                   add_interval_to_gis(workspace, root->name_test2, gis);
diff --git a/src/getfem_generic_assembly_semantic.cc 
b/src/getfem_generic_assembly_semantic.cc
index 956f801..f5ff172 100644
--- a/src/getfem_generic_assembly_semantic.cc
+++ b/src/getfem_generic_assembly_semantic.cc
@@ -1764,16 +1764,17 @@ namespace getfem {
             if (n == 1) {
               if (test) {
                 pnode->init_vector_tensor(1);
-                pnode->tensor()[0]=scalar_type(1);
+                pnode->tensor()[0] = scalar_type(1);
               }
-              else pnode->init_scalar_tensor(workspace.value(name)[0]);
+              else
+                pnode->init_scalar_tensor(workspace.value(name)[0]);
             } else {
               if (test) {
                 pnode->init_matrix_tensor(n,n);
                 for (size_type i = 0; i < n; ++i)
                   for (size_type j = 0; j < n; ++j)
-                    pnode->tensor()(i,j)
-                      = ((i == j) ? scalar_type(1) : scalar_type(0));
+                    pnode->tensor()(i,j) = (i == j) ? scalar_type(1)
+                                                    : scalar_type(0);
               } else {
                 pnode->t.adjust_sizes(workspace.qdims(name));
                 gmm::copy(workspace.value(name), pnode->tensor().as_vector());
@@ -1794,7 +1795,7 @@ namespace getfem {
                                    "Invalid null size of variable " << name);
             if (mii.size() > 6)
               ga_throw_error(pnode->expr, pnode->pos,
-                            "Tensor with too much dimensions. Limited to 6");
+                            "Tensor with too many dimensions. Limited to 6");
 
             switch (prefix_id) {
             case 0: // value
@@ -3891,7 +3892,7 @@ namespace getfem {
                       "transformations having an explicit expression");
 
         ga_tree trans_tree;
-        ga_read_string(expr_trans, trans_tree, workspace.macro_dictionnary());
+        ga_read_string(expr_trans, trans_tree, workspace.macro_dictionary());
         ga_semantic_analysis(trans_tree, workspace, m,
                              ref_elt_dim_of_mesh(m), false, false, 1);
         if (trans_tree.root) {
diff --git a/src/getfem_generic_assembly_tree.cc 
b/src/getfem_generic_assembly_tree.cc
index caffa4f..d3c5130 100644
--- a/src/getfem_generic_assembly_tree.cc
+++ b/src/getfem_generic_assembly_tree.cc
@@ -186,7 +186,7 @@ namespace getfem {
   //=========================================================================
 
   void ga_throw_error_msg(pstring expr, size_type pos,
-                         const std::string &msg) {
+                          const std::string &msg) {
     int length_before = 70, length_after = 70;
     if (expr && expr->size()) {
       int first = std::max(0, int(pos)-length_before);
@@ -218,12 +218,11 @@ namespace getfem {
   void ga_tree_node::mult_test(const pga_tree_node n0, const pga_tree_node n1) 
{
     
     size_type test0 = n0->test_function_type, test1 = n1->test_function_type;
-    if (test0 && test1 && (test0 == test1 ||
-                          test0 >= 3 || test1 >= 3))
+    if (test0 && test1 && (test0 == test1 || test0 >= 3 || test1 >= 3))
       ga_throw_error(expr, pos,
-                    "Incompatibility of test functions in product.");
+                     "Incompatibility of test functions in product.");
     GMM_ASSERT1(test0 != size_type(-1) && test1 != size_type(-1),
-               "internal error");
+                "internal error");
     
     test_function_type = test0 + test1;
     
@@ -288,7 +287,7 @@ namespace getfem {
   }
   
   void ga_tree::add_name(const char *name, size_type length, size_type pos,
-                        pstring expr) {
+                         pstring expr) {
     while (current_node && current_node->node_type != GA_NODE_OP)
       current_node = current_node->parent;
     if (current_node) {
@@ -304,23 +303,23 @@ namespace getfem {
   
   void ga_tree::add_sub_tree(ga_tree &sub_tree) {
     if (current_node &&
-       (current_node->node_type == GA_NODE_PARAMS ||
-        current_node->node_type == GA_NODE_INTERPOLATE_FILTER ||
-        current_node->node_type == GA_NODE_C_MATRIX)) {
+        (current_node->node_type == GA_NODE_PARAMS ||
+         current_node->node_type == GA_NODE_INTERPOLATE_FILTER ||
+         current_node->node_type == GA_NODE_C_MATRIX)) {
       GMM_ASSERT1(sub_tree.root, "Invalid tree operation");
       current_node->adopt_child(sub_tree.root);
     } else {
       GMM_ASSERT1(sub_tree.root, "Invalid tree operation");
       while (current_node && current_node->node_type != GA_NODE_OP)
-       current_node = current_node->parent;
+        current_node = current_node->parent;
       if (current_node) {
-       current_node->adopt_child(sub_tree.root);
-       current_node = sub_tree.root;
+        current_node->adopt_child(sub_tree.root);
+        current_node = sub_tree.root;
       }
       else {
-       GMM_ASSERT1(root == nullptr, "Invalid tree operation");
-       current_node = root = sub_tree.root;
-       root->parent = nullptr;
+        GMM_ASSERT1(root == nullptr, "Invalid tree operation");
+        current_node = root = sub_tree.root;
+        root->parent = nullptr;
       }
     }
     sub_tree.root = sub_tree.current_node = nullptr;
@@ -329,8 +328,8 @@ namespace getfem {
   void ga_tree::add_params(size_type pos, pstring expr) {
     GMM_ASSERT1(current_node, "internal error");
     while (current_node && current_node->parent &&
-          current_node->parent->node_type == GA_NODE_OP &&
-          ga_operator_priorities[current_node->parent->op_type] >= 4)
+           current_node->parent->node_type == GA_NODE_OP &&
+           ga_operator_priorities[current_node->parent->op_type] >= 4)
       current_node = current_node->parent;
     pga_tree_node new_node = new ga_tree_node(GA_NODE_PARAMS, pos, expr);
     new_node->parent = current_node->parent;
@@ -358,26 +357,26 @@ namespace getfem {
   }
   
   void ga_tree::add_op(GA_TOKEN_TYPE op_type, size_type pos,
-                      pstring expr) {
+                       pstring expr) {
     while (current_node && current_node->parent &&
-          current_node->parent->node_type == GA_NODE_OP &&
-          ga_operator_priorities[current_node->parent->op_type]
-          >= ga_operator_priorities[op_type])
+           current_node->parent->node_type == GA_NODE_OP &&
+           ga_operator_priorities[current_node->parent->op_type]
+           >= ga_operator_priorities[op_type])
       current_node = current_node->parent;
     pga_tree_node new_node = new ga_tree_node(op_type, pos, expr);
     if (current_node) {
       if (op_type == GA_UNARY_MINUS
-         || op_type == GA_SYM || op_type == GA_SKEW
-         || op_type == GA_TRACE || op_type == GA_DEVIATOR
-         || op_type == GA_PRINT) {
-       current_node->adopt_child(new_node);
+          || op_type == GA_SYM || op_type == GA_SKEW
+          || op_type == GA_TRACE || op_type == GA_DEVIATOR
+          || op_type == GA_PRINT) {
+        current_node->adopt_child(new_node);
       } else {
-       new_node->parent = current_node->parent;
-       if (current_node->parent)
-         current_node->parent->replace_child(current_node, new_node);
-       else
-         root = new_node;
-       new_node->adopt_child(current_node);
+        new_node->parent = current_node->parent;
+        if (current_node->parent)
+          current_node->parent->replace_child(current_node, new_node);
+        else
+          root = new_node;
+        new_node->adopt_child(current_node);
       }
     } else {
       if (root) new_node->adopt_child(root);
@@ -390,7 +389,7 @@ namespace getfem {
   void ga_tree::clear_node_rec(pga_tree_node pnode) {
     if (pnode) {
       for (pga_tree_node &child : pnode->children)
-       clear_node_rec(child);
+        clear_node_rec(child);
       delete pnode;
       current_node = nullptr;
     }
@@ -400,13 +399,13 @@ namespace getfem {
     if (pnode) {
       pga_tree_node parent = pnode->parent;
       if (parent) { // keep all siblings of pnode
-       size_type j = 0;
-       for (pga_tree_node &sibling : parent->children)
-         if (sibling != pnode)
-           parent->children[j++] = sibling;
-       parent->children.resize(j);
+        size_type j = 0;
+        for (pga_tree_node &sibling : parent->children)
+          if (sibling != pnode)
+            parent->children[j++] = sibling;
+        parent->children.resize(j);
       } else
-       root = nullptr;
+        root = nullptr;
     }
     clear_node_rec(pnode);
   }
@@ -432,7 +431,7 @@ namespace getfem {
   }
   
   void ga_tree::copy_node(pga_tree_node pnode, pga_tree_node parent,
-                         pga_tree_node &child) {
+                          pga_tree_node &child) {
     GMM_ASSERT1(child == nullptr, "Internal error");
     child = new ga_tree_node();
     *child = *pnode;
@@ -444,7 +443,7 @@ namespace getfem {
   }
   
   void ga_tree::duplicate_with_operation(pga_tree_node pnode,
-                                        GA_TOKEN_TYPE op_type) {
+                                         GA_TOKEN_TYPE op_type) {
     pga_tree_node newop = new ga_tree_node(op_type, pnode->pos, pnode->expr);
     newop->children.resize(2, nullptr);
     newop->children[0] = pnode;
@@ -1091,65 +1090,65 @@ namespace getfem {
       GMM_ASSERT1(pnode->nbc1 == pnode->tensor_order(), "Invalid C_MATRIX");
       switch (pnode->tensor_order()) {
       case 0:
-       ga_print_node(pnode->children[0], str);
-       break;
-       
+        ga_print_node(pnode->children[0], str);
+        break;
+
       case 1:
-       str << "[";
-       for (size_type i = 0; i < pnode->tensor_proper_size(0); ++i) {
-         if (i != 0) str << ",";
-         ga_print_node(pnode->children[i], str);
-       }
-       str << "]";
-       break;
-       
+        str << "[";
+        for (size_type i = 0; i < pnode->tensor_proper_size(0); ++i) {
+          if (i != 0) str << ",";
+          ga_print_node(pnode->children[i], str);
+        }
+        str << "]";
+        break;
+
       case 2: case 3: case 4:
-       {
-         size_type ii(0);
-         size_type n0 = pnode->tensor_proper_size(0);
-         size_type n1 = pnode->tensor_proper_size(1);
-         size_type n2 = ((pnode->tensor_order() > 2) ?
-                         pnode->tensor_proper_size(2) : 1);
-         size_type n3 = ((pnode->tensor_order() > 3) ?
-                         pnode->tensor_proper_size(3) : 1);
-         if (n3 > 1) str << "[";
-         for (size_type l = 0; l < n3; ++l) {
-           if (l != 0) str << ",";
-           if (n2 > 1) str << "[";
-           for (size_type k = 0; k < n2; ++k) {
-             if (k != 0) str << ",";
-             if (n1 > 1) str << "[";
-             for (size_type j = 0; j < n1; ++j) {
-               if (j != 0) str << ",";
-               if (n0 > 1) str << "[";
-               for (size_type i = 0; i < n0; ++i) {
-                 if (i != 0) str << ",";
-                 ga_print_node(pnode->children[ii++], str);
-               }
-               if (n0 > 1) str << "]";
-             }
-             if (n1 > 1) str << "]";
-           }
-           if (n2 > 1) str << "]";
-         }
-         if (n3 > 1) str << "]";
-       }
-       break;
-       
+        {
+          size_type ii(0);
+          size_type n0 = pnode->tensor_proper_size(0);
+          size_type n1 = pnode->tensor_proper_size(1);
+          size_type n2 = ((pnode->tensor_order() > 2) ?
+                          pnode->tensor_proper_size(2) : 1);
+          size_type n3 = ((pnode->tensor_order() > 3) ?
+                          pnode->tensor_proper_size(3) : 1);
+          if (n3 > 1) str << "[";
+          for (size_type l = 0; l < n3; ++l) {
+            if (l != 0) str << ",";
+            if (n2 > 1) str << "[";
+            for (size_type k = 0; k < n2; ++k) {
+              if (k != 0) str << ",";
+              if (n1 > 1) str << "[";
+              for (size_type j = 0; j < n1; ++j) {
+                if (j != 0) str << ",";
+                if (n0 > 1) str << "[";
+                for (size_type i = 0; i < n0; ++i) {
+                  if (i != 0) str << ",";
+                  ga_print_node(pnode->children[ii++], str);
+                }
+                if (n0 > 1) str << "]";
+              }
+              if (n1 > 1) str << "]";
+            }
+            if (n2 > 1) str << "]";
+          }
+          if (n3 > 1) str << "]";
+        }
+        break;
+
       case 5: case 6: case 7: case 8:
-       str << "Reshape([";
-       for (size_type i = 0; i < pnode->tensor_proper_size(); ++i) {
-         if (i != 0) str << ";";
-         ga_print_node(pnode->children[i], str);
-       }
-       str << "]";
-       for (size_type i = 0; i < pnode->tensor_order(); ++i) {
-         if (i != 0) str << ", ";
-         str << pnode->tensor_proper_size(i);
-       }
-       str << ")";
-       break;
-       
+        str << "Reshape([";
+        for (size_type i = 0; i < pnode->tensor_proper_size(); ++i) {
+          if (i != 0) str << ";";
+          ga_print_node(pnode->children[i], str);
+        }
+        str << "]";
+        for (size_type i = 0; i < pnode->tensor_order(); ++i) {
+          if (i != 0) str << ", ";
+          str << pnode->tensor_proper_size(i);
+        }
+        str << ")";
+        break;
+
       default: GMM_ASSERT1(false, "Invalid tensor dimension");
       }
       break;
@@ -1291,156 +1290,156 @@ namespace getfem {
       pga_tree_node pchild = children[pnode->nbc1+1];
 
       if (po || pt) {
-       if (!(pchild->children.empty()) || pchild->node_type != GA_NODE_NAME)
-         ga_throw_error(pchild->expr, pchild->pos, "Error in macro "
-                        "expansion. Only variable name are allowed for macro "
-                        "parameter preceded by Grad_ Hess_ Test_ or Test2_ "
-                        "prefixes.");
-       switch(pnode->op_type) {
-       case GA_NAME : pnode->node_type = GA_NODE_NAME; break;
-       case GA_INTERPOLATE : pnode->node_type = GA_NODE_INTERPOLATE; break;
-       case GA_ELEMENTARY : pnode->node_type = GA_NODE_ELEMENTARY; break;
-       case GA_SECONDARY_DOMAIN :
+        if (!(pchild->children.empty()) || pchild->node_type != GA_NODE_NAME)
+          ga_throw_error(pchild->expr, pchild->pos, "Error in macro "
+                         "expansion. Only variable name are allowed for macro "
+                         "parameter preceded by Grad_ Hess_ Test_ or Test2_ "
+                         "prefixes.");
+        switch(pnode->op_type) {
+        case GA_NAME : pnode->node_type = GA_NODE_NAME; break;
+        case GA_INTERPOLATE : pnode->node_type = GA_NODE_INTERPOLATE; break;
+        case GA_ELEMENTARY : pnode->node_type = GA_NODE_ELEMENTARY; break;
+        case GA_SECONDARY_DOMAIN :
           pnode->node_type = GA_NODE_SECONDARY_DOMAIN; break;
-       case GA_XFEM_PLUS : pnode->node_type = GA_NODE_XFEM_PLUS; break;
-       case GA_XFEM_MINUS: pnode->node_type = GA_NODE_XFEM_MINUS; break;
-       default:break;
-       }
-       pnode->name = pchild->name;
-       if (pt == 1) pnode->name = "Test_" + pnode->name;
-       if (pt == 2) pnode->name = "Test2_" + pnode->name;
-       if (po == 1) pnode->name = "Grad_" + pnode->name;
-       if (po == 2) pnode->name = "Hess_" + pnode->name;
-       if (po == 3) pnode->name = "Div_" + pnode->name;
+        case GA_XFEM_PLUS : pnode->node_type = GA_NODE_XFEM_PLUS; break;
+        case GA_XFEM_MINUS: pnode->node_type = GA_NODE_XFEM_MINUS; break;
+        default:break;
+        }
+        pnode->name = pchild->name;
+        if (pt == 1) pnode->name = "Test_" + pnode->name;
+        if (pt == 2) pnode->name = "Test2_" + pnode->name;
+        if (po == 1) pnode->name = "Grad_" + pnode->name;
+        if (po == 2) pnode->name = "Hess_" + pnode->name;
+        if (po == 3) pnode->name = "Div_" + pnode->name;
       } else {   
-       pga_tree_node pnode_old = pnode;
-       pnode = nullptr;
-       tree.copy_node(pchild, pnode_old->parent, pnode);
-       if (pnode_old->parent)
-         pnode_old->parent->replace_child(pnode_old, pnode);
-       else
-         tree.root = pnode;
-       GMM_ASSERT1(pnode_old->children.empty(), "Internal error");
-       delete pnode_old;
+        pga_tree_node pnode_old = pnode;
+        pnode = nullptr;
+        tree.copy_node(pchild, pnode_old->parent, pnode);
+        if (pnode_old->parent)
+          pnode_old->parent->replace_child(pnode_old, pnode);
+        else
+          tree.root = pnode;
+        GMM_ASSERT1(pnode_old->children.empty(), "Internal error");
+        delete pnode_old;
       }
     }
   }
   
   static void ga_expand_macro(ga_tree &tree, pga_tree_node pnode,
-                             const ga_macro_dictionnary &macro_dict) {
+                              const ga_macro_dictionary &macro_dict) {
     if (!pnode) return;
     
     if (pnode->node_type == GA_NODE_PARAMS) {
       
       for (size_type i = 1; i < pnode->children.size(); ++i)
-       ga_expand_macro(tree, pnode->children[i], macro_dict);
+        ga_expand_macro(tree, pnode->children[i], macro_dict);
 
       if (pnode->children[0]->node_type != GA_NODE_NAME) {
-       ga_expand_macro(tree, pnode->children[0], macro_dict);
+        ga_expand_macro(tree, pnode->children[0], macro_dict);
       } else {
-       
-       if (macro_dict.macro_exists(pnode->children[0]->name)) {
-         
-         const ga_macro &gam = macro_dict.get_macro(pnode->children[0]->name);
-         
-         if (gam.nb_params()==0) { // Macro without parameters
-           pga_tree_node pnode_old = pnode->children[0];
-           pnode->children[0] = nullptr;
-           tree.copy_node(gam.tree().root,
-                          pnode_old->parent,pnode->children[0]);
-           GMM_ASSERT1(pnode_old->children.empty(), "Internal error");
-           delete pnode_old;
-           
-         } else { // Macro with parameters
-           
-           if (gam.nb_params()+1 != pnode->children.size())
-             ga_throw_error(pnode->expr, pnode->pos,
-                            "Bad number of parameters in the use of macro '"
-                            << gam.name() << "'. Expected " << gam.nb_params()
-                            << " found " << pnode->children.size()-1 << ".");
-           
-           pga_tree_node pnode_old = pnode;
-           pnode = nullptr;
-           tree.copy_node(gam.tree().root, pnode_old->parent, pnode);
-           if (pnode_old->parent)
-             pnode_old->parent->replace_child(pnode_old, pnode);
-           else
-             tree.root = pnode;
-           ga_replace_macro_params(tree, pnode, pnode_old->children);
-         }
-       }
+
+        if (macro_dict.macro_exists(pnode->children[0]->name)) {
+
+          const ga_macro &gam = macro_dict.get_macro(pnode->children[0]->name);
+
+          if (gam.nb_params()==0) { // Macro without parameters
+            pga_tree_node pnode_old = pnode->children[0];
+            pnode->children[0] = nullptr;
+            tree.copy_node(gam.tree().root,
+                           pnode_old->parent,pnode->children[0]);
+            GMM_ASSERT1(pnode_old->children.empty(), "Internal error");
+            delete pnode_old;
+
+          } else { // Macro with parameters
+
+            if (gam.nb_params()+1 != pnode->children.size())
+              ga_throw_error(pnode->expr, pnode->pos,
+                             "Bad number of parameters in the use of macro '"
+                             << gam.name() << "'. Expected " << gam.nb_params()
+                             << " found " << pnode->children.size()-1 << ".");
+
+            pga_tree_node pnode_old = pnode;
+            pnode = nullptr;
+            tree.copy_node(gam.tree().root, pnode_old->parent, pnode);
+            if (pnode_old->parent)
+              pnode_old->parent->replace_child(pnode_old, pnode);
+            else
+              tree.root = pnode;
+            ga_replace_macro_params(tree, pnode, pnode_old->children);
+          }
+        }
       }
 
     } else if (pnode->node_type == GA_NODE_NAME &&
-              macro_dict.macro_exists(pnode->name)) {
+               macro_dict.macro_exists(pnode->name)) {
       // Macro without parameters
       const ga_macro &gam = macro_dict.get_macro(pnode->name);
       if (gam.nb_params() != 0)
-       ga_throw_error(pnode->expr, pnode->pos,
-                        "Bad number of parameters in the use of macro '"
-                        << gam.name() << "'. Expected " << gam.nb_params()
-                        << " none found.");
+        ga_throw_error(pnode->expr, pnode->pos,
+                       "Bad number of parameters in the use of macro '"
+                       << gam.name() << "'. Expected " << gam.nb_params()
+                       << " none found.");
 
       pga_tree_node pnode_old = pnode;
       pnode = nullptr;
       tree.copy_node(gam.tree().root, pnode_old->parent, pnode);
       if (pnode_old->parent)
-       pnode_old->parent->replace_child(pnode_old, pnode);
+        pnode_old->parent->replace_child(pnode_old, pnode);
       else
-       tree.root = pnode;
+        tree.root = pnode;
       GMM_ASSERT1(pnode_old->children.empty(), "Internal error");
       delete pnode_old;
     } else {
       for (size_type i = 0; i < pnode->children.size(); ++i)
-       ga_expand_macro(tree, pnode->children[i], macro_dict);
+        ga_expand_macro(tree, pnode->children[i], macro_dict);
     }
   }
 
   static void ga_mark_macro_params_rec(const pga_tree_node pnode,
-                                      const std::vector<std::string> &params) {
+                                       const std::vector<std::string> &params) 
{
     if (!pnode) return;
     for (size_type i = 0; i < pnode->children.size(); ++i)
       ga_mark_macro_params_rec(pnode->children[i], params);
     
     if (pnode->node_type == GA_NODE_NAME ||
-       pnode->node_type == GA_NODE_INTERPOLATE ||
-       pnode->node_type == GA_NODE_ELEMENTARY ||
-       pnode->node_type == GA_NODE_SECONDARY_DOMAIN ||
-       pnode->node_type == GA_NODE_XFEM_PLUS ||
-       pnode->node_type == GA_NODE_XFEM_MINUS) {
+        pnode->node_type == GA_NODE_INTERPOLATE ||
+        pnode->node_type == GA_NODE_ELEMENTARY ||
+        pnode->node_type == GA_NODE_SECONDARY_DOMAIN ||
+        pnode->node_type == GA_NODE_XFEM_PLUS ||
+        pnode->node_type == GA_NODE_XFEM_MINUS) {
       std::string name = pnode->name;
       size_type po = ga_parse_prefix_operator(name);
       size_type pt = ga_parse_prefix_test(name);
 
       for (size_type i = 0; i < params.size(); ++i)
-       if (name.compare(params[i]) == 0) {
-         pnode->name = name;
-         switch(pnode->node_type) {
-         case GA_NODE_NAME : pnode->op_type = GA_NAME; break;
-         case GA_NODE_INTERPOLATE : pnode->op_type = GA_INTERPOLATE; break;
-         case GA_NODE_ELEMENTARY : pnode->op_type = GA_ELEMENTARY; break;
-         case GA_NODE_SECONDARY_DOMAIN :
+        if (name.compare(params[i]) == 0) {
+          pnode->name = name;
+          switch(pnode->node_type) {
+          case GA_NODE_NAME : pnode->op_type = GA_NAME; break;
+          case GA_NODE_INTERPOLATE : pnode->op_type = GA_INTERPOLATE; break;
+          case GA_NODE_ELEMENTARY : pnode->op_type = GA_ELEMENTARY; break;
+          case GA_NODE_SECONDARY_DOMAIN :
             pnode->op_type = GA_SECONDARY_DOMAIN; break;
-         case GA_NODE_XFEM_PLUS : pnode->op_type = GA_XFEM_PLUS; break;
-         case GA_NODE_XFEM_MINUS: pnode->op_type = GA_XFEM_MINUS; break;
-         default:break;
-         }
-         pnode->node_type = GA_NODE_MACRO_PARAM;
-         pnode->nbc1 = i; pnode->nbc2 = po; pnode->nbc3 = pt;
-       }
+          case GA_NODE_XFEM_PLUS : pnode->op_type = GA_XFEM_PLUS; break;
+          case GA_NODE_XFEM_MINUS: pnode->op_type = GA_XFEM_MINUS; break;
+          default:break;
+          }
+          pnode->node_type = GA_NODE_MACRO_PARAM;
+          pnode->nbc1 = i; pnode->nbc2 = po; pnode->nbc3 = pt;
+        }
     }
   }
   
   static void ga_mark_macro_params(ga_macro &gam,
-                                  const std::vector<std::string> &params,
-                                  const ga_macro_dictionnary &macro_dict) {
+                                   const std::vector<std::string> &params,
+                                   const ga_macro_dictionary &macro_dict) {
     if (gam.tree().root) {
       ga_mark_macro_params_rec(gam.tree().root, params);
       ga_expand_macro(gam.tree(), gam.tree().root, macro_dict);
     }
   }
 
-  bool ga_macro_dictionnary::macro_exists(const std::string &name) const {
+  bool ga_macro_dictionary::macro_exists(const std::string &name) const {
     if (macros.find(name) != macros.end()) return true;
     if (parent && parent->macro_exists(name)) return true;
     return false;
@@ -1448,21 +1447,21 @@ namespace getfem {
 
 
   const ga_macro &
-  ga_macro_dictionnary::get_macro(const std::string &name) const {
+  ga_macro_dictionary::get_macro(const std::string &name) const {
     auto it = macros.find(name);
     if (it != macros.end()) return it->second;
     if (parent) return parent->get_macro(name);
     GMM_ASSERT1(false, "Undefined macro");
   }
 
-  void ga_macro_dictionnary::add_macro(const ga_macro &gam)
+  void ga_macro_dictionary::add_macro(const ga_macro &gam)
   { macros[gam.name()] = gam; }
 
-  void ga_macro_dictionnary::add_macro(const std::string &name,
-                                      const std::string &expr)
+  void ga_macro_dictionary::add_macro(const std::string &name,
+                                      const std::string &expr)
   { ga_tree tree; ga_read_string_reg("Def "+name+":="+expr, tree, *this); }
 
-  void ga_macro_dictionnary::del_macro(const std::string &name) {
+  void ga_macro_dictionary::del_macro(const std::string &name) {
     auto it = macros.find(name);
     GMM_ASSERT1(it != macros.end(), "Undefined macro (at this level)");
     macros.erase(it);
@@ -1476,7 +1475,7 @@ namespace getfem {
   // Read a term with an (implicit) pushdown automaton.
   static GA_TOKEN_TYPE ga_read_term(pstring expr, size_type &pos,
                                     ga_tree &tree,
-                                   ga_macro_dictionnary &macro_dict) {
+                                    ga_macro_dictionary &macro_dict) {
     size_type token_pos, token_length;
     GA_TOKEN_TYPE t_type;
     int state = 1; // 1 = reading term, 2 = reading after term
@@ -1528,60 +1527,60 @@ namespace getfem {
           tree.add_op(GA_DEVIATOR, token_pos, expr);
           state = 1; break;
 
-       case GA_DEF:
-         {
-           ga_macro gam;
-           t_type = ga_get_token(*expr, pos, token_pos, token_length);
-           if (t_type != GA_NAME)
+        case GA_DEF:
+          {
+            ga_macro gam;
+            t_type = ga_get_token(*expr, pos, token_pos, token_length);
+            if (t_type != GA_NAME)
               ga_throw_error(expr, pos,
                              "Macro definition should begin with macro name");
-           gam.name() = std::string(&((*expr)[token_pos]), token_length);
-           if (ga_check_name_validity(gam.name()))
-             ga_throw_error(expr, pos-1, "Invalid macro name.")
-           t_type = ga_get_token(*expr, pos, token_pos, token_length);
-           std::vector<std::string> params;
+            gam.name() = std::string(&((*expr)[token_pos]), token_length);
+            if (ga_check_name_validity(gam.name()))
+              ga_throw_error(expr, pos-1, "Invalid macro name.")
+            t_type = ga_get_token(*expr, pos, token_pos, token_length);
+            std::vector<std::string> params;
             if (t_type == GA_LPAR) {
-             t_type = ga_get_token(*expr, pos, token_pos, token_length);
-             while (t_type == GA_NAME) {
-               params.push_back(std::string(&((*expr)[token_pos]),
-                                            token_length));
-               if (ga_check_name_validity(params.back()))
-                 ga_throw_error(expr, pos-1, "Invalid macro parameter name.");
-               for (size_type i = 0; i+1 < params.size(); ++i)
-                 if (params.back().compare(params[i]) == 0)
-                   ga_throw_error(expr, pos-1,
-                                  "Invalid repeated macro parameter name.");
-               t_type = ga_get_token(*expr, pos, token_pos, token_length);
-               if (t_type == GA_COMMA)
-                 t_type = ga_get_token(*expr, pos, token_pos, token_length);
-             }
-             if (t_type != GA_RPAR)
-               ga_throw_error(expr, pos-1,
-                             "Missing right parenthesis in macro definition.");
-             t_type = ga_get_token(*expr, pos, token_pos, token_length);
-           }
-           if (t_type != GA_COLON_EQ)
-             ga_throw_error(expr, pos-1, "Missing := for macro definition.");
-
-           t_type = ga_read_term(expr, pos, gam.tree(), macro_dict);
-           if (gam.tree().root)
-             ga_expand_macro(gam.tree(), gam.tree().root, macro_dict);
-           gam.nb_params() = params.size();
-           if (params.size())
-             ga_mark_macro_params(gam, params, macro_dict);
-           macro_dict.add_macro(gam);
-
-           // cout << "macro \"" << gam.name() << "\" registered with "
-           //   << gam.nb_params() << " params  := "
-           //   << ga_tree_to_string(gam.tree()) << endl;
-           
-           if (t_type == GA_END) return t_type;
+              t_type = ga_get_token(*expr, pos, token_pos, token_length);
+              while (t_type == GA_NAME) {
+                params.push_back(std::string(&((*expr)[token_pos]),
+                                             token_length));
+                if (ga_check_name_validity(params.back()))
+                  ga_throw_error(expr, pos-1, "Invalid macro parameter name.");
+                for (size_type i = 0; i+1 < params.size(); ++i)
+                  if (params.back().compare(params[i]) == 0)
+                    ga_throw_error(expr, pos-1,
+                                   "Invalid repeated macro parameter name.");
+                t_type = ga_get_token(*expr, pos, token_pos, token_length);
+                if (t_type == GA_COMMA)
+                  t_type = ga_get_token(*expr, pos, token_pos, token_length);
+              }
+              if (t_type != GA_RPAR)
+                ga_throw_error(expr, pos-1,
+                              "Missing right parenthesis in macro 
definition.");
+              t_type = ga_get_token(*expr, pos, token_pos, token_length);
+            }
+            if (t_type != GA_COLON_EQ)
+              ga_throw_error(expr, pos-1, "Missing := for macro definition.");
+
+            t_type = ga_read_term(expr, pos, gam.tree(), macro_dict);
+            if (gam.tree().root)
+              ga_expand_macro(gam.tree(), gam.tree().root, macro_dict);
+            gam.nb_params() = params.size();
+            if (params.size())
+              ga_mark_macro_params(gam, params, macro_dict);
+            macro_dict.add_macro(gam);
+
+            // cout << "macro \"" << gam.name() << "\" registered with "
+            //      << gam.nb_params() << " params  := "
+            //      << ga_tree_to_string(gam.tree()) << endl;
+
+            if (t_type == GA_END) return t_type;
             else if (t_type != GA_SEMICOLON)
               ga_throw_error(expr, pos-1,
-                            "Syntax error at the end of macro definition.");
-           state = 1;
-         }
-         break;
+                             "Syntax error at the end of macro definition.");
+            state = 1;
+          }
+          break;
 
         case GA_INTERPOLATE:
           {
@@ -1780,66 +1779,66 @@ namespace getfem {
             size_type tensor_order(1);
             bool foundcomma(false), foundsemi(false);
 
-           r_type = ga_read_term(expr, pos, sub_tree, macro_dict);
-           size_type nb_comp = 0;
-           tree.add_matrix(token_pos, expr);
-           
-           if (sub_tree.root->node_type == GA_NODE_C_MATRIX) { // nested format
-             bgeot::multi_index mii;
-             do {
-               if (nb_comp) {
-                 sub_tree.clear();
-                 r_type = ga_read_term(expr, pos, sub_tree, macro_dict);
-               }
-               // in the nested format only "," and "]" are expected
-               if (sub_tree.root->node_type != GA_NODE_C_MATRIX || 
-                   (r_type != GA_COMMA && r_type != GA_RBRACKET))
+            r_type = ga_read_term(expr, pos, sub_tree, macro_dict);
+            size_type nb_comp = 0;
+            tree.add_matrix(token_pos, expr);
+
+            if (sub_tree.root->node_type == GA_NODE_C_MATRIX) { // nested 
format
+              bgeot::multi_index mii;
+              do {
+                if (nb_comp) {
+                  sub_tree.clear();
+                  r_type = ga_read_term(expr, pos, sub_tree, macro_dict);
+                }
+                // in the nested format only "," and "]" are expected
+                if (sub_tree.root->node_type != GA_NODE_C_MATRIX ||
+                    (r_type != GA_COMMA && r_type != GA_RBRACKET))
                   ga_throw_error(expr, pos-1, "Bad explicit "
                                  "vector/matrix/tensor format.");
 
-               // convert a row vector [a,b] to a column vector [a;b]
+                // convert a row vector [a,b] to a column vector [a;b]
                 if (sub_tree.root->marked &&
-                   sub_tree.root->tensor().sizes()[0] == 1 &&
-                   sub_tree.root->tensor().size() != 1) {
-                 bgeot::multi_index mi = sub_tree.root->tensor().sizes();
-                 for (size_type i = mi.size()-1; i > 0; i--)
-                   mi[i-1] = mi[i];
-                 mi.pop_back();
-                 sub_tree.root->tensor().adjust_sizes(mi);
-               }
-               if (!nb_comp) mii = sub_tree.root->tensor().sizes();
-               else {
-                 const bgeot::multi_index &mi=sub_tree.root->tensor().sizes();
-                 bool cmp = true;
-                 if (mii.size() == mi.size()) {
-                    for (size_type i = 0; i < mi.size(); ++i)
-                      if (mi[i] != mii[i]) cmp = false;
-                 } else cmp = false;
-                 if (!cmp)
-                   ga_throw_error(expr, pos-1, "Bad explicit "
+                    sub_tree.root->tensor().sizes()[0] == 1 &&
+                    sub_tree.root->tensor().size() != 1) {
+                  bgeot::multi_index mi = sub_tree.root->tensor().sizes();
+                  for (size_type i = mi.size()-1; i > 0; i--)
+                    mi[i-1] = mi[i];
+                  mi.pop_back();
+                  sub_tree.root->tensor().adjust_sizes(mi);
+                }
+                if (!nb_comp) mii = sub_tree.root->tensor().sizes();
+                else {
+                  const bgeot::multi_index &mi=sub_tree.root->tensor().sizes();
+                  bool cmp = true;
+                  if (mii.size() == mi.size()) {
+                     for (size_type i = 0; i < mi.size(); ++i)
+                       if (mi[i] != mii[i]) cmp = false;
+                  } else cmp = false;
+                  if (!cmp)
+                    ga_throw_error(expr, pos-1, "Bad explicit "
                                    "vector/matrix/tensor format.");
-               }
-               for (size_type i = 0; i < sub_tree.root->children.size(); ++i) {
-                 sub_tree.root->children[i]->parent = tree.current_node;
-                 tree.current_node->children.push_back
-                   (sub_tree.root->children[i]);
-               }
-               sub_tree.root->children.resize(0);
-               nb_comp++;
-             } while (r_type != GA_RBRACKET);
-             tree.current_node->marked = false;
-             mii.push_back(nb_comp);
-             tree.current_node->tensor().adjust_sizes(mii);
-           } else { // non nested format
-             do {
-               if (nb_comp) {
-                 sub_tree.clear();
-                 r_type = ga_read_term(expr, pos, sub_tree, macro_dict);
-               }
-               nb_comp++;
-               
-               tree.add_sub_tree(sub_tree);
-               
+                }
+                for (size_type i = 0; i < sub_tree.root->children.size(); ++i) 
{
+                  sub_tree.root->children[i]->parent = tree.current_node;
+                  tree.current_node->children.push_back
+                    (sub_tree.root->children[i]);
+                }
+                sub_tree.root->children.resize(0);
+                nb_comp++;
+              } while (r_type != GA_RBRACKET);
+              tree.current_node->marked = false;
+              mii.push_back(nb_comp);
+              tree.current_node->tensor().adjust_sizes(mii);
+            } else { // non nested format
+              do {
+                if (nb_comp) {
+                  sub_tree.clear();
+                  r_type = ga_read_term(expr, pos, sub_tree, macro_dict);
+                }
+                nb_comp++;
+
+                tree.add_sub_tree(sub_tree);
+
                 ++n1; ++n2; ++n3;
                 if (tensor_order < 2) ++nbc1;
                 if (tensor_order < 3) ++nbc2;
@@ -1888,42 +1887,42 @@ namespace getfem {
                                  "separated by ',', ';', ',,' and ';;' and "
                                  "be ended by ']'.");
                 }
-               
-             } while (r_type != GA_RBRACKET);
-             bgeot::multi_index mi;
-             nbc1 = tree.current_node->nbc1; nbc2 = tree.current_node->nbc2;
-             nbc3 = tree.current_node->nbc3;
-             
-             size_type nbl = tree.current_node->children.size()
-               / (nbc2 * nbc1 * nbc3);
-             switch(tensor_order) {
-             case 1:
-               /* mi.push_back(1); */ mi.push_back(nbc1); break;
-             case 2:
-               mi.push_back(nbl); if (nbc1 > 1) mi.push_back(nbc1); break; 
-             case 3:
-               mi.push_back(nbl); mi.push_back(nbc2);
-               mi.push_back(nbc1);
-               break;
-             case 4:
-               mi.push_back(nbl); mi.push_back(nbc3);
-               mi.push_back(nbc2); mi.push_back(nbc1);
-               break;
-             default: GMM_ASSERT1(false, "Internal error");
-             }
-             tree.current_node->tensor().adjust_sizes(mi);
-             std::vector<pga_tree_node> children = tree.current_node->children;
-             auto it = tree.current_node->children.begin();
-             for (size_type i = 0; i < nbc1; ++i)
-               for (size_type j = 0; j < nbc2; ++j)
-                 for (size_type k = 0; k < nbc3; ++k)
-                   for (size_type l = 0; l < nbl; ++l, ++it)
-                     *it = children[i+nbc1*(j+nbc2*(k+nbc3*l))];
-             tree.current_node->marked = true;
-           }
+
+              } while (r_type != GA_RBRACKET);
+              bgeot::multi_index mi;
+              nbc1 = tree.current_node->nbc1; nbc2 = tree.current_node->nbc2;
+              nbc3 = tree.current_node->nbc3;
+
+              size_type nbl = tree.current_node->children.size()
+                / (nbc2 * nbc1 * nbc3);
+              switch(tensor_order) {
+              case 1:
+                /* mi.push_back(1); */ mi.push_back(nbc1); break;
+              case 2:
+                mi.push_back(nbl); if (nbc1 > 1) mi.push_back(nbc1); break;
+              case 3:
+                mi.push_back(nbl); mi.push_back(nbc2);
+                mi.push_back(nbc1);
+                break;
+              case 4:
+                mi.push_back(nbl); mi.push_back(nbc3);
+                mi.push_back(nbc2); mi.push_back(nbc1);
+                break;
+              default: GMM_ASSERT1(false, "Internal error");
+              }
+              tree.current_node->tensor().adjust_sizes(mi);
+              std::vector<pga_tree_node> children = 
tree.current_node->children;
+              auto it = tree.current_node->children.begin();
+              for (size_type i = 0; i < nbc1; ++i)
+                for (size_type j = 0; j < nbc2; ++j)
+                  for (size_type k = 0; k < nbc3; ++k)
+                    for (size_type l = 0; l < nbl; ++l, ++it)
+                      *it = children[i+nbc1*(j+nbc2*(k+nbc3*l))];
+              tree.current_node->marked = true;
+            }
           }
-         tree.current_node->nbc1 = tree.current_node->tensor().sizes().size();
-         state = 2;
+          tree.current_node->nbc1 = tree.current_node->tensor().sizes().size();
+          state = 2;
           break;
 
         default:
@@ -1973,7 +1972,7 @@ namespace getfem {
 
   // Syntax analysis of a string. Conversion to a tree. register the macros.
   void ga_read_string_reg(const std::string &expr, ga_tree &tree,
-                         ga_macro_dictionnary &macro_dict) {
+                          ga_macro_dictionary &macro_dict) {
     size_type pos = 0, token_pos, token_length;
     tree.clear();
     GA_TOKEN_TYPE t = ga_get_token(expr, pos, token_pos, token_length);
@@ -1995,8 +1994,8 @@ namespace getfem {
   // Syntax analysis of a string. Conversion to a tree.
   // Do not register the macros (but expand them).
   void ga_read_string(const std::string &expr, ga_tree &tree,
-                     const ga_macro_dictionnary &macro_dict) {
-    ga_macro_dictionnary macro_dict_loc(true, macro_dict);
+                      const ga_macro_dictionary &macro_dict) {
+    ga_macro_dictionary macro_dict_loc(true, macro_dict);
     ga_read_string_reg(expr, tree, macro_dict_loc);
   }
 
diff --git a/src/getfem_generic_assembly_workspace.cc 
b/src/getfem_generic_assembly_workspace.cc
index af872ec..8d719f4 100644
--- a/src/getfem_generic_assembly_workspace.cc
+++ b/src/getfem_generic_assembly_workspace.cc
@@ -80,7 +80,6 @@ namespace getfem {
        gmm::vect_size(VV)/(imd.nb_filtered_index() * imd.nb_tensor_elem()));
   }
 
-
   bool ga_workspace::variable_exists(const std::string &name) const {
     return (md && md->variable_exists(name)) ||
       (parent_workspace && parent_workspace->variable_exists(name)) ||
@@ -216,7 +215,7 @@ namespace getfem {
   size_type ga_workspace::qdim(const std::string &name) const {
     VAR_SET::const_iterator it = variables.find(name);
     if (it != variables.end()) {
-      const mesh_fem *mf =  it->second.is_fem_dofs ? it->second.mf : 0;
+      const mesh_fem *mf = it->second.is_fem_dofs ? it->second.mf : 0;
       const im_data *imd = it->second.imd;
       size_type n = it->second.qdim();
       if (mf) {
@@ -475,8 +474,8 @@ namespace getfem {
         for (const var_trans_pair &var : expr_variables) {
           if (!(is_constant(var.varname))) {
             ga_tree dtree = (remain ? tree : *(trees[ind_tree].ptree));
-            // cout << "Derivation with respect to " << var.first << " : "
-            //     << var.second << " of " << ga_tree_to_string(dtree) << endl;
+            // cout << "Derivation with respect to " << var.varname << " : "
+            //     << var.transname << " of " << ga_tree_to_string(dtree) << 
endl;
             GA_TIC;
             ga_derivative(dtree, *this, m, var.varname, var.transname, 
1+order);
             // cout << "Result : " << ga_tree_to_string(dtree) << endl;
@@ -514,7 +513,7 @@ namespace getfem {
     GA_TIC;
     size_type max_order = 0;
     std::vector<ga_tree> ltrees(1);
-    ga_read_string(expr, ltrees[0], macro_dictionnary());
+    ga_read_string(expr, ltrees[0], macro_dictionary());
     if (secondary_dom.size()) {
       GMM_ASSERT1(secondary_domain_exists(secondary_dom),
                  "Unknow secondary domain " << secondary_dom);
@@ -563,7 +562,7 @@ namespace getfem {
 
   void ga_workspace::add_function_expression(const std::string &expr) {
     ga_tree tree;
-    ga_read_string(expr, tree, macro_dictionnary());
+    ga_read_string(expr, tree, macro_dictionary());
     ga_semantic_analysis(tree, *this, dummy_mesh(), 1, false, true);
     if (tree.root) {
       // GMM_ASSERT1(tree.root->nb_test_functions() == 0,
@@ -578,7 +577,7 @@ namespace getfem {
                                                   const mesh_region &rg_) {
     const mesh_region &rg = register_region(m, rg_);
     ga_tree tree;
-    ga_read_string(expr, tree, macro_dictionnary());
+    ga_read_string(expr, tree, macro_dictionary());
     ga_semantic_analysis(tree, *this, m, ref_elt_dim_of_mesh(m),
                          false, false);
     if (tree.root) {
@@ -594,7 +593,7 @@ namespace getfem {
     const mesh &m = mim.linked_mesh();
     const mesh_region &rg = register_region(m, rg_);
     ga_tree tree;
-    ga_read_string(expr, tree, macro_dictionnary());
+    ga_read_string(expr, tree, macro_dictionary());
     ga_semantic_analysis(tree, *this, m, ref_elt_dim_of_mesh(m),
                          false, false);
     if (tree.root) {
@@ -613,7 +612,7 @@ namespace getfem {
     const mesh &m = mim.linked_mesh();
     const mesh_region &rg = register_region(m, rg_);
     ga_tree tree;
-    ga_read_string(expr, tree, macro_dictionnary());
+    ga_read_string(expr, tree, macro_dictionary());
     ga_semantic_analysis(tree, *this, m, ref_elt_dim_of_mesh(m), false, false);
     if (tree.root) {
       GMM_ASSERT1(tree.root->nb_test_functions() == 0,
@@ -914,11 +913,11 @@ namespace getfem {
                             bool enable_all_variables)
     : md(&md_), parent_workspace(0),
       enable_all_md_variables(enable_all_variables),
-      macro_dict(md_.macro_dictionnary())
+      macro_dict(md_.macro_dictionary())
   { init(); }
   ga_workspace::ga_workspace(bool, const ga_workspace &gaw)
     : md(0), parent_workspace(&gaw), enable_all_md_variables(false),
-      macro_dict(gaw.macro_dictionnary())
+      macro_dict(gaw.macro_dictionary())
   { init(); }
   ga_workspace::ga_workspace()
     : md(0), parent_workspace(0), enable_all_md_variables(false)
diff --git a/src/getfem_models.cc b/src/getfem_models.cc
index b3a7831..9d0f275 100644
--- a/src/getfem_models.cc
+++ b/src/getfem_models.cc
@@ -338,7 +338,7 @@ namespace getfem {
 
 
     // In case of change in fems or mims, linear terms have to be recomputed
-    // We couls select which brick is to be recomputed if we would be able
+    // We could select which brick is to be recomputed if we would be able
     // to know which fem or mim is changed.
     for (dal::bv_visitor ib(valid_bricks); !ib.finished(); ++ib)
       bricks[ib].terms_to_be_computed = true;
@@ -346,7 +346,7 @@ namespace getfem {
     for (auto &&v : variables) {
       const std::string &vname = v.first;
       var_description &vdescr = v.second;
-      if (vdescr.is_fem_dofs && !(vdescr.is_affine_dependent)) {
+      if (vdescr.is_fem_dofs && !vdescr.is_affine_dependent) {
         if ((vdescr.filter & VDESCRFILTER_CTERM)
             || (vdescr.filter & VDESCRFILTER_INFSUP)) {
           VAR_SET::iterator vfilt = variables.find(vdescr.filter_var);
@@ -470,17 +470,14 @@ namespace getfem {
                         workspace.assembly(2);
                     );
                   } //accumulated_distro scope
-                  gmm::add
-                    (gmm::sub_matrix(rTM, vdescr.I, multdescr.I), MM);
+                  gmm::add(gmm::sub_matrix(rTM, vdescr.I, multdescr.I), MM);
                   gmm::add(gmm::transposed
-                           (gmm::sub_matrix(rTM, multdescr.I,
-                                            vdescr.I)), MM);
+                           (gmm::sub_matrix(rTM, multdescr.I, vdescr.I)), MM);
                   bupd = false;
                 }
               }
             }
 
-
             for (size_type j = 0; j < brick.tlist.size(); ++j) {
               const term_description &term = brick.tlist[j];
 
@@ -516,8 +513,8 @@ namespace getfem {
                     bupd = true;
                   }
                   if (cplx)
-                    gmm::add
-                      (gmm::transposed(gmm::real_part(brick.cmatlist[j])), MM);
+                    
gmm::add(gmm::transposed(gmm::real_part(brick.cmatlist[j])),
+                             MM);
                   else
                     gmm::add(gmm::transposed(brick.rmatlist[j]), MM);
                   termadded = true;
@@ -581,16 +578,13 @@ namespace getfem {
 
 //         #if GETFEM_PARA_LEVEL > 1
 //         if (!rk) cout << "Range basis for  multipliers for " << vname << " 
time : " << MPI_Wtime()-tt_ref << endl;
-
 //         #endif
 
       if (mults.size() > 1) {
         gmm::range_basis(MGLOB, glob_columns, 1E-12, gmm::col_major(), true);
 
-
 //         #if GETFEM_PARA_LEVEL > 1
 //         if (!rk) cout << "Producing partial mf for  multipliers for " << 
vname << " time : " << MPI_Wtime()-tt_ref << endl;
-
 //         #endif
 
         s = 0;
@@ -610,7 +604,6 @@ namespace getfem {
       }
 //       #if GETFEM_PARA_LEVEL > 1
 //       if (!rk) cout << "End compute size of  multipliers for " << vname << 
" time : " << MPI_Wtime()-tt_ref << endl;
-
 //       #endif
     }
 
@@ -653,15 +646,15 @@ namespace getfem {
         ost << endl;
       }
       for (auto it = variable_groups.begin();
-          it != variable_groups.end(); ++it) {
-       ost << "Variable group " << std::setw(30) << std::left
-           << it->first;
-       if (it->second.size()) {
-         auto it2 = it->second.begin();
-         ost << " " << *it2; ++it2;
-         for (; it2 != it->second.end(); ++it2) ost << ", " << *it2;
-         ost << endl;
-       } else ost << " empty" << endl;
+           it != variable_groups.end(); ++it) {
+        ost << "Variable group " << std::setw(30) << std::left
+            << it->first;
+        if (it->second.size()) {
+          auto it2 = it->second.begin();
+          ost << " " << *it2; ++it2;
+          for (; it2 != it->second.end(); ++it2) ost << ", " << *it2;
+          ost << endl;
+        } else ost << " empty" << endl;
       }
     }
   }
@@ -706,7 +699,6 @@ namespace getfem {
     variables[name].set_size();
   }
 
-
   void model::resize_fixed_size_variable(const std::string &name,
                                          size_type size) {
     GMM_ASSERT1(!(variables[name].is_fem_dofs),
@@ -729,7 +721,6 @@ namespace getfem {
     variables[name].set_size();
   }
 
-
   void model::add_fixed_size_data(const std::string &name, size_type size,
                                   size_type niter) {
     check_name_validity(name);
@@ -935,15 +926,15 @@ namespace getfem {
      valid_bricks.del(ib);
      active_bricks.del(ib);
 
-     for  (size_type i = 0; i < bricks[ib].mims.size(); ++i) {
+     for (size_type i = 0; i < bricks[ib].mims.size(); ++i) {
        const mesh_im *mim = bricks[ib].mims[i];
        bool found = false;
        for (dal::bv_visitor ibb(valid_bricks); !ibb.finished(); ++ibb) {
          for  (size_type j = 0; j < bricks[ibb].mims.size(); ++j)
            if (bricks[ibb].mims[j] == mim) found = true;
        }
-       for(VAR_SET::iterator it2 = variables.begin();
-           it2 != variables.end(); ++it2) {
+       for (VAR_SET::iterator it2 = variables.begin();
+            it2 != variables.end(); ++it2) {
          if (it2->second.is_fem_dofs &&
              (it2->second.filter & VDESCRFILTER_INFSUP) &&
              mim == it2->second.mim) found = true;
@@ -1740,22 +1731,16 @@ namespace getfem {
 
       if (brick.is_update_brick) //contributions of pure update bricks must be 
deleted
       {
-        for (size_type i = 0; i < brick.cmatlist.size(); ++i)
-        {
-          gmm::clear(brick.cmatlist[i]);
-        }
+        for (auto &&mat : brick.cmatlist)
+          gmm::clear(mat);
 
-        for (size_type i = 0; i < brick.cveclist.size(); ++i)
-          for (size_type j = 0; j < brick.cveclist[i].size(); ++j)
-          {
-            gmm::clear(brick.cveclist[i][j]);
-          }
+        for (auto &&vecs : brick.cveclist)
+          for (auto &&vec : vecs)
+            gmm::clear(vec);
 
-        for (size_type i = 0; i < brick.cveclist_sym.size(); ++i)
-          for (size_type j = 0; j < brick.cveclist_sym[i].size(); ++j)
-          {
-            gmm::clear(brick.cveclist_sym[i][j]);
-          }
+        for (auto &&vecs : brick.cveclist_sym)
+          for (auto &&vec : vecs)
+            gmm::clear(vec);
       }
     }
     else //not cplx
@@ -1792,22 +1777,16 @@ namespace getfem {
 
       if (brick.is_update_brick) //contributions of pure update bricks must be 
deleted
       {
-        for (size_type i = 0; i < brick.rmatlist.size(); ++i)
-        {
-          gmm::clear(brick.rmatlist[i]);
-        }
+        for (auto &&mat : brick.rmatlist)
+          gmm::clear(mat);
 
-        for (size_type i = 0; i < brick.rveclist.size(); ++i)
-          for (size_type j = 0; j < brick.rveclist[i].size(); ++j)
-          {
-            gmm::clear(brick.rveclist[i][j]);
-          }
+        for (auto &&vecs : brick.rveclist)
+          for (auto &&vec : vecs)
+            gmm::clear(vec);
 
-        for (size_type i = 0; i < brick.rveclist_sym.size(); ++i)
-          for (size_type j = 0; j < brick.rveclist_sym[i].size(); ++j)
-          {
-            gmm::clear(brick.rveclist_sym[i][j]);
-          }
+        for (auto &&vecs : brick.rveclist_sym)
+          for (auto &&vec : vecs)
+            gmm::clear(vec);
       }
     }
   }
@@ -2179,7 +2158,7 @@ namespace getfem {
                               .region_is_faces_of(m, brick.region,
                                                   pmim->linked_mesh()));
         GMM_ASSERT1(ifo >= 0,
-                   "The given region is only partially covered by "
+                    "The given region is only partially covered by "
                     "region of brick \"" << brick.pbr->brick_name()
                     << "\". Please subdivise the region");
         if (ifo == 1) {
@@ -2188,7 +2167,7 @@ namespace getfem {
 
           ga_workspace workspace(*this);
           size_type order = workspace.add_expression
-           (expr, dummy_mim, region);
+            (expr, dummy_mim, region);
           GMM_ASSERT1(order <= 1, "Wrong order for a Neumann term");
           expr = workspace.extract_Neumann_term(varname);
           if (expr.size()) {
@@ -2409,7 +2388,7 @@ namespace getfem {
               if (brick.pbr->is_linear()
                   && (!is_linear() || (version & BUILD_WITH_COMPLETE_RHS))) {
                 gmm::mult_add(gmm::transposed(brick.rmatlist[j]),
-                             gmm::scaled(it1->second.complex_value[0],
+                              gmm::scaled(it1->second.complex_value[0],
                                           complex_type(-alpha2)),
                               gmm::sub_vector(crhs, I2));
               }
@@ -2499,7 +2478,6 @@ namespace getfem {
 //           brick.rveclist[0] = real_veclist(brick.tlist.size());
 //         }
 
-
       if (version & BUILD_RHS) approx_external_load_ += brick.external_load;
     }
 
@@ -2549,8 +2527,8 @@ namespace getfem {
         )
       } //end of distro scope
 
-      if (version & BUILD_RHS) gmm::add(gmm::scaled(residual, 
scalar_type(-1)), rrhs);
-
+      if (version & BUILD_RHS)
+        gmm::add(gmm::scaled(residual, scalar_type(-1)), rrhs);
     }
 
     // Post simplification for dof constraints
@@ -2606,7 +2584,7 @@ namespace getfem {
             }
           }
         }
-      } else {
+      } else { // !is_complex()
         std::vector<size_type> dof_indices;
         std::vector<scalar_type> dof_pr_values;
         std::vector<scalar_type> dof_go_values;
@@ -2678,7 +2656,6 @@ namespace getfem {
       approx_external_load_ = MPI_SUM_SCALAR(approx_external_load_);
     }
 
-
     #if GETFEM_PARA_LEVEL > 1
     // int rk; MPI_Comm_rank(MPI_COMM_WORLD, &rk);
     if (MPI_IS_MASTER()) cout << "Assembly time " << MPI_Wtime()-t_ref << endl;
@@ -3164,7 +3141,7 @@ model_complex_plain_vector &
                                    const model::varnamelist &vl_test1_,
                                    const std::string &directvarname_,
                                    const std::string &directdataname_,
-                                  const std::string &secdom)
+                                   const std::string &secdom)
       : vl_test1(vl_test1_), secondary_domain(secdom) {
       if (brickname.size() == 0)
         brickname = "Generic source term assembly brick";
@@ -3199,7 +3176,7 @@ model_complex_plain_vector &
 
     ga_workspace workspace(md);
     size_type order = workspace.add_expression(expr, mim, region, 1,
-                                              secondary_domain);
+                                               secondary_domain);
     GMM_ASSERT1(order <= 1, "Wrong order for a source term");
     model::varnamelist vl, vl_test1, vl_test2, dl;
     bool is_lin = workspace.used_variables(vl, vl_test1, vl_test2, dl, 1);
@@ -3313,7 +3290,7 @@ model_complex_plain_vector &
                               bool is_coer, std::string brickname,
                               const model::varnamelist &vl_test1_,
                               const model::varnamelist &vl_test2_,
-                             const std::string &secdom)
+                              const std::string &secdom)
       : vl_test1(vl_test1_), vl_test2(vl_test2_), secondary_domain(secdom) {
       if (brickname.size() == 0) brickname = "Generic linear assembly brick";
       expr = expr_;
@@ -3359,7 +3336,7 @@ model_complex_plain_vector &
 
     ga_workspace workspace(md, true);
     size_type order = workspace.add_expression(expr, mim, region,
-                                              2, secondary_domain);
+                                               2, secondary_domain);
     model::varnamelist vl, vl_test1, vl_test2, dl;
     bool is_lin = workspace.used_variables(vl, vl_test1, vl_test2, dl, 2);
 
@@ -3371,7 +3348,7 @@ model_complex_plain_vector &
     if (const_expr.size()) {
       add_source_term_
         (md, mim, const_expr, region, brickname+" (source term)",
-        "", "", false, secondary_domain);
+         "", "", false, secondary_domain);
     }
 
     // GMM_ASSERT1(order <= 1,
@@ -3384,7 +3361,7 @@ model_complex_plain_vector &
     if (vl_test1.size()) {
       pbrick pbr = std::make_shared<gen_linear_assembly_brick>
         (expr, mim, is_sym, is_coercive, brickname, vl_test1, vl_test2,
-        secondary_domain);
+         secondary_domain);
       model::termlist tl;
       for (size_type i = 0; i < vl_test1.size(); ++i)
         tl.push_back(model::term_description(vl_test1[i], vl_test2[i], false));
@@ -3399,7 +3376,7 @@ model_complex_plain_vector &
    bool is_sym, bool is_coercive, const std::string &brickname,
    bool return_if_nonlin) {
     return add_linear_term_(md, mim, expr, region, is_sym, is_coercive,
-                           brickname, return_if_nonlin, "");
+                            brickname, return_if_nonlin, "");
   }
 
   size_type add_linear_twodomain_term
@@ -3407,7 +3384,7 @@ model_complex_plain_vector &
    const std::string &secondary_domain, bool is_sym, bool is_coercive,
    const std::string &brickname, bool return_if_nonlin) {
     return add_linear_term_(md, mim, expr, region, is_sym, is_coercive,
-                           brickname, return_if_nonlin, secondary_domain);
+                            brickname, return_if_nonlin, secondary_domain);
   }
 
 
@@ -3448,8 +3425,8 @@ model_complex_plain_vector &
     gen_nonlinear_assembly_brick(const std::string &expr_, const mesh_im &mim,
                                  bool is_sym,
                                  bool is_coer,
-                                std::string brickname,
-                                const std::string &secdom) {
+                                 std::string brickname,
+                                 const std::string &secdom) {
       if (brickname.size() == 0) brickname = "Generic linear assembly brick";
       expr = expr_;
       secondary_domain = secdom;
@@ -3468,7 +3445,7 @@ model_complex_plain_vector &
 
     ga_workspace workspace(md);
     size_type order = workspace.add_expression(expr, mim, region, 2,
-                                              secondary_domain);
+                                               secondary_domain);
     GMM_ASSERT1(order < 2, "Order two test functions (Test2) are not allowed"
                 " in assembly string for nonlinear terms");
     model::varnamelist vl, vl_test1, vl_test2, ddl, dl;
@@ -3490,7 +3467,7 @@ model_complex_plain_vector &
   (model &md, const mesh_im &mim, const std::string &expr, size_type region,
    bool is_sym, bool is_coercive, const std::string &brickname) {
     return add_nonlinear_term_(md, mim, expr, region, is_sym, is_coercive,
-                       brickname, "");
+                               brickname, "");
   }
 
   size_type add_nonlinear_twodomain_term
@@ -3498,7 +3475,7 @@ model_complex_plain_vector &
    const std::string &secondary_domain, bool is_sym, bool is_coercive,
    const std::string &brickname) {
     return add_nonlinear_term_(md, mim, expr, region, is_sym, is_coercive,
-                       brickname, secondary_domain);
+                               brickname, secondary_domain);
   }
 
 
@@ -3751,7 +3728,7 @@ model_complex_plain_vector &
       else
         expr = "Grad_"+varname+":Grad_"+test_varname;
       return add_linear_term(md, mim, expr, region, true, true,
-                            "Laplacian", false);
+                             "Laplacian", false);
     }
   }
 



reply via email to

[Prev in Thread] Current Thread [Next in Thread]