From: wolf Date: Mon, 4 Oct 2004 21:29:55 +0000 (+0000) Subject: Introduce a field that stores which constraints exist already. This is desirable... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=1bda7ae894d4a299917a9c2331a76b3663e78507;p=dealii-svn.git Introduce a field that stores which constraints exist already. This is desirable since we would otherwise have to walk through the list of existing constraints in add_line, which is an O(N) operation that has to be repeated O(N) times. This quadratic algorithm is now gone, making the SPEC testcase a full 5% faster than previously. The main reason is that it used to walk through a lot of memory and create a horrific number of cache misses, which are now all gone. git-svn-id: https://svn.dealii.org/trunk@9689 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/dof_constraints.h b/deal.II/deal.II/include/dofs/dof_constraints.h index 8e46abc763..77c5561045 100644 --- a/deal.II/deal.II/include/dofs/dof_constraints.h +++ b/deal.II/deal.II/include/dofs/dof_constraints.h @@ -898,6 +898,53 @@ class ConstraintMatrix : public Subscriptor * much slower. */ std::vector lines; + + /** + * A list of flags that indicate + * whether there is a constraint + * line for a given degree of + * freedom index. Note that this + * class has no notion of how + * many degrees of freedom there + * really are, so if we check + * whether there is a constraint + * line for a given degree of + * freedom, then this vector may + * actually be shorter than the + * index of the DoF we check for. + * + * This field exists since when adding a + * new constraint line we have to figure + * out whether it already + * exists. Previously, we would simply + * walk the unsorted list of constraint + * lines until we either hit the end or + * found it. This algorithm is O(N) if N + * is the number of constraints, which + * makes it O(N^2) when inserting all + * constraints. For large problems with + * many constraints, this could easily + * take 5-10 per cent of the total run + * time. With this field, we can at least + * save this time when checking whether a + * new constraint line already exists. + * + * To make things worse, traversing the + * list of existing constraints requires + * reads from many different places in + * memory. Thus, in large 3d + * applications, the add_line() function + * showed up very prominently in the + * overall compute time, mainly because + * it generated a lot of cache + * misses. This should also be fixed by + * using the O(1) algorithm to access the + * fields of this array. + * + * The field is useful in a number of + * other contexts as well, though. + */ + std::vector constraint_line_exists; /** * Store whether the arrays are sorted. @@ -917,4 +964,101 @@ class ConstraintMatrix : public Subscriptor }; + +/* ---------------- template and inline functions ----------------- */ + +inline +void +ConstraintMatrix::add_line (const unsigned int line) +{ + Assert (sorted==false, ExcMatrixIsClosed()); + + + + // check whether line already + // exists; it may, but then we need + // can just quit + if ((line < constraint_line_exists.size()) + && + (constraint_line_exists[line] == true)) + return; + + // if necessary enlarge vector of + // existing entries + if (line >= constraint_line_exists.size()) + constraint_line_exists.resize (line+1, false); + constraint_line_exists[line] = true; + + // push a new line to the end of the + // list + lines.push_back (ConstraintLine()); + lines.back().line = line; +} + + + +inline +void +ConstraintMatrix::add_entry (const unsigned int line, + const unsigned int column, + const double value) +{ + Assert (sorted==false, ExcMatrixIsClosed()); + + std::vector::iterator line_ptr; + const std::vector::const_iterator start=lines.begin(); + + // the usual case is that the line where + // a value is entered is the one we + // added last, so we search backward + for (line_ptr=(lines.end()-1); line_ptr!=start; --line_ptr) + if (line_ptr->line == line) + break; + + // if the loop didn't break, then + // line_ptr must be begin(). + // we have an error if that doesn't + // point to 'line' then + Assert (line_ptr->line==line, ExcLineInexistant(line)); + // second check: the respective + // flag in the + // constraint_line_exists field + // must exist + Assert (line < constraint_line_exists.size(), ExcInternalError()); + Assert (constraint_line_exists[line] == true, ExcInternalError()); + + + // if in debug mode, check whether an + // entry for this column already + // exists and if its the same as + // the one entered at present + // + // in any case: exit the function if an + // entry for this column already exists, + // since we don't want to enter it twice + for (std::vector >::const_iterator + p=line_ptr->entries.begin(); + p != line_ptr->entries.end(); ++p) + if (p->first == column) + { + Assert (p->second == value, + ExcEntryAlreadyExists(line, column, p->second, value)); + return; + } + + line_ptr->entries.push_back (std::make_pair(column,value)); +} + + + +inline +bool +ConstraintMatrix::is_constrained (const unsigned int index) const +{ + return ((index < constraint_line_exists.size()) + && + (constraint_line_exists[index] == true)); +} + + #endif diff --git a/deal.II/deal.II/source/dofs/dof_constraints.cc b/deal.II/deal.II/source/dofs/dof_constraints.cc index d35ca48ff7..51d0598b0b 100644 --- a/deal.II/deal.II/source/dofs/dof_constraints.cc +++ b/deal.II/deal.II/source/dofs/dof_constraints.cc @@ -80,75 +80,13 @@ ConstraintMatrix::ConstraintLine::memory_consumption () const -ConstraintMatrix::ConstraintMatrix () : - lines(), - sorted(false) +ConstraintMatrix::ConstraintMatrix () + : + lines (), + sorted (false) {} -void ConstraintMatrix::add_line (const unsigned int line) -{ - Assert (sorted==false, ExcMatrixIsClosed()); - - // check whether line already exists; - // it may, but then we need to quit - for (unsigned int i=0; i!=lines.size(); ++i) - if (lines[i].line == line) - return; - - // push a new line to the end of the - // list - lines.push_back (ConstraintLine()); - lines.back().line = line; -} - - - -void -ConstraintMatrix::add_entry (const unsigned int line, - const unsigned int column, - const double value) -{ - Assert (sorted==false, ExcMatrixIsClosed()); - - std::vector::iterator line_ptr; - const std::vector::const_iterator start=lines.begin(); - // the usual case is that the line where - // a value is entered is the one we - // added last, so we search backward - for (line_ptr=(lines.end()-1); line_ptr!=start; --line_ptr) - if (line_ptr->line == line) - break; - - // if the loop didn't break, then - // line_ptr must be begin(). - // we have an error if that doesn't - // point to 'line' then - Assert (line_ptr->line==line, ExcLineInexistant(line)); - - // if in debug mode, check whether an - // entry for this column already - // exists and if its the same as - // the one entered at present - // - // in any case: exit the function if an - // entry for this column already exists, - // since we don't want to enter it twice - for (std::vector >::const_iterator - p=line_ptr->entries.begin(); - p != line_ptr->entries.end(); ++p) - if (p->first == column) - { - Assert (p->second == value, - ExcEntryAlreadyExists(line, column, p->second, value)); - return; - }; - - line_ptr->entries.push_back (std::make_pair(column,value)); -} - - - void ConstraintMatrix::add_entries (const unsigned int line, const std::vector > &col_val_pairs) @@ -341,6 +279,15 @@ void ConstraintMatrix::merge (const ConstraintMatrix &other_constraints) const bool object_was_sorted = sorted; sorted = false; + // before we even start: merge the + // two flag arrays + if (other_constraints.constraint_line_exists.size() > + constraint_line_exists.size()) + constraint_line_exists.resize (other_constraints.constraint_line_exists.size(), + false); + for (unsigned int i=0; i::iterator i = lines.begin(); i != lines.end(); i++) { @@ -514,8 +464,16 @@ void ConstraintMatrix::shift (const unsigned int offset) void ConstraintMatrix::clear () { - std::vector tmp; - lines.swap (tmp); + { + std::vector tmp; + lines.swap (tmp); + } + + { + std::vector tmp; + constraint_line_exists.swap (tmp); + } + sorted = false; } @@ -1139,36 +1097,9 @@ unsigned int ConstraintMatrix::n_constraints () const -bool ConstraintMatrix::is_constrained (const unsigned int index) const -{ - if (lines.size() == 0) - return false; - - if (sorted == true) - { - ConstraintLine index_comparison; - index_comparison.line = index; - - return std::binary_search (lines.begin (), - lines.end (), - index_comparison); - } - else - { - for (std::vector::const_iterator i=lines.begin(); - i!=lines.end(); ++i) - if (i->line == index) - return true; - - return false; - }; -} - - - bool ConstraintMatrix::is_identity_constrained (const unsigned int index) const { - if (lines.size() == 0) + if (is_constrained(index) == false) return false; if (sorted == true) @@ -1241,6 +1172,7 @@ unsigned int ConstraintMatrix::memory_consumption () const { return (MemoryConsumption::memory_consumption (lines) + + MemoryConsumption::memory_consumption (constraint_line_exists) + MemoryConsumption::memory_consumption (sorted)); } diff --git a/deal.II/doc/news/c-5.0.html b/deal.II/doc/news/c-5.0.html index 9a1094dac1..81a95a6bf9 100644 --- a/deal.II/doc/news/c-5.0.html +++ b/deal.II/doc/news/c-5.0.html @@ -343,6 +343,19 @@ inconvenience this causes.

deal.II

    +
  1. + Fixed: The ConstraintMatrix class had some + algorithms that were linear in the number of constraints. Since these + functions had to be called for each constraint, this resulted in a + quadratic behavior. To make things worse, these algorithms traversed + large memory blocks leading to a vast number of cache misses which made + them really slow. This is now fixed: the algorithm is O(1) and should + only access single elements in memory, and one 3d testcase is now a + full 5 per cent faster on about 10 minutes runtime. +
    + (WB 2004/10/04) +

    +
  2. New: The GridOut::write_gnuplot function now supports curved boundaries also for dim==3.