### Matrix Class

137
views
0
10 weeks ago by
Hi all,
I'm Valerio and I'm using FEniCS  fenics/2016.1.0 in C++. I have a problem with the manipulation of an object of the class Matrix. In my main I have on object defined:
Matrix vss_fenics;​
How can I set the values and the dimensions of the Matrix vss_fenics?
Community: FEniCS Project

0
10 weeks ago by
Hi,

Just a few hints. Hope this will get you started.
To initialze the dimensions and the sparsity pattern of a Matrix object, use the TensorLayout class (to see how this is used, consult the AssemblerBase class). In fact, if your matrix has a layout similar to a known Form, let FEniCS do the work and use the AssemblerBase::init_global_tensor() method.

To set the values in a matrix, use Matrix::set_local(). Don't forget to finalize the manual assemblage of your Matrix object using Matrix::apply("insert") or Matrix::apply("add")
I tried to use Matrix::setrow becouse I don't understant what it is dolfin::la_index. My test code was:
Matrix A;

std::vector<double> values1(3);
values1[0]=10;
values1[1]=11;
values1[2]=12;
std::vector<double> values2(3);
values2[0]=20;
values2[1]=21;
values2[2]=22;
std::vector<std::size_t> columns(3);
columns[0]=1;
columns[1]=2;
columns[2]=3;
A.setrow(0,columns,values1);
A.setrow(1,columns,values2);​

but I had a segmentation fault. My aim is to copy an Eigen::Matrix in a dolfin::Matrix. Likely the segmentation fault is due to the fact that the size of the Matrix A isn't set. How do I have to correct the test code? thank you for your help

written 10 weeks ago by Valerio
Indeed, the segfault is because size and layout of your matrix are not initialized. What info do you have a priori about the Matrix? Is it somehow related to a certain Form?

By the way, dolfin::la_index is just an index type being compatible with the linear algebra backend.
written 9 weeks ago by Jakob Maljaars
The info that I know about the Matrix is that it is double and it isn't related to a certain Form. The dolfin::Matrix in fact have to store the values of an Eigen::Matrix of double, with Dynamic rows and Dynamic columns
written 9 weeks ago by Valerio
Then I think the only option is to study the TensorLayout class to initialize your matrix properly. Although you might have a look at this post on the old forum.

By the way, if your Matrix is not related to a certain Form why can't you use Eigen as the backend to do the matrix manipulation and why do you really want to cast it into a dolfin::Matrix?
written 9 weeks ago by Jakob Maljaars
In my code I have
Matrix A;
Vector f;
std::vector<std::shared_ptr<const DirichletBC>> dirichlet_matrici({dirichlet});
assemble_system(A,f,a,L,dirichlet_matrici);
solve(A, *(u.vector()), f);​

and the line

solve(A, *(u.vector()), f);

resolves the linear system
$Au=f$Au=ƒ

I don't want to resolve that system, but I want to resolve the system
$A_ru_r=f_r$Arur=ƒ r

With
$A_r=Y^TAY$Ar=YTAY   $f_r=Y^Tf$ƒ r=YTƒ
The problem is that A is a dolfin::Matrix, f is a dolfin::Vector but Y is a Eigen::Matrix. Then I need to convert or all in dolfin::Matrix (as i tried to do) or all in Eigen::Matrix (I don't know how, i searched and i didn't find anything)

written 9 weeks ago by Valerio
Clear. So in which FunctionSpace does $u_r$ur live?
written 9 weeks ago by Jakob Maljaars
My idea was to calculate  $u_r$ur just as a simple vector, then calculate the new vector
$u_2=Yu_r$u2=Yur
and then convert the vector  $u_2$u2 in a function of the same space of  $u$u of the system
solve(A, *(u.vector()), f);​
written 9 weeks ago by Valerio