pxfm/gamma/
digamma_coeffs.rs

1/*
2 * // Copyright (c) Radzivon Bartoshyk 8/2025. All rights reserved.
3 * //
4 * // Redistribution and use in source and binary forms, with or without modification,
5 * // are permitted provided that the following conditions are met:
6 * //
7 * // 1.  Redistributions of source code must retain the above copyright notice, this
8 * // list of conditions and the following disclaimer.
9 * //
10 * // 2.  Redistributions in binary form must reproduce the above copyright notice,
11 * // this list of conditions and the following disclaimer in the documentation
12 * // and/or other materials provided with the distribution.
13 * //
14 * // 3.  Neither the name of the copyright holder nor the names of its
15 * // contributors may be used to endorse or promote products derived from
16 * // this software without specific prior written permission.
17 * //
18 * // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
19 * // AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
20 * // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
21 * // DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
22 * // FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
23 * // DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
24 * // SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
25 * // CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
26 * // OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
27 * // OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
28 */
29
30/*
31Poly numerators for digamma with step 0.25 from 0 to 1: [-0.99999999,-0.75], [-0.75, -0.5] etc.
32Generated by Wolfram:
33<<FunctionApproximations`
34ClearAll["Global`*"]
35f[x_]:=PolyGamma[x+1]
36{err0,approx,err1}=MiniMaxApproximation[f[z],{z,{-0.99999999,-0.75},11,11},WorkingPrecision->75,MaxIterations->100]
37num=Numerator[approx];
38den=Denominator[approx];
39poly=num;
40coeffs=CoefficientList[poly,z];
41TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
42poly=den;
43coeffs=CoefficientList[poly,z];
44TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
45*/
46static DIGAMMA_P: [[(u64, u64); 12]; 5] = [
47    [
48        (0x3c56cb8ed2e8f89d, 0xbfe2788cfc6fb619),
49        (0xbc3f4a887143e1a0, 0x3face2294e40a23d),
50        (0x3c66cbeeb663e02b, 0x3ff93da24c048557),
51        (0xbc99dce69443c3a1, 0x3ffb9ca4e7e3b448),
52        (0x3c86642ccc010b7c, 0x3febdebcbbc26740),
53        (0xbc6b0b2f44dae75c, 0x3fcfec9b12222c41),
54        (0xbbf4498c8a1b0940, 0x3fa62080e59486a8),
55        (0x3c0b9dc350a33198, 0x3f72bc74facca682),
56        (0xbbcd99184e051c10, 0x3f32c94b7b2fb002),
57        (0xbb8738e809ec0cee, 0x3ee48ad4af286d1f),
58        (0x3b130f69e7502da5, 0x3e84788fd8c782e7),
59        (0xbaabf9ce6f6a39c0, 0x3e075df82b55c605),
60    ],
61    [
62        (0x3c56cb907008b543, 0xbfe2788cfc6fb619),
63        (0xbc100842dc7f9e7f, 0x3fb01301a3e98b73),
64        (0xbc823b1b4a0eb7aa, 0x3ff938b063dae04f),
65        (0xbc9ad21ad2ab96c5, 0x3ffb55c29a135ea7),
66        (0xbc7d5ba993767236, 0x3feb4f810b04c753),
67        (0x3c64772eb3244da3, 0x3fcee4afd0031184),
68        (0x3c31c1810f303099, 0x3fa513f613ed0729),
69        (0x3c0bd0845bf7f471, 0x3f717f05f5d10c47),
70        (0xbbd5725001e2660d, 0x3f311a7e461bfdbc),
71        (0xbb8598e1d6bf14c5, 0x3ee21971cbf5769f),
72        (0x3b2c987b53229a3f, 0x3e81461991174455),
73        (0xbaa40e6fc1c3cdbd, 0x3e029575e1cf2d2e),
74    ],
75    [
76        (0x3c56cb90701fbfa6, 0xbfe2788cfc6fb619),
77        (0xbc5fccfe13017bdd, 0x3fb1afc092258787),
78        (0xbc5e1697dc2cc4ab, 0x3ff9337b6bb944de),
79        (0x3c7ec613451b31b9, 0x3ffb0fd93c22113b),
80        (0x3c69294dd477ad47, 0x3feac3c7b64adea4),
81        (0xbc69ac8bff5565b5, 0x3fcde66e558d148f),
82        (0x3c41acf6837637f5, 0x3fa4154b581a555a),
83        (0xbc000c1483343db8, 0x3f70580604407f51),
84        (0x3bbc673a5ee2174f, 0x3f2f290c93113075),
85        (0x3b7ea2e3314e056e, 0x3edfed2db5039158),
86        (0xbaf8d7013fff21f7, 0x3e7d35fdf2716965),
87        (0xba739bf5ff69075f, 0x3dfdb070ae4d9275),
88    ],
89    [
90        (0x3c56cb90701fbfab, 0xbfe2788cfc6fb619),
91        (0x3c32d7bac835354c, 0x3fb34c4468e4cca5),
92        (0x3c83dfd7ab2663ee, 0x3ff92df4d81668e8),
93        (0x3c99f27ceb6a6462, 0x3ffaca12fb44387d),
94        (0xbc873d42ae8bf0ee, 0x3fea39db2aa6755c),
95        (0x3c6f0709b773a576, 0x3fccee948fe7d9a3),
96        (0xbc319c78a343f2ba, 0x3fa320f0d611b675),
97        (0xbba5feb85ab5acbc, 0x3f6e855445cb1dae),
98        (0x3bba9b56af9b20c5, 0x3f2c5ecc73cbcdd8),
99        (0xbb79a32258828ace, 0x3edc257fc6ffd42d),
100        (0x3b1efad7acd7d7f1, 0x3e78b3710a7cc9b4),
101        (0xba9b941fb3ff6b0c, 0x3df7c1965f2d7194),
102    ],
103    [
104        (0x3c56cb90701fbfab, 0xbfe2788cfc6fb619),
105        (0xbc0f99d5dc5f8b06, 0x3fb4e913ab76f9ab),
106        (0x3c3191e63046dc0e, 0x3ff9281aa8398c82),
107        (0xbc8dca8def824c6e, 0x3ffa845a59ce1115),
108        (0xbc84dbb9117bcea5, 0x3fe9b18e7fe7829a),
109        (0xbc5e38525daac915, 0x3fcbfcbc6601a1fe),
110        (0x3c4df8dcb4769fe3, 0x3fa2364b5ab136dd),
111        (0x3c03a220e74edede, 0x3f6c7b78bf3817d2),
112        (0xbb831144f7e90d73, 0x3f29d06f3db25586),
113        (0x3b7b4bf8f2ec756c, 0x3ed8cd0278bee187),
114        (0x3b0c2bf25bbf39b5, 0x3e74e2f6a044ec34),
115        (0x3a942086a4928038, 0x3df3096d2199706d),
116    ],
117];
118
119static DIGAMMA_Q: [[(u64, u64); 12]; 5] = [
120    [
121        (0x0000000000000000, 0x3ff0000000000000),
122        (0xbca6adb0ceca2bb2, 0x4006042dff16f40e),
123        (0xbc85b67e7a730640, 0x4008379d3bb3e8e0),
124        (0xbc9f004eb95ee438, 0x3ffc7e68e904512d),
125        (0xbc86ee1922dc4c80, 0x3fe407ea97357260),
126        (0xbc4a51002b1b8870, 0x3fc1a70e08dee4dd),
127        (0x3c3c858930786915, 0x3f93c1e4eef4b20b),
128        (0x3be9a73c7d428d88, 0x3f5baa2f04f498b1),
129        (0x3bbad1b8117f7756, 0x3f171c5667c8dffc),
130        (0x3b173ce742b4d895, 0x3ec4e29a45581499),
131        (0xbadceb3243cdae8f, 0x3e60a2d5d1c5f8e4),
132        (0xba72b2562de7b14b, 0x3ddb303e7518d74b),
133    ],
134    [
135        (0x0000000000000000, 0x3ff0000000000000),
136        (0x3c91779e812448f8, 0x4005ed8db27be60e),
137        (0x3ca08f641dca5f12, 0x4007fb6b0dcbcebe),
138        (0xbc9b5b059968427b, 0x3ffc005c8adea86d),
139        (0x3c7a6c4119840b60, 0x3fe37d5ff5b97171),
140        (0xbc530ac1661d05de, 0x3fc0f56578da17ee),
141        (0xbc0418ed2bc636c5, 0x3f92ad4e7387f851),
142        (0xbbc1f659b88045fd, 0x3f599f72920bc12d),
143        (0x3bba548ec376cb95, 0x3f14dae97a930698),
144        (0x3b440fec06d88466, 0x3ec239fbc98052c7),
145        (0xbae8bd7909c83aaf, 0x3e5bcab555b96b0c),
146        (0x3a7da912eade2739, 0x3dd5684408d505b7),
147    ],
148    [
149        (0x0000000000000000, 0x3ff0000000000000),
150        (0xbca91be54288a0e0, 0x4005d735308bb4dd),
151        (0xbcad30b7bbf2b783, 0x4007c03f8e4c4fe6),
152        (0xbc9eeb30365815fb, 0x3ffb854f0db08baa),
153        (0x3c782d2becaf6172, 0x3fe2f7622383d14d),
154        (0xbc5e010ee89a3662, 0x3fc04badf9bdbed9),
155        (0x3c2af3eb8a7d2bde, 0x3f91a94aee0050b0),
156        (0xbbd2450502e41cf9, 0x3f57be27fec5860c),
157        (0xbb4d401f4b8a5f96, 0x3f12d5723c46d383),
158        (0xbb56834644f0201f, 0x3ebfda938f459b86),
159        (0x3ada0645f44a4bcd, 0x3e574563cf974833),
160        (0xba7147ae6edcbc3d, 0x3dd0f0bfcad5f047),
161    ],
162    [
163        (0x0000000000000000, 0x3ff0000000000000),
164        (0xbc9045b43a44ffbd, 0x4005c0dfe199018b),
165        (0x3cadc1c9c8ce4d30, 0x40078563dd87aa1b),
166        (0xbc7ee1018a57916d, 0x3ffb0bbe462659fc),
167        (0x3c84b1e3cb3c10f4, 0x3fe2743f18457656),
168        (0x3c5ce0b8927bff8a, 0x3fbf4f3f5d1d2c88),
169        (0x3c227cc874c2781a, 0x3f90b1f8cc772eb7),
170        (0x3bea9648da0df3bd, 0x3f55fdcf409bac26),
171        (0x3bbc08910e142370, 0x3f11003e3494b94d),
172        (0x3b5a6e0762b57f19, 0x3ebbd3ea01684109),
173        (0xbae425ecbc8c7993, 0x3e537e51126d11b7),
174        (0xba656f95223fc44f, 0x3dcadec5fd0f622c),
175    ],
176    [
177        (0x0000000000000000, 0x3ff0000000000000),
178        (0x3caa0a86221a0755, 0x4005aa867d56dca4),
179        (0xbcac770719607c75, 0x40074ac4f6ca5a47),
180        (0x3c90d03a36128d7f, 0x3ffa93815f190569),
181        (0x3c4d48dbdc6c9a7c, 0x3fe1f3c3b4044619),
182        (0x3c55f2ca2ba72f08, 0x3fbe11c76c6789f6),
183        (0x3c07ca0eacaa9a9c, 0x3f8f8d146c709524),
184        (0xbbd3ae3ef7bb9792, 0x3f545be615cbebf2),
185        (0xbb86e82ec7f13b7f, 0x3f0ead31c1a5b093),
186        (0xbb3275f4a13e5259, 0x3eb84d337c38d867),
187        (0xbaf557dfebeee67e, 0x3e505544e8f1eb85),
188        (0xba664f9845fcddf2, 0x3dc559cc156a72a5),
189    ],
190];
191
192#[inline]
193#[allow(clippy::type_complexity)]
194pub(crate) fn digamma_coeefs(x: f64) -> (&'static [(u64, u64); 12], &'static [(u64, u64); 12]) {
195    if x < 1.25 {
196        const FOUR: u64 = 2 << 52;
197        let idx = f64::from_bits(x.to_bits().wrapping_add(FOUR)) as usize;
198        let p = &DIGAMMA_P[idx];
199        let q = &DIGAMMA_Q[idx];
200        (p, q)
201    } else if x <= 1.461632144 {
202        // Polynomial generated by Wolfram Mathematica:
203        // <<FunctionApproximations`
204        // ClearAll["Global`*"]
205        // f[x_]:=PolyGamma[x+1]
206        // {err0,approx,err1}=MiniMaxApproximation[f[z],{z,{0.25,0.4616321445},11,11},WorkingPrecision->75,MaxIterations->100]
207        // num=Numerator[approx];
208        // poly=den;
209        // coeffs=CoefficientList[poly,z];
210        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
211        static P: [(u64, u64); 12] = [
212            (0x3c56cb90701fbfac, 0xbfe2788cfc6fb619),
213            (0x3c519fc5fe329992, 0x3fb66b6f7ca77b95),
214            (0x3c934d7e8434b86d, 0x3ff92255627f42b7),
215            (0x3c637005ae2e1665, 0x3ffa43326e6dda36),
216            (0x3c3cadafce6428f4, 0x3fe9338d9b2539f1),
217            (0xbc5934a30d152b54, 0x3fcb1febac693a97),
218            (0xbc3d4ec7808fa04f, 0x3fa1636286a89dfe),
219            (0xbbf577e9f832bec3, 0x3f6aafbc9ad0d5ad),
220            (0xbbcd6cc934d6ba4e, 0x3f279eaf8168ee05),
221            (0xbb774ef4099bad0c, 0x3ed605efaec77016),
222            (0xbb13d911fcebafac, 0x3e71dae45db651bb),
223            (0x3a8f943200c1fd92, 0x3deefe67af5b7c69),
224        ];
225
226        // Polynomial generated by Wolfram Mathematica:
227        // <<FunctionApproximations`
228        // ClearAll["Global`*"]
229        // f[x_]:=PolyGamma[x+1]
230        // {err0,approx,err1}=MiniMaxApproximation[f[z],{z,{0.25,0.4616321445},11,11},WorkingPrecision->75,MaxIterations->100]
231        // den=Denominator[approx];
232        // poly=den;
233        // coeffs=CoefficientList[poly,z];
234        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
235        static Q: [(u64, u64); 12] = [
236            (0x0000000000000000, 0x3ff0000000000000),
237            (0x3ca0fb55209ad44d, 0x4005959bb2d38810),
238            (0xbc4a7dccb43a2a97, 0x40071428acb23be5),
239            (0xbc8b9d316e4d2eac, 0x3ffa24400e09040c),
240            (0xbc82963c4efb8cd1, 0x3fe17dfa742912bc),
241            (0x3c57e765f5582745, 0x3fbcf24965915c1f),
242            (0xbc2bfe1e01bc34f0, 0x3f8de96be6a73955),
243            (0x3bff6a6cc1679110, 0x3f52ef034cdc300f),
244            (0xbb95e707835497ae, 0x3f0bd94abae4b483),
245            (0xbb435b31a77f1eac, 0x3eb5673711efb89a),
246            (0xbab7879ef73e896a, 0x3e4bb1a59a808117),
247            (0xba50ff72e575b8bd, 0x3dc13fef9cd36cd7),
248        ];
249        (&P, &Q)
250    } else if x <= 1.75 {
251        // Polynomial generated by Wolfram Mathematica:
252        // <<FunctionApproximations`
253        // ClearAll["Global`*"]
254        // f[x_]:=PolyGamma[x+1]
255        // {err0,approx,err1}=MiniMaxApproximation[f[z],{z,{0.4616321445,0.75},11,11},WorkingPrecision->75,MaxIterations->100]
256        // num=Numerator[approx];
257        // den=Denominator[approx];
258        // poly=num;
259        // coeffs=CoefficientList[poly,z];
260        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
261        static P: [(u64, u64); 12] = [
262            (0x3c56cb907020a112, 0xbfe2788cfc6fb619),
263            (0xbc3cc4d6d2ddccf8, 0x3fb8010bfcb83141),
264            (0x3c953d47923ee59b, 0x3ff91bf87a8009ef),
265            (0x3c856cf1c819abcd, 0x3ff9fee67ab017a1),
266            (0xbc7e5e24a77dd53e, 0x3fe8b0e3f74e8969),
267            (0xbc60cf6e08c2898c, 0x3fca3dc2ac288997),
268            (0xbc4d6745fb4dd03b, 0x3fa08eb928fceb5d),
269            (0x3c0153615a07e895, 0x3f68e95d5669dc9c),
270            (0x3bcc104646590a9e, 0x3f25818c7facb551),
271            (0xbb6ce130dd405fa2, 0x3ed36f27c8982809),
272            (0xbae5a7d3777ea48d, 0x3e6e4b3d604133f4),
273            (0xba86b002676d830d, 0x3de9027ae9347c77),
274        ];
275
276        // Generated by Wolfram Mathematica:
277        // <<FunctionApproximations`
278        // ClearAll["Global`*"]
279        // f[x_]:=PolyGamma[x+1]
280        // {err0,approx,err1}=MiniMaxApproximation[f[z],{z,{0.4616321445,0.75},11,11},WorkingPrecision->75,MaxIterations->100]
281        // den=Denominator[approx];
282        // poly=den;
283        // coeffs=CoefficientList[poly,z];
284        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
285        static Q: [(u64, u64); 12] = [
286            (0x0000000000000000, 0x3ff0000000000000),
287            (0x3ca13038b3e9b2fe, 0x40057fa6128d0d08),
288            (0xbcacf21951ec09a9, 0x4006db1750ade888),
289            (0x3c9fd268ab672291, 0x3ff9b0c5a589896f),
290            (0xbc4b818d4493b88a, 0x3fe104d8aeed391e),
291            (0xbc5531e244a9d256, 0x3fbbce38cef1dcd3),
292            (0x3c14e017cdbf9c7b, 0x3f8c45ceb3fbe0a3),
293            (0x3bff3d4bd6748d1c, 0x3f51898acc32fa08),
294            (0xbb6c0ebdb9714f4b, 0x3f09266424ad60b1),
295            (0x3b30abd2021dc5ea, 0x3eb2ba49a890927f),
296            (0xba9629f3ce206d76, 0x3e474b9f374cf052),
297            (0x39f3beca069e2142, 0x3dbba0c4466d2dcd),
298        ];
299        (&P, &Q)
300    } else if x <= 2. {
301        // Polynomial generated by Wolfram Mathematica:
302        // <<FunctionApproximations`
303        // ClearAll["Global`*"]
304        // f[x_]:=PolyGamma[x+1]
305        // {err0,approx}=MiniMaxApproximation[f[z],{z,{0.75,1},11,11},WorkingPrecision->75,MaxIterations->100]
306        // num=Numerator[approx][[1]];
307        // poly=num;
308        // coeffs=CoefficientList[poly,z];
309        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
310        static P: [(u64, u64); 12] = [
311            (0x3c56cb907bb52809, 0xbfe2788cfc6fb619),
312            (0xbc420cb940ff1bad, 0x3fb9b8ce99af6686),
313            (0xbc9554bebef623c7, 0x3ff914b779852ffb),
314            (0x3c81daceacffe906, 0x3ff9b4fb845f87b4),
315            (0xbc589d02262e231e, 0x3fe8251c504ee4e7),
316            (0xbc48b053e0ac1f18, 0x3fc94f0c5b394803),
317            (0xbc12b44611d0cc3d, 0x3f9f6411022b8fee),
318            (0x3bf71626c3a722d6, 0x3f671bef7af39602),
319            (0x3bca802abbe6c032, 0x3f236afc88364acd),
320            (0xbb73ba0bffc04b11, 0x3ed0f71a4569fdd2),
321            (0x3ae6866974b3745b, 0x3e695b19eb1a11d0),
322            (0x3a53fcd00ba0f0cd, 0x3de3dbee8f81e7d1),
323        ];
324        // Generated by Wolfram Mathematica:
325        //<<FunctionApproximations`
326        // ClearAll["Global`*"]
327        // f[x_]:=PolyGamma[x+1]
328        // {err0,approx}=MiniMaxApproximation[f[z],{z,{0.75,1},11,11},WorkingPrecision->75,MaxIterations->100]
329        // den=Denominator[approx][[1]];
330        // poly=den;
331        // coeffs=CoefficientList[poly,z];
332        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
333        static Q: [(u64, u64); 12] = [
334            (0x0000000000000000, 0x3ff0000000000000),
335            (0x3ca79b06fde08bf7, 0x400567d727c212e8),
336            (0xbc9b0efe8b146dee, 0x40069d86c47f0240),
337            (0xbc8299a0687d2bc4, 0x3ff9351a45c9504d),
338            (0xbc7146fcb08e6211, 0x3fe0846a56e8c617),
339            (0xbc482c27fd2b0de3, 0x3fba9ca29930f546),
340            (0x3beaea4e7ed3058f, 0x3f8a96406d35271a),
341            (0x3bf776396bc4adbd, 0x3f5021eb0152afad),
342            (0xbb96d179aedf8a7a, 0x3f0682935c1f4e53),
343            (0xbb469533a5f262d4, 0x3eb033bbf3243fc1),
344            (0x3ad33e062bbd7ac2, 0x3e43535be5d463ef),
345            (0x3a5a7d127f1db626, 0x3db5c4427025e959),
346        ];
347        (&P, &Q)
348    } else if x <= 3. {
349        // x <= 3
350        // Polynomial generated by Wolfram Mathematica:
351        // <<FunctionApproximations`
352        // ClearAll["Global`*"]
353        // f[x_]:=PolyGamma[x+1]
354        // {err0,approx}=MiniMaxApproximation[f[z],{z,{1,2},11,11},WorkingPrecision->75,MaxIterations->100]
355        // num=Numerator[approx][[1]];
356        // poly=num;
357        // coeffs=CoefficientList[poly,z];
358        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
359        static P: [(u64, u64); 12] = [
360            (0x3c56d1cc8bcd2834, 0xbfe2788cfc6fb619),
361            (0x3c50701f78e64a63, 0x3fbd8bb4f2220f15),
362            (0xbc6f91b861a7c3f1, 0x3ff903416798b207),
363            (0x3c8c53da657116b3, 0x3ff910f28e38ea3f),
364            (0xbc8720408446640c, 0x3fe6f4e66e069dae),
365            (0xbc34e63b5e59b9bd, 0x3fc75316d7a07060),
366            (0xbc22cfa33a314334, 0x3f9bd312bdf809fb),
367            (0xbbfd6f2d0a941d18, 0x3f638321ec193179),
368            (0xbbbd1cabd7d7d204, 0x3f1edf0eeec2ed8f),
369            (0x3b62b3641fcc9055, 0x3ec908c824d24ebd),
370            (0x3b09f21b5bfbaac3, 0x3e61119012a56c3c),
371            (0xba7ccbf11ca33a44, 0x3dd7e364822d6038),
372        ];
373        // Generated by Wolfram Mathematica:
374        //<<FunctionApproximations`
375        // ClearAll["Global`*"]
376        // f[x_]:=PolyGamma[x+1]
377        // {err0,approx}=MiniMaxApproximation[f[z],{z,{1,2},11,11},WorkingPrecision->75,MaxIterations->100]
378        // den=Denominator[approx][[1]];
379        // poly=den;
380        // coeffs=CoefficientList[poly,z];
381        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
382        static Q: [(u64, u64); 12] = [
383            (0x0000000000000000, 0x3ff0000000000000),
384            (0x3ca8ae32e83e7971, 0x400532d7f57a1fb8),
385            (0x3ca7e06fe036f034, 0x4006159f6c8559a6),
386            (0x3c696aeb79408ada, 0x3ff8276de15b9949),
387            (0x3c6b4b6c4100ee75, 0x3fdee20a3a9510d1),
388            (0xbc5366197ca045c2, 0x3fb81bb8b46e3ac9),
389            (0x3c2f178bb9cecf1c, 0x3f8726ac2011b61b),
390            (0x3be51ce661a77196, 0x3f4abe5166a64e0d),
391            (0xbb84cb13a26b099d, 0x3f018de8d20c9ae7),
392            (0xbb4f127e563d4fce, 0x3ea772bc233a7f3d),
393            (0xbad8391ca6522a23, 0x3e398817e8c8a3f5),
394            (0x3a46dbfc8fc0a31f, 0x3da9c1eb2bb8a09d),
395        ];
396        (&P, &Q)
397    } else if x <= 4. {
398        // x <= 4
399        // Polynomial generated by Wolfram Mathematica:
400        // <<FunctionApproximations`
401        // ClearAll["Global`*"]
402        // f[x_]:=PolyGamma[x+1]
403        // {err0,approx}=MiniMaxApproximation[f[z],{z,{2,3},11,11},WorkingPrecision->75,MaxIterations->100]
404        // num=Numerator[approx][[1]];
405        // poly=num;
406        // coeffs=CoefficientList[poly,z];
407        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
408        static P: [(u64, u64); 12] = [
409            (0xbc1138539f5170ce, 0xbfe2788cfc6fb617),
410            (0x3c5c28276512f306, 0x3fc1fc2364555134),
411            (0xbc6a09be4de1611c, 0x3ff8e1befdb6bc25),
412            (0x3c74ca0e29a4c28f, 0x3ff7ff202fdee96b),
413            (0xbc81af890e1880a0, 0x3fe50bb10d8d3f83),
414            (0xbc63ea18b51ff6b6, 0x3fc444e2fde41a82),
415            (0x3c3778d316f51d1e, 0x3f969f9868302eb9),
416            (0x3be358e50a96a230, 0x3f5d3a1baf51425c),
417            (0x3bbc38a174d2733f, 0x3f14e8fe4a7dddfb),
418            (0x3b5f5ed1876db6a7, 0x3ebe07d4dbcc1ad4),
419            (0x3af907af35fa32d0, 0x3e51b3f54fae106f),
420            (0x3a62af4ca3fd0f35, 0x3dc4d88e2dc21ca8),
421        ];
422        // Generated by Wolfram Mathematica:
423        //<<FunctionApproximations`
424        // ClearAll["Global`*"]
425        // f[x_]:=PolyGamma[x+1]
426        // {err0,approx}=MiniMaxApproximation[f[z],{z,{2,3},11,11},WorkingPrecision->75,MaxIterations->100]
427        // den=Denominator[approx][[1]];
428        // poly=den;
429        // coeffs=CoefficientList[poly,z];
430        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
431        static Q: [(u64, u64); 12] = [
432            (0x0000000000000000, 0x3ff0000000000000),
433            (0x3c9648605b16dd34, 0x4004d9ced1a59873),
434            (0x3c7dc6b76d368927, 0x400534eafce36582),
435            (0xbc3ce3eb41ed98ed, 0x3ff673eeb71c9f83),
436            (0xbc66f13de6a48d20, 0x3fdb84d7f538c20d),
437            (0x3c5035f4fefd130a, 0x3fb45b63620b1db4),
438            (0x3c1923173b36285f, 0x3f8246e776c12954),
439            (0x3be4c733ad174a92, 0x3f436ec81a059c76),
440            (0x3b8012652ce0bc96, 0x3ef70fe1dc481482),
441            (0x3b2d0670cd6004bc, 0x3e9b497296f5f46e),
442            (0xba85189a8dab3039, 0x3e29bb6d630a41ac),
443            (0x3a3cb2bee18cf77a, 0x3d95f146af7f3514),
444        ];
445        (&P, &Q)
446    } else if x <= 5. {
447        // x <= 5
448        // Polynomial generated by Wolfram Mathematica:
449        // <<FunctionApproximations`
450        // ClearAll["Global`*"]
451        // f[x_]:=PolyGamma[x+1]
452        // {err0,approx}=MiniMaxApproximation[f[z],{z,{3,4},11,11},WorkingPrecision->75,MaxIterations->100]
453        // num=Numerator[approx][[1]];
454        // poly=num;
455        // coeffs=CoefficientList[poly,z];
456        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
457        static P: [(u64, u64); 12] = [
458            (0x3c82777617606035, 0xbfe2788cfc6fb22a),
459            (0x3c6a89168e92a2bd, 0x3fc5286b343ecbe5),
460            (0xbc7562c3f1327b91, 0x3ff8bb7dfd30d653),
461            (0xbc8167bf7208ef96, 0x3ff6f31495ad9560),
462            (0x3c82d2177e694dd6, 0x3fe3431bda031b31),
463            (0xbc6be5da70840171, 0x3fc19291a81dd8a2),
464            (0x3c3bf690d21eb357, 0x3f9258218b3ef7bf),
465            (0xbbffd8d4b0692fbb, 0x3f55d7a5623a2dcf),
466            (0x3babb591eb87dd6c, 0x3f0c56ac18bd4d6b),
467            (0xbb5c8b806eef4ae7, 0x3eb2209eb7acd22b),
468            (0xbad36fc4d27a4168, 0x3e42af070578d594),
469            (0xb9dfe0085849fd5a, 0x3db2d9b6d2913a2a),
470        ];
471        // Generated by Wolfram Mathematica:
472        //<<FunctionApproximations`
473        // ClearAll["Global`*"]
474        // f[x_]:=PolyGamma[x+1]
475        // {err0,approx}=MiniMaxApproximation[f[z],{z,{3,4},11,11},WorkingPrecision->75,MaxIterations->100]
476        // den=Denominator[approx][[1]];
477        // poly=den;
478        // coeffs=CoefficientList[poly,z];
479        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
480        static Q: [(u64, u64); 12] = [
481            (0x0000000000000000, 0x3ff0000000000000),
482            (0xbca9a57355580b15, 0x400481dafdd7bdb2),
483            (0x3c8c54a5f9c29d95, 0x40045b68f0902217),
484            (0xbc795347ec01fbd2, 0x3ff4daef6b912517),
485            (0x3c699ae8f19cf3ec, 0x3fd87d05ac8eb769),
486            (0x3c59ea5ce4904986, 0x3fb12877d62826ca),
487            (0x3be6bd7425bdd967, 0x3f7cce33e2dd0530),
488            (0xbbc460e353f4c8bc, 0x3f3c38bcf76c94b5),
489            (0xbb712ad762f9ec9d, 0x3eee601d38ab6d9b),
490            (0xbb23ce6868791dc5, 0x3e9006213ba059aa),
491            (0xbab95abdc33a7321, 0x3e1a7ac297a0a53d),
492            (0xba174f07c41d8cfd, 0x3d836fa7694a651d),
493        ];
494        (&P, &Q)
495    } else if x <= 6. {
496        // x <= 6
497        // Polynomial generated by Wolfram Mathematica:
498        // <<FunctionApproximations`
499        // ClearAll["Global`*"]
500        // f[x_]:=PolyGamma[x+1]
501        // {err0,approx}=MiniMaxApproximation[f[z],{z,{4,5},11,11},WorkingPrecision->75,MaxIterations->100]
502        // num=Numerator[approx][[1]];
503        // poly=num;
504        // coeffs=CoefficientList[poly,z];
505        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
506        static P: [(u64, u64); 12] = [
507            (0xbc8032a1d2c0ad13, 0xbfe2788cfc6e9dbe),
508            (0x3c38696772d75917, 0x3fc849d273d05d75),
509            (0x3c954ca48d09b1aa, 0x3ff890aef7dc27e1),
510            (0xbc963b4b7c9a58af, 0x3ff5ed5797be10f8),
511            (0x3c87390c0470f3c6, 0x3fe19a75e3a5dbe3),
512            (0xbc4a8fca8acf4f1f, 0x3fbe68afab720ce8),
513            (0x3c1ed20f546b9104, 0x3f8db0bf34ed5e72),
514            (0xbbd6b67d3f6d936e, 0x3f504f9138773832),
515            (0x3bae68dee81a333d, 0x3f033ff5436fb4e3),
516            (0xbb1fad498922c758, 0x3ea61338f9837209),
517            (0x3ad4cdce5341705f, 0x3e341863900b259b),
518            (0xba1368ab344aa52f, 0x3da1a36dc66478ff),
519        ];
520        // Generated by Wolfram Mathematica:
521        //<<FunctionApproximations`
522        // ClearAll["Global`*"]
523        // f[x_]:=PolyGamma[x+1]
524        // {err0,approx}=MiniMaxApproximation[f[z],{z,{4,5},11,11},WorkingPrecision->75,MaxIterations->100]
525        // den=Denominator[approx][[1]];
526        // poly=den;
527        // coeffs=CoefficientList[poly,z];
528        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
529        static Q: [(u64, u64); 12] = [
530            (0x0000000000000000, 0x3ff0000000000000),
531            (0xbc8415c937927b80, 0x40042b14ac1de31a),
532            (0xbca02e650cfc91cb, 0x4003893428b832f0),
533            (0x3c85d3d357beac28, 0x3ff35bb913ab979e),
534            (0x3c763ef88de8bbf4, 0x3fd5c47e7c28af74),
535            (0x3c4624d0ce63afac, 0x3face2a04d07592f),
536            (0x3befca3faad26d95, 0x3f76adbb8f0c9e71),
537            (0x3bd6bb21845ab51f, 0x3f3482259071c038),
538            (0xbb888fad8984adb0, 0x3ee4182fa08710cd),
539            (0xbb2ec0c609c2589a, 0x3e83082e53b39cd4),
540            (0xbaa2954cc8fe5ef2, 0x3e0bd8dfedf43944),
541            (0x39e30166e0cac476, 0x3d71dc88357652c3),
542        ];
543        (&P, &Q)
544    } else if x <= 7. {
545        // x <= 7
546        // Polynomial generated by Wolfram Mathematica:
547        // <<FunctionApproximations`
548        // ClearAll["Global`*"]
549        // f[x_]:=PolyGamma[x+1]
550        // {err0,approx}=MiniMaxApproximation[f[z],{z,{5,6},11,11},WorkingPrecision->75,MaxIterations->100]
551        // num=Numerator[approx][[1]];
552        // poly=num;
553        // coeffs=CoefficientList[poly,z];
554        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
555        static P: [(u64, u64); 12] = [
556            (0xbc81dfd5a5534592, 0xbfe2788cfc565c6c),
557            (0x3c5ee519bbeafb4a, 0x3fcb5e9b249cc2ae),
558            (0x3c8df433a8ab926f, 0x3ff86195ebe61e55),
559            (0xbc9c26e874cedb94, 0x3ff4eeb2d59b9e46),
560            (0xbc8b629f64d30eb4, 0x3fe011547d3b4deb),
561            (0x3c4877eafbf1d0bf, 0x3fba461689f2d09b),
562            (0xbc02cb677bd729e7, 0x3f8801390c64edcd),
563            (0x3be6bb411d2772bf, 0x3f4861722dbdb621),
564            (0x3b95c312dfadf0f8, 0x3efa4589f0866235),
565            (0x3b2e65576c0b0642, 0x3e9b2d2386424cdf),
566            (0x3acbcd0ab8139043, 0x3e260e89319914dd),
567            (0x3a2bde33fa1ac8f3, 0x3d9110236b17c522),
568        ];
569        // Generated by Wolfram Mathematica:
570        //<<FunctionApproximations`
571        // ClearAll["Global`*"]
572        // f[x_]:=PolyGamma[x+1]
573        // {err0,approx}=MiniMaxApproximation[f[z],{z,{5,6},11,11},WorkingPrecision->75,MaxIterations->100]
574        // den=Denominator[approx][[1]];
575        // poly=den;
576        // coeffs=CoefficientList[poly,z];
577        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
578        static Q: [(u64, u64); 12] = [
579            (0x0000000000000000, 0x3ff0000000000000),
580            (0x3c8ff6dd498c91b2, 0x4003d5ac26d888b6),
581            (0x3c9acb90be620fba, 0x4002be9b5e78f13f),
582            (0x3c80944df7c8a780, 0x3ff1f5e5ec0c3d71),
583            (0x3c50f7c451bbdd2e, 0x3fd355e0f87c5218),
584            (0xbc4cd8f8f78b44a1, 0x3fa84c8863aae54b),
585            (0xbc13df2d02c4d3c4, 0x3f71dc00f5ff8e9c),
586            (0x3bc942f7fa51682f, 0x3f2de1134f7b2de6),
587            (0xbb5bcf5da39fb215, 0x3edac2e18af5b8d1),
588            (0xbb11b483e6820b4e, 0x3e76e765628dffe0),
589            (0x3a89bf8734fc1f0e, 0x3dfdf503240dffa2),
590            (0xb9fcaf5858297b35, 0x3d610182e7d328a5),
591        ];
592        (&P, &Q)
593    } else if x <= 8. {
594        // x <= 8
595        // Polynomial generated by Wolfram Mathematica:
596        // <<FunctionApproximations`
597        // ClearAll["Global`*"]
598        // f[x_]:=PolyGamma[x+1]
599        // {err0,approx}=MiniMaxApproximation[f[z],{z,{6,7},11,11},WorkingPrecision->75,MaxIterations->100]
600        // num=Numerator[approx][[1]];
601        // poly=num;
602        // coeffs=CoefficientList[poly,z];
603        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
604        static P: [(u64, u64); 12] = [
605            (0xbc8162f86e7bad36, 0xbfe2788cfb513f92),
606            (0xbc6ddc13d58b297e, 0x3fce652fd3cec73c),
607            (0x3c701fd75f0e2153, 0x3ff82e7ee851190e),
608            (0x3c86494ccd407812, 0x3ff3f7d8f2cc7221),
609            (0x3c71d5a4a3016802, 0x3fdd4df807c7f6a2),
610            (0xbc55bd847484e7bd, 0x3fb6aeaa992b7e5e),
611            (0x3c1d99552c28a292, 0x3f8368db68048c4f),
612            (0xbbecf51acbe60763, 0x3f42437beb7b56ec),
613            (0x3b7ec2e94f896a2d, 0x3ef20a7a9ea0accc),
614            (0x3a9a0a1b60fffd51, 0x3e90f02feda221d2),
615            (0x3ab8a5da0c600246, 0x3e18b9818d3b86a1),
616            (0xba1dc40c982b9669, 0x3d810bda1d0bde80),
617        ];
618        // Generated by Wolfram Mathematica:
619        //<<FunctionApproximations`
620        // ClearAll["Global`*"]
621        // f[x_]:=PolyGamma[x+1]
622        // {err0,approx}=MiniMaxApproximation[f[z],{z,{6,7},11,11},WorkingPrecision->75,MaxIterations->100]
623        // den=Denominator[approx][[1]];
624        // poly=den;
625        // coeffs=CoefficientList[poly,z];
626        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
627        static Q: [(u64, u64); 12] = [
628            (0x0000000000000000, 0x3ff0000000000000),
629            (0x3c95f4382d9f9dbc, 0x400381cd52a27389),
630            (0x3ca36a2a63377659, 0x4001fbd9cbbdde5b),
631            (0xbc83cb5796f6bba0, 0x3ff0a8db1d651315),
632            (0x3c779f72905d08e8, 0x3fd12b8a9f02c16c),
633            (0xbc4dacaf20322afe, 0x3fa47160278b25ee),
634            (0xbc0262fab0a08a3c, 0x3f6c2a94e4e5cf4f),
635            (0x3bc553f072d7f874, 0x3f25d9d7c995af31),
636            (0x3b62f92e30dbacb9, 0x3ed1f75f5d275707),
637            (0x3b02b661ecc2ba69, 0x3e6bf75615fe0e2f),
638            (0x3a571b7290df94bc, 0x3df07d3db0951aef),
639            (0xb9f18f3d33a27e32, 0x3d50bfea8f77e0c3),
640        ];
641        (&P, &Q)
642    } else if x <= 9. {
643        // x <= 9
644        // Polynomial generated by Wolfram Mathematica:
645        // <<FunctionApproximations`
646        // ClearAll["Global`*"]
647        // f[x_]:=PolyGamma[x+1]
648        // {err0,approx}=MiniMaxApproximation[f[z],{z,{7,8},11,11},WorkingPrecision->75,MaxIterations->100]
649        // num=Numerator[approx][[1]];
650        // poly=num;
651        // coeffs=CoefficientList[poly,z];
652        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
653        static P: [(u64, u64); 12] = [
654            (0x3c8efd7fe9db6e79, 0xbfe2788cf4b3267f),
655            (0x3c6a2e8811472031, 0x3fd0ae258336bb2b),
656            (0x3c7399e3844d1e5f, 0x3ff7f7b99c5c43be),
657            (0x3c822ee42e5fe80d, 0x3ff30958a7e8d80e),
658            (0x3c6848d741adf8c4, 0x3fdab4b23acdeb18),
659            (0x3c5aa072e08707aa, 0x3fb393aa5016e08a),
660            (0xbc18b1985dc2c38c, 0x3f7f6b771fa6e163),
661            (0x3bdc86329c365200, 0x3f3b75a61ceee8c8),
662            (0x3b8e8f8ff597b2e6, 0x3ee8f7bf1b272a81),
663            (0x3b19b561962f8b4e, 0x3e8566c7b18e4a14),
664            (0x3aadf125c1ad7adb, 0x3e0c4f51c4444e9a),
665            (0x39f1e1df7a4c3a8b, 0x3d7190106017b110),
666        ];
667        // Generated by Wolfram Mathematica:
668        //<<FunctionApproximations`
669        // ClearAll["Global`*"]
670        // f[x_]:=PolyGamma[x+1]
671        // {err0,approx}=MiniMaxApproximation[f[z],{z,{7,8},11,11},WorkingPrecision->75,MaxIterations->100]
672        // den=Denominator[approx][[1]];
673        // poly=den;
674        // coeffs=CoefficientList[poly,z];
675        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
676        static Q: [(u64, u64); 12] = [
677            (0x0000000000000000, 0x3ff0000000000000),
678            (0x3ca684bdc1b74b83, 0x40032f9b62c83c86),
679            (0x3ca5d8d28ceedbb4, 0x4001410ed7fd4814),
680            (0x3c7686b33443b989, 0x3feee782a9ab5027),
681            (0x3c53620c43b8821f, 0x3fce7f4d1d2ffe3c),
682            (0xbc33e2abe9a40984, 0x3fa1362ff4c8e4cb),
683            (0xbbef2639d9860de0, 0x3f66430cca12747e),
684            (0xbbb8084e63a55a7b, 0x3f200faaa577da43),
685            (0x3b6b0dea347bf9f3, 0x3ec859b423a711e8),
686            (0x3b04cb8df32e9178, 0x3e61563b86bed1a2),
687            (0xba801506c55d0eb4, 0x3de292d18eb8fd2f),
688            (0xb9d06300f4b9226a, 0x3d4109bff966bb51),
689        ];
690        (&P, &Q)
691    } else if x <= 10. {
692        // x <= 10
693        // Polynomial generated by Wolfram Mathematica:
694        // <<FunctionApproximations`
695        // ClearAll["Global`*"]
696        // f[x_]:=PolyGamma[x+1]
697        // {err0,approx}=MiniMaxApproximation[f[z],{z,{8,9},11,11},WorkingPrecision->75,MaxIterations->100]
698        // num=Numerator[approx][[1]];
699        // poly=num;
700        // coeffs=CoefficientList[poly,z];
701        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
702        static P: [(u64, u64); 12] = [
703            (0xbc8593c1df8ef6c7, 0xbfe2788cd6e083bb),
704            (0x3c6415b1331fd1e0, 0x3fd2217df51ae592),
705            (0x3c9c29eeb5fe4731, 0x3ff7bd965beb1b87),
706            (0xbc85fe7985564d6e, 0x3ff2239b9feaa8a7),
707            (0xbc5cc9b9441046be, 0x3fd8541d3ac20fdc),
708            (0xbc243a1d71c1b7e2, 0x3fb0e6b0237a2f28),
709            (0xbc1534c74442766c, 0x3f797a05781bf99b),
710            (0xbbd8ba51ad9881fe, 0x3f34bc78d4108de1),
711            (0xbb81efd29f25726f, 0x3ee16d09b7ddc32f),
712            (0xba96ec9993a947f2, 0x3e7b6c967312dff2),
713            (0xbaa2a0b4acd633ea, 0x3e008d27d632d14c),
714            (0xb9ff596e6d378d22, 0x3d62a21a99d77977),
715        ];
716        // Generated by Wolfram Mathematica:
717        //<<FunctionApproximations`
718        // ClearAll["Global`*"]
719        // f[x_]:=PolyGamma[x+1]
720        // {err0,approx}=MiniMaxApproximation[f[z],{z,{8,9},11,11},WorkingPrecision->75,MaxIterations->100]
721        // den=Denominator[approx][[1]];
722        // poly=den;
723        // coeffs=CoefficientList[poly,z];
724        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
725        static Q: [(u64, u64); 12] = [
726            (0x0000000000000000, 0x3ff0000000000000),
727            (0x3c6f1ef07b1b5dc9, 0x4002df30458e626a),
728            (0xbc960127cb9e2eab, 0x40008e3f3217718d),
729            (0x3c8412d2aabe8bf8, 0x3fecab1dc0e9bcac),
730            (0x3c64edc0f7754efe, 0x3fcb18b1038e24f6),
731            (0xbc3abc16284838a3, 0x3f9d05522b423e2f),
732            (0xbc0e44998e55a65f, 0x3f61a64c0a74479e),
733            (0xbbb190403d2a5a61, 0x3f17c0a8b1842890),
734            (0x3b598fe19ed917ac, 0x3ec0ab7613d6b90e),
735            (0xbaf7c33460b01c25, 0x3e55d581be9c2732),
736            (0xba7278196b41c6f6, 0x3dd565843701949d),
737            (0x39912fb123ad4ce0, 0x3d31de54105a6ae4),
738        ];
739        (&P, &Q)
740    } else if x <= 11. {
741        // x <= 11
742        // Polynomial generated by Wolfram Mathematica:
743        // <<FunctionApproximations`
744        // ClearAll["Global`*"]
745        // f[x_]:=PolyGamma[x+1]
746        // {err0,approx}=MiniMaxApproximation[f[z],{z,{9,10},11,11},WorkingPrecision->75,MaxIterations->100]
747        // num=Numerator[approx][[1]];
748        // poly=num;
749        // coeffs=CoefficientList[poly,z];
750        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
751        static P: [(u64, u64); 12] = [
752            (0xbc7002de428656f8, 0xbfe2788c706086ea),
753            (0xbc7e5dea366e1224, 0x3fd38c4fa29daad3),
754            (0x3c95297f4120d6e1, 0x3ff780642930e686),
755            (0xbc7cdcec6276d962, 0x3ff146e905a85e7c),
756            (0xbc4dde11481dc9d7, 0x3fd629105bd7a86b),
757            (0x3c29f43835822d7e, 0x3fad342a0507d56a),
758            (0x3c1f98d2402273b2, 0x3f74b55e06fd0100),
759            (0x3bc61bf252df59de, 0x3f2f7b95f117d541),
760            (0x3b4f7f3c66a58c10, 0x3ed88c81b5a7bffc),
761            (0xbb1ba95a86f7955a, 0x3e71d2d46e31e04c),
762            (0xb9feae47109db197, 0x3df3c1020acf06fb),
763            (0xb9ed26ed15eff6d5, 0x3d5453226c1a47c8),
764        ];
765        // Generated by Wolfram Mathematica:
766        // <<FunctionApproximations`
767        // ClearAll["Global`*"]
768        // f[x_]:=PolyGamma[x+1]
769        // {err0,approx}=MiniMaxApproximation[f[z],{z,{9,10},11,11},WorkingPrecision->75,MaxIterations->100]
770        // den=Denominator[approx][[1]];
771        // poly=den;
772        // coeffs=CoefficientList[poly,z];
773        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
774        static Q: [(u64, u64); 12] = [
775            (0x0000000000000000, 0x3ff0000000000000),
776            (0xbc8f67d2f15a73c9, 0x4002909d472a6258),
777            (0x3c94a8e3edcab93b, 0x3fffc6b132050914),
778            (0xbc5b982f33a20a73, 0x3fea9a2aba70b06e),
779            (0x3c619d024f5cdba4, 0x3fc817cd43625de8),
780            (0xbc3dc40df7f1ecbb, 0x3f9882a708f4756f),
781            (0x3b88f0ba3ff4834f, 0x3f5c16f3fb9fe3e9),
782            (0xbb8ed7db88e71ba0, 0x3f11ae1da151fd8f),
783            (0x3b3437f8aa73a246, 0x3eb711185a5568ff),
784            (0x3ab876889e88ff72, 0x3e4bee89a4622cb4),
785            (0xba6ac2419dbe0e2f, 0x3dc930ba3560eb64),
786            (0xb9b2124aaf27b632, 0x3d234925cfe9d5d8),
787        ];
788        (&P, &Q)
789    } else {
790        // x <= 12
791        // Polynomial generated by Wolfram Mathematica:
792        // <<FunctionApproximations`
793        // ClearAll["Global`*"]
794        // f[x_]:=PolyGamma[x+1]
795        // {err0,approx}=MiniMaxApproximation[f[z],{z,{10,11},11,11},WorkingPrecision->75,MaxIterations->100]
796        // num=Numerator[approx][[1]];
797        // poly=num;
798        // coeffs=CoefficientList[poly,z];
799        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
800        static P: [(u64, u64); 12] = [
801            (0xbc88aa87989ef5a1, 0xbfe2788b51df7906),
802            (0xbc7af1e1932b591c, 0x3fd4ee69b56504f9),
803            (0x3c740a2c04d19bfc, 0x3ff7406f9cf46140),
804            (0xbc917572eb9329ef, 0x3ff07369547283f5),
805            (0x3c54e46eaf6429e3, 0x3fd4301a4c091a22),
806            (0x3c3c3ba8aa15ca20, 0x3fa9426e448ee7db),
807            (0xbc1f041777b9f035, 0x3f70e2032a465789),
808            (0xbbb783cb47e19f83, 0x3f28099ad5f73d16),
809            (0xbb6932c40b32e5e0, 0x3ed175058b1e67e6),
810            (0x3b024c00ae6427d0, 0x3e677fd917aaa847),
811            (0x3a50b5029ff4eadf, 0x3de80bf02bf0e221),
812            (0x39e405e48a7d9326, 0x3d46c1a62ef66417),
813        ];
814        // Generated by Wolfram Mathematica:
815        // <<FunctionApproximations`
816        // ClearAll["Global`*"]
817        // f[x_]:=PolyGamma[x+1]
818        // {err0,approx}=MiniMaxApproximation[f[z],{z,{10,11},11,11},WorkingPrecision->75,MaxIterations->100]
819        // den=Denominator[approx][[1]];
820        // poly=den;
821        // coeffs=CoefficientList[poly,z];
822        // TableForm[Table[Row[{"'",NumberForm[coeffs[[i+1]],{50,50}, ExponentFunction->(Null&)],"',"}],{i,0,Length[coeffs]-1}]]
823        static Q: [(u64, u64); 12] = [
824            (0x0000000000000000, 0x3ff0000000000000),
825            (0x3c67856052188c26, 0x400243ec3166f349),
826            (0x3c2b34f0f808ddac, 0x3ffe806cb15a71b3),
827            (0xbc823b5d5bc54e38, 0x3fe8b21e296d1320),
828            (0xbc52f7423180111c, 0x3fc571b5d26833f9),
829            (0xbc33e19defee24d6, 0x3f94bf1ad662d38a),
830            (0x3beb180e8d55583a, 0x3f56723080074543),
831            (0xbb7d192d116623ad, 0x3f0a818d001404be),
832            (0x3b194b9889dea9d1, 0x3eb021f1a8b7f978),
833            (0x3a8ca7cbf623a001, 0x3e4224e2c71fae88),
834            (0xba4d8314f1a91e06, 0x3dbe48212e08d6f2),
835            (0xb9a210925f26ebb8, 0x3d15627ecb7e02ac),
836        ];
837        (&P, &Q)
838    }
839}
pxfm/gamma/digamma_coeffs.rs

pxfm/gamma/
digamma_coeffs.rs
