← Back to branch summary

~tsarev/boostdc/cmake

« back to all changes in this revision

Viewing changes to boost/boost/math/tools/detail/rational_horner3_9.hpp

  • Committer: bigmuscle
  • Date: 2010-05-08 08:47:15 UTC
  • Revision ID: svn-v4:5fb55d53-692c-0410-a46a-e90ab66e00ee:trunk:497
removed old boost version

Show diffs side-by-side

added added

removed removed

Lines of Context:
boost/boost/math/tools/detail/rational_horner3_9.hpp
1
 
//  (C) Copyright John Maddock 2007.
2
 
//  Use, modification and distribution are subject to the
3
 
//  Boost Software License, Version 1.0. (See accompanying file
4
 
//  LICENSE_1_0.txt or copy at http://www.boost.org/LICENSE_1_0.txt)
5
 
//
6
 
//  This file is machine generated, do not edit by hand
7
 
 
8
 
// Polynomial evaluation using second order Horners rule
9
 
#ifndef BOOST_MATH_TOOLS_RAT_EVAL_9_HPP
10
 
#define BOOST_MATH_TOOLS_RAT_EVAL_9_HPP
11
 
 
12
 
namespace boost{ namespace math{ namespace tools{ namespace detail{
13
 
 
14
 
template <class T, class U, class V>
15
 
inline V evaluate_rational_c_imp(const T* a, const U* b, const V&, const mpl::int_<0>*)
16
 
{
17
 
   return static_cast<V>(0);
18
 
}
19
 
 
20
 
template <class T, class U, class V>
21
 
inline V evaluate_rational_c_imp(const T* a, const U* b, const V&, const mpl::int_<1>*)
22
 
{
23
 
   return static_cast<V>(a[0]) / static_cast<V>(b[0]);
24
 
}
25
 
 
26
 
template <class T, class U, class V>
27
 
inline V evaluate_rational_c_imp(const T* a, const U* b, const V& x, const mpl::int_<2>*)
28
 
{
29
 
   return static_cast<V>((a[1] * x + a[0]) / (b[1] * x + b[0]));
30
 
}
31
 
 
32
 
template <class T, class U, class V>
33
 
inline V evaluate_rational_c_imp(const T* a, const U* b, const V& x, const mpl::int_<3>*)
34
 
{
35
 
   return static_cast<V>(((a[2] * x + a[1]) * x + a[0]) / ((b[2] * x + b[1]) * x + b[0]));
36
 
}
37
 
 
38
 
template <class T, class U, class V>
39
 
inline V evaluate_rational_c_imp(const T* a, const U* b, const V& x, const mpl::int_<4>*)
40
 
{
41
 
   return static_cast<V>((((a[3] * x + a[2]) * x + a[1]) * x + a[0]) / (((b[3] * x + b[2]) * x + b[1]) * x + b[0]));
42
 
}
43
 
 
44
 
template <class T, class U, class V>
45
 
inline V evaluate_rational_c_imp(const T* a, const U* b, const V& x, const mpl::int_<5>*)
46
 
{
47
 
   if(x <= 1)
48
 
   {
49
 
      V x2 = x * x;
50
 
      V t[4];
51
 
      t[0] = a[4] * x2 + a[2];
52
 
      t[1] = a[3] * x2 + a[1];
53
 
      t[2] = b[4] * x2 + b[2];
54
 
      t[3] = b[3] * x2 + b[1];
55
 
      t[0] *= x2;
56
 
      t[2] *= x2;
57
 
      t[0] += static_cast<V>(a[0]);
58
 
      t[2] += static_cast<V>(b[0]);
59
 
      t[1] *= x;
60
 
      t[3] *= x;
61
 
      return (t[0] + t[1]) / (t[2] + t[3]);
62
 
   }
63
 
   else
64
 
   {
65
 
      V z = 1 / x;
66
 
      V z2 = 1 / (x * x);
67
 
      V t[4];
68
 
      t[0] = a[0] * z2 + a[2];
69
 
      t[1] = a[1] * z2 + a[3];
70
 
      t[2] = b[0] * z2 + b[2];
71
 
      t[3] = b[1] * z2 + b[3];
72
 
      t[0] *= z2;
73
 
      t[2] *= z2;
74
 
      t[0] += static_cast<V>(a[4]);
75
 
      t[2] += static_cast<V>(b[4]);
76
 
      t[1] *= z;
77
 
      t[3] *= z;
78
 
      return (t[0] + t[1]) / (t[2] + t[3]);
79
 
   }
80
 
}
81
 
 
82
 
template <class T, class U, class V>
83
 
inline V evaluate_rational_c_imp(const T* a, const U* b, const V& x, const mpl::int_<6>*)
84
 
{
85
 
   if(x <= 1)
86
 
   {
87
 
      V x2 = x * x;
88
 
      V t[4];
89
 
      t[0] = a[5] * x2 + a[3];
90
 
      t[1] = a[4] * x2 + a[2];
91
 
      t[2] = b[5] * x2 + b[3];
92
 
      t[3] = b[4] * x2 + b[2];
93
 
      t[0] *= x2;
94
 
      t[1] *= x2;
95
 
      t[2] *= x2;
96
 
      t[3] *= x2;
97
 
      t[0] += static_cast<V>(a[1]);
98
 
      t[1] += static_cast<V>(a[0]);
99
 
      t[2] += static_cast<V>(b[1]);
100
 
      t[3] += static_cast<V>(b[0]);
101
 
      t[0] *= x;
102
 
      t[2] *= x;
103
 
      return (t[0] + t[1]) / (t[2] + t[3]);
104
 
   }
105
 
   else
106
 
   {
107
 
      V z = 1 / x;
108
 
      V z2 = 1 / (x * x);
109
 
      V t[4];
110
 
      t[0] = a[0] * z2 + a[2];
111
 
      t[1] = a[1] * z2 + a[3];
112
 
      t[2] = b[0] * z2 + b[2];
113
 
      t[3] = b[1] * z2 + b[3];
114
 
      t[0] *= z2;
115
 
      t[1] *= z2;
116
 
      t[2] *= z2;
117
 
      t[3] *= z2;
118
 
      t[0] += static_cast<V>(a[4]);
119
 
      t[1] += static_cast<V>(a[5]);
120
 
      t[2] += static_cast<V>(b[4]);
121
 
      t[3] += static_cast<V>(b[5]);
122
 
      t[0] *= z;
123
 
      t[2] *= z;
124
 
      return (t[0] + t[1]) / (t[2] + t[3]);
125
 
   }
126
 
}
127
 
 
128
 
template <class T, class U, class V>
129
 
inline V evaluate_rational_c_imp(const T* a, const U* b, const V& x, const mpl::int_<7>*)
130
 
{
131
 
   if(x <= 1)
132
 
   {
133
 
      V x2 = x * x;
134
 
      V t[4];
135
 
      t[0] = a[6] * x2 + a[4];
136
 
      t[1] = a[5] * x2 + a[3];
137
 
      t[2] = b[6] * x2 + b[4];
138
 
      t[3] = b[5] * x2 + b[3];
139
 
      t[0] *= x2;
140
 
      t[1] *= x2;
141
 
      t[2] *= x2;
142
 
      t[3] *= x2;
143
 
      t[0] += static_cast<V>(a[2]);
144
 
      t[1] += static_cast<V>(a[1]);
145
 
      t[2] += static_cast<V>(b[2]);
146
 
      t[3] += static_cast<V>(b[1]);
147
 
      t[0] *= x2;
148
 
      t[2] *= x2;
149
 
      t[0] += static_cast<V>(a[0]);
150
 
      t[2] += static_cast<V>(b[0]);
151
 
      t[1] *= x;
152
 
      t[3] *= x;
153
 
      return (t[0] + t[1]) / (t[2] + t[3]);
154
 
   }
155
 
   else
156
 
   {
157
 
      V z = 1 / x;
158
 
      V z2 = 1 / (x * x);
159
 
      V t[4];
160
 
      t[0] = a[0] * z2 + a[2];
161
 
      t[1] = a[1] * z2 + a[3];
162
 
      t[2] = b[0] * z2 + b[2];
163
 
      t[3] = b[1] * z2 + b[3];
164
 
      t[0] *= z2;
165
 
      t[1] *= z2;
166
 
      t[2] *= z2;
167
 
      t[3] *= z2;
168
 
      t[0] += static_cast<V>(a[4]);
169
 
      t[1] += static_cast<V>(a[5]);
170
 
      t[2] += static_cast<V>(b[4]);
171
 
      t[3] += static_cast<V>(b[5]);
172
 
      t[0] *= z2;
173
 
      t[2] *= z2;
174
 
      t[0] += static_cast<V>(a[6]);
175
 
      t[2] += static_cast<V>(b[6]);
176
 
      t[1] *= z;
177
 
      t[3] *= z;
178
 
      return (t[0] + t[1]) / (t[2] + t[3]);
179
 
   }
180
 
}
181
 
 
182
 
template <class T, class U, class V>
183
 
inline V evaluate_rational_c_imp(const T* a, const U* b, const V& x, const mpl::int_<8>*)
184
 
{
185
 
   if(x <= 1)
186
 
   {
187
 
      V x2 = x * x;
188
 
      V t[4];
189
 
      t[0] = a[7] * x2 + a[5];
190
 
      t[1] = a[6] * x2 + a[4];
191
 
      t[2] = b[7] * x2 + b[5];
192
 
      t[3] = b[6] * x2 + b[4];
193
 
      t[0] *= x2;
194
 
      t[1] *= x2;
195
 
      t[2] *= x2;
196
 
      t[3] *= x2;
197
 
      t[0] += static_cast<V>(a[3]);
198
 
      t[1] += static_cast<V>(a[2]);
199
 
      t[2] += static_cast<V>(b[3]);
200
 
      t[3] += static_cast<V>(b[2]);
201
 
      t[0] *= x2;
202
 
      t[1] *= x2;
203
 
      t[2] *= x2;
204
 
      t[3] *= x2;
205
 
      t[0] += static_cast<V>(a[1]);
206
 
      t[1] += static_cast<V>(a[0]);
207
 
      t[2] += static_cast<V>(b[1]);
208
 
      t[3] += static_cast<V>(b[0]);
209
 
      t[0] *= x;
210
 
      t[2] *= x;
211
 
      return (t[0] + t[1]) / (t[2] + t[3]);
212
 
   }
213
 
   else
214
 
   {
215
 
      V z = 1 / x;
216
 
      V z2 = 1 / (x * x);
217
 
      V t[4];
218
 
      t[0] = a[0] * z2 + a[2];
219
 
      t[1] = a[1] * z2 + a[3];
220
 
      t[2] = b[0] * z2 + b[2];
221
 
      t[3] = b[1] * z2 + b[3];
222
 
      t[0] *= z2;
223
 
      t[1] *= z2;
224
 
      t[2] *= z2;
225
 
      t[3] *= z2;
226
 
      t[0] += static_cast<V>(a[4]);
227
 
      t[1] += static_cast<V>(a[5]);
228
 
      t[2] += static_cast<V>(b[4]);
229
 
      t[3] += static_cast<V>(b[5]);
230
 
      t[0] *= z2;
231
 
      t[1] *= z2;
232
 
      t[2] *= z2;
233
 
      t[3] *= z2;
234
 
      t[0] += static_cast<V>(a[6]);
235
 
      t[1] += static_cast<V>(a[7]);
236
 
      t[2] += static_cast<V>(b[6]);
237
 
      t[3] += static_cast<V>(b[7]);
238
 
      t[0] *= z;
239
 
      t[2] *= z;
240
 
      return (t[0] + t[1]) / (t[2] + t[3]);
241
 
   }
242
 
}
243
 
 
244
 
template <class T, class U, class V>
245
 
inline V evaluate_rational_c_imp(const T* a, const U* b, const V& x, const mpl::int_<9>*)
246
 
{
247
 
   if(x <= 1)
248
 
   {
249
 
      V x2 = x * x;
250
 
      V t[4];
251
 
      t[0] = a[8] * x2 + a[6];
252
 
      t[1] = a[7] * x2 + a[5];
253
 
      t[2] = b[8] * x2 + b[6];
254
 
      t[3] = b[7] * x2 + b[5];
255
 
      t[0] *= x2;
256
 
      t[1] *= x2;
257
 
      t[2] *= x2;
258
 
      t[3] *= x2;
259
 
      t[0] += static_cast<V>(a[4]);
260
 
      t[1] += static_cast<V>(a[3]);
261
 
      t[2] += static_cast<V>(b[4]);
262
 
      t[3] += static_cast<V>(b[3]);
263
 
      t[0] *= x2;
264
 
      t[1] *= x2;
265
 
      t[2] *= x2;
266
 
      t[3] *= x2;
267
 
      t[0] += static_cast<V>(a[2]);
268
 
      t[1] += static_cast<V>(a[1]);
269
 
      t[2] += static_cast<V>(b[2]);
270
 
      t[3] += static_cast<V>(b[1]);
271
 
      t[0] *= x2;
272
 
      t[2] *= x2;
273
 
      t[0] += static_cast<V>(a[0]);
274
 
      t[2] += static_cast<V>(b[0]);
275
 
      t[1] *= x;
276
 
      t[3] *= x;
277
 
      return (t[0] + t[1]) / (t[2] + t[3]);
278
 
   }
279
 
   else
280
 
   {
281
 
      V z = 1 / x;
282
 
      V z2 = 1 / (x * x);
283
 
      V t[4];
284
 
      t[0] = a[0] * z2 + a[2];
285
 
      t[1] = a[1] * z2 + a[3];
286
 
      t[2] = b[0] * z2 + b[2];
287
 
      t[3] = b[1] * z2 + b[3];
288
 
      t[0] *= z2;
289
 
      t[1] *= z2;
290
 
      t[2] *= z2;
291
 
      t[3] *= z2;
292
 
      t[0] += static_cast<V>(a[4]);
293
 
      t[1] += static_cast<V>(a[5]);
294
 
      t[2] += static_cast<V>(b[4]);
295
 
      t[3] += static_cast<V>(b[5]);
296
 
      t[0] *= z2;
297
 
      t[1] *= z2;
298
 
      t[2] *= z2;
299
 
      t[3] *= z2;
300
 
      t[0] += static_cast<V>(a[6]);
301
 
      t[1] += static_cast<V>(a[7]);
302
 
      t[2] += static_cast<V>(b[6]);
303
 
      t[3] += static_cast<V>(b[7]);
304
 
      t[0] *= z2;
305
 
      t[2] *= z2;
306
 
      t[0] += static_cast<V>(a[8]);
307
 
      t[2] += static_cast<V>(b[8]);
308
 
      t[1] *= z;
309
 
      t[3] *= z;
310
 
      return (t[0] + t[1]) / (t[2] + t[3]);
311
 
   }
312
 
}
313
 
 
314
 
 
315
 
}}}} // namespaces
316
 
 
317
 
#endif // include guard
318