← Back to branch summary

~ubuntu-branches/ubuntu/precise/eglibc/precise-201308281639

« back to all changes in this revision

Viewing changes to sysdeps/ieee754/dbl-64/mpexp.c

  • Committer: Package Import Robot
  • Author(s): Matthias Klose
  • Date: 2012-02-08 01:58:09 UTC
  • mfrom: (1.5.3) (288.1.12 precise)
  • Revision ID: package-import@ubuntu.com-20120208015809-ulscst7uteq3e22z
Tags: 2.15~pre6-0ubuntu10
Merge from Debian (r5151, 2.13-26).

Show diffs side-by-side

added added

removed removed

Lines of Context:
sysdeps/ieee754/dbl-64/mpexp.c
1
 
 
2
1
/*
3
2
 * IBM Accurate Mathematical Library
4
3
 * written by International Business Machines Corp.
5
 
 * Copyright (C) 2001 Free Software Foundation
 
4
 * Copyright (C) 2001, 2011 Free Software Foundation
6
5
 *
7
6
 * This program is free software; you can redistribute it and/or modify
8
7
 * it under the terms of the GNU  Lesser General Public License as published by
34
33
#include "mpa.h"
35
34
#include "mpexp.h"
36
35
 
 
36
#ifndef SECTION
 
37
# define SECTION
 
38
#endif
 
39
 
37
40
/* Multi-Precision exponential function subroutine (for p >= 4,          */
38
41
/* 2**(-55) <= abs(x) <= 1024).                                          */
39
 
void __mpexp(mp_no *x, mp_no *y, int p) {
 
42
void
 
43
SECTION
 
44
__mpexp(mp_no *x, mp_no *y, int p) {
40
45
 
41
46
  int i,j,k,m,m1,m2,n;
42
47
  double a,b;
43
48
  static const int np[33] = {0,0,0,0,3,3,4,4,5,4,4,5,5,5,6,6,6,6,6,6,
44
 
                             6,6,6,6,7,7,7,7,8,8,8,8,8};
 
49
                             6,6,6,6,7,7,7,7,8,8,8,8,8};
45
50
  static const int m1p[33]= {0,0,0,0,17,23,23,28,27,38,42,39,43,47,43,47,50,54,
46
 
                               57,60,64,67,71,74,68,71,74,77,70,73,76,78,81};
 
51
                               57,60,64,67,71,74,68,71,74,77,70,73,76,78,81};
47
52
  static const int m1np[7][18] = {
48
 
                 { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
49
 
                 { 0, 0, 0, 0,36,48,60,72, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
50
 
                 { 0, 0, 0, 0,24,32,40,48,56,64,72, 0, 0, 0, 0, 0, 0, 0},
51
 
                 { 0, 0, 0, 0,17,23,29,35,41,47,53,59,65, 0, 0, 0, 0, 0},
52
 
                 { 0, 0, 0, 0, 0, 0,23,28,33,38,42,47,52,57,62,66, 0, 0},
53
 
                 { 0, 0, 0, 0, 0, 0, 0, 0,27, 0, 0,39,43,47,51,55,59,63},
54
 
                 { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,43,47,50,54}};
 
53
                 { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
 
54
                 { 0, 0, 0, 0,36,48,60,72, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},
 
55
                 { 0, 0, 0, 0,24,32,40,48,56,64,72, 0, 0, 0, 0, 0, 0, 0},
 
56
                 { 0, 0, 0, 0,17,23,29,35,41,47,53,59,65, 0, 0, 0, 0, 0},
 
57
                 { 0, 0, 0, 0, 0, 0,23,28,33,38,42,47,52,57,62,66, 0, 0},
 
58
                 { 0, 0, 0, 0, 0, 0, 0, 0,27, 0, 0,39,43,47,51,55,59,63},
 
59
                 { 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,43,47,50,54}};
55
60
  mp_no mpone = {0,{0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
56
 
                    0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
57
 
                    0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0}};
 
61
                    0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
 
62
                    0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0}};
58
63
  mp_no mpk   = {0,{0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
59
 
                    0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
60
 
                    0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0}};
 
64
                    0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,
 
65
                    0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0,0.0}};
61
66
  mp_no mps,mpak,mpt1,mpt2;
62
67
 
63
68
  /* Choose m,n and compute a=2**(-m) */
64
 
  n = np[p];    m1 = m1p[p];    a = twomm1[p].d;
 
69
  n = np[p];    m1 = m1p[p];    a = __mpexp_twomm1[p].d;
65
70
  for (i=0; i<EX; i++)  a *= RADIXI;
66
71
  for (   ; i>EX; i--)  a *= RADIX;
67
72
  b = X[1]*RADIXI;   m2 = 24*EX;
81
86
 
82
87
  /* Evaluate the polynomial. Put result in mpt2 */
83
88
  mpone.e=1;   mpone.d[0]=ONE;   mpone.d[1]=ONE;
84
 
  mpk.e = 1;   mpk.d[0] = ONE;   mpk.d[1]=nn[n].d;
 
89
  mpk.e = 1;   mpk.d[0] = ONE;   mpk.d[1]=__mpexp_nn[n].d;
85
90
  __dvd(&mps,&mpk,&mpt1,p);
86
91
  __add(&mpone,&mpt1,&mpak,p);
87
92
  for (k=n-1; k>1; k--) {
88
93
    __mul(&mps,&mpak,&mpt1,p);
89
 
    mpk.d[1]=nn[k].d;
 
94
    mpk.d[1]=__mpexp_nn[k].d;
90
95
    __dvd(&mpt1,&mpk,&mpt2,p);
91
96
    __add(&mpone,&mpt2,&mpak,p);
92
97
  }