File: | libinterp/corefcn/quadcc.cc |
Location: | line 1469, column 28 |
Description: | The left operand of '+' is a garbage value |
1 | /* | |||
2 | ||||
3 | Copyright (C) 2010-2013 Pedro Gonnet | |||
4 | ||||
5 | This file is part of Octave. | |||
6 | ||||
7 | Octave is free software; you can redistribute it and/or modify it | |||
8 | under the terms of the GNU General Public License as published by the | |||
9 | Free Software Foundation; either version 3 of the License, or (at your | |||
10 | option) any later version. | |||
11 | ||||
12 | Octave is distributed in the hope that it will be useful, but WITHOUT | |||
13 | ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or | |||
14 | FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License | |||
15 | for more details. | |||
16 | ||||
17 | You should have received a copy of the GNU General Public License | |||
18 | along with Octave; see the file COPYING. If not, see | |||
19 | <http://www.gnu.org/licenses/>. | |||
20 | ||||
21 | */ | |||
22 | ||||
23 | #ifdef HAVE_CONFIG_H1 | |||
24 | #include <config.h> | |||
25 | #endif | |||
26 | ||||
27 | #include "lo-ieee.h" | |||
28 | #include "parse.h" | |||
29 | #include "variables.h" | |||
30 | ||||
31 | #include "defun.h" | |||
32 | #include "error.h" | |||
33 | #include "oct-obj.h" | |||
34 | #include "utils.h" | |||
35 | ||||
36 | ||||
37 | /* Extended debugging */ | |||
38 | #define DEBUG_QUADCC0 0 | |||
39 | ||||
40 | /* Define the minimum size of the interval heap. */ | |||
41 | #define min_cquad_heapsize200 200 | |||
42 | ||||
43 | ||||
44 | /* Data of a single interval */ | |||
45 | typedef struct | |||
46 | { | |||
47 | double a, b; | |||
48 | double c[64]; | |||
49 | double fx[33]; | |||
50 | double igral, err; | |||
51 | int depth, rdepth, ndiv; | |||
52 | } cquad_ival; | |||
53 | ||||
54 | /* Some constants and matrices that we'll need. */ | |||
55 | ||||
56 | static const double xi[33] = | |||
57 | { | |||
58 | -1., -0.99518472667219688624, -0.98078528040323044912, | |||
59 | -0.95694033573220886493, -0.92387953251128675612, | |||
60 | -0.88192126434835502970, -0.83146961230254523708, | |||
61 | -0.77301045336273696082, -0.70710678118654752440, | |||
62 | -0.63439328416364549822, -0.55557023301960222475, | |||
63 | -0.47139673682599764857, -0.38268343236508977173, | |||
64 | -0.29028467725446236764, -0.19509032201612826785, | |||
65 | -0.098017140329560601995, 0., 0.098017140329560601995, | |||
66 | 0.19509032201612826785, 0.29028467725446236764, 0.38268343236508977173, | |||
67 | 0.47139673682599764857, 0.55557023301960222475, 0.63439328416364549822, | |||
68 | 0.70710678118654752440, 0.77301045336273696082, 0.83146961230254523708, | |||
69 | 0.88192126434835502970, 0.92387953251128675612, 0.95694033573220886493, | |||
70 | 0.98078528040323044912, 0.99518472667219688624, 1. | |||
71 | }; | |||
72 | ||||
73 | static const double bee[68] = | |||
74 | { | |||
75 | 0.00000000000000e+00, 2.28868854108532e-01, 0.00000000000000e+00, | |||
76 | -8.15740215243451e-01, 0.00000000000000e+00, 5.31212715259731e-01, | |||
77 | 0.00000000000000e+00, 1.38538036812454e-02, 0.00000000000000e+00, | |||
78 | 3.74405228908818e-02, 0.00000000000000e+00, 2.12224115039342e-01, | |||
79 | 0.00000000000000e+00, -8.16362644507898e-01, 0.00000000000000e+00, | |||
80 | 5.35648426691481e-01, 0.00000000000000e+00, 1.52417902753662e-03, | |||
81 | 0.00000000000000e+00, 2.63058840550873e-03, 0.00000000000000e+00, | |||
82 | 4.15292106318904e-03, 0.00000000000000e+00, 6.97106011119775e-03, | |||
83 | 0.00000000000000e+00, 1.35535708431058e-02, 0.00000000000000e+00, | |||
84 | 3.52132898424856e-02, 0.00000000000000e+00, 2.06946714741884e-01, | |||
85 | 0.00000000000000e+00, -8.15674251283876e-01, 0.00000000000000e+00, | |||
86 | 5.38841175520580e-01, 0.00000000000000e+00, 1.84909689577590e-04, | |||
87 | 0.00000000000000e+00, 2.90936325007499e-04, 0.00000000000000e+00, | |||
88 | 3.84877750950089e-04, 0.00000000000000e+00, 4.86436656735046e-04, | |||
89 | 0.00000000000000e+00, 6.08688640346879e-04, 0.00000000000000e+00, | |||
90 | 7.66732830740331e-04, 0.00000000000000e+00, 9.82753336104205e-04, | |||
91 | 0.00000000000000e+00, 1.29359957505615e-03, 0.00000000000000e+00, | |||
92 | 1.76616363801885e-03, 0.00000000000000e+00, 2.53323433039089e-03, | |||
93 | 0.00000000000000e+00, 3.88872172121956e-03, 0.00000000000000e+00, | |||
94 | 6.58635106468291e-03, 0.00000000000000e+00, 1.30326736343254e-02, | |||
95 | 0.00000000000000e+00, 3.44353850696714e-02, 0.00000000000000e+00, | |||
96 | 2.05025409531915e-01, 0.00000000000000e+00, -8.14985893995401e-01, | |||
97 | 0.00000000000000e+00, 5.40679930965238e-01 | |||
98 | }; | |||
99 | ||||
100 | static const double Lalpha[33] = | |||
101 | { | |||
102 | 5.77350269189626e-01, 5.16397779494322e-01, 5.07092552837110e-01, | |||
103 | 5.03952630678970e-01, 5.02518907629606e-01, 5.01745206004255e-01, | |||
104 | 5.01280411827603e-01, 5.00979432868120e-01, 5.00773395667191e-01, | |||
105 | 5.00626174321759e-01, 5.00517330712619e-01, 5.00434593736979e-01, | |||
106 | 5.00370233297676e-01, 5.00319182924304e-01, 5.00278009473803e-01, | |||
107 | 5.00244319584578e-01, 5.00216403386025e-01, 5.00193012939056e-01, | |||
108 | 5.00173220168024e-01, 5.00156323280355e-01, 5.00141783641018e-01, | |||
109 | 5.00129182278347e-01, 5.00118189340972e-01, 5.00108542278496e-01, | |||
110 | 5.00100030010004e-01, 5.00092481273333e-01, 5.00085755939229e-01, | |||
111 | 5.00079738458365e-01, 5.00074332862969e-01, 5.00069458915387e-01, | |||
112 | 5.00065049112355e-01, 5.00061046334395e-01, 5.00057401986298e-01 | |||
113 | }; | |||
114 | ||||
115 | static const double Lgamma[33] = | |||
116 | { | |||
117 | 0.0, 0.0, 5.16397779494322e-01, 5.07092552837110e-01, 5.03952630678970e-01, | |||
118 | 5.02518907629606e-01, 5.01745206004255e-01, 5.01280411827603e-01, | |||
119 | 5.00979432868120e-01, 5.00773395667191e-01, 5.00626174321759e-01, | |||
120 | 5.00517330712619e-01, 5.00434593736979e-01, 5.00370233297676e-01, | |||
121 | 5.00319182924304e-01, 5.00278009473803e-01, 5.00244319584578e-01, | |||
122 | 5.00216403386025e-01, 5.00193012939056e-01, 5.00173220168024e-01, | |||
123 | 5.00156323280355e-01, 5.00141783641018e-01, 5.00129182278347e-01, | |||
124 | 5.00118189340972e-01, 5.00108542278496e-01, 5.00100030010003e-01, | |||
125 | 5.00092481273333e-01, 5.00085755939229e-01, 5.00079738458365e-01, | |||
126 | 5.00074332862969e-01, 5.00069458915387e-01, 5.00065049112355e-01, | |||
127 | 5.00061046334395e-01 | |||
128 | }; | |||
129 | ||||
130 | static const double V1inv[5 * 5] = | |||
131 | { | |||
132 | .47140452079103168293e-1, .37712361663282534635, .56568542494923801952, | |||
133 | .37712361663282534635, .47140452079103168293e-1, | |||
134 | -.81649658092772603273e-1, -.46188021535170061160, 0, | |||
135 | .46188021535170061160, .81649658092772603273e-1, .15058465048420853962, | |||
136 | .12046772038736683169, -.54210474174315074262, .12046772038736683169, | |||
137 | .15058465048420853962, -.21380899352993950775, .30237157840738178177, -0., | |||
138 | -.30237157840738178177, .21380899352993950775, .10774960475223581324, | |||
139 | -.21549920950447162648, .21549920950447162648, -.21549920950447162648, | |||
140 | .10774960475223581324 | |||
141 | }; | |||
142 | ||||
143 | static const double V2inv[9 * 9] = | |||
144 | { | |||
145 | .11223917161691230546e-1, .10339219839658349826, .19754094204576565761, | |||
146 | .25577315077753587922, .27835314560994251755, .25577315077753587922, | |||
147 | .19754094204576565761, .10339219839658349826, .11223917161691230546e-1, | |||
148 | -.19440394783993476970e-1, -.16544884625069155470, -.24193725566041460608, | |||
149 | -.16953338808305493604, 0.0, .16953338808305493604, .24193725566041460608, | |||
150 | .16544884625069155470, .19440394783993476970e-1, .26466393115406349388e-1, | |||
151 | .17766815796285469394, .11316664642449611462, -.16306601003711325980, | |||
152 | -.30847037493128779631, -.16306601003711325980, .11316664642449611462, | |||
153 | .17766815796285469394, .26466393115406349388e-1, | |||
154 | -.32395302049990834508e-1, -.15521142532414866547, | |||
155 | .88573492664788602740e-1, .29570405784974857322, 0.0, | |||
156 | -.29570405784974857322, -.88573492664788602740e-1, .15521142532414866547, | |||
157 | .32395302049990834508e-1, .41442155673936851246e-1, | |||
158 | .98186757907405608245e-1, -.23056908429499411784, | |||
159 | -.68047008326360625520e-1, .31797435808002456774, | |||
160 | -.68047008326360625520e-1, -.23056908429499411784, | |||
161 | .98186757907405608245e-1, .41442155673936851246e-1, | |||
162 | -.49981120317798783134e-1, -.24861810572835756217e-1, | |||
163 | .23561326072010832539, -.24472785656448415351, 0.0, .24472785656448415351, | |||
164 | -.23561326072010832539, .24861810572835756217e-1, | |||
165 | .49981120317798783134e-1, .79691635865674781228e-1, | |||
166 | -.95725617891693941833e-1, -.57957553356854386344e-1, | |||
167 | .21164072460540271452, -.27529837844505833514, .21164072460540271452, | |||
168 | -.57957553356854386344e-1, -.95725617891693941833e-1, | |||
169 | .79691635865674781228e-1, | |||
170 | -.10894869830716590913, .20131094491947531782, -.15407672674888869038, | |||
171 | .83385723639789791384e-1, 0.0, -.83385723639789791384e-1, | |||
172 | .15407672674888869038, -.20131094491947531782, .10894869830716590913, | |||
173 | .54581057089643838221e-1, -.10916211417928767644, .10916211417928767644, | |||
174 | -.10916211417928767644, .10916211417928767644, -.10916211417928767644, | |||
175 | .10916211417928767644, -.10916211417928767644, .54581057089643838221e-1 | |||
176 | }; | |||
177 | ||||
178 | static const double V3inv[17 * 17] = | |||
179 | { | |||
180 | .27729677693590098996e-2, .26423663180333065153e-1, | |||
181 | .53374068493933898312e-1, .77007854739523195947e-1, | |||
182 | .98257061072911596869e-1, .11538049741786835604, .12832134344120884559, | |||
183 | .13612785914022865001, .13888293186236181317, .13612785914022865001, | |||
184 | .12832134344120884559, .11538049741786835604, .98257061072911596869e-1, | |||
185 | .77007854739523195947e-1, .53374068493933898312e-1, | |||
186 | .26423663180333065153e-1, .27729677693590098996e-2, | |||
187 | -.48029210642807413690e-2, -.44887724635478800254e-1, | |||
188 | -.85409520147301089416e-1, -.11090267822061423050, -.12033983162705862441, | |||
189 | -.11102786862182788886, -.85054870109799336515e-1, | |||
190 | -.45998467987742225160e-1, 0.0, .45998467987742225160e-1, | |||
191 | .85054870109799336515e-1, .11102786862182788886, .12033983162705862441, | |||
192 | .11090267822061423050, .85409520147301089416e-1, .44887724635478800254e-1, | |||
193 | .48029210642807413690e-2, .62758546879582030087e-2, | |||
194 | .55561297093529155869e-1, | |||
195 | .93281491021051539742e-1, .92320151237493695139e-1, | |||
196 | .55077987469605684531e-1, | |||
197 | -.96998141716497488255e-2, -.80285961895427405567e-1, | |||
198 | -.13496839655913850224, | |||
199 | -.15512521776684524331, -.13496839655913850224, -.80285961895427405567e-1, | |||
200 | -.96998141716497488255e-2, .55077987469605684531e-1, | |||
201 | .92320151237493695139e-1, .93281491021051539742e-1, | |||
202 | .55561297093529155869e-1, .62758546879582030087e-2, | |||
203 | -.74850969394858555939e-2, -.61751608943839234096e-1, | |||
204 | -.82974150437304275958e-1, -.38437763431942633378e-1, | |||
205 | .45745502025779701366e-1, .12369235652734542162, .14720439712852868239, | |||
206 | .98768034347019704401e-1, 0.0, | |||
207 | -.98768034347019704401e-1, -.14720439712852868239, -.12369235652734542162, | |||
208 | -.45745502025779701366e-1, .38437763431942633378e-1, | |||
209 | .82974150437304275958e-1, .61751608943839234096e-1, | |||
210 | .74850969394858555939e-2, .86710099994384056338e-2, | |||
211 | .64006230103659573344e-1, .58517426396091675690e-1, | |||
212 | -.29743410528985802680e-1, | |||
213 | -.11934127779157114754, -.12686773515361299409, -.30729137153877447035e-1, | |||
214 | .97307836256600731568e-1, .15635811574451401023, .97307836256600731568e-1, | |||
215 | -.30729137153877447035e-1, -.12686773515361299409, -.11934127779157114754, | |||
216 | -.29743410528985802680e-1, .58517426396091675690e-1, | |||
217 | .64006230103659573344e-1, .86710099994384056338e-2, | |||
218 | -.97486395666294840165e-2, -.62995604908060224672e-1, | |||
219 | -.24373234450275529219e-1, .87760984413626872730e-1, | |||
220 | .12205204576993351394, | |||
221 | .16216004196864002088e-1, -.12422320942156845775, -.13682714580929614678, | |||
222 | 0.0, .13682714580929614678, .12422320942156845775, | |||
223 | -.16216004196864002088e-1, -.12205204576993351394, | |||
224 | -.87760984413626872730e-1, .24373234450275529219e-1, | |||
225 | .62995604908060224672e-1, .97486395666294840165e-2, | |||
226 | .10956271233750488468e-1, .58613204255294358939e-1, | |||
227 | -.13306063940736618859e-1, -.11606666444978454399, | |||
228 | -.52059598001115805639e-1, .10868540217796151849, .12594452879014618005, | |||
229 | -.44678658254872910434e-1, -.15617684362128533405, | |||
230 | -.44678658254872910434e-1, .12594452879014618005, .10868540217796151849, | |||
231 | -.52059598001115805639e-1, -.11606666444978454399, | |||
232 | -.13306063940736618859e-1, .58613204255294358939e-1, | |||
233 | .10956271233750488468e-1, -.12098893000863087230e-1, | |||
234 | -.51626244709126208453e-1, .48919433304746979330e-1, | |||
235 | .10467644465949427090, | |||
236 | -.48729879523084673782e-1, -.13668732103524749234, .28190838706814496438e-1, | |||
237 | .15434223333238741600, 0.0, -.15434223333238741600, | |||
238 | -.28190838706814496438e-1, .13668732103524749234, | |||
239 | .48729879523084673782e-1, -.10467644465949427090, | |||
240 | -.48919433304746979330e-1, .51626244709126208453e-1, | |||
241 | .12098893000863087230e-1, .13542668300437944822e-1, | |||
242 | .41712033418258689308e-1, | |||
243 | -.76190463272803434388e-1, -.58303943170068132010e-1, .12158068748245606853, | |||
244 | .42121099930651007882e-1, -.14684425840766337756, | |||
245 | -.16108203535058647043e-1, .15698075850757976092, | |||
246 | -.16108203535058647043e-1, -.14684425840766337756, | |||
247 | .42121099930651007882e-1, .12158068748245606853, | |||
248 | -.58303943170068132010e-1, -.76190463272803434388e-1, | |||
249 | .41712033418258689308e-1, .13542668300437944822e-1, | |||
250 | -.14939634995117694417e-1, -.30047246373341564039e-1, | |||
251 | .91624635082546425678e-1, -.79133374319110026377e-2, | |||
252 | -.12292558212072233355, .90013382617762643524e-1, | |||
253 | .84013717196539593395e-1, -.14813033309980695856, 0.0, | |||
254 | .14813033309980695856, -.84013717196539593395e-1, | |||
255 | -.90013382617762643524e-1, | |||
256 | .12292558212072233355, .79133374319110026377e-2, -.91624635082546425678e-1, | |||
257 | .30047246373341564039e-1, .14939634995117694417e-1, | |||
258 | .16986031342807474208e-1, | |||
259 | .15760203882617033601e-1, -.91494054040950941996e-1, | |||
260 | .70082459207876130806e-1, | |||
261 | .53390713710144539104e-1, -.14340746778352039430, .84048122493418898508e-1, | |||
262 | .72456667788091316868e-1, -.15564535320096811360, | |||
263 | .72456667788091316868e-1, .84048122493418898508e-1, | |||
264 | -.14340746778352039430, .53390713710144539104e-1, | |||
265 | .70082459207876130806e-1, -.91494054040950941996e-1, | |||
266 | .15760203882617033601e-1, | |||
267 | .16986031342807474208e-1, -.18994065631858742028e-1, | |||
268 | -.82901821370405592927e-3, .77239669773015192888e-1, | |||
269 | -.10850735431039424680, .47524484622086496464e-1, | |||
270 | .69148184871588737021e-1, -.14829314646228194928, .11992057742398672066, | |||
271 | 0.0, -.11992057742398672066, .14829314646228194928, | |||
272 | -.69148184871588737021e-1, -.47524484622086496464e-1, | |||
273 | .10850735431039424680, -.77239669773015192888e-1, | |||
274 | .82901821370405592927e-3, .18994065631858742028e-1, | |||
275 | .22761703826371535132e-1, -.17728848711449643358e-1, | |||
276 | -.47496371572480503788e-1, .10659958402328690063, -.11696013966166296514, | |||
277 | .63073750910894244526e-1, .32928881123602721303e-1, | |||
278 | -.12280950532497593683, .15926189077282729505, -.12280950532497593683, | |||
279 | .32928881123602721303e-1, .63073750910894244526e-1, | |||
280 | -.11696013966166296514, .10659958402328690063, -.47496371572480503788e-1, | |||
281 | -.17728848711449643358e-1, .22761703826371535132e-1, | |||
282 | -.26493215276042203434e-1, .35579780856128386192e-1, | |||
283 | .10447309718398935122e-1, -.68616154085314996709e-1, | |||
284 | .11775363082763954214, -.13918901977011837274, .12312819418827395690, | |||
285 | -.72053565748259077905e-1, 0.0, .72053565748259077905e-1, | |||
286 | -.12312819418827395690, .13918901977011837274, -.11775363082763954214, | |||
287 | .68616154085314996709e-1, -.10447309718398935122e-1, | |||
288 | -.35579780856128386192e-1, | |||
289 | .26493215276042203434e-1, .40742523354399706918e-1, | |||
290 | -.73124912999529117195e-1, .49317266444153837821e-1, | |||
291 | -.13686605413876015320e-1, -.28342624942191100464e-1, | |||
292 | .70371855298258216249e-1, -.10600251632853603875, .12981016288391131812, | |||
293 | -.13817029659318161476, .12981016288391131812, -.10600251632853603875, | |||
294 | .70371855298258216249e-1, -.28342624942191100464e-1, | |||
295 | -.13686605413876015320e-1, | |||
296 | .49317266444153837821e-1, -.73124912999529117195e-1, | |||
297 | .40742523354399706918e-1, -.54944368958699908688e-1, | |||
298 | .10777725663147408190, -.10152395581538265428, .91369146312596428468e-1, | |||
299 | -.77703071757424700773e-1, .61050911730999815031e-1, | |||
300 | -.42052599404498348871e-1, .21438229266251454773e-1, 0.0, | |||
301 | -.21438229266251454773e-1, .42052599404498348871e-1, | |||
302 | -.61050911730999815031e-1, .77703071757424700773e-1, | |||
303 | -.91369146312596428468e-1, | |||
304 | .10152395581538265428, -.10777725663147408190, .54944368958699908688e-1, | |||
305 | .27485608464748840573e-1, -.54971216929497681146e-1, | |||
306 | .54971216929497681146e-1, | |||
307 | -.54971216929497681146e-1, .54971216929497681146e-1, | |||
308 | -.54971216929497681146e-1, .54971216929497681146e-1, | |||
309 | -.54971216929497681146e-1, .54971216929497681146e-1, | |||
310 | -.54971216929497681146e-1, .54971216929497681146e-1, | |||
311 | -.54971216929497681146e-1, .54971216929497681146e-1, | |||
312 | -.54971216929497681146e-1, .54971216929497681146e-1, | |||
313 | -.54971216929497681146e-1, .27485608464748840573e-1 | |||
314 | }; | |||
315 | ||||
316 | static const double V4inv[33 * 33] = | |||
317 | { | |||
318 | .69120897476690862600e-3, .66419939766331555194e-2, | |||
319 | .13600665164323186111e-1, .20122785860913684493e-1, | |||
320 | .26583214101668429944e-1, .32712713318999268739e-1, | |||
321 | .38576221976287138036e-1, .44033030938268925133e-1, | |||
322 | .49092709529622799673e-1, .53657949874312515646e-1, | |||
323 | .57724533144734311859e-1, .61219564530655179096e-1, | |||
324 | .64138907503837875026e-1, .66427905189318792009e-1, | |||
325 | .68088956652280022887e-1, .69083051391555695878e-1, | |||
326 | .69422738116739271449e-1, .69083051391555695878e-1, | |||
327 | .68088956652280022887e-1, .66427905189318792009e-1, | |||
328 | .64138907503837875026e-1, .61219564530655179096e-1, | |||
329 | .57724533144734311859e-1, .53657949874312515646e-1, | |||
330 | .49092709529622799673e-1, .44033030938268925133e-1, | |||
331 | .38576221976287138036e-1, .32712713318999268739e-1, | |||
332 | .26583214101668429944e-1, .20122785860913684493e-1, | |||
333 | .13600665164323186111e-1, .66419939766331555194e-2, | |||
334 | .69120897476690862600e-3, -.11972090629438798134e-2, | |||
335 | -.11448874821643225573e-1, -.23104401104002905904e-1, | |||
336 | -.33352899418646530133e-1, -.42538626424075425908e-1, | |||
337 | -.49969730733911825941e-1, -.55555454015360728353e-1, | |||
338 | -.58955533624852604918e-1, -.60126044219122513907e-1, | |||
339 | -.58959430451175833624e-1, -.55546925396227130606e-1, | |||
340 | -.49984739749347973762e-1, -.42513009141170294365e-1, | |||
341 | -.33399140950669746346e-1, -.23007690803851790829e-1, | |||
342 | -.11728275717520066169e-1, 0.0, .11728275717520066169e-1, | |||
343 | .23007690803851790829e-1, .33399140950669746346e-1, | |||
344 | .42513009141170294365e-1, .49984739749347973762e-1, | |||
345 | .55546925396227130606e-1, .58959430451175833624e-1, | |||
346 | .60126044219122513907e-1, .58955533624852604918e-1, | |||
347 | .55555454015360728353e-1, .49969730733911825941e-1, | |||
348 | .42538626424075425908e-1, .33352899418646530133e-1, | |||
349 | .23104401104002905904e-1, .11448874821643225573e-1, | |||
350 | .11972090629438798134e-2, .15501585012936019146e-2, | |||
351 | .14628781502199620482e-1, .28684915921474815271e-1, | |||
352 | .39299396074628048026e-1, .46393418975496284204e-1, | |||
353 | .48756902531094699526e-1, .46331333488337494692e-1, | |||
354 | .39012645376980228775e-1, .27452795421085791153e-1, | |||
355 | .12430953621169863781e-1, -.47682978056024928800e-2, | |||
356 | -.22825828045428973853e-1, | |||
357 | -.40195512090720278312e-1, -.55503004262826221955e-1, | |||
358 | -.67424537752827046308e-1, -.75020199300113606452e-1, | |||
359 | -.77607844312483656131e-1, -.75020199300113606452e-1, | |||
360 | -.67424537752827046308e-1, -.55503004262826221955e-1, | |||
361 | -.40195512090720278312e-1, -.22825828045428973853e-1, | |||
362 | -.47682978056024928800e-2, .12430953621169863781e-1, | |||
363 | .27452795421085791153e-1, .39012645376980228775e-1, | |||
364 | .46331333488337494692e-1, .48756902531094699526e-1, | |||
365 | .46393418975496284204e-1, .39299396074628048026e-1, | |||
366 | .28684915921474815271e-1, .14628781502199620482e-1, | |||
367 | .15501585012936019146e-2, -.18377757558949194214e-2, | |||
368 | -.17050470050949761565e-1, -.31952119564923250836e-1, | |||
369 | -.40197423449026348155e-1, | |||
370 | -.41205649520281371624e-1, -.33909965817492272248e-1, | |||
371 | -.19393664422115332144e-1, .56661049630886784692e-3, | |||
372 | .22948272173686561721e-1, .44489719570904738207e-1, | |||
373 | .61790363672287920596e-1, .72121014727028013894e-1, | |||
374 | .73627151185287858579e-1, .65784665375961398923e-1, | |||
375 | .49369676372333667559e-1, .26444326317059715065e-1, 0.0, | |||
376 | -.26444326317059715065e-1, -.49369676372333667559e-1, | |||
377 | -.65784665375961398923e-1, -.73627151185287858579e-1, | |||
378 | -.72121014727028013894e-1, -.61790363672287920596e-1, | |||
379 | -.44489719570904738207e-1, -.22948272173686561721e-1, | |||
380 | -.56661049630886784692e-3, .19393664422115332144e-1, | |||
381 | .33909965817492272248e-1, .41205649520281371624e-1, | |||
382 | .40197423449026348155e-1, .31952119564923250836e-1, | |||
383 | .17050470050949761565e-1, .18377757558949194214e-2, | |||
384 | .20942714740729767769e-2, .18935902405146518232e-1, | |||
385 | .33335840852491735126e-1, .36770680999102286065e-1, | |||
386 | .28873194534132768509e-1, .10267303017729535513e-1, | |||
387 | -.14607738306201572890e-1, -.40139568545572305818e-1, | |||
388 | -.59808326733858291561e-1, -.68528358823372627506e-1, | |||
389 | -.63306535387619244879e-1, -.44508601817574921056e-1, | |||
390 | -.15449116105605395357e-1, .17941083795006546367e-1, | |||
391 | .48747356011657242123e-1, .70329553984201665523e-1, | |||
392 | .78106117292526169663e-1, .70329553984201665523e-1, | |||
393 | .48747356011657242123e-1, .17941083795006546367e-1, | |||
394 | -.15449116105605395357e-1, -.44508601817574921056e-1, | |||
395 | -.63306535387619244879e-1, -.68528358823372627506e-1, | |||
396 | -.59808326733858291561e-1, | |||
397 | -.40139568545572305818e-1, -.14607738306201572890e-1, | |||
398 | .10267303017729535513e-1, .28873194534132768509e-1, | |||
399 | .36770680999102286065e-1, .33335840852491735126e-1, | |||
400 | .18935902405146518232e-1, .20942714740729767769e-2, | |||
401 | -.23245285491878278419e-2, -.20401404737639389919e-1, | |||
402 | -.33019548231022514097e-1, -.29709828426463720091e-1, | |||
403 | -.11760070922697422156e-1, .15987584743850393793e-1, | |||
404 | .43619012891472813485e-1, .61177322409671487721e-1, | |||
405 | .61144030218486655594e-1, | |||
406 | .41895377620089086167e-1, .80232011820644308033e-2, | |||
407 | -.30574701186675900915e-1, | |||
408 | -.62072243008844865848e-1, -.76336186183574765586e-1, | |||
409 | -.68435466095345537115e-1, -.40237669208466966207e-1, 0.0, | |||
410 | .40237669208466966207e-1, .68435466095345537115e-1, | |||
411 | .76336186183574765586e-1, .62072243008844865848e-1, | |||
412 | .30574701186675900915e-1, -.80232011820644308033e-2, | |||
413 | -.41895377620089086167e-1, -.61144030218486655594e-1, | |||
414 | -.61177322409671487721e-1, -.43619012891472813485e-1, | |||
415 | -.15987584743850393793e-1, .11760070922697422156e-1, | |||
416 | .29709828426463720091e-1, .33019548231022514097e-1, | |||
417 | .20401404737639389919e-1, .23245285491878278419e-2, | |||
418 | .25451717261579269307e-2, .21480418595666878775e-1, | |||
419 | .31177212469293007998e-1, .19816333607013379373e-1, | |||
420 | -.72439496274458793681e-2, -.38404203906598342397e-1, | |||
421 | -.57633632255322221046e-1, -.54070547403585392952e-1, | |||
422 | -.26249823354368866005e-1, .15643058212336881516e-1, | |||
423 | .54539832735118677194e-1, .73283028002473989724e-1, | |||
424 | .62835303524135936213e-1, .26175977027801048141e-1, | |||
425 | -.22193636309998606610e-1, -.62597049956093311234e-1, | |||
426 | -.78206986173170212505e-1, -.62597049956093311234e-1, | |||
427 | -.22193636309998606610e-1, .26175977027801048141e-1, | |||
428 | .62835303524135936213e-1, | |||
429 | .73283028002473989724e-1, .54539832735118677194e-1, | |||
430 | .15643058212336881516e-1, | |||
431 | -.26249823354368866005e-1, -.54070547403585392952e-1, | |||
432 | -.57633632255322221046e-1, -.38404203906598342397e-1, | |||
433 | -.72439496274458793681e-2, .19816333607013379373e-1, | |||
434 | .31177212469293007998e-1, .21480418595666878775e-1, | |||
435 | .25451717261579269307e-2, -.27506573922483820005e-2, | |||
436 | -.22224442095099251870e-1, -.27949927254215773020e-1, | |||
437 | -.80918481053370034987e-2, .25121859354449306916e-1, | |||
438 | .51563535009373061074e-1, .51936965107145960512e-1, | |||
439 | .22146626648171527753e-1, | |||
440 | -.24172689882103382748e-1, -.61731229104853568296e-1, | |||
441 | -.68477262429344201201e-1, -.38311232728303704742e-1, | |||
442 | .14160578713659552679e-1, .61248813427564184033e-1, | |||
443 | .77136328841293031805e-1, .52514801765183697988e-1, 0.0, | |||
444 | -.52514801765183697988e-1, -.77136328841293031805e-1, | |||
445 | -.61248813427564184033e-1, -.14160578713659552679e-1, | |||
446 | .38311232728303704742e-1, | |||
447 | .68477262429344201201e-1, .61731229104853568296e-1, | |||
448 | .24172689882103382748e-1, | |||
449 | -.22146626648171527753e-1, -.51936965107145960512e-1, | |||
450 | -.51563535009373061074e-1, -.25121859354449306916e-1, | |||
451 | .80918481053370034987e-2, .27949927254215773020e-1, | |||
452 | .22224442095099251870e-1, .27506573922483820005e-2, | |||
453 | .29562461131654311467e-2, .22630271480554450613e-1, | |||
454 | .23547399831373800971e-1, -.43964593440902476642e-2, | |||
455 | -.39055315767504970597e-1, -.52369643937940066804e-1, | |||
456 | -.28506131614971613422e-1, .19906048093338832322e-1, | |||
457 | .60408880866392420279e-1, .62493397473656883090e-1, | |||
458 | .21391278377641297859e-1, -.37302864786623254746e-1, | |||
459 | -.73665127933539496872e-1, -.61706142476854010202e-1, | |||
460 | -.78065168882546327888e-2, .52335307373945544428e-1, | |||
461 | .78278746279419264777e-1, .52335307373945544428e-1, | |||
462 | -.78065168882546327888e-2, -.61706142476854010202e-1, | |||
463 | -.73665127933539496872e-1, -.37302864786623254746e-1, | |||
464 | .21391278377641297859e-1, .62493397473656883090e-1, | |||
465 | .60408880866392420279e-1, .19906048093338832322e-1, | |||
466 | -.28506131614971613422e-1, -.52369643937940066804e-1, | |||
467 | -.39055315767504970597e-1, -.43964593440902476642e-2, | |||
468 | .23547399831373800971e-1, .22630271480554450613e-1, | |||
469 | .29562461131654311467e-2, -.31515718415504761303e-2, | |||
470 | -.22739451096655080673e-1, -.18157123602272119779e-1, | |||
471 | .16496480897167303621e-1, .46921166788569301124e-1, | |||
472 | .40644395739978416354e-1, -.46275803430732216900e-2, | |||
473 | -.52883375891308909486e-1, -.61116483226324111734e-1, | |||
474 | -.17411698764545629853e-1, .44773430013166822765e-1, | |||
475 | .73441577962383869198e-1, .42127368371995472815e-1, | |||
476 | -.25504645957196772465e-1, -.74126818045972742488e-1, | |||
477 | -.62780077864719287317e-1, 0.0, .62780077864719287317e-1, | |||
478 | .74126818045972742488e-1, .25504645957196772465e-1, | |||
479 | -.42127368371995472815e-1, -.73441577962383869198e-1, | |||
480 | -.44773430013166822765e-1, .17411698764545629853e-1, | |||
481 | .61116483226324111734e-1, .52883375891308909486e-1, | |||
482 | .46275803430732216900e-2, -.40644395739978416354e-1, | |||
483 | -.46921166788569301124e-1, -.16496480897167303621e-1, | |||
484 | .18157123602272119779e-1, .22739451096655080673e-1, | |||
485 | .31515718415504761303e-2, .33536559294882188208e-2, | |||
486 | .22535348942792006185e-1, | |||
487 | .12048629300953560767e-1, -.27166076791299493403e-1, | |||
488 | -.47492745604230978367e-1, -.19246623430993153174e-1, | |||
489 | .36231297307556299322e-1, .61713617181636122004e-1, | |||
490 | .25928029734266134490e-1, -.40478700752883602818e-1, | |||
491 | -.71053889866326412049e-1, -.31870824482961751482e-1, | |||
492 | .41515251100219081281e-1, .76481960760098381651e-1, | |||
493 | .36726509155999912440e-1, -.40090067032627055969e-1, | |||
494 | -.78270742903374539397e-1, -.40090067032627055969e-1, | |||
495 | .36726509155999912440e-1, .76481960760098381651e-1, | |||
496 | .41515251100219081281e-1, -.31870824482961751482e-1, | |||
497 | -.71053889866326412049e-1, -.40478700752883602818e-1, | |||
498 | .25928029734266134490e-1, .61713617181636122004e-1, | |||
499 | .36231297307556299322e-1, -.19246623430993153174e-1, | |||
500 | -.47492745604230978367e-1, -.27166076791299493403e-1, | |||
501 | .12048629300953560767e-1, .22535348942792006185e-1, | |||
502 | .33536559294882188208e-2, | |||
503 | -.35481220456925318865e-2, -.22062913693073191150e-1, | |||
504 | -.54487362861834144999e-2, .35438821865804087489e-1, | |||
505 | .40733077820527411302e-1, -.67403098138950720914e-2, | |||
506 | -.55559584405239171054e-1, -.42417050790865158745e-1, | |||
507 | .24499901971884704925e-1, .68721232891705409302e-1, | |||
508 | .34086082787461126592e-1, -.43441000373118474002e-1, | |||
509 | -.73878085292669148950e-1, -.18846995664706657127e-1, | |||
510 | .59827776178286834498e-1, .70644634584085901794e-1, 0.0, | |||
511 | -.70644634584085901794e-1, -.59827776178286834498e-1, | |||
512 | .18846995664706657127e-1, .73878085292669148950e-1, | |||
513 | .43441000373118474002e-1, -.34086082787461126592e-1, | |||
514 | -.68721232891705409302e-1, -.24499901971884704925e-1, | |||
515 | .42417050790865158745e-1, .55559584405239171054e-1, | |||
516 | .67403098138950720914e-2, -.40733077820527411302e-1, | |||
517 | -.35438821865804087489e-1, .54487362861834144999e-2, | |||
518 | .22062913693073191150e-1, .35481220456925318865e-2, | |||
519 | .37554176816665075631e-2, .21297045781589919482e-1, | |||
520 | -.13327293083183431816e-2, | |||
521 | -.40635299172764596484e-1, -.27659860508374175359e-1, | |||
522 | .31089232744083445986e-1, .56113781541334176109e-1, | |||
523 | .37577840643257763400e-2, -.60511227350664590865e-1, | |||
524 | -.46670556446129053853e-1, .33263195878575888247e-1, | |||
525 | .72757324720645228775e-1, .15011712351692283635e-1, | |||
526 | -.65601212994924119078e-1, -.60016855838843789772e-1, | |||
527 | .26220858553188665966e-1, .78322776605833552980e-1, | |||
528 | .26220858553188665966e-1, -.60016855838843789772e-1, | |||
529 | -.65601212994924119078e-1, | |||
530 | .15011712351692283635e-1, .72757324720645228775e-1, | |||
531 | .33263195878575888247e-1, | |||
532 | -.46670556446129053853e-1, -.60511227350664590865e-1, | |||
533 | .37577840643257763400e-2, .56113781541334176109e-1, | |||
534 | .31089232744083445986e-1, -.27659860508374175359e-1, | |||
535 | -.40635299172764596484e-1, -.13327293083183431816e-2, | |||
536 | .21297045781589919482e-1, .37554176816665075631e-2, | |||
537 | -.39566995305720591229e-2, -.20291873414438919995e-1, | |||
538 | .80617453830770930551e-2, .42270189157016547906e-1, | |||
539 | .10332624526759093004e-1, -.48054759547616142024e-1, | |||
540 | -.37678032941171643972e-1, | |||
541 | .36617192625732482394e-1, .61009425973424865714e-1, | |||
542 | -.95589113168026591466e-2, | |||
543 | -.71023202645076922361e-1, -.25097788086808784456e-1, | |||
544 | .62406621963267050244e-1, .56907293171100693511e-1, | |||
545 | -.36435383083882206257e-1, -.75790105119208756348e-1, 0.0, | |||
546 | .75790105119208756348e-1, .36435383083882206257e-1, | |||
547 | -.56907293171100693511e-1, -.62406621963267050244e-1, | |||
548 | .25097788086808784456e-1, .71023202645076922361e-1, | |||
549 | .95589113168026591466e-2, | |||
550 | -.61009425973424865714e-1, -.36617192625732482394e-1, | |||
551 | .37678032941171643972e-1, .48054759547616142024e-1, | |||
552 | -.10332624526759093004e-1, -.42270189157016547906e-1, | |||
553 | -.80617453830770930551e-2, .20291873414438919995e-1, | |||
554 | .39566995305720591229e-2, .41776092289182138591e-2, | |||
555 | .19013221163904414395e-1, -.14420609729849899876e-1, | |||
556 | -.40259160586844441220e-1, .86327811113710831649e-2, | |||
557 | .53564430703021034399e-1, .65469185402150431933e-2, | |||
558 | -.60383116311280629856e-1, | |||
559 | -.25657793784058876939e-1, .58745680576829226900e-1, | |||
560 | .45649937869034420296e-1, | |||
561 | -.49167932056844167772e-1, -.62696614328552187977e-1, | |||
562 | .32540234556426699997e-1, .74280410383464269758e-1, | |||
563 | -.11425672633410999870e-1, -.78280649404686404903e-1, | |||
564 | -.11425672633410999870e-1, .74280410383464269758e-1, | |||
565 | .32540234556426699997e-1, -.62696614328552187977e-1, | |||
566 | -.49167932056844167772e-1, .45649937869034420296e-1, | |||
567 | .58745680576829226900e-1, -.25657793784058876939e-1, | |||
568 | -.60383116311280629856e-1, .65469185402150431933e-2, | |||
569 | .53564430703021034399e-1, | |||
570 | .86327811113710831649e-2, -.40259160586844441220e-1, | |||
571 | -.14420609729849899876e-1, .19013221163904414395e-1, | |||
572 | .41776092289182138591e-2, -.43935502082478059199e-2, | |||
573 | -.17528761237509401631e-1, .20208915249153872535e-1, | |||
574 | .34734743119040669109e-1, -.26275910172353637955e-1, | |||
575 | -.46368003346018878786e-1, | |||
576 | .26800056330709381025e-1, .56681476464606609921e-1, | |||
577 | -.24749011438127255898e-1, | |||
578 | -.64934612189056658992e-1, .20333742247679279535e-1, | |||
579 | .71429299070059318651e-1, | |||
580 | -.14452513210428671266e-1, -.75793341281736586582e-1, | |||
581 | .74717094137184935270e-2, .78034921554757317374e-1, 0.0, | |||
582 | -.78034921554757317374e-1, -.74717094137184935270e-2, | |||
583 | .75793341281736586582e-1, .14452513210428671266e-1, | |||
584 | -.71429299070059318651e-1, -.20333742247679279535e-1, | |||
585 | .64934612189056658992e-1, .24749011438127255898e-1, | |||
586 | -.56681476464606609921e-1, | |||
587 | -.26800056330709381025e-1, .46368003346018878786e-1, | |||
588 | .26275910172353637955e-1, | |||
589 | -.34734743119040669109e-1, -.20208915249153872535e-1, | |||
590 | .17528761237509401631e-1, .43935502082478059199e-2, | |||
591 | .46379089482818671473e-2, .15791188144791287229e-1, | |||
592 | -.25134290048737455284e-1, -.26249795071946841205e-1, | |||
593 | .39960457575789924651e-1, .28111892450146525404e-1, | |||
594 | -.51026476400767918226e-1, | |||
595 | -.27266747278681831364e-1, .60708796647861610865e-1, | |||
596 | .23532306960642115854e-1, | |||
597 | -.68169639871532441111e-1, -.18204924701958312032e-1, | |||
598 | .73822890510656128485e-1, .11373392486424717019e-1, | |||
599 | -.77133324017644609416e-1, -.39295877480342619961e-2, | |||
600 | .78351902829418987960e-1, -.39295877480342619961e-2, | |||
601 | -.77133324017644609416e-1, .11373392486424717019e-1, | |||
602 | .73822890510656128485e-1, -.18204924701958312032e-1, | |||
603 | -.68169639871532441111e-1, .23532306960642115854e-1, | |||
604 | .60708796647861610865e-1, -.27266747278681831364e-1, | |||
605 | -.51026476400767918226e-1, .28111892450146525404e-1, | |||
606 | .39960457575789924651e-1, -.26249795071946841205e-1, | |||
607 | -.25134290048737455284e-1, .15791188144791287229e-1, | |||
608 | .46379089482818671473e-2, -.48780095920069827068e-2, | |||
609 | -.13886961667516983541e-1, .29071311049368895844e-1, | |||
610 | .15480559452075811600e-1, -.47527977686242313065e-1, | |||
611 | -.31929089844361042178e-2, .58015667638415922967e-1, | |||
612 | -.14547915466597622925e-1, -.61067668299848923244e-1, | |||
613 | .35093678009090186851e-1, .55378399159800654657e-1, | |||
614 | -.54277226474891610385e-1, -.42023830782434076509e-1, | |||
615 | .69197384645944912066e-1, .22610783557709586445e-1, | |||
616 | -.77269275900637030185e-1, 0.0, .77269275900637030185e-1, | |||
617 | -.22610783557709586445e-1, | |||
618 | -.69197384645944912066e-1, .42023830782434076509e-1, | |||
619 | .54277226474891610385e-1, | |||
620 | -.55378399159800654657e-1, -.35093678009090186851e-1, | |||
621 | .61067668299848923244e-1, .14547915466597622925e-1, | |||
622 | -.58015667638415922967e-1, .31929089844361042178e-2, | |||
623 | .47527977686242313065e-1, -.15480559452075811600e-1, | |||
624 | -.29071311049368895844e-1, .13886961667516983541e-1, | |||
625 | .48780095920069827068e-2, .51591759101720291381e-2, | |||
626 | .11747497650231330965e-1, -.31777863364694653331e-1, | |||
627 | -.34555825499804605557e-2, .47914131921157015198e-1, | |||
628 | -.22573685920142225247e-1, -.45320344390022666738e-1, | |||
629 | .49660630547172186418e-1, .25707858143963615736e-1, | |||
630 | -.68132707341917233933e-1, .67534860185243140399e-2, | |||
631 | .69268150370037450063e-1, -.41585011920451477177e-1, | |||
632 | -.51622397460510041271e-1, .68408139576363036148e-1, | |||
633 | .18981259024768933323e-1, -.78265472429342305554e-1, | |||
634 | .18981259024768933323e-1, .68408139576363036148e-1, | |||
635 | -.51622397460510041271e-1, | |||
636 | -.41585011920451477177e-1, .69268150370037450063e-1, | |||
637 | .67534860185243140399e-2, | |||
638 | -.68132707341917233933e-1, .25707858143963615736e-1, | |||
639 | .49660630547172186418e-1, | |||
640 | -.45320344390022666738e-1, -.22573685920142225247e-1, | |||
641 | .47914131921157015198e-1, -.34555825499804605557e-2, | |||
642 | -.31777863364694653331e-1, .11747497650231330965e-1, | |||
643 | .51591759101720291381e-2, -.54365757412741340377e-2, | |||
644 | -.94862516619529080191e-2, .33240472093448190877e-1, | |||
645 | -.88698898099681552229e-2, | |||
646 | -.40973252097216337576e-1, .42995673349795657065e-1, | |||
647 | .17320914507876958783e-1, | |||
648 | -.62201292691914856803e-1, .24726274174637346693e-1, | |||
649 | .51320859246515407288e-1, | |||
650 | -.62882063373810501763e-1, -.11003569131725622672e-1, | |||
651 | .73842261324108943465e-1, -.39240120294802923208e-1, | |||
652 | -.49293966443941122807e-1, .73552644778818223475e-1, 0.0, | |||
653 | -.73552644778818223475e-1, .49293966443941122807e-1, | |||
654 | .39240120294802923208e-1, -.73842261324108943465e-1, | |||
655 | .11003569131725622672e-1, .62882063373810501763e-1, | |||
656 | -.51320859246515407288e-1, | |||
657 | -.24726274174637346693e-1, .62201292691914856803e-1, | |||
658 | -.17320914507876958783e-1, -.42995673349795657065e-1, | |||
659 | .40973252097216337576e-1, .88698898099681552229e-2, | |||
660 | -.33240472093448190877e-1, .94862516619529080191e-2, | |||
661 | .54365757412741340377e-2, .57750194549356126240e-2, | |||
662 | .69981166020044116791e-2, -.33274982140403110792e-1, | |||
663 | .20297071020698356116e-1, .27898517839646066582e-1, | |||
664 | -.53368678853282030262e-1, .16656482990394548343e-1, | |||
665 | .46342901447260614255e-1, | |||
666 | -.60536796508149003365e-1, .29109107483842596340e-2, | |||
667 | .63224486124385124504e-1, | |||
668 | -.59028872851312033411e-1, -.14783105962696191734e-1, | |||
669 | .74269399241069253865e-1, -.49053677339382384625e-1, | |||
670 | -.33525466624811186739e-1, .78397349622515386647e-1, | |||
671 | -.33525466624811186739e-1, -.49053677339382384625e-1, | |||
672 | .74269399241069253865e-1, -.14783105962696191734e-1, | |||
673 | -.59028872851312033411e-1, | |||
674 | .63224486124385124504e-1, .29109107483842596340e-2, | |||
675 | -.60536796508149003365e-1, | |||
676 | .46342901447260614255e-1, .16656482990394548343e-1, | |||
677 | -.53368678853282030262e-1, | |||
678 | .27898517839646066582e-1, .20297071020698356116e-1, | |||
679 | -.33274982140403110792e-1, | |||
680 | .69981166020044116791e-2, .57750194549356126240e-2, | |||
681 | -.61100308370519200637e-2, -.44383614355738148616e-2, | |||
682 | .32011283412619094811e-1, -.29965011866372897633e-1, | |||
683 | -.10560682331349193348e-1, .51110336443392506342e-1, | |||
684 | -.45012284729681775492e-1, -.94236825555873320102e-2, | |||
685 | .60860695783141264746e-1, | |||
686 | -.55014628647083368926e-1, -.73474782382499482121e-2, | |||
687 | .66640148475243034781e-1, -.62533116045749887988e-1, | |||
688 | -.38650525912400102585e-2, .68429769005837003777e-1, | |||
689 | -.66984505412544901945e-1, 0.0, .66984505412544901945e-1, | |||
690 | -.68429769005837003777e-1, .38650525912400102585e-2, | |||
691 | .62533116045749887988e-1, -.66640148475243034781e-1, | |||
692 | .73474782382499482121e-2, | |||
693 | .55014628647083368926e-1, -.60860695783141264746e-1, | |||
694 | .94236825555873320102e-2, | |||
695 | .45012284729681775492e-1, -.51110336443392506342e-1, | |||
696 | .10560682331349193348e-1, | |||
697 | .29965011866372897633e-1, -.32011283412619094811e-1, | |||
698 | .44383614355738148616e-2, | |||
699 | .61100308370519200637e-2, .65409373892036191538e-2, | |||
700 | .16350101107071157065e-2, -.29301957285983144319e-1, | |||
701 | .36838667173388832579e-1, -.81922703976491586393e-2, | |||
702 | -.36955670021050133434e-1, .58374851095540469865e-1, | |||
703 | -.31977016246946181856e-1, -.25311073698658094646e-1, | |||
704 | .66674413950106952577e-1, | |||
705 | -.54865713324521039571e-1, -.39797027891537985440e-2, | |||
706 | .62830285264808449064e-1, -.72226313251296100676e-1, | |||
707 | .22560232697133353980e-1, .46455784709904033738e-1, | |||
708 | -.78200930751070349956e-1, .46455784709904033738e-1, | |||
709 | .22560232697133353980e-1, -.72226313251296100676e-1, | |||
710 | .62830285264808449064e-1, -.39797027891537985440e-2, | |||
711 | -.54865713324521039571e-1, .66674413950106952577e-1, | |||
712 | -.25311073698658094646e-1, -.31977016246946181856e-1, | |||
713 | .58374851095540469865e-1, -.36955670021050133434e-1, | |||
714 | -.81922703976491586393e-2, .36838667173388832579e-1, | |||
715 | -.29301957285983144319e-1, .16350101107071157065e-2, | |||
716 | .65409373892036191538e-2, -.69686180931868703196e-2, | |||
717 | .11849538727632789870e-2, .25452286414610537766e-1, | |||
718 | -.40522480651713943230e-1, .25694679053362813183e-1, | |||
719 | .14057118113748390637e-1, -.52037614725803488893e-1, | |||
720 | .58849342223684035589e-1, | |||
721 | -.25075229077361409271e-1, -.29559771094034181083e-1, | |||
722 | .68296746944165720199e-1, -.62890462146423984955e-1, | |||
723 | .14457636466274596445e-1, .45787612031322361496e-1, | |||
724 | -.77231759014655809742e-1, .57881203613910543657e-1, 0.0, | |||
725 | -.57881203613910543657e-1, .77231759014655809742e-1, | |||
726 | -.45787612031322361496e-1, -.14457636466274596445e-1, | |||
727 | .62890462146423984955e-1, | |||
728 | -.68296746944165720199e-1, .29559771094034181083e-1, | |||
729 | .25075229077361409271e-1, | |||
730 | -.58849342223684035589e-1, .52037614725803488893e-1, | |||
731 | -.14057118113748390637e-1, -.25694679053362813183e-1, | |||
732 | .40522480651713943230e-1, -.25452286414610537766e-1, | |||
733 | -.11849538727632789870e-2, .69686180931868703196e-2, | |||
734 | .75611653617520254845e-2, -.43290610418608409141e-2, | |||
735 | -.20277062025115566914e-1, | |||
736 | .40362947027704828926e-1, -.38938808024132120254e-1, | |||
737 | .11831186195916702262e-1, | |||
738 | .28476667401744525357e-1, -.59320969056617684621e-1, | |||
739 | .61101629747436200186e-1, | |||
740 | -.29514834848355389223e-1, -.20668001885001084821e-1, | |||
741 | .62923592802445122793e-1, -.73558456263588833115e-1, | |||
742 | .45314556330160999776e-1, .79031645918426015574e-2, | |||
743 | -.58136953576334689357e-1, .78538474524006405758e-1, | |||
744 | -.58136953576334689357e-1, .79031645918426015574e-2, | |||
745 | .45314556330160999776e-1, -.73558456263588833115e-1, | |||
746 | .62923592802445122793e-1, -.20668001885001084821e-1, | |||
747 | -.29514834848355389223e-1, .61101629747436200186e-1, | |||
748 | -.59320969056617684621e-1, .28476667401744525357e-1, | |||
749 | .11831186195916702262e-1, -.38938808024132120254e-1, | |||
750 | .40362947027704828926e-1, -.20277062025115566914e-1, | |||
751 | -.43290610418608409141e-2, .75611653617520254845e-2, | |||
752 | -.81505692478987769484e-2, .74297333588288568430e-2, | |||
753 | .14314212513540223314e-1, -.36711242251332751607e-1, | |||
754 | .46240027755503814626e-1, -.34921532671769023773e-1, | |||
755 | .46930051972353714773e-2, | |||
756 | .32842770336385381562e-1, -.61317813706529588466e-1, | |||
757 | .67000809902468893103e-1, | |||
758 | -.45337449655535622885e-1, .35794459576271920867e-2, | |||
759 | .41830061526027213385e-1, | |||
760 | -.72091371931944711708e-1, .74150028530317793195e-1, | |||
761 | -.46487632538609942002e-1, 0.0, .46487632538609942002e-1, | |||
762 | -.74150028530317793195e-1, .72091371931944711708e-1, | |||
763 | -.41830061526027213385e-1, -.35794459576271920867e-2, | |||
764 | .45337449655535622885e-1, -.67000809902468893103e-1, | |||
765 | .61317813706529588466e-1, -.32842770336385381562e-1, | |||
766 | -.46930051972353714773e-2, .34921532671769023773e-1, | |||
767 | -.46240027755503814626e-1, .36711242251332751607e-1, | |||
768 | -.14314212513540223314e-1, -.74297333588288568430e-2, | |||
769 | .81505692478987769484e-2, .90693182942442189743e-2, | |||
770 | -.11121000903959576737e-1, -.71308296141317458546e-2, | |||
771 | .29219439765986671645e-1, -.45820286629778129593e-1, | |||
772 | .49088381175879124421e-1, -.35614888785023038938e-1, | |||
773 | .78906970900092777895e-2, | |||
774 | .26262843038404929480e-1, -.56143674270125757857e-1, | |||
775 | .71700220472378350694e-1, | |||
776 | -.66963544500697307945e-1, .42215091779892228883e-1, | |||
777 | -.41338867413966866997e-2, -.36164891772995367321e-1, | |||
778 | .66584367783847858225e-1, -.77874712365070098328e-1, | |||
779 | .66584367783847858225e-1, -.36164891772995367321e-1, | |||
780 | -.41338867413966866997e-2, .42215091779892228883e-1, | |||
781 | -.66963544500697307945e-1, | |||
782 | .71700220472378350694e-1, -.56143674270125757857e-1, | |||
783 | .26262843038404929480e-1, | |||
784 | .78906970900092777895e-2, -.35614888785023038938e-1, | |||
785 | .49088381175879124421e-1, | |||
786 | -.45820286629778129593e-1, .29219439765986671645e-1, | |||
787 | -.71308296141317458546e-2, -.11121000903959576737e-1, | |||
788 | .90693182942442189743e-2, -.99848472706332791043e-2, | |||
789 | .14701271465939718856e-1, -.32917820356048383366e-3, | |||
790 | -.19201195309873585230e-1, .38409681836626963278e-1, | |||
791 | -.51647324405878909521e-1, .54522171113149311354e-1, | |||
792 | -.45040302741689006270e-1, .24183738595685990149e-1, | |||
793 | .42204134165479735097e-2, -.34317295181348742251e-1, | |||
794 | .59542472465494579941e-1, -.74135115907618101263e-1, | |||
795 | .74491937840566532596e-1, -.60042604725161994304e-1, | |||
796 | .33437677409000083169e-1, 0.0, | |||
797 | -.33437677409000083169e-1, .60042604725161994304e-1, | |||
798 | -.74491937840566532596e-1, .74135115907618101263e-1, | |||
799 | -.59542472465494579941e-1, .34317295181348742251e-1, | |||
800 | -.42204134165479735097e-2, -.24183738595685990149e-1, | |||
801 | .45040302741689006270e-1, -.54522171113149311354e-1, | |||
802 | .51647324405878909521e-1, -.38409681836626963278e-1, | |||
803 | .19201195309873585230e-1, .32917820356048383366e-3, | |||
804 | -.14701271465939718856e-1, .99848472706332791043e-2, | |||
805 | .11775579274769383373e-1, -.19892153937316935880e-1, | |||
806 | .95335114477449041055e-2, .57661528440359081617e-2, | |||
807 | -.23382690532380910781e-1, .40237257037170725321e-1, | |||
808 | -.53280289903551636474e-1, .59974361806023689068e-1, | |||
809 | -.58701684061992853224e-1, .49033407111597129616e-1, | |||
810 | -.31818835267847249219e-1, .90800541261162098886e-2, | |||
811 | .16272906819312603838e-1, -.40863896581186229487e-1, | |||
812 | .61346046297517367703e-1, | |||
813 | -.74896047554167268919e-1, .79632642148310325817e-1, | |||
814 | -.74896047554167268919e-1, .61346046297517367703e-1, | |||
815 | -.40863896581186229487e-1, .16272906819312603838e-1, | |||
816 | .90800541261162098886e-2, -.31818835267847249219e-1, | |||
817 | .49033407111597129616e-1, -.58701684061992853224e-1, | |||
818 | .59974361806023689068e-1, -.53280289903551636474e-1, | |||
819 | .40237257037170725321e-1, -.23382690532380910781e-1, | |||
820 | .57661528440359081617e-2, .95335114477449041055e-2, | |||
821 | -.19892153937316935880e-1, | |||
822 | .11775579274769383373e-1, -.13562702617218467450e-1, | |||
823 | .24885419969649845849e-1, -.18368693901908875583e-1, | |||
824 | .81673147806084084638e-2, .47890591326129587131e-2, | |||
825 | -.19313752945227974024e-1, .34065953398362954708e-1, | |||
826 | -.47667045133463415672e-1, .58820377816690514309e-1, | |||
827 | -.66424139824618415970e-1, | |||
828 | .69667606260856092515e-1, -.68102459384364543253e-1, | |||
829 | .61683024923302547971e-1, | |||
830 | -.50771943476441639136e-1, .36110771847327189215e-1, | |||
831 | -.18758028464284563358e-1, 0.0, .18758028464284563358e-1, | |||
832 | -.36110771847327189215e-1, .50771943476441639136e-1, | |||
833 | -.61683024923302547971e-1, .68102459384364543253e-1, | |||
834 | -.69667606260856092515e-1, .66424139824618415970e-1, | |||
835 | -.58820377816690514309e-1, .47667045133463415672e-1, | |||
836 | -.34065953398362954708e-1, .19313752945227974024e-1, | |||
837 | -.47890591326129587131e-2, -.81673147806084084638e-2, | |||
838 | .18368693901908875583e-1, -.24885419969649845849e-1, | |||
839 | .13562702617218467450e-1, .20576545037980523979e-1, | |||
840 | -.40093155172981004337e-1, .36954083167944054826e-1, | |||
841 | -.31856506837591907746e-1, .24996323181546255126e-1, | |||
842 | -.16637165210473614136e-1, .71002706773325085237e-2, | |||
843 | .32478629093205201133e-2, | |||
844 | -.14009562579050569518e-1, .24771262248780618922e-1, | |||
845 | -.35119395835433647559e-1, .44656290368574753171e-1, | |||
846 | -.53015448339647394161e-1, .59875631995693046782e-1, | |||
847 | -.64973208326045193862e-1, .68112280331082143373e-1, | |||
848 | -.69172215234062186994e-1, .68112280331082143373e-1, | |||
849 | -.64973208326045193862e-1, .59875631995693046782e-1, | |||
850 | -.53015448339647394161e-1, .44656290368574753171e-1, | |||
851 | -.35119395835433647559e-1, .24771262248780618922e-1, | |||
852 | -.14009562579050569518e-1, .32478629093205201133e-2, | |||
853 | .71002706773325085237e-2, -.16637165210473614136e-1, | |||
854 | .24996323181546255126e-1, -.31856506837591907746e-1, | |||
855 | .36954083167944054826e-1, -.40093155172981004337e-1, | |||
856 | .20576545037980523979e-1, -.27584914609096156163e-1, | |||
857 | .54904171411058497973e-1, -.54109756419563083153e-1, | |||
858 | .52794234894345577483e-1, -.50970276026831042415e-1, | |||
859 | .48655445537990983379e-1, | |||
860 | -.45872036510847994332e-1, .42646854695899611372e-1, | |||
861 | -.39010960357087507670e-1, .34999369144476467749e-1, | |||
862 | -.30650714874402762189e-1, .26006877464703437057e-1, | |||
863 | -.21112579608213651273e-1, .16014956068786763273e-1, | |||
864 | -.10763099747751940252e-1, .54075888924374485533e-2, 0.0, | |||
865 | -.54075888924374485533e-2, .10763099747751940252e-1, | |||
866 | -.16014956068786763273e-1, | |||
867 | .21112579608213651273e-1, -.26006877464703437057e-1, | |||
868 | .30650714874402762189e-1, | |||
869 | -.34999369144476467749e-1, .39010960357087507670e-1, | |||
870 | -.42646854695899611372e-1, .45872036510847994332e-1, | |||
871 | -.48655445537990983379e-1, .50970276026831042415e-1, | |||
872 | -.52794234894345577483e-1, .54109756419563083153e-1, | |||
873 | -.54904171411058497973e-1, .27584914609096156163e-1, | |||
874 | .13794141262469565740e-1, -.27588282524939131481e-1, | |||
875 | .27588282524939131481e-1, -.27588282524939131481e-1, | |||
876 | .27588282524939131481e-1, -.27588282524939131481e-1, | |||
877 | .27588282524939131481e-1, | |||
878 | -.27588282524939131481e-1, .27588282524939131481e-1, | |||
879 | -.27588282524939131481e-1, .27588282524939131481e-1, | |||
880 | -.27588282524939131481e-1, .27588282524939131481e-1, | |||
881 | -.27588282524939131481e-1, .27588282524939131481e-1, | |||
882 | -.27588282524939131481e-1, .27588282524939131481e-1, | |||
883 | -.27588282524939131481e-1, .27588282524939131481e-1, | |||
884 | -.27588282524939131481e-1, .27588282524939131481e-1, | |||
885 | -.27588282524939131481e-1, .27588282524939131481e-1, | |||
886 | -.27588282524939131481e-1, .27588282524939131481e-1, | |||
887 | -.27588282524939131481e-1, .27588282524939131481e-1, | |||
888 | -.27588282524939131481e-1, .27588282524939131481e-1, | |||
889 | -.27588282524939131481e-1, .27588282524939131481e-1, | |||
890 | -.27588282524939131481e-1, .13794141262469565740e-1 | |||
891 | }; | |||
892 | ||||
893 | static const double Tleft[33 * 33] = | |||
894 | { | |||
895 | 1., -.86602540378443864678, 0., .33071891388307382381, 0., | |||
896 | -.20728904939721249057, 0., .15128841196122722208, 0., | |||
897 | -.11918864298744029244, 0., .98352013661686631224e-1, 0., | |||
898 | -.83727065404940845733e-1, 0., .72893399403505841203e-1, 0., | |||
899 | -.64544632643375022436e-1, 0., .57913170372415565639e-1, 0., | |||
900 | -.52518242575729562263e-1, 0., .48043311993977520457e-1, 0., | |||
901 | -.44271433659733990243e-1, 0., .41048928022856771981e-1, 0., | |||
902 | -.38263878662008271459e-1, 0., .35832844026365304501e-1, 0., 0., | |||
903 | .50000000000000000000, -.96824583655185422130, .57282196186948000082, | |||
904 | .21650635094610966169, -.35903516540862679125, -.97578093724974971969e-1, | |||
905 | .26203921611325660506, .55792409597991015609e-1, -.20644078533943456204, | |||
906 | -.36172381205961199479e-1, .17035068468874958194, | |||
907 | .25371838001497225980e-1, -.14501953125000000000, | |||
908 | -.18786835250972344757e-1, .12625507130328301066, | |||
909 | .14473795929590520582e-1, -.11179458309419422675, | |||
910 | -.11494434254897626155e-1, .10030855351241635862, | |||
911 | .93498556820544479096e-2, -.90964264465390582629e-1, | |||
912 | -.77546391824364392762e-2, .83213457337452292745e-1, | |||
913 | .65358085945588638605e-2, -.76680372422574234569e-1, | |||
914 | -.55835321940047427169e-2, .71098828931825789428e-1, | |||
915 | .48253327982967591019e-2, -.66274981937248958553e-1, | |||
916 | -.42118078245337801387e-2, .62064306433355646267e-1, | |||
917 | .37083386598903548973e-2, 0., 0., .25000000000000000000, | |||
918 | -.73950997288745200531, .83852549156242113615, -.23175620272173946716, | |||
919 | -.37791833195149451496, .25710129174850522325, .21608307321780204633, | |||
920 | -.22844049245646009157, -.14009503000335388415, .19897685605518413847, | |||
921 | .98264706042471226893e-1, -.17445445004279014046, | |||
922 | -.72761100054958328401e-1, .15463589893742108388, | |||
923 | .56056770591708784481e-1, -.13855313872640495158, | |||
924 | -.44517752443294564781e-1, .12534277657695128850, | |||
925 | .36211835346039665762e-1, -.11434398255136139683, | |||
926 | -.30033588409423828125e-1, .10506705408753910481, | |||
927 | .25313077840725783008e-1, -.97149327637744872155e-1, | |||
928 | -.21624927200393328444e-1, .90319582367202122625e-1, | |||
929 | .18688433567711780666e-1, -.84372291635345108584e-1, | |||
930 | -.16312261561845420752e-1, .79149526894804751586e-1, | |||
931 | .14362333871852474757e-1, 0., 0., 0., .12500000000000000000, | |||
932 | -.49607837082461073572, .82265291131801144317, -.59621200088559103072, | |||
933 | -.80054302859059362371e-1, .42612156697795759420, | |||
934 | -.90098145270865592887e-1, -.29769623255090078484, .13630307904779758221, | |||
935 | .21638835185708931831, -.14600247270306082052, -.16348801804014290453, | |||
936 | .14340708728599057249, .12755243353979286190, -.13661523715071346961, | |||
937 | -.10215585947881057394, .12864248070157166547, .83592528025348693602e-1, | |||
938 | -.12066728689302565222, -.69633728678718053052e-1, .11314245177331919532, | |||
939 | .58882939251410088028e-1, -.10621835858758221487, | |||
940 | -.50432266865187597572e-1, .99916834723527771581e-1, | |||
941 | .43672094283057258509e-1, -.94206380251950852413e-1, | |||
942 | -.38181356812697746418e-1, .89035739656537771225e-1, | |||
943 | .33661934598216332678e-1, 0., 0., 0., 0., .62500000000000000000e-1, | |||
944 | -.31093357409581873586, .67604086414949799246, -.75644205980613611039, | |||
945 | .28990586430124175741, .30648508196770360914, -.35801372616842500052, | |||
946 | -.91326869828709014708e-1, .31127929687500000000, | |||
947 | -.90915752838698393094e-2, -.25637381283965534330, | |||
948 | .57601077850322797594e-1, .21019685709225757945, | |||
949 | -.81244992138514014256e-1, -.17375078516720988858, | |||
950 | .92289437277967051125e-1, .14527351914265391374, | |||
951 | -.96675340792832019889e-1, -.12289485697108543415, | |||
952 | .97448175340011084006e-1, .10511755943298339844, | |||
953 | -.96242247086378239657e-1, -.90822942272780513537e-1, | |||
954 | .93966350452322132384e-1, .79189411876493712558e-1, | |||
955 | -.91139307067989309325e-1, -.69613039934383197265e-1, | |||
956 | .88062491671135767870e-1, .61646331729340817494e-1, 0., 0., 0., 0., 0., | |||
957 | .31250000000000000000e-1, -.18684782411095934408, .50176689760410660236, | |||
958 | -.74784031498626095398, .56472001151566251186, .14842464993721351203e-1, | |||
959 | -.41162920273003120936, .20243071230196532282, .23772054897172750436, | |||
960 | -.24963810923972235950, -.12116179938394678936, .24330535483519110663, | |||
961 | .47903849781124471359e-1, -.22133299683101224293, | |||
962 | -.20542915138527200983e-2, .19653465717678146728, | |||
963 | -.26818172626509178444e-1, -.17319122357631210944, | |||
964 | .45065391411065545445e-1, .15253391395444065941, | |||
965 | -.56543897711725408302e-1, -.13469154928743585367, | |||
966 | .63632471400208840155e-1, .11941684923913523817, | |||
967 | -.67828850207933293098e-1, -.10636309084510652670, | |||
968 | .70095786922999181504e-1, .95187373095150709082e-1, 0., 0., 0., 0., 0., | |||
969 | 0., .15625000000000000000e-1, -.10909562534194485289, | |||
970 | .34842348626527747318, -.64461114561628111443, .69382480527334683659, | |||
971 | -.29551102358528827763, -.25527584713978439819, .38878771718544715394, | |||
972 | -.82956185835347407489e-2, -.31183177761966943912, .12831420840372374767, | |||
973 | .22067618205599434368, -.17569196937129496961, -.14598057000132284135, | |||
974 | .18864406621763419484, .89921002550386645767e-1, -.18571835020187122114, | |||
975 | -.48967672227195481777e-1, .17584685670380332798, | |||
976 | .19267984545067426324e-1, -.16335437520503462738, | |||
977 | .22598055455032407594e-2, .15032800884170631129, | |||
978 | -.17883358353754640871e-1, -.13774837869432209951, | |||
979 | .29227555960587143675e-1, .12604194747513151053, 0., 0., 0., 0., 0., 0., | |||
980 | 0., .78125000000000000000e-2, -.62377810244809812496e-1, | |||
981 | .23080781467370883845, -.50841310636012325368, .69834547012574056043, | |||
982 | -.52572723156526459672, .11464215704954976471e-1, .38698869011491210342, | |||
983 | -.26125646622255207507, -.16951698812361607510, .29773875898928782269, | |||
984 | .20130501202570367491e-1, -.26332493149159310198, | |||
985 | .67734613690401207009e-1, .21207315477103762715, -.11541543390889415193, | |||
986 | -.16249634759782417533, .13885887405041735068, .11996491328010275427, | |||
987 | -.14810432001630926895, -.85177658352556243411e-1, .14918860659904380587, | |||
988 | .57317789510444151564e-1, -.14569827645586660151, | |||
989 | -.35213090145965327390e-1, .13975998126844578198, 0., 0., 0., 0., 0., 0., | |||
990 | 0., 0., .39062500000000000000e-2, -.35101954600803571207e-1, | |||
991 | .14761284084133737720, -.37655033076080192966, .62410290231517322776, | |||
992 | -.64335622317683389875, .28188168266139524244, .22488495672137010675, | |||
993 | -.39393811089283576186, .75184777995770096714e-1, .28472023119398293003, | |||
994 | -.20410910833705899572, -.15590046962908511750, .23814567544617953125, | |||
995 | .54442805556829031204e-1, -.22855930338589720954, | |||
996 | .16303223615756629897e-1, .20172722433875559213, | |||
997 | -.62723406421217419404e-1, -.17012230831020922010, | |||
998 | .91754642766136561612e-1, .13927644821381121197, -.10886600968068418181, | |||
999 | -.11139075654373395292, .11797455976331702879, 0., 0., 0., 0., 0., 0., 0., | |||
1000 | 0., 0., .19531250000000000000e-2, -.19506820659607596598e-1, | |||
1001 | .91865676095362231937e-1, -.26604607809696493849, .51425874205091288223, | |||
1002 | -.66047561132505329292, .48660109511591303851, -.17575661168678285615e-1, | |||
1003 | -.36594333408055703366, .29088854695378694533, .11318677346656537927, | |||
1004 | -.31110645235730182168, .60733219161008787341e-1, .24333848233620420826, | |||
1005 | -.15254312332655419708, -.15995968483455388613, .19010344455215289289, | |||
1006 | .86040636766440260000e-1, -.19652589954665259945, | |||
1007 | -.27633388517205837713e-1, .18660848552712880387, | |||
1008 | -.15942583868416775867e-1, -.16902042462382064786, | |||
1009 | .47278526495327740646e-1, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1010 | .97656250000000000000e-3, -.10731084460857378207e-1, | |||
1011 | .55939644713816406331e-1, -.18118487371914493668, .39914857299829864263, | |||
1012 | -.60812322949933902435, .60011887183061967583, -.26002695805835928795, | |||
1013 | -.20883922404786010096, .38988130966114638081, -.11797833550782589082, | |||
1014 | -.25231824756239520077, .24817859972953934712, .90516417677868996417e-1, | |||
1015 | -.26079073291293066798, .30259468817169480161e-1, .22178195264114178432, | |||
1016 | -.10569877864302048175, -.16679648389266977455, .14637718550245050850, | |||
1017 | .11219272032739559870, -.16359363640525750353, -.64358194509092101393e-1, | |||
1018 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., .48828125000000000000e-3, | |||
1019 | -.58542865274813470967e-2, .33461741635290096452e-1, | |||
1020 | -.11979993155896201271, .29580223766987206958, -.51874761979436016742, | |||
1021 | .62861483498014306968, -.44868895761051453296, .12567502628371529386e-1, | |||
1022 | .35040366183235474275, -.30466868455569500886, -.70903913601490112666e-1, | |||
1023 | .30822791893032512740, -.11969443264190207736, -.20764760317621313946, | |||
1024 | .20629838355452128532, .95269702915334718507e-1, -.22432624768705133300, | |||
1025 | -.33103381593477797101e-2, .20570036048155716333, | |||
1026 | -.62208282720094518964e-1, -.17095309330441436348, 0., 0., 0., 0., 0., 0., | |||
1027 | 0., 0., 0., 0., 0., 0., .24414062500000000000e-3, | |||
1028 | -.31714797501871532475e-2, .19721062526127334100e-1, | |||
1029 | -.77311181185536498246e-1, .21124871792841566575, -.41777980401893650886, | |||
1030 | .59401977834943551650, -.56132417807488349048, .23433675061367565951, | |||
1031 | .20222775295220942126, -.38280372496506190127, .14443804214023095767, | |||
1032 | .22268950939178466797, -.27211314150777981984, -.34184876506180717313e-1, | |||
1033 | .26006498895669734842, -.97650425186005090107e-1, -.19024527660129101293, | |||
1034 | .16789164198044635671, .10875811641651905252, -.19276785058805921298, 0., | |||
1035 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., .12207031250000000000e-3, | |||
1036 | -.17078941137247586143e-2, .11477733754843910060e-1, | |||
1037 | -.48887017020924625462e-1, .14634927241421789683, -.32156282683019547854, | |||
1038 | .52165811920227223937, -.60001958466396926460, .41208501541480733755, | |||
1039 | -.11366945503190350975e-2, -.33968093962672089159, .30955190935923386766, | |||
1040 | .40657421856578262210e-1, -.29873400409871531764, .16094481791768257440, | |||
1041 | .16876122436206497694, -.23650217045022161255, -.33070260090574765012e-1, | |||
1042 | .22985258456375907796, -.68645651043827097771e-1, 0., 0., 0., 0., 0., 0., | |||
1043 | 0., 0., 0., 0., 0., 0., 0., 0., .61035156250000000000e-4, | |||
1044 | -.91501857608428649078e-3, .66085179496951987952e-2, | |||
1045 | -.30383171695850355404e-1, .98840838845366876117e-1, | |||
1046 | -.23855447246420318989, .43322017468145613917, -.58049033744876107191, | |||
1047 | .52533893203742699346, -.20681056202371946180, -.20180000924562504384, | |||
1048 | .37503922291962681797, -.15988102869837429062, -.19823558102762374094, | |||
1049 | .28393023878803799622, -.11188133439357510403e-1, -.24730368377168229255, | |||
1050 | .14731529061377942839, .14878558042884266021, 0., 0., 0., 0., 0., 0., 0., | |||
1051 | 0., 0., 0., 0., 0., 0., 0., 0., .30517578125000000000e-4, | |||
1052 | -.48804277318479845551e-3, .37696080990601968396e-2, | |||
1053 | -.18603912108994738255e-1, .65325006755649582964e-1, | |||
1054 | -.17162960707938819795, .34411527956476971322, -.52289350347082497959, | |||
1055 | .57319653625674910592, -.37662253421045430413, -.14099055105384663902e-1, | |||
1056 | .33265570610216904208, -.30921265572647566661, -.19911390594166455281e-1, | |||
1057 | .28738590811031797718, -.18912130469738472647, -.13235936203215819193, | |||
1058 | .25076406142356675279, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1059 | 0., 0., 0., .15258789062500000000e-4, -.25928719280954633249e-3, | |||
1060 | .21327398937568540428e-2, -.11244626133630732010e-1, | |||
1061 | .42375605740664331966e-1, -.12031130345907846211, .26352562258934426830, | |||
1062 | -.44590628258512682078, .56682835613700749379, -.49116715128261660395, | |||
1063 | .17845943097110339078, .20541650677432497477, -.36739803642257458221, | |||
1064 | .16776034069210108273, .17920950989905112908, -.28867732805385066532, | |||
1065 | .46473465543376206337e-1, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1066 | 0., 0., 0., 0., 0., .76293945312500000000e-5, -.13727610943181290891e-3, | |||
1067 | .11979683091449349286e-2, -.67195313034570709806e-2, | |||
1068 | .27044920779931968175e-1, -.82472196498517457862e-1, | |||
1069 | .19570475044896150093, -.36391620788543817693, .52241392782736588032, | |||
1070 | -.54727504974907879912, .34211551468813581183, .31580472732719957762e-1, | |||
1071 | -.32830006549176759667, .30563797665254420769, .64905014620683140120e-2, | |||
1072 | -.27642986248995073032, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1073 | 0., 0., 0., 0., 0., 0., .38146972656250000000e-5, | |||
1074 | -.72454147007837596854e-4, .66859847582761390285e-3, | |||
1075 | -.39751311980366118437e-2, .17015198650201528366e-1, | |||
1076 | -.55443621868993855715e-1, .14157060481641692131, -.28641242619559616836, | |||
1077 | .45610665490966615415, -.55262786406029265394, .45818352706035500108, | |||
1078 | -.14984403004611673047, -.21163807462970713245, .36007252928843413718, | |||
1079 | -.17030961385712954159, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1080 | 0., 0., 0., 0., 0., 0., 0., .19073486328125000000e-5, | |||
1081 | -.38135049864067468562e-4, .37101393638555730015e-3, | |||
1082 | -.23305339886279723213e-2, .10569913448297127219e-1, | |||
1083 | -.36640175162216897547e-1, .10010476414320235508, -.21860074212675559892, | |||
1084 | .38124757096345313719, -.52020999209879669177, .52172632730659212045, | |||
1085 | -.30841620620308814614, -.50322546186721500184e-1, .32577618885114899053, | |||
1086 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1087 | 0., 0., .95367431640625000000e-6, -.20021483206955925244e-4, | |||
1088 | .20481807322420625431e-3, -.13553476938058909882e-2, | |||
1089 | .64919676350791905019e-2, -.23848725425069251903e-1, | |||
1090 | .69384632678886421292e-1, -.16249711393618776934, .30736618106830314788, | |||
1091 | -.46399909601971539157, .53765031034002467225, -.42598991476520183929, | |||
1092 | .12130445348350215652, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1093 | 0., 0., 0., 0., 0., 0., 0., 0., .47683715820312500000e-6, | |||
1094 | -.10487707828484902486e-4, .11254146162337528943e-3, | |||
1095 | -.78248929534271987118e-3, .39468337145306794566e-2, | |||
1096 | -.15313546659475671763e-1, .47249070825218564146e-1, | |||
1097 | -.11804374107101480543, .24031796927792491122, -.39629215049166341285, | |||
1098 | .51629108968402548545, -.49622372075429782915, 0., 0., 0., 0., 0., 0., 0., | |||
1099 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1100 | .23841857910156250000e-6, -.54823314130625337326e-5, | |||
1101 | .61575377321535518154e-4, -.44877834366497538134e-3, | |||
1102 | .23774612048621955857e-2, -.97136347645161687796e-2, | |||
1103 | .31671599547606636717e-1, -.84028665767000747480e-1, | |||
1104 | .18298487576742964949, -.32647878537696945218, .46970971486488895077, 0., | |||
1105 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1106 | 0., 0., 0., 0., .11920928955078125000e-6, -.28604020001177375838e-5, | |||
1107 | .33559227978295551013e-4, -.25583821662860610560e-3, | |||
1108 | .14201552747787302339e-2, -.60938046986874414969e-2, | |||
1109 | .20930869247951926793e-1, -.58745021125678072911e-1, | |||
1110 | .13613725780285953720, -.26083988356030237586, 0., 0., 0., 0., 0., 0., 0., | |||
1111 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1112 | .59604644775390625000e-7, -.14898180663526043291e-5, | |||
1113 | .18224991282807693921e-4, -.14504433444608833821e-3, | |||
1114 | .84184722720281809548e-3, -.37846965430000478789e-2, | |||
1115 | .13656355548211376864e-1, -.40409541997718853934e-1, | |||
1116 | .99226988101858325902e-1, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1117 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1118 | .29802322387695312500e-7, -.77471708843445529468e-6, | |||
1119 | .98649879372606876995e-5, -.81814934772838523887e-4, | |||
1120 | .49554483992403011328e-3, -.23290922072351413938e-2, | |||
1121 | .88068134250844034186e-2, -.27393666952485719070e-1, 0., 0., 0., 0., 0., | |||
1122 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1123 | 0., 0., 0., .14901161193847656250e-7, -.40226235946098233685e-6, | |||
1124 | .53236418690561306700e-5, -.45933829691164002269e-4, | |||
1125 | .28982005232838857913e-3, -.14212974043211018374e-2, | |||
1126 | .56192363087488842264e-2, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1127 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1128 | .74505805969238281250e-8, -.20858299254133430408e-6, | |||
1129 | .28648457300134381744e-5, -.25677535898258910850e-4, | |||
1130 | .16849420429491355445e-3, -.86062824010315834002e-3, 0., 0., 0., 0., 0., | |||
1131 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1132 | 0., 0., 0., 0., 0., .37252902984619140625e-8, -.10801736017613096861e-6, | |||
1133 | .15376606719887104015e-5, -.14296523739727437959e-4, | |||
1134 | .97419023656050887203e-4, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1135 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1136 | .18626451492309570312e-8, -.55871592916438890146e-7, | |||
1137 | .82331193828137454068e-6, -.79302250528382787666e-5, 0., 0., 0., 0., 0., | |||
1138 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1139 | 0., 0., 0., 0., 0., 0., 0., .93132257461547851562e-9, | |||
1140 | -.28867244235852488244e-7, .43982811713864556957e-6, 0., 0., 0., 0., 0., | |||
1141 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1142 | 0., 0., 0., 0., 0., 0., 0., 0., .46566128730773925781e-9, | |||
1143 | -.14899342093408253335e-7, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1144 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1145 | 0., 0., .23283064365386962891e-9 | |||
1146 | }; | |||
1147 | ||||
1148 | static const double Tright[33 * 33] = | |||
1149 | { | |||
1150 | 1., .86602540378443864678, 0., -.33071891388307382381, 0., | |||
1151 | .20728904939721249057, 0., -.15128841196122722208, 0., | |||
1152 | .11918864298744029244, 0., -.98352013661686631224e-1, 0., | |||
1153 | .83727065404940845733e-1, 0., -.72893399403505841203e-1, 0., | |||
1154 | .64544632643375022436e-1, 0., -.57913170372415565639e-1, 0., | |||
1155 | .52518242575729562263e-1, 0., -.48043311993977520457e-1, 0., | |||
1156 | .44271433659733990243e-1, 0., -.41048928022856771981e-1, 0., | |||
1157 | .38263878662008271459e-1, 0., -.35832844026365304501e-1, 0., 0., | |||
1158 | .50000000000000000000, .96824583655185422130, .57282196186948000082, | |||
1159 | -.21650635094610966169, -.35903516540862679125, .97578093724974971969e-1, | |||
1160 | .26203921611325660506, -.55792409597991015609e-1, -.20644078533943456204, | |||
1161 | .36172381205961199479e-1, .17035068468874958194, | |||
1162 | -.25371838001497225980e-1, -.14501953125000000000, | |||
1163 | .18786835250972344757e-1, .12625507130328301066, | |||
1164 | -.14473795929590520582e-1, -.11179458309419422675, | |||
1165 | .11494434254897626155e-1, .10030855351241635862, | |||
1166 | -.93498556820544479096e-2, -.90964264465390582629e-1, | |||
1167 | .77546391824364392762e-2, .83213457337452292745e-1, | |||
1168 | -.65358085945588638605e-2, -.76680372422574234569e-1, | |||
1169 | .55835321940047427169e-2, .71098828931825789428e-1, | |||
1170 | -.48253327982967591019e-2, -.66274981937248958553e-1, | |||
1171 | .42118078245337801387e-2, .62064306433355646267e-1, | |||
1172 | -.37083386598903548973e-2, 0., 0., .25000000000000000000, | |||
1173 | .73950997288745200531, .83852549156242113615, .23175620272173946716, | |||
1174 | -.37791833195149451496, -.25710129174850522325, .21608307321780204633, | |||
1175 | .22844049245646009157, -.14009503000335388415, -.19897685605518413847, | |||
1176 | .98264706042471226893e-1, .17445445004279014046, | |||
1177 | -.72761100054958328401e-1, -.15463589893742108388, | |||
1178 | .56056770591708784481e-1, .13855313872640495158, | |||
1179 | -.44517752443294564781e-1, -.12534277657695128850, | |||
1180 | .36211835346039665762e-1, .11434398255136139683, | |||
1181 | -.30033588409423828125e-1, -.10506705408753910481, | |||
1182 | .25313077840725783008e-1, .97149327637744872155e-1, | |||
1183 | -.21624927200393328444e-1, -.90319582367202122625e-1, | |||
1184 | .18688433567711780666e-1, .84372291635345108584e-1, | |||
1185 | -.16312261561845420752e-1, -.79149526894804751586e-1, | |||
1186 | .14362333871852474757e-1, 0., 0., 0., .12500000000000000000, | |||
1187 | .49607837082461073572, .82265291131801144317, .59621200088559103072, | |||
1188 | -.80054302859059362371e-1, -.42612156697795759420, | |||
1189 | -.90098145270865592887e-1, .29769623255090078484, .13630307904779758221, | |||
1190 | -.21638835185708931831, -.14600247270306082052, .16348801804014290453, | |||
1191 | .14340708728599057249, -.12755243353979286190, -.13661523715071346961, | |||
1192 | .10215585947881057394, .12864248070157166547, -.83592528025348693602e-1, | |||
1193 | -.12066728689302565222, .69633728678718053052e-1, .11314245177331919532, | |||
1194 | -.58882939251410088028e-1, -.10621835858758221487, | |||
1195 | .50432266865187597572e-1, .99916834723527771581e-1, | |||
1196 | -.43672094283057258509e-1, -.94206380251950852413e-1, | |||
1197 | .38181356812697746418e-1, .89035739656537771225e-1, | |||
1198 | -.33661934598216332678e-1, 0., 0., 0., 0., .62500000000000000000e-1, | |||
1199 | .31093357409581873586, .67604086414949799246, .75644205980613611039, | |||
1200 | .28990586430124175741, -.30648508196770360914, -.35801372616842500052, | |||
1201 | .91326869828709014708e-1, .31127929687500000000, .90915752838698393094e-2, | |||
1202 | -.25637381283965534330, -.57601077850322797594e-1, .21019685709225757945, | |||
1203 | .81244992138514014256e-1, -.17375078516720988858, | |||
1204 | -.92289437277967051125e-1, .14527351914265391374, | |||
1205 | .96675340792832019889e-1, -.12289485697108543415, | |||
1206 | -.97448175340011084006e-1, .10511755943298339844, | |||
1207 | .96242247086378239657e-1, -.90822942272780513537e-1, | |||
1208 | -.93966350452322132384e-1, .79189411876493712558e-1, | |||
1209 | .91139307067989309325e-1, -.69613039934383197265e-1, | |||
1210 | -.88062491671135767870e-1, .61646331729340817494e-1, 0., 0., 0., 0., 0., | |||
1211 | .31250000000000000000e-1, .18684782411095934408, .50176689760410660236, | |||
1212 | .74784031498626095398, .56472001151566251186, -.14842464993721351203e-1, | |||
1213 | -.41162920273003120936, -.20243071230196532282, .23772054897172750436, | |||
1214 | .24963810923972235950, -.12116179938394678936, -.24330535483519110663, | |||
1215 | .47903849781124471359e-1, .22133299683101224293, | |||
1216 | -.20542915138527200983e-2, -.19653465717678146728, | |||
1217 | -.26818172626509178444e-1, .17319122357631210944, | |||
1218 | .45065391411065545445e-1, -.15253391395444065941, | |||
1219 | -.56543897711725408302e-1, .13469154928743585367, | |||
1220 | .63632471400208840155e-1, -.11941684923913523817, | |||
1221 | -.67828850207933293098e-1, .10636309084510652670, | |||
1222 | .70095786922999181504e-1, -.95187373095150709082e-1, 0., 0., 0., 0., 0., | |||
1223 | 0., .15625000000000000000e-1, .10909562534194485289, | |||
1224 | .34842348626527747318, .64461114561628111443, .69382480527334683659, | |||
1225 | .29551102358528827763, -.25527584713978439819, -.38878771718544715394, | |||
1226 | -.82956185835347407489e-2, .31183177761966943912, .12831420840372374767, | |||
1227 | -.22067618205599434368, -.17569196937129496961, .14598057000132284135, | |||
1228 | .18864406621763419484, -.89921002550386645767e-1, -.18571835020187122114, | |||
1229 | .48967672227195481777e-1, .17584685670380332798, | |||
1230 | -.19267984545067426324e-1, -.16335437520503462738, | |||
1231 | -.22598055455032407594e-2, .15032800884170631129, | |||
1232 | .17883358353754640871e-1, -.13774837869432209951, | |||
1233 | -.29227555960587143675e-1, .12604194747513151053, 0., 0., 0., 0., 0., 0., | |||
1234 | 0., .78125000000000000000e-2, .62377810244809812496e-1, | |||
1235 | .23080781467370883845, .50841310636012325368, .69834547012574056043, | |||
1236 | .52572723156526459672, .11464215704954976471e-1, -.38698869011491210342, | |||
1237 | -.26125646622255207507, .16951698812361607510, .29773875898928782269, | |||
1238 | -.20130501202570367491e-1, -.26332493149159310198, | |||
1239 | -.67734613690401207009e-1, .21207315477103762715, .11541543390889415193, | |||
1240 | -.16249634759782417533, -.13885887405041735068, .11996491328010275427, | |||
1241 | .14810432001630926895, -.85177658352556243411e-1, -.14918860659904380587, | |||
1242 | .57317789510444151564e-1, .14569827645586660151, | |||
1243 | -.35213090145965327390e-1, -.13975998126844578198, 0., 0., 0., 0., 0., 0., | |||
1244 | 0., 0., .39062500000000000000e-2, .35101954600803571207e-1, | |||
1245 | .14761284084133737720, .37655033076080192966, .62410290231517322776, | |||
1246 | .64335622317683389875, .28188168266139524244, -.22488495672137010675, | |||
1247 | -.39393811089283576186, -.75184777995770096714e-1, .28472023119398293003, | |||
1248 | .20410910833705899572, -.15590046962908511750, -.23814567544617953125, | |||
1249 | .54442805556829031204e-1, .22855930338589720954, .16303223615756629897e-1, | |||
1250 | -.20172722433875559213, -.62723406421217419404e-1, .17012230831020922010, | |||
1251 | .91754642766136561612e-1, -.13927644821381121197, -.10886600968068418181, | |||
1252 | .11139075654373395292, .11797455976331702879, 0., 0., 0., 0., 0., 0., 0., | |||
1253 | 0., 0., .19531250000000000000e-2, .19506820659607596598e-1, | |||
1254 | .91865676095362231937e-1, .26604607809696493849, .51425874205091288223, | |||
1255 | .66047561132505329292, .48660109511591303851, .17575661168678285615e-1, | |||
1256 | -.36594333408055703366, -.29088854695378694533, .11318677346656537927, | |||
1257 | .31110645235730182168, .60733219161008787341e-1, -.24333848233620420826, | |||
1258 | -.15254312332655419708, .15995968483455388613, .19010344455215289289, | |||
1259 | -.86040636766440260000e-1, -.19652589954665259945, | |||
1260 | .27633388517205837713e-1, .18660848552712880387, .15942583868416775867e-1, | |||
1261 | -.16902042462382064786, -.47278526495327740646e-1, 0., 0., 0., 0., 0., 0., | |||
1262 | 0., 0., 0., 0., .97656250000000000000e-3, .10731084460857378207e-1, | |||
1263 | .55939644713816406331e-1, .18118487371914493668, .39914857299829864263, | |||
1264 | .60812322949933902435, .60011887183061967583, .26002695805835928795, | |||
1265 | -.20883922404786010096, -.38988130966114638081, -.11797833550782589082, | |||
1266 | .25231824756239520077, .24817859972953934712, -.90516417677868996417e-1, | |||
1267 | -.26079073291293066798, -.30259468817169480161e-1, .22178195264114178432, | |||
1268 | .10569877864302048175, -.16679648389266977455, -.14637718550245050850, | |||
1269 | .11219272032739559870, .16359363640525750353, -.64358194509092101393e-1, | |||
1270 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., .48828125000000000000e-3, | |||
1271 | .58542865274813470967e-2, .33461741635290096452e-1, .11979993155896201271, | |||
1272 | .29580223766987206958, .51874761979436016742, .62861483498014306968, | |||
1273 | .44868895761051453296, .12567502628371529386e-1, -.35040366183235474275, | |||
1274 | -.30466868455569500886, .70903913601490112666e-1, .30822791893032512740, | |||
1275 | .11969443264190207736, -.20764760317621313946, -.20629838355452128532, | |||
1276 | .95269702915334718507e-1, .22432624768705133300, | |||
1277 | -.33103381593477797101e-2, -.20570036048155716333, | |||
1278 | -.62208282720094518964e-1, .17095309330441436348, 0., 0., 0., 0., 0., 0., | |||
1279 | 0., 0., 0., 0., 0., 0., .24414062500000000000e-3, | |||
1280 | .31714797501871532475e-2, .19721062526127334100e-1, | |||
1281 | .77311181185536498246e-1, .21124871792841566575, .41777980401893650886, | |||
1282 | .59401977834943551650, .56132417807488349048, .23433675061367565951, | |||
1283 | -.20222775295220942126, -.38280372496506190127, -.14443804214023095767, | |||
1284 | .22268950939178466797, .27211314150777981984, -.34184876506180717313e-1, | |||
1285 | -.26006498895669734842, -.97650425186005090107e-1, .19024527660129101293, | |||
1286 | .16789164198044635671, -.10875811641651905252, -.19276785058805921298, 0., | |||
1287 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., .12207031250000000000e-3, | |||
1288 | .17078941137247586143e-2, .11477733754843910060e-1, | |||
1289 | .48887017020924625462e-1, .14634927241421789683, .32156282683019547854, | |||
1290 | .52165811920227223937, .60001958466396926460, .41208501541480733755, | |||
1291 | .11366945503190350975e-2, -.33968093962672089159, -.30955190935923386766, | |||
1292 | .40657421856578262210e-1, .29873400409871531764, .16094481791768257440, | |||
1293 | -.16876122436206497694, -.23650217045022161255, .33070260090574765012e-1, | |||
1294 | .22985258456375907796, .68645651043827097771e-1, 0., 0., 0., 0., 0., 0., | |||
1295 | 0., 0., 0., 0., 0., 0., 0., 0., .61035156250000000000e-4, | |||
1296 | .91501857608428649078e-3, .66085179496951987952e-2, | |||
1297 | .30383171695850355404e-1, .98840838845366876117e-1, .23855447246420318989, | |||
1298 | .43322017468145613917, .58049033744876107191, .52533893203742699346, | |||
1299 | .20681056202371946180, -.20180000924562504384, -.37503922291962681797, | |||
1300 | -.15988102869837429062, .19823558102762374094, .28393023878803799622, | |||
1301 | .11188133439357510403e-1, -.24730368377168229255, -.14731529061377942839, | |||
1302 | .14878558042884266021, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1303 | 0., 0., .30517578125000000000e-4, .48804277318479845551e-3, | |||
1304 | .37696080990601968396e-2, .18603912108994738255e-1, | |||
1305 | .65325006755649582964e-1, .17162960707938819795, .34411527956476971322, | |||
1306 | .52289350347082497959, .57319653625674910592, .37662253421045430413, | |||
1307 | -.14099055105384663902e-1, -.33265570610216904208, -.30921265572647566661, | |||
1308 | .19911390594166455281e-1, .28738590811031797718, .18912130469738472647, | |||
1309 | -.13235936203215819193, -.25076406142356675279, 0., 0., 0., 0., 0., 0., | |||
1310 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., .15258789062500000000e-4, | |||
1311 | .25928719280954633249e-3, .21327398937568540428e-2, | |||
1312 | .11244626133630732010e-1, .42375605740664331966e-1, .12031130345907846211, | |||
1313 | .26352562258934426830, .44590628258512682078, .56682835613700749379, | |||
1314 | .49116715128261660395, .17845943097110339078, -.20541650677432497477, | |||
1315 | -.36739803642257458221, -.16776034069210108273, .17920950989905112908, | |||
1316 | .28867732805385066532, .46473465543376206337e-1, 0., 0., 0., 0., 0., 0., | |||
1317 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., .76293945312500000000e-5, | |||
1318 | .13727610943181290891e-3, .11979683091449349286e-2, | |||
1319 | .67195313034570709806e-2, .27044920779931968175e-1, | |||
1320 | .82472196498517457862e-1, .19570475044896150093, .36391620788543817693, | |||
1321 | .52241392782736588032, .54727504974907879912, .34211551468813581183, | |||
1322 | -.31580472732719957762e-1, -.32830006549176759667, -.30563797665254420769, | |||
1323 | .64905014620683140120e-2, .27642986248995073032, 0., 0., 0., 0., 0., 0., | |||
1324 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., .38146972656250000000e-5, | |||
1325 | .72454147007837596854e-4, .66859847582761390285e-3, | |||
1326 | .39751311980366118437e-2, .17015198650201528366e-1, | |||
1327 | .55443621868993855715e-1, .14157060481641692131, .28641242619559616836, | |||
1328 | .45610665490966615415, .55262786406029265394, .45818352706035500108, | |||
1329 | .14984403004611673047, -.21163807462970713245, -.36007252928843413718, | |||
1330 | -.17030961385712954159, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1331 | 0., 0., 0., 0., 0., 0., 0., .19073486328125000000e-5, | |||
1332 | .38135049864067468562e-4, .37101393638555730015e-3, | |||
1333 | .23305339886279723213e-2, .10569913448297127219e-1, | |||
1334 | .36640175162216897547e-1, .10010476414320235508, .21860074212675559892, | |||
1335 | .38124757096345313719, .52020999209879669177, .52172632730659212045, | |||
1336 | .30841620620308814614, -.50322546186721500184e-1, -.32577618885114899053, | |||
1337 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1338 | 0., 0., .95367431640625000000e-6, .20021483206955925244e-4, | |||
1339 | .20481807322420625431e-3, .13553476938058909882e-2, | |||
1340 | .64919676350791905019e-2, .23848725425069251903e-1, | |||
1341 | .69384632678886421292e-1, .16249711393618776934, .30736618106830314788, | |||
1342 | .46399909601971539157, .53765031034002467225, .42598991476520183929, | |||
1343 | .12130445348350215652, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1344 | 0., 0., 0., 0., 0., 0., 0., 0., .47683715820312500000e-6, | |||
1345 | .10487707828484902486e-4, .11254146162337528943e-3, | |||
1346 | .78248929534271987118e-3, .39468337145306794566e-2, | |||
1347 | .15313546659475671763e-1, .47249070825218564146e-1, .11804374107101480543, | |||
1348 | .24031796927792491122, .39629215049166341285, .51629108968402548545, | |||
1349 | .49622372075429782915, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1350 | 0., 0., 0., 0., 0., 0., 0., 0., 0., .23841857910156250000e-6, | |||
1351 | .54823314130625337326e-5, .61575377321535518154e-4, | |||
1352 | .44877834366497538134e-3, .23774612048621955857e-2, | |||
1353 | .97136347645161687796e-2, .31671599547606636717e-1, | |||
1354 | .84028665767000747480e-1, .18298487576742964949, .32647878537696945218, | |||
1355 | .46970971486488895077, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1356 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., .11920928955078125000e-6, | |||
1357 | .28604020001177375838e-5, .33559227978295551013e-4, | |||
1358 | .25583821662860610560e-3, .14201552747787302339e-2, | |||
1359 | .60938046986874414969e-2, .20930869247951926793e-1, | |||
1360 | .58745021125678072911e-1, .13613725780285953720, .26083988356030237586, | |||
1361 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1362 | 0., 0., 0., 0., 0., 0., .59604644775390625000e-7, | |||
1363 | .14898180663526043291e-5, .18224991282807693921e-4, | |||
1364 | .14504433444608833821e-3, .84184722720281809548e-3, | |||
1365 | .37846965430000478789e-2, .13656355548211376864e-1, | |||
1366 | .40409541997718853934e-1, .99226988101858325902e-1, 0., 0., 0., 0., 0., | |||
1367 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1368 | 0., 0., .29802322387695312500e-7, .77471708843445529468e-6, | |||
1369 | .98649879372606876995e-5, .81814934772838523887e-4, | |||
1370 | .49554483992403011328e-3, .23290922072351413938e-2, | |||
1371 | .88068134250844034186e-2, .27393666952485719070e-1, 0., 0., 0., 0., 0., | |||
1372 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1373 | 0., 0., 0., .14901161193847656250e-7, .40226235946098233685e-6, | |||
1374 | .53236418690561306700e-5, .45933829691164002269e-4, | |||
1375 | .28982005232838857913e-3, .14212974043211018374e-2, | |||
1376 | .56192363087488842264e-2, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1377 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1378 | .74505805969238281250e-8, .20858299254133430408e-6, | |||
1379 | .28648457300134381744e-5, .25677535898258910850e-4, | |||
1380 | .16849420429491355445e-3, .86062824010315834002e-3, 0., 0., 0., 0., 0., | |||
1381 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1382 | 0., 0., 0., 0., 0., .37252902984619140625e-8, .10801736017613096861e-6, | |||
1383 | .15376606719887104015e-5, .14296523739727437959e-4, | |||
1384 | .97419023656050887203e-4, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1385 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1386 | .18626451492309570312e-8, .55871592916438890146e-7, | |||
1387 | .82331193828137454068e-6, .79302250528382787666e-5, 0., 0., 0., 0., 0., | |||
1388 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1389 | 0., 0., 0., 0., 0., 0., 0., .93132257461547851562e-9, | |||
1390 | .28867244235852488244e-7, .43982811713864556957e-6, 0., 0., 0., 0., 0., | |||
1391 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1392 | 0., 0., 0., 0., 0., 0., 0., 0., .46566128730773925781e-9, | |||
1393 | .14899342093408253335e-7, 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1394 | 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., | |||
1395 | 0., 0., .23283064365386962891e-9 | |||
1396 | }; | |||
1397 | ||||
1398 | /* Allocates a workspace for the given maximum number of intervals. | |||
1399 | Note that if the workspace gets filled, the intervals with the | |||
1400 | lowest error estimates are dropped. The maximum number of | |||
1401 | intervals is therefore not the maximum number of intervals | |||
1402 | that will be computed, but merely the size of the buffer. | |||
1403 | */ | |||
1404 | ||||
1405 | /* Compute the product of the fx with one of the inverse | |||
1406 | Vandermonde-like matrices. */ | |||
1407 | ||||
1408 | void | |||
1409 | Vinvfx (const double *fx, double *c, const int d) | |||
1410 | { | |||
1411 | ||||
1412 | int i, j; | |||
1413 | ||||
1414 | switch (d) | |||
1415 | { | |||
1416 | case 0: | |||
1417 | for (i = 0; i <= 4; i++) | |||
1418 | { | |||
1419 | c[i] = 0.0; | |||
1420 | for (j = 0; j <= 4; j++) | |||
1421 | c[i] += V1inv[i * 5 + j] * fx[j * 8]; | |||
1422 | } | |||
1423 | break; | |||
1424 | case 1: | |||
1425 | for (i = 0; i <= 8; i++) | |||
1426 | { | |||
1427 | c[i] = 0.0; | |||
1428 | for (j = 0; j <= 8; j++) | |||
1429 | c[i] += V2inv[i * 9 + j] * fx[j * 4]; | |||
1430 | } | |||
1431 | break; | |||
1432 | case 2: | |||
1433 | for (i = 0; i <= 16; i++) | |||
1434 | { | |||
1435 | c[i] = 0.0; | |||
1436 | for (j = 0; j <= 16; j++) | |||
1437 | c[i] += V3inv[i * 17 + j] * fx[j * 2]; | |||
1438 | } | |||
1439 | break; | |||
1440 | case 3: | |||
1441 | for (i = 0; i <= 32; i++) | |||
1442 | { | |||
1443 | c[i] = 0.0; | |||
1444 | for (j = 0; j <= 32; j++) | |||
1445 | c[i] += V4inv[i * 33 + j] * fx[j]; | |||
1446 | } | |||
1447 | break; | |||
1448 | } | |||
1449 | ||||
1450 | } | |||
1451 | ||||
1452 | ||||
1453 | /* Downdate the interpolation given by the n coefficients c | |||
1454 | by removing the nodes with indices in nans. */ | |||
1455 | ||||
1456 | void | |||
1457 | downdate (double *c, int n, int d, int *nans, int nnans) | |||
1458 | { | |||
1459 | ||||
1460 | static const int bidx[4] = { 0, 6, 16, 34 }; | |||
1461 | double b_new[34], alpha; | |||
1462 | int i, j; | |||
1463 | ||||
1464 | for (i = 0; i <= n + 1; i++) | |||
| ||||
1465 | b_new[i] = bee[bidx[d] + i]; | |||
1466 | for (i = 0; i < nnans; i++) | |||
1467 | { | |||
1468 | b_new[n + 1] = b_new[n + 1] / Lalpha[n]; | |||
1469 | b_new[n] = (b_new[n] + xi[nans[i]] * b_new[n + 1]) / Lalpha[n - 1]; | |||
| ||||
1470 | for (j = n - 1; j > 0; j--) | |||
1471 | b_new[j] = | |||
1472 | (b_new[j] + xi[nans[i]] * b_new[j + 1] - | |||
1473 | Lgamma[j + 1] * b_new[j + 2]) / Lalpha[j - 1]; | |||
1474 | for (j = 0; j <= n; j++) | |||
1475 | b_new[j] = b_new[j + 1]; | |||
1476 | alpha = c[n] / b_new[n]; | |||
1477 | for (j = 0; j < n; j++) | |||
1478 | c[j] -= alpha * b_new[j]; | |||
1479 | c[n] = 0; | |||
1480 | n--; | |||
1481 | } | |||
1482 | ||||
1483 | } | |||
1484 | ||||
1485 | ||||
1486 | /* The actual integration routine. */ | |||
1487 | ||||
1488 | DEFUN (quadcc, args, nargout,octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1489 | "-*- texinfo -*-\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1490 | @deftypefn {Function File} {@var{q} =} quadcc (@var{f}, @var{a}, @var{b})\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1491 | @deftypefnx {Function File} {@var{q} =} quadcc (@var{f}, @var{a}, @var{b}, @var{tol})\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1492 | @deftypefnx {Function File} {@var{q} =} quadcc (@var{f}, @var{a}, @var{b}, @var{tol}, @var{sing})\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1493 | @deftypefnx {Function File} {[@var{q}, @var{err}, @var{nr_points}] =} quadcc (@dots{})\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1494 | Numerically evaluate the integral of @var{f} from @var{a} to @var{b}\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1495 | using the doubly-adaptive Clenshaw-Curtis quadrature described by P. Gonnet\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1496 | in @cite{Increasing the Reliability of Adaptive Quadrature Using Explicit\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1497 | Interpolants}.\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1498 | @var{f} is a function handle, inline function, or string\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1499 | containing the name of the function to evaluate.\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1500 | The function @var{f} must be vectorized and must return a vector of output\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1501 | values if given a vector of input values. For example,\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1502 | \n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1503 | @example\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1504 | f = @@(x) x .* sin (1./x) .* sqrt (abs (1 - x));\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1505 | @end example\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1506 | \n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1507 | @noindent\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1508 | which uses the element-by-element ``dot'' form for all operators.\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1509 | \n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1510 | @var{a} and @var{b} are the lower and upper limits of integration. Either\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1511 | or both limits may be infinite. @code{quadcc} handles an inifinite limit\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1512 | by substituting the variable of integration with @code{x = tan (pi/2*u)}.\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1513 | \n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1514 | The optional argument @var{tol} defines the relative tolerance used to stop\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1515 | the integration procedure. The default value is @math{1e^{-6}}.\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1516 | \n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1517 | The optional argument @var{sing} contains a list of points where the\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1518 | integrand has known singularities, or discontinuities\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1519 | in any of its derivatives, inside the integration interval.\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1520 | For the example above, which has a discontinuity at x=1, the call to\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1521 | @code{quadcc} would be as follows\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1522 | \n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1523 | @example\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1524 | int = quadcc (f, a, b, 1.0e-6, [ 1 ]);\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1525 | @end example\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1526 | \n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1527 | The result of the integration is returned in @var{q}.\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1528 | @var{err} is an estimate of the absolute integration error and\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1529 | @var{nr_points} is the number of points at which the integrand was evaluated.\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1530 | If the adaptive integration did not converge, the value of\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1531 | @var{err} will be larger than the requested tolerance. Therefore, it is\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1532 | recommended to verify this value for difficult integrands.\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1533 | \n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1534 | @code{quadcc} is capable of dealing with non-numeric\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1535 | values of the integrand such as @code{NaN} or @code{Inf}.\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1536 | If the integral diverges, and @code{quadcc} detects this,\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1537 | then a warning is issued and @code{Inf} or @code{-Inf} is returned.\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1538 | \n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1539 | Note: @code{quadcc} is a general purpose quadrature algorithm\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1540 | and, as such, may be less efficient for a smooth or otherwise\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1541 | well-behaved integrand than other methods such as @code{quadgk}.\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1542 | \n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1543 | The algorithm uses Clenshaw-Curtis quadrature rules of increasing\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1544 | degree in each interval and bisects the interval if either the\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1545 | function does not appear to be smooth or a rule of maximum\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1546 | degree has been reached. The error estimate is computed from the\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1547 | L2-norm of the difference between two successive interpolations\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1548 | of the integrand over the nodes of the respective quadrature rules.\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1549 | \n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1550 | Reference: P. Gonnet, @cite{Increasing the Reliability of Adaptive\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1551 | Quadrature Using Explicit Interpolants}, ACM Transactions on\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1552 | Mathematical Software, Vol. 37, Issue 3, Article No. 3, 2010.\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1553 | @seealso{quad, quadv, quadl, quadgk, trapz, dblquad, triplequad}\n\octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1554 | @end deftypefn")octave_value_list Fquadcc (const octave_value_list& args, int nargout) | |||
1555 | { | |||
1556 | octave_value_list retval; | |||
1557 | ||||
1558 | /* Some constants that we will need. */ | |||
1559 | static const int n[4] = { 4, 8, 16, 32 }; | |||
1560 | static const int skip[4] = { 8, 4, 2, 1 }; | |||
1561 | static const int idx[4] = { 0, 5, 14, 31 }; | |||
1562 | static const double w = M_SQRT21.41421356237309504880 / 2; | |||
1563 | static const int ndiv_max = 20; | |||
1564 | ||||
1565 | /* Arguments left and right */ | |||
1566 | int nargin = args.length (); | |||
1567 | octave_function *fcn; | |||
1568 | double a, b, tol, *sing; | |||
1569 | ||||
1570 | /* Variables needed for transforming the integrand. */ | |||
1571 | bool wrap = false; | |||
1572 | double xw; | |||
1573 | ||||
1574 | /* Stuff we will need to call the integrand. */ | |||
1575 | octave_value_list fargs, fvals; | |||
1576 | ||||
1577 | /* Actual variables (as opposed to constants above). */ | |||
1578 | double m, h, ml, hl, mr, hr, temp; | |||
1579 | double igral, err, igral_final, err_final; | |||
1580 | int nivals, neval = 0; | |||
1581 | int i, j, d, split, t; | |||
1582 | int nnans, nans[33]; | |||
1583 | cquad_ival *iv, *ivl, *ivr; | |||
1584 | double nc, ncdiff; | |||
1585 | ||||
1586 | ||||
1587 | /* Parse the input arguments. */ | |||
1588 | if (nargin < 3) | |||
1589 | { | |||
1590 | print_usage (); | |||
1591 | return retval; | |||
1592 | } | |||
1593 | ||||
1594 | if (args(0).is_function_handle () || args(0).is_inline_function ()) | |||
1595 | fcn = args(0).function_value (); | |||
1596 | else | |||
1597 | { | |||
1598 | std::string fcn_name = unique_symbol_name ("__quadcc_fcn__"); | |||
1599 | std::string fname = "function y = "; | |||
1600 | fname.append (fcn_name); | |||
1601 | fname.append ("(x) y = "); | |||
1602 | fcn = extract_function (args(0), "quadcc", fcn_name, fname, | |||
1603 | "; endfunction"); | |||
1604 | } | |||
1605 | ||||
1606 | if (! args(1).is_real_scalar ()) | |||
1607 | { | |||
1608 | error ("quadcc: lower limit of integration (A) must be a single real scalar"); | |||
1609 | return retval; | |||
1610 | } | |||
1611 | else | |||
1612 | a = args(1).double_value (); | |||
1613 | ||||
1614 | if (! args(2).is_real_scalar ()) | |||
1615 | { | |||
1616 | error ("quadcc: upper limit of integration (B) must be a single real scalar"); | |||
1617 | return retval; | |||
1618 | } | |||
1619 | else | |||
1620 | b = args(2).double_value (); | |||
1621 | ||||
1622 | if (nargin < 4 || args(3).is_empty ()) | |||
1623 | tol = 1.0e-6; | |||
1624 | else if (! args(3).is_real_scalar () || args(3).double_value () <= 0) | |||
1625 | { | |||
1626 | error ("quadcc: tolerance (TOL) must be a single real scalar > 0"); | |||
1627 | return retval; | |||
1628 | } | |||
1629 | else | |||
1630 | tol = args(3).double_value (); | |||
1631 | ||||
1632 | if (nargin < 5) | |||
1633 | { | |||
1634 | nivals = 1; | |||
1635 | } | |||
1636 | else if (!(args(4).is_real_scalar () || args(4).is_real_matrix ())) | |||
1637 | { | |||
1638 | error ("quadcc: list of singularities (SING) must be a vector of real values"); | |||
1639 | return retval; | |||
1640 | } | |||
1641 | else | |||
1642 | { | |||
1643 | nivals = 1 + args(4).numel (); | |||
1644 | } | |||
1645 | ||||
1646 | int cquad_heapsize = (nivals >= min_cquad_heapsize200 ? nivals + 1 | |||
1647 | : min_cquad_heapsize200); | |||
1648 | /* The interval heap. */ | |||
1649 | OCTAVE_LOCAL_BUFFER (cquad_ival, ivals, cquad_heapsize)octave_local_buffer<cquad_ival> _buffer_ivals (cquad_heapsize ); cquad_ival *ivals = _buffer_ivals; | |||
1650 | OCTAVE_LOCAL_BUFFER (double, iivals, cquad_heapsize)octave_local_buffer<double> _buffer_iivals (cquad_heapsize ); double *iivals = _buffer_iivals; | |||
1651 | OCTAVE_LOCAL_BUFFER (int, heap, cquad_heapsize)octave_local_buffer<int> _buffer_heap (cquad_heapsize); int *heap = _buffer_heap; | |||
1652 | ||||
1653 | if (nivals == 1) | |||
1654 | { | |||
1655 | iivals[0] = a; | |||
1656 | iivals[1] = b; | |||
1657 | } | |||
1658 | else | |||
1659 | { | |||
1660 | // Intervals around singularities | |||
1661 | sing = args(4).array_value ().fortran_vec (); | |||
1662 | iivals[0] = a; | |||
1663 | for (i = 0; i < nivals - 1; i++) | |||
1664 | iivals[i + 1] = sing[i]; | |||
1665 | iivals[nivals] = b; | |||
1666 | } | |||
1667 | ||||
1668 | /* If a or b are +/-Inf, transform the integral. */ | |||
1669 | if (xisinf (a) || xisinf (b)) | |||
1670 | { | |||
1671 | wrap = true; | |||
1672 | for (i = 0; i < nivals + 1; i++) | |||
1673 | if (xisinf (iivals[i])) | |||
1674 | iivals[i] = gnulib::copysign (1.0, iivals[i]); | |||
1675 | else | |||
1676 | iivals[i] = 2.0 * atan (iivals[i]) / M_PI3.14159265358979323846; | |||
1677 | } | |||
1678 | ||||
1679 | ||||
1680 | /* Initialize the heaps. */ | |||
1681 | for (i = 0; i < cquad_heapsize; i++) | |||
1682 | heap[i] = i; | |||
1683 | ||||
1684 | /* Create the first interval(s). */ | |||
1685 | igral = 0.0; | |||
1686 | err = 0.0; | |||
1687 | for (j = 0; j < nivals; j++) | |||
1688 | { | |||
1689 | ||||
1690 | /* Initialize the interval. */ | |||
1691 | iv = &(ivals[heap[j]]); | |||
1692 | m = (iivals[j] + iivals[j + 1]) / 2; | |||
1693 | h = (iivals[j + 1] - iivals[j]) / 2; | |||
1694 | nnans = 0; | |||
1695 | ColumnVector ex (33); | |||
1696 | if (wrap) | |||
1697 | { | |||
1698 | for (i = 0; i <= n[3]; i++) | |||
1699 | ex(i) = tan (M_PI3.14159265358979323846 / 2 * (m + xi[i] * h)); | |||
1700 | } | |||
1701 | else | |||
1702 | { | |||
1703 | for (i = 0; i <= n[3]; i++) | |||
1704 | ex(i) = m + xi[i] * h; | |||
1705 | } | |||
1706 | fargs(0) = ex; | |||
1707 | fvals = feval (fcn, fargs, 1); | |||
1708 | if (fvals.length () != 1 || ! fvals(0).is_real_matrix ()) | |||
1709 | { | |||
1710 | error ("quadcc: integrand F must return a single, real-valued vector"); | |||
1711 | return retval; | |||
1712 | } | |||
1713 | Matrix effex = fvals(0).matrix_value (); | |||
1714 | if (effex.length () != ex.length ()) | |||
1715 | { | |||
1716 | error ("quadcc: integrand F must return a single, real-valued vector of the same size as the input"); | |||
1717 | return retval; | |||
1718 | } | |||
1719 | for (i = 0; i <= n[3]; i++) | |||
1720 | { | |||
1721 | iv->fx[i] = effex(i); | |||
1722 | if (wrap) | |||
1723 | { | |||
1724 | xw = ex(i); | |||
1725 | iv->fx[i] *= (1.0 + xw * xw) * M_PI3.14159265358979323846 / 2; | |||
1726 | } | |||
1727 | neval++; | |||
1728 | if (! xfinite (iv->fx[i])) | |||
1729 | { | |||
1730 | nans[nnans++] = i; | |||
1731 | iv->fx[i] = 0.0; | |||
1732 | } | |||
1733 | } | |||
1734 | Vinvfx (iv->fx, &(iv->c[idx[3]]), 3); | |||
1735 | Vinvfx (iv->fx, &(iv->c[idx[2]]), 2); | |||
1736 | Vinvfx (iv->fx, &(iv->c[0]), 0); | |||
1737 | for (i = 0; i < nnans; i++) | |||
1738 | iv->fx[nans[i]] = octave_NaN; | |||
1739 | iv->a = iivals[j]; | |||
1740 | iv->b = iivals[j + 1]; | |||
1741 | iv->depth = 3; | |||
1742 | iv->rdepth = 1; | |||
1743 | iv->ndiv = 0; | |||
1744 | iv->igral = 2 * h * iv->c[idx[3]] * w; | |||
1745 | nc = 0.0; | |||
1746 | for (i = n[2] + 1; i <= n[3]; i++) | |||
1747 | { | |||
1748 | temp = iv->c[idx[3] + i]; | |||
1749 | nc += temp * temp; | |||
1750 | } | |||
1751 | ncdiff = nc; | |||
1752 | for (i = 0; i <= n[2]; i++) | |||
1753 | { | |||
1754 | temp = iv->c[idx[2] + i] - iv->c[idx[3] + i]; | |||
1755 | ncdiff += temp * temp; | |||
1756 | nc += iv->c[idx[3] + i] * iv->c[idx[3] + i]; | |||
1757 | } | |||
1758 | ncdiff = sqrt (ncdiff); | |||
1759 | nc = sqrt (nc); | |||
1760 | iv->err = ncdiff * 2 * h; | |||
1761 | if (ncdiff / nc > 0.1 && iv->err < 2 * h * nc) | |||
1762 | iv->err = 2 * h * nc; | |||
1763 | ||||
1764 | /* Tabulate this interval's data. */ | |||
1765 | igral += iv->igral; | |||
1766 | err += iv->err; | |||
1767 | ||||
1768 | /* Sift it up the heap. */ | |||
1769 | i = j; | |||
1770 | while (i > 0 && ivals[heap[i / 2]].err < ivals[heap[i]].err) | |||
1771 | { | |||
1772 | temp = heap[i]; | |||
1773 | heap[i] = heap[i / 2]; | |||
1774 | heap[i / 2] = temp; | |||
1775 | i /= 2; | |||
1776 | } | |||
1777 | ||||
1778 | } | |||
1779 | ||||
1780 | ||||
1781 | /* Initialize some global values. */ | |||
1782 | igral_final = 0.0; | |||
1783 | err_final = 0.0; | |||
1784 | ||||
1785 | ||||
1786 | /* Main loop. */ | |||
1787 | while (nivals > 0 && err > 0.0 && err > fabs (igral) * tol | |||
1788 | && !(err_final > fabs (igral) * tol | |||
1789 | && err - err_final < fabs (igral) * tol)) | |||
1790 | { | |||
1791 | ||||
1792 | /* Allow the user to interrupt. */ | |||
1793 | OCTAVE_QUIToctave_quit (); | |||
1794 | ||||
1795 | /* Put our finger on the interval with the largest error. */ | |||
1796 | iv = &(ivals[heap[0]]); | |||
1797 | m = (iv->a + iv->b) / 2; | |||
1798 | h = (iv->b - iv->a) / 2; | |||
1799 | ||||
1800 | #if (DEBUG_QUADCC0) | |||
1801 | printf ("quadcc: processing ival %i (of %i) with [%e,%e] int=%e, err=%e, depth=%i\n", | |||
1802 | heap[0], nivals, iv->a, iv->b, iv->igral, iv->err, iv->depth); | |||
1803 | #endif | |||
1804 | ||||
1805 | /* Should we try to increase the degree? */ | |||
1806 | if (iv->depth < 3) | |||
1807 | { | |||
1808 | ||||
1809 | /* Keep tabs on some variables. */ | |||
1810 | d = ++iv->depth; | |||
1811 | ||||
1812 | /* Get the new (missing) function values */ | |||
1813 | { | |||
1814 | ColumnVector ex (n[d] / 2); | |||
1815 | if (wrap) | |||
1816 | { | |||
1817 | for (i = 0; i < n[d] / 2; i++) | |||
1818 | ex(i) = tan (M_PI3.14159265358979323846 / 2 * (m + xi[(2 * i + 1) * skip[d]] * h)); | |||
1819 | } | |||
1820 | else | |||
1821 | { | |||
1822 | for (i = 0; i < n[d] / 2; i++) | |||
1823 | ex(i) = m + xi[(2 * i + 1) * skip[d]] * h; | |||
1824 | } | |||
1825 | fargs(0) = ex; | |||
1826 | fvals = feval (fcn, fargs, 1); | |||
1827 | if (fvals.length () != 1 || ! fvals(0).is_real_matrix ()) | |||
1828 | { | |||
1829 | error ("quadcc: integrand F must return a single, real-valued vector"); | |||
1830 | return retval; | |||
1831 | } | |||
1832 | Matrix effex = fvals(0).matrix_value (); | |||
1833 | if (effex.length () != ex.length ()) | |||
1834 | { | |||
1835 | error ("quadcc: integrand F must return a single, real-valued vector of the same size as the input"); | |||
1836 | return retval; | |||
1837 | } | |||
1838 | neval += effex.length (); | |||
1839 | for (i = 0; i < n[d] / 2; i++) | |||
1840 | { | |||
1841 | j = (2 * i + 1) * skip[d]; | |||
1842 | iv->fx[j] = effex(i); | |||
1843 | if (wrap) | |||
1844 | { | |||
1845 | xw = ex(i); | |||
1846 | iv->fx[j] *= (1.0 + xw * xw) * M_PI3.14159265358979323846 / 2; | |||
1847 | } | |||
1848 | } | |||
1849 | } | |||
1850 | nnans = 0; | |||
1851 | for (i = 0; i <= 32; i += skip[d]) | |||
1852 | { | |||
1853 | if (! xfinite (iv->fx[i])) | |||
1854 | { | |||
1855 | nans[nnans++] = i; | |||
1856 | iv->fx[i] = 0.0; | |||
1857 | } | |||
1858 | } | |||
1859 | ||||
1860 | /* Compute the new coefficients. */ | |||
1861 | Vinvfx (iv->fx, &(iv->c[idx[d]]), d); | |||
1862 | /* Downdate any NaNs. */ | |||
1863 | if (nnans > 0) | |||
1864 | { | |||
1865 | downdate (&(iv->c[idx[d]]), n[d], d, nans, nnans); | |||
1866 | for (i = 0; i < nnans; i++) | |||
1867 | iv->fx[nans[i]] = octave_NaN; | |||
1868 | } | |||
1869 | ||||
1870 | /* Compute the error estimate. */ | |||
1871 | nc = 0.0; | |||
1872 | for (i = n[d - 1] + 1; i <= n[d]; i++) | |||
1873 | { | |||
1874 | temp = iv->c[idx[d] + i]; | |||
1875 | nc += temp * temp; | |||
1876 | } | |||
1877 | ncdiff = nc; | |||
1878 | for (i = 0; i <= n[d - 1]; i++) | |||
1879 | { | |||
1880 | temp = iv->c[idx[d - 1] + i] - iv->c[idx[d] + i]; | |||
1881 | ncdiff += temp * temp; | |||
1882 | nc += iv->c[idx[d] + i] * iv->c[idx[d] + i]; | |||
1883 | } | |||
1884 | ncdiff = sqrt (ncdiff); | |||
1885 | nc = sqrt (nc); | |||
1886 | iv->err = ncdiff * 2 * h; | |||
1887 | /* Compute the local integral. */ | |||
1888 | iv->igral = 2 * h * w * iv->c[idx[d]]; | |||
1889 | /* Split the interval prematurely? */ | |||
1890 | split = (nc > 0 && ncdiff / nc > 0.1); | |||
1891 | } | |||
1892 | ||||
1893 | /* Maximum degree reached, just split. */ | |||
1894 | else | |||
1895 | { | |||
1896 | split = 1; | |||
1897 | } | |||
1898 | ||||
1899 | ||||
1900 | /* Should we drop this interval? */ | |||
1901 | if ((m + h * xi[0]) >= (m + h * xi[1]) | |||
1902 | || (m + h * xi[31]) >= (m + h * xi[32]) | |||
1903 | || iv->err < fabs (iv->igral) | |||
1904 | * std::numeric_limits<double>::epsilon () * 10) | |||
1905 | { | |||
1906 | ||||
1907 | #if (DEBUG_QUADCC0) | |||
1908 | printf ("quadcc: dropping ival %i (of %i) with [%e,%e] int=%e, err=%e, depth=%i\n", | |||
1909 | heap[0], nivals, iv->a, iv->b, iv->igral, iv->err, iv->depth); | |||
1910 | #endif | |||
1911 | ||||
1912 | /* Keep this interval's contribution */ | |||
1913 | err_final += iv->err; | |||
1914 | igral_final += iv->igral; | |||
1915 | /* Swap with the last element on the heap */ | |||
1916 | t = heap[nivals - 1]; | |||
1917 | heap[nivals - 1] = heap[0]; | |||
1918 | heap[0] = t; | |||
1919 | nivals--; | |||
1920 | /* Fix up the heap */ | |||
1921 | i = 0; | |||
1922 | while (2 * i + 1 < nivals) | |||
1923 | { | |||
1924 | /* Get the kids */ | |||
1925 | j = 2 * i + 1; | |||
1926 | /* If the j+1st entry exists and is larger than the jth, | |||
1927 | use it instead. */ | |||
1928 | if (j + 1 < nivals | |||
1929 | && ivals[heap[j + 1]].err >= ivals[heap[j]].err) | |||
1930 | j++; | |||
1931 | /* Do we need to move the ith entry up? */ | |||
1932 | if (ivals[heap[j]].err <= ivals[heap[i]].err) | |||
1933 | break; | |||
1934 | else | |||
1935 | { | |||
1936 | t = heap[j]; | |||
1937 | heap[j] = heap[i]; | |||
1938 | heap[i] = t; | |||
1939 | i = j; | |||
1940 | } | |||
1941 | } | |||
1942 | ||||
1943 | } | |||
1944 | ||||
1945 | /* Do we need to split this interval? */ | |||
1946 | else if (split) | |||
1947 | { | |||
1948 | ||||
1949 | /* Some values we will need often... */ | |||
1950 | d = iv->depth; | |||
1951 | /* Generate the interval on the left */ | |||
1952 | ivl = &(ivals[heap[nivals++]]); | |||
1953 | ivl->a = iv->a; | |||
1954 | ivl->b = m; | |||
1955 | ml = (ivl->a + ivl->b) / 2; | |||
1956 | hl = h / 2; | |||
1957 | ivl->depth = 0; | |||
1958 | ivl->rdepth = iv->rdepth + 1; | |||
1959 | ivl->fx[0] = iv->fx[0]; | |||
1960 | ivl->fx[32] = iv->fx[16]; | |||
1961 | { | |||
1962 | ColumnVector ex (n[0] - 1); | |||
1963 | if (wrap) | |||
1964 | { | |||
1965 | for (i = 0; i < n[0] - 1; i++) | |||
1966 | ex(i) = tan (M_PI3.14159265358979323846 / 2 * (ml + xi[(i + 1) * skip[0]] * hl)); | |||
1967 | } | |||
1968 | else | |||
1969 | { | |||
1970 | for (i = 0; i < n[0] - 1; i++) | |||
1971 | ex(i) = ml + xi[(i + 1) * skip[0]] * hl; | |||
1972 | } | |||
1973 | fargs(0) = ex; | |||
1974 | fvals = feval (fcn, fargs, 1); | |||
1975 | if (fvals.length () != 1 || ! fvals(0).is_real_matrix ()) | |||
1976 | { | |||
1977 | error ("quadcc: integrand F must return a single, real-valued vector"); | |||
1978 | return retval; | |||
1979 | } | |||
1980 | Matrix effex = fvals(0).matrix_value (); | |||
1981 | if (effex.length () != ex.length ()) | |||
1982 | { | |||
1983 | error ("quadcc: integrand F must return a single, real-valued vector of the same size as the input"); | |||
1984 | return retval; | |||
1985 | } | |||
1986 | neval += effex.length (); | |||
1987 | for (i = 0; i < n[0] - 1; i++) | |||
1988 | { | |||
1989 | j = (i + 1) * skip[0]; | |||
1990 | ivl->fx[j] = effex(i); | |||
1991 | if (wrap) | |||
1992 | { | |||
1993 | xw = ex(i); | |||
1994 | ivl->fx[j] *= (1.0 + xw * xw) * M_PI3.14159265358979323846 / 2; | |||
1995 | } | |||
1996 | } | |||
1997 | } | |||
1998 | nnans = 0; | |||
1999 | for (i = 0; i <= 32; i += skip[0]) | |||
2000 | { | |||
2001 | if (! xfinite (ivl->fx[i])) | |||
2002 | { | |||
2003 | nans[nnans++] = i; | |||
2004 | ivl->fx[i] = 0.0; | |||
2005 | } | |||
2006 | } | |||
2007 | Vinvfx (ivl->fx, ivl->c, 0); | |||
2008 | if (nnans > 0) | |||
2009 | { | |||
2010 | downdate (ivl->c, n[0], 0, nans, nnans); | |||
2011 | for (i = 0; i < nnans; i++) | |||
2012 | ivl->fx[nans[i]] = octave_NaN; | |||
2013 | } | |||
2014 | for (i = 0; i <= n[d]; i++) | |||
2015 | { | |||
2016 | ivl->c[idx[d] + i] = 0.0; | |||
2017 | for (j = i; j <= n[d]; j++) | |||
2018 | ivl->c[idx[d] + i] += Tleft[i * 33 + j] * iv->c[idx[d] + j]; | |||
2019 | } | |||
2020 | ncdiff = 0.0; | |||
2021 | for (i = 0; i <= n[0]; i++) | |||
2022 | { | |||
2023 | temp = ivl->c[i] - ivl->c[idx[d] + i]; | |||
2024 | ncdiff += temp * temp; | |||
2025 | } | |||
2026 | for (i = n[0] + 1; i <= n[d]; i++) | |||
2027 | { | |||
2028 | temp = ivl->c[idx[d] + i]; | |||
2029 | ncdiff += temp * temp; | |||
2030 | } | |||
2031 | ncdiff = sqrt (ncdiff); | |||
2032 | ivl->err = ncdiff * h; | |||
2033 | /* Check for divergence. */ | |||
2034 | ivl->ndiv = iv->ndiv + (fabs (iv->c[0]) > 0 | |||
2035 | && ivl->c[0] / iv->c[0] > 2); | |||
2036 | if (ivl->ndiv > ndiv_max && 2 * ivl->ndiv > ivl->rdepth) | |||
2037 | { | |||
2038 | igral = gnulib::copysign (octave_Inf, igral); | |||
2039 | warning ("quadcc: divergent integral detected"); | |||
2040 | break; | |||
2041 | } | |||
2042 | ||||
2043 | /* Compute the local integral. */ | |||
2044 | ivl->igral = h * w * ivl->c[0]; | |||
2045 | ||||
2046 | ||||
2047 | /* Generate the interval on the right */ | |||
2048 | ivr = &(ivals[heap[nivals++]]); | |||
2049 | ivr->a = m; | |||
2050 | ivr->b = iv->b; | |||
2051 | mr = (ivr->a + ivr->b) / 2; | |||
2052 | hr = h / 2; | |||
2053 | ivr->depth = 0; | |||
2054 | ivr->rdepth = iv->rdepth + 1; | |||
2055 | ivr->fx[0] = iv->fx[16]; | |||
2056 | ivr->fx[32] = iv->fx[32]; | |||
2057 | { | |||
2058 | ColumnVector ex (n[0] - 1); | |||
2059 | if (wrap) | |||
2060 | { | |||
2061 | for (i = 0; i < n[0] - 1; i++) | |||
2062 | ex(i) = tan (M_PI3.14159265358979323846 / 2 * (mr + xi[(i + 1) * skip[0]] * hr)); | |||
2063 | } | |||
2064 | else | |||
2065 | { | |||
2066 | for (i = 0; i < n[0] - 1; i++) | |||
2067 | ex(i) = mr + xi[(i + 1) * skip[0]] * hr; | |||
2068 | } | |||
2069 | fargs(0) = ex; | |||
2070 | fvals = feval (fcn, fargs, 1); | |||
2071 | if (fvals.length () != 1 || ! fvals(0).is_real_matrix ()) | |||
2072 | { | |||
2073 | error ("quadcc: integrand F must return a single, real-valued vector"); | |||
2074 | return retval; | |||
2075 | } | |||
2076 | Matrix effex = fvals(0).matrix_value (); | |||
2077 | if (effex.length () != ex.length ()) | |||
2078 | { | |||
2079 | error ("quadcc: integrand F must return a single, real-valued vector of the same size as the input"); | |||
2080 | return retval; | |||
2081 | } | |||
2082 | neval += effex.length (); | |||
2083 | for (i = 0; i < n[0] - 1; i++) | |||
2084 | { | |||
2085 | j = (i + 1) * skip[0]; | |||
2086 | ivr->fx[j] = effex(i); | |||
2087 | if (wrap) | |||
2088 | { | |||
2089 | xw = ex(i); | |||
2090 | ivr->fx[j] *= (1.0 + xw * xw) * M_PI3.14159265358979323846 / 2; | |||
2091 | } | |||
2092 | } | |||
2093 | } | |||
2094 | nnans = 0; | |||
2095 | for (i = 0; i <= 32; i += skip[0]) | |||
2096 | { | |||
2097 | if (! xfinite (ivr->fx[i])) | |||
2098 | { | |||
2099 | nans[nnans++] = i; | |||
2100 | ivr->fx[i] = 0.0; | |||
2101 | } | |||
2102 | } | |||
2103 | Vinvfx (ivr->fx, ivr->c, 0); | |||
2104 | if (nnans > 0) | |||
2105 | { | |||
2106 | downdate (ivr->c, n[0], 0, nans, nnans); | |||
2107 | for (i = 0; i < nnans; i++) | |||
2108 | ivr->fx[nans[i]] = octave_NaN; | |||
2109 | } | |||
2110 | for (i = 0; i <= n[d]; i++) | |||
2111 | { | |||
2112 | ivr->c[idx[d] + i] = 0.0; | |||
2113 | for (j = i; j <= n[d]; j++) | |||
2114 | ivr->c[idx[d] + i] += Tright[i * 33 + j] * iv->c[idx[d] + j]; | |||
2115 | } | |||
2116 | ncdiff = 0.0; | |||
2117 | for (i = 0; i <= n[0]; i++) | |||
2118 | { | |||
2119 | temp = ivr->c[i] - ivr->c[idx[d] + i]; | |||
2120 | ncdiff += temp * temp; | |||
2121 | } | |||
2122 | for (i = n[0] + 1; i <= n[d]; i++) | |||
2123 | { | |||
2124 | temp = ivr->c[idx[d] + i]; | |||
2125 | ncdiff += temp * temp; | |||
2126 | } | |||
2127 | ncdiff = sqrt (ncdiff); | |||
2128 | ivr->err = ncdiff * h; | |||
2129 | /* Check for divergence. */ | |||
2130 | ivr->ndiv = iv->ndiv + (fabs (iv->c[0]) > 0 | |||
2131 | && ivr->c[0] / iv->c[0] > 2); | |||
2132 | if (ivr->ndiv > ndiv_max && 2 * ivr->ndiv > ivr->rdepth) | |||
2133 | { | |||
2134 | igral = gnulib::copysign (octave_Inf, igral); | |||
2135 | warning ("quadcc: divergent integral detected"); | |||
2136 | break; | |||
2137 | } | |||
2138 | ||||
2139 | /* Compute the local integral. */ | |||
2140 | ivr->igral = h * w * ivr->c[0]; | |||
2141 | ||||
2142 | ||||
2143 | /* Fix-up the heap: we now have one interval on top | |||
2144 | that we don't need any more and two new, unsorted | |||
2145 | ones at the bottom. */ | |||
2146 | /* Flip the last interval to the top of the heap and | |||
2147 | sift down. */ | |||
2148 | t = heap[nivals - 1]; | |||
2149 | heap[nivals - 1] = heap[0]; | |||
2150 | heap[0] = t; | |||
2151 | nivals--; | |||
2152 | /* Sift this interval back down the heap. */ | |||
2153 | i = 0; | |||
2154 | while (2 * i + 1 < nivals - 1) | |||
2155 | { | |||
2156 | j = 2 * i + 1; | |||
2157 | if (j + 1 < nivals - 1 | |||
2158 | && ivals[heap[j + 1]].err >= ivals[heap[j]].err) | |||
2159 | j++; | |||
2160 | if (ivals[heap[j]].err <= ivals[heap[i]].err) | |||
2161 | break; | |||
2162 | else | |||
2163 | { | |||
2164 | t = heap[j]; | |||
2165 | heap[j] = heap[i]; | |||
2166 | heap[i] = t; | |||
2167 | i = j; | |||
2168 | } | |||
2169 | } | |||
2170 | ||||
2171 | /* Now grab the last interval and sift it up the heap. */ | |||
2172 | i = nivals - 1; | |||
2173 | while (i > 0) | |||
2174 | { | |||
2175 | j = (i - 1) / 2; | |||
2176 | if (ivals[heap[j]].err < ivals[heap[i]].err) | |||
2177 | { | |||
2178 | t = heap[j]; | |||
2179 | heap[j] = heap[i]; | |||
2180 | heap[i] = t; | |||
2181 | i = j; | |||
2182 | } | |||
2183 | else | |||
2184 | break; | |||
2185 | } | |||
2186 | ||||
2187 | ||||
2188 | } | |||
2189 | ||||
2190 | /* Otherwise, just fix-up the heap. */ | |||
2191 | else | |||
2192 | { | |||
2193 | i = 0; | |||
2194 | while (2 * i + 1 < nivals) | |||
2195 | { | |||
2196 | j = 2 * i + 1; | |||
2197 | if (j + 1 < nivals | |||
2198 | && ivals[heap[j + 1]].err >= ivals[heap[j]].err) | |||
2199 | j++; | |||
2200 | if (ivals[heap[j]].err <= ivals[heap[i]].err) | |||
2201 | break; | |||
2202 | else | |||
2203 | { | |||
2204 | t = heap[j]; | |||
2205 | heap[j] = heap[i]; | |||
2206 | heap[i] = t; | |||
2207 | i = j; | |||
2208 | } | |||
2209 | } | |||
2210 | ||||
2211 | } | |||
2212 | ||||
2213 | /* If the heap is about to overflow, remove the last two intervals. */ | |||
2214 | while (nivals > cquad_heapsize - 2) | |||
2215 | { | |||
2216 | iv = &(ivals[heap[nivals - 1]]); | |||
2217 | #if (DEBUG_QUADCC0) | |||
2218 | printf ("quadcc: dropping ival %i (of %i) with [%e,%e] int=%e, err=%e, depth=%i\n", | |||
2219 | heap[0], nivals, iv->a, iv->b, iv->igral, iv->err, iv->depth); | |||
2220 | #endif | |||
2221 | err_final += iv->err; | |||
2222 | igral_final += iv->igral; | |||
2223 | nivals--; | |||
2224 | } | |||
2225 | ||||
2226 | /* Collect the value of the integral and error. */ | |||
2227 | igral = igral_final; | |||
2228 | err = err_final; | |||
2229 | for (i = 0; i < nivals; i++) | |||
2230 | { | |||
2231 | igral += ivals[heap[i]].igral; | |||
2232 | err += ivals[heap[i]].err; | |||
2233 | } | |||
2234 | ||||
2235 | } | |||
2236 | ||||
2237 | /* Dump the contents of the heap. */ | |||
2238 | #if (DEBUG_QUADCC0) | |||
2239 | for (i = 0; i < nivals; i++) | |||
2240 | { | |||
2241 | iv = &(ivals[heap[i]]); | |||
2242 | printf ("quadcc: ival %i (%i) with [%e,%e], int=%e, err=%e, depth=%i, rdepth=%i, ndiv=%i\n", | |||
2243 | i, heap[i], iv->a, iv->b, iv->igral, iv->err, iv->depth, | |||
2244 | iv->rdepth, iv->ndiv); | |||
2245 | } | |||
2246 | #endif | |||
2247 | ||||
2248 | /* Clean up and present the results. */ | |||
2249 | if (nargout > 2) | |||
2250 | retval(2) = neval; | |||
2251 | if (nargout > 1) | |||
2252 | retval(1) = err; | |||
2253 | retval(0) = igral; | |||
2254 | /* All is well that ends well. */ | |||
2255 | return retval; | |||
2256 | } | |||
2257 | ||||
2258 | ||||
2259 | /* | |||
2260 | %!assert (quadcc (@sin, -pi, pi), 0, 1e-6) | |||
2261 | %!assert (quadcc (inline ("sin"),- pi, pi), 0, 1e-6) | |||
2262 | %!assert (quadcc ("sin", -pi, pi), 0, 1e-6) | |||
2263 | ||||
2264 | %!assert (quadcc (@sin, -pi, 0), -2, 1e-6) | |||
2265 | %!assert (quadcc (@sin, 0, pi), 2, 1e-6) | |||
2266 | %!assert (quadcc (@(x) 1./sqrt (x), 0, 1), 2, 1e-6) | |||
2267 | %!assert (quadcc (@(x) 1./(sqrt (x).*(x+1)), 0, Inf), pi, 1e-6) | |||
2268 | ||||
2269 | %!assert (quadcc (@(x) exp (-x .^ 2), -Inf, Inf), sqrt (pi), 1e-6) | |||
2270 | %!assert (quadcc (@(x) exp (-x .^ 2), -Inf, 0), sqrt (pi)/2, 1e-6) | |||
2271 | ||||
2272 | ## Test function with NaNs in interval | |||
2273 | %!function y = __nansin (x) | |||
2274 | %! nan_locs = [-3*pi/4, -pi/4, 0, pi/3, pi/2, pi]; | |||
2275 | %! y = sin (x); | |||
2276 | %! idx = min (abs (bsxfun (@minus, x(:), nan_locs)), [], 2); | |||
2277 | %! y(idx < 1e-10) = NaN; | |||
2278 | %!endfunction | |||
2279 | ||||
2280 | %!test | |||
2281 | %! [q, err, npoints] = quadcc ("__nansin", -pi, pi); | |||
2282 | %! assert (q, 0, 1e-6); | |||
2283 | %! assert (err, 0, 15*eps); | |||
2284 | ||||
2285 | %% Test input validation | |||
2286 | %!error (quadcc ()) | |||
2287 | %!error (quadcc (@sin)) | |||
2288 | %!error (quadcc (@sin, 0)) | |||
2289 | %!error (quadcc (@sin, ones (2), pi)) | |||
2290 | %!error (quadcc (@sin, -i, pi)) | |||
2291 | %!error (quadcc (@sin, 0, ones (2))) | |||
2292 | %!error (quadcc (@sin, 0, i)) | |||
2293 | %!error (quadcc (@sin, 0, pi, 0)) | |||
2294 | %!error (quadcc (@sin, 0, pi, 1e-6, [ i ])) | |||
2295 | */ |