Bug Summary

File:root/firefox-clang/modules/fdlibm/src/k_rem_pio2.cpp
Warning:line 370, column 19
The left operand of '==' is a garbage value due to array index out of bounds

Annotated Source Code

Press '?' to see keyboard shortcuts

clang -cc1 -cc1 -triple x86_64-pc-linux-gnu -O2 -analyze -disable-free -clear-ast-before-backend -disable-llvm-verifier -discard-value-names -main-file-name k_rem_pio2.cpp -analyzer-checker=core -analyzer-checker=apiModeling -analyzer-checker=unix -analyzer-checker=deadcode -analyzer-checker=cplusplus -analyzer-checker=security.insecureAPI.UncheckedReturn -analyzer-checker=security.insecureAPI.getpw -analyzer-checker=security.insecureAPI.gets -analyzer-checker=security.insecureAPI.mktemp -analyzer-checker=security.insecureAPI.mkstemp -analyzer-checker=security.insecureAPI.vfork -analyzer-checker=nullability.NullPassedToNonnull -analyzer-checker=nullability.NullReturnedFromNonnull -analyzer-output plist -w -setup-static-analyzer -analyzer-config-compatibility-mode=true -mrelocation-model pic -pic-level 2 -fhalf-no-semantic-interposition -mframe-pointer=all -relaxed-aliasing -ffp-contract=off -fno-rounding-math -mconstructor-aliases -funwind-tables=2 -target-cpu x86-64 -tune-cpu generic -debugger-tuning=gdb -fdebug-compilation-dir=/root/firefox-clang/obj-x86_64-pc-linux-gnu/modules/fdlibm/src -fcoverage-compilation-dir=/root/firefox-clang/obj-x86_64-pc-linux-gnu/modules/fdlibm/src -resource-dir /usr/lib/llvm-23/lib/clang/23 -include /root/firefox-clang/config/gcc_hidden.h -include /root/firefox-clang/obj-x86_64-pc-linux-gnu/mozilla-config.h -I /root/firefox-clang/obj-x86_64-pc-linux-gnu/dist/stl_wrappers -D _GLIBCXX_ASSERTIONS=1 -I /root/firefox-clang/obj-x86_64-pc-linux-gnu/dist/system_wrappers -U _FORTIFY_SOURCE -D _FORTIFY_SOURCE=2 -D DEBUG=1 -I /root/firefox-clang/modules/fdlibm/src -I /root/firefox-clang/obj-x86_64-pc-linux-gnu/modules/fdlibm/src -I /root/firefox-clang/obj-x86_64-pc-linux-gnu/dist/include -I /root/firefox-clang/obj-x86_64-pc-linux-gnu/dist/include/nspr -I /root/firefox-clang/obj-x86_64-pc-linux-gnu/dist/include/nss -D MOZILLA_CLIENT -internal-isystem /usr/lib/gcc/x86_64-linux-gnu/16/../../../../include/c++/16 -internal-isystem /usr/lib/gcc/x86_64-linux-gnu/16/../../../../include/x86_64-linux-gnu/c++/16 -internal-isystem /usr/lib/gcc/x86_64-linux-gnu/16/../../../../include/c++/16/backward -internal-isystem /usr/lib/llvm-23/lib/clang/23/include -internal-isystem /usr/local/include -internal-isystem /usr/lib/gcc/x86_64-linux-gnu/16/../../../../x86_64-linux-gnu/include -internal-externc-isystem /usr/include/x86_64-linux-gnu -internal-externc-isystem /include -internal-externc-isystem /usr/include -Wno-error=pessimizing-move -Wno-error=large-by-value-copy=128 -Wno-error=implicit-int-float-conversion -Wno-error=thread-safety-analysis -Wno-error=tautological-type-limit-compare -Wno-invalid-offsetof -Wno-range-loop-analysis -Wno-deprecated-anon-enum-enum-conversion -Wno-deprecated-enum-enum-conversion -Wno-inline-new-delete -Wno-error=deprecated-declarations -Wno-error=array-bounds -Wno-error=free-nonheap-object -Wno-error=atomic-alignment -Wno-error=deprecated-builtins -Wno-psabi -Wno-error=builtin-macro-redefined -Wno-vla-cxx-extension -Wno-unknown-warning-option -Wno-character-conversion -Wno-parentheses -Wno-sign-compare -Wno-dangling-else -std=gnu++20 -fdeprecated-macro -ferror-limit 19 -fstrict-flex-arrays=1 -stack-protector 2 -fstack-clash-protection -ftrivial-auto-var-init=pattern -fno-rtti -fgnuc-version=4.2.1 -fno-implicit-modules -fskip-odr-check-in-gmf -fno-sized-deallocation -fno-aligned-allocation -fdiagnostics-absolute-paths -vectorize-loops -vectorize-slp -analyzer-checker optin.performance.Padding -analyzer-output=html -analyzer-config stable-report-filename=true -mllvm -dwarf-linkage-names=Abstract -faddrsig -fdwarf2-cfi-asm -o /tmp/scan-build-2026-09-01-224014-2642839-1 -x c++ /root/firefox-clang/modules/fdlibm/src/k_rem_pio2.cpp
1
2/* @(#)k_rem_pio2.c 1.3 95/01/18 */
3/*
4 * ====================================================
5 * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
6 *
7 * Developed at SunSoft, a Sun Microsystems, Inc. business.
8 * Permission to use, copy, modify, and distribute this
9 * software is freely granted, provided that this notice
10 * is preserved.
11 * ====================================================
12 */
13
14//#include <sys/cdefs.h>
15//__FBSDID("$FreeBSD$");
16
17/*
18 * __kernel_rem_pio2(x,y,e0,nx,prec)
19 * double x[],y[]; int e0,nx,prec;
20 *
21 * __kernel_rem_pio2 return the last three digits of N with
22 * y = x - N*pi/2
23 * so that |y| < pi/2.
24 *
25 * The method is to compute the integer (mod 8) and fraction parts of
26 * (2/pi)*x without doing the full multiplication. In general we
27 * skip the part of the product that are known to be a huge integer (
28 * more accurately, = 0 mod 8 ). Thus the number of operations are
29 * independent of the exponent of the input.
30 *
31 * (2/pi) is represented by an array of 24-bit integers in ipio2[].
32 *
33 * Input parameters:
34 * x[] The input value (must be positive) is broken into nx
35 * pieces of 24-bit integers in double precision format.
36 * x[i] will be the i-th 24 bit of x. The scaled exponent
37 * of x[0] is given in input parameter e0 (i.e., x[0]*2^e0
38 * match x's up to 24 bits.
39 *
40 * Example of breaking a double positive z into x[0]+x[1]+x[2]:
41 * e0 = ilogb(z)-23
42 * z = scalbn(z,-e0)
43 * for i = 0,1,2
44 * x[i] = floor(z)
45 * z = (z-x[i])*2**24
46 *
47 *
48 * y[] output result in an array of double precision numbers.
49 * The dimension of y[] is:
50 * 24-bit precision 1
51 * 53-bit precision 2
52 * 64-bit precision 2
53 * 113-bit precision 3
54 * The actual value is the sum of them. Thus for 113-bit
55 * precision, one may have to do something like:
56 *
57 * long double t,w,r_head, r_tail;
58 * t = (long double)y[2] + (long double)y[1];
59 * w = (long double)y[0];
60 * r_head = t+w;
61 * r_tail = w - (r_head - t);
62 *
63 * e0 The exponent of x[0]. Must be <= 16360 or you need to
64 * expand the ipio2 table.
65 *
66 * nx dimension of x[]
67 *
68 * prec an integer indicating the precision:
69 * 0 24 bits (single)
70 * 1 53 bits (double)
71 * 2 64 bits (extended)
72 * 3 113 bits (quad)
73 *
74 * External function:
75 * double scalbn(), floor();
76 *
77 *
78 * Here is the description of some local variables:
79 *
80 * jk jk+1 is the initial number of terms of ipio2[] needed
81 * in the computation. The minimum and recommended value
82 * for jk is 3,4,4,6 for single, double, extended, and quad.
83 * jk+1 must be 2 larger than you might expect so that our
84 * recomputation test works. (Up to 24 bits in the integer
85 * part (the 24 bits of it that we compute) and 23 bits in
86 * the fraction part may be lost to cancellation before we
87 * recompute.)
88 *
89 * jz local integer variable indicating the number of
90 * terms of ipio2[] used.
91 *
92 * jx nx - 1
93 *
94 * jv index for pointing to the suitable ipio2[] for the
95 * computation. In general, we want
96 * ( 2^e0*x[0] * ipio2[jv-1]*2^(-24jv) )/8
97 * is an integer. Thus
98 * e0-3-24*jv >= 0 or (e0-3)/24 >= jv
99 * Hence jv = max(0,(e0-3)/24).
100 *
101 * jp jp+1 is the number of terms in PIo2[] needed, jp = jk.
102 *
103 * q[] double array with integral value, representing the
104 * 24-bits chunk of the product of x and 2/pi.
105 *
106 * q0 the corresponding exponent of q[0]. Note that the
107 * exponent for q[i] would be q0-24*i.
108 *
109 * PIo2[] double precision array, obtained by cutting pi/2
110 * into 24 bits chunks.
111 *
112 * f[] ipio2[] in floating point
113 *
114 * iq[] integer array by breaking up q[] in 24-bits chunk.
115 *
116 * fq[] final product of x*(2/pi) in fq[0],..,fq[jk]
117 *
118 * ih integer. If >0 it indicates q[] is >= 0.5, hence
119 * it also indicates the *sign* of the result.
120 *
121 */
122
123
124/*
125 * Constants:
126 * The hexadecimal values are the intended ones for the following
127 * constants. The decimal values may be used, provided that the
128 * compiler will convert from decimal to binary accurately enough
129 * to produce the hexadecimal values shown.
130 */
131
132#include <float.h>
133#include <math.h>
134
135#include "math_private.h"
136
137static const int init_jk[] = {3,4,4,6}; /* initial value for jk */
138
139/*
140 * Table of constants for 2/pi, 396 Hex digits (476 decimal) of 2/pi
141 *
142 * integer array, contains the (24*i)-th to (24*i+23)-th
143 * bit of 2/pi after binary point. The corresponding
144 * floating value is
145 *
146 * ipio2[i] * 2^(-24(i+1)).
147 *
148 * NB: This table must have at least (e0-3)/24 + jk terms.
149 * For quad precision (e0 <= 16360, jk = 6), this is 686.
150 */
151static const int32_t ipio2[] = {
1520xA2F983, 0x6E4E44, 0x1529FC, 0x2757D1, 0xF534DD, 0xC0DB62,
1530x95993C, 0x439041, 0xFE5163, 0xABDEBB, 0xC561B7, 0x246E3A,
1540x424DD2, 0xE00649, 0x2EEA09, 0xD1921C, 0xFE1DEB, 0x1CB129,
1550xA73EE8, 0x8235F5, 0x2EBB44, 0x84E99C, 0x7026B4, 0x5F7E41,
1560x3991D6, 0x398353, 0x39F49C, 0x845F8B, 0xBDF928, 0x3B1FF8,
1570x97FFDE, 0x05980F, 0xEF2F11, 0x8B5A0A, 0x6D1F6D, 0x367ECF,
1580x27CB09, 0xB74F46, 0x3F669E, 0x5FEA2D, 0x7527BA, 0xC7EBE5,
1590xF17B3D, 0x0739F7, 0x8A5292, 0xEA6BFB, 0x5FB11F, 0x8D5D08,
1600x560330, 0x46FC7B, 0x6BABF0, 0xCFBC20, 0x9AF436, 0x1DA9E3,
1610x91615E, 0xE61B08, 0x659985, 0x5F14A0, 0x68408D, 0xFFD880,
1620x4D7327, 0x310606, 0x1556CA, 0x73A8C9, 0x60E27B, 0xC08C6B,
163
164#if LDBL_MAX_EXP16384 > 1024
165#if LDBL_MAX_EXP16384 > 16384
166#error "ipio2 table needs to be expanded"
167#endif
1680x47C419, 0xC367CD, 0xDCE809, 0x2A8359, 0xC4768B, 0x961CA6,
1690xDDAF44, 0xD15719, 0x053EA5, 0xFF0705, 0x3F7E33, 0xE832C2,
1700xDE4F98, 0x327DBB, 0xC33D26, 0xEF6B1E, 0x5EF89F, 0x3A1F35,
1710xCAF27F, 0x1D87F1, 0x21907C, 0x7C246A, 0xFA6ED5, 0x772D30,
1720x433B15, 0xC614B5, 0x9D19C3, 0xC2C4AD, 0x414D2C, 0x5D000C,
1730x467D86, 0x2D71E3, 0x9AC69B, 0x006233, 0x7CD2B4, 0x97A7B4,
1740xD55537, 0xF63ED7, 0x1810A3, 0xFC764D, 0x2A9D64, 0xABD770,
1750xF87C63, 0x57B07A, 0xE71517, 0x5649C0, 0xD9D63B, 0x3884A7,
1760xCB2324, 0x778AD6, 0x23545A, 0xB91F00, 0x1B0AF1, 0xDFCE19,
1770xFF319F, 0x6A1E66, 0x615799, 0x47FBAC, 0xD87F7E, 0xB76522,
1780x89E832, 0x60BFE6, 0xCDC4EF, 0x09366C, 0xD43F5D, 0xD7DE16,
1790xDE3B58, 0x929BDE, 0x2822D2, 0xE88628, 0x4D58E2, 0x32CAC6,
1800x16E308, 0xCB7DE0, 0x50C017, 0xA71DF3, 0x5BE018, 0x34132E,
1810x621283, 0x014883, 0x5B8EF5, 0x7FB0AD, 0xF2E91E, 0x434A48,
1820xD36710, 0xD8DDAA, 0x425FAE, 0xCE616A, 0xA4280A, 0xB499D3,
1830xF2A606, 0x7F775C, 0x83C2A3, 0x883C61, 0x78738A, 0x5A8CAF,
1840xBDD76F, 0x63A62D, 0xCBBFF4, 0xEF818D, 0x67C126, 0x45CA55,
1850x36D9CA, 0xD2A828, 0x8D61C2, 0x77C912, 0x142604, 0x9B4612,
1860xC459C4, 0x44C5C8, 0x91B24D, 0xF31700, 0xAD43D4, 0xE54929,
1870x10D5FD, 0xFCBE00, 0xCC941E, 0xEECE70, 0xF53E13, 0x80F1EC,
1880xC3E7B3, 0x28F8C7, 0x940593, 0x3E71C1, 0xB3092E, 0xF3450B,
1890x9C1288, 0x7B20AB, 0x9FB52E, 0xC29247, 0x2F327B, 0x6D550C,
1900x90A772, 0x1FE76B, 0x96CB31, 0x4A1679, 0xE27941, 0x89DFF4,
1910x9794E8, 0x84E6E2, 0x973199, 0x6BED88, 0x365F5F, 0x0EFDBB,
1920xB49A48, 0x6CA467, 0x427271, 0x325D8D, 0xB8159F, 0x09E5BC,
1930x25318D, 0x3974F7, 0x1C0530, 0x010C0D, 0x68084B, 0x58EE2C,
1940x90AA47, 0x02E774, 0x24D6BD, 0xA67DF7, 0x72486E, 0xEF169F,
1950xA6948E, 0xF691B4, 0x5153D1, 0xF20ACF, 0x339820, 0x7E4BF5,
1960x6863B2, 0x5F3EDD, 0x035D40, 0x7F8985, 0x295255, 0xC06437,
1970x10D86D, 0x324832, 0x754C5B, 0xD4714E, 0x6E5445, 0xC1090B,
1980x69F52A, 0xD56614, 0x9D0727, 0x50045D, 0xDB3BB4, 0xC576EA,
1990x17F987, 0x7D6B49, 0xBA271D, 0x296996, 0xACCCC6, 0x5414AD,
2000x6AE290, 0x89D988, 0x50722C, 0xBEA404, 0x940777, 0x7030F3,
2010x27FC00, 0xA871EA, 0x49C266, 0x3DE064, 0x83DD97, 0x973FA3,
2020xFD9443, 0x8C860D, 0xDE4131, 0x9D3992, 0x8C70DD, 0xE7B717,
2030x3BDF08, 0x2B3715, 0xA0805C, 0x93805A, 0x921110, 0xD8E80F,
2040xAF806C, 0x4BFFDB, 0x0F9038, 0x761859, 0x15A562, 0xBBCB61,
2050xB989C7, 0xBD4010, 0x04F2D2, 0x277549, 0xF6B6EB, 0xBB22DB,
2060xAA140A, 0x2F2689, 0x768364, 0x333B09, 0x1A940E, 0xAA3A51,
2070xC2A31D, 0xAEEDAF, 0x12265C, 0x4DC26D, 0x9C7A2D, 0x9756C0,
2080x833F03, 0xF6F009, 0x8C402B, 0x99316D, 0x07B439, 0x15200C,
2090x5BC3D8, 0xC492F5, 0x4BADC6, 0xA5CA4E, 0xCD37A7, 0x36A9E6,
2100x9492AB, 0x6842DD, 0xDE6319, 0xEF8C76, 0x528B68, 0x37DBFC,
2110xABA1AE, 0x3115DF, 0xA1AE00, 0xDAFB0C, 0x664D64, 0xB705ED,
2120x306529, 0xBF5657, 0x3AFF47, 0xB9F96A, 0xF3BE75, 0xDF9328,
2130x3080AB, 0xF68C66, 0x15CB04, 0x0622FA, 0x1DE4D9, 0xA4B33D,
2140x8F1B57, 0x09CD36, 0xE9424E, 0xA4BE13, 0xB52333, 0x1AAAF0,
2150xA8654F, 0xA5C1D2, 0x0F3F0B, 0xCD785B, 0x76F923, 0x048B7B,
2160x721789, 0x53A6C6, 0xE26E6F, 0x00EBEF, 0x584A9B, 0xB7DAC4,
2170xBA66AA, 0xCFCF76, 0x1D02D1, 0x2DF1B1, 0xC1998C, 0x77ADC3,
2180xDA4886, 0xA05DF7, 0xF480C6, 0x2FF0AC, 0x9AECDD, 0xBC5C3F,
2190x6DDED0, 0x1FC790, 0xB6DB2A, 0x3A25A3, 0x9AAF00, 0x9353AD,
2200x0457B6, 0xB42D29, 0x7E804B, 0xA707DA, 0x0EAA76, 0xA1597B,
2210x2A1216, 0x2DB7DC, 0xFDE5FA, 0xFEDB89, 0xFDBE89, 0x6C76E4,
2220xFCA906, 0x70803E, 0x156E85, 0xFF87FD, 0x073E28, 0x336761,
2230x86182A, 0xEABD4D, 0xAFE7B3, 0x6E6D8F, 0x396795, 0x5BBF31,
2240x48D784, 0x16DF30, 0x432DC7, 0x356125, 0xCE70C9, 0xB8CB30,
2250xFD6CBF, 0xA200A4, 0xE46C05, 0xA0DD5A, 0x476F21, 0xD21262,
2260x845CB9, 0x496170, 0xE0566B, 0x015299, 0x375550, 0xB7D51E,
2270xC4F133, 0x5F6E13, 0xE4305D, 0xA92E85, 0xC3B21D, 0x3632A1,
2280xA4B708, 0xD4B1EA, 0x21F716, 0xE4698F, 0x77FF27, 0x80030C,
2290x2D408D, 0xA0CD4F, 0x99A520, 0xD3A2B3, 0x0A5D2F, 0x42F9B4,
2300xCBDA11, 0xD0BE7D, 0xC1DB9B, 0xBD17AB, 0x81A2CA, 0x5C6A08,
2310x17552E, 0x550027, 0xF0147F, 0x8607E1, 0x640B14, 0x8D4196,
2320xDEBE87, 0x2AFDDA, 0xB6256B, 0x34897B, 0xFEF305, 0x9EBFB9,
2330x4F6A68, 0xA82A4A, 0x5AC44F, 0xBCF82D, 0x985AD7, 0x95C7F4,
2340x8D4D0D, 0xA63A20, 0x5F57A4, 0xB13F14, 0x953880, 0x0120CC,
2350x86DD71, 0xB6DEC9, 0xF560BF, 0x11654D, 0x6B0701, 0xACB08C,
2360xD0C0B2, 0x485551, 0x0EFB1E, 0xC37295, 0x3B06A3, 0x3540C0,
2370x7BDC06, 0xCC45E0, 0xFA294E, 0xC8CAD6, 0x41F3E8, 0xDE647C,
2380xD8649B, 0x31BED9, 0xC397A4, 0xD45877, 0xC5E369, 0x13DAF0,
2390x3C3ABA, 0x461846, 0x5F7555, 0xF5BDD2, 0xC6926E, 0x5D2EAC,
2400xED440E, 0x423E1C, 0x87C461, 0xE9FD29, 0xF3D6E7, 0xCA7C22,
2410x35916F, 0xC5E008, 0x8DD7FF, 0xE26A6E, 0xC6FDB0, 0xC10893,
2420x745D7C, 0xB2AD6B, 0x9D6ECD, 0x7B723E, 0x6A11C6, 0xA9CFF7,
2430xDF7329, 0xBAC9B5, 0x5100B7, 0x0DB2E2, 0x24BA74, 0x607DE5,
2440x8AD874, 0x2C150D, 0x0C1881, 0x94667E, 0x162901, 0x767A9F,
2450xBEFDFD, 0xEF4556, 0x367ED9, 0x13D9EC, 0xB9BA8B, 0xFC97C4,
2460x27A831, 0xC36EF1, 0x36C594, 0x56A8D8, 0xB5A8B4, 0x0ECCCF,
2470x2D8912, 0x34576F, 0x89562C, 0xE3CE99, 0xB920D6, 0xAA5E6B,
2480x9C2A3E, 0xCC5F11, 0x4A0BFD, 0xFBF4E1, 0x6D3B8E, 0x2C86E2,
2490x84D4E9, 0xA9B4FC, 0xD1EEEF, 0xC9352E, 0x61392F, 0x442138,
2500xC8D91B, 0x0AFC81, 0x6A4AFB, 0xD81C2F, 0x84B453, 0x8C994E,
2510xCC2254, 0xDC552A, 0xD6C6C0, 0x96190B, 0xB8701A, 0x649569,
2520x605A26, 0xEE523F, 0x0F117F, 0x11B5F4, 0xF5CBFC, 0x2DBC34,
2530xEEBC34, 0xCC5DE8, 0x605EDD, 0x9B8E67, 0xEF3392, 0xB817C9,
2540x9B5861, 0xBC57E1, 0xC68351, 0x103ED8, 0x4871DD, 0xDD1C2D,
2550xA118AF, 0x462C21, 0xD7F359, 0x987AD9, 0xC0549E, 0xFA864F,
2560xFC0656, 0xAE79E5, 0x362289, 0x22AD38, 0xDC9367, 0xAAE855,
2570x382682, 0x9BE7CA, 0xA40D51, 0xB13399, 0x0ED7A9, 0x480569,
2580xF0B265, 0xA7887F, 0x974C88, 0x36D1F9, 0xB39221, 0x4A827B,
2590x21CF98, 0xDC9F40, 0x5547DC, 0x3A74E1, 0x42EB67, 0xDF9DFE,
2600x5FD45E, 0xA4677B, 0x7AACBA, 0xA2F655, 0x23882B, 0x55BA41,
2610x086E59, 0x862A21, 0x834739, 0xE6E389, 0xD49EE5, 0x40FB49,
2620xE956FF, 0xCA0F1C, 0x8A59C5, 0x2BFA94, 0xC5C1D3, 0xCFC50F,
2630xAE5ADB, 0x86C547, 0x624385, 0x3B8621, 0x94792C, 0x876110,
2640x7B4C2A, 0x1A2C80, 0x12BF43, 0x902688, 0x893C78, 0xE4C4A8,
2650x7BDBE5, 0xC23AC4, 0xEAF426, 0x8A67F7, 0xBF920D, 0x2BA365,
2660xB1933D, 0x0B7CBD, 0xDC51A4, 0x63DD27, 0xDDE169, 0x19949A,
2670x9529A8, 0x28CE68, 0xB4ED09, 0x209F44, 0xCA984E, 0x638270,
2680x237C7E, 0x32B90F, 0x8EF5A7, 0xE75614, 0x08F121, 0x2A9DB5,
2690x4D7E6F, 0x5119A5, 0xABF9B5, 0xD6DF82, 0x61DD96, 0x023616,
2700x9F3AC4, 0xA1A283, 0x6DED72, 0x7A8D39, 0xA9B882, 0x5C326B,
2710x5B2746, 0xED3400, 0x7700D2, 0x55F4FC, 0x4D5901, 0x8071E0,
272#endif
273
274};
275
276static const double PIo2[] = {
277 1.57079625129699707031e+00, /* 0x3FF921FB, 0x40000000 */
278 7.54978941586159635335e-08, /* 0x3E74442D, 0x00000000 */
279 5.39030252995776476554e-15, /* 0x3CF84698, 0x80000000 */
280 3.28200341580791294123e-22, /* 0x3B78CC51, 0x60000000 */
281 1.27065575308067607349e-29, /* 0x39F01B83, 0x80000000 */
282 1.22933308981111328932e-36, /* 0x387A2520, 0x40000000 */
283 2.73370053816464559624e-44, /* 0x36E38222, 0x80000000 */
284 2.16741683877804819444e-51, /* 0x3569F31D, 0x00000000 */
285};
286
287static const double
288zero = 0.0,
289one = 1.0,
290two24 = 1.67772160000000000000e+07, /* 0x41700000, 0x00000000 */
291twon24 = 5.96046447753906250000e-08; /* 0x3E700000, 0x00000000 */
292
293int
294__kernel_rem_pio2(double *x, double *y, int e0, int nx, int prec)
295{
296 int32_t jz,jx,jv,jp,jk,carry,n,iq[20],i,j,k,m,q0,ih;
297 double z,fw,f[20],fq[20],q[20];
298
299 /* initialize jk*/
300 jk = init_jk[prec];
301 jp = jk;
302
303 /* determine jx,jv,q0, note that 3>q0 */
304 jx = nx-1;
305 jv = (e0-3)/24; if(jv<0) jv=0;
1
Assuming 'jv' is >= 0
→
2
←
Taking false branch
→
306 q0 = e0-24*(jv+1);
307
308 /* set up f[0] to f[jx+jk] where f[jx+jk] = ipio2[jv+jk] */
309 j = jv-jx; m = jx+jk;
310 for(i=0;i<=m;i++,j++) f[i] = (j<0)? zero : (double) ipio2[j];
3
←
Assuming 'i' is > 'm'
→
4
←
Loop condition is false. Execution continues on line 313
→
311
312 /* compute q[0],q[1],...q[jk] */
313 for (i=0;i<=jk;i++) {
5
←
Assuming 'i' is <= 'jk'
→
6
←
Loop condition is true. Entering loop body
→
9
←
Assuming 'i' is > 'jk'
→
10
←
Loop condition is false. Execution continues on line 318
→
314 for(j=0,fw=0.0;j<=jx;j++) fw += x[j]*f[jx+i-j];
7
←
Assuming 'j' is > 'jx'
→
8
←
Loop condition is false. Execution continues on line 315
→
315 q[i] = fw;
316 }
317
318 jz = jk;
319recompute:
320 /* distill q[] into iq[] reversingly */
321 for(i=0,j=jz,z=q[jz];j>0;i++,j--) {
11
←
Loop condition is false. Execution continues on line 328
→
322 fw = (double)((int32_t)(twon24* z));
323 iq[i] = (int32_t)(z-two24*fw);
324 z = q[j-1]+fw;
325 }
326
327 /* compute n */
328 z = scalbn(z,q0); /* actual value of z */
329 z -= 8.0*floor(z*0.125); /* trim off integer >= 8 */
330 n = (int32_t) z;
331 z -= (double)n;
332 ih = 0;
333 if(q0>0) { /* need iq[jz-1] to determine n */
12
←
Assuming 'q0' is <= 0
→
13
←
Taking false branch
→
334 i = (iq[jz-1]>>(24-q0)); n += i;
335 iq[jz-1] -= i<<(24-q0);
336 ih = iq[jz-1]>>(23-q0);
337 }
338 else if(q0==0) ih = iq[jz-1]>>23;
14
←
Assuming 'q0' is not equal to 0
→
15
←
Taking false branch
→
339 else if(z>=0.5) ih=2;
16
←
Assuming the condition is false
→
17
←
Taking false branch
→
340
341 if(ih
17.1
'ih' is <= 0
>0) { /* q > 0.5 */
18
←
Taking false branch
→
342 n += 1; carry = 0;
343 for(i=0;i<jz ;i++) { /* compute 1-q */
344 j = iq[i];
345 if(carry==0) {
346 if(j!=0) {
347 carry = 1; iq[i] = 0x1000000- j;
348 }
349 } else iq[i] = 0xffffff - j;
350 }
351 if(q0>0) { /* rare case: chance is 1 in 12 */
352 switch(q0) {
353 case 1:
354 iq[jz-1] &= 0x7fffff; break;
355 case 2:
356 iq[jz-1] &= 0x3fffff; break;
357 }
358 }
359 if(ih==2) {
360 z = one - z;
361 if(carry!=0) z -= scalbn(one,q0);
362 }
363 }
364
365 /* check if recomputation is needed */
366 if(z==zero) {
19
←
Assuming 'z' is equal to 'zero'
→
20
←
Taking true branch
→
367 j = 0;
368 for (i=jz-1;i>=jk;i--) j |= iq[i];
21
←
Loop condition is false. Execution continues on line 369
→
369 if(j
21.1
'j' is equal to 0
==0) { /* need recomputation */
22
←
Taking true branch
→
370 for(k=1;iq[jk-k]==0;k++); /* k = no. of terms needed */
23
←
The left operand of '==' is a garbage value due to array index out of bounds
371
372 for(i=jz+1;i<=jz+k;i++) { /* add q[jz+1] to q[jz+k] */
373 f[jx+i] = (double) ipio2[jv+i];
374 for(j=0,fw=0.0;j<=jx;j++) fw += x[j]*f[jx+i-j];
375 q[i] = fw;
376 }
377 jz += k;
378 goto recompute;
379 }
380 }
381
382 /* chop off zero terms */
383 if(z==0.0) {
384 jz -= 1; q0 -= 24;
385 while(iq[jz]==0) { jz--; q0-=24;}
386 } else { /* break z into 24-bit if necessary */
387 z = scalbn(z,-q0);
388 if(z>=two24) {
389 fw = (double)((int32_t)(twon24*z));
390 iq[jz] = (int32_t)(z-two24*fw);
391 jz += 1; q0 += 24;
392 iq[jz] = (int32_t) fw;
393 } else iq[jz] = (int32_t) z ;
394 }
395
396 /* convert integer "bit" chunk to floating-point value */
397 fw = scalbn(one,q0);
398 for(i=jz;i>=0;i--) {
399 q[i] = fw*(double)iq[i]; fw*=twon24;
400 }
401
402 /* compute PIo2[0,...,jp]*q[jz,...,0] */
403 for(i=jz;i>=0;i--) {
404 for(fw=0.0,k=0;k<=jp&&k<=jz-i;k++) fw += PIo2[k]*q[i+k];
405 fq[jz-i] = fw;
406 }
407
408 /* compress fq[] into y[] */
409 switch(prec) {
410 case 0:
411 fw = 0.0;
412 for (i=jz;i>=0;i--) fw += fq[i];
413 y[0] = (ih==0)? fw: -fw;
414 break;
415 case 1:
416 case 2:
417 fw = 0.0;
418 for (i=jz;i>=0;i--) fw += fq[i];
419 STRICT_ASSIGN(double,fw,fw)((fw) = (fw));
420 y[0] = (ih==0)? fw: -fw;
421 fw = fq[0]-fw;
422 for (i=1;i<=jz;i++) fw += fq[i];
423 y[1] = (ih==0)? fw: -fw;
424 break;
425 case 3: /* painful */
426 for (i=jz;i>0;i--) {
427 fw = fq[i-1]+fq[i];
428 fq[i] += fq[i-1]-fw;
429 fq[i-1] = fw;
430 }
431 for (i=jz;i>1;i--) {
432 fw = fq[i-1]+fq[i];
433 fq[i] += fq[i-1]-fw;
434 fq[i-1] = fw;
435 }
436 for (fw=0.0,i=jz;i>=2;i--) fw += fq[i];
437 if(ih==0) {
438 y[0] = fq[0]; y[1] = fq[1]; y[2] = fw;
439 } else {
440 y[0] = -fq[0]; y[1] = -fq[1]; y[2] = -fw;
441 }
442 }
443 return n&7;
444}