Nice catch! I forgot to mention it is bilinear. But thanks for all the
explanation.
在 2018年3月28日星期三 UTC-7下午12:42:27，Wolfgang Bangerth写道：
>
>
> > |
> > // The vector I want to output is the solution for vector-valued
> problem.
> > // I want to output the cell-average value for the 0th component.
> > Vector<double>modi (solu);
> > Quadrature<dim>qq(fe.get_unit_support_points());
> > FEValues<dim>fv(fe,qq,update_values |update_quadrature_points);
> > std::vector<bool>selected_dofs (dof_handler.n_dofs());
> >
> > // comp[0] is the 0th component of
> > std::vector<FEValuesExtractors::Scalar> for the system
> >
> DoFTools::extract_dofs(dof_handler,fe.component_mask(comp[0]),selected_dofs);
>
> > std::vector<types::global_dof_index>local_to_global_dof
> (fe.dofs_per_cell);
> >
> >
> > // modify the 0th component
> > for(typenameDoFHandler<dim>::active_cell_iterator
> > cell=dof_handler.begin_active();cell!=dof_handler.end();++cell){
> > fv.reinit(cell);
> > std::vector<double>local_sol(qq.size());
> > fv[comp[0]].get_function_values(solu_out,local_sol);
> > doubleavg
> > =std::accumulate(local_sol.begin(),local_sol.end(),0.)/local_sol.size();
> > cell->get_dof_indices (local_to_global_dof);
> > for(autoi=0;i<fe.dofs_per_cell;++i){
> > if(selected_dofs[local_to_global_dof[i]]){
> > solu_out[local_to_global_dof[i]]=avg;
> > }
> > }
> > }
> > |
>
> This code gives an approximation in that it takes the average of the
> values of the degrees of freedom on each cell, but that is not equal to
> the spatial average
> 1/|K| \int_K u_h(x) dx
> that you probably want to compute. That's because if you expand u_h on
> cell K into
> u_h(x) = \sum_i U_i \phi_i(x)
> then the *real* average is given by
>
> \sum_i U_i (1/|K| \int_K phi_i(x) dx)
>
> but in general you do not have that
>
> 1/|K| \int_K phi_i(x) dx == const = 1/dofs_per_cell
>
> (This is true for Q1 elements, but not for higher order elements.)
> Consequently, you probably want to use a proper integration formula.
> Maybe something like this:
>
>
> QGauss<dim> quadrature(degree + 2);
> FEValues<dim> fe_values (...);
> Vector<double> cellwise_averages (triangulation.n_active_cells());
>
> for (auto cell : ...)
> {
> fe_values.reinit (cell);
> fe_values.get_function_values (solution, local_values);
> double average = 0;
> double cell_volume = 0;
> for (q=0.....)
> {
> average += local_values[q] * fe_values.JxW(q);
> cell_volume += fe_values.JxW(q);
> }
> average /= cell_volume;
>
> cellwise_averages(cell->active_cell_index()) = average;
> }
>
> You can then attach this vector with as many entries as there are cells
> to a DataOut and output it, or do whatever else you need to do with it!
>
> Best
> W.
>
>
>
>
>
> --
> ------------------------------------------------------------------------
> Wolfgang Bangerth email: bang...@colostate.edu
> <javascript:>
> www: http://www.math.colostate.edu/~bangerth/
>
--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see
https://groups.google.com/d/forum/dealii?hl=en
---
You received this message because you are subscribed to the Google Groups
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to dealii+unsubscr...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.
