dolfin team mailing list archive

["Martin Sandve Alnæs" <martinal@xxxxxxxxx>]

2007/3/9, Garth N. Wells <g.n.wells@xxxxxxxxxx>:
> Anders Logg wrote:
> > On Fri, Mar 09, 2007 at 07:10:31PM +0100, Garth N. Wells wrote:
> >> Anders Logg wrote:
> >>> On Fri, Mar 09, 2007 at 03:07:43PM +0100, Garth N. Wells wrote:
> >>>> Why do we need a class AssemblyMatrix? GenericMatrix/GenericVector
> >>>> should provide the necessary  virtual functions for assembly.
> >>>>
> >>>> Garth
> >>> Yes, definitely. The assembler assumes it gets a GenericTensor that
> >>> has a function
> >>>
> >>>     void add(real* block, uint* size, uint** entries)
> >>>
> >>> (entries is probably a bad choice of name here, let's change it)
> >>> This makes it possible to implement just one assembly function and
> >>> keep the Assembler class small.
> >>>
> >> Then GenericVector and GenericMatrix should then derive from GenericTensor?
> >
> > Yes. So far it is just two functions that are needed, init() and add().
> >
> >>> Anyway, the class AssemblyMatrix is just a temporary implementation of
> >>> something that implements the GenericTensor interface. My suggestion
> >>> is that we overload the add() function above in both GenericMatrix and
> >>> GenericVector. The overloaded functions will map to the correct add()
> >>> function in the respective classes. For GenericMatrix, we'll have
> >>>
> >>>     void add(real* block, uint* size, uint** entries)
> >>>     {
> >>>         add(block, size[0], entries[0], size[1], entries[1]);
> >>>     }
> >>>
> >>> So, even if AssemblyMatrix could be of some use (it seems to be a fast
> >>> option for assembly, the STL vector of maps beat PETSc in at least one
> >>> benchmark I made)
> >> Don't draw conclusions too quickly. STL vector of maps is fast to
> >> assemble, but is close to useless for linear algebra. The current
> >> implementation for the assembly uBLAS matrices use something like STL
> >> vector of maps  (fast to assemble but slow for linear algebra) and
> >> converts it to compressed row storage (fast for linear algebra). The
> >> conversion can take some time. A STL vector of integer maps should be
> >> created and then used to initialise a particular sparse storage format
> >> (before first use only). Then we wouldn't have to convert matrix formats.
> >
> > Do you mean this should be done in DofMap?
> >
>
> Yes. The STL vector of maps is cheap to construct, and can be used to
> initialise the non-zero layout of a sparse matrix.

You should use vector<set<int>>, not vector<map<int>>, if the goal is
only to build the sparsity pattern, as this gives lower memory
overhead since the set is a map with no value.

Anders: I have an implementation of this in
pycc/MatSparse/CRSGraph.cpp and pycc/FEMAssembler/CRSAssembler.cpp,
feel free to use any suitable parts of it.

In my brief tests, using vector<set<int>> and using it for allocating
the pycc::CRS matrix is about as fast as Epetra_CRSGraph/Epetra_CRS.

martin