[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
[Getfem-commits] r5450 - in /trunk/getfem: contrib/opt_assembly/ src/ sr
|
From: |
Yves . Renard |
|
Subject: |
[Getfem-commits] r5450 - in /trunk/getfem: contrib/opt_assembly/ src/ src/getfem/ |
|
Date: |
Fri, 28 Oct 2016 08:15:43 -0000 |
Author: renard
Date: Fri Oct 28 10:15:42 2016
New Revision: 5450
URL: http://svn.gna.org/viewcvs/getfem?rev=5450&view=rev
Log:
minor optimizations
Modified:
trunk/getfem/contrib/opt_assembly/opt_assembly.cc
trunk/getfem/src/getfem/bgeot_geometric_trans.h
trunk/getfem/src/getfem/getfem_mesh_fem.h
trunk/getfem/src/getfem_generic_assembly.cc
Modified: trunk/getfem/contrib/opt_assembly/opt_assembly.cc
URL:
http://svn.gna.org/viewcvs/getfem/trunk/getfem/contrib/opt_assembly/opt_assembly.cc?rev=5450&r1=5449&r2=5450&view=diff
==============================================================================
--- trunk/getfem/contrib/opt_assembly/opt_assembly.cc (original)
+++ trunk/getfem/contrib/opt_assembly/opt_assembly.cc Fri Oct 28 10:15:42 2016
@@ -255,7 +255,7 @@
chrono ch;
bool all = false;
- bool select = true;
+ bool select = false;
int only_one = 1;
if (all || select || only_one == 1) {
@@ -424,7 +424,7 @@
MAT_TEST_1("Test for linear non homogeneous elasticity stiffness matrix",
ndofu, ndofu, "(Div_Test_u*(lambda2*Id(meshdim)) "
- "+ mu2*(Grad_Test_u'+Grad_Test_u)):Grad_Test2_u",
+ "+ mu2*Sym(Grad_Test_u)):Grad_Test2_u",
mim2, Iu, Iu,
getfem::asm_stiffness_matrix_for_linear_elasticity
(K, mim2, mf_u, mf_p, lambda2, mu2));
@@ -438,7 +438,7 @@
GMM_SET_EXCEPTION_DEBUG; // Exceptions make a memory fault, to debug.
FE_ENABLE_EXCEPT; // Enable floating point exception for Nan.
- bool all = true;
+ bool all = false;
int only_one = 2;
// Mesured times for
@@ -453,49 +453,49 @@
// new | old | sto | asse | exec | Ins |
if (all || only_one == 1) // ndofu = 321602 ndofp = 160801 ndofchi = 1201
test_new_assembly(2, 400, 1);
- // Vector source term : 0.25 | 0.68 |
- // Nonlinear residual : 0.34 | |
- // Mass (scalar) : 0.18 | 0.58 | 0.04 | 0.06 | 0.06 | 0.06 |
- // Mass (vector) : 0.30 | 0.82 | 0.09 | 0.15 | 0.06 | 0.09 |
+ // Vector source term : 0.20 | 0.66 |
+ // Nonlinear residual : 0.31 | |
+ // Mass (scalar) : 0.18 | 0.56 | 0.04 | 0.06 | 0.06 | 0.06 |
+ // Mass (vector) : 0.30 | 0.80 | 0.09 | 0.15 | 0.06 | 0.09 |
// Laplacian : 0.16 | 0.80 | 0.04 | 0.05 | 0.06 | 0.05 |
- // Homogeneous elas : 0.31 | 1.88 | 0.08 | 0.13 | 0.06 | 0.10 |
- // Non-homogeneous elast: 0.36 | 2.26 | 0.09 | 0.16 | 0.06 | 0.14 |
+ // Homogeneous elas : 0.30 | 1.88 | 0.08 | 0.13 | 0.06 | 0.09 |
+ // Non-homogeneous elast: 0.35 | 2.26 | 0.09 | 0.16 | 0.06 | 0.13 |
if (all || only_one == 2) // ndofu = 151959 ndofp = 50653 ndofchi = 6553
test_new_assembly(3, 36, 1);
- // Vector source term : 0.31 | 0.82 |
- // Nonlinear residual : 0.59 | |
+ // Vector source term : 0.26 | 0.79 |
+ // Nonlinear residual : 0.54 | |
// Mass (scalar) : 0.21 | 0.58 | 0.05 | 0.08 | 0.08 | 0.06 |
- // Mass (vector) : 0.69 | 1.37 | 0.17 | 0.27 | 0.08 | 0.34 |
- // Laplacian : 0.18 | 1.15 | 0.03 | 0.07 | 0.08 | 0.03 |
- // Homogeneous elas : 0.79 | 4.29 | 0.26 | 0.33 | 0.08 | 0.38 |
- // Non-homogeneous elast: 0.86 | 6.35 | 0.26 | 0.33 | 0.08 | 0.45 |
+ // Mass (vector) : 0.68 | 1.37 | 0.17 | 0.27 | 0.08 | 0.33 |
+ // Laplacian : 0.17 | 1.15 | 0.03 | 0.06 | 0.08 | 0.03 |
+ // Homogeneous elas : 0.79 | 4.25 | 0.26 | 0.33 | 0.08 | 0.38 |
+ // Non-homogeneous elast: 0.83 | 6.29 | 0.26 | 0.33 | 0.08 | 0.42 |
if (all || only_one == 3) // ndofu = 321602 ndofp = 160801 ndofchi = 1201
test_new_assembly(2, 200, 2);
- // Vector source term : 0.11 | 0.24 |
- // Nonlinear residual : 0.17 | |
+ // Vector source term : 0.09 | 0.23 |
+ // Nonlinear residual : 0.15 | |
// Mass (scalar) : 0.09 | 0.25 | 0.02 | 0.03 | 0.03 | 0.03 |
// Mass (vector) : 0.26 | 0.44 | 0.05 | 0.09 | 0.03 | 0.14 |
// Laplacian : 0.09 | 0.37 | 0.02 | 0.03 | 0.03 | 0.03 |
// Homogeneous elas : 0.26 | 1.28 | 0.06 | 0.09 | 0.03 | 0.14 |
- // Non-homogeneous elast: 0.30 | 2.38 | 0.07 | 0.10 | 0.03 | 0.17 |
+ // Non-homogeneous elast: 0.28 | 2.38 | 0.06 | 0.10 | 0.03 | 0.16 |
if (all || only_one == 4) // ndofu = 151959 ndofp = 50653 ndofchi = 6553
test_new_assembly(3, 18, 2);
- // Vector source term : 0.16 | 0.24 |
- // Nonlinear residual : 0.39 | |
+ // Vector source term : 0.13 | 0.23 |
+ // Nonlinear residual : 0.37 | |
// Mass (scalar) : 0.11 | 0.25 | 0.05 | 0.05 | 0.03 | 0.03 |
// Mass (vector) : 1.13 | 0.89 | 0.21 | 0.35 | 0.03 | 0.75 |
// Laplacian : 0.10 | 0.53 | 0.03 | 0.04 | 0.03 | 0.03 |
// Homogeneous elas : 1.66 | 3.35 | 0.59 | 0.73 | 0.03 | 0.90 |
- // Non-homogeneous elast: 1.72 | 9.15 | 0.59 | 0.73 | 0.03 | 0.96 |
+ // Non-homogeneous elast: 1.70 | 9.08 | 0.59 | 0.73 | 0.03 | 0.94 |
if (all || only_one == 5) // ndofu = 151959 ndofp = 50653 ndofchi = 6553
test_new_assembly(3, 9, 4);
- // Vector source term : 0.13 | 0.20 |
- // Nonlinear residual : 0.38 | |
+ // Vector source term : 0.11 | 0.19 |
+ // Nonlinear residual : 0.36 | |
// Mass (scalar) : 0.51 | 0.34 | 0.09 | 0.16 | 0.01 | 0.33 |
// Mass (vector) : 4.37 | 1.31 | 0.41 | 1.27 | 0.01 | 3.10 |
- // Laplacian : 0.40 | 0.77 | 0.09 | 0.14 | 0.01 | 0.25 |
- // Homogeneous elas : 8.93 | 5.23 | 0.88 | 1.43 | 0.01 | 7.49 |
- // Non-homogeneous elast: 9.01 | 47.4 | 0.76 | 1.35 | 0.01 | 7.65 |
+ // Laplacian : 0.40 | 0.76 | 0.09 | 0.14 | 0.01 | 0.25 |
+ // Homogeneous elas : 8.93 | 5.23 | 0.82 | 1.41 | 0.01 | 7.49 |
+ // Non-homogeneous elast: 8.94 | 47.4 | 0.82 | 1.41 | 0.01 | 7.50 |
// Conclusions :
// - Deactivation of debug test has no sensible effect.
Modified: trunk/getfem/src/getfem/bgeot_geometric_trans.h
URL:
http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/getfem/bgeot_geometric_trans.h?rev=5450&r1=5449&r2=5450&view=diff
==============================================================================
--- trunk/getfem/src/getfem/bgeot_geometric_trans.h (original)
+++ trunk/getfem/src/getfem/bgeot_geometric_trans.h Fri Oct 28 10:15:42 2016
@@ -511,6 +511,10 @@
/* Optimized operation for small matrices */
scalar_type lu_det(const scalar_type *A, size_type N);
scalar_type lu_inverse(scalar_type *A, size_type N, bool doassert = true);
+ inline scalar_type lu_det(const base_matrix &A)
+ { return lu_det(&(*(A.begin())), A.nrows()); }
+ inline scalar_type lu_inverse(base_matrix &A, bool doassert = true)
+ { return lu_inverse(&(*(A.begin())), A.nrows(), doassert); }
} /* end of namespace bgeot. */
Modified: trunk/getfem/src/getfem/getfem_mesh_fem.h
URL:
http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/getfem/getfem_mesh_fem.h?rev=5450&r1=5449&r2=5450&view=diff
==============================================================================
--- trunk/getfem/src/getfem/getfem_mesh_fem.h (original)
+++ trunk/getfem/src/getfem/getfem_mesh_fem.h Fri Oct 28 10:15:42 2016
@@ -637,47 +637,29 @@
template <typename VEC1, typename VEC2>
void slice_vector_on_basic_dof_of_element(const mesh_fem &mf,
const VEC1 &vec,
- size_type cv, VEC2 &coeff) {
- size_type nbdof = mf.nb_basic_dof();
- size_type qmult = gmm::vect_size(vec) / nbdof;
- GMM_ASSERT1(gmm::vect_size(vec) == qmult * nbdof, "Bad dof vector size");
-
+ size_type cv, VEC2 &coeff,
+ size_type qmult1 = size_type(-1),
+ size_type qmult2 = size_type(-1)) {
+ if (qmult1 == size_type(-1)) {
+ size_type nbdof = mf.nb_basic_dof();
+ qmult1 = gmm::vect_size(vec) / nbdof;
+ GMM_ASSERT1(gmm::vect_size(vec) == qmult1 * nbdof, "Bad dof vector
size");
+ }
+ if (qmult2 == size_type(-1)) {
+ qmult2 = mf.get_qdim();
+ if (qmult2 > 1) qmult2 /= mf.fem_of_element(cv)->target_dim();
+ }
+ size_type qmultot = qmult1*qmult2;
auto &ct = mf.ind_scalar_basic_dof_of_element(cv);
- size_type qmult2 = mf.get_qdim();
- if (qmult2 > 1) qmult2 /= mf.fem_of_element(cv)->target_dim();
- size_type qmultot = qmult*qmult2;
gmm::resize(coeff, ct.size()*qmultot);
-
+
auto it = ct.begin();
auto itc = coeff.begin();
if (qmultot == 1) {
for (; it != ct.end(); ++it) *itc++ = vec[*it];
} else {
for (; it != ct.end(); ++it) {
- auto itv = vec.begin()+(*it)*qmult;
- for (size_type m = 0; m < qmultot; ++m) *itc++ = *itv++;
- }
- }
- }
-
- template <typename VEC1, typename VEC2>
- void slice_vector_on_basic_dof_of_element(const mesh_fem &mf,
- const VEC1 &vec,
- size_type cv, VEC2 &coeff,
- size_type qmult) {
- auto &ct = mf.ind_scalar_basic_dof_of_element(cv);
- size_type qmult2 = mf.get_qdim();
- if (qmult2 > 1) qmult2 /= mf.fem_of_element(cv)->target_dim();
- size_type qmultot = qmult*qmult2;
- gmm::resize(coeff, ct.size()*qmultot);
-
- auto it = ct.begin();
- auto itc = coeff.begin();
- if (qmultot == 1) {
- for (; it != ct.end(); ++it) *itc++ = vec[*it];
- } else {
- for (; it != ct.end(); ++it) {
- auto itv = vec.begin()+(*it)*qmult;
+ auto itv = vec.begin()+(*it)*qmult1;
for (size_type m = 0; m < qmultot; ++m) *itc++ = *itv++;
}
}
Modified: trunk/getfem/src/getfem_generic_assembly.cc
URL:
http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/getfem_generic_assembly.cc?rev=5450&r1=5449&r2=5450&view=diff
==============================================================================
--- trunk/getfem/src/getfem_generic_assembly.cc (original)
+++ trunk/getfem/src/getfem_generic_assembly.cc Fri Oct 28 10:15:42 2016
@@ -2340,41 +2340,48 @@
size_type N = args[0]->sizes()[0];
__mat_aux1().base_resize(N, N);
gmm::copy(args[0]->as_vector(), __mat_aux1().as_vector());
- scalar_type det = bgeot::lu_inverse(&(*(__mat_aux1().begin())), N);
+ scalar_type det = bgeot::lu_inverse(__mat_aux1());
if (det == scalar_type(0))
gmm::clear(result.as_vector());
else {
- base_tensor::iterator it = result.begin();
- for (size_type j = 0; j < N; ++j)
- for (size_type i = 0; i < N; ++i, ++it)
- *it = __mat_aux1()(j, i) * det;
+ auto it = result.begin();
+ auto ita = __mat_aux1().begin();
+ for (size_type j = 0; j < N; ++j, ++ita) {
+ auto itaa = ita;
+ for (size_type i = 0; i < N; ++i, ++it, itaa += N)
+ *it = (*itaa) * det;
+ }
GA_DEBUG_ASSERT(it == result.end(), "Internal error");
}
}
// Second derivative : det(M)(address@hidden - M^{-T}_{jk}M^{-T}_{li})
void second_derivative(const arg_list &args, size_type, size_type,
- base_tensor &result) const { // To be verified
+ base_tensor &result) const { // to be verified
size_type N = args[0]->sizes()[0];
__mat_aux1().base_resize(N, N);
gmm::copy(args[0]->as_vector(), __mat_aux1().as_vector());
- scalar_type det = bgeot::lu_inverse(&(*(__mat_aux1().begin())), N);
+ scalar_type det = bgeot::lu_inverse(__mat_aux1());
if (det == scalar_type(0))
gmm::clear(result.as_vector());
else {
- base_tensor::iterator it = result.begin();
- for (size_type l = 0; l < N; ++l)
- for (size_type k = 0; k < N; ++k)
- for (size_type j = 0; j < N; ++j)
- for (size_type i = 0; i < N; ++i, ++it)
- *it = (__mat_aux1()(j, i) * __mat_aux1()(l,k)
- - __mat_aux1()(j,k) * __mat_aux1()(l, i)) * det;
+ auto it = result.begin();
+ auto ita = __mat_aux1().begin(), ita_l = ita;
+ for (size_type l = 0; l < N; ++l, ++ita_l) {
+ auto ita_lk = ita_l, ita_jk = ita;
+ for (size_type k = 0; k < N; ++k, ita_lk += N, ita_jk += N) {
+ auto ita_j = ita;
+ for (size_type j = 0; j < N; ++j, ++ita_j, ++ita_jk) {
+ auto ita_ji = ita_j, ita_li = ita_l;
+ for (size_type i = 0; i < N; ++i, ++it, ita_ji += N, ita_li += N)
+ *it = ((*ita_ji) * (*ita_lk) - (*ita_jk) * (*ita_li)) * det;
+ }
+ }
+ }
GA_DEBUG_ASSERT(it == result.end(), "Internal error");
}
}
};
-
-
// Inverse Operator (for square matrices)
struct inverse_operator : public ga_nonlinear_operator {
@@ -2390,7 +2397,7 @@
size_type N = args[0]->sizes()[0];
__mat_aux1().base_resize(N, N);
gmm::copy(args[0]->as_vector(), __mat_aux1().as_vector());
- bgeot::lu_inverse(&(*(__mat_aux1().begin())), N);
+ bgeot::lu_inverse(__mat_aux1());
gmm::copy(__mat_aux1().as_vector(), result.as_vector());
}
@@ -2400,13 +2407,20 @@
size_type N = args[0]->sizes()[0];
__mat_aux1().base_resize(N, N);
gmm::copy(args[0]->as_vector(), __mat_aux1().as_vector());
- bgeot::lu_inverse(&(*(__mat_aux1().begin())), N);
- base_tensor::iterator it = result.begin();
- for (size_type l = 0; l < N; ++l)
- for (size_type k = 0; k < N; ++k)
- for (size_type j = 0; j < N; ++j)
- for (size_type i = 0; i < N; ++i, ++it)
- *it = -__mat_aux1()(i,k)*__mat_aux1()(l,j);
+ bgeot::lu_inverse(__mat_aux1());
+ auto it = result.begin();
+ auto ita = __mat_aux1().begin(), ita_l = ita;
+ for (size_type l = 0; l < N; ++l, ++ita_l) {
+ auto ita_k = ita;
+ for (size_type k = 0; k < N; ++k, ita_k += N) {
+ auto ita_lj = ita_l;
+ for (size_type j = 0; j < N; ++j, ++ita_lj) {
+ auto ita_ik = ita_k;
+ for (size_type i = 0; i < N; ++i, ++it, ++ita_ik)
+ *it = -(*ita_ik) * (*ita_lj);
+ }
+ }
+ }
GA_DEBUG_ASSERT(it == result.end(), "Internal error");
}
@@ -2418,7 +2432,7 @@
size_type N = args[0]->sizes()[0];
__mat_aux1().base_resize(N, N);
gmm::copy(args[0]->as_vector(), __mat_aux1().as_vector());
- bgeot::lu_inverse(&(*(__mat_aux1().begin())), N);
+ bgeot::lu_inverse(__mat_aux1());
base_tensor::iterator it = result.begin();
for (size_type n = 0; n < N; ++n)
for (size_type m = 0; m < N; ++m)
@@ -2431,9 +2445,6 @@
GA_DEBUG_ASSERT(it == result.end(), "Internal error");
}
};
-
-
-
//=========================================================================
// Initialization of predefined functions and operators.
@@ -2708,26 +2719,27 @@
cv_old(-1) {}
};
-
struct ga_instruction_slice_local_dofs : public ga_instruction {
const mesh_fem &mf;
const base_vector &U;
const fem_interpolation_context &ctx;
base_vector &coeff;
const size_type &ipt;
- size_type qmult;
+ size_type qmult1, qmult2;
virtual int exec() {
if (ipt == 0) {
GA_DEBUG_INFO("Instruction: Slice local dofs");
- slice_vector_on_basic_dof_of_element(mf,U,ctx.convex_num(),coeff,qmult);
+ slice_vector_on_basic_dof_of_element(mf, U, ctx.convex_num(),
+ coeff, qmult1, qmult2);
}
return 0;
}
ga_instruction_slice_local_dofs(const mesh_fem &mf_, const base_vector &U_,
const fem_interpolation_context &ctx_,
base_vector &coeff_, const size_type &ipt_,
- size_type qmult_)
- : mf(mf_), U(U_), ctx(ctx_), coeff(coeff_), ipt(ipt_), qmult(qmult_) {}
+ size_type qmult1_, size_type qmult2_)
+ : mf(mf_), U(U_), ctx(ctx_), coeff(coeff_), ipt(ipt_),
+ qmult1(qmult1_), qmult2(qmult2_) {}
};
struct ga_instruction_update_pfp : public ga_instruction {
@@ -10598,10 +10610,14 @@
rmi.local_dofs_hierarchy[pnode->name].push_back(if_hierarchy);
extend_variable_in_gis(workspace, pnode->name, gis);
// cout << "local dof of " << pnode->name << endl;
+ size_type qmult2 = mf->get_qdim();
+ if (qmult2 > 1 && !(mf->is_uniformly_vectorized()))
+ qmult2 = size_type(-1);
pgai = std::make_shared<ga_instruction_slice_local_dofs>
(*mf, *(gis.extended_vars[pnode->name]), gis.ctx,
rmi.local_dofs[pnode->name], gis.ipt,
- gis.extended_vars[pnode->name]->size() / mf->nb_basic_dof());
+ gis.extended_vars[pnode->name]->size() / mf->nb_basic_dof(),
+ qmult2);
rmi.instructions.push_back(std::move(pgai));
}
| [Prev in Thread] |
Current Thread |
[Next in Thread] |
- [Getfem-commits] r5450 - in /trunk/getfem: contrib/opt_assembly/ src/ src/getfem/,
Yves . Renard <=