OpenCores
URL https://opencores.org/ocsvn/openrisc_me/openrisc_me/trunk

Subversion Repositories openrisc_me

[/] [openrisc/] [trunk/] [gnu-src/] [gcc-4.5.1/] [gcc/] [testsuite/] [gcc.dg/] [pr36584.c] - Blame information for rev 298

Details | Compare with Previous | View Log

Line No. Rev Author Line
1 298 jeremybenn
/* { dg-do run } */
2
/* { dg-options "-O2 -lm" } */
3
/* { dg-options "-O2 -msse2 -mfpmath=sse" { target { { i?86-*-* x86_64-*-* } && ilp32 } } } */
4
/* { dg-require-effective-target sse2 { target { { i?86-*-* x86_64-*-* } && ilp32 } } } */
5
/* { dg-require-effective-target sse2_runtime { target { { i?86-*-* x86_64-*-* } && ilp32 } } } */
6
 
7
extern double fabs (double);
8
extern void abort (void);
9
 
10
const int MAX_ITERATIONS = 50;
11
const double SMALL_ENOUGH = 1.0e-10;
12
const double RELERROR = 1.0e-12;
13
 
14
typedef struct p
15
{
16
  int ord;
17
  double coef[7];
18
}
19
polynomial;
20
 
21
static double
22
polyeval (double x, int n, double *Coeffs)
23
{
24
  register int i;
25
  double val;
26
 
27
  val = Coeffs[n];
28
  for (i = n - 1; i >= 0; i--)
29
    val = val * x + Coeffs[i];
30
 
31
  return (val);
32
}
33
 
34
static int
35
regula_falsa (int order, double *coef, double a, double b, double *val)
36
{
37
  int its;
38
  double fa, fb, x, fx, lfx;
39
 
40
  fa = polyeval (a, order, coef);
41
  fb = polyeval (b, order, coef);
42
 
43
  if (fa * fb > 0.0)
44
    return 0;
45
 
46
  if (fabs (fa) < SMALL_ENOUGH)
47
    {
48
      *val = a;
49
      return 1;
50
    }
51
 
52
  if (fabs (fb) < SMALL_ENOUGH)
53
    {
54
      *val = b;
55
      return 1;
56
    }
57
 
58
  lfx = fa;
59
 
60
  for (its = 0; its < MAX_ITERATIONS; its++)
61
    {
62
      x = (fb * a - fa * b) / (fb - fa);
63
      fx = polyeval (x, order, coef);
64
      if (fabs (x) > RELERROR)
65
        {
66
          if (fabs (fx / x) < RELERROR)
67
            {
68
              *val = x;
69
              return 1;
70
            }
71
        }
72
      else
73
        {
74
          if (fabs (fx) < RELERROR)
75
            {
76
              *val = x;
77
              return 1;
78
            }
79
        }
80
 
81
      if (fa < 0)
82
        {
83
          if (fx < 0)
84
            {
85
              a = x;
86
              fa = fx;
87
              if ((lfx * fx) > 0)
88
                fb /= 2;
89
            }
90
          else
91
            {
92
              b = x;
93
              fb = fx;
94
              if ((lfx * fx) > 0)
95
                fa /= 2;
96
            }
97
        }
98
      else
99
        {
100
          if (fx < 0)
101
            {
102
              b = x;
103
              fb = fx;
104
              if ((lfx * fx) > 0)
105
                fa /= 2;
106
            }
107
          else
108
            {
109
              a = x;
110
              fa = fx;
111
              if ((lfx * fx) > 0)
112
                fb /= 2;
113
            }
114
        }
115
 
116
      if (fabs (b - a) < RELERROR)
117
        {
118
          *val = x;
119
          return 1;
120
        }
121
 
122
      lfx = fx;
123
    }
124
 
125
  return 0;
126
}
127
 
128
static int
129
numchanges (int np, polynomial * sseq, double a)
130
{
131
  int changes;
132
  double f, lf;
133
  polynomial *s;
134
  changes = 0;
135
 
136
  lf = polyeval (a, sseq[0].ord, sseq[0].coef);
137
 
138
  for (s = sseq + 1; s <= sseq + np; s++)
139
    {
140
      f = polyeval (a, s->ord, s->coef);
141
      if (lf == 0.0 || lf * f < 0)
142
        changes++;
143
 
144
      lf = f;
145
    }
146
 
147
  return changes;
148
}
149
 
150
int
151
sbisect (int np, polynomial * sseq, double min_value, double max_value,
152
         int atmin, int atmax, double *roots)
153
{
154
  double mid;
155
  int n1, n2, its, atmid;
156
 
157
  if ((atmin - atmax) == 1)
158
    {
159
      if (regula_falsa (sseq->ord, sseq->coef, min_value, max_value, roots))
160
        return 1;
161
      else
162
        {
163
          for (its = 0; its < MAX_ITERATIONS; its++)
164
            {
165
              mid = (min_value + max_value) / 2;
166
              atmid = numchanges (np, sseq, mid);
167
              if ((atmid < atmax) || (atmid > atmin))
168
                return 0;
169
 
170
              if (fabs (mid) > RELERROR)
171
                {
172
                  if (fabs ((max_value - min_value) / mid) < RELERROR)
173
                    {
174
                      roots[0] = mid;
175
                      return 1;
176
                    }
177
                }
178
              else
179
                {
180
                  if (fabs (max_value - min_value) < RELERROR)
181
                    {
182
                      roots[0] = mid;
183
                      return 1;
184
                    }
185
                }
186
 
187
              if ((atmin - atmid) == 0)
188
                min_value = mid;
189
              else
190
                max_value = mid;
191
            }
192
 
193
          roots[0] = mid;
194
          return 1;
195
        }
196
    }
197
 
198
  for (its = 0; its < MAX_ITERATIONS; its++)
199
    {
200
      mid = (min_value + max_value) / 2;
201
      atmid = numchanges (np, sseq, mid);
202
      if ((atmid < atmax) || (atmid > atmin))
203
        return 0;
204
 
205
      if (fabs (mid) > RELERROR)
206
        {
207
          if (fabs ((max_value - min_value) / mid) < RELERROR)
208
            {
209
              roots[0] = mid;
210
              return 1;
211
            }
212
        }
213
      else
214
        {
215
          if (fabs (max_value - min_value) < RELERROR)
216
            {
217
              roots[0] = mid;
218
              return 1;
219
            }
220
        }
221
 
222
      n1 = atmin - atmid;
223
      n2 = atmid - atmax;
224
 
225
      if ((n1 != 0) && (n2 != 0))
226
        {
227
          n1 = sbisect (np, sseq, min_value, mid, atmin, atmid, roots);
228
          n2 = sbisect (np, sseq, mid, max_value, atmid, atmax, &roots[n1]);
229
 
230
          return (n1 + n2);
231
        }
232
 
233
      if (n1 == 0)
234
        min_value = mid;
235
      else
236
        max_value = mid;
237
    }
238
 
239
  roots[0] = mid;
240
  return 1;
241
}
242
 
243
int
244
main ()
245
{
246
  polynomial sseq[7] = {
247
    {6, {0.15735259075109281, -5.1185263411378736, 1.8516070705868664,
248
         7.348009172322695, -2.2152395279161343, -2.7543325329350692, 1.0}},
249
    {5, {-0.8530877235229789, 0.61720235686228875, 3.6740045861613475,
250
         -1.4768263519440896, -2.2952771107792245, 1.0}},
251
    {4, {0.13072124257049417, 2.2220687798791126, -1.6299431586726509,
252
         -1.6718404582408546, 1.0}},
253
    {3, {0.86776597575462633, -2.1051099695282511, -0.49008580100694688,
254
         1.0}},
255
    {2, {-11.117984175064155, 10.89886635045883, 1.0}},
256
    {1, {0.94453099602191237, -1.0}},
257
    {0, {-0.068471716890574186}}
258
  };
259
 
260
  double roots[7];
261
  int nroots;
262
 
263
  nroots = sbisect (6, sseq, 0.0, 10000000.0, 5, 1, roots);
264
  if (nroots != 4)
265
    abort ();
266
 
267
  return 0;
268
}

powered by: WebSVN 2.1.0

© copyright 1999-2024 OpenCores.org, equivalent to Oliscience, all rights reserved. OpenCores®, registered trademark.