template <int dim>
void FEValuesBase<dim>::get_function_values (const Vector<double> &fe_function,
- vector< vector<double> > &values) const
+ vector< Vector<double> > &values) const
{
Assert (n_quadrature_points == values.size(),
ExcWrongNoOfComponents());
// functions
for (unsigned int point=0; point<n_quadrature_points; ++point)
for (unsigned int shape_func=0; shape_func<total_dofs; ++shape_func)
- values[point][fe->system_to_component_index(shape_func).first]
+ values[point](fe->system_to_component_index(shape_func).first)
+= (dof_values(shape_func) * shape_values[selected_dataset](shape_func, point));
};
FEValues<dim> fe(dofs->get_fe(), points, UpdateFlags(update_q_points),
dofs->get_tria().get_boundary());
- vector< vector <vector<double> > >
+ vector< vector <Vector<double> > >
values (dof_data.size(),
- vector< vector<double> >(points.n_quadrature_points,
- vector<double>(dofs->get_fe().n_components
+ vector< Vector<double> >(points.n_quadrature_points,
+ Vector<double>(dofs->get_fe().n_components
)));
unsigned int cell_index=0;
out << pt << " ";
for (unsigned int i=0; i!=dof_data.size(); ++i)
for (unsigned int j=0; j < dofs->get_fe().n_components; ++j)
- out << values[i][supp_pt][j]
+ out << values[i][supp_pt](j)
<< ' ';
for (unsigned int i=0; i<cell_data.size(); ++i)
out << (*cell_data[i].data)(cell_index)
for (unsigned int i=0; i!=dof_data.size(); ++i)
for (unsigned int j=0; j < dofs->get_fe().n_components; ++j)
- out << values[i][supp_pt][j]
+ out << values[i][supp_pt](j)
<< ' ';
for (unsigned int i=0; i<cell_data.size(); ++i)
out << (*cell_data[i].data)(cell_index)
for (unsigned int i=0; i!=dof_data.size(); ++i)
for (unsigned int j=0; j < dofs->get_fe().n_components; ++j)
- out << values[i][supp_pt][j]
+ out << values[i][supp_pt](j)
<< ' ';
for (unsigned int i=0; i<cell_data.size(); ++i)
out << (*cell_data[i].data)(cell_index)
// fe_function
if (true)
{
- vector< vector<double> > function_values (n_q_points,
- vector<double>(fe.n_components));
+ vector< Vector<double> > function_values (n_q_points,
+ Vector<double>(fe.n_components));
fe_values.get_function_values (fe_function, function_values);
/* transform (psi.begin(), psi.end(),