Strange issue with KellyErrorEstimator functionality

[Deepak GUpta]

Dear All,

I am trying to use KellyErrorEstimator of deal.II and I saw various
tutorials which use it. I ran those examples and they work well.
Recently, I tried it with my work and when I print the
estimated_error_per_cell vector, all of the values in that vector are
zero. I have been trying to find the bug for past two days but all in
vain. May be someone can comment (if they also encountered ever) on why
this error can be zero, based on the way it is implemented in deal.II.
The piece of code is as follows: (this is an hp-adaptive framework):

//Quadrature collection for FE
hp::QCollection<dim-1> quad_collection;
for (unsigned int qrule = 1; qrule <= (fem->fe_collection).size();
++qrule){
quad_collection.push_back(QGauss<dim-1>(qrule));
}

std::cout<<"No. of dofs : "<<fem->dof_handler->n_dofs()<<std::endl;
Vector<float> estimated_error_per_cell
((*(fem->triangulation)).n_active_cells());
std::cout<<"Size of error vector :
"<<estimated_error_per_cell.size()<<std::endl;
KellyErrorEstimator<dim>::estimate (*(fem->dof_handler),
quad_collection,
typename
FunctionMap<dim>::type(),
fem->solution,
estimated_error_per_cell);


dof_hander->n_dofs() gives correct output. Size of
estimated_error_per_cell is correct. fem-> solution vector is correct
(cross-checked). A minimal working example for my problem would be
similar to step-27 but then step-27 works fine.


I hope someone has seen such an issue and can help.

Best

Deepak

--
Deepak K. Gupta

PhD Candidate/Researcher
Structural Optimization & Mechanics
TU Delft / Precision & Microsystems Engineering
Faculty of Mechanical, Maritime & Materials Engineering (3mE)
Mekelweg 2, 2628 CD Delft, The Netherlands

T +31 (0)15 27 86818
D.K....@tudelft.nl

Room 3mE.G.1.150

[Jean-Paul Pelteret]

Dear Deepak,

Have you run your code in debug mode to confirm that no errors are thrown? I have used the Kelly error estimator in the context of an hp problem with no issues. Have you tried increasing the number of face quadrature points (I notice that qrule starts at 1, but without knowing which type and degree FE you are using, its hard to judge if the way you choose the quadrature order is correct). 

Regards,

J-P

[Deepak GUpta]

Dear Jean-Paul

Thanks for the prompt reply. Indeed it does not look like throwing any error in debug mode as well. I agree that the Kelly estimator works fine for hp-problems since the tutorial example on hp-adaptivity works fine for me as well. One thing that bugs  me how we can relate which element should use which face_quadrature rule. I mean we can specify which quadrature rule to use m=by passing the index of that location in the quadrature_collection when using hp_fe_values.reinit(). However, I did not see how to specify the face_quadrature for every cell. Thus, I assumed that face_quadrature_collection and quadrature_collection should be arranged in the same order. 

However, this should not be causing problems for now, since I am using linear elements right now  and Q<1>(1) should also be ok for the same. Anyways I checked by putting the face_quadrature rules starting from 5 onwards and higher. That does not help.

Best,

Deepak

--
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+un...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

[Wolfgang Bangerth]

On 08/26/2016 05:23 AM, Deepak GUpta wrote:
>
> I am trying to use KellyErrorEstimator of deal.II and I saw various tutorials
> which use it. I ran those examples and they work well. Recently, I tried it
> with my work and when I print the estimated_error_per_cell vector, all of the
> values in that vector are zero.

If you use it for mesh refinement, is the mesh refined? If it is, then your
routine for printing the elements has a bug. If the mesh is not refined, then
what is the solution? If the (finite element) solution is a linear function,
then the Kelly estimator will compute zeros.

I'm afraid there is nothing we can help with from a distance without seeing
the solution and/or the whole code...

Cheers
W.

--
------------------------------------------------------------------------
Wolfgang Bangerth email: bang...@colostate.edu
www: http://www.math.colostate.edu/~bangerth/

[Deepak Gupta]

Dear All,

In my problem case, I do not have any distributed load at the boundary which means no boundary Neumann boundary condition. However, I have a point load at one of the nodes (found the node and change the value) which I do not state in the KellyEstimator module and I am blind-guessing that this might be problem. Do I need to pass this information to the estimate() function, if yes, could you help me on how to do this?

Best

Deepak

--
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+unsubscribe@googlegroups.com.

[Wolfgang Bangerth]

On 08/30/2016 03:28 AM, Deepak Gupta wrote:
>
> In my problem case, I do not have any distributed load at the boundary
> which means no boundary Neumann boundary condition. However, I have a
> point load at one of the nodes (found the node and change the value)
> which I do not state in the KellyEstimator module and I am
> blind-guessing that this might be problem. Do I need to pass this
> information to the estimate() function, if yes, could you help me on how
> to do this?

No. The Kelly estimator does not need to know about the right hand side
of your problem.

Best