On Mon, 02 Sep 2013 11:52:55 +0200 Johan Hake <[email protected]> wrote:
> On Monday September 2 2013 11:45:21 Jan Blechta wrote: > > On Mon, 02 Sep 2013 11:36:28 +0200 > > > > Johan Hake <[email protected]> wrote: > > > On Monday September 2 2013 10:09:56 Garth N. Wells wrote: > > > > On 2 September 2013 09:58, Johan Hake <[email protected]> > > > > wrote: > > > > > On Monday September 2 2013 09:42:49 Garth N. Wells wrote: > > > > >> On 2 September 2013 09:30, Johan Hake <[email protected]> > > > > >> wrote: > > > > >> > Seems like updating the doc string wont help as enough > > > > >> > people have tried > > > > >> > to > > > > >> > use the vertex_to_dof_map and failed. > > > > >> > > > > > >> > I agree that the left to right reading does not apply to > > > > >> > the example Garth > > > > >> > presented. If that is the expected behavior, and I guess > > > > >> > it is given the > > > > >> > comments in this treahd, we should just rename the methods. > > > > >> > That would generalize the methods and probably fit better > > > > >> > into the general interface > > > > >> > of DofMap. > > > > >> > > > > > >> > However that would limit the scope of the map and remove > > > > >> > one important motivation for adding the map in the first > > > > >> > place, namely to turn general > > > > >> > > > > >> > vector function values ordered as: > > > > >> How can renaming limit scope? The functionality remains the > > > > >> same. > > > > > > > > > > Because it is not enough to just rename it. We also need to > > > > > remove the functionality for vector function spaces. Your > > > > > example does not make sense > > > > > if> > > > > > > > > > > you change: > > > > > V = FunctionSpace(mesh, "Lagrange", 1) > > > > > > > > > > to > > > > > > > > > > V = VectorFunctionSpace(mesh, "Lagrange", 1) > > > > > > > > My example was deliberately simple. The functionality can be > > > > retained. > > > > > > How? With a map of maps, where the second is the local dof index? > > > > > > > >> > vertex_index*dofs_per_vertex+local_dof > > > > >> > > > > > >> > to an array which could be feed directly into a vector of a > > > > >> > Function in a > > > > >> > VectorFunctionSpace (or similar mixed CG1 function spaces). > > > > >> > The present functionality also works for parallel runs, as > > > > >> > seen by the following example: > > > > >> > > > > > >> > mpirunt -np 2 python vertex_to_dofs.py > > > > >> > > > > > >> > # vertex_to_dofs.py > > > > >> > > > > > >> > from dolfin import * > > > > >> > import numpy as np > > > > >> > > > > > >> > mesh = UnitSquareMesh(20,20) > > > > >> > V = VectorFunctionSpace(mesh, "CG", 1) > > > > >> > u = Function(V) > > > > >> > vertex_to_dof_map = V.dofmap().vertex_to_dof_map(mesh) > > > > >> > > > > > >> > data = np.reshape(mesh.coordinates()[:], > > > > >> > (mesh.num_vertices()*2)) > > > > >> > > > > >> This is problematic - it makes an assumption of the ordering > > > > >> in mesh.coordinates(). > > > > > > > > > > The only assumption is that you have some data (possible > > > > > vector or tensor data) which are ordered based on the mesh > > > > > (vertices). > > > > > > > > > >> I have seen that a good re-ordering of mesh data > > > > >> can give up to a 50% speed up for assembly, and which will be > > > > >> added in the future. We should not be exposing low-level > > > > >> storage. > > > > > > > > > > Not sure what you mean. This has nothing to do with assemble. > > > > > Only transferring vertex based data into a Function. > > > > > > > > Exposing low level storage (e.g. (mesh.coordinates()), violates > > > > data hiding, which then can affect all parts of a code. If the > > > > mesh data ordering is changed, say to make assembler faster, > > > > your example code will likely break. > > > > > > Using mesh.coordinates() was just an example on some data which > > > comes together with the mesh. I just high-jacked coordinates to > > > represent some vector field aligned with the mesh. Real data > > > often comes aligned with the mesh and we need a simple and > > > intuitive way to read such data into a Function. This is basic > > > functionality alot of users need for their applications. > > > > > > If you intend to include some mesh reordering, I suggest that you > > > also include some mapping that brings mesh data to reordered mesh > > > data, and then we need a map to bring reordered mesh data to dof > > > ordering. > > > > Why? This is just changing vertex indices. Currently they're already > > ordered "somehow irregularly". So they will ordered in the other > > way... > > Because the average Joe have a mesh from somewhere together with some > data. If the mesh is reordered by some algorithm we need to provide a > way for Joe to map his data to the new mesh ordering, so he can call > > u.assign_vertex_data(reordered_vretex_based_data) > > or what not. If I get it correctly, mesh will be reordered once when loaded by DOLFIN. This is already being done and the method how to reorder some data accordingly is not available. For example, there is an open bug on launchpad, that when unordered DOLFIN XML mesh with mesh function is loaded and reordered, mesh function breaks. Loading of mesh functions is currently supported and safe to the extent supported in `meshconvert.py`. Note also that changing vertex indices will require to renumber cell-local vertices and facets according to UFC numbering scheme (=mesh.order()). Hence in addition you will have to provide a way how to keep up-to-date facet based data (FacetFunction), as this is not currently available as I said above. Jan > > Johan > > > Jan > > > > > > >> > u.vector().set_local(data[vertex_to_dof_map]) > > > > >> > plot(u, interactive=True) > > > > >> > > > > >> Why not just use Function::compute_vertex_values(...) (plus > > > > >> any necessary generalisation)? > > > > > > > > > > The comparison with compute_vertex_values is appropriate. It > > > > > was raised when we discussed the inclusion of the map in the > > > > > first place. However the (present) vertex_to_dof_map give the > > > > > mapping from vertex based data to a Function, where > > > > > compute_vertex_values does the opposite. > > > > > > > > Yes, but two functions were added to GenericDofMap. One seems to > > > > duplicate existing functionality. > > > > > > True. But the two maps, can only be used on data defined on > > > vertices (CG1). compute_vertex_values are more general as it > > > works for Functions on a lot more FunctionSpaces (CG2, DG0, aso)? > > > > > > > > The map is also just > > > > > computed once and can therefore be reused by the user if that > > > > > is needed. > > > > > > > > I don't see the benefit if one can use > > > > Function::compute_vertex_values. > > > > > > See Martin's answer. > > > > > > Johan > > > > > > > Garth > > > > > > > > > Johan > > > > > > > > > >> Garth > > > > >> > > > > >> > Johan > > > > >> > > > > > >> > On Saturday August 31 2013 10:20:21 Simone Pezzuto wrote: > > > > >> >> Hi, > > > > >> >> > > > > >> >> I'm familiar with these two maps since I use > > > > >> >> them > > > > >> >> > > > > >> >> for a gradient > > > > >> >> > > > > >> >> recovery technique. > > > > >> >> > > > > >> >> I can assure you that first time I used vertex_to_dof_map > > > > >> >> I was a bit confused, > > > > >> >> since the convention should be left to right (as Garth > > > > >> >> pointed out). > > > > >> >> > > > > >> >> Example: eps2pdf fig.eps ---> fig.pdf > > > > >> >> > > > > >> >> vertex2dof vertex_id --> dof_id > > > > >> >> dof2vertex dof_id --> vertex_id > > > > >> >> > > > > >> >> So at the moment is really confusing. Maybe we can > > > > >> >> introduce new functions > > > > >> >> {vertex2dof,dof2vertex}_map > > > > >> >> (no name collision) and deprecate the old one, so the > > > > >> >> user is aware of the > > > > >> >> change but its code doesn't brake. > > > > >> >> > > > > >> >> Simone > > > > >> >> > > > > >> >> 2013/8/31 Jan Blechta <[email protected]> > > > > >> >> > > > > >> >> > On Fri, 30 Aug 2013 23:47:35 +0100 > > > > >> >> > > > > > >> >> > "Garth N. Wells" <[email protected]> wrote: > > > > >> >> > > On 30 August 2013 23:37, Johan Hake > > > > >> >> > > > > > > >> >> > > <[email protected]> wrote: > > > > >> >> > > > On Friday August 30 2013 23:19:09 Garth N. Wells > > > > >> >> > > > wrote: > > > > >> >> > > >> On 30 August 2013 22:50, Johan Hake > > > > >> >> > > >> <[email protected]> > > > > > > wrote: > > > > >> >> > > >> > On Friday August 30 2013 15:47:28 Garth N. Wells > > > > >> >> > > >> > > > > > >> >> > > >> > wrote: > > > > >> >> > > >> >> The functions GenericDofmap::vertex_to_dof_map > > > > >> >> > > >> >> and GenericDofMap::dof_to_vertex_map are not > > > > >> >> > > >> >> properly documented (the doc string is the same > > > > >> >> > > >> >> for both), and I think that they are back to > > > > >> >> > > >> >> front. The docstring in DofMap has > > > > >> >> > > >> >> inconsistencies. I would expect that > > > > >> >> > > >> >> > > > > >> >> > > >> >> map0 = > > > > >> >> > > >> >> GenericDofmap::vertex_to_dof_map(...) > > > > >> >> > > >> >> > > > > >> >> > > >> >> would mean a map from vertex to dof, i.e. > > > > >> >> > > >> >> > > > > >> >> > > >> >> map0[vertex_index] -> dof index > > > > >> >> > > >> >> > > > > >> >> > > >> >> and that > > > > >> >> > > >> >> > > > > >> >> > > >> >> map1 = > > > > >> >> > > >> >> GenericDofmap::dof_to_vertex_map(...) > > > > >> >> > > >> >> > > > > >> >> > > >> >> would mean a map from dof index to > > > > >> >> > > >> >> > > > > >> >> > > >> >> map1[dof_index] -> vertex index > > > > >> >> > > >> >> > > > > >> >> > > >> >> Tests (see below code) and the return types also > > > > >> >> > > >> >> indicate that > > > > >> >> > > >> >> things are back to front. Can someone clarify > > > > >> >> > > >> >> the situation? > > > > >> >> > > >> > > > > > >> >> > > >> > The map was introduced to help a user map vertex > > > > >> >> > > >> > based data onto > > > > >> >> > > >> > a Function.> > > > > >> >> > > >> > > > > > >> >> > > >> > from dolfin import * > > > > >> >> > > >> > import numpy as np > > > > >> >> > > >> > > > > > >> >> > > >> > mesh = UnitSquareMesh(20,20) > > > > >> >> > > >> > V = VectorFunctionSpace(mesh, "CG", 1) > > > > >> >> > > >> > u = Function(V) > > > > >> >> > > >> > vertex_to_dof_map = > > > > >> >> > > >> > > > > > >> >> > > >> > V.dofmap().vertex_to_dof_map(mesh) > > > > >> >> > > >> > > > > > >> >> > > >> > data = np.reshape(mesh.coordinates()[:], > > > > >> >> > > >> > > > > > >> >> > > >> > (mesh.num_vertices()*2)) u.vector()[:] = > > > > >> >> > > >> > data[vertex_to_dof_map] > > > > >> >> > > >> > > > > > >> >> > > >> > plot(u, interactive=True) > > > > >> >> > > >> > > > > > >> >> > > >> > The size of the data array should be: > > > > >> >> > > >> > mesh.num_vertices()*u.value_size() > > > > >> >> > > >> > > > > > >> >> > > >> > The documentation should be improved, and not > > > > >> >> > > >> > least properly mapped from C++ to Python. > > > > >> >> > > >> > > > > > >> >> > > >> > The name refer to the mapping that turn vertex > > > > >> >> > > >> > based data to dof > > > > >> >> > > >> > based and reads quite well when used as above. I > > > > >> >> > > >> > can see that the word map can be missleading. It > > > > >> >> > > >> > is not a "map" data structure. It is an index > > > > >> >> > > >> > set that "maps values". > > > > >> >> > > >> > > > > > >> >> > > >> > Still confused? > > > > >> >> > > >> > > > > >> >> > > >> I'm not confused. It's clear that the function > > > > >> >> > > >> names are back-to-front. It doesn't matter what > > > > >> >> > > >> they were included for - they > > > > >> >> > > >> are members of GenericDofMap and must make sense in > > > > >> >> > > >> that context. > > > > >> >> > > >> > > > > >> >> > > >> Since reading from left to right is a well > > > > >> >> > > >> established convention, > > > > >> >> > > >> I propose that (a) the function names be fixed by > > > > >> >> > > >> reversing them; > > > > >> >> > > >> and (b) the doc strings be fixed. > > > > >> >> > > > > > > > >> >> > > > Agree on (b). I am not fully convinced by (a). > > > > >> >> > > > > > > > >> >> > > > I am not sure what your example tries to show. You > > > > >> >> > > > are not using the mapping the intended way and I am > > > > >> >> > > > therefore confused about the > > > > >> >> > > > whole back-to-front, front-to-back discussion. > > > > >> >> > > > > > > >> >> > > Just read the function names aloud from left to right > > > > >> >> > > - 'vertex_to_dof_map' should be a 'vertex to dof > > > > >> >> > > map', i.e. a map from > > > > >> >> > > a > > > > >> >> > > vertex *to* a dof. > > > > >> >> > > > > > >> >> > Just read from left to right - 'vertex_to_dof_map' > > > > >> >> > stands for a map which turns a vertex map into a dof > > > > >> >> > map (when used as a right composition). > > > > >> >> > > > > > >> >> > Yes, I was confused at first when I saw this and agree > > > > >> >> > with Garth it should be 'left to right'. But does it > > > > >> >> > worth switching it? Is the whole > > > > >> >> > concept of indexing by > > > > >> >> > > > > > >> >> > vertex_index*dofs_per_vertex+local_dof > > > > >> >> > > > > > >> >> > sustainable? Or should it be replaced by some more > > > > >> >> > robust types which > > > > >> >> > would handle non-injective map (and its inversion)? > > > > >> >> > > > > > >> >> > There were some user codes using these functions as > > > > >> >> > seen in discussions. > > > > >> >> > > > > > >> >> > Jan > > > > >> >> > > > > > >> >> > > Garth > > > > >> >> > > > > > > >> >> > > > Johan > > > > >> >> > > > > > > > >> >> > > >> Garth > > > > >> >> > > >> > > > > >> >> > > >> > Johan > > > > >> >> > > >> > > > > > >> >> > > >> >> Garth > > > > >> >> > > >> >> > > > > >> >> > > >> >> > > > > >> >> > > >> >> > > > > >> >> > > >> >> from dolfin import * > > > > >> >> > > >> >> > > > > >> >> > > >> >> mesh = UnitSquareMesh(4, 4) > > > > >> >> > > >> >> V = FunctionSpace(mesh, "Lagrange", 1) > > > > >> >> > > >> >> > > > > >> >> > > >> >> dof_to_vertex = > > > > >> >> > > >> >> V.dofmap().dof_to_vertex_map(mesh) > > > > >> >> > > >> >> vertex_to_dof = > > > > >> >> > > >> >> V.dofmap().vertex_to_dof_map(mesh) > > > > >> >> > > >> >> > > > > >> >> > > >> >> for c in cells(mesh): > > > > >> >> > > >> >> print "Cell index:", c.index() > > > > >> >> > > >> >> > > > > >> >> > > >> >> # Get cell dofs > > > > >> >> > > >> >> dofs = V.dofmap().cell_dofs(c.index()) > > > > >> >> > > >> >> print " Cell dofs:", dofs > > > > >> >> > > >> >> > > > > >> >> > > >> >> # Get vertices from cell > > > > >> >> > > >> >> cell_vertices0 = sorted([v.index() for v in > > > > >> >> > > >> >> vertices(c)]) > > > > >> >> > > >> >> print " Cell vertex indices (from cell):", > > > > >> >> > > >> >> cell_vertices0 > > > > >> >> > > >> >> > > > > >> >> > > >> >> # Get vertices from dof_to_vertex > > > > >> >> > > >> >> cell_vertices1 = sorted([dof_to_vertex[dof] > > > > >> >> > > >> >> for > > > > >> >> > > >> >> > > > > >> >> > > >> >> dof in > > > > >> >> > > >> >> > > > > >> >> > > >> >> dofs]) print " Cell vertex indices (from > > > > >> >> > > >> >> dof_to_vertex_map):", > > > > >> >> > > >> >> > > > > >> >> > > >> >> cell_vertices1 > > > > >> >> > > >> >> > > > > >> >> > > >> >> # Get vertices from vertex_to_dof_map > > > > >> >> > > >> >> cell_vertices2 = sorted([vertex_to_dof[dof] > > > > >> >> > > >> >> for > > > > >> >> > > >> >> > > > > >> >> > > >> >> dof in > > > > >> >> > > >> >> > > > > >> >> > > >> >> dofs]) print " Cell vertex indices (from > > > > >> >> > > >> >> vertex_to_dof_map):", > > > > >> >> > > >> >> > > > > >> >> > > >> >> cell_vertices2 > > > > >> >> > > >> >> > > > > >> >> > > >> >> _______________________________________________ > > > > >> >> > > >> >> fenics mailing list > > > > >> >> > > >> >> [email protected] > > > > >> >> > > >> >> http://fenicsproject.org/mailman/listinfo/fenics > > > > >> >> > > > > > > >> >> > > _______________________________________________ > > > > >> >> > > fenics mailing list > > > > >> >> > > [email protected] > > > > >> >> > > http://fenicsproject.org/mailman/listinfo/fenics > > > > >> >> > > > > > >> >> > _______________________________________________ > > > > >> >> > fenics mailing list > > > > >> >> > [email protected] > > > > >> >> > http://fenicsproject.org/mailman/listinfo/fenics > > > > >> > > > > > >> > _______________________________________________ > > > > >> > fenics mailing list > > > > >> > [email protected] > > > > >> > http://fenicsproject.org/mailman/listinfo/fenics > > > > >> > > > > >> _______________________________________________ > > > > >> fenics mailing list > > > > >> [email protected] > > > > >> http://fenicsproject.org/mailman/listinfo/fenics > > > > > > > > _______________________________________________ > > > > fenics mailing list > > > > [email protected] > > > > http://fenicsproject.org/mailman/listinfo/fenics > > > > > > _______________________________________________ > > > fenics mailing list > > > [email protected] > > > http://fenicsproject.org/mailman/listinfo/fenics > > > > _______________________________________________ > > fenics mailing list > > [email protected] > > http://fenicsproject.org/mailman/listinfo/fenics > _______________________________________________ > fenics mailing list > [email protected] > http://fenicsproject.org/mailman/listinfo/fenics _______________________________________________ fenics mailing list [email protected] http://fenicsproject.org/mailman/listinfo/fenics
