Ada binding to GMP and MPFR

["vincent.diemunsch@gmail.com" <vincent.diemunsch@gmail.com>]

Hello everybody,

As the libraries GMP ( http://gmplib.org/ ) and MPFR ( http://www.mpfr.org/
) becomes part of the GNAT
free compiler (since they are part of the new GCC) and since these
libraries have excellent performances,
I thought it could be interesting to create an Ada binding for them.

I searched on Internet and found two old binding but non convinced me
and they seemed to have been droped
out. Therefore, I have undertaken to write my own binding for GMP and
MPFR in Ada 2005, and I have just managed to successfuly test basic
operations for them such, as adding or multiplying MPFR_Floats of
different precision...

Now I wonder if someone could be interested in hosting the sources on
an Internet page, so that anybody could
download them and test them and eventually make improvements...

The Binding is as follows :

* A THIN BINDING composed of two files :
- gmp.ads
- mpfr.ads

It basically translates in Ada most functions of gmp.h and mpfr.h but
it's not exhaustive.
There are also some specific files for target dependent types :
- mp_x86_32bits.ads
- mp_x86_64bits.ads...

* A THICK BINDING with the following specification files :
- gmp.Integers.ads;
- gmp.Rationals;ads;
- mpfr.Floats.ads;

These files declare the following types : Unbounded_Integer,
Unbounded_Fraction, MPFR_Float
that can be seen as extensions of the typical Ada Integer and Float,
with operator overloading
"+", "*", ... and elementary functions.


I know that some work needs to be done to :
- complete the gmp.ads and mpfr.ads files
- fully test the binding
- have a good integration in the Ada "spirit".

Therfore, I require your help :
- Where may i found a place to put the files (or a SVN repository ;-)
- Would you be interested in giving comments on implementation
especially for MPFR_Types ?

Thanks a lot for your help.
Regards,

Vincent

EXTRACT OF MPFR-FLOAT.ADS :

-- Copyright (c) 2009 Vincent DIEMUNSCH.
-- GNU General Public License.
--
-- This file is distributed in the hope that it will be useful, but
WITHOUT ANY
-- WARRANTY; without even the implied warranty of MERCHANTABILITY or
FITNESS FOR
-- A PARTICULAR PURPOSE.  See the license for more details.
--
--------------------------------------------------------------------------------
-- This package is a thick Ada binding to the Multiple Precision
Floating Point
-- Reliable Library. See : http://www.mpfr.org/
--
-- It declares a new private type : MPFR_Float, that is an abstraction
of a
-- floating point value of arbitrary precision.
--
-- This packages uses a new approach to dynamically adjust the
precision of a
-- MPFR_Float, everytime a value is computed (and not, as mpfr does,
at the
-- initialization of the mpfr_t variable). Before any operation,
-- the precision required for the result is set according to the
following rules:
--    1) if the result is exact, then precision is the default
precision
--       for instance : , "1/3", "pi/4", "10.5", "1.45e-3".
--    2) for a unary operation, the precision of the result is the
precision
--       of the operand.
--       for instance : "abs(x)", "1/x", "x+n" where n is an integer.
--    3) for a binary multipliying operation the precision of the
result is the
--       smallest precision of the two operands.
--       for instance : "x*y", "x mod y", etc.
--    4) for a binary adding operation, the precision of the result is
the
--       smallest precision of :
--        * the greatest operand in absolute value,
--        * the smallest operand in absolute value, adjusted to the
difference
--          of exponants.
--       for instance : "1.111111 + 1.0e-5" (resp. a precision of 7
digits and
--                       2 digits) will give "1.111121" (7 digits = 2
- (0-5))
--                       and not "1.1" (2 digits).
--
-- If, an adjustement of precision is required, in a computation, for
instance
-- when the result is computed through the use of a convergent
sequence, one
-- can use the functions "Unchecked_Pecision_Adjustement", with 2
variants :
--   1) an absolute precision expressed in digits
--   2) a relative precision, throug the comparison to another value :
--      new precision = exponant(x) - exponant(x - value) + n digits
(n > 3).
--      for instance :
--          u(n) = 1/u(n-1) + 1   (Computing the golden number)
--          u(n) = 1.615 (4 digits) and u(n-1) = 1.625, error : 1e-2,
--          then the precision will be set to 0 - (-2) + n digits, at
least
--          5 digits ("1.6150")


with Ada.Finalization;
with GMP.Integers;

package MPFR.Floats is

   type MPFR_Float is private;

   -- Private type so that all the following functions are NOT
primitive
   -- operations.

   -- PRECISION :
--------------------------------------------------------------

   procedure Set_Default_Precision (Decimal_Digits : Positive);
   -- Set the default precision whenever it is used.
   -- Note : if you want to precisely set this precision in bits
simply use the
   --        procedure mpfr_set_default_prec (prec : mpfr_prec_t) from
the
   --        MPFR package.

   procedure Unchecked_Precision_Adjustement
     (x              : in out MPFR_Float;
      Decimal_Digits : Positive);

   procedure Unchecked_Precision_Adjustement
     (x              : in out MPFR_Float;
      Error          : in     MPFR_Float;
      Decimal_Digits : Positive := 3);

   -- Return the value of x, either rounded or extended with trailing
zeros.
   -- Adjust the precision either to the absolute precision given in
digits, or
   -- by adjusting it to the given digits plus the order of magnitude
of the
   -- the difference between x and the reference value.

   type Rounding_Mode_Kind is (Round_to_Nearest,
                               Round_Toward_Zero,
                               Round_Toward_Plus_Infinity,
                               Round_Toward_Minus_Infinity);
   -- we declare here exactly the same type as the mpfr_rnd_t from the
   -- MPFR package, but for the sake of clarity and simplicity of the
user.

   procedure Set_Default_Rounding_Mode (Rounding :
Rounding_Mode_Kind);
   -- Set the default rounding used by this package, every time an
operation
   -- is done on a MPFR_Float. The default rounding mode is
   -- to nearest initially.

   procedure Set_Precision_Check (Decimal_Digits : Positive);
   procedure Inhibate_Precision_Check;
   -- These functions enable or inhibate the checking during
assignement to a
   -- variable, that the value has a precision above the curret
minimal required
   -- precision.
   -- If needed, the Not_Accurate exception is raised.


   -- INITIALISATION :
---------------------------------------------------------
   function To_MPFR_Float (Image          : String;
                           Decimal_Digits : Natural := 0 )
                           return MPFR_Float;

   -- This function creates Floating Point from a value given by a
string
   -- that must conform to the Ada Standard for Float literals
(decimal literals
   -- and base literals), and whose precision is at least the number
of
   -- decimal_digits given, or if this value is nul the
Default_Precision.
   -- The base can be above 16, but it must be under (or equal to) 36
(And in
   -- that case we use digits from 0 .. 9, A .. Z).
   -- We tolerate a negative value (i.e. the '-' at the beginning),
even if in
   -- Ada the '-' is seen as a unary operator and not part of the
literal.
   -- Exemples :
   --      Decimal literals : 12.0   0.0   0.456  3.14159_26
   --      Based   literals : 16#F.FF#E+2   2#1.1111_1111_1110#E11  =
4095.0

   Data_Error : Exception;

   subtype Base_Range is Integer range 2 .. 36;

   function Image (Item : MPFR_Float;
                   Base : Base_Range := 10)   return String;
   -- This function return a Float literal (decimal literal or base
literal)
   -- that conform to the Ada Standard except the '-' at the beginning
and the
   -- fact that base can be any integer between 2 and 36.


   function Not_A_Number   return MPFR_Float;
   function Plus_Infinity  return MPFR_Float;
   function Minus_Infinity return MPFR_Float;

   function Pi return MPFR_Float;
   -- Pi at the current default precision.

   ----------------------------
   -- "Standard" description --
   ----------------------------

   -- (Mimics Standard package)
   -- The predefined operators for this type are as follows:

   function "="   (Left, Right : MPFR_Float) return Boolean;
   function "<"   (Left, Right : MPFR_Float) return Boolean;
   function "<="  (Left, Right : MPFR_Float) return Boolean;
   function ">"   (Left, Right : MPFR_Float) return Boolean;
   function ">="  (Left, Right : MPFR_Float) return Boolean;
   pragma Inline ("=", "<", "<=", ">", ">=");


   function "abs" (Right : MPFR_Float) return MPFR_Float;
   function "-"   (Right : MPFR_Float) return MPFR_Float;
   pragma Inline ("abs", "-");

   function "+"   (Left, Right : MPFR_Float) return MPFR_Float;
   function "-"   (Left, Right : MPFR_Float) return MPFR_Float;
   function "*"   (Left, Right : MPFR_Float) return MPFR_Float;
   function "/"   (Left, Right : MPFR_Float) return MPFR_Float;

   function "**"  (Left : MPFR_Float; Right : Integer'Base) return
MPFR_Float;

   -- The specification of each operator for the type root_real, or
for
   -- any additional predefined floating point type, is obtained by
   -- replacing Float by the name of the type in the specification of
the
   -- corresponding operator of the type Float.

   -- In addition, the following operators are predefined for the root
   -- numeric types:
   function "*" (Left : MPFR_Float;   Right : Integer'Base) return
MPFR_Float;
   function "*" (Left : Integer'Base; Right : MPFR_Float)   return
MPFR_Float;
   function "/" (Left : MPFR_Float;   Right : Integer'Base) return
MPFR_Float;

   -- and we add :
   function "+" (Left : MPFR_Float;   Right : Integer'Base) return
MPFR_Float;
   function "+" (Left : Integer'Base; Right : MPFR_Float)   return
MPFR_Float;
   function "-" (Left : MPFR_Float;   Right : Integer'Base) return
MPFR_Float;
   function "-" (Left : Integer'Base; Right : MPFR_Float)   return
MPFR_Float;
   function "/" (Left : Integer'Base; Right : MPFR_Float)   return
MPFR_Float;

   --------------------------
   -- with Standard Floats --
   --------------------------
   -- use Float'Image and Float'Digits to create the MPFR_Float

   -----------------------------
   -- with Unbounded_Integers --
   -----------------------------

   subtype Unbounded_Integer is GMP.Integers.Unbounded_Integer;
   function To_MPFR_Float (Value : Unbounded_Integer) return
MPFR_Float;

   function Round_To_Unbounded_Integer(Item     : MPFR_Float)
                                       return Unbounded_Integer;
   function Round_To_Unbounded_Integer(Item     : MPFR_Float;
                                       Rounding : Rounding_Mode_Kind)
                                       return Unbounded_Integer;
   -- Round to an unbounded integer, using either the default rounding
mode
   -- or the specified mode.


   --------------------------
   -- Elementary functions --
   --------------------------

   function Sqrt (X    : MPFR_Float)        return MPFR_Float;
   function Log  (X    : MPFR_Float)        return MPFR_Float;
   function Log  (X    : MPFR_Float;
                  Base : Integer'Base)      return MPFR_Float;
   function Exp  (X    : MPFR_Float)        return MPFR_Float;
   function "**" (Left : MPFR_Float;
                  Right: Unbounded_Integer) return MPFR_Float;
   function "**" (Left, Right : MPFR_Float) return MPFR_Float;

   function Sin    (X : MPFR_Float) return MPFR_Float;
   function Cos    (X : MPFR_Float) return MPFR_Float;
   function Tan    (X : MPFR_Float) return MPFR_Float;
   function Sec    (X : MPFR_Float) return MPFR_Float;
   function Csc    (X : MPFR_Float) return MPFR_Float;
   function Cot    (X : MPFR_Float) return MPFR_Float;

   function Arcsin (X : MPFR_Float) return MPFR_Float;
   function Arccos (X : MPFR_Float) return MPFR_Float;
   function Arctan (X : MPFR_Float) return MPFR_Float;

   function Sinh   (X : MPFR_Float) return MPFR_Float;
   function Cosh   (X : MPFR_Float) return MPFR_Float;
   function Tanh   (X : MPFR_Float) return MPFR_Float;
   function Sech   (X : MPFR_Float) return MPFR_Float;
   function Csch   (X : MPFR_Float) return MPFR_Float;
   function Coth   (X : MPFR_Float) return MPFR_Float;
   function Arcsinh(X : MPFR_Float) return MPFR_Float;
   function Arccosh(X : MPFR_Float) return MPFR_Float;
   function Arctanh(X : MPFR_Float) return MPFR_Float;

   function Factorial (N : Natural'Base) return MPFR_Float;
   function Gamma     (X : MPFR_Float)   return MPFR_Float;

   ------------------
   -- Exceptions : --
   ------------------

   MPFR_Underflow    : exception;
   -- occurs when the exact result of a function is a non-zero real
number and
   -- the result obtained after the rounding, assuming an unbounded
exponent
   -- range (for the rounding), has an exponent smaller than the
minimum exponent
   -- of the current range. In the round-to-nearest mode, the halfway
case is
   -- rounded toward zero.
   MPFR_Overflow     : exception;
   -- occurs when the exact result of a function is a non-zero real
number and
   -- the result obtained after the rounding, assuming an unbounded
exponent
   -- range (for the rounding), has an exponent larger than the
maximum exponent
   -- of the current range. In the round-to-nearest mode, the result
is infinite.
   MPFR_Not_a_Number : exception;
   -- occurs when the result of a function is a NaN.
   MPFR_Inexact      : exception;
   -- occurs when the result of a function cannot be represented
exactly and must
   -- be rounded.
   MPFR_Range_Error  : exception;
   -- occurs when a function that does not return a MPFR number (such
as
   -- comparisons and conversions to an integer) has an invalid
result.

   Not_Accurate : exception;
   -- This exception might be raised by the "Precision_Check".
   -- It is not a MPFR exception in fact !


private

   type MPFR_Float is new Ada.Finalization.Controlled with record
      Value    : MPFR.mpfr_t;
   end record;

   overriding procedure Initialize (Object: in out MPFR_Float);
   overriding procedure Adjust     (Object: in out MPFR_Float);
   overriding procedure Finalize   (Object: in out MPFR_Float);

end MPFR.Floats;