John Maddock wrote:
For me the rational data type is very important and should wrap the gmp mpq_type. It would be nice if the library could use mp_int with Boost.Rational to implement its own standalone mp_rational datatype, but I would prefer to use the gmp type whenever possible.
Support for mpq_t is on my (long) todo list.
I'm glad to hear. In theory a fixed precision floating point type with 128 bits or more in the mantissa would also satisfy my needs, but I have only tested with multi-precison rational.
By the way, I have put extensive thought into how to build a good expression template system for dealing with multi-precision arithmetic and would like to be involved in some of the discussion of the design of the system. Particularly, the most important thing about providing expression templates for multi-precsision is eliminating allocations and deallocations. The allocations can easily dominate the runtime of a multiprecision algorithm.
I suspect this may depend on the data type - for real-valued types, I've done some experimenting with a port of the LINPACK benchmark to C++ that shows *almost* no difference between expression template enabled code (whether mine, or existing wrappers such as mpf_class) and more traditional class abstractions. I need to do some experimenting and profiling, but it seems like the runtime is completely dominated by the cost of multiplication/division, and even though my code generates fewer temporaries that either mpf_class or traditional non-ET wrappers, that simply doesn't take much off the runtime. Of course LINPACK is an old Fortran program written in a very idiomatic style that probably doesn't really stretch the expression template code at all. So arguably we need better test cases - the special functions in Boost.Math would make good candidates I guess.
Gmp mpq_class expressions allocate for temporaries in compound expressions. If you declare your numerical types outside the loop and write only single operation expressions inside the loop with no implicit temporaries created then mpq_class won't allocate at all in the loop (except when the max size of the buffer in the numerical types needs to be occasionally increased.) The result can be allocation free execution in the common case. And I have seen significant speedups due to rewriting code to avoid allocations in both the line segment intersection and voronoi diagram use cases. I would like to get the same benefits without needing to rewrite arithmatic code so that replacing a builtin numerical datatype with a multiprecision in a templated context would produce optimal code even if the author followed the normal style of declaring variables in the inner most scope possible and large compound expressions with many temporaries.
There are many options and considerations for how to optimize the allocations, such as object pools, thread safety etc. Anyone who uses multiprecision numerical data types in the same way they would use built-in numerical data types is writing code that runs significantly slower than the equivalent C code that doesn't abstract away allocation and deallocation and forces the programmer to choose where and when to allocate so that they choose the proper place. I saw large speedups in my Polygon library by recycling gmpq_type variables rather than declaring them in the inner most scope in which they were used.
I take it you mean that within the "big number" class type initialized mpq_t variables get cached and recycled?
Yes. One simple way is the cache a temoporary along with the primary value in each wrapeed numerical object and return reference to the tempoarary as the result of an operator. That way if the numerical object stays in scope through many iterations of a loop the tempoarary gets recycled also. However, this still requires that numerical objects be declared at the outer rather than inner scope, which is not normal style.
That's an interesting idea, and not something I've investigated so far - not least because of the issues you highlight.
One thing I do want to investigate for fixed-precision real-valued types (with mpfr or mpf backends) is to eliminate the allocation altogether by placing the storage directly in the class (if it makes sense for not-too large precisions).
I like the idea of fixed sized array directly in the class. It would be a huge win and beat typical usage of gmp in C++ with even a slower arithmetic algorithm because allocation dominates arithmetic for modest sized values.
I have a couple of questions relating to the polygon library if that's OK?
* Is there a concept check program anywhere that I can plug my types into and verify they meet the libraries conceptual requirements?
Not exactly. I considered writing one, but basically what you need to do is make sure the traits for your type are defined and that it is registered as the right concept type then write a test that exercises the traits through the free function of the same name that expects the concept type. For example point concept has get and set free functions that exercise the get and set traits. If those two functions work correctly (and construct()) then your point is a model of the concept. If you copy-paste the example code for mapping a user defined type to the concept and replace my type with yours and compile the test is part of the example code because I exercise the traits.
* Is there a good program for testing "big number" performance within your library (for either large integers or rationals)?
If you look at the gmp_override.hpp you will see how I abstract away the gmpq_class as a high_precision_type<> metafunction lookup. You can replace gmpq_class with anything you want and benchmark. However, you will find that the performance of the high precision type is almost irrelevant to performance of the larger algorthm because lazy-exact arithmetic implies that it is very rarely if ever used. What you can do instead is dig into my polygon_aribtrary_formation.hpp in details and find the function that computes the intersection point of two line segments using the high precision type. You can then write a benchmark that generates random line segments and intersects them. Once you have that working you can swap in different high precision types by replacing gmp_override.hpp with your own override header file that you include just after polygon.hpp and before your own benchmark code. I have a "fitness" test for intersection points so you can also use the benchmark code as a stress test for the numerical data type chosen to see if it would provide robust execution of my algorithm. You will see the fitness test applied at the end of the lazy version of line segment intersection just before it conditionally calls the exact version.
With regard to expression templates, can I direct you to the "big_number" directory of the sandbox and the wildly out of date docs at http://svn.boost.org/svn/boost/sandbox/big_number/libs/math/doc/html/index.h... In particular the conceptual requirements required for backend types are constantly evolving - both as the needs of new backends are evaluated - and as I come up with new "good ideas" ;-) Nonetheless I would much welcome your input, and will try and add a rationale type backend ASAP for you to experiment with.
I'll try to look at is as time permits. Thanks, Luke