1 | /*************************************************************************
|
---|
2 | Cephes Math Library Release 2.8: June, 2000
|
---|
3 | Copyright by Stephen L. Moshier
|
---|
4 |
|
---|
5 | Contributors:
|
---|
6 | * Sergey Bochkanov (ALGLIB project). Translation from C to
|
---|
7 | pseudocode.
|
---|
8 |
|
---|
9 | See subroutines comments for additional copyrights.
|
---|
10 |
|
---|
11 | >>> SOURCE LICENSE >>>
|
---|
12 | This program is free software; you can redistribute it and/or modify
|
---|
13 | it under the terms of the GNU General Public License as published by
|
---|
14 | the Free Software Foundation (www.fsf.org); either version 2 of the
|
---|
15 | License, or (at your option) any later version.
|
---|
16 |
|
---|
17 | This program is distributed in the hope that it will be useful,
|
---|
18 | but WITHOUT ANY WARRANTY; without even the implied warranty of
|
---|
19 | MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
---|
20 | GNU General Public License for more details.
|
---|
21 |
|
---|
22 | A copy of the GNU General Public License is available at
|
---|
23 | http://www.fsf.org/licensing/licenses
|
---|
24 |
|
---|
25 | >>> END OF LICENSE >>>
|
---|
26 | *************************************************************************/
|
---|
27 |
|
---|
28 | using System;
|
---|
29 |
|
---|
30 | namespace alglib
|
---|
31 | {
|
---|
32 | public class trigintegrals
|
---|
33 | {
|
---|
34 | /*************************************************************************
|
---|
35 | Sine and cosine integrals
|
---|
36 |
|
---|
37 | Evaluates the integrals
|
---|
38 |
|
---|
39 | x
|
---|
40 | -
|
---|
41 | | cos t - 1
|
---|
42 | Ci(x) = eul + ln x + | --------- dt,
|
---|
43 | | t
|
---|
44 | -
|
---|
45 | 0
|
---|
46 | x
|
---|
47 | -
|
---|
48 | | sin t
|
---|
49 | Si(x) = | ----- dt
|
---|
50 | | t
|
---|
51 | -
|
---|
52 | 0
|
---|
53 |
|
---|
54 | where eul = 0.57721566490153286061 is Euler's constant.
|
---|
55 | The integrals are approximated by rational functions.
|
---|
56 | For x > 8 auxiliary functions f(x) and g(x) are employed
|
---|
57 | such that
|
---|
58 |
|
---|
59 | Ci(x) = f(x) sin(x) - g(x) cos(x)
|
---|
60 | Si(x) = pi/2 - f(x) cos(x) - g(x) sin(x)
|
---|
61 |
|
---|
62 |
|
---|
63 | ACCURACY:
|
---|
64 | Test interval = [0,50].
|
---|
65 | Absolute error, except relative when > 1:
|
---|
66 | arithmetic function # trials peak rms
|
---|
67 | IEEE Si 30000 4.4e-16 7.3e-17
|
---|
68 | IEEE Ci 30000 6.9e-16 5.1e-17
|
---|
69 |
|
---|
70 | Cephes Math Library Release 2.1: January, 1989
|
---|
71 | Copyright 1984, 1987, 1989 by Stephen L. Moshier
|
---|
72 | *************************************************************************/
|
---|
73 | public static void sinecosineintegrals(double x,
|
---|
74 | ref double si,
|
---|
75 | ref double ci)
|
---|
76 | {
|
---|
77 | double z = 0;
|
---|
78 | double c = 0;
|
---|
79 | double s = 0;
|
---|
80 | double f = 0;
|
---|
81 | double g = 0;
|
---|
82 | int sg = 0;
|
---|
83 | double sn = 0;
|
---|
84 | double sd = 0;
|
---|
85 | double cn = 0;
|
---|
86 | double cd = 0;
|
---|
87 | double fn = 0;
|
---|
88 | double fd = 0;
|
---|
89 | double gn = 0;
|
---|
90 | double gd = 0;
|
---|
91 |
|
---|
92 | if( (double)(x)<(double)(0) )
|
---|
93 | {
|
---|
94 | sg = -1;
|
---|
95 | x = -x;
|
---|
96 | }
|
---|
97 | else
|
---|
98 | {
|
---|
99 | sg = 0;
|
---|
100 | }
|
---|
101 | if( (double)(x)==(double)(0) )
|
---|
102 | {
|
---|
103 | si = 0;
|
---|
104 | ci = -AP.Math.MaxRealNumber;
|
---|
105 | return;
|
---|
106 | }
|
---|
107 | if( (double)(x)>(double)(1.0E9) )
|
---|
108 | {
|
---|
109 | si = 1.570796326794896619-Math.Cos(x)/x;
|
---|
110 | ci = Math.Sin(x)/x;
|
---|
111 | return;
|
---|
112 | }
|
---|
113 | if( (double)(x)<=(double)(4) )
|
---|
114 | {
|
---|
115 | z = x*x;
|
---|
116 | sn = -8.39167827910303881427E-11;
|
---|
117 | sn = sn*z+4.62591714427012837309E-8;
|
---|
118 | sn = sn*z-9.75759303843632795789E-6;
|
---|
119 | sn = sn*z+9.76945438170435310816E-4;
|
---|
120 | sn = sn*z-4.13470316229406538752E-2;
|
---|
121 | sn = sn*z+1.00000000000000000302E0;
|
---|
122 | sd = 2.03269266195951942049E-12;
|
---|
123 | sd = sd*z+1.27997891179943299903E-9;
|
---|
124 | sd = sd*z+4.41827842801218905784E-7;
|
---|
125 | sd = sd*z+9.96412122043875552487E-5;
|
---|
126 | sd = sd*z+1.42085239326149893930E-2;
|
---|
127 | sd = sd*z+9.99999999999999996984E-1;
|
---|
128 | s = x*sn/sd;
|
---|
129 | cn = 2.02524002389102268789E-11;
|
---|
130 | cn = cn*z-1.35249504915790756375E-8;
|
---|
131 | cn = cn*z+3.59325051419993077021E-6;
|
---|
132 | cn = cn*z-4.74007206873407909465E-4;
|
---|
133 | cn = cn*z+2.89159652607555242092E-2;
|
---|
134 | cn = cn*z-1.00000000000000000080E0;
|
---|
135 | cd = 4.07746040061880559506E-12;
|
---|
136 | cd = cd*z+3.06780997581887812692E-9;
|
---|
137 | cd = cd*z+1.23210355685883423679E-6;
|
---|
138 | cd = cd*z+3.17442024775032769882E-4;
|
---|
139 | cd = cd*z+5.10028056236446052392E-2;
|
---|
140 | cd = cd*z+4.00000000000000000080E0;
|
---|
141 | c = z*cn/cd;
|
---|
142 | if( sg!=0 )
|
---|
143 | {
|
---|
144 | s = -s;
|
---|
145 | }
|
---|
146 | si = s;
|
---|
147 | ci = 0.57721566490153286061+Math.Log(x)+c;
|
---|
148 | return;
|
---|
149 | }
|
---|
150 | s = Math.Sin(x);
|
---|
151 | c = Math.Cos(x);
|
---|
152 | z = 1.0/(x*x);
|
---|
153 | if( (double)(x)<(double)(8) )
|
---|
154 | {
|
---|
155 | fn = 4.23612862892216586994E0;
|
---|
156 | fn = fn*z+5.45937717161812843388E0;
|
---|
157 | fn = fn*z+1.62083287701538329132E0;
|
---|
158 | fn = fn*z+1.67006611831323023771E-1;
|
---|
159 | fn = fn*z+6.81020132472518137426E-3;
|
---|
160 | fn = fn*z+1.08936580650328664411E-4;
|
---|
161 | fn = fn*z+5.48900223421373614008E-7;
|
---|
162 | fd = 1.00000000000000000000E0;
|
---|
163 | fd = fd*z+8.16496634205391016773E0;
|
---|
164 | fd = fd*z+7.30828822505564552187E0;
|
---|
165 | fd = fd*z+1.86792257950184183883E0;
|
---|
166 | fd = fd*z+1.78792052963149907262E-1;
|
---|
167 | fd = fd*z+7.01710668322789753610E-3;
|
---|
168 | fd = fd*z+1.10034357153915731354E-4;
|
---|
169 | fd = fd*z+5.48900252756255700982E-7;
|
---|
170 | f = fn/(x*fd);
|
---|
171 | gn = 8.71001698973114191777E-2;
|
---|
172 | gn = gn*z+6.11379109952219284151E-1;
|
---|
173 | gn = gn*z+3.97180296392337498885E-1;
|
---|
174 | gn = gn*z+7.48527737628469092119E-2;
|
---|
175 | gn = gn*z+5.38868681462177273157E-3;
|
---|
176 | gn = gn*z+1.61999794598934024525E-4;
|
---|
177 | gn = gn*z+1.97963874140963632189E-6;
|
---|
178 | gn = gn*z+7.82579040744090311069E-9;
|
---|
179 | gd = 1.00000000000000000000E0;
|
---|
180 | gd = gd*z+1.64402202413355338886E0;
|
---|
181 | gd = gd*z+6.66296701268987968381E-1;
|
---|
182 | gd = gd*z+9.88771761277688796203E-2;
|
---|
183 | gd = gd*z+6.22396345441768420760E-3;
|
---|
184 | gd = gd*z+1.73221081474177119497E-4;
|
---|
185 | gd = gd*z+2.02659182086343991969E-6;
|
---|
186 | gd = gd*z+7.82579218933534490868E-9;
|
---|
187 | g = z*gn/gd;
|
---|
188 | }
|
---|
189 | else
|
---|
190 | {
|
---|
191 | fn = 4.55880873470465315206E-1;
|
---|
192 | fn = fn*z+7.13715274100146711374E-1;
|
---|
193 | fn = fn*z+1.60300158222319456320E-1;
|
---|
194 | fn = fn*z+1.16064229408124407915E-2;
|
---|
195 | fn = fn*z+3.49556442447859055605E-4;
|
---|
196 | fn = fn*z+4.86215430826454749482E-6;
|
---|
197 | fn = fn*z+3.20092790091004902806E-8;
|
---|
198 | fn = fn*z+9.41779576128512936592E-11;
|
---|
199 | fn = fn*z+9.70507110881952024631E-14;
|
---|
200 | fd = 1.00000000000000000000E0;
|
---|
201 | fd = fd*z+9.17463611873684053703E-1;
|
---|
202 | fd = fd*z+1.78685545332074536321E-1;
|
---|
203 | fd = fd*z+1.22253594771971293032E-2;
|
---|
204 | fd = fd*z+3.58696481881851580297E-4;
|
---|
205 | fd = fd*z+4.92435064317881464393E-6;
|
---|
206 | fd = fd*z+3.21956939101046018377E-8;
|
---|
207 | fd = fd*z+9.43720590350276732376E-11;
|
---|
208 | fd = fd*z+9.70507110881952025725E-14;
|
---|
209 | f = fn/(x*fd);
|
---|
210 | gn = 6.97359953443276214934E-1;
|
---|
211 | gn = gn*z+3.30410979305632063225E-1;
|
---|
212 | gn = gn*z+3.84878767649974295920E-2;
|
---|
213 | gn = gn*z+1.71718239052347903558E-3;
|
---|
214 | gn = gn*z+3.48941165502279436777E-5;
|
---|
215 | gn = gn*z+3.47131167084116673800E-7;
|
---|
216 | gn = gn*z+1.70404452782044526189E-9;
|
---|
217 | gn = gn*z+3.85945925430276600453E-12;
|
---|
218 | gn = gn*z+3.14040098946363334640E-15;
|
---|
219 | gd = 1.00000000000000000000E0;
|
---|
220 | gd = gd*z+1.68548898811011640017E0;
|
---|
221 | gd = gd*z+4.87852258695304967486E-1;
|
---|
222 | gd = gd*z+4.67913194259625806320E-2;
|
---|
223 | gd = gd*z+1.90284426674399523638E-3;
|
---|
224 | gd = gd*z+3.68475504442561108162E-5;
|
---|
225 | gd = gd*z+3.57043223443740838771E-7;
|
---|
226 | gd = gd*z+1.72693748966316146736E-9;
|
---|
227 | gd = gd*z+3.87830166023954706752E-12;
|
---|
228 | gd = gd*z+3.14040098946363335242E-15;
|
---|
229 | g = z*gn/gd;
|
---|
230 | }
|
---|
231 | si = 1.570796326794896619-f*c-g*s;
|
---|
232 | if( sg!=0 )
|
---|
233 | {
|
---|
234 | si = -si;
|
---|
235 | }
|
---|
236 | ci = f*s-g*c;
|
---|
237 | }
|
---|
238 |
|
---|
239 |
|
---|
240 | /*************************************************************************
|
---|
241 | Hyperbolic sine and cosine integrals
|
---|
242 |
|
---|
243 | Approximates the integrals
|
---|
244 |
|
---|
245 | x
|
---|
246 | -
|
---|
247 | | | cosh t - 1
|
---|
248 | Chi(x) = eul + ln x + | ----------- dt,
|
---|
249 | | | t
|
---|
250 | -
|
---|
251 | 0
|
---|
252 |
|
---|
253 | x
|
---|
254 | -
|
---|
255 | | | sinh t
|
---|
256 | Shi(x) = | ------ dt
|
---|
257 | | | t
|
---|
258 | -
|
---|
259 | 0
|
---|
260 |
|
---|
261 | where eul = 0.57721566490153286061 is Euler's constant.
|
---|
262 | The integrals are evaluated by power series for x < 8
|
---|
263 | and by Chebyshev expansions for x between 8 and 88.
|
---|
264 | For large x, both functions approach exp(x)/2x.
|
---|
265 | Arguments greater than 88 in magnitude return MAXNUM.
|
---|
266 |
|
---|
267 |
|
---|
268 | ACCURACY:
|
---|
269 |
|
---|
270 | Test interval 0 to 88.
|
---|
271 | Relative error:
|
---|
272 | arithmetic function # trials peak rms
|
---|
273 | IEEE Shi 30000 6.9e-16 1.6e-16
|
---|
274 | Absolute error, except relative when |Chi| > 1:
|
---|
275 | IEEE Chi 30000 8.4e-16 1.4e-16
|
---|
276 |
|
---|
277 | Cephes Math Library Release 2.8: June, 2000
|
---|
278 | Copyright 1984, 1987, 2000 by Stephen L. Moshier
|
---|
279 | *************************************************************************/
|
---|
280 | public static void hyperbolicsinecosineintegrals(double x,
|
---|
281 | ref double shi,
|
---|
282 | ref double chi)
|
---|
283 | {
|
---|
284 | double k = 0;
|
---|
285 | double z = 0;
|
---|
286 | double c = 0;
|
---|
287 | double s = 0;
|
---|
288 | double a = 0;
|
---|
289 | int sg = 0;
|
---|
290 | double b0 = 0;
|
---|
291 | double b1 = 0;
|
---|
292 | double b2 = 0;
|
---|
293 |
|
---|
294 | if( (double)(x)<(double)(0) )
|
---|
295 | {
|
---|
296 | sg = -1;
|
---|
297 | x = -x;
|
---|
298 | }
|
---|
299 | else
|
---|
300 | {
|
---|
301 | sg = 0;
|
---|
302 | }
|
---|
303 | if( (double)(x)==(double)(0) )
|
---|
304 | {
|
---|
305 | shi = 0;
|
---|
306 | chi = -AP.Math.MaxRealNumber;
|
---|
307 | return;
|
---|
308 | }
|
---|
309 | if( (double)(x)<(double)(8.0) )
|
---|
310 | {
|
---|
311 | z = x*x;
|
---|
312 | a = 1.0;
|
---|
313 | s = 1.0;
|
---|
314 | c = 0.0;
|
---|
315 | k = 2.0;
|
---|
316 | do
|
---|
317 | {
|
---|
318 | a = a*z/k;
|
---|
319 | c = c+a/k;
|
---|
320 | k = k+1.0;
|
---|
321 | a = a/k;
|
---|
322 | s = s+a/k;
|
---|
323 | k = k+1.0;
|
---|
324 | }
|
---|
325 | while( (double)(Math.Abs(a/s))>=(double)(AP.Math.MachineEpsilon) );
|
---|
326 | s = s*x;
|
---|
327 | }
|
---|
328 | else
|
---|
329 | {
|
---|
330 | if( (double)(x)<(double)(18.0) )
|
---|
331 | {
|
---|
332 | a = (576.0/x-52.0)/10.0;
|
---|
333 | k = Math.Exp(x)/x;
|
---|
334 | b0 = 1.83889230173399459482E-17;
|
---|
335 | b1 = 0.0;
|
---|
336 | chebiterationshichi(a, -9.55485532279655569575E-17, ref b0, ref b1, ref b2);
|
---|
337 | chebiterationshichi(a, 2.04326105980879882648E-16, ref b0, ref b1, ref b2);
|
---|
338 | chebiterationshichi(a, 1.09896949074905343022E-15, ref b0, ref b1, ref b2);
|
---|
339 | chebiterationshichi(a, -1.31313534344092599234E-14, ref b0, ref b1, ref b2);
|
---|
340 | chebiterationshichi(a, 5.93976226264314278932E-14, ref b0, ref b1, ref b2);
|
---|
341 | chebiterationshichi(a, -3.47197010497749154755E-14, ref b0, ref b1, ref b2);
|
---|
342 | chebiterationshichi(a, -1.40059764613117131000E-12, ref b0, ref b1, ref b2);
|
---|
343 | chebiterationshichi(a, 9.49044626224223543299E-12, ref b0, ref b1, ref b2);
|
---|
344 | chebiterationshichi(a, -1.61596181145435454033E-11, ref b0, ref b1, ref b2);
|
---|
345 | chebiterationshichi(a, -1.77899784436430310321E-10, ref b0, ref b1, ref b2);
|
---|
346 | chebiterationshichi(a, 1.35455469767246947469E-9, ref b0, ref b1, ref b2);
|
---|
347 | chebiterationshichi(a, -1.03257121792819495123E-9, ref b0, ref b1, ref b2);
|
---|
348 | chebiterationshichi(a, -3.56699611114982536845E-8, ref b0, ref b1, ref b2);
|
---|
349 | chebiterationshichi(a, 1.44818877384267342057E-7, ref b0, ref b1, ref b2);
|
---|
350 | chebiterationshichi(a, 7.82018215184051295296E-7, ref b0, ref b1, ref b2);
|
---|
351 | chebiterationshichi(a, -5.39919118403805073710E-6, ref b0, ref b1, ref b2);
|
---|
352 | chebiterationshichi(a, -3.12458202168959833422E-5, ref b0, ref b1, ref b2);
|
---|
353 | chebiterationshichi(a, 8.90136741950727517826E-5, ref b0, ref b1, ref b2);
|
---|
354 | chebiterationshichi(a, 2.02558474743846862168E-3, ref b0, ref b1, ref b2);
|
---|
355 | chebiterationshichi(a, 2.96064440855633256972E-2, ref b0, ref b1, ref b2);
|
---|
356 | chebiterationshichi(a, 1.11847751047257036625E0, ref b0, ref b1, ref b2);
|
---|
357 | s = k*0.5*(b0-b2);
|
---|
358 | b0 = -8.12435385225864036372E-18;
|
---|
359 | b1 = 0.0;
|
---|
360 | chebiterationshichi(a, 2.17586413290339214377E-17, ref b0, ref b1, ref b2);
|
---|
361 | chebiterationshichi(a, 5.22624394924072204667E-17, ref b0, ref b1, ref b2);
|
---|
362 | chebiterationshichi(a, -9.48812110591690559363E-16, ref b0, ref b1, ref b2);
|
---|
363 | chebiterationshichi(a, 5.35546311647465209166E-15, ref b0, ref b1, ref b2);
|
---|
364 | chebiterationshichi(a, -1.21009970113732918701E-14, ref b0, ref b1, ref b2);
|
---|
365 | chebiterationshichi(a, -6.00865178553447437951E-14, ref b0, ref b1, ref b2);
|
---|
366 | chebiterationshichi(a, 7.16339649156028587775E-13, ref b0, ref b1, ref b2);
|
---|
367 | chebiterationshichi(a, -2.93496072607599856104E-12, ref b0, ref b1, ref b2);
|
---|
368 | chebiterationshichi(a, -1.40359438136491256904E-12, ref b0, ref b1, ref b2);
|
---|
369 | chebiterationshichi(a, 8.76302288609054966081E-11, ref b0, ref b1, ref b2);
|
---|
370 | chebiterationshichi(a, -4.40092476213282340617E-10, ref b0, ref b1, ref b2);
|
---|
371 | chebiterationshichi(a, -1.87992075640569295479E-10, ref b0, ref b1, ref b2);
|
---|
372 | chebiterationshichi(a, 1.31458150989474594064E-8, ref b0, ref b1, ref b2);
|
---|
373 | chebiterationshichi(a, -4.75513930924765465590E-8, ref b0, ref b1, ref b2);
|
---|
374 | chebiterationshichi(a, -2.21775018801848880741E-7, ref b0, ref b1, ref b2);
|
---|
375 | chebiterationshichi(a, 1.94635531373272490962E-6, ref b0, ref b1, ref b2);
|
---|
376 | chebiterationshichi(a, 4.33505889257316408893E-6, ref b0, ref b1, ref b2);
|
---|
377 | chebiterationshichi(a, -6.13387001076494349496E-5, ref b0, ref b1, ref b2);
|
---|
378 | chebiterationshichi(a, -3.13085477492997465138E-4, ref b0, ref b1, ref b2);
|
---|
379 | chebiterationshichi(a, 4.97164789823116062801E-4, ref b0, ref b1, ref b2);
|
---|
380 | chebiterationshichi(a, 2.64347496031374526641E-2, ref b0, ref b1, ref b2);
|
---|
381 | chebiterationshichi(a, 1.11446150876699213025E0, ref b0, ref b1, ref b2);
|
---|
382 | c = k*0.5*(b0-b2);
|
---|
383 | }
|
---|
384 | else
|
---|
385 | {
|
---|
386 | if( (double)(x)<=(double)(88.0) )
|
---|
387 | {
|
---|
388 | a = (6336.0/x-212.0)/70.0;
|
---|
389 | k = Math.Exp(x)/x;
|
---|
390 | b0 = -1.05311574154850938805E-17;
|
---|
391 | b1 = 0.0;
|
---|
392 | chebiterationshichi(a, 2.62446095596355225821E-17, ref b0, ref b1, ref b2);
|
---|
393 | chebiterationshichi(a, 8.82090135625368160657E-17, ref b0, ref b1, ref b2);
|
---|
394 | chebiterationshichi(a, -3.38459811878103047136E-16, ref b0, ref b1, ref b2);
|
---|
395 | chebiterationshichi(a, -8.30608026366935789136E-16, ref b0, ref b1, ref b2);
|
---|
396 | chebiterationshichi(a, 3.93397875437050071776E-15, ref b0, ref b1, ref b2);
|
---|
397 | chebiterationshichi(a, 1.01765565969729044505E-14, ref b0, ref b1, ref b2);
|
---|
398 | chebiterationshichi(a, -4.21128170307640802703E-14, ref b0, ref b1, ref b2);
|
---|
399 | chebiterationshichi(a, -1.60818204519802480035E-13, ref b0, ref b1, ref b2);
|
---|
400 | chebiterationshichi(a, 3.34714954175994481761E-13, ref b0, ref b1, ref b2);
|
---|
401 | chebiterationshichi(a, 2.72600352129153073807E-12, ref b0, ref b1, ref b2);
|
---|
402 | chebiterationshichi(a, 1.66894954752839083608E-12, ref b0, ref b1, ref b2);
|
---|
403 | chebiterationshichi(a, -3.49278141024730899554E-11, ref b0, ref b1, ref b2);
|
---|
404 | chebiterationshichi(a, -1.58580661666482709598E-10, ref b0, ref b1, ref b2);
|
---|
405 | chebiterationshichi(a, -1.79289437183355633342E-10, ref b0, ref b1, ref b2);
|
---|
406 | chebiterationshichi(a, 1.76281629144264523277E-9, ref b0, ref b1, ref b2);
|
---|
407 | chebiterationshichi(a, 1.69050228879421288846E-8, ref b0, ref b1, ref b2);
|
---|
408 | chebiterationshichi(a, 1.25391771228487041649E-7, ref b0, ref b1, ref b2);
|
---|
409 | chebiterationshichi(a, 1.16229947068677338732E-6, ref b0, ref b1, ref b2);
|
---|
410 | chebiterationshichi(a, 1.61038260117376323993E-5, ref b0, ref b1, ref b2);
|
---|
411 | chebiterationshichi(a, 3.49810375601053973070E-4, ref b0, ref b1, ref b2);
|
---|
412 | chebiterationshichi(a, 1.28478065259647610779E-2, ref b0, ref b1, ref b2);
|
---|
413 | chebiterationshichi(a, 1.03665722588798326712E0, ref b0, ref b1, ref b2);
|
---|
414 | s = k*0.5*(b0-b2);
|
---|
415 | b0 = 8.06913408255155572081E-18;
|
---|
416 | b1 = 0.0;
|
---|
417 | chebiterationshichi(a, -2.08074168180148170312E-17, ref b0, ref b1, ref b2);
|
---|
418 | chebiterationshichi(a, -5.98111329658272336816E-17, ref b0, ref b1, ref b2);
|
---|
419 | chebiterationshichi(a, 2.68533951085945765591E-16, ref b0, ref b1, ref b2);
|
---|
420 | chebiterationshichi(a, 4.52313941698904694774E-16, ref b0, ref b1, ref b2);
|
---|
421 | chebiterationshichi(a, -3.10734917335299464535E-15, ref b0, ref b1, ref b2);
|
---|
422 | chebiterationshichi(a, -4.42823207332531972288E-15, ref b0, ref b1, ref b2);
|
---|
423 | chebiterationshichi(a, 3.49639695410806959872E-14, ref b0, ref b1, ref b2);
|
---|
424 | chebiterationshichi(a, 6.63406731718911586609E-14, ref b0, ref b1, ref b2);
|
---|
425 | chebiterationshichi(a, -3.71902448093119218395E-13, ref b0, ref b1, ref b2);
|
---|
426 | chebiterationshichi(a, -1.27135418132338309016E-12, ref b0, ref b1, ref b2);
|
---|
427 | chebiterationshichi(a, 2.74851141935315395333E-12, ref b0, ref b1, ref b2);
|
---|
428 | chebiterationshichi(a, 2.33781843985453438400E-11, ref b0, ref b1, ref b2);
|
---|
429 | chebiterationshichi(a, 2.71436006377612442764E-11, ref b0, ref b1, ref b2);
|
---|
430 | chebiterationshichi(a, -2.56600180000355990529E-10, ref b0, ref b1, ref b2);
|
---|
431 | chebiterationshichi(a, -1.61021375163803438552E-9, ref b0, ref b1, ref b2);
|
---|
432 | chebiterationshichi(a, -4.72543064876271773512E-9, ref b0, ref b1, ref b2);
|
---|
433 | chebiterationshichi(a, -3.00095178028681682282E-9, ref b0, ref b1, ref b2);
|
---|
434 | chebiterationshichi(a, 7.79387474390914922337E-8, ref b0, ref b1, ref b2);
|
---|
435 | chebiterationshichi(a, 1.06942765566401507066E-6, ref b0, ref b1, ref b2);
|
---|
436 | chebiterationshichi(a, 1.59503164802313196374E-5, ref b0, ref b1, ref b2);
|
---|
437 | chebiterationshichi(a, 3.49592575153777996871E-4, ref b0, ref b1, ref b2);
|
---|
438 | chebiterationshichi(a, 1.28475387530065247392E-2, ref b0, ref b1, ref b2);
|
---|
439 | chebiterationshichi(a, 1.03665693917934275131E0, ref b0, ref b1, ref b2);
|
---|
440 | c = k*0.5*(b0-b2);
|
---|
441 | }
|
---|
442 | else
|
---|
443 | {
|
---|
444 | if( sg!=0 )
|
---|
445 | {
|
---|
446 | shi = -AP.Math.MaxRealNumber;
|
---|
447 | }
|
---|
448 | else
|
---|
449 | {
|
---|
450 | shi = AP.Math.MaxRealNumber;
|
---|
451 | }
|
---|
452 | chi = AP.Math.MaxRealNumber;
|
---|
453 | return;
|
---|
454 | }
|
---|
455 | }
|
---|
456 | }
|
---|
457 | if( sg!=0 )
|
---|
458 | {
|
---|
459 | s = -s;
|
---|
460 | }
|
---|
461 | shi = s;
|
---|
462 | chi = 0.57721566490153286061+Math.Log(x)+c;
|
---|
463 | }
|
---|
464 |
|
---|
465 |
|
---|
466 | private static void chebiterationshichi(double x,
|
---|
467 | double c,
|
---|
468 | ref double b0,
|
---|
469 | ref double b1,
|
---|
470 | ref double b2)
|
---|
471 | {
|
---|
472 | b2 = b1;
|
---|
473 | b1 = b0;
|
---|
474 | b0 = x*b1-b2+c;
|
---|
475 | }
|
---|
476 | }
|
---|
477 | }
|
---|