Recurrence relations is the mathematical form used to express many important problems. Being able to rapidly and efficiently prototype a recurrence will aid research and engineering in diverse areas as heat flow, seismic processing, dynamic programming, etc. Unfortunately, the conventional techniques available is either to program the recurrence as a recursive program with exponential growth in time complexity, or unfold the recursion and hand code an efficient solution, with the associated high development costs. Here we present an approach that combines the best of these two methods.
Computational modelling using partial differential equations (PDEs) is becoming an important tool for industry. Currently this is hampered by the time and effort needed to develop good computational models, and the time and cost needed to port such a model onto a high performance computer (HPC).
Recurrence relations form the core of a computational model of a PDE. For many problems, such as those in computational fluid dynamics, formulating an appropriate recurrence relation is a non-trivial task. Often it will require several iterations of formulating a recurrence, creating a prototype implementation and trying it on selected data sets, before a usable strategy is found. The test problems may require the use of an HPC to increase the turn-around time during the computational model development cycle. Being able to rapidly prototype reasonably efficient implementation of recurrence relations will significantly reduce the time and effort needed to find the appropriate recurrence.
An order k linearly dependent recurrence r with natural numbers as index domain is given by the form
F(n) = if n > 1 then F(n-1) + F(n-2) else 1
Conventional recursive execution
Optimal execution strategy
The natural way of expressing a recurrence relation in a programming
language such as C++ would be as recursive function calls mimicing the
recursive equations from the mathematical formulation. Unfortunately, this
results in a program growing exponentially both in storage and time usage,
while a non-recursive implementation may need only a fixed amount of
storage and linear time (as in the case of the Fibonacci relation), or
linear growth of storage and n log n or n (log
n)^2 growth of time as in the FFT and bitonic sorting cases,
respectively. Thus naively writing a recurrence as a recursive function is
not a suitable path for prototyping recurrences.
The other alternative, hand coding a non-recursive version, gives efficient code, but is a difficult, hence costly, task if the stencil is complex. Thus this is not a viable approach for prototyping recurrences.
Our approach is to use the recursive formulation of the recurrence, but augment the definition of the stencil D, with the definition of the opposite of the stencil, Ð. This allows us to compile the recursively structured recurrence function into non-recursive code, giving the same time and space complexity as one would expect from a reasonably efficient implementation. One can also augment D and Ð with information on the distribution of the computation on the processors of a parallel computer. This information will also define the communication paths, thus yielding full control over the parallel execution of the code.
To improve the usability of this approach to prototyping recurrences, we will also develop a library of DDA's representing different standard stencils. For a given stencil, the library will also include instructions for the parallelisation of the recurrences using it, thus allowing a prototype implementation of a recurrence to be run efficiently on several HPC architectures without requiring the user to rewrite the recurrence program.