1
$\begingroup$

Any differential operator can be written as a finite difference matrix acting on a vector of values. I recently used it to solve $\Delta A = j$ by writing the Laplacian as a finite difference tridiagonal (second order approximation) matrix. I then used SciPy's spsolve, which solves matrix equations of the form $Ax=b$ by not explicitly inverting the sparse matrix (which blows up really fast).

It seems to me, that because of the direct accuracy estimate you get from the finite difference scheme, and the extreme parallellization that can be achieved with methods like spsolve, this would be an optimal way to solve this equation.

Is this the case (or is it at least a competitive method)? Why are all classes and lectures on numerical methods focusing on all other methods and I can only find sparse mention of the above described method, although it is both intuitive, simple to implement, and at least relatively fast.

I only saw finite difference matrices in light of eigenvalue calculations (in which I believe it is also quite optimal due to extreme optimization of the sparse eigenvalue calculation procedures), and thought of this "integration" technique myself. I might have just overlooked all books and other literature on the subject, hence this question directed to the experts.

Improve this question
$\endgroup$
  • 2
    $\begingroup$ Welcome to SciComp.SE! I'm not sure I follow you -- I would certainly call what you're describing a "normal method" (it's the first thing we teach in any undergraduate course on numerical differential equations). There are indeed other methods (e.g., finite element discretization and preconditioned conjugate gradient methods, respectively) which can perform even better in certain situations. $\endgroup$ – Christian Clason Dec 3 '15 at 13:14
  • $\begingroup$ @Christian Thanks! Well, when googling for methods to solve Schrödinger equations, I always get Numerov or odeint-style solutions presented, which are either quite illogical (Numerov, it's pretty much a guessing game) or black-box style things. I'm aware my view may be biased, I'm just wondering if the method is easily bested or not. For example, I also see FFT solutions to the Laplacian equation frequently, which means infinitely more than an approach with finite differences. $\endgroup$ – rubenvb Dec 3 '15 at 13:18
  • 1
    $\begingroup$ Oh, you're coming from quantum mechanics. Yes, there's a whole zoo of special purpose methods there. If you're looking at the time-dependent Schrödinger equation, you need time integration, which indeed requires some sort of ODEint routine (this is also the case if you go from Laplace to heat equation). Your approach is usually known as FDTD (finite difference time domain). The FFT approach (when applicable, i.e., special potentials and periodic boundary conditions) can be even more efficient, since you don't even need to solve any equation (just pointwise division). $\endgroup$ – Christian Clason Dec 3 '15 at 13:23
  • $\begingroup$ Yeah, time has special methods, but no, the equations and methods I'm talking about are stationary. Vector potential from current, 2D wave eigenfunctions of the Schrödinger equation, etc... The issue I have with the FFT method to calculate the vector potential is that is seems to result in a complex vector potential, which is probably a gauge transformation from a real one, but I'd like to avoid all that post-processing. $\endgroup$ – rubenvb Dec 3 '15 at 13:27
  • $\begingroup$ Then, yes, what you're doing is certainly a standard method. There might be faster/more accurate methods that take advantage of your specific situation (domain, boundary conditions, potential term), but as a default method it makes sense. $\endgroup$ – Christian Clason Dec 3 '15 at 13:32
2
$\begingroup$

If your question is why the word "sparse matrix" never appears in presentations on finite difference (or finite element, or finite volume) solvers for partial differential equations, then that's because that's always implied: it simply makes no sense to store the matrix in any way that doesn't exploit sparsity. That's true for sequential and parallel applications.

Because it's the only reasonable way to do it, it's not usually mentioned.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

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