[Trilinos-Users] [EXTERNAL] Block Vectors and Matrices

[Hoemmen, Mark]

On 4/26/16, 10:07 AM, "Trilinos-Users on behalf of trilinos-users-request at trilinos.org" <trilinos-users-bounces at trilinos.org on behalf of trilinos-users-request at trilinos.org> wrote:

>Date: Fri, 22 Apr 2016 11:07:56 +1000
>From: Cihan Altinay <c.altinay at uq.edu.au>
>To: "trilinos-users at trilinos.org" <trilinos-users at trilinos.org>
>Subject: [Trilinos-Users] Block Vectors and Matrices
>Message-ID: <5719796C.8060709 at uq.edu.au>
>Content-Type: text/plain; charset=utf-8; format=flowed
>
>Dear Trilinos Devs,
>
>Further to my initial queries in February I have now finally started in 
>earnest trying to use the (experimental) Tpetra block matrix and vector 
>classes and I am encountering some issues.
>
>To give some context, our project (escript [1]) is a FEM tool which 
>currently comes with its own solvers, however we'd like to incorporate 
>Trilinos as much as possible in the future to benefit from the wealth of 
>preconditioners and solvers.
>
>To date I have successfully integrated the Tpetra stack with Belos and 
>Amesos2 solvers, as well as Ifpack2 & MueLu preconditioners for single 
>DOF per node matrices (i.e. block size 1). Now I am trying to add 
>support for block sizes > 1 (always constant-sized square blocks) using 
>Tpetra's BlockCrsMatrix & BlockVector, while keeping the interface as 
>similar as possible.
>
>Here are the questions/issues:
>
>0) Is there an example for solving a block system using the experimental 
>classes somewhere?

There are tests in Tpetra and Ifpack2, but no public example codes.  Remember that BlockCrsMatrix implements Operator and RowMatrix, so you can give it to Belos.  Ifpack2::Relaxation does block Jacobi / Gauss-Seidel / SOR if you give it a BlockCrsMatrix.  There's also an experimental block ILU in Ifpack2.

We reserve the right to make backwards-compatible interface changes as long as these classes live in the "Experimental" namespace.

>1) I have been initializing the "regular" vectors using ArrayViews from
>local std::vectors. The BlockVector does not have a constructor that 
>takes an ArrayView - do I need to create a regular Vector first, then 
>pass that to the BlockVector constructor or is there a more direct way?

You can get the MultiVector out, then work with the Kokkos::View.  See MultiVector documentation.

>2) Unlike regular CrsMatrix/Vector the Block versions appear to take
>references to graph, map and vectors rather than RCPs. This is not 
>really a problem but is there a reason for this?

It reflects a gradual movement towards interface simplification.  Reference counting is handled inside, and operator= and copy constructor do shallow copy with reference counting.

>3) The block matrix does not appear to support resumeFill() or
>fillComplete(). Are these not required?

The Graph must be fill complete.  See public documentation of BlockCrsMatrix.  You may change values at any time.

>4) Is there a MatrixMarket interface for the block matrix yet? That
>would ease debugging since I could dump to file...

See tpetra/core/src/Tpetra_Experimental_BlockCrsMatrix_Helpers_decl.hpp.

mfh