Author: renard Date: Wed Oct 12 14:21:01 2016 New Revision: 5410 URL: http://svn.gna.org/viewcvs/getfem?rev=5410&view=rev Log: small fixes
Modified: trunk/getfem/contrib/opt_assembly/opt_assembly.cc trunk/getfem/src/bgeot_geometric_trans.cc trunk/getfem/src/getfem_generic_assembly.cc trunk/getfem/src/gmm/gmm_vector.h 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=5410&r1=5409&r2=5410&view=diff ============================================================================== --- trunk/getfem/contrib/opt_assembly/opt_assembly.cc (original) +++ trunk/getfem/contrib/opt_assembly/opt_assembly.cc Wed Oct 12 14:21:01 2016 @@ -433,29 +433,29 @@ // - Instructions execution except for assembly ones // new | old | sto | asse | exec | J | Ins | test_new_assembly(2, 400, 1); // ndofu = 321602 ndofp = 160801 ndofchi = 1201 - // Mass (scalar) : 0.37 | 0.61 | 0.07 | 0.10 | 0.21 | 0.06 | 0.06 | + // Mass (scalar) : 0.36 | 0.61 | 0.07 | 0.10 | 0.20 | 0.06 | 0.06 | // Mass (vector) : 0.46 | 0.82 | 0.12 | 0.22 | 0.15 | 0.05 | 0.09 | // Laplacian : 0.26 | 0.83 | 0.05 | 0.07 | 0.14 | 0.06 | 0.05 | - // Homogeneous elas : 0.47 | 1.88 | 0.14 | 0.22 | 0.14 | 0.05 | 0.11 | - // Non-homogeneous elast: 0.58 | 2.26 | 0.14 | 0.22 | 0.14 | 0.05 | 0.21 | + // Homogeneous elas : 0.46 | 1.88 | 0.14 | 0.22 | 0.13 | 0.05 | 0.11 | + // Non-homogeneous elast: 0.57 | 2.26 | 0.14 | 0.22 | 0.13 | 0.05 | 0.21 | test_new_assembly(3, 36, 1); // ndofu = 151959 ndofp = 50653 ndofchi = 6553 // Mass (scalar) : 0.44 | 0.77 | 0.08 | 0.12 | 0.24 | 0.11 | 0.08 | // Mass (vector) : 0.94 | 1.57 | 0.18 | 0.39 | 0.20 | 0.11 | 0.35 | - // Laplacian : 0.42 | 1.38 | 0.06 | 0.09 | 0.19 | 0.11 | 0.14 | + // Laplacian : 0.41 | 1.38 | 0.06 | 0.09 | 0.19 | 0.11 | 0.14 | // Homogeneous elas : 1.28 | 4.58 | 0.48 | 0.59 | 0.19 | 0.11 | 0.51 | - // Non-homogeneous elast: 1.44 | 6.81 | 0.49 | 0.63 | 0.19 | 0.11 | 0.62 | + // Non-homogeneous elast: 1.42 | 6.81 | 0.49 | 0.63 | 0.18 | 0.11 | 0.61 | test_new_assembly(2, 200, 2); // ndofu = 321602 ndofp = 160801 ndofchi = 1201 // Mass (scalar) : 0.15 | 0.25 | 0.04 | 0.07 | 0.05 | 0.02 | 0.03 | // Mass (vector) : 0.34 | 0.45 | 0.07 | 0.15 | 0.05 | 0.02 | 0.14 | // Laplacian : 0.13 | 0.37 | 0.03 | 0.05 | 0.04 | 0.02 | 0.04 | // Homogeneous elas : 0.37 | 1.28 | 0.13 | 0.19 | 0.04 | 0.01 | 0.14 | - // Non-homogeneous elast: 0.46 | 2.40 | 0.13 | 0.18 | 0.04 | 0.01 | 0.24 | + // Non-homogeneous elast: 0.45 | 2.40 | 0.13 | 0.18 | 0.04 | 0.01 | 0.23 | test_new_assembly(3, 18, 2); // ndofu = 151959 ndofp = 50653 ndofchi = 6553 // Mass (scalar) : 0.23 | 0.30 | 0.10 | 0.15 | 0.04 | 0.02 | 0.04 | // Mass (vector) : 1.46 | 0.90 | 0.33 | 0.64 | 0.04 | 0.02 | 0.78 | // Laplacian : 0.18 | 0.55 | 0.08 | 0.10 | 0.03 | 0.02 | 0.05 | // Homogeneous elas : 2.37 | 3.47 | 1.20 | 1.37 | 0.03 | 0.02 | 0.97 | - // Non-homogeneous elast: 2.46 | 9.25 | 1.20 | 1.37 | 0.03 | 0.02 | 1.06 | + // Non-homogeneous elast: 2.45 | 9.25 | 1.20 | 1.37 | 0.03 | 0.02 | 1.05 | test_new_assembly(3, 9, 4); // ndofu = 151959 ndofp = 50653 ndofchi = 6553 // Mass (scalar) : 0.75 | 0.38 | 0.25 | 0.39 | 0.01 | .005 | 0.35 | // Mass (vector) : 5.10 | 1.33 | 0.42 | 1.86 | 0.01 | .005 | 3.23 | Modified: trunk/getfem/src/bgeot_geometric_trans.cc URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/bgeot_geometric_trans.cc?rev=5410&r1=5409&r2=5410&view=diff ============================================================================== --- trunk/getfem/src/bgeot_geometric_trans.cc (original) +++ trunk/getfem/src/bgeot_geometric_trans.cc Wed Oct 12 14:21:01 2016 @@ -36,7 +36,7 @@ const base_node& geotrans_interpolation_context::xref() const { if (!have_xref()) { if (pspt_) xref_ = (*pspt_)[ii_]; - else GMM_ASSERT1(false, "missing xref"); + else GMM_ASSERT1(false, "Missing xref"); } return xref_; } @@ -51,7 +51,7 @@ } void geotrans_interpolation_context::compute_J(void) const { - GMM_ASSERT1(have_G() && have_pgt(), "unable to compute J\n"); + GMM_ASSERT1(have_G() && have_pgt(), "Unable to compute J\n"); size_type P = pgt_->structure()->dim(); if (P != N()) { B_factors.base_resize(P, P); @@ -80,15 +80,15 @@ const base_matrix& geotrans_interpolation_context::K() const { if (!have_K()) { - GMM_ASSERT1(have_G() && have_pgt(), "unable to compute K\n"); + GMM_ASSERT1(have_G() && have_pgt(), "Unable to compute K\n"); size_type P = pgt_->structure()->dim(); K_.base_resize(N(), P); if (have_pgp()) { - pgt_->compute_K_matrix(G(), pgp_->grad(ii_), K_); + pgt_->compute_K_matrix(G(), pgp_->grad(ii_), K_); } else { PC.base_resize(pgt()->nb_points(), P); pgt()->poly_vector_grad(xref(), PC); - pgt_->compute_K_matrix(G(), PC, K_); + pgt_->compute_K_matrix(G(), PC, K_); } have_K_ = true; } @@ -97,7 +97,7 @@ const base_matrix& geotrans_interpolation_context::B() const { if (!have_B()) { - GMM_ASSERT1(have_G() && have_pgt(), "unable to compute B\n"); + GMM_ASSERT1(have_G() && have_pgt(), "Unable to compute B\n"); size_type P = pgt_->structure()->dim(); const base_matrix &KK = K(); B_.base_resize(N(), P); @@ -108,20 +108,20 @@ gmm::mult(KK, PC, B_); } else if (P == 1) { scalar_type det = KK(0, 0); - GMM_ASSERT1(det != scalar_type(0), "non invertible matrix"); + GMM_ASSERT1(det != scalar_type(0), "Non invertible matrix"); B_(0, 0) = scalar_type(1)/det; J_ = gmm::abs(det); } else if (P == 2) { const scalar_type *p = &(KK(0,0)); scalar_type det = (*p) * (*(p+3)) - (*(p+1)) * (*(p+2)); - GMM_ASSERT1(det != scalar_type(0), "non invertible matrix"); + GMM_ASSERT1(det != scalar_type(0), "Non invertible matrix"); scalar_type *q = &(B_(0,0)); *q++ = (*(p+3)) / det; *q++ = -(*(p+2)) / det; *q++ = -(*(p+1)) / det; *q++ = (*p) / det; J_ = gmm::abs(det); } else { scalar_type det = J(); - GMM_ASSERT1(det != scalar_type(0), "non invertible matrix"); + GMM_ASSERT1(det != scalar_type(0), "Non invertible matrix"); gmm::lu_inverse(B_factors, ipvt, B_); } have_B_ = true; Modified: trunk/getfem/src/getfem_generic_assembly.cc URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/getfem_generic_assembly.cc?rev=5410&r1=5409&r2=5410&view=diff ============================================================================== --- trunk/getfem/src/getfem_generic_assembly.cc (original) +++ trunk/getfem/src/getfem_generic_assembly.cc Wed Oct 12 14:21:01 2016 @@ -2749,7 +2749,7 @@ size_type cv_1 = ctx.is_convex_num_valid() ? ctx.convex_num() : mf.convex_index().first_true(); pfem pf = mf.fem_of_element(cv_1); - GMM_ASSERT1(pf, "An element without finite element methode defined"); + GMM_ASSERT1(pf, "An element without finite element method defined"); size_type Qmult = qdim / pf->target_dim(); size_type s = pf->nb_dof(cv_1) * Qmult; if (t.sizes()[0] != s) @@ -6481,10 +6481,16 @@ 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() + GA_TIC; ga_instruction_set gis; ga_compile(*this, gis, order); - size_type ndof = gis.nb_dof, max_dof = gis.max_dof; + ndof = gis.nb_dof; + size_type max_dof = gis.max_dof; GA_TOCTIC("Compile time"); if (order == 2) { @@ -11721,7 +11727,7 @@ static void ga_exec(ga_instruction_set &gis, ga_workspace &workspace) { base_matrix G; - base_small_vector un, up; + base_small_vector un; scalar_type J(0); for (const std::string &t : gis.transformations) @@ -11739,10 +11745,11 @@ mesh_region rg(region); m.intersect_with_mpi_region(rg); size_type old_cv = size_type(-1); - bgeot::pgeometric_trans pgt = 0; + bgeot::pgeometric_trans pgt = 0, pgt_old = 0; pintegration_method pim = 0; papprox_integration pai = 0; bgeot::pstored_point_tab pspt = 0, old_pspt = 0; + bgeot::pgeotrans_precomp pgp = 0; bool first_gp = true; for (getfem::mr_visitor v(rg, m); !v.finished(); ++v) { if (mim.convex_index().is_in(v.cv())) { @@ -11750,27 +11757,27 @@ if (v.cv() != old_cv) { m.points_of_convex(v.cv(), G); pgt = m.trans_of_convex(v.cv()); - up.resize(G.nrows()); - un.resize(pgt->dim()); - pim = mim.int_method_of_element(v.cv()); + pim = mim.int_method_of_element(v.cv()); if (pim->type() == IM_NONE) continue; GMM_ASSERT1(pim->type() == IM_APPROX, "Sorry, exact methods cannot " "be used in high level generic assembly"); pai = pim->approx_method(); pspt = pai->pintegration_points(); if (pspt->size()) { - if (gis.ctx.have_pgp() && gis.pai == pai && gis.ctx.pgt()==pgt) { - gis.ctx.change(gis.ctx.pgp(), 0, 0, G, v.cv(), v.f()); + if (pgp && gis.pai == pai && pgt_old == pgt) { + gis.ctx.change(pgp, 0, 0, G, v.cv(), v.f()); } else { if (pai->is_built_on_the_fly()) { gis.ctx.change(pgt, 0, (*pspt)[0], G, v.cv(), v.f()); + pgp = 0; } else { - gis.ctx.change(gis.gp_pool(pgt, pspt), 0,0, G, v.cv(), v.f()); + pgp = gis.gp_pool(pgt, pspt); + gis.ctx.change(pgp, 0, 0, G, v.cv(), v.f()); } + pgt_old = pgt; gis.pai = pai; } - gis.pai = pai; - if (gis.need_elt_size) - gis.elt_size = m.convex_radius_estimate(v.cv())*scalar_type(2); + if (gis.need_elt_size) + gis.elt_size = convex_radius_estimate(pgt, G)*scalar_type(2); } old_cv = v.cv(); } else { @@ -11787,19 +11794,20 @@ first_ind = pai->ind_first_point_on_face(v.f()); } for (gis.ipt = 0; gis.ipt < gis.nbpt; ++(gis.ipt)) { - if (gis.ctx.have_pgp()) gis.ctx.set_ii(first_ind+gis.ipt); + if (pgp) gis.ctx.set_ii(first_ind+gis.ipt); else gis.ctx.set_xref((*pspt)[first_ind+gis.ipt]); if (gis.ipt == 0 || !(pgt->is_linear())) { J = gis.ctx.J(); // Computation of unit normal vector in case of a boundary if (v.f() != short_type(-1)) { - gmm::copy(pgt->normals()[v.f()], un); - gmm::mult(gis.ctx.B(), un, up); - scalar_type nup = gmm::vect_norm2(up); + gis.Normal.resize(G.nrows()); + un.resize(pgt->dim()); + gmm::copy(pgt->normals()[v.f()], un); + gmm::mult(gis.ctx.B(), un, gis.Normal); + scalar_type nup = gmm::vect_norm2(gis.Normal); J *= nup; - gmm::scale(up,1.0/nup); - gmm::clean(up, 1e-13); - gis.Normal = up; + gmm::scale(gis.Normal,1.0/nup); + gmm::clean(gis.Normal, 1e-13); } else gis.Normal.resize(0); } gis.coeff = J * pai->coeff(first_ind+gis.ipt); Modified: trunk/getfem/src/gmm/gmm_vector.h URL: http://svn.gna.org/viewcvs/getfem/trunk/getfem/src/gmm/gmm_vector.h?rev=5410&r1=5409&r2=5410&view=diff ============================================================================== --- trunk/getfem/src/gmm/gmm_vector.h (original) +++ trunk/getfem/src/gmm/gmm_vector.h Wed Oct 12 14:21:01 2016 @@ -1065,8 +1065,8 @@ base_type_::resize(this->nb_stored()+1, ev); if (ind != this->nb_stored() - 1) { it = this->begin() + ind; - iterator ite = this->end(); --ite; iterator itee = ite; --itee; - for (; ite != it; --ite, --itee) *ite = *itee; + iterator ite = this->end(); --ite; iterator itee = ite; + for (; ite != it; --ite) { --itee; *ite = *itee; } *it = ev; } } @@ -1092,8 +1092,8 @@ base_type_::resize(this->nb_stored()+1, ev); if (ind != this->nb_stored() - 1) { it = this->begin() + ind; - iterator ite = this->end(); --ite; iterator itee = ite; --itee; - for (; itee != it; --ite, --itee) *ite = *itee; + iterator ite = this->end(); --ite; iterator itee = ite; + for (; ite != it; --ite) { --itee; *ite = *itee; } *it = ev; } } _______________________________________________ Getfem-commits mailing list Getfem-commits@gna.org https://mail.gna.org/listinfo/getfem-commits