Sollya User's Manual

Name:

implementpoly implements a polynomial using double, double-double and triple-double arithmetic and generates a Gappa proof

Library names:

sollya_obj_t sollya_lib_implementpoly(sollya_obj_t, sollya_obj_t, sollya_obj_t, sollya_obj_t, sollya_obj_t, sollya_obj_t, ...) sollya_obj_t sollya_lib_v_implementpoly(sollya_obj_t, sollya_obj_t, sollya_obj_t, sollya_obj_t, sollya_obj_t, sollya_obj_t, va_list)

Usage:

Description:

The command implementpoly implements the polynomial polynomial in range range as a function called functionname in C code using double, double-double and triple-double arithmetic in a way that the rounding error (estimated at its first order) is bounded by error bound. The produced code is output in a file named filename. The argument format indicates the double, double-double or triple-double format of the variable in which the polynomial varies, influencing also in the signature of the C function.

If a seventh or eighth argument proof filename is given and if this argument evaluates to a variable of type string, the command implementpoly will produce a Gappa proof that the rounding error is less than the given bound. This proof will be output in Gappa syntax in a file name proof filename.

The command implementpoly returns the polynomial that has been implemented. As the command implementpoly tries to adapt the precision needed in each evaluation step to its strict minimum and as it applies renormalization to double-double and triple-double precision coefficients to bring them to a round-to-nearest expansion form, the returned polynomial may differ from the polynomial polynomial. Nevertheless the difference will be small enough that the rounding error bound with regard to the polynomial polynomial (estimated at its first order) will be less than the given error bound.

If a seventh argument honor coefficient precisions is given and evaluates to a variable honorcoeffprec of type honorcoeffprec, implementpoly will honor the precision of the given polynomial polynomials. This means if a coefficient needs a double-double or a triple-double to be exactly stored, implementpoly will allocate appropriate space and use a double-double or triple-double operation even if the automatic (heuristic) determination implemented in command implementpoly indicates that the coefficient could be stored on less precision or, respectively, the operation could be performed with less precision. The use of honorcoeffprec has advantages and disadvantages. If the polynomial polynomial given has not been determined by a process considering directly polynomials with floating-point coefficients, honorcoeffprec should not be indicated. The implementpoly command can then determine the needed precision using the same error estimation as used for the determination of the precisions of the operations. Generally, the coefficients will get rounded to double, double-double and triple-double precision in a way that minimizes their number and respects the rounding error bound error bound. Indicating honorcoeffprec may in this case short-circuit most precision estimations leading to sub-optimal code. On the other hand, if the polynomial polynomial has been determined with floating-point precisions in mind, honorcoeffprec should be indicated because such polynomials often are very sensitive in terms of error propagation with regard to their coefficients' values. Indicating honorcoeffprec prevents the implementpoly command from rounding the coefficients and altering by many orders of magnitude the approximation error of the polynomial with regard to the function it approximates.

The implementer behind the implementpoly command makes some assumptions on its input and verifies them. If some assumption cannot be verified, the implementation will not succeed and implementpoly will evaluate to a variable error of type error. The same behaviour is observed if some file is not writable or some other side-effect fails, e.g. if the implementer runs out of memory.

As error estimation is performed only on the first order, the code produced by the implementpoly command should be considered valid iff a Gappa proof has been produced and successfully run in Gappa.

Example 1:

   > implementpoly(1 - 1/6 * x^2 + 1/120 * x^4, [-1b-10;1b-10], 1b-30, D, "p","implementation.c");
   1 + x^2 * (-0.166666666666666657414808128123695496469736099243164 + x^2 * 8.3333333333333332176851016015461937058717012405395e-3)
   > readfile("implementation.c");
   /*
       This code was generated using non-trivial code generation commands of
       the Sollya software program.

       Before using, modifying and/or integrating this code into other
       software, review the copyright and license status of this generated
       code. In particular, see the exception below.

       Sollya is

       Copyright 2006-2013 by

       Laboratoire de l'Informatique du Parallelisme, UMR CNRS - ENS Lyon -
       UCB Lyon 1 - INRIA 5668,

       Laboratoire d'Informatique de Paris 6, equipe PEQUAN, UPMC Universite
       Paris 06 - CNRS - UMR 7606 - LIP6, Paris, France

       and by

       Centre de recherche INRIA Sophia-Antipolis Mediterranee, equipe APICS,
       Sophia Antipolis, France.

       Contributors Ch. Lauter, S. Chevillard, M. Joldes

       christoph.lauter@ens-lyon.org
       sylvain.chevillard@ens-lyon.org
       joldes@lass.fr

       The Sollya software is a computer program whose purpose is to provide
       an environment for safe floating-point code development. It is
       particularily targeted to the automatized implementation of
       mathematical floating-point libraries (libm). Amongst other features,
       it offers a certified infinity norm, an automatic polynomial
       implementer and a fast Remez algorithm.

       The Sollya software is governed by the CeCILL-C license under French
       law and abiding by the rules of distribution of free software.  You
       can use, modify and/ or redistribute the software under the terms of
       the CeCILL-C license as circulated by CEA, CNRS and INRIA at the
       following URL "http://www.cecill.info".

       As a counterpart to the access to the source code and rights to copy,
       modify and redistribute granted by the license, users are provided
       only with a limited warranty and the software's author, the holder of
       the economic rights, and the successive licensors have only limited
       liability.

       In this respect, the user's attention is drawn to the risks associated
       with loading, using, modifying and/or developing or reproducing the
       software by the user in light of its specific status of free software,
       that may mean that it is complicated to manipulate, and that also
       therefore means that it is reserved for developers and experienced
       professionals having in-depth computer knowledge. Users are therefore
       encouraged to load and test the software's suitability as regards
       their requirements in conditions enabling the security of their
       systems and/or data to be ensured and, more generally, to use and
       operate it in the same conditions as regards security.

       The fact that you are presently reading this means that you have had
       knowledge of the CeCILL-C license and that you accept its terms.

       The Sollya program is distributed WITHOUT ANY WARRANTY; without even
       the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
       PURPOSE.

       This generated program is distributed WITHOUT ANY WARRANTY; without
       even the implied warranty of MERCHANTABILITY or FITNESS FOR A
       PARTICULAR PURPOSE.

       As a special exception, you may create a larger work that contains
       part or all of this software generated using Sollya and distribute
       that work under terms of your choice, so long as that work isn't
       itself a numerical code generator using the skeleton of this code or a
       modified version thereof as a code skeleton.  Alternatively, if you
       modify or redistribute this generated code itself, or its skeleton,
       you may (at your option) remove this special exception, which will
       cause this generated code and its skeleton and the resulting Sollya
       output files to be licensed under the CeCILL-C licence without this
       special exception.

       This special exception was added by the Sollya copyright holders in
       version 4.1 of Sollya.

   */


   #define p_coeff_0h 1.00000000000000000000000000000000000000000000000000000000000000000000000000000000e+00
   #define p_coeff_2h -1.66666666666666657414808128123695496469736099243164062500000000000000000000000000e-01
   #define p_coeff_4h 8.33333333333333321768510160154619370587170124053955078125000000000000000000000000e-03


   void p(double *p_resh, double x) {
   double p_x_0_pow2h;


   p_x_0_pow2h = x * x;


   double p_t_1_0h;
   double p_t_2_0h;
   double p_t_3_0h;
   double p_t_4_0h;
   double p_t_5_0h;



   p_t_1_0h = p_coeff_4h;
   p_t_2_0h = p_t_1_0h * p_x_0_pow2h;
   p_t_3_0h = p_coeff_2h + p_t_2_0h;
   p_t_4_0h = p_t_3_0h * p_x_0_pow2h;
   p_t_5_0h = p_coeff_0h + p_t_4_0h;
   *p_resh = p_t_5_0h;


   }

Example 2:

> implementpoly(1 - 1/6 * x^2 + 1/120 * x^4, [-1b-10;1b-10], 1b-30, D, "p","implementation.c","implementation.gappa");
1 + x^2 * (-0.166666666666666657414808128123695496469736099243164 + x^2 * 8.3333333333333332176851016015461937058717012405395e-3)

Example 3:

   > verbosity = 1!;
   > q = implementpoly(1 - dirtysimplify(TD(1/6)) * x^2,[-1b-10;1b-10],1b-60,DD,"p","implementation.c");
   Warning: at least one of the coefficients of the given polynomial has been rounded in a way
   that the target precision can be achieved at lower cost. Nevertheless, the implemented polynomial
   is different from the given one.
   > printexpansion(q);
   0x3ff0000000000000 + x^2 * 0xbfc5555555555555
   > r = implementpoly(1 - dirtysimplify(TD(1/6)) * x^2,[-1b-10;1b-10],1b-60,DD,"p","implementation.c",honorcoeffprec);
   Warning: the infered precision of the 2th coefficient of the polynomial is greater than
   the necessary precision computed for this step. This may make the automatic determination
   of precisions useless.
   > printexpansion(r);
   0x3ff0000000000000 + x^2 * (0xbfc5555555555555 + 0xbc65555555555555 + 0xb905555555555555)

Example 4:

   > p = 0x3ff0000000000000 + x * (0x3ff0000000000000 + x * (0x3fe0000000000000 + x * (0x3fc5555555555559 + x * (0x3fa55555555555bd + x * (0x3f811111111106e2 + x * (0x3f56c16c16bf5eb7 + x * (0x3f2a01a01a292dcd + x * (0x3efa01a0218a016a + x * (0x3ec71de360331aad + x * (0x3e927e42e3823bf3 + x * (0x3e5ae6b2710c2c9a + x * (0x3e2203730c0a7c1d + x * 0x3de5da557e0781df))))))))))));
   > q = implementpoly(p,[-1/2;1/2],1b-60,D,"p","implementation.c",honorcoeffprec,"implementation.gappa");
   > if (q != p) then print("During implementation, rounding has happened.") else print("Polynomial implemented as given.");
   Polynomial implemented as given.

Go back to the list of commands