Issue with Facet/Cell tag

mesh
1

Hello everyone,
Greetings of the day.

I am using Dolfinx 0.7.1 and was trying to get calculated values on nodes represented by the facet tags.

I read the Gmsh mesh (which is a sqr as attached in .geo file) as below-

domain, cell_markers, facet_markers = gmshio.read_from_msh("try_sqr.msh", mesh_comm, gmsh_model_rank =2 , gdim=1)

Now as per my understanding, facet_markers should consider the tags and the corresponding nodal indices.

When I print the indices

facet_markers.indices

I get

array([ 0,  2,  5, 10, 18, 22, 28, 36, 44, 46, 62, 63, 64, 65, 69, 70],
      dtype=int32)

In the mesh, I have only 30 nodes. So I do not understand why it gives indices which are treater than number of nodes.

The way I plan to get the nodal values of the boundaries is-

u_surface = uh.x.array[ct.find(3)]

where uh is my calculated solution and boundary with tag 3 is where I am interested. Is this the correct way? Or there is a better way

When I do this at present, since node indices are above the actual, I get the out-of-bound error (which is obvious) since the maximum index for the solution is the number of nodes.

Looking forward to your kind support.

The .geo file-

lc = 0.25 ;

Point(1) = {0 , 0, 0, lc};
Point(2) = {1 , 0, 0, lc};
Point(3) = {1 , 1, 0, lc};
Point(4) = {0 , 1, 0, lc};


Line(1) = {1, 2};
Line(2) = {2, 3};
Line(3) = {3, 4};
Line(4) = {4, 1};

Curve Loop(1) = {1, 2, 3, 4};

Plane Surface(1) = {1};

Physical Curve("bottom") = {1}; 
Physical Curve("right") = {2}; 
Physical Curve("top") = {3}; 
Physical Curve("left") = {4}; 

Physical Surface("Copper") = {1};
2

You are actually meshing a 2D mesh here, so why are you sending in gdim=1?

3

Extremely sorry, that happened to be a typo when I was coping the stuff. Its 2 , I checked my code.

4

This is the actual

domain, cell_markers, facet_markers = gmshio.read_from_msh("try_sqr.msh", mesh_comm, gmsh_model_rank, gdim=2)
5

If your mesh is 2D, then the facets are line segments, not mesh nodes.

6

Ok, I understand now. That makes it clear why my boundary condition is working fine, but I am not able to get the nodal values using the same approach.
So how can I get the nodal values and coordinates for a particular facet tag? I tried using -

facet_geo = dolfinx.mesh.locate_entities(domain, domain.topology.dim -1, marker)

to get the nodal data, but I am not able to specify the marker properly. Should the marker be the tag number, like marker=3, or something like that?

7

Use entities_to_geometry, i.e.

    indices = dolfinx.cpp.mesh.entities_to_geometry(domain._cpp_object, mesh.topology.dim-1,
                   facet_markers.find(marker), False)
    cell_vertices = domain.geometry.x[indices]
    print(cell_vertices)

ref: https://github.com/search?q=repo%3AFEniCS%2Fdolfinx%20entities_to_geometry&type=code
and various other discourse posts (search for entities_to_geometry)

8

Thank you very much for your guidance.