[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: |
Wed, 1 Jan 2020 05:45:46 -0500 (EST) |
branch: master
commit cac4df66312b3b6155eab78c794b820a5e25ff57
Author: Konstantinos Poulios <address@hidden>
Date: Wed Jan 1 11:45:37 2020 +0100
Code readability and typo fixes
---
src/getfem/getfem_generic_assembly.h | 14 +-
src/getfem/getfem_models.h | 4 +-
src/getfem_generic_assembly_compile_and_exec.cc | 8 +-
src/getfem_generic_assembly_workspace.cc | 5 +
src/getfem_models.cc | 218 ++++++++++++------------
src/getfem_nonlinear_elasticity.cc | 2 +-
6 files changed, 131 insertions(+), 120 deletions(-)
diff --git a/src/getfem/getfem_generic_assembly.h
b/src/getfem/getfem_generic_assembly.h
index 68ea575..8d1f8d4 100644
--- a/src/getfem/getfem_generic_assembly.h
+++ b/src/getfem/getfem_generic_assembly.h
@@ -288,10 +288,10 @@ namespace getfem {
}
var_description(bool is_var, const mesh_fem *mf_, const im_data *imd_,
- gmm::sub_interval I_, const model_real_plain_vector *v,
+ gmm::sub_interval I_, const model_real_plain_vector *V_,
size_type Q)
: is_variable(is_var), is_fem_dofs(mf_ != 0), mf(mf_), imd(imd_),
- I(I_), V(v), qdims(1)
+ I(I_), V(V_), qdims(1)
{
GMM_ASSERT1(Q > 0, "Bad dimension");
qdims[0] = Q;
@@ -340,10 +340,12 @@ namespace getfem {
const mesh_region ®ister_region(const mesh &m, const mesh_region &rg);
// variables and variable groups
- mutable std::map<std::string, gmm::sub_interval> int_disabled_variables;
-
typedef std::map<std::string, var_description> VAR_SET;
VAR_SET variables;
+
+ mutable std::map<std::string, gmm::sub_interval> int_disabled_variables;
+ std::map<std::string, gmm::sub_interval> tmp_var_intervals;
+
std::map<std::string, pinterpolate_transformation> transformations;
std::map<std::string, pelementary_transformation> elem_transformations;
std::map<std::string, psecondary_domain> secondary_domains;
@@ -362,7 +364,6 @@ namespace getfem {
bool scalar_expr, operation_type op_type=ASSEMBLY,
const std::string varname_interpolation="");
-
std::shared_ptr<model_real_sparse_matrix> K;
std::shared_ptr<base_vector> V;
model_real_sparse_matrix col_unreduced_K,
@@ -372,8 +373,6 @@ namespace getfem {
base_tensor assemb_t;
bool include_empty_int_pts = false;
- std::map<std::string, gmm::sub_interval> tmp_var_intervals;
-
public:
// setter functions
void set_assembled_matrix(model_real_sparse_matrix &K_) {
@@ -460,6 +459,7 @@ namespace getfem {
std::vector<std::string> &vl_test2,
std::vector<std::string> &dl,
size_type order);
+ bool is_linear(size_type order);
bool variable_exists(const std::string &name) const;
diff --git a/src/getfem/getfem_models.h b/src/getfem/getfem_models.h
index 0d72462..54b62d7 100644
--- a/src/getfem/getfem_models.h
+++ b/src/getfem/getfem_models.h
@@ -2304,7 +2304,7 @@ namespace getfem {
/** Linear elasticity brick ( @f$ \int \sigma(u):\varepsilon(v) @f$ ).
for isotropic material. Parametrized by Young modulus and Poisson ratio
- For two-dimensional problems, corresponds to the plain strain
+ For two-dimensional problems, corresponds to the plane strain
approximation
( @f$ \lambda = E\nu/((1+\nu)(1-2\nu)), \mu = E/(2(1+\nu)) @f$ ).
Corresponds to the standard model for three-dimensional problems.
@@ -2317,7 +2317,7 @@ namespace getfem {
/**
Linear elasticity brick ( @f$ \int \sigma(u):\varepsilon(v) @f$ ).
for isotropic material. Parametrized by Young modulus and Poisson ratio.
- For two-dimensional problems, corresponds to the plain stress
+ For two-dimensional problems, corresponds to the plane stress
approximation
( @f$ \lambda^* = E\nu/(1-\nu^2), \mu = E/(2(1+\nu)) @f$ ).
Corresponds to the standard model for three-dimensional problems.
diff --git a/src/getfem_generic_assembly_compile_and_exec.cc
b/src/getfem_generic_assembly_compile_and_exec.cc
index dcf0921..f3d4ac9 100644
--- a/src/getfem_generic_assembly_compile_and_exec.cc
+++ b/src/getfem_generic_assembly_compile_and_exec.cc
@@ -5174,7 +5174,7 @@ namespace getfem {
(is_elementary ? pnode->elementary_target :
pnode->name)
<< " has to be defined on the same mesh than the "
<< "integration method or interpolation used");
-
+
// An instruction for extracting local dofs of the variable.
if (rmi.local_dofs.count(pnode->name) == 0) {
rmi.local_dofs[pnode->name] = base_vector(1);
@@ -5206,7 +5206,7 @@ namespace getfem {
(*mf, rmi.pfps[mf], gis.ctx, gis.fp_pool);
rmi.instructions.push_back(std::move(pgai));
}
-
+
// An instruction for the base value
pgai = pga_instruction();
switch (pnode->node_type) {
@@ -5289,7 +5289,7 @@ namespace getfem {
(rmi.xfem_minus_hess[mf], gis.ctx, *mf, rmi.pfps[mf]);
}
break;
-
+
default : GMM_ASSERT1(false, "Internal error");
}
if (pgai) rmi.instructions.push_back(std::move(pgai));
@@ -5419,7 +5419,7 @@ namespace getfem {
}
}
break;
-
+
case GA_NODE_SECONDARY_DOMAIN_VAL: case GA_NODE_SECONDARY_DOMAIN_GRAD:
case GA_NODE_SECONDARY_DOMAIN_HESS: case GA_NODE_SECONDARY_DOMAIN_DIVERG:
{
diff --git a/src/getfem_generic_assembly_workspace.cc
b/src/getfem_generic_assembly_workspace.cc
index 1430f7d..b8db69d 100644
--- a/src/getfem_generic_assembly_workspace.cc
+++ b/src/getfem_generic_assembly_workspace.cc
@@ -707,6 +707,11 @@ namespace getfem {
return islin;
}
+ bool ga_workspace::is_linear(size_type order) {
+ std::vector<std::string> vl, vl_test1, vl_test2, dl;
+ return used_variables(vl, vl_test1, vl_test2, dl, order);
+ }
+
void ga_workspace::define_variable_group(const std::string &group_name,
const std::vector<std::string> &nl)
{
GMM_ASSERT1(!(variable_exists(group_name)), "The name of a group of "
diff --git a/src/getfem_models.cc b/src/getfem_models.cc
index 8f6a42a..29754ab 100644
--- a/src/getfem_models.cc
+++ b/src/getfem_models.cc
@@ -946,8 +946,8 @@ namespace getfem {
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 (it != it2 && it2->second.is_fem_dofs &&
(it2->second.filter & VDESCRFILTER_INFSUP) &&
mim == it2->second.filter_mim) found = true;
@@ -981,7 +981,7 @@ namespace getfem {
active_bricks.add(ib);
valid_bricks.add(ib);
- // The brick itself already react to a mesh_im change in update_brick()
+ // The brick itself already reacts to a mesh_im change in update_brick()
// for (size_type i = 0; i < bricks[ib].mims.size(); ++i)
// add_dependency(*(bricks[ib].mims[i]));
@@ -2312,31 +2312,35 @@ namespace getfem {
size_type nbgdof = nb_dof();
scalar_type alpha = coeff0, alpha1 = coeff0, alpha2 = coeff0;
gmm::sub_interval I1(0,nbgdof), I2(0,nbgdof);
- VAR_SET::iterator it1, it2;
+ var_description *var1=nullptr, *var2=nullptr;
if (!isg) {
- it1 = variables.find(term.var1);
- GMM_ASSERT1(it1->second.is_variable, "Assembly of data not allowed");
- I1 = it1->second.I;
- }
- if (term.is_matrix_term && !isg) {
- it2 = variables.find(term.var2);
- I2 = it2->second.I;
- if (!(it2->second.is_variable)) {
- std::string vorgname = sup_previous_and_dot_to_varname(term.var2);
- VAR_SET::iterator it3 = variables.find(vorgname);
- GMM_ASSERT1(it3->second.is_variable,
- "Assembly of data not allowed");
- I2 = it3->second.I;
- isprevious = true;
+ VAR_SET::iterator it1 = variables.find(term.var1);
+ GMM_ASSERT1(it1 != variables.end(), "Internal error");
+ var1 = &(it1->second);
+ GMM_ASSERT1(var1->is_variable, "Assembly of data not allowed");
+ I1 = var1->I;
+ if (term.is_matrix_term) {
+ VAR_SET::iterator it2 = variables.find(term.var2);
+ GMM_ASSERT1(it2 != variables.end(), "Internal error");
+ var2 = &(it2->second);
+ I2 = var2->I;
+ if (!(var2->is_variable)) {
+ std::string vorgname =
sup_previous_and_dot_to_varname(term.var2);
+ VAR_SET::iterator it3 = variables.find(vorgname);
+ GMM_ASSERT1(it3->second.is_variable,
+ "Assembly of data not allowed");
+ I2 = it3->second.I;
+ isprevious = true;
+ }
+ alpha *= var1->alpha * var2->alpha;
+ alpha1 *= var1->alpha;
+ alpha2 *= var2->alpha;
}
- alpha *= it1->second.alpha * it2->second.alpha;
- alpha1 *= it1->second.alpha;
- alpha2 *= it2->second.alpha;
}
- if (cplx) {
+ if (cplx) { // complex term in complex model
if (term.is_matrix_term && (version & BUILD_MATRIX) && !isprevious
- && (isg || (!(it1->second.is_disabled) &&
!(it2->second.is_disabled)))) {
+ && (isg || (!(var1->is_disabled) && !(var2->is_disabled)))) {
gmm::add(gmm::scaled(brick.cmatlist[j], alpha),
gmm::sub_matrix(cTM, I1, I2));
if (term.is_symmetric && I1.first() != I2.first()) {
@@ -2345,7 +2349,7 @@ namespace getfem {
}
}
if (version & BUILD_RHS) {
- if (isg || !(it1->second.is_disabled)) {
+ if (isg || !(var1->is_disabled)) {
if (brick.pdispatch) {
for (size_type k = 0; k < brick.nbrhs; ++k)
gmm::add(gmm::scaled(brick.cveclist[k][j],
@@ -2359,31 +2363,31 @@ namespace getfem {
}
}
if (term.is_matrix_term && brick.pbr->is_linear() && is_linear()) {
- if (it2->second.is_affine_dependent
- && !(it1->second.is_disabled))
+ if (var2->is_affine_dependent
+ && !(var1->is_disabled))
gmm::mult_add(brick.cmatlist[j],
- gmm::scaled(it2->second.affine_complex_value,
+ gmm::scaled(var2->affine_complex_value,
complex_type(-alpha1)),
gmm::sub_vector(crhs, I1));
if (term.is_symmetric && I1.first() != I2.first()
- && it1->second.is_affine_dependent
- && !(it2->second.is_disabled)) {
+ && var1->is_affine_dependent
+ && !(var2->is_disabled)) {
gmm::mult_add(gmm::conjugated(brick.cmatlist[j]),
- gmm::scaled(it1->second.affine_complex_value,
+ gmm::scaled(var1->affine_complex_value,
complex_type(-alpha2)),
gmm::sub_vector(crhs, I2));
}
}
if (term.is_matrix_term && brick.pbr->is_linear()
&& (!is_linear() || (version & BUILD_WITH_COMPLETE_RHS))) {
- if (!(it1->second.is_disabled))
+ if (!(var1->is_disabled))
gmm::mult_add(brick.cmatlist[j],
- gmm::scaled(it2->second.complex_value[0],
+ gmm::scaled(var2->complex_value[0],
complex_type(-alpha1)),
gmm::sub_vector(crhs, I1));
}
if (term.is_symmetric && I1.first() != I2.first() &&
- !(it2->second.is_disabled)) {
+ !(var2->is_disabled)) {
if (brick.pdispatch) {
for (size_type k = 0; k < brick.nbrhs; ++k)
gmm::add(gmm::scaled(brick.cveclist_sym[k][j],
@@ -2397,15 +2401,15 @@ namespace getfem {
if (brick.pbr->is_linear()
&& (!is_linear() || (version & BUILD_WITH_COMPLETE_RHS))) {
gmm::mult_add(gmm::conjugated(brick.cmatlist[j]),
- gmm::scaled(it1->second.complex_value[0],
+ gmm::scaled(var1->complex_value[0],
complex_type(-alpha2)),
gmm::sub_vector(crhs, I2));
}
}
}
- } else if (is_complex()) {
+ } else if (is_complex()) { // real term in complex model
if (term.is_matrix_term && (version & BUILD_MATRIX) && !isprevious
- && (isg || (!(it1->second.is_disabled) &&
!(it2->second.is_disabled)))) {
+ && (isg || (!(var1->is_disabled) && !(var2->is_disabled)))) {
gmm::add(gmm::scaled(brick.rmatlist[j], alpha),
gmm::sub_matrix(cTM, I1, I2));
if (term.is_symmetric && I1.first() != I2.first()) {
@@ -2414,7 +2418,7 @@ namespace getfem {
}
}
if (version & BUILD_RHS) {
- if (isg || !(it1->second.is_disabled)) {
+ if (isg || !(var1->is_disabled)) {
if (brick.pdispatch) {
for (size_type k = 0; k < brick.nbrhs; ++k)
gmm::add(gmm::scaled(brick.rveclist[k][j],
@@ -2427,31 +2431,31 @@ namespace getfem {
}
}
if (term.is_matrix_term && brick.pbr->is_linear() && is_linear()) {
- if (it2->second.is_affine_dependent
- && !(it1->second.is_disabled))
+ if (var2->is_affine_dependent
+ && !(var1->is_disabled))
gmm::mult_add(brick.rmatlist[j],
- gmm::scaled(it2->second.affine_complex_value,
+ gmm::scaled(var2->affine_complex_value,
complex_type(-alpha1)),
gmm::sub_vector(crhs, I1));
if (term.is_symmetric && I1.first() != I2.first()
- && it1->second.is_affine_dependent
- && !(it2->second.is_disabled)) {
+ && var1->is_affine_dependent
+ && !(var2->is_disabled)) {
gmm::mult_add(gmm::transposed(brick.rmatlist[j]),
- gmm::scaled(it1->second.affine_complex_value,
+ gmm::scaled(var1->affine_complex_value,
complex_type(-alpha2)),
gmm::sub_vector(crhs, I2));
}
}
if (term.is_matrix_term && brick.pbr->is_linear()
&& (!is_linear() || (version & BUILD_WITH_COMPLETE_RHS))) {
- if (!(it1->second.is_disabled))
+ if (!(var1->is_disabled))
gmm::mult_add(brick.rmatlist[j],
- gmm::scaled(it2->second.complex_value[0],
+ gmm::scaled(var2->complex_value[0],
complex_type(-alpha1)),
gmm::sub_vector(crhs, I1));
}
if (term.is_symmetric && I1.first() != I2.first() &&
- !(it2->second.is_disabled)) {
+ !(var2->is_disabled)) {
if (brick.pdispatch) {
for (size_type k = 0; k < brick.nbrhs; ++k)
gmm::add(gmm::scaled(brick.rveclist_sym[k][j],
@@ -2465,16 +2469,16 @@ 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(var1->complex_value[0],
complex_type(-alpha2)),
gmm::sub_vector(crhs, I2));
}
}
}
- } else {
+ } else { // real term in real model
if (term.is_matrix_term && (version & BUILD_MATRIX) && !isprevious
- && (isg || (!(it1->second.is_disabled)
- && !(it2->second.is_disabled)))) {
+ && (isg || (!(var1->is_disabled)
+ && !(var2->is_disabled)))) {
gmm::add(gmm::scaled(brick.rmatlist[j], alpha),
gmm::sub_matrix(rTM, I1, I2));
if (term.is_symmetric && I1.first() != I2.first()) {
@@ -2483,60 +2487,56 @@ namespace getfem {
}
}
if (version & BUILD_RHS) {
- if (isg || !(it1->second.is_disabled)) {
+ // Contributions to interval I1 of var1
+ if (isg || !(var1->is_disabled)) {
if (brick.pdispatch) {
for (size_type k = 0; k < brick.nbrhs; ++k)
gmm::add(gmm::scaled(brick.rveclist[k][j],
brick.coeffs[k]),
gmm::sub_vector(rrhs, I1));
}
- else {
+ else
gmm::add(gmm::scaled(brick.rveclist[0][j], alpha1),
gmm::sub_vector(rrhs, I1));
- }
}
- if (term.is_matrix_term && brick.pbr->is_linear() && is_linear()) {
- if (it2->second.is_affine_dependent
- && !(it1->second.is_disabled))
+ if (!(var1->is_disabled)) {
+ // Contributions from affine dependent variables
+ if (term.is_matrix_term && brick.pbr->is_linear() && is_linear()
+ && var2->is_affine_dependent)
gmm::mult_add(brick.rmatlist[j],
- gmm::scaled(it2->second.affine_real_value,
- -alpha1),
+ gmm::scaled(var2->affine_real_value, -alpha1),
gmm::sub_vector(rrhs, I1));
- if (term.is_symmetric && I1.first() != I2.first()
- && it1->second.is_affine_dependent
- && !(it2->second.is_disabled)) {
- gmm::mult_add(gmm::transposed(brick.rmatlist[j]),
- gmm::scaled(it1->second.affine_real_value,
- -alpha2),
- gmm::sub_vector(rrhs, I2));
- }
- }
- if (term.is_matrix_term && brick.pbr->is_linear()
- && (!is_linear() || (version & BUILD_WITH_COMPLETE_RHS))) {
- if (!(it1->second.is_disabled))
+ // Contributions from linear terms
+ if (term.is_matrix_term && brick.pbr->is_linear()
+ && (!is_linear() || (version & BUILD_WITH_COMPLETE_RHS)))
gmm::mult_add(brick.rmatlist[j],
- gmm::scaled(it2->second.real_value[0],
- -alpha1),
+ gmm::scaled(var2->real_value[0], -alpha1),
gmm::sub_vector(rrhs, I1));
}
+ // Contributions to interval I2 of var2 due to symmetric terms
if (term.is_symmetric && I1.first() != I2.first() &&
- !(it2->second.is_disabled)) {
+ !(var2->is_disabled)) {
if (brick.pdispatch) {
for (size_type k = 0; k < brick.nbrhs; ++k)
gmm::add(gmm::scaled(brick.rveclist_sym[k][j],
brick.coeffs[k]),
gmm::sub_vector(rrhs, I2));
}
- else {
+ else
gmm::add(gmm::scaled(brick.rveclist_sym[0][j], alpha2),
gmm::sub_vector(rrhs, I2));
- }
+ // Contributions from affine dependent variables
+ if (term.is_matrix_term && brick.pbr->is_linear() && is_linear()
+ && var1->is_affine_dependent)
+ gmm::mult_add(gmm::transposed(brick.rmatlist[j]),
+ gmm::scaled(var1->affine_real_value, -alpha2),
+ gmm::sub_vector(rrhs, I2));
+ // Contributions from linear terms
if (brick.pbr->is_linear()
- && (!is_linear() || (version & BUILD_WITH_COMPLETE_RHS))) {
+ && (!is_linear() || (version & BUILD_WITH_COMPLETE_RHS)))
gmm::mult_add(gmm::transposed(brick.rmatlist[j]),
- gmm::scaled(it1->second.real_value[0], -alpha2),
+ gmm::scaled(var1->real_value[0], -alpha2),
gmm::sub_vector(rrhs, I2));
- }
}
}
}
@@ -3169,16 +3169,18 @@ namespace getfem {
gmm::clear(vecl[0]);
if (expr.size()) {
- size_type nbgdof = md.nb_dof();
ga_workspace workspace(md, true);
GMM_TRACE2(name << ": generic source term assembly");
workspace.add_expression(expr, *(mims[0]), region, 1,
secondary_domain);
- model::varnamelist vlmd; md.variable_list(vlmd);
- for (size_type i = 0; i < vlmd.size(); ++i)
- if (md.is_disabled_variable(vlmd[i]))
- nbgdof = std::max(nbgdof,
- workspace.interval_of_variable(vlmd[i]).last());
- GMM_TRACE2(name << ": generic matrix assembly");
+ size_type nbgdof = md.nb_dof();
+ {
+ model::varnamelist vlmd;
+ md.variable_list(vlmd);
+ for (const auto &varname : vlmd)
+ if (md.is_disabled_variable(varname))
+ nbgdof = std::max(nbgdof,
+
workspace.interval_of_variable(varname).last());
+ }
model_real_plain_vector V(nbgdof);
workspace.set_assembled_vector(V);
workspace.assembly(1);
@@ -4912,8 +4914,7 @@ namespace getfem {
ga_workspace workspace(md, true);
size_type order = workspace.add_expression(Neumannterm, mim, region, 1);
GMM_ASSERT1(order == 0, "Wrong expression of the Neumann term");
- model::varnamelist vl, vl_test1, vl_test2, dl;
- bool is_lin = workspace.used_variables(vl, vl_test1, vl_test2, dl, 1);
+ bool is_lin = workspace.is_linear(1);
std::string condition = "("+varname + (datag.size() ? "-("+datag+"))":")");
std::string gamma = "(("+datagamma0+")*element_size)";
@@ -4943,11 +4944,11 @@ namespace getfem {
const std::string &datagamma0, size_type region, scalar_type theta_,
const std::string &datag) {
std::string theta = std::to_string(theta_);
+
ga_workspace workspace(md, true);
size_type order = workspace.add_expression(Neumannterm, mim, region, 1);
GMM_ASSERT1(order == 0, "Wrong expression of the Neumann term");
- model::varnamelist vl, vl_test1, vl_test2, dl;
- bool is_lin = workspace.used_variables(vl, vl_test1, vl_test2, dl, 1);
+ bool is_lin = workspace.is_linear(1);
std::string condition = "("+varname+".Normal"
+ (datag.size() ? "-("+datag+"))":")");
@@ -4976,11 +4977,11 @@ namespace getfem {
const std::string &datagamma0, size_type region, scalar_type theta_,
const std::string &datag, const std::string &dataH) {
std::string theta = std::to_string(theta_);
+
ga_workspace workspace(md, true);
size_type order = workspace.add_expression(Neumannterm, mim, region, 1);
GMM_ASSERT1(order == 0, "Wrong expression of the Neumann term");
- model::varnamelist vl, vl_test1, vl_test2, dl;
- bool is_lin = workspace.used_variables(vl, vl_test1, vl_test2, dl, 1);
+ bool is_lin = workspace.is_linear(1);
std::string condition = "(("+dataH+")*"+varname
+ (datag.size() ? "-("+datag+"))":")");
@@ -5988,11 +5989,14 @@ namespace getfem {
std::string expr2 = "(Div_"+varname+"*(("+dataexpr1+")*Id(meshdim))"
+"+(2*("+dataexpr2+"))*Sym(Grad_"+varname+")):Grad_"+test_varname;
- ga_workspace workspace(md, true);
- workspace.add_expression(expr2, mim, region);
- model::varnamelist vl, vl_test1, vl_test2, dl;
- bool is_lin = workspace.used_variables(vl, vl_test1, vl_test2, dl, 2);
-
+ bool is_lin;
+ model::varnamelist vl, dl;
+ {
+ ga_workspace workspace(md, true);
+ workspace.add_expression(expr2, mim, region);
+ model::varnamelist vl_test1, vl_test2;
+ is_lin = workspace.used_variables(vl, vl_test1, vl_test2, dl, 2);
+ }
if (is_lin) {
pbrick pbr = std::make_shared<iso_lin_elasticity_new_brick>
(expr2, dataname3);
@@ -6021,11 +6025,12 @@ namespace getfem {
std::string expr = lambda+"*Div_"+varname+"*Div_"+test_varname
+ "+"+mu+"*(Grad_"+varname+"+Grad_"+varname+"'):Grad_"+test_varname;
- ga_workspace workspace(md, true);
- workspace.add_expression(expr, mim, region);
- model::varnamelist vl, vl_test1, vl_test2, dl;
- bool is_lin = workspace.used_variables(vl, vl_test1, vl_test2, dl, 2);
-
+ bool is_lin;
+ {
+ ga_workspace workspace(md, true);
+ workspace.add_expression(expr, mim, region);
+ is_lin = workspace.is_linear(2);
+ }
if (is_lin) {
return add_linear_term(md, mim, expr, region, false, false,
"Linearized isotropic elasticity");
@@ -6055,11 +6060,12 @@ namespace getfem {
std::string expr = lambda+"*Div_"+varname+"*Div_"+test_varname
+ "+"+mu+"*(Grad_"+varname+"+Grad_"+varname+"'):Grad_"+test_varname;
- ga_workspace workspace(md, true);
- workspace.add_expression(expr, mim, region);
- model::varnamelist vl, vl_test1, vl_test2, dl;
- bool is_lin = workspace.used_variables(vl, vl_test1, vl_test2, dl, 2);
-
+ bool is_lin;
+ {
+ ga_workspace workspace(md, true);
+ workspace.add_expression(expr, mim, region);
+ is_lin = workspace.is_linear(2);
+ }
if (is_lin) {
return add_linear_term(md, mim, expr, region, false, false,
"Linearized isotropic elasticity");
diff --git a/src/getfem_nonlinear_elasticity.cc
b/src/getfem_nonlinear_elasticity.cc
index 3542e57..7aeaa6e 100644
--- a/src/getfem_nonlinear_elasticity.cc
+++ b/src/getfem_nonlinear_elasticity.cc
@@ -2313,7 +2313,7 @@ namespace getfem {
"fem variables");
size_type Q = mf->get_qdim();
GMM_ASSERT1(Q == N, "Finite strain elasticity brick can only be applied "
- "on a fem variable having the same dimension than the mesh");
+ "on a fem variable having the same dimension as the mesh");
std::string adapted_lawname = adapt_law_name(lawname, N);