Trefethen & Bau & MATLAB & Julia: Iterative methods


Toby Driscoll


February 3, 2017

I’m going to wrap up the long-paused MATLAB versus Julia comparison on Trefethen & Bau by chugging through all the lectures on iterative methods in one post.

I’m back to using gists–not thrilled with any of the mechanisms for sharing this stuff.

These are remarkable mainly in that they have such striking similarity in both languages. Aside from square brackets and working around the 1x1/scalar distinction in Julia, little differs besides the syntax of the eigs command.

One frustration, though. I decided to try an interesting alternative to PyPlot in Julia, the Plots package. Actually Plots tries to be a generalization of and alternative route to using PyPlot/matplotlib. I decided to try the PlotlyJS backend instead, however. It makes lovely graphics with very responsive interaction. Since the rendering is in Javascript, I thought it would be perfectly portable, but you can’t see the output in the gist above, even though it should be embedded in the notebook.

I liked using Plots OK; for the most part it’s just different, not better or worse that I could see. I found it awkward to work with subplots. I ended up creating 4 plots individually and then displaying them in a table using another call to plot. I find MATLAB’s setup more convenient. I also could not figure out how to coax a contour plot with a contour at a specified value, which seems like a big lack.

Again the differences are minor. In sparse and iterative methods I found Julia to place a greater emphasis on keyword arguments. For example,

(xCG,~,~,~,resnorm) = cg(A,b,tol=1e-14,maxIter=100);

There are default values for tol and maxIter, but if you want to override them you must type the keyword. On the other hand, MATLAB’s arguments are purely positional:

[xCG,~,~,~,resnorm] = pcg(A,b,1e-14,100);

If I wanted to specify the maximum number of iterations without changing the default tolerance, then I would need to use an empty matrix in the third position. When one uses a command that does take named parameters as inputs, it’s typically done using 'propname',propval pairs. Except when it isn’t, such as for ODEs and optimization. Confusing! As a user I don’t love typing out the keywords, but Julia at least lets me skip the quote marks. I also know from experience that Julia’s version is a lot easier and clearer to implement on the other side.

So that’s that. I feel that I am at least ready to get off the bunny slopes with Julia. I haven’t found a compelling reason to switch to it, aside from supporting open source software for science (no small thing). Of course I’ve barely scratched the surface. On the flip side, MATLAB has a lot of well-designed and -maintained packages, and its environment still makes a smoother experience for newcomers. If you can afford it, it’s still a great option for interactive numerical computing.

I wonder about the future of Julia. Had Python not gotten a head start, I could see an outpouring of effort to make high-quality Julia packages and Julia being a complete MATLAB reboot. But numpy and scipy do exist, and despite their flaws, they have a huge first-mover advantage. It’s a snap to use Python packages in Julia, so there’s not a dichotomy here. But if the package you want to use a lot exists only in Python, the case for Julia weakens. Overall though, it’s a nice thing that we have several strong, expressive high-level environments for numerical computing. Happy coding!