[SciPy-dev] FE+sparse utilities

[Robert Cimrman]

Hi Lisandro,

Lisandro Dalcin wrote:
> On 12/11/06, Robert Cimrman <cimrman3 at ntc.zcu.cz> wrote:
>> I would greatly appreciate any feedback,
> 
> Hi, Robert... I will take a look tomorow to your code.  I am very
> interested in this. I am developing a big, complete wrapper to PETSc
> libraries (petsc4py) and this subject is a priority for me. I

I have always wanted to have PETSc accessible in Python so that I could 
try it easily, great!

> currently try to avoid assembles in the Python side. However, the most
> efficient way I've found is to write in the C side a new (Python) type
> object containing a pointer to the real matrix structure, and
> implement the mapping protocol. This way, you can almost skip Python
> bytecode.

Feutils is a very simple package good for:
1. given finite element connectivity of degrees of freedom (= mesh 
connectivity, if 1 DOF per node), prepare a global CSR matrix with 
nonzero pattern exactly corresponding to that connectivity. Thus, when 
assembling, no reallocations, data shifts etc. are needed;
2. given a bunch of element contributions, assemble them into the global 
matrix.

The actual matrix allocation as well as element contribution assembling 
are done in C (via swig), too. But I use standard numpy (array) and 
scipy (CSR sparse matrix) data types to store the data. The key point is 
to compute in one C call many (at least 1000) local element 
contributions, and them assemble all of them also in one function call - 
that way the time spent in Python is minimized.

The C code assumes 32-bit integers for indices and 64-bit floats for 
values. I can make a 64-bit version, too, if there is need.

I will add this e-mail's content to the web page 
(http://ui505p06-mbs.ntc.zcu.cz/mymoin/SFE).

Thanks for your interest,
r.