Extracting DoF indice from FEEvaluation / Distributing from FEEvaluation to row of a matrix

[Michał Wichrowski]

Dear all,

I need to construct a sparse matrix from MatrixFree operator. I want to do it in a similar way to present in step-37 for computing diagonal. The only difference will be that the computed values from  FEEValuation and FEFaceEvaluation have to be distributed to the specified row of a matrix. The code will look like, I need to get local_dof_indices somehow:

// Loop is executed without threads

for (unsigned int cell=0; cell<n_cells; ++cell)

{

phi.reinit (cell);

for (unsigned int i=0; i<phi.dofs_per_cell; ++i)

{

for (unsigned int j=0; j<phi.dofs_per_cell; ++j)

phi.submit_dof_value(VectorizedArray<number>(), j);

phi.submit_dof_value(make_vectorized_array<number>(1.), i);

// local integration operator:

integrate_cell(phi,cell);

phi.integrate (false, true);

// Here I have phi filled with a row of the matrix

// now I need to know global indices of dofs...

// And I have no clue how to obtain local_dof_indices needed below.

for (unsigned int j=0; j<phi.dofs_per_cell; ++j){

for (unsigned int v=0; v<this->data->n_components_filled(cell); ++v)

matrix(local_dof_indices[i][v],local_dof_indices[j][v] )+=phi.get_dof_value(j)[v];

}

}

}

The code for FEFaceVAlues will be very similar.

Moreover, I'm working with a vector-valued problem so local_dof_indices will need 3 indices ( local dof number, index in vectorization and component of the vector)

It would be great it FEEvaluation had a function for distributing to a matrix so that two loops are replaced with phi.distribute_local_to_global (matrix, i);

I only need to do that somehow, it does not have to be very efficient. There are three main purposes for having that function:

1) coarse matrix for direct solver

2) initialization of preconditioner based on sparse matrix

3) system matrix for validating code

Having this functionality would allow me to write a class for an arbitrary matrix-free operator that form local integration formula will be able to automatically:

1) compute matrix-vector product

2) compute diagonal

3) compute sparse matrix of operator.

The operator will be specyfied by defining a function :

void integrate_cell(FEEValuation & phi, unsigned int cell);

And similar ones for face and boundary in case of DG. I think it will be possible for both vector and scalar based 2nd order operators. Also non-square operator could be possible. And using this, the block operators may be constructed (eg. Stokes).

Best,

Michał

[Martin Kronbichler]

Dear Michal,

To answer the first request, we can provide an iterator from `MatrixFree` via `MatrixFree::get_cell_iterator(cell_index, v)`:

https://www.dealii.org/developer/doxygen/deal.II/classMatrixFree.html#adb12e97b3905c1cb12d6f90527604476

That would work for the cell part. However, we are currently lacking an equivalent functionality for the faces (it should be easy to add, though).

Regarding the second request to let `FEEvaluation` handle the transfer into a matrix: Peter Munch already did some work into that direction. I'm adding him in CC so we can discuss together how to best provide a solution for your problem. In the meantime, I'm adding an issue as well to get that information.

Best,
Martin

--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see https://groups.google.com/d/forum/dealii?hl=en
---
You received this message because you are subscribed to the Google Groups "deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email to dealii+un...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/e3f68a4b-4361-4167-9514-d89e746bfd8b%40googlegroups.com.