*/
void get_dof_indices (vector<int> &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
*/
void get_dof_indices (vector<int> &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#.
* 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<double> &dof_values) const;
+template <int dim, class BaseClass>
+double
+DoFLineAccessor<dim,BaseClass>::diameter () const {
+ return sqrt((vertex(1)-vertex(0)).square());
+};
+
+
+
template <int dim, class BaseClass>
TriaIterator<dim,DoFLineAccessor<dim,BaseClass> >
DoFLineAccessor<dim,BaseClass>::child (const unsigned int i) const {
+template <int dim, class BaseClass>
+double
+DoFQuadAccessor<dim,BaseClass>::diameter () const {
+ return sqrt(max((vertex(2)-vertex(0)).square(),
+ (vertex(3)-vertex(1)).square()));
+};
+
+
+
template <int dim, class BaseClass>
TriaIterator<dim,DoFLineAccessor<dim,LineAccessor<dim> > >
vector<double> &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<double>::iterator next_dof_value=dof_values.begin();
for (unsigned int vertex=0; vertex<2; ++vertex)
for (unsigned int d=0; d<dof_handler->get_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; d<dof_handler->get_selected_fe().dofs_per_line; ++d)
- dof_values.push_back (values(dof_index(d)));
+ *next_dof_value++ = values(dof_index(d));
};
vector<double> &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<double>::iterator next_dof_value=dof_values.begin();
for (unsigned int vertex=0; vertex<4; ++vertex)
for (unsigned int d=0; d<dof_handler->get_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; d<dof_handler->get_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; d<dof_handler->get_selected_fe().dofs_per_quad; ++d)
- dof_values.push_back (values(dof_index(d)));
+ *next_dof_value++ = values(dof_index(d));
};