interpolate_boundary_values does not work

[Zhuoran Wang]

Hi,

I have a question about the Dirichlet boundary conditions.

I'm using fe(FE_DGP<dim>(1),1, FE_FaceP<dim>(1),1) as the finite element space, FE_DGP is defined in interiors, FE_FaceP is defined on faces.

When I try to set the Dirichlet boundary conditions to the system by VectorTools::interpolate_boundary_values to constraints, I got the error message:

An error occurred in line <1080> of file </var/folders/8z/hlb6vc015qjggytkxn84m6_c0000gn/T/heltai/spack-stage/spack-stage-dealii-9.2.0-747mf54t7nef2iff4la5ab3qxlrpd2jb/spack-src/source/fe/fe.cc> in function

    virtual Point<dim - 1> dealii::FiniteElement<3, 3>::unit_face_support_point(const unsigned int) const

The violated condition was: 

    unit_face_support_points.size() == this->dofs_per_face

Additional information: 

    You are trying to access the support points of a finite element that either has no support points at all, or for which the corresponding tables have not been implemented.

I thought that FE_FaceP has support points on boundaries. Does it mean that I can't set the Dirichlet boundary conditions in this way? Or did I misunderstand anything?

Thank you.

Best,

Zhuoran

[Wolfgang Bangerth]

On 6/6/21 3:55 AM, Zhuoran Wang wrote:
>
> I thought that FE_FaceP has support points on boundaries. Does it mean that I
> can't set the Dirichlet boundary conditions in this way? Or did I
> misunderstand anything?

Hi Zhuoran,
nice to see that you're still using deal.II :-)

I believe that the FaceP spaces don't have support points in 3d whereas the
FaceQ spaces do. That's for the same reason as why a FE_DGP<2> can have no
support points: because they're not a complete tensor product polynomial.

Do you need to use FE_FaceP? Could you use FE_FaceQ?

Best
W.


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

[Zhuoran Wang]

Hi Dr.Bangerth,

Thank you for your quick reply! 

Yes, I'm still using deal.II to implement the WG.

I'm now trying to use WG(P1,P1) to solve Darcy. For WG(P0,P0), I used FE_DGQ and FE_FaceQ which are the same as FE_DGP and FE_FaceP. But for higher elements, I need to use FE_FaceP.

If the face dof is correctly figured out, I suppose I can assign the D.B.C to constraints by myself. Is that correct? 

Thank you.

Best,

Zhuoran

[Wolfgang Bangerth]

On 6/6/21 10:35 PM, Zhuoran Wang wrote:
> I'm now trying to use WG(P1,P1) to solve Darcy. For WG(P0,P0), I used FE_DGQ
> and FE_FaceQ which are the same as FE_DGP and FE_FaceP. But for higher
> elements, I need to use FE_FaceP.
> If the face dof is correctly figured out, I suppose I can assign the D.B.C to
> constraints by myself. Is that correct?

The issue is that FE_DGP is not interpolatory (i.e., the basis functions do
not have the phi_i(x_j)=delta_ij property), and consequently the
VectorTools::interpolate_boundary_values function does not work.

But you can do the equivalent work yourself:
for (cell=...)
for (face=...)
if (face->at_boundary())
{
get the DoF indices for that face, which in your case will only
correspond to the FE_FaceP component
determine which values these DoFs need to have
add the corresponding constraints to the AffineConstraints object
}

The second step in the inner loop might be the difficult one if you have
nontrivial boundary values since then you probably want to compute the
projection of these boundary values onto the finite element space on that
face. But if you have zero boundary values, for example, then you'd just set
all DoFs to zero in that loop.

[Zhuoran Wang]

Thanks for the clarification, Dr. Bangerth. 

I think I will use the zero boundary condition to test whether finite element spaces are appropriate.

Best,

Zhuoran

--
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 a topic in the Google Groups "deal.II User Group" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/dealii/TX_xy-KFDMQ/unsubscribe.
To unsubscribe from this group and all its topics, send an email to dealii+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/1ef2ff40-d865-dba5-847f-d5ad859ab04a%40colostate.edu.