Hello,

I've been having trouble evaluating a cell-centered field (stored in a
ExplicitSystem) at a point using *Number libMesh::System::point_value (
unsigned int var, const Point & p, const Elem * e). *More specifically, I
want to evaluate the field at the centroid of the elements. For most
points/elements, *libMesh::System::point_value* returns zero, but for some
elements it returns the expected non-zero value. Please see the toy example
below.

Does anyone have suggestions on how to deal with this?

Libmesh version: 1.5.1

Thanks,
Laryssa

 This is a toy example:
#include "libmesh/exodusII_io.h"
#include "libmesh/mesh.h"
#include "libmesh/elem.h"
#include "libmesh/mesh_base.h"
#include "libmesh/explicit_system.h"
#include "libmesh/equation_systems.h"
#include "libmesh/dof_map.h"
int main(int argc, char ** argv)
{
    using namespace libMesh;
    LibMeshInit init(argc, argv, MPI_COMM_WORLD);
    Mesh mesh_fiber(init.comm());
    ExodusII_IO importer(mesh_fiber);
    importer.read("file.e");
    mesh_fiber.prepare_for_use();

    libMesh::EquationSystems  es_fiber(mesh_fiber);
    typedef libMesh::ExplicitSystem  FiberSystem;
    FiberSystem& fiber_system = es_fiber.add_system<FiberSystem>("fibers");
    fiber_system.add_variable("fibersX", libMesh::CONSTANT, libMesh
::MONOMIAL);
    fiber_system.add_variable("fibersY", libMesh::CONSTANT, libMesh
::MONOMIAL);
    fiber_system.add_variable("fibersZ", libMesh::CONSTANT, libMesh
::MONOMIAL);
    fiber_system.init();

    unsigned int step=1;

    //Read fiber data into mesh_fiber
    importer.copy_elemental_solution(fiber_system, "fibersX", "fibersX",
step);
    importer.copy_elemental_solution(fiber_system, "fibersY", "fibersY",
step);
    importer.copy_elemental_solution(fiber_system, "fibersZ", "fibersZ",
step);
    fiber_system.print_info();

    libMesh::MeshBase::const_element_iterator el = mesh_fiber.
active_local_elements_begin();
    const libMesh::MeshBase::const_element_iterator end_el = mesh_fiber.
active_local_elements_end();

    std::vector<libMesh::Number> ff(3);
    for (; el != end_el; ++el)
    {
        Elem * elem= *el;
        libMesh::Point centroid = elem->centroid();

        ff[0] = fiber_system.point_value(fiber_system.variable_number(
"fibersX"), centroid, *elem);
        ff[1] = fiber_system.point_value(fiber_system.variable_number(
"fibersY"), centroid, *elem);
        ff[2] = fiber_system.point_value(fiber_system.variable_number(
"fibersZ"), centroid, *elem);
        std::cout << "ff = (" << ff[0] << ", " << ff[1] << ", " << ff[2] <<")
, norm=" << std::sqrt(ff[0]*ff[0]+ff[1]*ff[1]+ff[2]*ff[2])  << std::endl;

    }

    return 0;
}

_______________________________________________
Libmesh-users mailing list
Libmesh-users@lists.sourceforge.net
https://lists.sourceforge.net/lists/listinfo/libmesh-users

Reply via email to