[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: |
Mon, 25 Mar 2019 09:16:23 -0400 (EDT) |
branch: master
commit 4dec27b98690f70da0f20ca0ae59a80eb5e0431a
Author: Konstantinos Poulios <address@hidden>
Date: Mon Mar 25 14:16:13 2019 +0100
Handle temporary dof counting for unreduced variables in ga_workspace
instead of ga_instruction_set
---
src/getfem/getfem_generic_assembly.h | 44 +++++--
.../getfem_generic_assembly_compile_and_exec.h | 4 +-
src/getfem_generic_assembly_compile_and_exec.cc | 62 ++++-----
src/getfem_generic_assembly_workspace.cc | 143 +++++++++++----------
4 files changed, 135 insertions(+), 118 deletions(-)
diff --git a/src/getfem/getfem_generic_assembly.h
b/src/getfem/getfem_generic_assembly.h
index fc84c17..ac36a51 100644
--- a/src/getfem/getfem_generic_assembly.h
+++ b/src/getfem/getfem_generic_assembly.h
@@ -265,6 +265,7 @@ namespace getfem {
const model *md;
const ga_workspace *parent_workspace;
bool enable_all_md_variables;
+ size_type nb_prim_dof, nb_tmp_dof;
void init();
@@ -369,16 +370,9 @@ namespace getfem {
base_tensor assemb_t;
bool include_empty_int_pts = false;
- public:
+ std::map<std::string, gmm::sub_interval> tmp_var_intervals;
- const model_real_sparse_matrix &assembled_matrix() const { return *K;}
- model_real_sparse_matrix &assembled_matrix() { return *K; }
- scalar_type &assembled_potential()
- { GMM_ASSERT1(assemb_t.size() == 1, "Bad result size"); return
assemb_t[0]; }
- const scalar_type &assembled_potential() const
- { GMM_ASSERT1(assemb_t.size() == 1, "Bad result size"); return
assemb_t[0]; }
- const base_vector &assembled_vector() const { return *V; }
- base_vector &assembled_vector() { return *V; }
+ public:
// setter functions
void set_assembled_matrix(model_real_sparse_matrix &K_) {
K = std::shared_ptr<model_real_sparse_matrix>
@@ -389,9 +383,20 @@ namespace getfem {
(std::shared_ptr<base_vector>(), &V_); // alias
}
// getter functions
+ const model_real_sparse_matrix &assembled_matrix() const { return *K; }
+ model_real_sparse_matrix &assembled_matrix() { return *K; }
+ const base_vector &assembled_vector() const { return *V; }
+ base_vector &assembled_vector() { return *V; }
const base_tensor &assembled_tensor() const { return assemb_t; }
base_tensor &assembled_tensor() { return assemb_t; }
-
+ const scalar_type &assembled_potential() const {
+ GMM_ASSERT1(assemb_t.size() == 1, "Bad result size");
+ return assemb_t[0];
+ }
+ scalar_type &assembled_potential() {
+ GMM_ASSERT1(assemb_t.size() == 1, "Bad result size");
+ return assemb_t[0];
+ }
model_real_sparse_matrix &unreduced_matrix()
{ return unreduced_K; }
base_vector &unreduced_vector() { return unreduced_V; }
@@ -530,8 +535,6 @@ namespace getfem {
psecondary_domain secondary_domain(const std::string &name) const;
-
-
// extract terms
std::string extract_constant_term(const mesh &m);
std::string extract_order1_term(const std::string &varname);
@@ -544,6 +547,23 @@ namespace getfem {
void set_include_empty_int_points(bool include);
bool include_empty_int_points() const;
+ size_type nb_primary_dof() const { return nb_prim_dof; }
+ size_type nb_temporary_dof() const { return nb_tmp_dof; }
+
+ void add_temporary_interval_for_unreduced_variable(const std::string
&name);
+
+ void clear_temporary_variable_intervals() {
+ tmp_var_intervals.clear();
+ nb_tmp_dof = 0;
+ }
+
+ const gmm::sub_interval &
+ temporary_interval_of_variable(const std::string &name) const {
+ static const gmm::sub_interval empty_interval;
+ const auto it = tmp_var_intervals.find(name);
+ return (it != tmp_var_intervals.end()) ? it->second : empty_interval;
+ }
+
ga_workspace(const getfem::model &md_, bool enable_all_variables = false);
ga_workspace(bool, const ga_workspace &gaw);
ga_workspace();
diff --git a/src/getfem/getfem_generic_assembly_compile_and_exec.h
b/src/getfem/getfem_generic_assembly_compile_and_exec.h
index 888b2ce..b60bf0c 100644
--- a/src/getfem/getfem_generic_assembly_compile_and_exec.h
+++ b/src/getfem/getfem_generic_assembly_compile_and_exec.h
@@ -106,8 +106,6 @@ namespace getfem {
std::map<std::string, const base_vector *> extended_vars;
std::map<std::string, base_vector> really_extended_vars;
- std::map<std::string, gmm::sub_interval> var_intervals;
- size_type nb_dof, max_dof;
struct variable_group_info {
const mesh_fem *mf;
@@ -200,7 +198,7 @@ namespace getfem {
std::map<region_mim, region_mim_instructions> all_instructions;
- ga_instruction_set() { max_dof = nb_dof = 0; need_elt_size = false; ipt=0;
}
+ ga_instruction_set() : need_elt_size(false), nbpt(0), ipt(0) {}
};
diff --git a/src/getfem_generic_assembly_compile_and_exec.cc
b/src/getfem_generic_assembly_compile_and_exec.cc
index b25184c..fb14ee5 100644
--- a/src/getfem_generic_assembly_compile_and_exec.cc
+++ b/src/getfem_generic_assembly_compile_and_exec.cc
@@ -1228,16 +1228,13 @@ namespace getfem {
= inin.m ? workspace.variable_in_group(gname, *(inin.m))
: workspace.first_variable_of_group(gname);
vgi.mf = workspace.associated_mf(varname);
- const auto it1 = gis.var_intervals.find(varname);
- GA_DEBUG_ASSERT(it1 != gis.var_intervals.end(),
- "Variable " << varname << " not in gis variables");
- vgi.Ir = it1->second;
+ vgi.Ir = workspace.temporary_interval_of_variable(varname);
vgi.In = workspace.interval_of_variable(varname);
vgi.alpha = workspace.factor_of_variable(varname);
- const auto it2 = gis.extended_vars.find(varname);
- GA_DEBUG_ASSERT(it2 != gis.extended_vars.end(),
+ const auto it = gis.extended_vars.find(varname);
+ GA_DEBUG_ASSERT(it != gis.extended_vars.end(),
"Variable " << varname << " not in extended variables");
- vgi.U = it2->second;
+ vgi.U = it->second;
vgi.varname = &varname;
return 0;
}
@@ -4720,25 +4717,6 @@ namespace getfem {
// Compilation of assembly trees into a list of basic instructions
//=========================================================================
- static void add_interval_to_gis(const ga_workspace &workspace,
- const std::string &varname,
- ga_instruction_set &gis) {
- if (workspace.variable_group_exists(varname)) {
- for (const std::string &v : workspace.variable_group(varname))
- add_interval_to_gis(workspace, v, gis);
- } else {
- if (gis.var_intervals.count(varname) == 0) {
- 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);
- gis.nb_dof += nd;
- }
- gis.max_dof = std::max(gis.max_dof,
- workspace.interval_of_variable(varname).last());
- }
- }
-
static void extend_variable_in_gis(const ga_workspace &workspace,
const std::string &varname,
ga_instruction_set &gis) {
@@ -4771,8 +4749,10 @@ namespace getfem {
ga_clear_node_list(pnode->children[i], node_list);
}
+ // workspace argument is not const because of declaration of temporary
+ // unreduced variables
static void ga_compile_node(const pga_tree_node pnode,
- const ga_workspace &workspace,
+ ga_workspace &workspace,
ga_instruction_set &gis,
ga_instruction_set::region_mim_instructions &rmi,
const mesh &m, bool function_case,
@@ -5796,7 +5776,7 @@ namespace getfem {
}
if (pgai) rmi.instructions.push_back(std::move(pgai));
}
- add_interval_to_gis(workspace, pnode->name, gis);
+ workspace.add_temporary_interval_for_unreduced_variable(pnode->name);
}
break;
@@ -5910,7 +5890,7 @@ namespace getfem {
}
if (pgai) rmi.instructions.push_back(std::move(pgai));
}
- add_interval_to_gis(workspace, pnode->name, gis);
+ workspace.add_temporary_interval_for_unreduced_variable(pnode->name);
}
break;
@@ -5952,7 +5932,7 @@ namespace getfem {
rmi.interpolate_infos[intn], gis.fp_pool);
}
rmi.instructions.push_back(std::move(pgai));
- add_interval_to_gis(workspace, pnode->name, gis);
+ workspace.add_temporary_interval_for_unreduced_variable(pnode->name);
}
break;
@@ -6807,7 +6787,8 @@ namespace getfem {
= workspace.associated_mf(root->name_test1);
const im_data *imd
= workspace.associated_im_data(root->name_test1);
- add_interval_to_gis(workspace, root->name_test1, gis);
+ workspace.add_temporary_interval_for_unreduced_variable
+ (root->name_test1);
if (mf) {
const mesh_fem **mfg = 0;
@@ -6825,7 +6806,8 @@ namespace getfem {
mfg = &(vgi.mf);
mf = 0;
} else {
- Ir = &(gis.var_intervals[root->name_test1]);
+ Ir = &(workspace.temporary_interval_of_variable
+ (root->name_test1));
In = &(workspace.interval_of_variable(root->name_test1));
}
fem_interpolation_context
@@ -6883,8 +6865,10 @@ namespace getfem {
!(intn2.empty() || intn2 ==
"neighbour_elt"
|| secondary2);
- add_interval_to_gis(workspace, root->name_test1, gis);
- add_interval_to_gis(workspace, root->name_test2, gis);
+ workspace.add_temporary_interval_for_unreduced_variable
+ (root->name_test1);
+ workspace.add_temporary_interval_for_unreduced_variable
+ (root->name_test2);
const gmm::sub_interval *Ir1 = 0, *In1 = 0, *Ir2 = 0, *In2=0;
const scalar_type *alpha1 = 0, *alpha2 = 0;
@@ -6901,7 +6885,8 @@ namespace getfem {
alpha1 = &(vgi.alpha);
} else {
alpha1 = &(workspace.factor_of_variable(root->name_test1));
- Ir1 = &(gis.var_intervals[root->name_test1]);
+ Ir1 = &(workspace.temporary_interval_of_variable
+ (root->name_test1));
In1 = &(workspace.interval_of_variable(root->name_test1));
}
@@ -6917,7 +6902,8 @@ namespace getfem {
alpha2 = &(vgi.alpha);
} else {
alpha2 = &(workspace.factor_of_variable(root->name_test2));
- Ir2 = &(gis.var_intervals[root->name_test2]);
+ Ir2 = &(workspace.temporary_interval_of_variable
+ (root->name_test2));
In2 = &(workspace.interval_of_variable(root->name_test2));
}
@@ -6946,8 +6932,8 @@ namespace getfem {
else
pgai = std::make_shared
<ga_instruction_matrix_assembly_standard_vector>
- (root->tensor(),
workspace.assembled_matrix(),ctx1,ctx2,
- *In1, *In2, mf1, mf2,
+ (root->tensor(), workspace.assembled_matrix(),
+ ctx1, ctx2, *In1, *In2, mf1, mf2,
gis.coeff, *alpha1, *alpha2, gis.nbpt, gis.ipt);
} else {
diff --git a/src/getfem_generic_assembly_workspace.cc
b/src/getfem_generic_assembly_workspace.cc
index 6ad3bc0..f7b411a 100644
--- a/src/getfem_generic_assembly_workspace.cc
+++ b/src/getfem_generic_assembly_workspace.cc
@@ -45,18 +45,21 @@ namespace getfem {
void ga_workspace::add_fem_variable
(const std::string &name, const mesh_fem &mf,
const gmm::sub_interval &I, const model_real_plain_vector &VV) {
+ nb_prim_dof = std::max(nb_prim_dof, I.last());
variables.emplace(name, var_description(true, &mf, 0, I, &VV, 1));
}
void ga_workspace::add_im_variable
(const std::string &name, const im_data &imd,
const gmm::sub_interval &I, const model_real_plain_vector &VV) {
+ nb_prim_dof = std::max(nb_prim_dof, I.last());
variables.emplace(name, var_description(true, 0, &imd, I, &VV, 1));
}
void ga_workspace::add_fixed_size_variable
(const std::string &name,
const gmm::sub_interval &I, const model_real_plain_vector &VV) {
+ nb_prim_dof = std::max(nb_prim_dof, I.last());
variables.emplace(name, var_description(true, 0, 0, I, &VV,
dim_type(gmm::vect_size(VV))));
}
@@ -743,33 +746,30 @@ namespace getfem {
void ga_workspace::assembly(size_type order) {
- size_type ndof;
const ga_workspace *w = this;
while (w->parent_workspace) w = w->parent_workspace;
- if (w->md) ndof = w->md->nb_dof(); // To eventually call actualize_sizes()
+ if (w->md) w->md->nb_dof(); // To eventually call actualize_sizes()
GA_TIC;
ga_instruction_set gis;
ga_compile(*this, gis, order);
- ndof = gis.nb_dof;
- size_type max_dof = gis.max_dof;
GA_TOCTIC("Compile time");
if (order == 2) {
if (K.use_count()) {
gmm::clear(*K);
- gmm::resize(*K, max_dof, max_dof);
+ gmm::resize(*K, nb_prim_dof, nb_prim_dof);
}
gmm::clear(unreduced_K);
- gmm::resize(unreduced_K, ndof, ndof);
+ gmm::resize(unreduced_K, nb_tmp_dof, nb_tmp_dof);
}
if (order == 1) {
if (V.use_count()) {
gmm::clear(*V);
- gmm::resize(*V, max_dof);
+ gmm::resize(*V, nb_prim_dof);
}
gmm::clear(unreduced_V);
- gmm::resize(unreduced_V, ndof);
+ gmm::resize(unreduced_V, nb_tmp_dof);
}
gmm::clear(assembled_tensor().as_vector());
GA_TOCTIC("Init time");
@@ -790,68 +790,60 @@ namespace getfem {
std::set<std::pair<std::string, std::string> > vars_mat_done;
for (ga_tree &tree : gis.trees) {
if (tree.root) {
+ std::string &name1 = tree.root->name_test1;
+ std::string &name2 = tree.root->name_test2;
+ const std::vector<std::string> vnames1_(1,name1),
+ vnames2_(1,name2);
+ const std::vector<std::string>
+ &vnames1 = variable_group_exists(name1) ? variable_group(name1)
+ : vnames1_,
+ &vnames2 = variable_group_exists(name2) ? variable_group(name2)
+ : vnames2_;
if (order == 1) {
- const std::string &name = tree.root->name_test1;
- const std::vector<std::string> vnames_(1,name);
- const std::vector<std::string> &vnames
- = variable_group_exists(name) ? variable_group(name)
- : vnames_;
- for (const std::string &vname : vnames) {
- const mesh_fem *mf = associated_mf(vname);
- if (mf && mf->is_reduced() &&
- vars_vec_done.find(vname) == vars_vec_done.end()) {
- gmm::mult_add(gmm::transposed(mf->extension_matrix()),
- gmm::sub_vector(unreduced_V,
- gis.var_intervals[vname]),
- gmm::sub_vector(*V,
- interval_of_variable(vname)));
- vars_vec_done.insert(vname);
+ for (const std::string &vname1 : vnames1) {
+ const mesh_fem *mf1 = associated_mf(vname1);
+ if (mf1 && mf1->is_reduced() && vars_vec_done.count(vname1) == 0)
+ {
+ gmm::sub_interval uI1 = temporary_interval_of_variable(vname1),
+ I1 = interval_of_variable(vname1);
+ gmm::mult_add(gmm::transposed(mf1->extension_matrix()),
+ gmm::sub_vector(unreduced_V, uI1),
+ gmm::sub_vector(*V, I1));
+ vars_vec_done.insert(vname1);
}
}
} else {
- std::string &name1 = tree.root->name_test1;
- std::string &name2 = tree.root->name_test2;
- const std::vector<std::string> vnames1_(1,name1),
- vnames2_(1,name2);
- const std::vector<std::string> &vnames1
- = variable_group_exists(name1) ? variable_group(name1)
- : vnames1_;
- const std::vector<std::string> &vnames2
- = variable_group_exists(name2) ? variable_group(name2)
- : vnames2_;
for (const std::string &vname1 : vnames1) {
for (const std::string &vname2 : vnames2) {
- const mesh_fem *mf1 = associated_mf(vname1);
- const mesh_fem *mf2 = associated_mf(vname2);
- if (((mf1 && mf1->is_reduced())
- || (mf2 && mf2->is_reduced()))) {
- std::pair<std::string, std::string> p(vname1, vname2);
- if (vars_mat_done.find(p) == vars_mat_done.end()) {
- gmm::sub_interval uI1 = gis.var_intervals[vname1];
- gmm::sub_interval uI2 = gis.var_intervals[vname2];
- gmm::sub_interval I1 = interval_of_variable(vname1);
- gmm::sub_interval I2 = interval_of_variable(vname2);
- if ((mf1 && mf1->is_reduced()) &&
- (mf2 && mf2->is_reduced())) {
- model_real_sparse_matrix aux(I1.size(), uI2.size());
- model_real_row_sparse_matrix M(I1.size(), I2.size());
- gmm::mult(gmm::transposed(mf1->extension_matrix()),
- gmm::sub_matrix(unreduced_K, uI1, uI2), aux);
- gmm::mult(aux, mf2->extension_matrix(), M);
- gmm::add(M, gmm::sub_matrix(*K, I1, I2));
- } else if (mf1 && mf1->is_reduced()) {
- model_real_sparse_matrix M(I1.size(), I2.size());
- gmm::mult(gmm::transposed(mf1->extension_matrix()),
- gmm::sub_matrix(unreduced_K, uI1, uI2), M);
- gmm::add(M, gmm::sub_matrix(*K, I1, I2));
- } else {
- model_real_row_sparse_matrix M(I1.size(), I2.size());
- gmm::mult(gmm::sub_matrix(unreduced_K, uI1, uI2),
- mf2->extension_matrix(), M);
- gmm::add(M, gmm::sub_matrix(*K, I1, I2));
- }
- vars_mat_done.insert(p);
+ const mesh_fem *mf1 = associated_mf(vname1),
+ *mf2 = associated_mf(vname2);
+ if (((mf1 && mf1->is_reduced()) || (mf2 && mf2->is_reduced()))
+ && vars_mat_done.count(std::make_pair(vname1,vname2)) == 0)
+ {
+ gmm::sub_interval
+ uI1 = temporary_interval_of_variable(vname1),
+ uI2 = temporary_interval_of_variable(vname2),
+ I1 = interval_of_variable(vname1),
+ I2 = interval_of_variable(vname2);
+ if (mf1 && mf1->is_reduced() && mf2 && mf2->is_reduced()) {
+ model_real_sparse_matrix aux(I1.size(), uI2.size());
+ model_real_row_sparse_matrix M(I1.size(), I2.size());
+ gmm::mult(gmm::transposed(mf1->extension_matrix()),
+ gmm::sub_matrix(unreduced_K, uI1, uI2), aux);
+ gmm::mult(aux, mf2->extension_matrix(), M);
+ gmm::add(M, gmm::sub_matrix(*K, I1, I2));
+ } else if (mf1 && mf1->is_reduced()) {
+ model_real_sparse_matrix M(I1.size(), I2.size());
+ gmm::mult(gmm::transposed(mf1->extension_matrix()),
+ gmm::sub_matrix(unreduced_K, uI1, uI2), M);
+ gmm::add(M, gmm::sub_matrix(*K, I1, I2));
+ } else {
+ model_real_row_sparse_matrix M(I1.size(), I2.size());
+ gmm::mult(gmm::sub_matrix(unreduced_K, uI1, uI2),
+ mf2->extension_matrix(), M);
+ gmm::add(M, gmm::sub_matrix(*K, I1, I2));
}
+ vars_mat_done.insert(std::make_pair(vname1,vname2));
}
}
}
@@ -869,6 +861,21 @@ namespace getfem {
return include_empty_int_pts;
}
+ void ga_workspace::add_temporary_interval_for_unreduced_variable
+ (const std::string &name)
+ {
+ if (variable_group_exists(name)) {
+ for (const std::string &v : variable_group(name))
+ add_temporary_interval_for_unreduced_variable(v);
+ } else if (tmp_var_intervals.count(name) == 0) {
+ const mesh_fem *mf = associated_mf(name);
+ size_type nd = mf ? mf->nb_basic_dof()
+ : gmm::vect_size(value(name));
+ tmp_var_intervals[name] = gmm::sub_interval(nb_tmp_dof, nd);
+ nb_tmp_dof += nd;
+ }
+ }
+
void ga_workspace::clear_expressions() { trees.clear(); }
void ga_workspace::print(std::ostream &str) {
@@ -906,14 +913,20 @@ namespace getfem {
bool enable_all_variables)
: md(&md_), parent_workspace(0),
enable_all_md_variables(enable_all_variables),
+ nb_prim_dof(0), nb_tmp_dof(0),
macro_dict(md_.macro_dictionary())
- { init(); }
+ {
+ init();
+ nb_prim_dof = md->nb_dof();
+ }
ga_workspace::ga_workspace(bool, const ga_workspace &gaw)
: md(0), parent_workspace(&gaw), enable_all_md_variables(false),
+ nb_prim_dof(gaw.nb_primary_dof()), nb_tmp_dof(0),
macro_dict(gaw.macro_dictionary())
{ init(); }
ga_workspace::ga_workspace()
- : md(0), parent_workspace(0), enable_all_md_variables(false)
+ : md(0), parent_workspace(0), enable_all_md_variables(false),
+ nb_prim_dof(0), nb_tmp_dof(0)
{ init(); }
ga_workspace::~ga_workspace() { clear_expressions(); }