Getting the system matrix for a Poisson problem in C++
5 months ago by
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
to get the block of the matrix after assembling it from the bilinear form
A.get(double *block,size_t m,const int *row,size_t n,const int*col)
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.
Community: FEniCS Project
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']
Please login to add an answer/comment or follow this question.