Attaching mesh value collections to mesh domains fails

[Cian Wilson]

I'm trying to read in an xml file containing facet ids and attach it to the mesh domains class beneath a mesh. However this currently fails. The example main.cpp (compile with a demo CMakeLists.txt) below illustrates this problem (in a slightly contrived way - without requiring an extra xml file).

Several suggestions have been made to fix this in question 174556 ( https://answers.launchpad.net/dolfin/+question/174566 ). I've implemented the last one, an operator to copy the data from a MeshFunction to a MeshValueCollection:
const MeshValueCollection<T>& operator=(const MeshFunction<T>& mesh_function);
The implementation is essentially the same as the existing constructor. Am happy to propose a merge of this if this is the preferred solution.

Many thanks,
Cian

main.cpp:
-------------------------------------------------------------
#include <dolfin.h>

class Side : public dolfin::SubDomain
{
    public:
    Side(const double &side) : side_(side)
    {
    }

    bool inside(const dolfin::Array<double>& x, bool on_boundary) const
    {
      return (std::fabs(x[0] - side_) < DOLFIN_EPS && on_boundary);
    }

  private:

    double side_;
};

int main()
{
  dolfin::UnitInterval mesh(10);

  Side left(0.0);
  Side right(1.0);

  dolfin::MeshFunction<dolfin::uint> edgeids(mesh, 0, 0);
  left.mark(edgeids, 1);
  right.mark(edgeids, 2);
  std::cout << "edgeids:" << std::endl;
  std::cout << edgeids.str(true) << std::endl;

  boost::shared_ptr< dolfin::MeshValueCollection<dolfin::uint> >
   edgecoll( new dolfin::MeshValueCollection<dolfin::uint>(edgeids) );
  std::cout << "edgecoll:" << std::endl;
  std::cout << edgecoll->str(false) << std::endl;

  std::cout << "makers(0) before assignment:" << std::endl;
  std::cout << mesh.domains().markers(0).str(false) << std::endl;

  mesh.domains().markers_shared_ptr(0) = edgecoll;

  std::cout << "makers(0) after assignment:" << std::endl;
  std::cout << mesh.domains().markers(0).str(false) << std::endl;
  if(!mesh.domains().facet_domains(mesh))
    dolfin::error("No facet_domains found.");
  else
    std::cout << "facet_domains found." << std::endl;

  return 0;
}
------------------------------------------------------------------

[Johan Hake (johan-hake)]

Also see:
https://blueprints.launchpad.net/dolfin/+spec/extend-meshvaluecollection