dolfin team mailing list archive

[James Avery <question158315@xxxxxxxxxxxxxxxxxxxxx>]

Question #158315 on DOLFIN changed:
https://answers.launchpad.net/dolfin/+question/158315

James Avery posted a new comment:
1. The dolfin::Matrix uses the default matrix backend, which can be a
number of different systems. It can be a PETSc, Trilinos Epetra, MTL4,
or uBLAS. This means that you can't count on the backend being uBLAS
unless you specifically specify it to be.

If you really want an uBLAS matrix (this means your program won't be
able to use distributed memory for the matrix when PETSc or Trilinos is
enabled), you can just do so like this:

uBLASSparseMatrix T, M;

assemble(T, stiffness);
assemble(M,mass);

2. If you directly use the uBlasSparseMatrix type, you can simply call
M.invert() (This is a void operation, so don't try to assign the result
to anything). However, invert() is only defined for uBLAS matrices. Keep
in mind that the inverse of a sparse matrix is not sparse, so you can
only practically do this on small meshes. Do you really need the
inverse, or can you just solve the linear systems directly?

3. With uBLASMatrix, you simply have M_{ij} = M(i,j). However, the
general Matrix class has to assume a distributed memory sparse matrix
type like Trilinos Epetra or PETSc. This means that random access to
matrix element is an incredibly inefficient operation, and is not
implemented. Instead you have to extract a row of nonzero elements at a
time:

 vector<uint> index;
 vector<double> value;

 M.getrow(i, index, value); // index now contains the j for which M(i,j) != 0, and value[k] = M(i,index[k]).
 ...
 M.setrow(i, index, value); // commits the row to the distributed sparse matrix.

Hope this helps.

Cheers,
James

-- 
You received this question notification because you are a member of
DOLFIN Team, which is an answer contact for DOLFIN.