[Common] Use bindable grid function adapter

* The reason for this change is that although dune-grid has its own conversion to a bindable functions it chooses a wrong VTKInfo for tensors. With this change, we select the correct VTKInfo.
parent 3f8b940c
......@@ -16,6 +16,43 @@
namespace Dune{
namespace Dorie{
template<class GF>
class VTKGridFunctionAdapter
{
using GridFunction = typename std::decay<GF>::type;
using GridView = typename GF::Traits::GridViewType;
using Entity = typename GridView::template Codim<0>::Entity;
using Range = typename GF::Traits::RangeType;
using Domain = typename GF::Traits::DomainType;
public:
VTKGridFunctionAdapter(std::shared_ptr<const GF> gf)
: _gf(gf)
{}
void bind(const Entity& e)
{
_entity = &e;
}
void unbind()
{
_entity = nullptr;
}
Range operator()(const Domain& x) const
{
Range y;
_gf->evaluate(*_entity,x,y);
return y;
}
private:
const Entity* _entity;
std::shared_ptr<const GF> _gf;
};
/*-------------------------------------------------------------------------*//**
* @brief Class for grid fucntion vtk sequence writer. This class works
* exacly the same way as a VTKSequence writer but it receives
......@@ -44,8 +81,16 @@ public:
static_assert(std::is_same<typename GF::Traits::GridViewType,GV>::value,
"GridFunctionVTKSequenceWriter only works with only one GridView type.");
auto vtk_gf_ptr = std::make_shared<Dune::PDELab::VTKGridFunctionAdapter<GF> >(gf_ptr, name);
Dune::VTKSequenceWriter<GV>::addVertexData(vtk_gf_ptr);
Dune::Dorie::VTKGridFunctionAdapter<GF> vtk_gf(gf_ptr);
const int vtk_ncomps = GF::Traits::dimRange;
const int dim = GF::Traits::dimDomain;
auto vtk_type = (vtk_ncomps == dim) ? VTK::FieldInfo::Type::vector
: VTK::FieldInfo::Type::scalar;
auto vtk_pc = Dune::VTK::Precision::float32;
auto vtk_info = VTK::FieldInfo(name,vtk_type,vtk_ncomps,vtk_pc);
this->vtkWriter()->addVertexData(vtk_gf,vtk_info);
}
/*-----------------------------------------------------------------------*//**
......@@ -62,8 +107,16 @@ public:
static_assert(std::is_same<typename GF::Traits::GridViewType,GV>::value,
"GridFunctionVTKSequenceWriter only works with only one GridView type.");
auto vtk_gf_ptr = std::make_shared<Dune::PDELab::VTKGridFunctionAdapter<GF> >(gf_ptr, name);
Dune::VTKSequenceWriter<GV>::addCellData(vtk_gf_ptr);
Dune::Dorie::VTKGridFunctionAdapter<GF> vtk_gf(gf_ptr);
const int vtk_ncomps = GF::Traits::dimRange;
const int dim = GF::Traits::dimDomain;
auto vtk_type = (vtk_ncomps == dim) ? VTK::FieldInfo::Type::vector
: VTK::FieldInfo::Type::scalar;
auto vtk_pc = Dune::VTK::Precision::float32;
auto vtk_info = VTK::FieldInfo(name,vtk_type,vtk_ncomps,vtk_pc);
this->vtkWriter()->addCellData(vtk_gf,vtk_info);
}
/*-----------------------------------------------------------------------*//**
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment