How to use BDM element in a quadrilateral mesh?

1 month 3 days ago #3688 by Younghigh
Hi,

When running the following code:
from ngsolve import *
import ngsolve.meshes as ngm
mesh = ngm.MakeStructured2DMesh(nx=1,ny=1)
V = HDiv(mesh, order=2, RT=False)
gfu = GridFunction(V)
vtk = VTKOutput(ma=mesh, coefs=[gfu], names=["bdm-basis"], filename="bdm2basis", subdivision=3)
for i in range(gfu.vec.size):
    gfu.vec[:] = 0
    gfu.vec[i] = 1
    vtk.Do()
I find there are 24 vtk documents shown, which should be RT2 basis in a rectangle. But I have used
RT=False
, my question is: how to use BDM element in a quadrilateral mesh?


Best,
Di Yang

Please Log in or Create an account to join the conversation.

1 month 3 days ago #3689 by mneunteufel
Hi Di Yang,

to check the number of degrees of freedom on a single quad for different polynomial orders you can use the following code
from ngsolve import *
from ngsolve.meshes import MakeStructured2DMesh
mesh = MakeStructured2DMesh(quads=True,nx=1,ny=1)
Draw(mesh)

print("RT")
for p in range(6):
    fes = HDiv(mesh, order=p, RT=True)
    print(p,":",fes.ndof)

print("BDM")
for p in range(6):
    fes = HDiv(mesh, order=p+1)
    print(p,":",fes.ndof)
which gives the output
RT: 0 : 4, 1 : 12, 2 : 24, 3 : 40, 4 : 60, 5 : 84
BDM: 0 : 12, 1 : 24, 2 : 40, 3 : 60, 4 : 84, 5 : 112

As you can see they coincide (RT is "shifted"). The reason for this is that RT elements are available for triangles and tetrahedra, but not for quadrilaterals, hexahedra, or prisms. So there are only BDM elements for quads.

Best
Michael

Please Log in or Create an account to join the conversation.

1 month 3 days ago - 1 month 3 days ago #3690 by Younghigh
Thank u first for your specific replies. However, the dofs of BDMk on any rectangle K are devised as follows.
\[(\mathbf{w}\cdot\mathbf{n}, p_k)_{e_i},\quad \forall\,p_k\in P_k(e_i),\ i=1,2,3,4;\]
\[(\mathbf{w}, \mathbf{p}_{k-2})_K,\quad \forall\,\mathbf{p}_{k-2}\in (P_{k-2}(K))^2.\]

So, it should be
\[\mathrm{dim}(BDMk)=k^2+3k+4,\]

i.e., BDMk
k=1: 8;
k=2: 14;
k=3: 22;
k=4: 32;
...

But if running
V = HDiv(mesh, order=2, RT=False)
the number of dofs of V is 24. I think it is not a BDM2.

Best,
Di Yang

Please Log in or Create an account to join the conversation.

1 month 3 days ago #3691 by mneunteufel
Could be that it is exactly the other way round, sorry.

The Quad-HDiv basis in NGSolve is implemented as described in this thesis .
The basis span Q_{k+1,k} x Q_{k,k+1} if I remember correctly and thus has dimension 2*(k+2)*(k+1)=2k**2+6k+4.

Best
Michael

Please Log in or Create an account to join the conversation.

1 month 3 days ago #3692 by Younghigh
Ok, that means RTk elements indeed. Thank u very much.

Please Log in or Create an account to join the conversation.

© 2019 Netgen/NGSolve