Re: [boost] Math Constants Formal Review

From: Terje Slettebø (tslettebo_at_[hidden])
Date: 2003-06-21 20:32:52

I've finally managed to catch up with mailing list postings, and have now
read through the review postings of this library.

First, I'd like to thank Paul for the work he's done with the library,
having gone through no less two formal reviews (a dubious pleasure that
Robert Ramey may experience with his serialisation library, as well).

I've read through the docs, and also tested the library. I'd first like to
do a mini-recap of Paul's recap of the options available in the library
(present in the library in order to be able to discuss them in the review).
There are several ways the constants are represented:

1. Macros - #define BOOST_PI 3.14
2. Constant variables (for float, long and long double, a file for each) -
const type pi=3.14
3. Functions - type pi() { return 3.14; }
4. Kenniston's approach - an object with implicit conversions to any type
it's defined for. More of that below.
5. Interval constants. I won't go into these here. These appear independent
of the above constants, and may be implemented using the same ways.

Then we have Daniel's approach, as described in a posting. Also more of that
below.

Kenniston's approach
=================
It works as follows (quoting from the library):

real_cast.hpp
-------------------

template< typename Tag >
struct real_type
{
  template< typename Target >
  inline operator Target() const { return constant_value( static_cast<
Target * >( 0 ), *this ); }
};

pi_constant.hpp
----------------------

#include <real_cast.hpp>

struct pi_tag;

typedef real_type< pi_tag > pi_t;

extern pi_t pi; // Defined in pi_constants.cpp

inline float constant_value( const float *, const pi_t & )
{
  return 3.1F; // Simulating VERY low precision pi as float.
}

// Same for double and long double

That's it.

There are some components in the Kenniston version that appear to be not
needed, at least if you use other appropriate ways, such as real_cast:

template< typename Target, typename realType >
inline Target real_cast( const realType & real, const Target * target = 0 )
{
  return constant_value( target, real );
}

This is used in a couple of places, such as a few convenience overloaded
operators, like:

  template< typename Target, typename Tag >
   inline Target operator+( const Target & lhs, const real_type< Tag > &
rhs )
  { // For example: 2.F + pi
   return lhs + real_cast< Target >( rhs );
  }

However, at least when testing the library on Intel C++ 7.1, these operators
weren't needed. The reason is that the templated conversion operator of
real_type automatically performs the appropriate conversions, according to
the arguments used, just as if ordinary variables and constants are used in
expressions.

Thus these operators doesn't appear to be needed. Another place real_cast is
used is when the constants are referred to directly, without other
operators, so you may need to provide the type you want the constant value
as. The library uses:

cout << boost::math::real_cast<float>(pi) << endl;

However, a just as good way, which may well be clearer, as it requires no
new casts, is:

cout << (double) boost::math::pi << endl;

Note that, as has been pointed out in this thread, this doesn't perform any
conversion of the constant value - the cast is just used to select the
appropriate conversion operator instantiation, just like real_cast, above.

Now, let's cover Daniel's approach.

Daniel's approach
==============
This is similar to Kenniston's, except that it uses specialisation, rather
than overloading.

    // Generic base class for all constants
    template< typename T, template< class > class F > struct constant
    {
        // A cast-to-anything-operator :)
        template< typename U > operator U() const { return F< U >()(); }

#define ADD_OPERATOR( OP ) \
        template< typename U > friend U operator OP( const T& lhs, const U&
rhs ) \
        { U nrv( static_cast< U >( lhs ) ); nrv OP##= rhs; return nrv; } \
        template< typename U > friend U operator OP( const U& lhs, const T&
rhs ) \
        { U nrv( lhs ); nrv OP##= static_cast< U >( rhs ); return nrv; }

        ADD_OPERATOR( + );
        ADD_OPERATOR( - );
        ADD_OPERATOR( * );
        ADD_OPERATOR( / );
#undef ADD_OPERATOR
    };

    // Here's the definition for pi for all types (can be extended by UDTs):
    template< typename T > struct pi_value;
    template<> struct pi_value< float > { float operator()() const { return
3.14; } };
    template<> struct pi_value< double > { double operator()() const {
return 3.1416; } };
    template<> struct pi_value< long double > { long double operator()()
const { return 3.1415927; } };

    // Here's the single line to create a useful interface
    struct pi_t : constant< pi_t, pi_value > {} pi;

That's all there is to it.

There has also been a discussion of using expression template techniques for
optimising computations (and Daniel's library shows an example of that). A
couple of observations on that:

- This is possible using either Kenniston's or Daniel's approach, or
something similar. It relies on just overloading the appropriate operators.

- For built-in types, this kind of optimisation may well have no effect, as
the compiler may be able to perform the constant folding, etc., itself

- For UDTs, it should be up to the providers of the UDT to provide
expression templates, if wanted, as the constant library can't provide this
in the general case.

Another point is that this approach uses template template parameters, so it
won't work on MSVC 6. However, there are known workarounds for that (such as
class with member template).

Evaluation
========

Macros
-----------
- These are the only one with guaranteed no overhead (space or time) for the
executable.

- They don't respect scope, such as namespace scopes.

- There's no obvious way of extending them for new types. Even the current
ones only define constants for long double, relying on conversions to other
types. E.g.:
  #define BOOST_PI 3.14L

  To get alternatives for other types, you need to define other names, such
as:

  #define BOOST_FLOAT_PI 3.14F

  However, this works badly in generic code, like:

  template<class T>
  T area(T radius)
  {
    return math::pi*radius*radius;
  }

Constant variables
--------------------------
- Also these work badly in generic code like the above. You can't "overload"
variables.

- They may be extended to new types, by simply defining constant variable
for new types in the same way. However, as they can't be used generically,
and two definitions of the same name (and scope) can't co-exist in the same
program, this is of limited value.

- If a header, e.g. for the float constants are included, you can't include
any other corresponding header for another type. As has pointed out in the
discussion, this gives the possibility of "silent conversions", for example
from double to float, in this case.

Functions
-------------
- These require the pi() syntax, and it appears that this is not well
received.

- They may be extended to new types by overloading the functions.

Kenniston's approach
------------------------------
- Together with Daniel's approach, this is the only approach (besides
functions) that works well with generic code.

- They may be extended to new types, using overloading of the
"constant_value()" function.

- There is no danger in including several incompatible headers (unlike the
case for "constant variables"), as each constant is defined only once.
Specialisations are used to provide values.

- There are no "silent conversions" resulting from including a wrong header.
The constants will provide the best match available in a given expression,
or give an ambiguity if there's no best match. Moreover, it will do this
conversion automatically, with no need for user-specification (unless
there's ambiguity). This way, they work the same way as the built-in types,
which is a big advantage.

Daniel's approach
-------------------------
- As mentioned, this is similar to Kenniston's, with the same advantages,
and uses specialisation of the "<constant name>_value" class template.

Either overloading or specialisation may be fine. An advantage with
specialisation, albeit probably academic in this case, is that overloading
is order-dependent, whereas specialisations are not (as long as the
specialisation exists in the program).

Having considered Kenniston's and Daniel's approach, I wonder if it's
possible to simplify it. The following uses the same approach as Daniel's -
specialisation, as it allows it to be done with just one class for each
constant, and nothing else:

--- Start ---

#include <iostream>

#define BOOST_DEFINE_MATH_CONSTANT(name)\
const struct name##_type\
{\
  name##_type() {}\
  template<class T>\
  operator T() const;\
} name;

#define BOOST_MATH_CONSTANT_VALUE(name, type, value)\
template<>\
name##_type::operator type() const { return value; }

namespace boost {
namespace math {

BOOST_DEFINE_MATH_CONSTANT(pi)

BOOST_MATH_CONSTANT_VALUE(pi, float, 3.1)
BOOST_MATH_CONSTANT_VALUE(pi, double, 3.14)
BOOST_MATH_CONSTANT_VALUE(pi, long double, 3.141)

} // namespace math
} // namespace boost

using namespace boost;

template<class T>
T area(T radius)
{
  return math::pi*radius*radius;
}

int main()
{
  std::cout << "Area (float)=" << area(1.0f) << '\n'
            << "Area (double)=" << area(1.0) << '\n'
            << "Area (long double)=" <<  area(1.0l) << '\n'
            << "PI=" << (float) math::pi << ", " << (double) math::pi << '",
" << (long double) math::pi << '\n';
}

--- End ---

Output:

Area (float)=3.1
Area (double)=3.14
Area (long double)=3.141
PI=3.1, 3.14, 3.141

This is tested on Intel C++ 7.1, but it should work on most compilers,
especially if workarounds are applied, if needed.

As shown, this works well with generic code, like the "area()" function.

My question is: Is there any reason why it can't be done this simple?

Note that there are no extra classes or functions needed. There's only the
class representing the constant, and the templated conversion operator, and
any specialisations of it. This should be trivial for a compiler to inline,
and should also be easy to examine with a debugger (an issue that also has
come up).

This is very easy to extend with new constants and types. It's an open-ended
system.

It uses macros for easy definition of constants and their values for
different types, for easy definition (just like Kenniston's and Daniel's
approach).

The above addresses the presentation of the constants (their interface), as
identified by Paul as one of the primary objectives. What values that are
used in the definitions is orthogonal to this.

Regards,

Terje