Add `MeshFunction` constructor from `std::map`/`dict`

[Jan Blechta]

As MeshDomains interface changed recently, it would be useful for users to be able to construct MeshFunction from std::map / dict markers obtained for example by marker = mesh.domains().markers(i). It would enable to adapt old code with minimal changes.

[Prof Garth Wells]

Two possibilities are:

• A free function to build a MeshFunction from a Mesh and MeshData
• A MeshFunction constructor that takes a Mesh and a std::map (as in Issue title)

The former is better for hiding data, the latter is perhaps a bit cleaner. Both are memory robust and explicit (the old code was neither).

If we do not apply the SWIG typemap to the std::map, it should avoid the copying(?)

Preferences?

[Johan Hake]

I think it would be nice to keep the typemap, returning a usable dict instead of a std::map. If not only for debugging.

Maybe we could subclass MeshFunction (or FooFunction) to something like:

FooMarkers(mesh)

Where Foo is one of Cell, Facet, aso.

A MeshFunction of the proper dimension would then be created based on stored markers in the mesh. It could also raise an error* if no markers are registered for the particular dimension.

Johan

• Errors in C++ constructors is maybe a bad thing?

[Prof Garth Wells]

The first point (keeping the dict) points towards a free function to create the MeshFunction.

[Martin Sandve Alnæs]

So we'll be able to do

boundaries = MeshFunction("size_t", mesh, dim-1, mesh.domains())

?

It would be nice if FacetFunction and CellFunction are updated to mirror this. Then we can have just

mesh = Mesh(filename)
boundaries = FacetFunction(mesh, mesh.domains())

[Johan Hake]

Yes something like that. But passing MeshDomains as a second argument is superfluous as it is already contained in the mesh. That is why I suggested to subclass MeshFunction and/or FooFunction:

markers = FacetMarkers(mesh)

Maybe add a type argument.

[Prof Garth Wells]

I prefer Martin's suggestion. Once can dream of scenarios where one would want to use different MeshDomains:

boundaries0 = MeshFunction("size_t", mesh, dim-1, mesh_domains0)
boundaries1 = MeshFunction("size_t", mesh, dim-1, mesh_domains1)

More generally, perhaps we should allow a MeshFunction to be constructed from a std:

b = MeshFunction("size_t", mesh, dim-1, mesh.domains().markers(dim-1))

?

[Jan Blechta]

Second option is looking more transparently and seems more flexible. But I guess that first one is easier to handle in python.

[Prof Garth Wells]

I think (would need to test) that the second is easy in Python if we forgo the possibility of getting the domain marks as a dict.

[Johan Hake]

Agree. Also not sure proliferation in classes is worth it.

Not sure why we would need to wrap the second version to Python? To facilitate the above syntax in python without the dict typemap, we would need a C++ object to construct the std::map, for example a MeshDomain. But that is covered by the first syntax.

[Prof Garth Wells]

I don't quite follow the second paragraph. Which interface do you prefer?

Presently in C++ interface MeshDomains::markers(d) returns a std::map, and I guess in the absence of a typemap SWIG can transparently pass the std::map object into another function. Is this the case?

[Johan Hake]

Yes, but what would the point be?

In C++ one can alter the std::map, but in Python it would only be a stub, only functional to be passed to a constructor. What with including both constructors and ignore the last one in Python.

[Prof Garth Wells]

I don't follow your points. Could you illustrate with code examples?

[Johan Hake]

The only way we can do this effective is to remove the outtypemap of std::map so it returns a unusefull stubb.

markers = mesh.domains().markers(dim-1)

markers cannot be used to anything but be passed to a C++ function accepting a std::map, like the proposed MeshFunction constructor:

domains = MeshFunction("size_t", mesh, dim-1, markers)

Why not just pass the MeshDomains?

domains = MeshFunction("size_t", mesh, dim-1, mesh.domains())

In C++, however, it makes perfect sense as a std::map can be useful.

[Prof Garth Wells]

And your interface vote goes to . . . .

[Jan Blechta]

Does doing

domains = MeshFunction("size_t", mesh, dim-1, markers)

with markers being dict brings some overhead? Or is markers object only pointer to std::map and actual construction would be executed without converting std::map to dict?

[Johan Hake]

Well, I thought I was making an argument for including both constructors, but ignoring the std::map one in Python.

If that is a no-go I vote for the first one with MeshDomains

[Johan Hake]

@blechta: Converting a std::map to a dict brings piles of overhead, and yes to your second questions.

[Martin Sandve Alnæs]

The MeshDomains interface version will be more robust w.r.t. changes in the return type of markers(), and avoid adding the constraint that we cannot wrap std::map to dict. That could come back and bite us in an unrelated situation.

[Jan Blechta]

@JohanHake: So why typemap to dict needs to be canceled?

[Johan Hake]

If we land on using this interface in python:

b = MeshFunction("size_t", mesh, dim-1, mesh.domains().markers(dim-1))

We really need to remove the typemap, otherwise we might run into performance issues. By removing the dict typemap this syntax would work. Then we just pass the pointer. Which is something I would not prefer.

[Prof Garth Wells]

Let's do

domains = MeshFunction("size_t", mesh, dim-1, mesh.domains())

As Martin says, it hides the storage type (presently std::map) which will make any future changes less intrusive, and we can retain the dict typemap for debugging.

[Jan Blechta]

MeshFunction::_values type will be changed to pointer to mesh.domains().markers(i) or will data be copied?

[Johan Hake]

A MeshFunction stores data in contiguous arrays. A MeshDomain in std::map, so there will be coping.

[Prof Garth Wells]

Now in master

[Jan Blechta]

Thanks for very quick resolution, Garth, good job!