Computational cost of numerical methods for PDEs

Say I need to solve a PDE numerically. Depending on the problem and the numerical method chosen, I can usually see some issues coming: implementation issues (e.g. boundary conditions, parallelization), discretization issues (e.g. the stabilization of finite elements for hyperbolic problems is not funny while finite volumes "just work", the "opposite" for elliptic problems), time integration...

I usually base the initial choice of the discretization on these "intuitions" only. Still, I would love to base them on facts and proof regarding the question: how much computational effort would I have to put to solve my problem to a desired solution error (i.e. dollars/error)?

• Is there any resource I could use to weight the suitability of a method with this criterion?

Although one could compute the complexity of the methods one still would need a priori error estimators that, if they exist at all, would have different "sharpness", depending on the method. Is there a non-empirical way of finding this out? If testing all possible methods cannot be avoided, how would one go on about it? There are standard numerical testcases for most fields, would they be meaningful? Are there any public databases comparing implementations for these testcases?

Disclaimer: This should not be the only criterion, implementing finite difference for complex boundaries is "doable" but a pain. Still, is it worth it to go through that pain? If the cost of using FD would be much lower than the alternatives, for some people, it might.

• The answer may depend on very many things. Not just the nature of your PDE, but also boundary conditions, initial conditions, domain, coefficients. Then there is the computer architecture, compilers and so forth. Additionally, the quality of your implementation, from choice of language down to fine details. Each of these things can make a large difference under the right circumstances. – David Ketcheson Aug 30 '12 at 17:26
• I forgot to mention: the desired accuracy. Almost certainly, the optimal method for low accuracy will not be optimal for higher accuracy. – David Ketcheson Aug 30 '12 at 17:34
• @DavidKetcheson that's a nice remark! The optimal method depends on the desired accuracy indeed! Standard testcases exist in most fields, however the quality of the implementation has a huge impact. Still, one should prefer to chose a better algorithm than to optimize to the bare metal. It would be nice to have a way of measuring numerical methods to guide this decision. – gnzlbg Aug 30 '12 at 22:51

You are asking for something overly ambitious if you hope to be able to reliably get within a constant factor (e.g. a couple orders of magnitude) of an optimal method, even just on one type of computer. If someone claimed to have a business analysis technique that could tell you the value of any company in five years, they would be accused of selling snake oil. A uniform optimality analysis for PDE solution methods is similarly unattainable. You will find that as you construct a classification methodology, you will have to severely limit the scope of problems and methods that you can treat, or pay the price of such huge uncertainty factors as to make the method analysis useless. The lack of a simple answer is one reason why "optimal" PDE solution methods is still an active research area.

• There are a lot of standard numerical testcases in every field and the number of papers that show "cost/error" is very low. There are even less papers comparing "cost/error" against other methods. When people say that their 8th-accurate DG method is faster than their 2nd-order accurate DG method I cannot but wonder if it would also be faster than my second order FV/FE method. It would be good if there were some resource where the cost of solving standard testcases could be compared, but as I am learning from the answers this might just be too difficult. I hope it is doable as it would be nice. – gnzlbg Aug 30 '12 at 22:58

I don't think is possible given that the nature of solutions to PDEs varies so much. However, a general heuristic you can have in your head is something like:

1. Find out what type of equation (Eliptic, Parabolic, Hyperbolic).
2. Look at the coefficients in the problem / the Jacobian. Is the problem stiff? If so you'll want some kind of implicit method (i.e. Crank-Nicolson), otherwise you can get away with some explicit or predictor corrector method.
3. Look at your geometry. Is your geometry simple? Then try a finite difference method. When the geometry is complex, use a finite element method.

This process works for most problems. Parabolic stiff equation on a torus? Give an implicit FEM solver a try. However, this won't work for all problems, you might have to go back and try another method. Or worse: if the problems hard enough, getting a computational method to solve the PDE is a research topic in its own right. That's why computational PDEs is still an active field of research.

Also remember that any time spent coding counts towards your time spent solving. If you pick a method which takes 5 minutes to write but solves it in 2 days, it is still pragmatically "faster" than a method which takes 3 days to write but solves it in minutes (and you can get a lot more work done). So don't get stuck worrying about "the optimal method" for your problem unless you're going to be solving it for years or are a researcher in numerical methods.

The computational cost depends not just on the number of DOF, but also on how efficient you can solve the resulting (nonlinear?) equation system.

As an extreme case, the time harmonic Maxwell equations for a structure periodic in x- and y-direction can be expanded into Fourier modes (in x- and y-direction). The resulting system has very few degrees of freedom. It can be solved by various methods, but it ultimately boils down to solving a two point boundary value problem by some sort of shooting method. If you resort to matrix methods for the stable solution of this problem, the runtime scales with the third power of the area of the periodic unit cell. If you adapt a simple FDTD scheme to the solution of the time harmonic Maxwell equation instead, you end up with many degrees of freedom. You don't need to solve an equation system in this case, and the runtime scales linear with the area of the periodic unit cell. (Both approaches can be modified to overcome their main limitations, but these modifications are also not without issues.)

Note that the time harmonic Maxwell equations are elliptic, but time domain problem solved by FDTD is hyperbolic.

• So if I understood your answer correctly, what you are saying is that choosing a different numerical method can have a large impact on the computational cost, is that right? – gnzlbg Aug 30 '12 at 23:08
• @gnzlbg No, at least that's not the message I wanted to bring across. I had initially tried to describe some approaches to overcome the main limitations of the two described solution methods, but later decided against it. It was something along the line of: "If you use your physical intuition to devise an iterative scheme for the stabilization of the shooting procedure instead, you end up with a poorly understood, unbelievable fast and unreliable solver, which you basically have to tune manually for every new problem instance." The message is more: "There are unbelievable many tradeoffs." – Thomas Klimpel Aug 31 '12 at 8:10