I am wondering how Dirichlet boundary conditions in global sparse finite element matrices are actually implemented efficiently. For example lets say that our global finite element matrix was:
K=[5 2 0 -1; 2 4 1 0; 0 1 6 3; -1 0 3 7]
And right-hand side vector b:
b=[b1; b2; b3; b4]
Then to apply a Dirichlet condition on the first node (x1=c) we would zero out the first row, put a 1 at K11, and subtract the first column from the right-hand side. For example our system would become:
K=[1 0 0 0; 0 4 1 0; 0 1 6 3; 0 0 3 7]
Right-hand side:
b=[c; b2-2*c; b3-0*c; b4+1*c]
This is all well and good in theory, but if our K matrix is stored in compressed row format (CRS) then moving the columns to the right-hand side becomes expensive for large systems (with many nodes being dirichlet). An alternative would be to not move the columns corresponding to a Dirichlet condition to the right-hand side, i.e. our system would become:
K=[1 0 0 0; 2 4 1 0; 0 1 6 3; -1 0 3 7]
b=[c; b2;b3;b4]
This however has a major draw back in that the system is no longer symmetric and so we could no longer use preconditioned conjugate gradient (or other symmetric solvers). As far as I know there are a lot of commercial and non-commercial codes using csr format and solving boundary condition problems and
I am assuming that there is probably some standard efficient method used in these commercial or non-commercial codes that solves this problem (obviously I not expecting people to know all the inner workings of every commercial finite element solver, but this problem seems fundamental enough that someone likely has worked on such projects and could provide guidance). I would be grateful if you could help me with detail.
P.S
I am solving the poisson equation with finite element method using simplex triangular element. I stored the non_zero value of global stiffness matrix in compressed row storage format. Also I'm using preconditioned conjugate gradient method parallelized with MPI for solving the system of linear equations. The bottleneck of my work is how to apply the boundary condition in this format.
- Badges
- Science topic