Incorrect visualization of submesh transferred to mesh

[vladimir-yu-fedorov]

Hello,

In the MFEM example below I have a simple 3D box geometry (geometry.zip) where I extract one of the faces as a 2D submesh. On this submesh, I create a vector grid function u_submesh and project some constant vector onto it. I then transfer u_submesh to the vector grid function u defined on the original 3D mesh.

When I use glvis to visualizeu_submesh on a 2D submesh, I see the correct values given by the projected vector.
However, when ask glvis to visualize u on the original 3D mesh, I see some unexpected values (see the figure at the end).

Meanwhile, by manually checking the u values on the corresponding face, I can verify that the transfer from u_submesh was correct.

Please help me understand why glvis shows wrong values ​​in 3D case.

#include "mfem.hpp"

using namespace mfem;
using namespace std;


int main(int argc, char *argv[]) {

    Mpi::Init(argc, argv);
    Hypre::Init();

    int order = 1;   // Finite element polynomial degree


    // 3D mesh -----------------------------------------------------------------
    Mesh mesh_serial("netgen/geometry.nm");
    ParMesh mesh(MPI_COMM_WORLD, mesh_serial);
    mesh_serial.Clear();   // the serial mesh is no longer needed
    mesh.Save("mesh");   // save mesh to the external file

    ND_FECollection fec(order, mesh.Dimension());
    ParFiniteElementSpace fespace(&mesh, &fec);
    ParGridFunction u(&fespace);
    u = 0;


    // 2D submesh --------------------------------------------------------------
    int iface = 4;

    Array<int> marker(1);
    marker[0] = iface + 1;

    ParSubMesh submesh(ParSubMesh::CreateFromBoundary(mesh, marker));
    submesh.Save("submesh");

    ND_FECollection fec_submesh(order, submesh.Dimension());
    ParFiniteElementSpace fespace_submesh(&submesh, &fec_submesh);
    ParGridFunction u_submesh(&fespace_submesh);
    u_submesh = 0;


    // Initialize the submesh grid function ------------------------------------
    Vector v(3);
    v[0] = 1.0; v[1] = 2.0; v[2] = 3.0;
    VectorConstantCoefficient vcoef(v);

    u_submesh.ProjectCoefficient(vcoef);
    u_submesh.Save("u_submesh");


    // Transfer from submesh to mesh -------------------------------------------
    submesh.Transfer(u_submesh, u);
    u.Save("u");


    // Inspect the boundary values of the grid function ------------------------
    VectorGridFunctionCoefficient u_coef(&u);
    Geometry geometry;

    for(int bfi = 0; bfi < mesh.GetNBE(); ++bfi) {
        const int bface_attr = mesh.GetBdrAttribute(bfi);

        if(bface_attr != iface+1) {
        continue;
        }

        const FiniteElement *fe = fespace.GetBE(bfi);
        ElementTransformation *tr = mesh.GetBdrElementTransformation(bfi);

        // Apply the transformation:
        const auto ip = geometry.GetCenter(fe->GetGeomType());
        tr->SetIntPoint(&ip);

        Vector vec(3);
        u_coef.Eval(vec, *tr, ip);

        cout << "Boundary Vec: [" << vec(0) << ", "
                                  << vec(1) << ", "
                                  << vec(2) << "]" << endl;
    }


return 0;
}

The x component of u is expected to be 1 on the front face, but what I see is