C++ wrapper

[Gerard Jungman]

Dave Morrison wrote:
> http://www.fnal.gov/docs/working-groups/fpcltf/Pkg/SpecialFunctions/doc/html/0SpecialFunctions.html

Hmmm. Ok. This is not meant as a criticism (really, it's not),
but how come this is the first time I have heard about this?
I know a little about ZOOM, but I didn't know they were
wrapping GSL special functions (other stuff too...?).

I'm just interested in making sure that we keep open lines
of communication. I've been thinking alot about native
C++ implementations for some of the GSL functionality,
so I am interested in what people want to do. Anyway,
I see that this is a fairly straightforward wrapper
project, so maybe they thought it would not be
sufficiently whiz-bang to bring up in this forum.
But still, any sort of general discussion about
C++ is welcome here.

I am especially interested in new object-oriented
designs. As an example, consider the following.
You can get these situations where you want to calculate
a series of values, or you want both the function and its
derivative (say for the functions associated to second-order
equations, Bessel being the canonical example). Then there
are often ways to calculate the necessary extra values 
more optimally than just doing consecutive calculations
(this is clear for the case of a sequence of values
related by recursion, for example). In the procedural
world you deal with this by providing multiple interfaces,
such as
    o calculate f(x)
    o calculate both f(x) and f'(x)
    o calculate f_n(x), f_{n+1}(x), ... f_m(x)

The multiplicity of interfaces here reminds one
of the kind of multiplicity which occurs in the
standard procedural interfaces for things like
linear algebra (lapack or blas, with their exponential
explosion of interface for all the different types and
algorithms). It is not nearly so bad here, but we are still
talking about an increase in the size of the interface
specification by a factor of 3 or so, over the "naive"
notion of "Bessel function". That's annoying for various
reasons.

  o user perspective: more functions to understand and
    wade through when looking for what you want

  o advanced user perspective: may want to reuse the
    library algorithms in his own setting, but there
    is no way to look into a library procedure, other than
    a disagreeable copy-paste-edit-debug cycle

  o implementor perspective: more isolated but closely
    related functionalities to support, therefore more
    opportunity for problems like functional duplication

The abstract problem is that these more involved
calculations and composite calculations have
state, and the basic design problem is encapsulating
that state. In the procedural design you encapsulate
that state solely in the procedure local variables.
I would like to see an object-oriented design where
that state gets stored in the obvious place: some
object designed for it.

Note that there is an implicit assumption here
about layering the design, so that average end
users are not exposed to complexity that they
are not interested in. However, the layer will
be there to support users with special needs.
In particular, the user with the most demanding
special needs is the implementor himself who
would ideally be able to compose pieces of the
implementation from well-encapsulated and orthogonal
base parts.


This may seem like alot of discussion for something
as simple as "calculate J0(x)". But when you
start thinking about efficient and clean
implementations, you realize it is not so
simple to get what you want out of the obvious
mapping from functions of mathematics to functions
in a procedural language.

I have written codes which were naively dominated
by function evaluation and had to be structured
in non-naive ways in order to make use of properties
like recurrences to get the optimal performance
(we're talking factors of 10 here, not fiddling
at the 20% level). That sort of restructuring is
often not the way one thinks about the problem
from a mathematical point of view, and so you
get an extra layer of confusion. I wonder if it
would be possible to express what you mean in
these sorts of situations rather than how the
computer must do it.

Furthermore, I don't see anybody talking about this
at all. Most of the work I see in high-performance C++
is for vector/matrix/linear-algebra functionality.
Very important stuff. But the question of how to create
efficient clean and usable high-performance components
in other, perhaps more classical (and therefore boring?),
areas of numerical computing is not often addressed.

Here is a basic question: Is this going to be harder
or easier than the work being done on vector/matrix
libraries? I don't know. Maybe some new things
will have to be invented. Certainly we get new
abstractions to play with, along with their
associated implementation problems.


-- 
G. Jungman