1
$\begingroup$

I am new to solving PDEs, but have been looking at different implementations of finite difference and finite volume schemes. One thing I have noticed in different implementations is that some implementations use sparse matrices while other implementations--of the same algorithms--use basic loops.

I can understand that for some solvers, using a matrix is necessary, such as implicit solvers that have a linear solve step. Also, some of the confusion may stem from the fact that elliptic pdes are presented first. Since linear elliptic pdes are solved using differentiation matrices, this might explain some the attention on sparse matrices.

I was trying to understand which approach was more efficient--sparse matrices or basic loops? If the goal of practicing writing PDE solvers is to subsequently write codes for real problems, then I would like to learn best practices for writing efficient solvers.

Can anyone tell me which method, sparse matrices or loops, is more efficient method when dealing with larger problems with fine discretizations? Of larger problems we run into issues like breaking up the problem into separate pieces that can be parallelized. We also run into memory issues, etc. I have not seen any of these issues addressed in the simple demos I have found.

I am not sure if the specific language makes a difference, but I generally write in Julia.

Improve this question
$\endgroup$
2
  • 6
    $\begingroup$ It depends on what solvers you use and the problem size. For multidimensional problems with <250k unknowns, explicitly forming the matrix and using a sparse direct solver is often the most efficient approach. For larger problems of >1M unknowns and especially in 3D, you'll want iterative solvers for which the ability to perform matrix-vector multiplication is more important than having the matrix itself. (Search for "matrix-free".) You need to know about both and you should choose a software framework that makes it easy to switch. $\endgroup$
    – Daniel Shapero
    May 20, 2022 at 0:54
  • $\begingroup$ @DanielShapero YES, this is very helpful. This is the kind of information I was looking for. You are correct that I will always need to work through a framework, and that the way I write the solver will depend on how the framework handles solvers. But at least I know that when I write practice solvers, I should write a version for <250k and for >1M. This really helps to narrow down how I will write the basic scaffolding code. Thanks again. $\endgroup$
    – krishnab
    May 20, 2022 at 1:42

1 Answer 1

Reset to default
4
$\begingroup$

There are aspects of software engineering involved. Generally speaking a framework that separates the solvers from the particular physical model makes it really easy to take advantage of the full spectrum of mathematical research when it comes to solving large linear systems. If your software is set up this way, you can easily switch between iterative solvers, direct solvers, sequential/parallel, block based ... (you get the idea).

When it comes to pure performance the matrix-free methods and hand-written loops may in some cases be more efficient, but you loose that flexibility.

Improve this answer
$\endgroup$
1
  • $\begingroup$ Another important factor is the use of preconditioners for iterative methods. Some preconditioners (e.g. incomplete LU factorization) require the matrix in sparse form. $\endgroup$
    – Brian Borchers
    May 24, 2022 at 3:43

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.