[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: |
Fri, 9 Nov 2018 12:27:02 -0500 (EST) |
branch: master
commit a14b31c7f24724b7de18ba967f8bb129e18bec76
Author: Konstantinos Poulios <address@hidden>
Date: Fri Nov 9 18:26:52 2018 +0100
Code clean up
---
interface/src/gf_asm.cc | 2 +-
src/bgeot_geometric_trans.cc | 20 ++++++------
src/getfem/getfem_projected_fem.h | 8 ++---
src/getfem_contact_and_friction_common.cc | 30 ++++++++---------
src/getfem_generic_assembly_semantic.cc | 19 +++++------
src/getfem_projected_fem.cc | 2 +-
src/getfem_regular_meshes.cc | 53 ++++++++++++++-----------------
7 files changed, 64 insertions(+), 70 deletions(-)
diff --git a/interface/src/gf_asm.cc b/interface/src/gf_asm.cc
index d746be6..42416af 100644
--- a/interface/src/gf_asm.cc
+++ b/interface/src/gf_asm.cc
@@ -1254,7 +1254,7 @@ void gf_asm(getfemint::mexargs_in& m_in,
getfemint::mexargs_out& m_out) {
);
- /address@hidden ('expression analysis', @str expression [, address@hidden
mesh | @tmim mim}] [, der_order] [, @tmodel model] [, @str varname, @int
is_variable[, address@hidden mf | @tmimd mimd}], ...])
+ /address@hidden ('expression analysis', @str expression [, address@hidden
mesh | @tmim mim}] [, der_order] [, @tmodel model] [, @str varname, @int
is_variable[, address@hidden mf | @tmimd mimd}], ...])
Analyse a high-level generic assembly expression and print
information about the provided address@hidden/
sub_command
diff --git a/src/bgeot_geometric_trans.cc b/src/bgeot_geometric_trans.cc
index 43d1c7d..7d29aa1 100644
--- a/src/bgeot_geometric_trans.cc
+++ b/src/bgeot_geometric_trans.cc
@@ -242,7 +242,7 @@ namespace bgeot {
scalar_type a0 = A[4]*A[8] - A[5]*A[7], a1 = A[5]*A[6] - A[3]*A[8];
scalar_type a2 = A[3]*A[7] - A[4]*A[6];
scalar_type det = A[0] * a0 + A[1] * a1 + A[2] * a2;
- GMM_ASSERT1(det != scalar_type(0), "Non invertible matrix");
+ GMM_ASSERT1(det != scalar_type(0), "Non invertible matrix");
scalar_type a3 = (A[2]*A[7] - A[1]*A[8]), a6 = (A[1]*A[5] - A[2]*A[4]);
scalar_type a4 = (A[0]*A[8] - A[2]*A[6]), a7 = (A[2]*A[3] - A[0]*A[5]);
scalar_type a5 = (A[1]*A[6] - A[0]*A[7]), a8 = (A[0]*A[4] - A[1]*A[3]);
@@ -343,9 +343,9 @@ namespace bgeot {
J__ = it[0] * a0 + it[1] * a1 + it[2] * a2;
} break;
default:
- B_factors.base_resize(P, P); // store factorization for B computation
- gmm::copy(gmm::transposed(KK), B_factors);
- ipvt.resize(P);
+ B_factors.base_resize(P, P); // store factorization for B computation
+ gmm::copy(gmm::transposed(KK), B_factors);
+ ipvt.resize(P);
bgeot::lu_factor(&(*(B_factors.begin())), ipvt, P);
J__ = bgeot::lu_det(&(*(B_factors.begin())), ipvt, P);
break;
@@ -387,8 +387,10 @@ namespace bgeot {
case 2:
{
auto it = &(*(KK.begin())); auto itB = &(*(B_.begin()));
- *itB++ = it[3] / J__; *itB++ = -it[2] / J__;
- *itB++ = -it[1] / J__; *itB = (*it) / J__;
+ *itB++ = it[3] / J__;
+ *itB++ = -it[2] / J__;
+ *itB++ = -it[1] / J__;
+ *itB = (*it) / J__;
} break;
case 3:
{
@@ -1073,11 +1075,9 @@ namespace bgeot {
/* norm of returned vector is the ratio between the face surface on
the reference element and the face surface on the real element.
IT IS NOT UNITARY
-
- pt is the position of the evaluation point on the reference element
*/
- base_small_vector compute_normal(const geotrans_interpolation_context& c,
- size_type face) {
+ base_small_vector
+ compute_normal(const geotrans_interpolation_context& c, size_type face) {
GMM_ASSERT1(c.G().ncols() == c.pgt()->nb_points(), "dimensions mismatch");
base_small_vector un(c.N());
gmm::mult(c.B(), c.pgt()->normals()[face], un);
diff --git a/src/getfem/getfem_projected_fem.h
b/src/getfem/getfem_projected_fem.h
index dc06452..fb8b3d6 100644
--- a/src/getfem/getfem_projected_fem.h
+++ b/src/getfem/getfem_projected_fem.h
@@ -95,8 +95,8 @@ namespace getfem {
mutable bgeot::kdtree tree; // Tree containing the nodes of the
// projected mf_source dofs
mutable std::vector<size_type> ind_dof; /* all functions using this work
- array should keep it full of
- size_type(-1) */
+ array should keep it filled
+ with size_type(-1) */
mutable bgeot::geotrans_inv_convex gic;
mutable base_tensor taux;
mutable fem_interpolation_context fictx;
@@ -106,9 +106,9 @@ namespace getfem {
mutable bgeot::multi_index mi2, mi3;
mutable base_node ptref;
- void build_kdtree(void) const;
+ void build_kdtree() const;
- bool find_a_projected_point(base_node pt, base_node &ptr_proj,
+ bool find_a_projected_point(const base_node &pt, base_node &ptr_proj,
size_type &cv_proj, short_type &fc_proj) const;
virtual void update_from_context(void) const;
diff --git a/src/getfem_contact_and_friction_common.cc
b/src/getfem_contact_and_friction_common.cc
index fada408..e4d05a6 100644
--- a/src/getfem_contact_and_friction_common.cc
+++ b/src/getfem_contact_and_friction_common.cc
@@ -1933,7 +1933,7 @@ namespace getfem {
gmm::add(gmm::identity_matrix(), F_x);
gmm::copy(F_x, F_x_inv);
bgeot::lu_inverse(&(*(F_x_inv.begin())), N);
-
+
base_tensor base_ux;
base_matrix vbase_ux;
@@ -1941,7 +1941,7 @@ namespace getfem {
size_type qdim_ux = pfu_x->target_dim();
size_type ndof_ux = pfu_x->nb_dof(cv_x) * N / qdim_ux;
vectorize_base_tensor(base_ux, vbase_ux, ndof_ux, qdim_ux, N);
-
+
base_tensor base_uy;
base_matrix vbase_uy;
size_type ndof_uy = 0;
@@ -1951,7 +1951,7 @@ namespace getfem {
ndof_uy = pfu_y->nb_dof(cv_y) * N / qdim_uy;
vectorize_base_tensor(base_uy, vbase_uy, ndof_uy, qdim_uy, N);
}
-
+
base_tensor grad_base_ux, vgrad_base_ux;
ctx_x.grad_base_value(grad_base_ux);
vectorize_grad_base_tensor(grad_base_ux, vgrad_base_ux, ndof_ux,
@@ -1965,7 +1965,7 @@ namespace getfem {
gmm::mult(F_y_inv, I_nxny, M1);
base_matrix der_x(ndof_ux, N);
gmm::mult(vbase_ux, gmm::transposed(M1), der_x);
-
+
// -F_y^{-1}*I_nxny*Test_u(Y)
base_matrix der_y(ndof_uy, N);
if (ret_type == 1) {
@@ -2743,22 +2743,22 @@ namespace getfem {
scalar_type norm(0);
if (tau > scalar_type(0)) {
- gmm::add(lambda, gmm::scaled(Vs, -r), F);
- scalar_type mu = gmm::vect_sp(F, n)/nn;
- gmm::add(gmm::scaled(n, -mu/nn), F);
+ gmm::add(lambda, gmm::scaled(Vs, -r), F); // F <-- lambda -r*Vs
+ scalar_type mu = gmm::vect_sp(F, n)/nn; // mu <-- (lambda
-r*Vs).n/|n|
+ gmm::add(gmm::scaled(n, -mu/nn), F); // F <-- (lambda -r*Vs)*(I-n
x n / |n|²)
norm = gmm::vect_norm2(F);
- gmm::copy(gmm::identity_matrix(), dn);
- gmm::scale(dn, -mu/nn);
- gmm::rank_one_update(dn, gmm::scaled(n, mu/(nn*nn*nn)), n);
- gmm::rank_one_update(dn, gmm::scaled(n, scalar_type(-1)/(nn*nn)), F);
- gmm::copy(gmm::identity_matrix(), dVs);
- gmm::rank_one_update(dVs, n, gmm::scaled(n, scalar_type(-1)/(nn*nn)));
+ gmm::copy(gmm::identity_matrix(), dn); // dn
<-- I
+ gmm::scale(dn, -mu/nn); // dn
<-- -(lambda -r*Vs).n/|n|² I
+ gmm::rank_one_update(dn, gmm::scaled(n, mu/(nn*nn*nn)), n); // dn
<-- -(lambda -r*Vs).n/|n|² (I - n x n/|n|²)
+ gmm::rank_one_update(dn, gmm::scaled(n, scalar_type(-1)/(nn*nn)), F);
// dn <-- -(lambda -r*Vs).n/|n|² (I - n x n/|n|²) + n x ((lambda -r*Vs)*(I-n x
n / |n|²)) /|n|²
+ gmm::copy(gmm::identity_matrix(), dVs);
// dVs <-- I
+ gmm::rank_one_update(dVs, n, gmm::scaled(n, scalar_type(-1)/(nn*nn)));
// dVs <-- I - n x n/|n|²
- if (norm > tau) {
+ if (norm > tau) { // slip
gmm::rank_one_update(dVs, F,
gmm::scaled(F, scalar_type(-1)/(norm*norm)));
gmm::scale(dVs, tau / norm);
- gmm::copy(gmm::scaled(F, scalar_type(1)/norm), dg);
+ gmm::copy(gmm::scaled(F, scalar_type(1)/norm), dg);
// dg <-- Normalized((lambda -r*Vs)*(I-n x n / |n|²))
gmm::rank_one_update(dn, gmm::scaled(F, mu/(norm*norm*nn)), F);
gmm::scale(dn, tau / norm);
gmm::scale(F, tau / norm);
diff --git a/src/getfem_generic_assembly_semantic.cc
b/src/getfem_generic_assembly_semantic.cc
index 524e05e..956f801 100644
--- a/src/getfem_generic_assembly_semantic.cc
+++ b/src/getfem_generic_assembly_semantic.cc
@@ -375,13 +375,13 @@ namespace getfem {
}
# define ga_valid_operand(pnode) \
- { \
- if (pnode && (pnode->node_type == GA_NODE_PREDEF_FUNC || \
- pnode->node_type == GA_NODE_SPEC_FUNC || \
- pnode->node_type == GA_NODE_NAME || \
- pnode->node_type == GA_NODE_OPERATOR || \
- pnode->node_type == GA_NODE_ALLINDICES)) \
- ga_throw_error(pnode->expr, pnode->pos, "Invalid term"); \
+ { \
+ if (pnode && (pnode->node_type == GA_NODE_PREDEF_FUNC || \
+ pnode->node_type == GA_NODE_SPEC_FUNC || \
+ pnode->node_type == GA_NODE_NAME || \
+ pnode->node_type == GA_NODE_OPERATOR || \
+ pnode->node_type == GA_NODE_ALLINDICES)) \
+ ga_throw_error(pnode->expr, pnode->pos, "Invalid term"); \
}
static void ga_node_analysis(ga_tree &tree,
@@ -3175,7 +3175,7 @@ namespace getfem {
pnode_trans = pnode->parent->children[1];
}
- if (ivar) {
+ if (ivar) { // Derivative wrt the interpolated variable
mi.resize(1); mi[0] = 2;
for (size_type i = 0; i < pnode->tensor_order(); ++i)
mi.push_back(pnode->tensor_proper_size(i));
@@ -3191,7 +3191,8 @@ namespace getfem {
pnode->test_function_type = order;
}
- if (itrans) {
+ if (itrans) { // Derivative with respect to a variable that the
+ // interpolate transformation depends on
const mesh_fem *mf = workspace.associated_mf(pnode_trans->name);
size_type q = workspace.qdim(pnode_trans->name);
size_type n = mf->linked_mesh().dim();
diff --git a/src/getfem_projected_fem.cc b/src/getfem_projected_fem.cc
index 34adf8d..0050b0b 100644
--- a/src/getfem_projected_fem.cc
+++ b/src/getfem_projected_fem.cc
@@ -251,7 +251,7 @@ namespace getfem {
tree.add_point_with_id(mf_source.point_of_basic_dof(dof), dof);
}
- bool projected_fem::find_a_projected_point(base_node pt, base_node &ptr_proj,
+ bool projected_fem::find_a_projected_point(const base_node &pt, base_node
&ptr_proj,
size_type &cv_proj, short_type
&fc_proj) const {
bgeot::index_node_pair ipt;
diff --git a/src/getfem_regular_meshes.cc b/src/getfem_regular_meshes.cc
index c319277..d856c18 100644
--- a/src/getfem_regular_meshes.cc
+++ b/src/getfem_regular_meshes.cc
@@ -294,17 +294,15 @@ namespace getfem
std::string GT = PARAM.string_value("GT");
GMM_ASSERT1(!GT.empty(), "regular mesh : you have at least to "
"specify the geometric transformation");
- bgeot::pgeometric_trans pgt =
- bgeot::geometric_trans_descriptor(GT);
+ bgeot::pgeometric_trans pgt = bgeot::geometric_trans_descriptor(GT);
size_type N = pgt->dim();
base_small_vector org(N); gmm::clear(org);
- const std::vector<bgeot::md_param::param_value> &o
- = PARAM.array_value("ORG");
+ const auto &o = PARAM.array_value("ORG");
if (o.size() > 0) {
- GMM_ASSERT1(o.size() == N, "ORG parameter should be an array of size "
- << N);
+ GMM_ASSERT1(o.size() == N,
+ "ORG parameter should be an array of size " << N);
for (size_type i = 0; i < N; ++i) {
GMM_ASSERT1(o[i].type_of_param() == bgeot::md_param::REAL_VALUE,
"ORG should be a real array.");
@@ -316,14 +314,13 @@ namespace getfem
std::vector<size_type> nsubdiv(N);
gmm::fill(nsubdiv, 2);
- const std::vector<bgeot::md_param::param_value> &ns
- = PARAM.array_value("NSUBDIV");
+ const auto &ns = PARAM.array_value("NSUBDIV");
if (ns.size() > 0) {
GMM_ASSERT1(ns.size() == N,
"NSUBDIV parameter should be an array of size " << N);
for (size_type i = 0; i < N; ++i) {
GMM_ASSERT1(ns[i].type_of_param() == bgeot::md_param::REAL_VALUE,
- "NSUBDIV should be an integer array.");
+ "NSUBDIV should be an integer array");
nsubdiv[i] = size_type(ns[i].real()+0.5);
}
}
@@ -331,14 +328,13 @@ namespace getfem
base_small_vector sizes(N);
gmm::fill(sizes, 1.0);
- const std::vector<bgeot::md_param::param_value> &si
- = PARAM.array_value("SIZES");
+ const auto &si = PARAM.array_value("SIZES");
if (si.size() > 0) {
GMM_ASSERT1(si.size() == N,
"SIZES parameter should be an array of size " << N);
for (size_type i = 0; i < N; ++i) {
GMM_ASSERT1(si[i].type_of_param() == bgeot::md_param::REAL_VALUE,
- "SIZES should be a real array.");
+ "SIZES should be a real array");
sizes[i] = si[i].real();
}
}
@@ -361,20 +357,18 @@ namespace getfem
std::string GT = PARAM.string_value("GT");
GMM_ASSERT1(!GT.empty(), "regular ball mesh : you have at least to "
"specify the geometric transformation");
- bgeot::pgeometric_trans pgt =
- bgeot::geometric_trans_descriptor(GT);
+ bgeot::pgeometric_trans pgt = bgeot::geometric_trans_descriptor(GT);
size_type N = pgt->dim();
base_small_vector org(N);
- const std::vector<bgeot::md_param::param_value> &o
- = PARAM.array_value("ORG");
+ const auto &o = PARAM.array_value("ORG");
if (o.size() > 0) {
- GMM_ASSERT1(o.size() == N, "ORG parameter should be an array of size "
- << N);
+ GMM_ASSERT1(o.size() == N,
+ "ORG parameter should be an array of size " << N);
for (size_type i = 0; i < N; ++i) {
GMM_ASSERT1(o[i].type_of_param() == bgeot::md_param::REAL_VALUE,
- "ORG should be a real array.");
+ "ORG should be a real array");
org[i] = o[i].real();
}
}
@@ -383,29 +377,26 @@ namespace getfem
bool noised = (PARAM.int_value("NOISED") != 0);
size_type nsubdiv0(2), nsubdiv1(2);
- const std::vector<bgeot::md_param::param_value> &ns
- = PARAM.array_value("NSUBDIV");
+ const auto &ns = PARAM.array_value("NSUBDIV");
if (ns.size() > 0) {
GMM_ASSERT1(ns.size() == 2,
- "NSUBDIV parameter should be an array of size " << 2);
+ "NSUBDIV parameter should be an array of size 2");
for (size_type i = 0; i < 2; ++i)
GMM_ASSERT1(ns[i].type_of_param() == bgeot::md_param::REAL_VALUE,
- "NSUBDIV should be an integer array.");
+ "NSUBDIV should be an integer array");
nsubdiv0 = size_type(ns[0].real()+0.5);
nsubdiv1 = size_type(ns[1].real()+0.5);
}
scalar_type radius(1), core_ratio(M_SQRT1_2);
- const std::vector<bgeot::md_param::param_value> &si
- = PARAM.array_value("SIZES");
+ const auto &si = PARAM.array_value("SIZES");
if (si.size() > 0) {
GMM_ASSERT1(si.size() == 1,
- "SIZES parameter should be an array of size " << 1);
+ "SIZES parameter should be an array of size 1");
GMM_ASSERT1(si[0].type_of_param() == bgeot::md_param::REAL_VALUE,
- "SIZES should be a real array.");
+ "SIZES should be a real array");
radius = si[0].real();
}
-
std::vector<size_type> nsubdiv(N);
gmm::fill(nsubdiv, nsubdiv0);
@@ -442,7 +433,8 @@ namespace getfem
trsl[i] = core_ratio;
mm[i].translation(trsl);
for (dal::bv_visitor cv(mm[i].convex_index()); !cv.finished(); ++cv)
- m.add_convex_by_points(mm[i].trans_of_convex(cv),
mm[i].points_of_convex(cv).begin());
+ m.add_convex_by_points(mm[i].trans_of_convex(cv),
+ mm[i].points_of_convex(cv).begin());
}
std::vector<base_node> pts(m.points().card(), base_node(N));
@@ -503,7 +495,8 @@ namespace getfem
m0.copy_from(m);
m0.transformation(M);
for (dal::bv_visitor cv(m0.convex_index()); !cv.finished(); ++cv)
- m.add_convex_by_points(m0.trans_of_convex(cv),
m0.points_of_convex(cv).begin());
+ m.add_convex_by_points(m0.trans_of_convex(cv),
+ m0.points_of_convex(cv).begin());
}
m.translation(org);