dolfin team mailing list archive

[Johan Hake <hake@xxxxxxxxx>]

On Thursday 19 February 2009 13:57:09 A Navaei wrote:
> 2009/2/19 Johan Hake <hake@xxxxxxxxx>:
> > The previous email was sent a bit premature...
> >
> > [snip]
> >
> >> > I have also thought of this :)
> >> >
> >> > Then you can take an ordered mesh, iterate over the vertices, and use
> >> > the present implementation of the eval function (which in its present
> >> > shape is a piecewise constant interpolator isn't it?), to fill the new
> >> > _vector, all done in the constructor. Then you have a discrete
> >> > function in the first place :)
> >>
> >> Is it possible to set the pointer to _vector data directly where we
> >> can pass the pointer to image data? If this works, no loop will be
> >> involve in data conversion and it will be super fast.
> >
> > Yes but you do not know how the mesh is ordered.
>
> Consider a 2D image data mapped on a structured grid, with size
> (number of nodes in each direction) of (sx, sy), the value of image
> brightness at coordinate (i, j) is given by I(k), where:
>
> k = i + sx * j
>
> Is this not enough information to create the mesh?

For structured grids yes. A print out of the coordinates of a UnitSquare shows 
that this is actually the case. If you trust dolfin to not change the 
algorithm for how it creates a UnitSquare you can probably assume this.

> >> It seems _vector.set()
> >> (http://www.fenics.org/pub/documents/dolfin/dolfin-progr-reference/d6/da
> >>4/c lassdolfin_1_1GenericVector.html#bc4ecede55f2846ecadae02426df8f63)
> >> accepts the pointer to each row, is that right?
> >
> > There is no such thing as a row in a vector. It stores the data in a
> > contigouos way (if not distributed).
>
> No idea why it says rows: virtual void dolfin::GenericVector::set
> (const double *block, const uint *num_rows, const uint *const *
> 	rows).

GenericTensor interface is for Vectors, Matrices and Tensors for rank > 2.

If you take a look at GenericVector.set(double* values) you probably have what 
you want. 

> >> Is there a better way of assigning the whole data at once?
> >
> > Maybe but you then have to be sure that your data is arrange in exact
> > same manner as the mesh.
>
> As long as we work with complete rectangular images, which is the
> frequently used case, it should be fine. Given this, how is it
> possible to directly assign the image data pointer to _vector?

See above.

Johan

> >> > Note that this will only work for a piecewise linear lagrange
> >> > FunctionSpace.
> >>
> >> It should be ok, it would be equivalent of having the image function
> >> initialised with Lagrange element, but much more faster.
> >>
> >> > But haveing this in place you can always interpolate to other
> >> > FunctionSpaces. However these interpolations will only be linear of
> >> > course.
> >>
> >> It is not possible to interpolate a function to another space using
> >> non-linear interpolators?
> >
> > If your raw image data is represented by a linear FunctionSpace any
> > interpolation using this space would be linear.
>
> Right, that's yet another problem.
>
>
> -Ali
>
> > Johan
> >
> >> To Anders: Do you think if this is more practical and efficient than
> >> the other discussed approach? Would this be general enough to work
> >> with an arbitrary PDE of any finite element type?
> >>
> >>
> >> -Ali
> >>
> >> > Johan
> >> >
> >> >> -Ali
> >> >>
> >> >> > -Ali
> >> >> >
> >> >> >>> >> Another way of doing this could be by the use of an existing
> >> >> >>> >> FunctionSpace:
> >> >> >>> >>
> >> >> >>> >>   UnitSquare dummy_mesh(1, 1);
> >> >> >>> >>   PoissonFunctionSpace V(dummy_mesh);
> >> >> >>> >>   ImageFunction v(image, V);
> >> >> >>> >>
> >> >> >>> >> Then, in the constructor of ImageFunction, V.element and
> >> >> >>> >> V.dofmap can be used to create another FunctionSpace which has
> >> >> >>> >> a mesh created using image:
> >> >> >>> >>
> >> >> >>> >>   ImageToFunction(Image image, const FunctionSpace& V)
> >> >> >>> >>   {
> >> >> >>> >>     // Create the function space
> >> >> >>> >>     UnitSquare mesh(image.get_size()[0] - 1,
> >> >> >>> >> image.get_size()[1] - 1); FunctionSpace IV(mesh, V.element(),
> >> >> >>> >> V.dofmap());
> >> >> >>> >>
> >> >> >>> >>     // ...
> >> >> >>> >>   };
> >> >> >>> >>
> >> >> >>> >> The problem with this approach is that it involves the use of
> >> >> >>> >> a dummy mesh.
> >> >> >>> >>
> >> >> >>> >> A mesh-independent constructor added to FunctionSpace could
> >> >> >>> >> help. Alternatively, if a (protected) default (empty)
> >> >> >>> >> constructor is added to FunctionSpace,
> >> >> >>> >> ImageFunctionSpace:FunctionSpace can have a mesh-independent
> >> >> >>> >> constructor. However, the FFC-generated function spaces, eg
> >> >> >>> >> PoissonFunctionSpace, still need a mesh.
> >> >> >>> >>
> >> >> >>> >> Hope this makes the problem more clear now.
> >> >> >>> >
> >> >> >>> > Create the mesh and the FunctionSpace subclass inside the
> >> >> >>> > ImageFunction constructor. Neither the mesh nor the function
> >> >> >>> > space need to be visible outside.
> >> >> >>>
> >> >> >>> Again, there is no mesh-free ctor for FunctionSpace, and it
> >> >> >>> doesn't come with a default ctor so that the subclass can
> >> >> >>> implement a mesh-free ctor. I quote the above again:
> >> >> >>
> >> >> >> You don't need a mesh-free constructor for FunctionSpace, see
> >> >> >> above.
> >> >> >>
> >> >> >>> ' A mesh-independent constructor added to FunctionSpace could
> >> >> >>> help. Alternatively, if a (protected) default (empty) constructor
> >> >> >>> is added to FunctionSpace, ImageFunctionSpace:FunctionSpace can
> >> >> >>> have a mesh-independent constructor.'
> >> >> >>
> >> >> >> It's not needed if you do as I suggest.
> >> >> >>
> >> >> >> --
> >> >> >> Anders
> >> >> >>
> >> >> >> -----BEGIN PGP SIGNATURE-----
> >> >> >> Version: GnuPG v1.4.9 (GNU/Linux)
> >> >> >>
> >> >> >> iEYEARECAAYFAkmdMN8ACgkQTuwUCDsYZdGnOQCdFBJSD4FymLnVPbheRt63aJJa
> >> >> >> yyoAn3KDuOmwd8ZX5YR1KucbafvieNBc
> >> >> >> =lpyI
> >> >> >> -----END PGP SIGNATURE-----
> >> >> >>
> >> >> >> _______________________________________________
> >> >> >> DOLFIN-dev mailing list
> >> >> >> DOLFIN-dev@xxxxxxxxxx
> >> >> >> http://www.fenics.org/mailman/listinfo/dolfin-dev
> >> >>
> >> >> _______________________________________________
> >> >> DOLFIN-dev mailing list
> >> >> DOLFIN-dev@xxxxxxxxxx
> >> >> http://www.fenics.org/mailman/listinfo/dolfin-dev