[Getfem-users] Unexpected mesh indexing and fem_dof

[Yousef Elkurdi]

Hi,

I am trying to build a simple mesh but I am seeing inconsistencies between the mesh point indexes of each convex and their FEM DOF numbering.  The mesh I am building is made of four triangles making up a square such as.  When I output the mesh node indexes for each convex and their corresponding FEM DOF, the last convex produces unexpected results.  The FEM DOF numbering does not correspond to the last convex added to the mesh, however, the mesh indexing is correct.

Cnvx: 0

Mesh indx [ 0 ]: 0 fem_dof: 0

Mesh indx [ 1 ]: 1 fem_dof: 1

Mesh indx [ 2 ]: 4 fem_dof: 2

Cnvx: 1

Mesh indx [ 0 ]: 1 fem_dof: 1

Mesh indx [ 1 ]: 4 fem_dof: 2

Mesh indx [ 2 ]: 2 fem_dof: 4

Cnvx: 2

Mesh indx [ 0 ]: 2 fem_dof: 4

Mesh indx [ 1 ]: 4 fem_dof: 2

Mesh indx [ 2 ]: 3 fem_dof: 3

Cnvx: 3

Mesh indx [ 0 ]: 3 fem_dof: 3

Mesh indx [ 1 ]: 4 fem_dof: 2  <--- ERROR: this should be (4)

Mesh indx [ 2 ]: 0 fem_dof: 0

The code that produces this output is below, I would appreciate it if you can tell if there is anything wrong I am doing,

Thank you,

Yousef

---------------- code ----------------------------------

    getfem::mesh mesh; /* the mesh */

    getfem::mesh_im mim(mesh); /* the integration methods. */

    getfem::mesh_fem mf_u(mesh); /* the main mesh_fem, for the Laplacian solution */

    mf_u.set_auto_add(getfem::fem_descriptor("FEM_PK(2,1)"));

    mim.set_auto_add(getfem::int_method_descriptor("IM_EXACT_SIMPLEX(2)"));

    // Adding points to mesh and obtaining their indexes.

    std::vector<bgeot::size_type> ind(5);

    ind[0] = mesh.add_point(bgeot::base_node(0.0, 0.0, 0.0));

    ind[1] = mesh.add_point(bgeot::base_node(1.0, 0.0, 0.0));

    ind[2] = mesh.add_point(bgeot::base_node(1.0, 1.0, 0.0));

    ind[3] = mesh.add_point(bgeot::base_node(0.0, 1.0, 0.0));

    ind[4] = mesh.add_point(bgeot::base_node(0.5, 0.5, 0.0));

    std::vector<bgeot::size_type> cp(3);

    // add convex 1

    cp[0] = ind[0];

    cp[1] = ind[1];

    cp[2] = ind[4];

    mesh.add_convex(bgeot::simplex_geotrans(2, 1), cp.begin());

    // add convex 2;

    cp[0] = ind[1];

    cp[1] = ind[4];

    cp[2] = ind[2];

    mesh.add_convex(bgeot::simplex_geotrans(2, 1), cp.begin());

    // add convex 3;

    cp[0] = ind[2];

    cp[1] = ind[4];

    cp[2] = ind[3];

    mesh.add_convex(bgeot::simplex_geotrans(2, 1), cp.begin());

    // add convex 4;

    cp[0] = ind[3];

    cp[1] = ind[0];

    cp[2] = ind[4];

    mesh.add_convex(bgeot::simplex_geotrans(2, 1), cp.begin());

    std::cout << "Hello..." << std::endl;

    bgeot::size_type ip_n;

    dal::bit_vector nn = mesh.convex_index();

    unsigned int msh_indx, dof_fem;

    for (ip_n << nn; ip_n != bgeot::size_type(-1); ip_n << nn) {

        std::cout << "Cnvx: " << ip_n << std::endl;

        for (int i = 0; i < mesh.nb_points_of_convex(ip_n); i++) { // global dof enumeration

            dof_fem = mf_u.ind_basic_dof_of_element(ip_n)[i];

            msh_indx = mesh.ind_points_of_convex(ip_n)[i];

            std::cout << "Mesh indx [ " << i << " ]: " << msh_indx << " fem_dof: " << dof_fem << std::endl;

        }

        std::cout << std::endl;

    }

------------ code end ----------------