dolfin team mailing list archive

[Johan Hake <hake@xxxxxxxxx>]

On Wednesday 18 February 2009 04:47:13 N wrote:
> ----- Original Message ----
>
> > From: Johan Hake <hake@xxxxxxxxx>
> > To: dolfin-dev@xxxxxxxxxx
> > Cc: N <tenpoundwalleye@xxxxxxxxx>
> > Sent: Tuesday, February 17, 2009 4:42:46 PM
> > Subject: Re: [DOLFIN-dev] access to CSR data in python
> >
> > On Monday 16 February 2009 23:29:48 N wrote:
> > > Take a simple example like
> > >
> > > from dolfin import *
> > > mesh = UnitSquare(4, 4)
> > > V = FunctionSpace(mesh, "CG", 1)
> > > v = TestFunction(V)
> > > u = TrialFunction(V)
> > > a = dot(grad(v), grad(u))*dx
> > > A = assemble(a)
> > >
> > > The matrix A (with boost as the linear algebra base) gives the
> > > following for attribute data: data(self) ->
> > > std::tr1::tuple<(p.q(const).std::size_t,p.q(const).std::size_t,p.q(cons
> > >t).d ouble,int)>
> > >
> > >     Return pointers to underlying compressed storage data. See
> > >     GenericMatrix for documentation.
> > >
> > >
> > > How can this be used to gain access to the column data, row pointer,
> > > and data elements in the matrix?  I would like to have access to this
> > > data in Numpy arrays without copying the data.  Is this possible?
> >
> > You now can. Sorry that I did not made it to the release ;)
> >
> > You can now do:
> >   >>> rows, cols, values = A.data()
> >
> > Here rows, cols and values are numpy arrays based on pointers to the
> > underlaying data structures in c++.
> >
> > It works for the linear algebra backends that supports the feature in
> > C++, uBLAS and MTL4. Error is raised for the other backends.
> >
> > Also added direct access to Vector data through the same data() function:
> >   >>> values = v.data()
> >
> > This is similare to v.array(), but with no copying. This feature is also
> > only supported for the uBLAS and MTL4 backends.
>
> Perfect.  Works well.

Good, I have to sort out a unit test that fails on 64bits(!?) but not on the 
other buildbot...

> I'm curious how implemented this.  I worked long and hard with both swig
> and direct python extensions to no avail.  Can you point me to your
> changes?

Sure, anyone interesting in swigs magical alleys I would gladely help! ;) 

Take a look at the web interface of our hg repository:

  <http://www.fenics.org/hg/dolfin>

There you see that the changes are done in dolfin/swig/dolfin_la{post,pre}.i.

The "standard" procedure when you want to wrap a function to python when it is 
not done automaticly with swig is to %rename/%ignore it first. This is often 
done in the XXX_pre.i files. Then we %extend the swig generated c++ and/or 
python proxy class in XXX_pre.i, with the actuall function that will do the 
correct wrapping for us. In this case I do this by defining a macro, so the 
readability goes down a bit, but I hope you get the points.

Btw: The distinction between pre and post is that pre files are read by swig 
before the class declarations and post files are read after.

Cheers!

Johan