* Exception
*/
DeclException0 (ExcVectorNotEmpty);
+ /**
+ * Exception
+ */
+ DeclException0 (ExcVectorDoesNotMatch);
protected:
/**
* on vertex 0, dofs on vertex 1,
* dofs on line.
*
- * The dof indices 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_indices (vector<int> &dof_indices) const;
* dofs on line 0, dofs on line 1, etc,
* dofs on quad 0, etc.
*
- * The dof indices 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_indices (vector<int> &dof_indices) const;
void get_dof_values (const dVector &values,
vector<double> &dof_values) const;
- /**
- * Exception
- */
- DeclException0 (ExcVectorDoesNotMatch);
/**
* Exception
*/
DoFLineAccessor<dim,BaseClass>::get_dof_indices (vector<int> &dof_indices) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (dof_indices.size() == 0, ExcVectorNotEmpty());
-
- dof_indices.reserve (dof_handler->get_selected_fe().total_dofs);
+ Assert (dof_indices.size() == (2*dof_handler->get_selected_fe().dofs_per_vertex +
+ dof_handler->get_selected_fe().dofs_per_line),
+ ExcVectorDoesNotMatch());
+
+ vector<int>::iterator next = dof_indices.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_indices.push_back (vertex_dof_index(vertex,d));
+ *next++ = vertex_dof_index(vertex,d);
for (unsigned int d=0; d<dof_handler->get_selected_fe().dofs_per_line; ++d)
- dof_indices.push_back (dof_index(d));
+ *next++ = dof_index(d);
};
DoFQuadAccessor<dim,BaseClass>::get_dof_indices (vector<int> &dof_indices) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (dof_indices.size() == 0, ExcVectorNotEmpty());
-
- dof_indices.reserve (dof_handler->get_selected_fe().total_dofs);
+ Assert (dof_indices.size() == (4*dof_handler->get_selected_fe().dofs_per_vertex +
+ 4*dof_handler->get_selected_fe().dofs_per_line +
+ dof_handler->get_selected_fe().dofs_per_quad),
+ ExcVectorDoesNotMatch());
+
+ vector<int>::iterator next = dof_indices.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_indices.push_back (vertex_dof_index(vertex,d));
+ *next++ = 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_indices.push_back (this->line(line)->dof_index(d));
+ *next++ = this->line(line)->dof_index(d);
for (unsigned int d=0; d<dof_handler->get_selected_fe().dofs_per_quad; ++d)
- dof_indices.push_back (dof_index(d));
+ *next++ = dof_index(d);
};
if (new_cell->active() && old_cell->active())
// both cells have no children
{
- vector<int> old_dofs, new_dofs;
+ vector<int> old_dofs(selected_fe->total_dofs);
+ vector<int> new_dofs(selected_fe->total_dofs);
old_cell->get_dof_indices (old_dofs);
new_cell->get_dof_indices (new_dofs);
- Assert (old_dofs.size() == selected_fe->total_dofs,
- ExcInternalError ());
- Assert (new_dofs.size() == selected_fe->total_dofs,
- ExcInternalError ());
// copy dofs one-by-one
for (unsigned int j=0; j<old_dofs.size(); ++j)
};
// numbers of old dofs
- vector<int> old_dof_indices;
+ vector<int> old_dof_indices (selected_fe->total_dofs);
old_cell->get_dof_indices (old_dof_indices);
- Assert (old_dof_indices.size() == selected_fe->total_dofs,
- ExcInternalError ());
-
for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
{
// numbers of child dofs
- vector<int> child_dof_indices;
+ vector<int> child_dof_indices(selected_fe->total_dofs);
child[c]->get_dof_indices (child_dof_indices);
- Assert (child_dof_indices.size() == selected_fe->total_dofs,
- ExcInternalError ());
-
for (unsigned int k=0; k<selected_fe->total_dofs; ++k)
for (unsigned int j=0; j<selected_fe->total_dofs; ++j)
if (selected_fe->prolongate(c)(k,j) != 0.)
};
// numbers of new dofs
- vector<int> new_dof_indices;
+ vector<int> new_dof_indices(selected_fe->total_dofs);
new_cell->get_dof_indices(new_dof_indices);
-
- Assert (new_dof_indices.size() == selected_fe->total_dofs,
- ExcInternalError ());
for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
{
// numbers of child dofs
- vector<int> child_dof_indices;
+ vector<int> child_dof_indices (selected_fe->total_dofs);
child[c]->get_dof_indices (child_dof_indices);
- Assert (child_dof_indices.size() == selected_fe->total_dofs,
- ExcInternalError ());
-
for (unsigned int k=0; k<selected_fe->total_dofs; ++k)
for (unsigned int j=0; j<selected_fe->total_dofs; ++j)
if (selected_fe->restrict(c)(k,j) != 0.)
if (new_cell->active() && old_cell->active())
// both cells have no children
{
- vector<int> old_dofs, new_dofs;
+ vector<int> old_dofs (selected_fe->total_dofs);
+ vector<int> new_dofs (selected_fe->total_dofs);
old_cell->get_dof_indices (old_dofs);
new_cell->get_dof_indices (new_dofs);
- Assert (old_dofs.size() == selected_fe->total_dofs,
- ExcInternalError ());
- Assert (new_dofs.size() == selected_fe->total_dofs,
- ExcInternalError ());
// copy dofs one-by-one
for (unsigned int j=0; j<old_dofs.size(); ++j)
};
// numbers of old dofs
- vector<int> old_dof_indices;
+ vector<int> old_dof_indices (selected_fe->total_dofs);
old_cell->get_dof_indices (old_dof_indices);
- Assert (old_dof_indices.size() == selected_fe->total_dofs,
- ExcInternalError ());
-
for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
{
// numbers of child dofs
- vector<int> child_dof_indices;
+ vector<int> child_dof_indices (selected_fe->total_dofs);
child[c]->get_dof_indices (child_dof_indices);
- Assert (child_dof_indices.size() == selected_fe->total_dofs,
- ExcInternalError ());
-
for (unsigned int k=0; k<selected_fe->total_dofs; ++k)
for (unsigned int j=0; j<selected_fe->total_dofs; ++j)
if (selected_fe->prolongate(c)(k,j) != 0.)
};
// numbers of new dofs
- vector<int> new_dof_indices;
+ vector<int> new_dof_indices (selected_fe->total_dofs);
new_cell->get_dof_indices(new_dof_indices);
-
- Assert (new_dof_indices.size() == selected_fe->total_dofs,
- ExcInternalError ());
for (unsigned int c=0; c<GeometryInfo<dim>::children_per_cell; ++c)
{
// numbers of child dofs
- vector<int> child_dof_indices;
+ vector<int> child_dof_indices (selected_fe->total_dofs);
child[c]->get_dof_indices (child_dof_indices);
- Assert (child_dof_indices.size() == selected_fe->total_dofs,
- ExcInternalError ());
-
for (unsigned int k=0; k<selected_fe->total_dofs; ++k)
for (unsigned int j=0; j<selected_fe->total_dofs; ++j)
if (selected_fe->restrict(c)(k,j) != 0.)
active_cell_iterator cell = begin_active(),
endc = end();
unsigned int present_cell = 0;
- vector<int> dof_indices;
const unsigned int total_dofs = selected_fe->total_dofs;
+ vector<int> dof_indices (total_dofs);
for (; cell!=endc; ++cell, ++present_cell)
{
- dof_indices.resize (0);
cell->get_dof_indices (dof_indices);
for (unsigned int i=0; i<total_dofs; ++i)
{
// get indices of dofs
- vector<int> dofs;
+ vector<int> dofs (n_dofs);
get_dof_indices (dofs);
// distribute cell matrix
void ProblemBase<dim>::clear () {
tria = 0;
dof_handler = 0;
- system_sparsity.reinit (0,0,0);
+ system_sparsity.reinit (1,1,1);
system_matrix.clear ();
- right_hand_side.reinit (0);
- solution.reinit (0);
+ right_hand_side.reinit (1);
+ solution.reinit (1);
constraints.clear ();
};
// and assign that to the vector
// \psi.
const unsigned int n_q_points = q.n_quadrature_points;
- vector<double> psi;
+ vector<double> psi (n_q_points);
// in praxi: first compute
// exact solution vector
// same procedure as above, but now
// psi is a vector of gradients
const unsigned int n_q_points = q.n_quadrature_points;
- vector<Point<dim> > psi;
+ vector<Point<dim> > psi (n_q_points);
// in praxi: first compute
// exact solution vector
FunctionMap::const_iterator function_ptr;
// field to store the indices of dofs
- // initialize once to get the size right
- // for the following fields.
- vector<int> face_dofs;
- face->get_dof_indices (face_dofs);
+ vector<int> face_dofs (fe.dofs_per_face);
vector<Point<dim> > dof_locations (face_dofs.size(), Point<dim>());
- vector<double> dof_values;
+ vector<double> dof_values (fe.dofs_per_face);
for (; face!=endf; ++face)
if ((function_ptr = dirichlet_bc.find(face->boundary_indicator())) !=
// get indices, physical location and
// boundary values of dofs on this
// face
- face_dofs.erase (face_dofs.begin(), face_dofs.end());
- dof_values.erase (dof_values.begin(), dof_values.end());
face->get_dof_indices (face_dofs);
fe.get_face_ansatz_points (face, boundary, dof_locations);
function_ptr->second->value_list (dof_locations, dof_values);
if (coefficient != 0)
{
- vector<double> coefficient_values;
+ vector<double> coefficient_values (fe_values.n_quadrature_points);
coefficient->value_list (fe_values.get_quadrature_points(),
coefficient_values);
for (unsigned int point=0; point<fe_values.n_quadrature_points; ++point)
const dFMatrix &values = fe_values.get_shape_values ();
const vector<double> &weights = fe_values.get_JxW_values ();
- vector<double> rhs_values;
+ vector<double> rhs_values (fe_values.n_quadrature_points);
right_hand_side->value_list (fe_values.get_quadrature_points(), rhs_values);
if (coefficient != 0)
{
- vector<double> coefficient_values;
+ vector<double> coefficient_values (fe_values.n_quadrature_points);
coefficient->value_list (fe_values.get_quadrature_points(),
coefficient_values);
for (unsigned int point=0; point<fe_values.n_quadrature_points; ++point)
const dFMatrix &values = fe_values.get_shape_values ();
const vector<double> &weights = fe_values.get_JxW_values ();
- vector<double> rhs_values;
+ vector<double> rhs_values(fe_values.n_quadrature_points);
right_hand_side->value_list (fe_values.get_quadrature_points(), rhs_values);
for (unsigned int point=0; point<fe_values.n_quadrature_points; ++point)
const vector<vector<Point<dim> > >&gradients = fe_values.get_shape_grads ();
const dFMatrix &values = fe_values.get_shape_values ();
- vector<double> rhs_values;
+ vector<double> rhs_values(fe_values.n_quadrature_points);
const vector<double> &weights = fe_values.get_JxW_values ();
right_hand_side->value_list (fe_values.get_quadrature_points(), rhs_values);
if (coefficient != 0)
{
- vector<double> coefficient_values;
+ vector<double> coefficient_values(fe_values.n_quadrature_points);
coefficient->value_list (fe_values.get_quadrature_points(),
coefficient_values);
for (unsigned int point=0; point<fe_values.n_quadrature_points; ++point)
if (coefficient != 0)
{
- vector<double> coefficient_values;
+ vector<double> coefficient_values(fe_values.n_quadrature_points);
coefficient->value_list (fe_values.get_quadrature_points(),
coefficient_values);
for (unsigned int point=0; point<fe_values.n_quadrature_points; ++point)
const dFMatrix &values = fe_values.get_shape_values ();
const vector<double> &weights = fe_values.get_JxW_values ();
- vector<double> rhs_values;
+ vector<double> rhs_values(fe_values.n_quadrature_points);
right_hand_side->value_list (fe_values.get_quadrature_points(), rhs_values);
for (unsigned int point=0; point<fe_values.n_quadrature_points; ++point)