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>
<https://dealii.org/9.0.0/doxygen/deal.II/classVectorizedArray.html>(), 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ł
--
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 [email protected].
To view this discussion on the web visit
https://groups.google.com/d/msgid/dealii/e3f68a4b-4361-4167-9514-d89e746bfd8b%40googlegroups.com.
