Use unsigned integers for the colnums array in dSMatrixStruct. This
would enhance safety since colnum=-1 would no longer point to a
- valid address.
+ valid address. How do you mark non-used columns? (gk)
Let dSMatrixStruct::compress free the memory of colnums which is
no longer needed.
-FeValues: Replace 'ansatz', add flexibility for update flags and put
-the element into the constructor
+FeValues: add flexibility for update flags
+Remove all fe& in vectors.h
+
+Collection of transfer matrices: how to design?
#include <vector>
#include <map>
#include <base/exceptions.h>
-
-
-
+#include <base/smartpointer.h>
template <int dim> class TriaAccessor;
template <int dim> class LineAccessor;
DoFHandler (Triangulation<dim> *tria);
/**
- * Destructor
+ * Destructor.
*/
virtual ~DoFHandler ();
* function is inline, so it should
* be reasonably fast.
*/
- const FiniteElementBase<dim> & get_selected_fe () const;
+ const FiniteElementBase<dim> & get_fe () const;
/**
* Return a constant reference to the
* dofs per line, but also restriction
* and prolongation matrices, etc.
*/
- const FiniteElementBase<dim> *selected_fe;
+ SmartPointer<const FiniteElementBase<dim> > selected_fe;
private:
/**
template <int dim>
inline
-const FiniteElementBase<dim> & DoFHandler<dim>::get_selected_fe () const {
+const FiniteElementBase<dim> & DoFHandler<dim>::get_fe () const {
return *selected_fe;
};
update_JxW_values = 8,
/**
* Compute the points on the real cell
- * on which the ansatz functions are
+ * on which the trial functions are
* located.
*
* Giving this flag to the
* #FESubfaceValues# class will result
- * in an error, since ansatz points are
+ * in an error, since support points are
* not useful in that case.
*/
- update_ansatz_points = 16,
+ update_support_points = 16,
/**
* Update the outward normal vectors
* to the face relative to this cell.
* created the value of the function given as argument. This is identical
* to saying that the resulting finite element function (which is isomorphic
* to the output vector) has exact function values in all off-points of
- * ansatz functions. The off-point of an ansatz function is the point where
- * it assumes its nominal value, e.g. for linear ansatz functions the
+ * trial functions. The off-point of an trial function is the point where
+ * it assumes its nominal value, e.g. for linear trial functions the
* off-points are th corners of an element. This function therefore relies
* on the assumption that a finite element is used for which the degrees
* of freedom are function values (Lagrange elements) rather than gradients,
*
* It seems inevitable that some values of the vector to be created are set
* twice or even more than that. The reason is that we have to loop over
- * all cells and get the function values for each of the ansatz functions
+ * all cells and get the function values for each of the trial functions
* located thereon. This applies also to the functions located on faces and
* corners which we thus visit more than once. While setting the value
* in the vector is not an expensive operation, the evaluation of the
/**
* Compute the interpolation of
- * #function# at the ansatz points to
+ * #function# at the support points to
* the finite element space.
*
* See the general documentation of this
const Function<dim> &function,
dVector &vec);
+ /**
+ * Interpolate different finite element
+ * spaces. The interpolation is
+ * executed from the higher order
+ * space represented by #high# to
+ the lower order space
+ #low#. The interpolation on each
+ cell is represented by the matrix
+ #transfer#. Curved boundaries are
+ neglected so far.
+ */
+ static void interpolate(const DoFHandler<dim> &high_dof,
+ const DoFHandler<dim> &low_dof,
+ const dFMatrix& transfer,
+ const dVector& high,
+ dVector& low);
+
+
/**
* Compute the projection of
* #function# to the finite element space.
DoFLineAccessor<dim,BaseClass>::get_dof_indices (vector<int> &dof_indices) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (dof_handler->selected_fe != 0, ExcInvalidObject());
- Assert (dof_indices.size() == (2*dof_handler->get_selected_fe().dofs_per_vertex +
- dof_handler->get_selected_fe().dofs_per_line),
+ Assert (dof_indices.size() == (2*dof_handler->get_fe().dofs_per_vertex +
+ dof_handler->get_fe().dofs_per_line),
ExcVectorDoesNotMatch());
- const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_selected_fe().dofs_per_line;
+ const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
+ dofs_per_line = dof_handler->get_fe().dofs_per_line;
vector<int>::iterator next = dof_indices.begin();
for (unsigned int vertex=0; vertex<2; ++vertex)
for (unsigned int d=0; d<dofs_per_vertex; ++d)
dVector &global_destination) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (dof_handler->selected_fe != 0, ExcInvalidObject());
- Assert (local_source.size() == (2*dof_handler->get_selected_fe().dofs_per_vertex +
- dof_handler->get_selected_fe().dofs_per_line),
+ Assert (local_source.size() == (2*dof_handler->get_fe().dofs_per_vertex +
+ dof_handler->get_fe().dofs_per_line),
ExcVectorDoesNotMatch());
Assert (dof_handler->n_dofs() == global_destination.size(),
ExcVectorDoesNotMatch());
dSMatrix &global_destination) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (dof_handler->selected_fe != 0, ExcInvalidObject());
- Assert (local_source.m() == (2*dof_handler->get_selected_fe().dofs_per_vertex +
- dof_handler->get_selected_fe().dofs_per_line),
+ Assert (local_source.m() == (2*dof_handler->get_fe().dofs_per_vertex +
+ dof_handler->get_fe().dofs_per_line),
ExcVectorDoesNotMatch());
Assert (local_source.m() == local_source.n(),
ExcMatrixDoesNotMatch());
DoFQuadAccessor<dim,BaseClass>::get_dof_indices (vector<int> &dof_indices) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (dof_handler->selected_fe != 0, ExcInvalidObject());
- 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),
+ Assert (dof_indices.size() == (4*dof_handler->get_fe().dofs_per_vertex +
+ 4*dof_handler->get_fe().dofs_per_line +
+ dof_handler->get_fe().dofs_per_quad),
ExcVectorDoesNotMatch());
- const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_selected_fe().dofs_per_line,
- dofs_per_quad = dof_handler->get_selected_fe().dofs_per_quad;
+ const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
+ dofs_per_line = dof_handler->get_fe().dofs_per_line,
+ dofs_per_quad = dof_handler->get_fe().dofs_per_quad;
vector<int>::iterator next = dof_indices.begin();
for (unsigned int vertex=0; vertex<4; ++vertex)
for (unsigned int d=0; d<dofs_per_vertex; ++d)
dVector &global_destination) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (dof_handler->selected_fe != 0, ExcInvalidObject());
- Assert (local_source.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),
+ Assert (local_source.size() == (4*dof_handler->get_fe().dofs_per_vertex +
+ 4*dof_handler->get_fe().dofs_per_line +
+ dof_handler->get_fe().dofs_per_quad),
ExcVectorDoesNotMatch());
Assert (dof_handler->n_dofs() == global_destination.size(),
ExcVectorDoesNotMatch());
dSMatrix &global_destination) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (dof_handler->selected_fe != 0, ExcInvalidObject());
- Assert (local_source.m() == (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),
+ Assert (local_source.m() == (4*dof_handler->get_fe().dofs_per_vertex +
+ 4*dof_handler->get_fe().dofs_per_line +
+ dof_handler->get_fe().dofs_per_quad),
ExcMatrixDoesNotMatch());
Assert (local_source.m() == local_source.n(),
ExcMatrixDoesNotMatch());
DoFCellAccessor<1>::get_dof_values (const dVector &values,
dVector &dof_values) const {
Assert (dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (dof_values.size() == dof_handler->get_selected_fe().total_dofs,
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (dof_values.size() == dof_handler->get_fe().total_dofs,
ExcVectorDoesNotMatch());
Assert (values.size() == dof_handler->n_dofs(),
ExcVectorDoesNotMatch());
- const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_selected_fe().dofs_per_line;
+ const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
+ dofs_per_line = dof_handler->get_fe().dofs_per_line;
vector<double>::iterator next_dof_value=dof_values.begin();
for (unsigned int vertex=0; vertex<2; ++vertex)
for (unsigned int d=0; d<dofs_per_vertex; ++d)
DoFCellAccessor<2>::get_dof_values (const dVector &values,
dVector &dof_values) const {
Assert (dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (dof_values.size() == dof_handler->get_selected_fe().total_dofs,
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (dof_values.size() == dof_handler->get_fe().total_dofs,
ExcVectorDoesNotMatch());
Assert (values.size() == dof_handler->n_dofs(),
ExcVectorDoesNotMatch());
- const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_selected_fe().dofs_per_line,
- dofs_per_quad = dof_handler->get_selected_fe().dofs_per_quad;
+ const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
+ dofs_per_line = dof_handler->get_fe().dofs_per_line,
+ dofs_per_quad = dof_handler->get_fe().dofs_per_quad;
vector<double>::iterator next_dof_value=dof_values.begin();
for (unsigned int vertex=0; vertex<4; ++vertex)
for (unsigned int d=0; d<dofs_per_vertex; ++d)
template <int dim>
DoFHandler<dim>::DoFHandler (Triangulation<dim> *tria) :
tria(tria),
- selected_fe (0),
used_dofs (0) {};
template <int dim>
-DoFHandler<dim>::~DoFHandler () {
- if (selected_fe != 0)
- delete selected_fe;
-};
+DoFHandler<dim>::~DoFHandler ()
+{}
template <int dim>
-void DoFHandler<dim>::distribute_dofs (const FiniteElementBase<dim> &fe) {
+void DoFHandler<dim>::distribute_dofs (const FiniteElementBase<dim> &ff) {
Assert (tria->n_levels() > 0, ExcInvalidTriangulation());
- if (selected_fe != 0) delete selected_fe;
- selected_fe = new FiniteElementBase<dim>(fe);
+ selected_fe = &ff;
+
reserve_space ();
// clear user flags because we will
Assert (mg_dof_handler != 0, ExcInvalidObject());
// make sure a FE has been selected
// and enough room was reserved
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (i<dof_handler->get_selected_fe().dofs_per_line,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_line));
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (i<dof_handler->get_fe().dofs_per_line,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_line));
return mg_dof_handler->mg_levels[present_level]
- ->line_dofs[present_index*dof_handler->get_selected_fe().dofs_per_line+i];
+ ->line_dofs[present_index*dof_handler->get_fe().dofs_per_line+i];
};
Assert (mg_dof_handler != 0, ExcInvalidObject());
// make sure a FE has been selected
// and enough room was reserved
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (i<dof_handler->get_selected_fe().dofs_per_line,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_line));
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (i<dof_handler->get_fe().dofs_per_line,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_line));
mg_dof_handler->mg_levels[present_level]
- ->line_dofs[present_index*dof_handler->get_selected_fe().dofs_per_line+i] = index;
+ ->line_dofs[present_index*dof_handler->get_fe().dofs_per_line+i] = index;
};
const unsigned int i) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
Assert (vertex<2, ExcInvalidIndex (i,0,2));
- Assert (i<dof_handler->get_selected_fe().dofs_per_vertex,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_vertex));
+ Assert (i<dof_handler->get_fe().dofs_per_vertex,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_vertex));
return (mg_dof_handler->mg_vertex_dofs[vertex_index(vertex)]
- .get_index (present_level, i, dof_handler->get_selected_fe().dofs_per_vertex));
+ .get_index (present_level, i, dof_handler->get_fe().dofs_per_vertex));
};
const int index) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
Assert (vertex<2, ExcInvalidIndex (i,0,2));
- Assert (i<dof_handler->get_selected_fe().dofs_per_vertex,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_vertex));
+ Assert (i<dof_handler->get_fe().dofs_per_vertex,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_vertex));
mg_dof_handler->mg_vertex_dofs[vertex_index(vertex)]
- .set_index (present_level, i, dof_handler->get_selected_fe().dofs_per_vertex, index);
+ .set_index (present_level, i, dof_handler->get_fe().dofs_per_vertex, index);
};
MGDoFLineAccessor<dim,BaseClass>::get_mg_dof_indices (vector<int> &dof_indices) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (dof_indices.size() == (2*dof_handler->get_selected_fe().dofs_per_vertex +
- dof_handler->get_selected_fe().dofs_per_line),
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (dof_indices.size() == (2*dof_handler->get_fe().dofs_per_vertex +
+ dof_handler->get_fe().dofs_per_line),
ExcVectorDoesNotMatch());
- const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_selected_fe().dofs_per_line;
+ const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
+ dofs_per_line = dof_handler->get_fe().dofs_per_line;
vector<int>::iterator next = dof_indices.begin();
for (unsigned int vertex=0; vertex<2; ++vertex)
for (unsigned int d=0; d<dofs_per_vertex; ++d)
Assert (mg_dof_handler != 0, ExcInvalidObject());
// make sure a FE has been selected
// and enough room was reserved
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (i<dof_handler->get_selected_fe().dofs_per_quad,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_quad));
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (i<dof_handler->get_fe().dofs_per_quad,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_quad));
return mg_dof_handler->mg_levels[present_level]
- ->quad_dofs[present_index*dof_handler->get_selected_fe().dofs_per_quad+i];
+ ->quad_dofs[present_index*dof_handler->get_fe().dofs_per_quad+i];
};
Assert (mg_dof_handler != 0, ExcInvalidObject());
// make sure a FE has been selected
// and enough room was reserved
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (i<dof_handler->get_selected_fe().dofs_per_quad,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_quad));
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (i<dof_handler->get_fe().dofs_per_quad,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_quad));
mg_dof_handler->mg_levels[present_level]
- ->quad_dofs[present_index*dof_handler->get_selected_fe().dofs_per_quad+i] = index;
+ ->quad_dofs[present_index*dof_handler->get_fe().dofs_per_quad+i] = index;
};
const unsigned int i) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
Assert (vertex<4, ExcInvalidIndex (i,0,4));
- Assert (i<dof_handler->get_selected_fe().dofs_per_vertex,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_vertex));
+ Assert (i<dof_handler->get_fe().dofs_per_vertex,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_vertex));
return (mg_dof_handler->mg_vertex_dofs[vertex_index(vertex)]
- .get_index (present_level, i, dof_handler->get_selected_fe().dofs_per_vertex));
+ .get_index (present_level, i, dof_handler->get_fe().dofs_per_vertex));
};
const int index) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
Assert (vertex<4, ExcInvalidIndex (i,0,4));
- Assert (i<dof_handler->get_selected_fe().dofs_per_vertex,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_vertex));
+ Assert (i<dof_handler->get_fe().dofs_per_vertex,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_vertex));
mg_dof_handler->mg_vertex_dofs[vertex_index(vertex)]
- .set_index (present_level, i, dof_handler->get_selected_fe().dofs_per_vertex, index);
+ .set_index (present_level, i, dof_handler->get_fe().dofs_per_vertex, index);
};
MGDoFQuadAccessor<dim,BaseClass>::get_mg_dof_indices (vector<int> &dof_indices) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- 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),
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (dof_indices.size() == (4*dof_handler->get_fe().dofs_per_vertex +
+ 4*dof_handler->get_fe().dofs_per_line +
+ dof_handler->get_fe().dofs_per_quad),
ExcVectorDoesNotMatch());
- const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_selected_fe().dofs_per_line,
- dofs_per_quad = dof_handler->get_selected_fe().dofs_per_quad;
+ const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
+ dofs_per_line = dof_handler->get_fe().dofs_per_line,
+ dofs_per_quad = dof_handler->get_fe().dofs_per_quad;
vector<int>::iterator next = dof_indices.begin();
for (unsigned int vertex=0; vertex<4; ++vertex)
for (unsigned int d=0; d<dofs_per_vertex; ++d)
dVector &dof_values) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (dof_values.size() == dof_handler->get_selected_fe().total_dofs,
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (dof_values.size() == dof_handler->get_fe().total_dofs,
ExcVectorDoesNotMatch());
Assert (values.size() == dof_handler->n_dofs(),
ExcVectorDoesNotMatch());
- const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_selected_fe().dofs_per_line;
+ const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
+ dofs_per_line = dof_handler->get_fe().dofs_per_line;
vector<double>::iterator next_dof_value=dof_values.begin();
for (unsigned int vertex=0; vertex<2; ++vertex)
for (unsigned int d=0; d<dofs_per_vertex; ++d)
dVector &dof_values) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (dof_values.size() == dof_handler->get_selected_fe().total_dofs,
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (dof_values.size() == dof_handler->get_fe().total_dofs,
ExcVectorDoesNotMatch());
Assert (values.size() == dof_handler->n_dofs(),
ExcVectorDoesNotMatch());
- const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_selected_fe().dofs_per_line,
- dofs_per_quad = dof_handler->get_selected_fe().dofs_per_quad;
+ const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
+ dofs_per_line = dof_handler->get_fe().dofs_per_line,
+ dofs_per_quad = dof_handler->get_fe().dofs_per_quad;
vector<double>::iterator next_dof_value=dof_values.begin();
for (unsigned int vertex=0; vertex<4; ++vertex)
for (unsigned int d=0; d<dofs_per_vertex; ++d)
template <int dim>
FEValuesBase<dim>::FEValuesBase (const unsigned int n_q_points,
- const unsigned int n_ansatz_points,
+ const unsigned int n_support_points,
const unsigned int n_dofs,
const unsigned int n_transform_functions,
const unsigned int n_values_arrays,
weights (n_q_points, 0),
JxW_values (n_q_points, 0),
quadrature_points (n_q_points, Point<dim>()),
- ansatz_points (n_ansatz_points, Point<dim>()),
+ support_points (n_support_points, Point<dim>()),
jacobi_matrices (n_q_points, dFMatrix(dim,dim)),
shape_values_transform (n_values_arrays,
dFMatrix(n_transform_functions,
// initialize with zero
fill_n (values.begin(), n_quadrature_points, 0);
- // add up contributions of ansatz
+ // add up contributions of trial
// functions
for (unsigned int point=0; point<n_quadrature_points; ++point)
for (unsigned int shape_func=0; shape_func<total_dofs; ++shape_func)
// initialize with zero
fill_n (gradients.begin(), n_quadrature_points, Point<dim>());
- // add up contributions of ansatz
+ // add up contributions of trial
// functions
for (unsigned int point=0; point<n_quadrature_points; ++point)
for (unsigned int shape_func=0; shape_func<total_dofs; ++shape_func)
template <int dim>
-const Point<dim> & FEValuesBase<dim>::ansatz_point (const unsigned int i) const {
- Assert (i<ansatz_points.size(), ExcInvalidIndex(i, ansatz_points.size()));
- Assert (update_flags & update_ansatz_points, ExcAccessToUninitializedField());
+const Point<dim> & FEValuesBase<dim>::support_point (const unsigned int i) const {
+ Assert (i<support_points.size(), ExcInvalidIndex(i, support_points.size()));
+ Assert (update_flags & update_support_points, ExcAccessToUninitializedField());
- return ansatz_points[i];
+ return support_points[i];
};
(update_flags & update_JxW_values)||
(update_flags & update_q_points) ||
(update_flags & update_gradients) ||
- (update_flags & update_ansatz_points))
+ (update_flags & update_support_points))
fe->fill_fe_values (cell,
unit_quadrature_points,
jacobi_matrices,
update_flags & (update_jacobians |
update_JxW_values |
update_gradients),
- ansatz_points,
- update_flags & update_ansatz_points,
+ support_points,
+ update_flags & update_support_points,
quadrature_points,
update_flags & update_q_points,
shape_values_transform[0], unit_shape_gradients_transform,
template <int dim>
FEFaceValuesBase<dim>::FEFaceValuesBase (const unsigned int n_q_points,
- const unsigned int n_ansatz_points,
+ const unsigned int n_support_points,
const unsigned int n_dofs,
const unsigned int n_transform_functions,
const unsigned int n_faces_or_subfaces,
const UpdateFlags update_flags) :
FEValuesBase<dim> (n_q_points,
- n_ansatz_points,
+ n_support_points,
n_dofs,
n_transform_functions,
n_faces_or_subfaces,
(update_flags & update_JxW_values)||
(update_flags & update_q_points) ||
(update_flags & update_gradients) ||
- (update_flags & update_ansatz_points) ||
+ (update_flags & update_support_points) ||
(update_flags & update_JxW_values))
fe->fill_fe_face_values (cell,
face_no,
update_flags & (update_jacobians |
update_gradients |
update_JxW_values),
- ansatz_points,
- update_flags & update_ansatz_points,
+ support_points,
+ update_flags & update_support_points,
quadrature_points,
update_flags & update_q_points,
face_jacobi_determinants,
update_flags),
fe(&fe)
{
- Assert ((update_flags & update_ansatz_points) == false,
+ Assert ((update_flags & update_support_points) == false,
ExcInvalidUpdateFlag());
unit_face_quadrature_points = quadrature.get_quad_points();
Assert (mg_dof_handler != 0, ExcInvalidObject());
// make sure a FE has been selected
// and enough room was reserved
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (i<dof_handler->get_selected_fe().dofs_per_line,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_line));
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (i<dof_handler->get_fe().dofs_per_line,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_line));
return mg_dof_handler->mg_levels[present_level]
- ->line_dofs[present_index*dof_handler->get_selected_fe().dofs_per_line+i];
+ ->line_dofs[present_index*dof_handler->get_fe().dofs_per_line+i];
};
Assert (mg_dof_handler != 0, ExcInvalidObject());
// make sure a FE has been selected
// and enough room was reserved
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (i<dof_handler->get_selected_fe().dofs_per_line,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_line));
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (i<dof_handler->get_fe().dofs_per_line,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_line));
mg_dof_handler->mg_levels[present_level]
- ->line_dofs[present_index*dof_handler->get_selected_fe().dofs_per_line+i] = index;
+ ->line_dofs[present_index*dof_handler->get_fe().dofs_per_line+i] = index;
};
const unsigned int i) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
Assert (vertex<2, ExcInvalidIndex (i,0,2));
- Assert (i<dof_handler->get_selected_fe().dofs_per_vertex,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_vertex));
+ Assert (i<dof_handler->get_fe().dofs_per_vertex,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_vertex));
return (mg_dof_handler->mg_vertex_dofs[vertex_index(vertex)]
- .get_index (present_level, i, dof_handler->get_selected_fe().dofs_per_vertex));
+ .get_index (present_level, i, dof_handler->get_fe().dofs_per_vertex));
};
const int index) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
Assert (vertex<2, ExcInvalidIndex (i,0,2));
- Assert (i<dof_handler->get_selected_fe().dofs_per_vertex,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_vertex));
+ Assert (i<dof_handler->get_fe().dofs_per_vertex,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_vertex));
mg_dof_handler->mg_vertex_dofs[vertex_index(vertex)]
- .set_index (present_level, i, dof_handler->get_selected_fe().dofs_per_vertex, index);
+ .set_index (present_level, i, dof_handler->get_fe().dofs_per_vertex, index);
};
MGDoFLineAccessor<dim,BaseClass>::get_mg_dof_indices (vector<int> &dof_indices) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (dof_indices.size() == (2*dof_handler->get_selected_fe().dofs_per_vertex +
- dof_handler->get_selected_fe().dofs_per_line),
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (dof_indices.size() == (2*dof_handler->get_fe().dofs_per_vertex +
+ dof_handler->get_fe().dofs_per_line),
ExcVectorDoesNotMatch());
- const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_selected_fe().dofs_per_line;
+ const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
+ dofs_per_line = dof_handler->get_fe().dofs_per_line;
vector<int>::iterator next = dof_indices.begin();
for (unsigned int vertex=0; vertex<2; ++vertex)
for (unsigned int d=0; d<dofs_per_vertex; ++d)
Assert (mg_dof_handler != 0, ExcInvalidObject());
// make sure a FE has been selected
// and enough room was reserved
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (i<dof_handler->get_selected_fe().dofs_per_quad,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_quad));
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (i<dof_handler->get_fe().dofs_per_quad,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_quad));
return mg_dof_handler->mg_levels[present_level]
- ->quad_dofs[present_index*dof_handler->get_selected_fe().dofs_per_quad+i];
+ ->quad_dofs[present_index*dof_handler->get_fe().dofs_per_quad+i];
};
Assert (mg_dof_handler != 0, ExcInvalidObject());
// make sure a FE has been selected
// and enough room was reserved
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (i<dof_handler->get_selected_fe().dofs_per_quad,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_quad));
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (i<dof_handler->get_fe().dofs_per_quad,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_quad));
mg_dof_handler->mg_levels[present_level]
- ->quad_dofs[present_index*dof_handler->get_selected_fe().dofs_per_quad+i] = index;
+ ->quad_dofs[present_index*dof_handler->get_fe().dofs_per_quad+i] = index;
};
const unsigned int i) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
Assert (vertex<4, ExcInvalidIndex (i,0,4));
- Assert (i<dof_handler->get_selected_fe().dofs_per_vertex,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_vertex));
+ Assert (i<dof_handler->get_fe().dofs_per_vertex,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_vertex));
return (mg_dof_handler->mg_vertex_dofs[vertex_index(vertex)]
- .get_index (present_level, i, dof_handler->get_selected_fe().dofs_per_vertex));
+ .get_index (present_level, i, dof_handler->get_fe().dofs_per_vertex));
};
const int index) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
Assert (vertex<4, ExcInvalidIndex (i,0,4));
- Assert (i<dof_handler->get_selected_fe().dofs_per_vertex,
- ExcInvalidIndex (i, 0, dof_handler->get_selected_fe().dofs_per_vertex));
+ Assert (i<dof_handler->get_fe().dofs_per_vertex,
+ ExcInvalidIndex (i, 0, dof_handler->get_fe().dofs_per_vertex));
mg_dof_handler->mg_vertex_dofs[vertex_index(vertex)]
- .set_index (present_level, i, dof_handler->get_selected_fe().dofs_per_vertex, index);
+ .set_index (present_level, i, dof_handler->get_fe().dofs_per_vertex, index);
};
MGDoFQuadAccessor<dim,BaseClass>::get_mg_dof_indices (vector<int> &dof_indices) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- 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),
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (dof_indices.size() == (4*dof_handler->get_fe().dofs_per_vertex +
+ 4*dof_handler->get_fe().dofs_per_line +
+ dof_handler->get_fe().dofs_per_quad),
ExcVectorDoesNotMatch());
- const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_selected_fe().dofs_per_line,
- dofs_per_quad = dof_handler->get_selected_fe().dofs_per_quad;
+ const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
+ dofs_per_line = dof_handler->get_fe().dofs_per_line,
+ dofs_per_quad = dof_handler->get_fe().dofs_per_quad;
vector<int>::iterator next = dof_indices.begin();
for (unsigned int vertex=0; vertex<4; ++vertex)
for (unsigned int d=0; d<dofs_per_vertex; ++d)
dVector &dof_values) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (dof_values.size() == dof_handler->get_selected_fe().total_dofs,
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (dof_values.size() == dof_handler->get_fe().total_dofs,
ExcVectorDoesNotMatch());
Assert (values.size() == dof_handler->n_dofs(),
ExcVectorDoesNotMatch());
- const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_selected_fe().dofs_per_line;
+ const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
+ dofs_per_line = dof_handler->get_fe().dofs_per_line;
vector<double>::iterator next_dof_value=dof_values.begin();
for (unsigned int vertex=0; vertex<2; ++vertex)
for (unsigned int d=0; d<dofs_per_vertex; ++d)
dVector &dof_values) const {
Assert (dof_handler != 0, ExcInvalidObject());
Assert (mg_dof_handler != 0, ExcInvalidObject());
- Assert (&dof_handler->get_selected_fe() != 0, ExcInvalidObject());
- Assert (dof_values.size() == dof_handler->get_selected_fe().total_dofs,
+ Assert (&dof_handler->get_fe() != 0, ExcInvalidObject());
+ Assert (dof_values.size() == dof_handler->get_fe().total_dofs,
ExcVectorDoesNotMatch());
Assert (values.size() == dof_handler->n_dofs(),
ExcVectorDoesNotMatch());
- const unsigned int dofs_per_vertex = dof_handler->get_selected_fe().dofs_per_vertex,
- dofs_per_line = dof_handler->get_selected_fe().dofs_per_line,
- dofs_per_quad = dof_handler->get_selected_fe().dofs_per_quad;
+ const unsigned int dofs_per_vertex = dof_handler->get_fe().dofs_per_vertex,
+ dofs_per_line = dof_handler->get_fe().dofs_per_line,
+ dofs_per_quad = dof_handler->get_fe().dofs_per_quad;
vector<double>::iterator next_dof_value=dof_values.begin();
for (unsigned int vertex=0; vertex<4; ++vertex)
for (unsigned int d=0; d<dofs_per_vertex; ++d)
const int index,
const AssemblerData<dim> *local_data) :
DoFCellAccessor<dim> (tria,level,index, &local_data->dof),
- cell_matrix (dof_handler->get_selected_fe().total_dofs),
- cell_vector (dVector(dof_handler->get_selected_fe().total_dofs)),
+ cell_matrix (dof_handler->get_fe().total_dofs),
+ cell_vector (dVector(dof_handler->get_fe().total_dofs)),
assemble_matrix (local_data->assemble_matrix),
assemble_rhs (local_data->assemble_rhs),
matrix(local_data->matrix),
ExcInvalidData());
Assert (!assemble_matrix || (matrix.n() == dof_handler->n_dofs()),
ExcInvalidData());
- Assert (((AssemblerData<dim>*)local_data)->fe == dof_handler->get_selected_fe(),
+ Assert (((AssemblerData<dim>*)local_data)->fe == dof_handler->get_fe(),
ExcInvalidData());
Assert (!assemble_rhs || (rhs_vector.size()==dof_handler->n_dofs()),
ExcInvalidData());
present_index,
dof_handler),
boundary);
- const unsigned int n_dofs = dof_handler->get_selected_fe().total_dofs;
+ const unsigned int n_dofs = dof_handler->get_fe().total_dofs;
if (assemble_matrix)
cell_matrix.clear ();
template <int dim>
void DataOut<dim>::write_ucd (ostream &out) const {
Assert (dofs != 0, ExcNoDoFHandlerSelected());
- Assert (dofs->get_selected_fe().dofs_per_vertex==1,
+ Assert (dofs->get_fe().dofs_per_vertex==1,
ExcIncorrectDofsPerVertex());
Assert ((1<=dim) && (dim<=3), ExcNotImplemented());
Assert ((1<=dim) && (dim<=3), ExcNotImplemented());
if (dim==1)
- Assert (dofs->get_selected_fe().dofs_per_vertex==1,
+ Assert (dofs->get_fe().dofs_per_vertex==1,
ExcIncorrectDofsPerVertex());
// write preamble
Assert (matrix.n() == matrix.m(), ExcInternalError());
Assert (matrix.n() == rhs_vector.size(), ExcInternalError());
Assert (rhs.size() != 0, ExcInternalError());
- Assert (dof.get_selected_fe() == fe, ExcInternalError());
+ Assert (dof.get_fe() == fe, ExcInternalError());
Assert (dof_to_boundary_mapping.size() == dof.n_dofs(), ExcInternalError());
Assert (*max_element(dof_to_boundary_mapping.begin(),dof_to_boundary_mapping.end()) ==
(signed int)matrix.n()-1,
};
+template <int dim> void
+VectorTools<dim>::interpolate(const DoFHandler<dim> &high_dof,
+ const DoFHandler<dim> &low_dof,
+ const dFMatrix &transfer,
+ const dVector &high,
+ dVector &low)
+{
+
+}
#if deal_II_dimension == 1
-template <>
+
void VectorTools<1>::project (const DoFHandler<1> &,
const ConstraintMatrix &,
const FiniteElement<1> &,
const FiniteElement<dim> &fe,
const NormType &norm,
const Boundary<dim> &boundary) {
- Assert (fe == dof.get_selected_fe(), ExcInvalidFE());
+ Assert (fe == dof.get_fe(), ExcInvalidFE());
difference.reinit (dof.get_tria().n_active_cells());