My name is Reid Atcheson. I have a PhD in Computational and Applied Mathematics (CAAM) from Rice University. I currently
work for Numerical Algorithms Group (NAG). I periodically write here about my interests. My interests are mainly

- Scientific computation and numerical analysis
- High performance computing
- Functional programming

One of my recent tasks working at NAG was to establish a collaboration with Julia Computing. We decided to write a Julia interface to NAG’s nearest correlation matrix routines and demonstrate their usage. The post may be found either on NAG’s blog or on Julia Computing’s.

Read full post.

Not surprisingly from its name, Linear Algebra gives tools to analyze functions which
are linear. *Computational* linear algebra however customarily requires a
matrix representation of the linear functions under study, usually with floating point number entries.
A number of algorithms explicitly require this, such as LU and Cholesky factorizations which are important
tools for solving systems of linear equations. This
somewhat breaks away from the “interface” of linearity however, where we are no longer
only using the fact that the function is linear and are instead using only a *particular*
kind of linear function, namely a matrix. Algorithms do exist which only use the underlying
vector space structure and linearity of the function on that vector space. These are often
called “matrix-free,” these are such algorithms as Richardson iteration, GMRES, BiCGSTAB,
Arnoldi iteration, and many more; These algorithms are especially constructed to work
without knowledge of a matrix representation. This results in a gap between theory and
computation. On the theoretical side, usually little more information than linearity
is required to prove interesting and useful results about a linear function. There is
a very powerful spectral theory for bounded operators, self adjoint operators, and
normal operators for example. Yet when entering computation more information is demanded.
I decided to try and bridge this gap using a progrmaming language called OCaml. I want
to see if by following arguments made on the theoretical side I could reproduce some
elements of spectral theory for linear operators in a way that is also computational -
i.e. gives us an algorithm for computing these quantities. This will
open the algorithms of computational linear algebra to linear functions which can not
be easily represented by a matrix. I have been moderately successful already, so I will
document the early results here.

Read full post.

I just recently sorted out how to run MPI jobs on google compute cloud instances. For my own memory’s sake and also in the remote possibility that this could be useful to anyone I give an overview of the steps I took.

Read full post.

I have recently encountered several significant bugs in C codes with a common source: variable length arrays (VLAs). The bug ultimately arises from assumptions made years ago about how big a VLA could get, which would later be violated - leading to undefined behavior in C. The tricky part about this kind of bug is that there is no way within C to check for successful allocation of a VLA; incorrect VLA construction simply is undefined behavior leading usually to crashing code.

Read full post.

I came across an International Mathematical Olympiad problem from 1978 which I had the feeling I would be able to solve quickly. The problem is stated as follows:

Read full post.

I have successfully defended my PhD thesis “Accelerated Plane-wave Discontinuous Galerkin Methods for Heterogeneous Scattering Problems,” and today submitted my committee’s signatures. Hopefully now I will have a little more time to pursue some personal projects :)

Read full post.

In a previous post I gave a nice proof that there are infinitely many prime numbers. The proof however is somewhat nonconstructive because by way of a contradiction argument I sneak around issues of actually constructing a new prime not in a prespecified finite list, which would be the “constructive” manner of proving this same fact.

Read full post.

Here is the YouTube video of my JuliaCon presentation on the implementation of finite element methods in Julia. I had a great time at the conference and was really impressed with all of the presentations. I’m really glad so many people have come together to make this a reality.

Read full post.

In a recent post I showed that there are
infinitely many prime numbers. A little stronger fact than this is that the sum
of the reciprocals of prime numbers diverges, which I attempted to concisely explain in
this Quora response.
I have always wanted to find a simple proof of this fact, but so far I have not been successful. Euler’s proof to date
appears to me to be the most straightforward argument. However at the time of writing this post,
this Wikipedia entry
lists Euler’s proof and calls its steps *a sequence of audacious leaps of logic.*

Read full post.

*(Note: I have updated this post to reflect that the origin
of this argument was Euler, although it was unknown to me
when I first posted this)*

Read full post.

I decided to try my hand at various technologies of the web. Whether I have been a failure I will let you decide. The experience has however left me fairly confident that I should keep the majority of my focus on mathematics and scientific computation.

Read full post.