← Back to team overview

dolfin team mailing list archive

Re: [HG] Add GenericMatrix class. This templated class acts as a base class for dense and sparse matrices. It provides polymorphism without having to use virtual functions.

 

I found a nice compromise solution for dense and sparse matrices. It's
often called the Barton-Nackman trick. 

I've added a base class GenericMatrix which is a template. The classes
DenseMatrix and Matrix (which is sparse) are then derived from
GenericMatrix (this has only been done for DenseMatrix at this stage).
The user can create matrices as usual,

DenseMatrix A;
Matrix B;

and use any of the public member functions of the respective class
(which don't have to be the same). Functions which accept either a dense
or sparse matrix then accept GenericMatrix<class T> as an argument, so
they have to be templated. GenericMatrix only needs to have member
functions that are used by the templated functions. At this stage,
GenericMatrix contains only the functions required for assembly.


For example, it is possible to assemble a dense or sparse matrix by

Mesh mesh;
BilinearForm a;
FEM::assemble(a, A, mesh)
FEM::assemble(a, B, mesh)

Take a look in src/test/main.cpp where a sparse and a dense matrix are
assembled. I've only templated one function in FEM.h for testing. 

The advantages of this approach are:

1) Avoids virtual functions which are called from within loops.
2) No changes at the users end when declaring or assembling matrices.
3) GenericMatrix is simple as it needs to contain only member functions
which are called by functions which accept both dense and sparse
matrices.

The disadvantage is that functions which accept dense and sparse
matrices need to be templated. I don't think that this is such a big
deal, especially since FEM (which is where the function template will be
needed) has been cleaned up now. 

Let me know what you think. If everyone thinks that it's OK, I'll go
ahead an complete the implementation. Also, any preference for "Matrix"
or "SparseMatrix" as names for a sparse matrix? Most problems involve
sparse matrices only, so I have a very slight preference for Matrix.  

Garth

On Mon, 2006-04-24 at 16:02 +0200, DOLFIN wrote:
> One or more new changesets pushed to the primary DOLFIN repository.
> A short summary of the last three changesets is included below.
> 
> changeset:   1893:545d61886439a362a658b0a48e5ef2808441320c
> tag:         tip
> user:        "Garth N. Wells <g.n.wells@xxxxxxxxxx>"
> date:        Mon Apr 24 15:54:19 2006 +0200
> files:       src/kernel/fem/dolfin/FEM.h src/kernel/la/DenseMatrix.cpp src/kernel/la/Makefile.am src/kernel/la/Makefile.in src/kernel/la/dolfin/DenseMatrix.h src/kernel/la/dolfin/Makefile.am src/kernel/la/dolfin/Makefile.in src/test/main.cpp
> description:
> Add GenericMatrix class. This templated class acts as a base class for dense and sparse matrices. It provides polymorphism without having to use virtual functions.
> 
> 
> changeset:   1892:779ec524d4b0c72ec15c7f77f09015645d696f62
> user:        "Anders Logg <logg@xxxxxxxxx>"
> date:        Sun Apr 23 07:44:47 2006 -0500
> files:       src/bench/ode/reaction/main.cpp src/bench/ode/wave/bench.log src/bench/ode/wave/bench.py src/bench/ode/wave/main.cpp
> description:
> Fixes for second multi-adaptive benchmark (wave equation).
> Seems to work fine now:
> 
> CPU time cG(1):  188.093 s
> CPU time mcG(1): 45.4958 s
> 
> 
> changeset:   1891:abcc54c4a3f446f20c29cf8eb032bbf5f9fbaa19
> user:        "Anders Logg <logg@xxxxxxxxx>"
> date:        Sat Apr 22 15:41:44 2006 -0500
> files:       src/bench/ode/reaction/bench-domain.log src/bench/ode/reaction/bench-domain.py src/bench/ode/reaction/bench-tol.log src/bench/ode/reaction/bench-tol.py src/bench/ode/reaction/benchutil.py src/bench/ode/reaction/main.cpp src/bench/ode/reaction/parameters-bench.xml
> description:
> Recompute benchmarks with improved adaptivity.
> 
> Results look very good:
> 
> bench-tol: cg/fixed-point
> ---------------------------------
> 
> tol     N       Error           CPU time        Steps   Iterations      Index
> --------------------------------------------------------------------------------------
> 1.0e-06 1000    2.345e-05       28.050          117089 (1)      4.000 (0.000)   1.0
> 5.0e-07 1000    1.173e-05       39.520          165586 (1)      4.000 (0.000)   1.0
> 1.0e-07 1000    2.344e-06       71.940          370254 (1)      3.000 (0.000)   1.0
> 5.0e-08 1000    1.172e-06       101.730         523615 (1)      3.000 (0.000)   1.0
> 
> bench-tol: mcg/fixed-point
> ---------------------------------
> 
> tol     N       Error           CPU time        Steps   Iterations      Index
> --------------------------------------------------------------------------------------
> 1.0e-06 1000    1.771e-05       14.240          1922 (5)        3.990 (1.498)   95.266
> 5.0e-07 1000    1.091e-05       23.300          1912 (9)        4.822 (1.544)   138.240
> 1.0e-07 1000    1.921e-06       48.050          1929 (7)        4.905 (1.594)   142.619
> 5.0e-08 1000    8.979e-07       49.770          1917 (7)        4.131 (1.680)   172.388
> 
> 
> bench-domain: cg/fixed-point
> ---------------------------------
> 
> tol     N       Error           CPU time        Steps   Iterations      Index
> --------------------------------------------------------------------------------------
> 1.0e-06 1000    2.345e-05       28.100          117089 (1)      4.000 (0.000)   1.0
> 1.0e-06 2000    2.233e-05       64.840          117091 (1)      4.000 (0.000)   1.0
> 1.0e-06 4000    2.166e-05       101.290         117090 (1)      4.000 (0.000)   1.0
> 1.0e-06 8000    2.219e-05       175.110         117089 (1)      4.000 (0.000)   1.0
> 1.0e-06 16000   2.242e-05       327.680         117089 (1)      4.000 (0.000)   1.0
> 
> bench-domain: mcg/fixed-point
> ---------------------------------
> 
> tol     N       Error           CPU time        Steps   Iterations      Index
> --------------------------------------------------------------------------------------
> 1.0e-06 1000    1.771e-05       13.640          1922 (5)        3.990 (1.498)   95.266
> 1.0e-06 2000    1.680e-05       17.310          1923 (5)        3.991 (1.154)   140.461
> 1.0e-06 4000    1.643e-05       24.000          1920 (6)        3.991 (1.004)   184.969
> 1.0e-06 8000    1.691e-05       33.670          1918 (5)        3.991 (1.004)   218.782
> 1.0e-06 16000   1.697e-05       57.850          1919 (5)        3.990 (1.004)   239.992
> 
> 
> Newton benchmarks are temporarily disabled: takes longer
> time anyway and there seems to be a bug in the multi-adaptive
> Newton solver (takes more iterations than fixed point). Also
> changed kmax to 1e-3 in bench-domain (works better with new
> choice of threshold). kmax is now same for both benchmars and
> set in parameters-bench.xml.
> 
> 
> -------------------------------------------------------
> For more details, visit http://www.fenics.org/hg/dolfin
> 
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/cgi-bin/mailman/listinfo/dolfin-dev




Follow ups

References