dolfin team mailing list archive

[Johan Hake <hake.dev@xxxxxxxxx>]

On 01/28/2013 11:36 PM, Anders Logg wrote:
> On Mon, Jan 28, 2013 at 10:33:15PM +0000, Garth N. Wells wrote:
>> On 28 January 2013 22:27, Anders Logg <logg@xxxxxxxxx> wrote:
>>> On Mon, Jan 28, 2013 at 11:19:13PM +0100, Johan Hake wrote:
>>>> On 01/28/2013 11:08 PM, Anders Logg wrote:
>>>>> On Mon, Jan 28, 2013 at 10:37:46PM +0100, Johan Hake wrote:
>>>>>> On 01/28/2013 10:22 PM, Anders Logg wrote:
>>>>>>> On Mon, Jan 28, 2013 at 10:10:40PM +0100, Johan Hake wrote:
>>>>>>>> Hello!
>>>>>>>>
>>>>>>>> I have now added a method to DofMap, which tabulate a map between
>>>>>>>> vertices and dofs. It will only work for DofMaps with dofs on vertices.
>>>>>>>> So basically for any CG 1 FunctionSpaces and will raise error else wise.
>>>>>>>>
>>>>>>>> Usage:
>>>>>>>>
>>>>>>>>   from dolfin import *
>>>>>>>>   import numpy as np
>>>>>>>>   mesh = UnitSquareMesh(10,10)
>>>>>>>>   V = VectorFunctionSpace(mesh, "CG", 1)
>>>>>>>>   u = Function(V)
>>>>>>>>   u.interpolate(Constant((1,2)))
>>>>>>>>   vert_values = np.zeros(mesh.num_vertices()*2)
>>>>>>>>   vert_values[V.dofmap().tabulate_vertex_map(mesh)] = \
>>>>>>>>   u.vector().array()
>>>>>>>>   print vert_values
>>>>>>>>
>>>>>>>> In parallel the map will follow the local dofs. This means that some
>>>>>>>> values in vert_values above will still be 0. The 0 values will then
>>>>>>>> correspond to vertex values which are owned by another process.
>>>>>>>>
>>>>>>>> In C++ the method returns a std::vector<std::size_t>.
>>>>>>>>
>>>>>>>> Questions:
>>>>>>>> Does the name "tabulate_vertex_map" work?
>>>>>>>
>>>>>>> Sounds ok.
>>>>>>>
>>>>>>>> Should a version of this method exist in other classes. FunctionSpace,
>>>>>>>> Function (without the mesh argument of course)?
>>>>>>>
>>>>>>> Yes, would be good to have for convenience.
>>>>>>>
>>>>>>> What is the relation to Function::compute_vertex_values?
>>>>>>
>>>>>> compute_vertex_values works on any GenericFunction. However, for all
>>>>>> Functions from CG 1 families I guess one could use this map to just
>>>>>> tabulate the entries.
>>>>>
>>>>> What does the tabulate_vertex_map method do which
>>>>> compute_vertex_values doesn't? Enable writing to dofs?
>>>>
>>>> RTFC! ;)
>>>>
>>>> tabulate_vertex_map only computes a map based on the prior knowledge of
>>>> the dofs coinciding with the vertices. compute_vertex_values are much
>>>> more general and used
>>>>
>>>>   element.interpolate_vertex_values
>>>>
>>>> to do the heavy lifting.
>>>
>>> Yes, so why is tabulate_vertex_map needed, if compute_vertex_values
>>> already computes the vertex values and is much more general?
>>>
>>>> I see tabulate_vertex_map as a convenience mapping for users who want
>>>> their data on a vertex based format, or users who want data which
>>>> resides on vertices into dolfin la structures.
>>>
>>> What is the use case? If it is for plotting or postprocessing, then
>>> compute_vertex_values should be enough.
>>>
>>> Or is the intention to use it for modifying values of CG1 functions
>>> from the outside?
>>>
>>
>> I don't think that it's required, but I guess Johan is adding it by
>> popular demand.

;)

Well I have a similar mapping in my own tools, and instead of giving
that code away all the freaking time I pushed it to DOLFIN. My code was
also Python only and recently we got some C++ request.

> I'm not objecting to it. I just came to think about the
> compute_vertex_values function and it should be enough to cover most
> of the use cases that are the background for the popular demand: users
> of Hans Petter's tutorial who want to get the vertex values on a
> regular grid for postprocessing in external software like MATLAB.

I think the compute_vertex_values could fill the purpose of most use
cases, but I think having the mapping available could be nice and
convenient.

Even if we discourage indexing into vectors, people tend to ignore it,
because it is often just the natural way to do stuff. Especially for
Functions in CG1 family spaces. By forcing a user to use local process
data for all things numpy, we can be useful for both numpy and HPC folks.

Johan

> --
> Anders
>