Getting the system matrix for a Poisson problem in C++


187
views
0
5 months ago by
Hello,

I would like to use FEniCS as the initial interface to get my system matrix from a physical problem. Let's say I have a Poisson problem. The form file would be the same as here and the the main file would be here . with the difference that I use a quadrilateral element instead of a triangle element. But how do I get the values of the matrix into three arrays (the row pointer, column index and the value array) for the CSR format. I understand that I can use
A.get(double *block,size_t m,const int *row,size_t n,const int*col)​
to get the block of the matrix after assembling it from the bilinear form
assemble(A,a);​
but that doesnt seem to give me the values that I want, which leads to my next question.

which is, I would expect that if I have a square mesh with quadrilateral elements instead of triangle elements, I should get the system matrix same as that of the finite difference with a 5 point stencil. Am I wrong here ? Sorry, I am not that adept at Finite element method.

Thank you.
Community: FEniCS Project
1
Some possibly-helpful answers for importing matrices into Matlab are given in this post:

https://www.allanswered.com/post/nqzoo/
written 5 months ago by David Kamensky  

1 Answer


2
5 months ago by

Here's python code to create assemble the stiffness matrix and save it in CSR format to a file.

import numpy as np
from fenics import *
parameters['linear_algebra_backend'] = 'Eigen'
mesh = UnitSquareMesh(10, 10)
V = FunctionSpace(mesh, 'Lagrange', 1)
u = TrialFunction(V)
v = TestFunction(V)
a = (inner(grad(u), grad(v)) + u*v)*dx
A = assemble(a)
Aspa = as_backend_type(A).sparray()
np.savez_compressed('stiffness-matrix.npz', data=Aspa.data, indices=Aspa.indices, indptr=Aspa.indptr)
​


Later you can read the file with

stiffmat = np.load(file1)
data = stiffmat['data']
indices = stiffmat['indices']
indptr = stiffmat['indptr']​
Thanks, but is there any way to do it in C++ , using the available Matrix class?
written 5 months ago by dommath  
Please login to add an answer/comment or follow this question.

Similar posts:
Search »