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: Tue, 19 Feb 2019 17:15:14 -0500 (EST)

branch: master
commit f374fd75eabd2596d3a2fd766ec08932c3307e40
Author: Konstantinos Poulios <address@hidden>
Date:   Tue Feb 19 23:15:04 2019 +0100

    Whitespace, typos, file encodings
---
 interface/src/gf_asm.cc                       |  80 +++----
 src/getfem/getfem_assembling.h                |   8 +-
 src/getfem/getfem_generic_assembly_semantic.h |  34 +--
 src/getfem/getfem_generic_assembly_tree.h     |  16 +-
 src/getfem/getfem_models.h                    |  16 +-
 src/getfem_generic_assembly_tree.cc           |   2 +-
 src/getfem_generic_assembly_workspace.cc      |  14 +-
 src/gmm/gmm_matrix.h                          | 330 +++++++++++++-------------
 8 files changed, 248 insertions(+), 252 deletions(-)

diff --git a/interface/src/gf_asm.cc b/interface/src/gf_asm.cc
index 42416af..30c2010 100644
--- a/interface/src/gf_asm.cc
+++ b/interface/src/gf_asm.cc
@@ -70,7 +70,7 @@ public:
   }
   const bgeot::multi_index &sizes(getfem::size_type) const {  return sizes_; }
   virtual void compute(getfem::fem_interpolation_context& ctx,
-                      bgeot::base_tensor &t) {
+                       bgeot::base_tensor &t) {
     bgeot::size_type cv = ctx.convex_num();
     coeff.resize(mf.nb_basic_dof_of_element(cv));
     gmm::copy
@@ -95,7 +95,7 @@ void asm_lsneuman_matrix
     nterm(ls.get_mesh_fem(), ls.values());
 
   getfem::generic_assembly assem("t=comp(Base(#2).Grad(#1).NonLin(#3));"
-                                "M(#2, #1)+= t(:,:,i,i)");
+                                 "M(#2, #1)+= t(:,:,i,i)");
   assem.push_mi(mim);
   assem.push_mf(mf);
   assem.push_mf(mf_mult);
@@ -107,7 +107,7 @@ void asm_lsneuman_matrix
 
 
     /**
-       generic normal grad level set matrix (on the whole mesh level set or
+        generic normal grad level set matrix (on the whole mesh level set or
         on the specified element set level set or boundary level set)
     */
 
@@ -123,7 +123,7 @@ template<typename MAT>  void asm_nlsgrad_matrix
 
   getfem::generic_assembly
     assem("t=comp(Grad(#1).NonLin(#3).Grad(#2).NonLin(#3));"
-         "M(#1, #2)+= sym(t(:,i,i,:,j,j))");
+          "M(#1, #2)+= sym(t(:,i,i,:,j,j))");
   assem.push_mi(mim);
   assem.push_mf(mf1);
     assem.push_mf(mf2);
@@ -228,11 +228,11 @@ void asm_stabilization_patch_matrix
   int options[5] = {0,0,0,0,0};
   // float ubvec[1] = {1.03f};
   //METIS_mCPartGraphKway(&ne, &ncon, &(xadj[0]), &(adjncy[0]), &(vwgt[0]), 
&(adjwgt[0]), &wgtflag,
-  //               &numflag, &nparts, &(ubvec[0]),  options, &edgecut, 
&(part[0]));
+  //                      &numflag, &nparts, &(ubvec[0]),  options, &edgecut, 
&(part[0]));
   //METIS_mCPartGraphRecursive(&ne, &ncon, &(xadj[0]), &(adjncy[0]), 
&(vwgt[0]), &(adjwgt[0]), &wgtflag,
-  //                    &numflag, &nparts,  options, &edgecut, &(part[0]));
+  //                           &numflag, &nparts,  options, &edgecut, 
&(part[0]));
   //METIS_PartGraphKway(&ne, &(xadj[0]), &(adjncy[0]), &(vwgt[0]), 
&(adjwgt[0]), &wgtflag,
-  //     &numflag, &nparts, options, &edgecut, &(part[0]));
+  //                    &numflag, &nparts, options, &edgecut, &(part[0]));
   METIS_PartGraphRecursive(&ne, &(xadj[0]), &(adjncy[0]), &(vwgt[0]), 
&(adjwgt[0]), &wgtflag,
                            &numflag, &nparts, options, &edgecut, &(part[0]));
 # else
@@ -280,7 +280,7 @@ void asm_stabilization_patch_matrix
     bgeot::size_type cv = mf_P0.first_convex_of_basic_dof(r);
     int p=part[indelt[cv]];
     gmm::copy(gmm::scaled(gmm::mat_row(MAT_aux, p), 1./size_patch[p]),
-             gmm::mat_row(MAT_proj, r));
+              gmm::mat_row(MAT_proj, r));
   }
 
   gmm::mult(M0, MAT_proj, M1);
@@ -446,16 +446,16 @@ static void do_high_level_generic_assembly(mexargs_in& 
in, mexargs_out& out) {
   while (in.remaining()) {
     std::string varname = in.pop().to_string();
     if (varname.compare("select_output") == 0 ||
-       varname.compare("select output") == 0) {
+        varname.compare("select output") == 0) {
       GMM_ASSERT1(order > 0, "select_output option is for order 1 or 2"
                   "assemblies only");
       with_select_output = true;
       select_var1 = in.pop().to_string();
       if (order == 2) select_var2 = in.pop().to_string();
     } else if (varname.compare("Secondary_domain") == 0 ||
-       varname.compare("Secondary_Domain") == 0) {
+        varname.compare("Secondary_Domain") == 0) {
       GMM_ASSERT1(!with_secondary,
-                 "Only one secondary domain can be specified");
+                  "Only one secondary domain can be specified");
       secondary_domain = in.pop().to_string();
       with_secondary = true;
     } else {
@@ -463,37 +463,37 @@ static void do_high_level_generic_assembly(mexargs_in& 
in, mexargs_out& out) {
       const getfem::mesh_fem *mf(0);
       const getfem::im_data *mimd(0);
       if (is_meshfem_object(in.front()))
-       mf = to_meshfem_object(in.pop());
+        mf = to_meshfem_object(in.pop());
       else if (is_meshimdata_object(in.front()))
-       mimd = to_meshimdata_object(in.pop());
+        mimd = to_meshimdata_object(in.pop());
       darray U = in.pop().to_darray();
-      GMM_ASSERT1(vectors.find(varname) == vectors.end(),
-                 "The same variable/constant name is repeated twice: "
-                 << varname);
+      GMM_ASSERT1(vectors.count(varname) == 0,
+                  "The same variable/constant name is repeated twice: "
+                  << varname);
       GMM_ASSERT1(!with_model || !md.variable_exists(varname),
-                 "The same variable/constant name is already defined in "
-                 "the model: " << varname);
+                  "The same variable/constant name is already defined in "
+                  "the model: " << varname);
       gmm::resize(vectors[varname], U.size());
       gmm::copy(U, vectors[varname]);
       if (is_cte) {
-       if (mf)
-         workspace.add_fem_constant(varname, *mf, vectors[varname]);
-       else if (mimd)
-         workspace.add_im_data(varname, *mimd, vectors[varname]);
-       else
-         workspace.add_fixed_size_constant(varname, vectors[varname]);
+        if (mf)
+          workspace.add_fem_constant(varname, *mf, vectors[varname]);
+        else if (mimd)
+          workspace.add_im_data(varname, *mimd, vectors[varname]);
+        else
+          workspace.add_fixed_size_constant(varname, vectors[varname]);
       } else {
-       if (mf) {
-         gmm::sub_interval I(nbdof, mf->nb_dof());
-         nbdof += mf->nb_dof();
-         workspace.add_fem_variable(varname, *mf, I, vectors[varname]);
-       }  else if (mimd) {
-         THROW_BADARG("Data defined on integration points can not be a 
variable");
-       }  else {
-         gmm::sub_interval I(nbdof, U.size());
-         nbdof += U.size();
-         workspace.add_fixed_size_variable(varname, I, vectors[varname]);
-       }
+        if (mf) {
+          gmm::sub_interval I(nbdof, mf->nb_dof());
+          nbdof += mf->nb_dof();
+          workspace.add_fem_variable(varname, *mf, I, vectors[varname]);
+        }  else if (mimd) {
+          THROW_BADARG("Data defined on integration points can not be a 
variable");
+        }  else {
+          gmm::sub_interval I(nbdof, U.size());
+          nbdof += U.size();
+          workspace.add_fixed_size_variable(varname, I, vectors[varname]);
+        }
       }
     }
   }
@@ -748,7 +748,7 @@ template <typename T> static inline void dummy_func(T &) {}
                        getfemint::mexargs_out& out)                     \
       { dummy_func(in); dummy_func(out); code }                         \
     };                                                                  \
-    psub_command psubc = std::make_shared<subc>();                     \
+    psub_command psubc = std::make_shared<subc>();                      \
     psubc->arg_in_min = arginmin; psubc->arg_in_max = arginmax;         \
     psubc->arg_out_min = argoutmin; psubc->arg_out_max = argoutmax;     \
     subc_tab[cmd_normalize(name)] = psubc;                              \
@@ -799,9 +799,9 @@ void gf_asm(getfemint::mexargs_in& m_in, 
getfemint::mexargs_out& m_out) {
       functions are only available for variables, not for constants.
 
       `select_output` is an optional parameter which allows to reduce the
-      output vecotr (for `order` equal to 1) or the matrix (for `order`
+      output vector (for `order` equal to 1) or the matrix (for `order`
       equal to 2) to the degrees of freedom of the specified variables.
-      One variable has to be specified for a vector ouptut and two for a
+      One variable has to be specified for a vector output and two for a
       matrix output.
 
       Note that if several variables are given, the assembly of the
@@ -963,7 +963,7 @@ void gf_asm(getfemint::mexargs_in& m_in, 
getfemint::mexargs_out& m_out) {
            getfem::asm_nonlinear_elasticity_rhs(B, *mim, *mf_u, U, mf_d,
                                                 param,*law, rg);
          } else if (cmd_strmatch(what, "incompressible tangent matrix")) {
-          const getfem::mesh_fem *mf_p = to_meshfem_object(in.pop());
+           const getfem::mesh_fem *mf_p = to_meshfem_object(in.pop());
            darray P = in.pop().to_darray(int(mf_p->nb_dof()));
            gf_real_sparse_by_col  K(mf_u->nb_dof(), mf_u->nb_dof());
            gf_real_sparse_by_col  B(mf_u->nb_dof(), mf_p->nb_dof());
@@ -972,7 +972,7 @@ void gf_asm(getfemint::mexargs_in& m_in, 
getfemint::mexargs_out& m_out) {
            out.pop().from_sparse(K);
            out.pop().from_sparse(B);
          } else if (cmd_strmatch(what, "incompressible rhs")) {
-          const getfem::mesh_fem *mf_p = to_meshfem_object(in.pop());
+           const getfem::mesh_fem *mf_p = to_meshfem_object(in.pop());
            darray P = in.pop().to_darray(int(mf_p->nb_dof()));
            darray RU = out.pop().create_darray_v(unsigned(mf_u->nb_dof()));
            darray RB = out.pop().create_darray_v(unsigned(mf_p->nb_dof()));
diff --git a/src/getfem/getfem_assembling.h b/src/getfem/getfem_assembling.h
index d958318..885cbdf 100644
--- a/src/getfem/getfem_assembling.h
+++ b/src/getfem/getfem_assembling.h
@@ -944,7 +944,7 @@ namespace getfem {
   }
   
   /** 
-      Stiffness matrix for linear elasticity, with Lam� coefficients
+      Stiffness matrix for linear elasticity, with Lamé coefficients
       @ingroup asm
   */
   template<class MAT, class VECT>
@@ -986,7 +986,7 @@ namespace getfem {
 
 
   /** 
-      Stiffness matrix for linear elasticity, with constant Lam� coefficients
+      Stiffness matrix for linear elasticity, with constant Lamé coefficients
       @ingroup asm
   */
   template<class MAT, class VECT>
@@ -1589,7 +1589,7 @@ namespace getfem {
              -> the constraint is simplified:
              we replace \int{(H_j.psi_j)*phi_i}=\int{R_j.psi_j} (sum over j)
              with             H_j*phi_i = R_j     
-             --> Le principe peut �tre faux : non identique � la projection
+             --> Le principe peut être faux : non identique à la projection
              L^2 et peut entrer en conccurence avec les autres ddl -> a revoir
            */
            if (tdof_u == tdof_rh &&
@@ -1684,7 +1684,7 @@ namespace getfem {
        gmm::mult(H, e, aux);
        for (size_type j = 0; j < nb_bimg; ++j) { 
          T c = gmm::vect_sp(aux, base_img[j]);
-         // if (gmm::abs(c > 1.0E-6) { // � scaler sur l'ensemble de H ...
+         // if (gmm::abs(c > 1.0E-6) { // à scaler sur l'ensemble de H ...
          if (c != T(0)) {
            gmm::add(gmm::scaled(base_img[j], -c), aux);
            gmm::add(gmm::scaled(base_img_inv[j], -c), f);
diff --git a/src/getfem/getfem_generic_assembly_semantic.h 
b/src/getfem/getfem_generic_assembly_semantic.h
index 98ae7b4..f129065 100644
--- a/src/getfem/getfem_generic_assembly_semantic.h
+++ b/src/getfem/getfem_generic_assembly_semantic.h
@@ -50,33 +50,33 @@ namespace getfem {
      option = 0 : strict analysis,
               1 : do not complain about incompatible test functions but
                   store them,
-             2 : cut incompatible test function branches with respect to the
+              2 : cut incompatible test function branches with respect to the
                   one in workspace.selected_pair
               3 : do not complain about incompatible test functions neither
                   store them.
   */
   void ga_semantic_analysis(ga_tree &tree,
-                           const ga_workspace &workspace,
-                           const mesh &m,
-                           size_type ref_elt_dim,
-                           bool eval_fixed_size,
-                           bool ignore_X, int option = 0);
+                            const ga_workspace &workspace,
+                            const mesh &m,
+                            size_type ref_elt_dim,
+                            bool eval_fixed_size,
+                            bool ignore_X, int option = 0);
 
   /* Extract the variables used in a sub-tree, ignoring or not the data.
      Variable groups are taken into account. Return if at least one variable
      as been detected or not */
   bool ga_extract_variables(const pga_tree_node pnode,
-                           const ga_workspace &workspace,
-                           const mesh &m,
-                           std::set<var_trans_pair> &vars,
-                           bool ignore_data);
+                            const ga_workspace &workspace,
+                            const mesh &m,
+                            std::set<var_trans_pair> &vars,
+                            bool ignore_data);
   
   /* Extract a sub tree which consists of the corresponding node and of
      the terms multiplying this term, but not the term in addition.
      The aim is to expand an expression is a sum of elementary factors.
      Complains if a nonlinear term is encountered. */
   void ga_extract_factor(ga_tree &result_tree, pga_tree_node pnode,
-                        pga_tree_node &new_pnode);
+                         pga_tree_node &new_pnode);
 
   /* Extract the constant term of degree 1 expressions. */
   bool ga_node_extract_constant_term
@@ -93,16 +93,16 @@ namespace getfem {
      The tree is modified and should be copied first and passed to
      ga_semantic_analysis after for enrichment. */
   void ga_derivative(ga_tree &tree, const ga_workspace &workspace,
-                    const mesh &m, const std::string &varname,
-                    const std::string &interpolatename,
-                    size_type order);
+                     const mesh &m, const std::string &varname,
+                     const std::string &interpolatename,
+                     size_type order);
 
   std::string ga_derivative_scalar_function(const std::string &expr,
-                                           const std::string &var);
+                                            const std::string &var);
 
   bool ga_is_affine(const ga_tree &tree, const ga_workspace &workspace,
-                   const std::string &varname,
-                   const std::string &interpolatename);
+                    const std::string &varname,
+                    const std::string &interpolatename);
   
   // Function of internal use
   inline size_type ref_elt_dim_of_mesh(const mesh &m) {
diff --git a/src/getfem/getfem_generic_assembly_tree.h 
b/src/getfem/getfem_generic_assembly_tree.h
index 84b0dcb..f7628d6 100644
--- a/src/getfem/getfem_generic_assembly_tree.h
+++ b/src/getfem/getfem_generic_assembly_tree.h
@@ -201,11 +201,11 @@ namespace getfem {
   typedef std::shared_ptr<std::string> pstring;
   // Print error message indicating the position in the assembly string
   void ga_throw_error_msg(pstring expr, size_type pos,
-                         const std::string &msg);
+                          const std::string &msg);
 
 # define ga_throw_error(expr, pos, msg)               \
-  { std::stringstream ss; ss << msg;                 \
-    ga_throw_error_msg(expr, pos, ss.str());         \
+  { std::stringstream ss; ss << msg;                  \
+    ga_throw_error_msg(expr, pos, ss.str());          \
     GMM_ASSERT1(false, "Error in assembly string" );  \
   }
   
@@ -391,7 +391,7 @@ namespace getfem {
       : node_type(GA_NODE_VOID), test_function_type(-1), qdim1(0), qdim2(0),
       nbc1(0), nbc2(0), nbc3(0), pos(0), expr(0), der1(0), der2(0),
         symmetric_op(false), hash_value(0) {}
-  ga_tree_node(GA_NODE_TYPE ty, size_type p, pstring expr_)
+    ga_tree_node(GA_NODE_TYPE ty, size_type p, pstring expr_)
       : node_type(ty), test_function_type(-1),
         qdim1(0), qdim2(0), nbc1(0), nbc2(0), nbc3(0),
       pos(p), expr(expr_), der1(0), der2(0), symmetric_op(false),
@@ -421,7 +421,7 @@ namespace getfem {
     void add_scalar(scalar_type val, size_type pos, pstring expr);
     void add_allindices(size_type pos, pstring expr);
     void add_name(const char *name, size_type length, size_type pos,
-                 pstring expr);
+                  pstring expr);
     void add_sub_tree(ga_tree &sub_tree);
     void add_params(size_type pos, pstring expr);
     void add_matrix(size_type pos, pstring expr);
@@ -446,9 +446,9 @@ namespace getfem {
       std::swap(secondary_domain, tree.secondary_domain);
     }
 
-  ga_tree() : root(nullptr), current_node(nullptr), secondary_domain() {}
+    ga_tree() : root(nullptr), current_node(nullptr), secondary_domain() {}
 
-  ga_tree(const ga_tree &tree) : root(nullptr), current_node(nullptr),
+    ga_tree(const ga_tree &tree) : root(nullptr), current_node(nullptr),
       secondary_domain(tree.secondary_domain)
     { if (tree.root) copy_node(tree.root, nullptr, root); }
 
@@ -473,7 +473,7 @@ namespace getfem {
   // Transform the expression of a node and its sub-nodes in the equivalent
   // assembly string sent to ostream str
   void ga_print_node(const pga_tree_node pnode,
-                    std::ostream &str);
+                     std::ostream &str);
   // The same for the whole tree, the result is a std::string
   std::string ga_tree_to_string(const ga_tree &tree);
 
diff --git a/src/getfem/getfem_models.h b/src/getfem/getfem_models.h
index 3dd3b4c..b8b697c 100644
--- a/src/getfem/getfem_models.h
+++ b/src/getfem/getfem_models.h
@@ -185,7 +185,7 @@ namespace getfem {
       model_complex_plain_vector affine_complex_value;
       scalar_type alpha;    // Factor for the affine dependent variables
       std::string org_name; // Name of the original variable for affine
-                            //  dependent variables
+                            // dependent variables
 
       // im data description
       const im_data *pim_data;
@@ -239,7 +239,7 @@ namespace getfem {
       { return is_complex ? complex_value[0].size() : real_value[0].size(); }
 
       void set_size();
-    };
+    }; // struct var_description
 
   public:
 
@@ -506,21 +506,21 @@ namespace getfem {
       active_bricks.add(ib);
     }
 
-    /** Disable a variable (and its attached mutlipliers).  */
+    /** Disable a variable (and its attached mutlipliers). */
     void disable_variable(const std::string &name);
 
-    /** Enable a variable (and its attached mutlipliers).  */
+    /** Enable a variable (and its attached mutlipliers). */
     void enable_variable(const std::string &name);
 
-    /** Says if a name corresponds to a declared variable.  */
+    /** States if a name corresponds to a declared variable. */
     bool variable_exists(const std::string &name) const;
 
     bool is_disabled_variable(const std::string &name) const;
 
-    /** Says if a name corresponds to a declared data or disabled variable.  */
+    /** States if a name corresponds to a declared data or disabled variable. 
*/
     bool is_data(const std::string &name) const;
 
-    /** Says if a name corresponds to a declared data.  */
+    /** States if a name corresponds to a declared data. */
     bool is_true_data(const std::string &name) const;
 
     bool is_affine_dependent_variable(const std::string &name) const;
@@ -882,7 +882,7 @@ namespace getfem {
       return cTM;
     }
 
-    /** Gives the access to the right hand side of the tangent linear system.
+    /** Gives access to the right hand side of the tangent linear system.
         For the real version. An assembly of the rhs has to be done first. */
     const model_real_plain_vector &real_rhs() const {
       GMM_ASSERT1(!complex_version, "This model is a complex one");
diff --git a/src/getfem_generic_assembly_tree.cc 
b/src/getfem_generic_assembly_tree.cc
index c553596..2d673c7 100644
--- a/src/getfem_generic_assembly_tree.cc
+++ b/src/getfem_generic_assembly_tree.cc
@@ -1559,7 +1559,7 @@ namespace getfem {
         case GA_MINUS: // unary -
           tree.add_op(GA_UNARY_MINUS, token_pos, expr);
           state = 1; break;
-          
+
         case GA_PLUS:  // unary +
           state = 1; break;
 
diff --git a/src/getfem_generic_assembly_workspace.cc 
b/src/getfem_generic_assembly_workspace.cc
index 15cb98e..506ad16 100644
--- a/src/getfem_generic_assembly_workspace.cc
+++ b/src/getfem_generic_assembly_workspace.cc
@@ -59,7 +59,7 @@ namespace getfem {
   (const std::string &name, const mesh_fem &mf,
    const model_real_plain_vector &VV) {
     GMM_ASSERT1(mf.nb_dof(), "The provided mesh_fem of variable" << name
-                             << "has zero degrees of freedom.");
+                             << "has zero degrees of freedom");
     size_type Q = gmm::vect_size(VV)/mf.nb_dof();
     if (Q == 0) Q = size_type(1);
     variables[name] = var_description(false, true, &mf,
@@ -374,9 +374,6 @@ namespace getfem {
   }
 
 
-
-  
-
   const mesh_region &
   ga_workspace::register_region(const mesh &m, const mesh_region &region) {
     if (&m == &dummy_mesh())
@@ -389,7 +386,6 @@ namespace getfem {
     return lmr.back();
   }
 
-
   void ga_workspace::add_tree(ga_tree &tree, const mesh &m,
                               const mesh_im &mim, const mesh_region &rg,
                               const std::string &expr,
@@ -454,7 +450,7 @@ namespace getfem {
         trees.push_back(tree_description());
         trees.back().mim = &mim; trees.back().m = &m;
         trees.back().rg = &rg;
-       trees.back().secondary_domain = tree.secondary_domain;
+        trees.back().secondary_domain = tree.secondary_domain;
         trees.back().ptree = new ga_tree;
         trees.back().ptree->swap(tree);
         pga_tree_node root = trees.back().ptree->root;
@@ -507,7 +503,7 @@ namespace getfem {
                                          const mesh_im &mim,
                                          const mesh_region &rg_,
                                          size_type add_derivative_order,
-                                        const std::string &secondary_dom) {
+                                         const std::string &secondary_dom) {
     const mesh_region &rg = register_region(mim.linked_mesh(), rg_);
     // cout << "adding expression " << expr << endl;
     GA_TIC;
@@ -516,7 +512,7 @@ namespace getfem {
     ga_read_string(expr, ltrees[0], macro_dictionary());
     if (secondary_dom.size()) {
       GMM_ASSERT1(secondary_domain_exists(secondary_dom),
-                 "Unknow secondary domain " << secondary_dom);
+                  "Unknown secondary domain " << secondary_dom);
       ltrees[0].secondary_domain = secondary_dom;
     }
     // cout << "read : " << ga_tree_to_string(ltrees[0])  << endl;
@@ -910,7 +906,7 @@ namespace getfem {
   { if (ptree) delete ptree; ptree = 0; }
 
   ga_workspace::ga_workspace(const getfem::model &md_,
-                            bool enable_all_variables)
+                             bool enable_all_variables)
     : md(&md_), parent_workspace(0),
       enable_all_md_variables(enable_all_variables),
       macro_dict(md_.macro_dictionary())
diff --git a/src/gmm/gmm_matrix.h b/src/gmm/gmm_matrix.h
index 23fb9d2..3b133bb 100644
--- a/src/gmm/gmm_matrix.h
+++ b/src/gmm/gmm_matrix.h
@@ -48,9 +48,9 @@ namespace gmm
 {
 
   /* ******************************************************************** */
-  /*                                                                     */
-  /*           Identity matrix                                           */
-  /*                                                                     */
+  /*                                                                      */
+  /*            Identity matrix                                           */
+  /*                                                                      */
   /* ******************************************************************** */
 
   struct identity_matrix {
@@ -71,7 +71,7 @@ namespace gmm
   void mult(const identity_matrix&, const V1 &v1, V2 &v2)
   { copy(v1, v2); }
   template <typename V1, typename V2> inline
-  void mult(const identity_matrix&, const V1 &v1, const V2 &v2) 
+  void mult(const identity_matrix&, const V1 &v1, const V2 &v2)
   { copy(v1, v2); }
   template <typename V1, typename V2, typename V3> inline
   void mult(const identity_matrix&, const V1 &v1, const V2 &v2, V3 &v3)
@@ -83,25 +83,25 @@ namespace gmm
   void left_mult(const identity_matrix&, const V1 &v1, V2 &v2)
   { copy(v1, v2); }
   template <typename V1, typename V2> inline
-  void left_mult(const identity_matrix&, const V1 &v1, const V2 &v2) 
+  void left_mult(const identity_matrix&, const V1 &v1, const V2 &v2)
   { copy(v1, v2); }
   template <typename V1, typename V2> inline
   void right_mult(const identity_matrix&, const V1 &v1, V2 &v2)
   { copy(v1, v2); }
   template <typename V1, typename V2> inline
-  void right_mult(const identity_matrix&, const V1 &v1, const V2 &v2) 
+  void right_mult(const identity_matrix&, const V1 &v1, const V2 &v2)
   { copy(v1, v2); }
   template <typename V1, typename V2> inline
   void transposed_left_mult(const identity_matrix&, const V1 &v1, V2 &v2)
   { copy(v1, v2); }
   template <typename V1, typename V2> inline
-  void transposed_left_mult(const identity_matrix&, const V1 &v1,const V2 &v2) 
+  void transposed_left_mult(const identity_matrix&, const V1 &v1,const V2 &v2)
   { copy(v1, v2); }
   template <typename V1, typename V2> inline
   void transposed_right_mult(const identity_matrix&, const V1 &v1, V2 &v2)
   { copy(v1, v2); }
   template <typename V1, typename V2> inline
-  void transposed_right_mult(const identity_matrix&,const V1 &v1,const V2 &v2) 
+  void transposed_right_mult(const identity_matrix&,const V1 &v1,const V2 &v2)
   { copy(v1, v2); }
   template <typename M> void copy_ident(const identity_matrix&, M &m) {
     size_type i = 0, n = std::min(mat_nrows(m), mat_ncols(m));
@@ -109,7 +109,7 @@ namespace gmm
     for (; i < n; ++i) m(i,i) = typename linalg_traits<M>::value_type(1);
   }
   template <typename M> inline void copy(const identity_matrix&, M &m)
-  { copy_ident(identity_matrix(), m); } 
+  { copy_ident(identity_matrix(), m); }
   template <typename M> inline void copy(const identity_matrix &, const M &m)
   { copy_ident(identity_matrix(), linalg_const_cast(m)); }
   template <typename V1, typename V2> inline
@@ -124,24 +124,24 @@ namespace gmm
   inline bool is_identity(const identity_matrix&) { return true; }
 
   /* ******************************************************************** */
-  /*                                                                     */
-  /*           Row matrix                                                */
-  /*                                                                     */
+  /*                                                                      */
+  /*            Row matrix                                                */
+  /*                                                                      */
   /* ******************************************************************** */
 
   template<typename V> class row_matrix {
   protected :
     std::vector<V> li; /* array of rows.                                  */
     size_type nc;
-    
+
   public :
-    
+
     typedef typename linalg_traits<V>::reference reference;
     typedef typename linalg_traits<V>::value_type value_type;
-    
+
     row_matrix(size_type r, size_type c) : li(r, V(c)), nc(c) {}
     row_matrix(void) : nc(0) {}
-    reference operator ()(size_type l, size_type c) 
+    reference operator ()(size_type l, size_type c)
     { return li[l][c]; }
     value_type operator ()(size_type l, size_type c) const
     { return li[l][c]; }
@@ -151,19 +151,19 @@ namespace gmm
 
     typename std::vector<V>::iterator begin(void)
     { return li.begin(); }
-    typename std::vector<V>::iterator end(void)  
+    typename std::vector<V>::iterator end(void)
     { return li.end(); }
     typename std::vector<V>::const_iterator begin(void) const
     { return li.begin(); }
     typename std::vector<V>::const_iterator end(void) const
     { return li.end(); }
-    
-    
+
+
     V& row(size_type i) { return li[i]; }
     const V& row(size_type i) const { return li[i]; }
     V& operator[](size_type i) { return li[i]; }
     const V& operator[](size_type i) const { return li[i]; }
-    
+
     inline size_type nrows(void) const { return li.size(); }
     inline size_type ncols(void) const { return nc;        }
 
@@ -176,7 +176,7 @@ namespace gmm
     li.resize(m);
     for (size_type i=nr; i < m; ++i) gmm::resize(li[i], n);
     if (n != nc) {
-      for (size_type i=0; i < nr; ++i) gmm::resize(li[i], n);    
+      for (size_type i=0; i < nr; ++i) gmm::resize(li[i], n);
       nc = n;
     }
   }
@@ -213,7 +213,7 @@ namespace gmm
     { return m.end(); }
     static const_sub_row_type row(const const_row_iterator &it)
     { return const_sub_row_type(*it); }
-    static sub_row_type row(const row_iterator &it) 
+    static sub_row_type row(const row_iterator &it)
     { return sub_row_type(*it); }
     static origin_type* origin(this_type &m) { return &m; }
     static const origin_type* origin(const this_type &m) { return &m; }
@@ -232,21 +232,21 @@ namespace gmm
     (std::ostream &o, const row_matrix<V>& m) { gmm::write(o,m); return o; }
 
   /* ******************************************************************** */
-  /*                                                                     */
-  /*           Column matrix                                             */
-  /*                                                                     */
+  /*                                                                      */
+  /*            Column matrix                                             */
+  /*                                                                      */
   /* ******************************************************************** */
 
   template<typename V> class col_matrix {
   protected :
     std::vector<V> li; /* array of columns.                               */
     size_type nr;
-    
+
   public :
-    
+
     typedef typename linalg_traits<V>::reference reference;
     typedef typename linalg_traits<V>::value_type value_type;
-    
+
     col_matrix(size_type r, size_type c) : li(c, V(r)), nr(r) { }
     col_matrix(void) : nr(0) {}
     reference operator ()(size_type l, size_type c)
@@ -264,13 +264,13 @@ namespace gmm
 
     typename std::vector<V>::iterator begin(void)
     { return li.begin(); }
-    typename std::vector<V>::iterator end(void)  
+    typename std::vector<V>::iterator end(void)
     { return li.end(); }
     typename std::vector<V>::const_iterator begin(void) const
     { return li.begin(); }
     typename std::vector<V>::const_iterator end(void) const
     { return li.end(); }
-    
+
     inline size_type ncols(void) const { return li.size(); }
     inline size_type nrows(void) const { return nr; }
 
@@ -283,7 +283,7 @@ namespace gmm
     li.resize(n);
     for (size_type i=nc; i < n; ++i) gmm::resize(li[i], m);
     if (m != nr) {
-      for (size_type i=0; i < nc; ++i) gmm::resize(li[i], m);    
+      for (size_type i=0; i < nc; ++i) gmm::resize(li[i], m);
       nr = m;
     }
   }
@@ -319,7 +319,7 @@ namespace gmm
     { return m.end(); }
     static const_sub_col_type col(const const_col_iterator &it)
     { return *it; }
-    static sub_col_type col(const col_iterator &it) 
+    static sub_col_type col(const col_iterator &it)
     { return *it; }
     static origin_type* origin(this_type &m) { return &m; }
     static const origin_type* origin(const this_type &m) { return &m; }
@@ -338,9 +338,9 @@ namespace gmm
     (std::ostream &o, const col_matrix<V>& m) { gmm::write(o,m); return o; }
 
   /* ******************************************************************** */
-  /*                                                                     */
-  /*           Dense matrix                                              */
-  /*                                                                     */
+  /*                                                                      */
+  /*            Dense matrix                                              */
+  /*                                                                      */
   /* ******************************************************************** */
 
   template<typename T> class dense_matrix : public std::vector<T> {
@@ -350,12 +350,12 @@ namespace gmm
     typedef typename std::vector<T>::const_iterator const_iterator;
     typedef typename std::vector<T>::reference reference;
     typedef typename std::vector<T>::const_reference const_reference;
-    
+
   protected:
     size_type nbc, nbl;
-    
+
   public:
-    
+
     inline const_reference operator ()(size_type l, size_type c) const {
       GMM_ASSERT2(l < nbl && c < nbc, "out of range");
       return *(this->begin() + c*nbl+l);
@@ -371,13 +371,13 @@ namespace gmm
     void resize(size_type, size_type);
     void base_resize(size_type, size_type);
     void reshape(size_type, size_type);
-    
+
     void fill(T a, T b = T(0));
     inline size_type nrows(void) const { return nbl; }
     inline size_type ncols(void) const { return nbc; }
     void swap(dense_matrix<T> &m)
     { std::vector<T>::swap(m); std::swap(nbc, m.nbc); std::swap(nbl, m.nbl); }
-    
+
     dense_matrix(size_type l, size_type c)
       : std::vector<T>(c*l), nbc(c), nbl(l)  {}
     dense_matrix(void) { nbl = nbc = 0; }
@@ -389,33 +389,33 @@ namespace gmm
   }
 
   template<typename T> void dense_matrix<T>::base_resize(size_type m,
-                                                        size_type n)
+                                                         size_type n)
   { std::vector<T>::resize(n*m); nbl = m; nbc = n; }
-  
+
   template<typename T> void dense_matrix<T>::resize(size_type m, size_type n) {
     if (n*m > nbc*nbl) std::vector<T>::resize(n*m);
     if (m < nbl) {
       for (size_type i = 1; i < std::min(nbc, n); ++i)
-       std::copy(this->begin()+i*nbl, this->begin()+(i*nbl+m),
-                 this->begin()+i*m);
+        std::copy(this->begin()+i*nbl, this->begin()+(i*nbl+m),
+                  this->begin()+i*m);
       for (size_type i = std::min(nbc, n); i < n; ++i)
-       std::fill(this->begin()+(i*m), this->begin()+(i+1)*m, T(0));
+        std::fill(this->begin()+(i*m), this->begin()+(i+1)*m, T(0));
       }
     else if (m > nbl) { /* do nothing when the nb of rows does not change */
       for (size_type i = std::min(nbc, n); i > 1; --i)
-       std::copy(this->begin()+(i-1)*nbl, this->begin()+i*nbl,
-                 this->begin()+(i-1)*m);
+        std::copy(this->begin()+(i-1)*nbl, this->begin()+i*nbl,
+                  this->begin()+(i-1)*m);
       for (size_type i = 0; i < std::min(nbc, n); ++i)
-       std::fill(this->begin()+(i*m+nbl), this->begin()+(i+1)*m, T(0));
+        std::fill(this->begin()+(i*m+nbl), this->begin()+(i+1)*m, T(0));
     }
     if (n*m < nbc*nbl) std::vector<T>::resize(n*m);
     nbl = m; nbc = n;
   }
-  
+
   template<typename T> void dense_matrix<T>::fill(T a, T b) {
     std::fill(this->begin(), this->end(), b);
     size_type n = std::min(nbl, nbc);
-    if (a != b) for (size_type i = 0; i < n; ++i) (*this)(i,i) = a; 
+    if (a != b) for (size_type i = 0; i < n; ++i) (*this)(i,i) = a;
   }
 
   template <typename T> struct linalg_traits<dense_matrix<T> > {
@@ -427,25 +427,25 @@ namespace gmm
     typedef T& reference;
     typedef abstract_dense storage_type;
     typedef tab_ref_reg_spaced_with_origin<typename this_type::iterator,
-                                          this_type> sub_row_type;
+                                           this_type> sub_row_type;
     typedef tab_ref_reg_spaced_with_origin<typename this_type::const_iterator,
-                                          this_type> const_sub_row_type;
+                                           this_type> const_sub_row_type;
     typedef dense_compressed_iterator<typename this_type::iterator,
-                                     typename this_type::iterator,
-                                     this_type *> row_iterator;
+                                      typename this_type::iterator,
+                                      this_type *> row_iterator;
     typedef dense_compressed_iterator<typename this_type::const_iterator,
-                                     typename this_type::iterator,
-                                     const this_type *> const_row_iterator;
-    typedef tab_ref_with_origin<typename this_type::iterator, 
-                               this_type> sub_col_type;
+                                      typename this_type::iterator,
+                                      const this_type *> const_row_iterator;
+    typedef tab_ref_with_origin<typename this_type::iterator,
+                                this_type> sub_col_type;
     typedef tab_ref_with_origin<typename this_type::const_iterator,
-                               this_type> const_sub_col_type;
+                                this_type> const_sub_col_type;
     typedef dense_compressed_iterator<typename this_type::iterator,
-                                     typename this_type::iterator,
-                                     this_type *> col_iterator;
+                                      typename this_type::iterator,
+                                      this_type *> col_iterator;
     typedef dense_compressed_iterator<typename this_type::const_iterator,
-                                     typename this_type::iterator,
-                                     const this_type *> const_col_iterator;
+                                      typename this_type::iterator,
+                                      const this_type *> const_col_iterator;
     typedef col_and_row sub_orientation;
     typedef linalg_true index_sorted;
     static size_type nrows(const this_type &m) { return m.nrows(); }
@@ -493,7 +493,7 @@ namespace gmm
 
   /* ******************************************************************** */
   /*                                                                      */
-  /*           Read only compressed sparse column matrix                 */
+  /*             Read only compressed sparse column matrix                */
   /*                                                                      */
   /* ******************************************************************** */
 
@@ -518,7 +518,7 @@ namespace gmm
     template <typename PT1, typename PT2, typename PT3, int cshift>
     void init_with(const csc_matrix_ref<PT1,PT2,PT3,cshift>& B)
     { init_with_good_format(B); }
-    template <typename U, int cshift>    
+    template <typename U, int cshift>
     void init_with(const csc_matrix<U, cshift>& B)
     { init_with_good_format(B); }
 
@@ -529,9 +529,9 @@ namespace gmm
 
     size_type nrows(void) const { return nr; }
     size_type ncols(void) const { return nc; }
-    void swap(csc_matrix<T, shift> &m) { 
-      std::swap(pr, m.pr); 
-      std::swap(ir, m.ir); std::swap(jc, m.jc); 
+    void swap(csc_matrix<T, shift> &m) {
+      std::swap(pr, m.pr);
+      std::swap(ir, m.ir); std::swap(jc, m.jc);
       std::swap(nc, m.nc); std::swap(nr, m.nr);
     }
     value_type operator()(size_type i, size_type j) const
@@ -552,30 +552,30 @@ namespace gmm
     for (size_type j = 0; j < nc; ++j) {
       col_type col = mat_const_col(B, j);
       typename linalg_traits<typename org_type<col_type>::t>::const_iterator
-       it = vect_const_begin(col), ite = vect_const_end(col);
+        it = vect_const_begin(col), ite = vect_const_end(col);
       for (size_type k = 0; it != ite; ++it, ++k) {
-       pr[jc[j]-shift+k] = *it;
-       ir[jc[j]-shift+k] = IND_TYPE(it.index() + shift);
+        pr[jc[j]-shift+k] = *it;
+        ir[jc[j]-shift+k] = IND_TYPE(it.index() + shift);
       }
     }
   }
-  
+
   template <typename T, int shift> template <typename Matrix>
   void csc_matrix<T, shift>::init_with(const Matrix &A) {
     col_matrix<wsvector<T> > B(mat_nrows(A), mat_ncols(A));
     copy(A, B);
     init_with_good_format(B);
   }
-  
+
   template <typename T, int shift>
   void csc_matrix<T, shift>::init_with_identity(size_type n) {
-    nc = nr = n; 
+    nc = nr = n;
     pr.resize(nc); ir.resize(nc); jc.resize(nc+1);
     for (size_type j = 0; j < nc; ++j)
       { ir[j] = jc[j] = shift + j; pr[j] = T(1); }
     jc[nc] = shift + nc;
   }
-  
+
   template <typename T, int shift>
   csc_matrix<T, shift>::csc_matrix(size_type nnr, size_type nnc)
     : nc(nnc), nr(nnr) {
@@ -601,7 +601,7 @@ namespace gmm
     typedef cs_vector_ref<const T *, const IND_TYPE *, shift>
     const_sub_col_type;
     typedef sparse_compressed_iterator<const T *, const IND_TYPE *,
-                                      const IND_TYPE *, shift>
+                                       const IND_TYPE *, shift>
     const_col_iterator;
     typedef abstract_null_type col_iterator;
     typedef col_major sub_orientation;
@@ -612,12 +612,12 @@ namespace gmm
     { return const_col_iterator(&m.pr[0],&m.ir[0],&m.jc[0], m.nr, &m.pr[0]); }
     static const_col_iterator col_end(const this_type &m) {
       return const_col_iterator(&m.pr[0],&m.ir[0],&m.jc[0]+m.nc,
-                               m.nr,&m.pr[0]);
+                                m.nr,&m.pr[0]);
     }
     static const_sub_col_type col(const const_col_iterator &it) {
       return const_sub_col_type(it.pr + *(it.jc) - shift,
-                               it.ir + *(it.jc) - shift,
-                               *(it.jc + 1) - *(it.jc), it.n);
+                                it.ir + *(it.jc) - shift,
+                                *(it.jc + 1) - *(it.jc), it.n);
     }
     static const origin_type* origin(const this_type &m) { return &m.pr[0]; }
     static void do_clear(this_type &m) { m.do_clear(); }
@@ -629,7 +629,7 @@ namespace gmm
   std::ostream &operator <<
     (std::ostream &o, const csc_matrix<T, shift>& m)
   { gmm::write(o,m); return o; }
-  
+
   template <typename T, int shift>
   inline void copy(const identity_matrix &, csc_matrix<T, shift>& M)
   { M.init_with_identity(mat_nrows(M)); }
@@ -640,7 +640,7 @@ namespace gmm
 
   /* ******************************************************************** */
   /*                                                                      */
-  /*           Read only compressed sparse row matrix                    */
+  /*             Read only compressed sparse row matrix                   */
   /*                                                                      */
   /* ******************************************************************** */
 
@@ -678,16 +678,16 @@ namespace gmm
 
     size_type nrows(void) const { return nr; }
     size_type ncols(void) const { return nc; }
-    void swap(csr_matrix<T, shift> &m) { 
-      std::swap(pr, m.pr); 
-      std::swap(ir,m.ir); std::swap(jc, m.jc); 
+    void swap(csr_matrix<T, shift> &m) {
+      std::swap(pr, m.pr);
+      std::swap(ir,m.ir); std::swap(jc, m.jc);
       std::swap(nc, m.nc); std::swap(nr,m.nr);
     }
-   
+
     value_type operator()(size_type i, size_type j) const
     { return mat_row(*this, i)[j]; }
   };
-  
+
   template <typename T, int shift> template <typename Matrix>
   void csr_matrix<T, shift>::init_with_good_format(const Matrix &B) {
     typedef typename linalg_traits<Matrix>::const_sub_row_type row_type;
@@ -702,24 +702,24 @@ namespace gmm
     for (size_type j = 0; j < nr; ++j) {
       row_type row = mat_const_row(B, j);
       typename linalg_traits<typename org_type<row_type>::t>::const_iterator
-       it = vect_const_begin(row), ite = vect_const_end(row);
+        it = vect_const_begin(row), ite = vect_const_end(row);
       for (size_type k = 0; it != ite; ++it, ++k) {
-       pr[jc[j]-shift+k] = *it;
-       ir[jc[j]-shift+k] = IND_TYPE(it.index()+shift);
+        pr[jc[j]-shift+k] = *it;
+        ir[jc[j]-shift+k] = IND_TYPE(it.index()+shift);
       }
     }
   }
 
-  template <typename T, int shift> template <typename Matrix> 
-  void csr_matrix<T, shift>::init_with(const Matrix &A) { 
-    row_matrix<wsvector<T> > B(mat_nrows(A), mat_ncols(A)); 
-    copy(A, B); 
+  template <typename T, int shift> template <typename Matrix>
+  void csr_matrix<T, shift>::init_with(const Matrix &A) {
+    row_matrix<wsvector<T> > B(mat_nrows(A), mat_ncols(A));
+    copy(A, B);
     init_with_good_format(B);
   }
 
-  template <typename T, int shift> 
+  template <typename T, int shift>
   void csr_matrix<T, shift>::init_with_identity(size_type n) {
-    nc = nr = n; 
+    nc = nr = n;
     pr.resize(nr); ir.resize(nr); jc.resize(nr+1);
     for (size_type j = 0; j < nr; ++j)
       { ir[j] = jc[j] = shift + j; pr[j] = T(1); }
@@ -753,7 +753,7 @@ namespace gmm
     typedef cs_vector_ref<const T *, const IND_TYPE *, shift>
     const_sub_row_type;
     typedef sparse_compressed_iterator<const T *, const IND_TYPE *,
-                                      const IND_TYPE *, shift>
+                                       const IND_TYPE *, shift>
     const_row_iterator;
     typedef abstract_null_type row_iterator;
     typedef row_major sub_orientation;
@@ -766,8 +766,8 @@ namespace gmm
     { return const_row_iterator(&m.pr[0], &m.ir[0], &m.jc[0] + m.nr, m.nc, 
&m.pr[0]); }
     static const_sub_row_type row(const const_row_iterator &it) {
       return const_sub_row_type(it.pr + *(it.jc) - shift,
-                               it.ir + *(it.jc) - shift,
-                               *(it.jc + 1) - *(it.jc), it.n);
+                                it.ir + *(it.jc) - shift,
+                                *(it.jc + 1) - *(it.jc), it.n);
     }
     static const origin_type* origin(const this_type &m) { return &m.pr[0]; }
     static void do_clear(this_type &m) { m.do_clear(); }
@@ -779,7 +779,7 @@ namespace gmm
   std::ostream &operator <<
     (std::ostream &o, const csr_matrix<T, shift>& m)
   { gmm::write(o,m); return o; }
-  
+
   template <typename T, int shift>
   inline void copy(const identity_matrix &, csr_matrix<T, shift>& M)
   { M.init_with_identity(mat_nrows(M)); }
@@ -789,9 +789,9 @@ namespace gmm
   { M.init_with(A); }
 
   /* ******************************************************************** */
-  /*                                                                     */
-  /*           Block matrix                                              */
-  /*                                                                     */
+  /*                                                                      */
+  /*             Block matrix                                             */
+  /*                                                                      */
   /* ******************************************************************** */
 
   template <typename MAT> class block_matrix {
@@ -811,7 +811,7 @@ namespace gmm
     size_type ncolblocks(void) const { return ncolblocks_; }
     const sub_interval &subrowinterval(size_type i) const { return introw[i]; }
     const sub_interval &subcolinterval(size_type i) const { return intcol[i]; }
-    const MAT &block(size_type i, size_type j) const 
+    const MAT &block(size_type i, size_type j) const
     { return blocks[j*ncolblocks_+i]; }
     MAT &block(size_type i, size_type j)
     { return blocks[j*ncolblocks_+i]; }
@@ -820,20 +820,20 @@ namespace gmm
     value_type operator() (size_type i, size_type j) const {
       size_type k, l;
       for (k = 0; k < nrowblocks_; ++k)
-       if (i >= introw[k].min && i <  introw[k].max) break;
+        if (i >= introw[k].min && i <  introw[k].max) break;
       for (l = 0; l < nrowblocks_; ++l)
-       if (j >= introw[l].min && j <  introw[l].max) break;
+        if (j >= introw[l].min && j <  introw[l].max) break;
       return (block(k, l))(i - introw[k].min, j - introw[l].min);
     }
     reference operator() (size_type i, size_type j) {
       size_type k, l;
       for (k = 0; k < nrowblocks_; ++k)
-       if (i >= introw[k].min && i <  introw[k].max) break;
+        if (i >= introw[k].min && i <  introw[k].max) break;
       for (l = 0; l < nrowblocks_; ++l)
-       if (j >= introw[l].min && j <  introw[l].max) break;
+        if (j >= introw[l].min && j <  introw[l].max) break;
       return (block(k, l))(i - introw[k].min, j - introw[l].min);
     }
-    
+
     template <typename CONT> void resize(const CONT &c1, const CONT &c2);
     template <typename CONT> block_matrix(const CONT &c1, const CONT &c2)
     { resize(c1, c2); }
@@ -864,17 +864,17 @@ namespace gmm
     static origin_type* origin(this_type &m) { return &m; }
     static const origin_type* origin(const this_type &m) { return &m; }
     static void do_clear(this_type &m) { m.do_clear(); }
-    // access to be done ...    
+    // access to be done ...
     static void resize(this_type &, size_type , size_type)
     { GMM_ASSERT1(false, "Sorry, to be done"); }
     static void reshape(this_type &, size_type , size_type)
     { GMM_ASSERT1(false, "Sorry, to be done"); }
   };
 
-  template <typename MAT> void block_matrix<MAT>::do_clear(void) { 
+  template <typename MAT> void block_matrix<MAT>::do_clear(void) {
     for (size_type j = 0, l = 0; j < ncolblocks_; ++j)
       for (size_type i = 0, k = 0; i < nrowblocks_; ++i)
-       clear(block(i,j));
+        clear(block(i,j));
   }
 
   template <typename MAT> template <typename CONT>
@@ -886,8 +886,8 @@ namespace gmm
     for (size_type j = 0, l = 0; j < ncolblocks_; ++j) {
       intcol[j] = sub_interval(l, c2[j]); l += c2[j];
       for (size_type i = 0, k = 0; i < nrowblocks_; ++i) {
-       if (j == 0) { introw[i] = sub_interval(k, c1[i]); k += c1[i]; }
-       block(i, j) = MAT(c1[i], c2[j]);
+        if (j == 0) { introw[i] = sub_interval(k, c1[i]); k += c1[i]; }
+        block(i, j) = MAT(c1[i], c2[j]);
       }
     }
   }
@@ -896,14 +896,14 @@ namespace gmm
   void copy(const block_matrix<M1> &m1, M2 &m2) {
     for (size_type j = 0; j < m1.ncolblocks(); ++j)
       for (size_type i = 0; i < m1.nrowblocks(); ++i)
-       copy(m1.block(i,j), sub_matrix(m2, m1.subrowinterval(i), 
-                                      m1.subcolinterval(j)));
+        copy(m1.block(i,j), sub_matrix(m2, m1.subrowinterval(i),
+                                       m1.subcolinterval(j)));
   }
 
   template <typename M1, typename M2>
   void copy(const block_matrix<M1> &m1, const M2 &m2)
   { copy(m1, linalg_const_cast(m2)); }
-  
+
 
   template <typename MAT, typename V1, typename V2>
   void mult(const block_matrix<MAT> &m, const V1 &v1, V2 &v2) {
@@ -911,9 +911,9 @@ namespace gmm
     typename sub_vector_type<V2 *, sub_interval>::vector_type sv;
     for (size_type i = 0; i < m.nrowblocks() ; ++i)
       for (size_type j = 0; j < m.ncolblocks() ; ++j) {
-       sv = sub_vector(v2, m.subrowinterval(i));
-       mult(m.block(i,j),
-            sub_vector(v1, m.subcolinterval(j)), sv, sv);
+        sv = sub_vector(v2, m.subrowinterval(i));
+        mult(m.block(i,j),
+             sub_vector(v1, m.subcolinterval(j)), sv, sv);
       }
   }
 
@@ -922,16 +922,16 @@ namespace gmm
     typename sub_vector_type<V3 *, sub_interval>::vector_type sv;
     for (size_type i = 0; i < m.nrowblocks() ; ++i)
       for (size_type j = 0; j < m.ncolblocks() ; ++j) {
-       sv = sub_vector(v3, m.subrowinterval(i));
-       if (j == 0)
-         mult(m.block(i,j),
-              sub_vector(v1, m.subcolinterval(j)),
-              sub_vector(v2, m.subrowinterval(i)), sv);
-       else
-         mult(m.block(i,j),
-              sub_vector(v1, m.subcolinterval(j)), sv, sv);
+        sv = sub_vector(v3, m.subrowinterval(i));
+        if (j == 0)
+          mult(m.block(i,j),
+               sub_vector(v1, m.subcolinterval(j)),
+               sub_vector(v2, m.subrowinterval(i)), sv);
+        else
+          mult(m.block(i,j),
+               sub_vector(v1, m.subcolinterval(j)), sv, sv);
       }
-    
+
   }
 
   template <typename MAT, typename V1, typename V2>
@@ -939,15 +939,15 @@ namespace gmm
   { mult(m, v1, linalg_const_cast(v2)); }
 
   template <typename MAT, typename V1, typename V2, typename V3>
-  void mult(const block_matrix<MAT> &m, const V1 &v1, const V2 &v2, 
-           const V3 &v3)
+  void mult(const block_matrix<MAT> &m, const V1 &v1, const V2 &v2,
+            const V3 &v3)
   { mult_const(m, v1, v2, linalg_const_cast(v3)); }
 
 }
   /* ******************************************************************** */
-  /*                                                                     */
-  /*           Distributed matrices                                      */
-  /*                                                                     */
+  /*                                                                      */
+  /*             Distributed matrices                                     */
+  /*                                                                      */
   /* ******************************************************************** */
 
 #ifdef GMM_USES_MPI
@@ -955,8 +955,8 @@ namespace gmm
 
 namespace gmm {
 
-  
-  
+
+
   template <typename T> inline MPI_Datatype mpi_type(T)
   { GMM_ASSERT1(false, "Sorry unsupported type"); return MPI_FLOAT; }
   inline MPI_Datatype mpi_type(double) { return MPI_DOUBLE; }
@@ -980,7 +980,7 @@ namespace gmm {
     const MAT &local_matrix(void) const { return M; }
     MAT &local_matrix(void) { return M; }
   };
-  
+
   template <typename MAT> inline MAT &eff_matrix(MAT &m) { return m; }
   template <typename MAT> inline
   const MAT &eff_matrix(const MAT &m) { return m; }
@@ -988,29 +988,29 @@ namespace gmm {
   MAT &eff_matrix(mpi_distributed_matrix<MAT> &m) { return m.M; }
   template <typename MAT> inline
   const MAT &eff_matrix(const mpi_distributed_matrix<MAT> &m) { return m.M; }
-  
+
 
   template <typename MAT1, typename MAT2>
   inline void copy(const mpi_distributed_matrix<MAT1> &m1,
-                  mpi_distributed_matrix<MAT2> &m2)
+                   mpi_distributed_matrix<MAT2> &m2)
   { copy(eff_matrix(m1), eff_matrix(m2)); }
   template <typename MAT1, typename MAT2>
   inline void copy(const mpi_distributed_matrix<MAT1> &m1,
-                  const mpi_distributed_matrix<MAT2> &m2)
+                   const mpi_distributed_matrix<MAT2> &m2)
   { copy(m1.M, m2.M); }
-  
+
   template <typename MAT1, typename MAT2>
   inline void copy(const mpi_distributed_matrix<MAT1> &m1, MAT2 &m2)
   { copy(m1.M, m2); }
   template <typename MAT1, typename MAT2>
   inline void copy(const mpi_distributed_matrix<MAT1> &m1, const MAT2 &m2)
   { copy(m1.M, m2); }
-  
+
 
   template <typename MATSP, typename V1, typename V2> inline
   typename strongest_value_type3<V1,V2,MATSP>::value_type
   vect_sp(const mpi_distributed_matrix<MATSP> &ps, const V1 &v1,
-         const V2 &v2) {
+          const V2 &v2) {
     typedef typename strongest_value_type3<V1,V2,MATSP>::value_type T;
     T res = vect_sp(ps.M, v1, v2), rest;
     MPI_Allreduce(&res, &rest, 1, mpi_type(T()), MPI_SUM,MPI_COMM_WORLD);
@@ -1019,7 +1019,7 @@ namespace gmm {
 
   template <typename MAT, typename V1, typename V2>
   inline void mult_add(const mpi_distributed_matrix<MAT> &m, const V1 &v1,
-                      V2 &v2) {
+                       V2 &v2) {
     typedef typename linalg_traits<V2>::value_type T;
     std::vector<T> v3(vect_size(v2)), v4(vect_size(v2));
     static double tmult_tot = 0.0;
@@ -1029,7 +1029,7 @@ namespace gmm {
     if (is_sparse(v2)) GMM_WARNING2("Using a plain temporary, here.");
     double t_ref2 = MPI_Wtime();
     MPI_Allreduce(&(v3[0]), &(v4[0]),gmm::vect_size(v2), mpi_type(T()),
-                 MPI_SUM,MPI_COMM_WORLD);
+                  MPI_SUM,MPI_COMM_WORLD);
     tmult_tot2 = MPI_Wtime()-t_ref2;
     cout << "reduce mult mpi = " << tmult_tot2 << endl;
     gmm::add(v4, v2);
@@ -1039,59 +1039,59 @@ namespace gmm {
 
   template <typename MAT, typename V1, typename V2>
   void mult_add(const mpi_distributed_matrix<MAT> &m, const V1 &v1,
-               const V2 &v2_)
+                const V2 &v2_)
   { mult_add(m, v1, const_cast<V2 &>(v2_)); }
 
   template <typename MAT, typename V1, typename V2>
   inline void mult(const mpi_distributed_matrix<MAT> &m, const V1 &v1,
-                  const V2 &v2_)
+                   const V2 &v2_)
   { V2 &v2 = const_cast<V2 &>(v2_); clear(v2); mult_add(m, v1, v2); }
 
   template <typename MAT, typename V1, typename V2>
   inline void mult(const mpi_distributed_matrix<MAT> &m, const V1 &v1,
-                  V2 &v2)
+                   V2 &v2)
   { clear(v2); mult_add(m, v1, v2); }
 
   template <typename MAT, typename V1, typename V2, typename V3>
   inline void mult(const mpi_distributed_matrix<MAT> &m, const V1 &v1,
-                  const V2 &v2, const V3 &v3_)
+                   const V2 &v2, const V3 &v3_)
   { V3 &v3 = const_cast<V3 &>(v3_); gmm::copy(v2, v3); mult_add(m, v1, v3); }
 
   template <typename MAT, typename V1, typename V2, typename V3>
   inline void mult(const mpi_distributed_matrix<MAT> &m, const V1 &v1,
-                  const V2 &v2, V3 &v3)
+                   const V2 &v2, V3 &v3)
   { gmm::copy(v2, v3); mult_add(m, v1, v3); }
-  
+
 
   template <typename MAT> inline
-  size_type mat_nrows(const mpi_distributed_matrix<MAT> &M) 
+  size_type mat_nrows(const mpi_distributed_matrix<MAT> &M)
   { return mat_nrows(M.M); }
   template <typename MAT> inline
-  size_type mat_ncols(const mpi_distributed_matrix<MAT> &M) 
+  size_type mat_ncols(const mpi_distributed_matrix<MAT> &M)
   { return mat_nrows(M.M); }
   template <typename MAT> inline
   void resize(mpi_distributed_matrix<MAT> &M, size_type m, size_type n)
   { resize(M.M, m, n); }
   template <typename MAT> inline void clear(mpi_distributed_matrix<MAT> &M)
   { clear(M.M); }
-  
+
 
   // For compute reduced system
   template <typename MAT1, typename MAT2> inline
   void mult(const MAT1 &M1, const mpi_distributed_matrix<MAT2> &M2,
-           mpi_distributed_matrix<MAT2> &M3)
+            mpi_distributed_matrix<MAT2> &M3)
   { mult(M1, M2.M, M3.M); }
   template <typename MAT1, typename MAT2> inline
   void mult(const mpi_distributed_matrix<MAT2> &M2,
-           const MAT1 &M1, mpi_distributed_matrix<MAT2> &M3)
+            const MAT1 &M1, mpi_distributed_matrix<MAT2> &M3)
   { mult(M2.M, M1, M3.M); }
   template <typename MAT1, typename MAT2, typename MAT3> inline
   void mult(const MAT1 &M1, const mpi_distributed_matrix<MAT2> &M2,
-                  MAT3 &M3)
+                   MAT3 &M3)
   { mult(M1, M2.M, M3); }
   template <typename MAT1, typename MAT2, typename MAT3> inline
   void mult(const MAT1 &M1, const mpi_distributed_matrix<MAT2> &M2,
-                  const MAT3 &M3)
+                   const MAT3 &M3)
   { mult(M1, M2.M, M3); }
 
   template <typename M, typename SUBI1, typename SUBI2>
@@ -1112,16 +1112,16 @@ namespace gmm {
   template <typename MAT, typename SUBI1, typename SUBI2>  inline
   typename select_return<typename sub_matrix_type<const MAT *, SUBI1, SUBI2>
   ::matrix_type, typename sub_matrix_type<MAT *, SUBI1, SUBI2>::matrix_type,
-                        const MAT *>::return_type
+                         const MAT *>::return_type
   sub_matrix(const mpi_distributed_matrix<MAT> &m, const SUBI1 &si1,
-            const SUBI2 &si2)
+             const SUBI2 &si2)
   { return sub_matrix(m.M, si1, si2);  }
 
   template <typename M, typename SUBI1>  inline
     typename select_return<typename sub_matrix_type<const M *, SUBI1, SUBI1>
     ::matrix_type, typename sub_matrix_type<M *, SUBI1, SUBI1>::matrix_type,
     M *>::return_type
-  sub_matrix(mpi_distributed_matrix<M> &m, const SUBI1 &si1) 
+  sub_matrix(mpi_distributed_matrix<M> &m, const SUBI1 &si1)
   { return sub_matrix(m.M, si1, si1); }
 
   template <typename M, typename SUBI1>  inline
@@ -1132,11 +1132,11 @@ namespace gmm {
   { return sub_matrix(m.M, si1, si1); }
 
 
-  template <typename L> struct transposed_return<const 
mpi_distributed_matrix<L> *> 
+  template <typename L> struct transposed_return<const 
mpi_distributed_matrix<L> *>
   { typedef abstract_null_type return_type; };
-  template <typename L> struct transposed_return<mpi_distributed_matrix<L> *> 
+  template <typename L> struct transposed_return<mpi_distributed_matrix<L> *>
   { typedef abstract_null_type return_type; };
-  
+
   template <typename L> inline typename transposed_return<const L 
*>::return_type
   transposed(const mpi_distributed_matrix<L> &l)
   { return transposed(l.M); }
@@ -1185,10 +1185,10 @@ namespace std {
   template <typename T>
   void swap(gmm::dense_matrix<T> &m1, gmm::dense_matrix<T> &m2)
   { m1.swap(m2); }
-  template <typename T, int shift> void 
+  template <typename T, int shift> void
   swap(gmm::csc_matrix<T,shift> &m1, gmm::csc_matrix<T,shift> &m2)
   { m1.swap(m2); }
-  template <typename T, int shift> void 
+  template <typename T, int shift> void
   swap(gmm::csr_matrix<T,shift> &m1, gmm::csr_matrix<T,shift> &m2)
   { m1.swap(m2); }
 }



reply via email to

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