From: Wolfgang Bangerth Date: Tue, 28 Apr 1998 12:30:43 +0000 (+0000) Subject: Small performance improvement. X-Git-Tag: v8.0.0~23042 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=14cfea790bcda473be1c78fb9a1bf994b19fd592;p=dealii.git Small performance improvement. git-svn-id: https://svn.dealii.org/trunk@212 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/dof_accessor.h b/deal.II/deal.II/include/dofs/dof_accessor.h index 0e0d39c042..ea259bf1e4 100644 --- a/deal.II/deal.II/include/dofs/dof_accessor.h +++ b/deal.II/deal.II/include/dofs/dof_accessor.h @@ -180,6 +180,15 @@ class DoFLineAccessor : public DoFAccessor, public BaseClass { */ void get_dof_indices (vector &dof_indices) const; + /** + * Return the length of the line. If the + * line describes part of the boundary + * and is not a straight one, ask the + * finite element class for the correct + * length! + */ + double diameter () const; + /** * Return the #i#th child as a DoF line * iterator. This function is needed since @@ -266,6 +275,18 @@ class DoFQuadAccessor : public DoFAccessor, public BaseClass { */ void get_dof_indices (vector &dof_indices) const; + /** + * Return the diameter of the quad. If the + * quad describes part of the boundary + * and is not a plane one, ask the + * finite element class for the correct + * length! + * + * The diameter of a quad is computed to + * be the larger of the two diagonals. + */ + double diameter () const; + /** * Return a pointer to the #i#th line * bounding this #Quad#. @@ -453,9 +474,8 @@ class DoFCellAccessor : public DoFSubstructAccessor { * dofs on line 0, dofs on line 1, etc, * dofs on quad 0, etc. * - * The values are appended to the end - * of the vector. It is assumed that - * the vector be empty beforehand. + * It is assumed that the vector already + * has the right size beforehand. */ void get_dof_values (const dVector &values, vector &dof_values) const; diff --git a/deal.II/deal.II/source/dofs/dof_accessor.cc b/deal.II/deal.II/source/dofs/dof_accessor.cc index 1aee0a8d59..ca220cae1a 100644 --- a/deal.II/deal.II/source/dofs/dof_accessor.cc +++ b/deal.II/deal.II/source/dofs/dof_accessor.cc @@ -101,6 +101,14 @@ DoFLineAccessor::get_dof_indices (vector &dof_indices) const +template +double +DoFLineAccessor::diameter () const { + return sqrt((vertex(1)-vertex(0)).square()); +}; + + + template TriaIterator > DoFLineAccessor::child (const unsigned int i) const { @@ -215,6 +223,15 @@ DoFQuadAccessor::get_dof_indices (vector &dof_indices) const +template +double +DoFQuadAccessor::diameter () const { + return sqrt(max((vertex(2)-vertex(0)).square(), + (vertex(3)-vertex(1)).square())); +}; + + + template TriaIterator > > @@ -316,16 +333,17 @@ DoFCellAccessor<1>::get_dof_values (const dVector &values, vector &dof_values) const { Assert (dof_handler != 0, ExcInvalidObject()); Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject()); - Assert (dof_values.size() == 0, ExcVectorNotEmpty()); + Assert (dof_values.size() == dof_handler->get_selected_fe().total_dofs, + ExcVectorDoesNotMatch()); Assert (values.n() == dof_handler->n_dofs(), ExcVectorDoesNotMatch()); - - dof_values.reserve (dof_handler->get_selected_fe().total_dofs); + + vector::iterator next_dof_value=dof_values.begin(); for (unsigned int vertex=0; vertex<2; ++vertex) for (unsigned int d=0; dget_selected_fe().dofs_per_vertex; ++d) - dof_values.push_back (values(vertex_dof_index(vertex,d))); + *next_dof_value++ = values(vertex_dof_index(vertex,d)); for (unsigned int d=0; dget_selected_fe().dofs_per_line; ++d) - dof_values.push_back (values(dof_index(d))); + *next_dof_value++ = values(dof_index(d)); }; @@ -336,19 +354,20 @@ DoFCellAccessor<2>::get_dof_values (const dVector &values, vector &dof_values) const { Assert (dof_handler != 0, ExcInvalidObject()); Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject()); - Assert (dof_values.size() == 0, ExcVectorNotEmpty()); + Assert (dof_values.size() == dof_handler->get_selected_fe().total_dofs, + ExcVectorDoesNotMatch()); Assert (values.n() == dof_handler->n_dofs(), ExcVectorDoesNotMatch()); - - dof_values.reserve (dof_handler->get_selected_fe().total_dofs); + + vector::iterator next_dof_value=dof_values.begin(); for (unsigned int vertex=0; vertex<4; ++vertex) for (unsigned int d=0; dget_selected_fe().dofs_per_vertex; ++d) - dof_values.push_back (values(vertex_dof_index(vertex,d))); + *next_dof_value++ = values(vertex_dof_index(vertex,d)); for (unsigned int line=0; line<4; ++line) for (unsigned int d=0; dget_selected_fe().dofs_per_line; ++d) - dof_values.push_back (values(this->line(line)->dof_index(d))); + *next_dof_value++ = values(this->line(line)->dof_index(d)); for (unsigned int d=0; dget_selected_fe().dofs_per_quad; ++d) - dof_values.push_back (values(dof_index(d))); + *next_dof_value++ = values(dof_index(d)); };