/src/gmp/mpz/cfdiv_q_2exp.c
Line | Count | Source (jump to first uncovered line) |
1 | | /* mpz_cdiv_q_2exp, mpz_fdiv_q_2exp -- quotient from mpz divided by 2^n. |
2 | | |
3 | | Copyright 1991, 1993, 1994, 1996, 1998, 1999, 2001, 2002, 2004, 2012, |
4 | | 2015 Free Software Foundation, Inc. |
5 | | |
6 | | This file is part of the GNU MP Library. |
7 | | |
8 | | The GNU MP Library is free software; you can redistribute it and/or modify |
9 | | it under the terms of either: |
10 | | |
11 | | * the GNU Lesser General Public License as published by the Free |
12 | | Software Foundation; either version 3 of the License, or (at your |
13 | | option) any later version. |
14 | | |
15 | | or |
16 | | |
17 | | * the GNU General Public License as published by the Free Software |
18 | | Foundation; either version 2 of the License, or (at your option) any |
19 | | later version. |
20 | | |
21 | | or both in parallel, as here. |
22 | | |
23 | | The GNU MP Library is distributed in the hope that it will be useful, but |
24 | | WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY |
25 | | or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License |
26 | | for more details. |
27 | | |
28 | | You should have received copies of the GNU General Public License and the |
29 | | GNU Lesser General Public License along with the GNU MP Library. If not, |
30 | | see https://www.gnu.org/licenses/. */ |
31 | | |
32 | | #include "gmp-impl.h" |
33 | | |
34 | | |
35 | | /* dir==1 for ceil, dir==-1 for floor */ |
36 | | |
37 | | static void __gmpz_cfdiv_q_2exp (REGPARM_3_1 (mpz_ptr, mpz_srcptr, mp_bitcnt_t, int)) REGPARM_ATTR (1); |
38 | 0 | #define cfdiv_q_2exp(w,u,cnt,dir) __gmpz_cfdiv_q_2exp (REGPARM_3_1 (w,u,cnt,dir)) |
39 | | |
40 | | REGPARM_ATTR (1) static void |
41 | | cfdiv_q_2exp (mpz_ptr w, mpz_srcptr u, mp_bitcnt_t cnt, int dir) |
42 | 0 | { |
43 | 0 | mp_size_t wsize, usize, abs_usize, limb_cnt, i; |
44 | 0 | mp_srcptr up; |
45 | 0 | mp_ptr wp; |
46 | 0 | mp_limb_t round, rmask; |
47 | |
|
48 | 0 | usize = SIZ (u); |
49 | 0 | abs_usize = ABS (usize); |
50 | 0 | limb_cnt = cnt / GMP_NUMB_BITS; |
51 | 0 | wsize = abs_usize - limb_cnt; |
52 | 0 | if (wsize <= 0) |
53 | 0 | { |
54 | | /* u < 2**cnt, so result 1, 0 or -1 according to rounding */ |
55 | 0 | MPZ_NEWALLOC (w, 1)[0] = 1; |
56 | 0 | SIZ(w) = (usize == 0 || (usize ^ dir) < 0 ? 0 : dir); |
57 | 0 | return; |
58 | 0 | } |
59 | | |
60 | | /* +1 limb to allow for mpn_add_1 below */ |
61 | 0 | wp = MPZ_REALLOC (w, wsize+1); |
62 | | |
63 | | /* Check for rounding if direction matches u sign. |
64 | | Set round if we're skipping non-zero limbs. */ |
65 | 0 | up = PTR(u); |
66 | 0 | round = 0; |
67 | 0 | rmask = ((usize ^ dir) >= 0 ? MP_LIMB_T_MAX : 0); |
68 | 0 | if (rmask != 0) |
69 | 0 | for (i = 0; i < limb_cnt && round == 0; i++) |
70 | 0 | round = up[i]; |
71 | |
|
72 | 0 | cnt %= GMP_NUMB_BITS; |
73 | 0 | if (cnt != 0) |
74 | 0 | { |
75 | 0 | round |= rmask & mpn_rshift (wp, up + limb_cnt, wsize, cnt); |
76 | 0 | wsize -= (wp[wsize - 1] == 0); |
77 | 0 | } |
78 | 0 | else |
79 | 0 | MPN_COPY_INCR (wp, up + limb_cnt, wsize); |
80 | |
|
81 | 0 | if (round != 0) |
82 | 0 | { |
83 | 0 | if (wsize != 0) |
84 | 0 | { |
85 | 0 | mp_limb_t cy; |
86 | 0 | cy = mpn_add_1 (wp, wp, wsize, CNST_LIMB(1)); |
87 | 0 | wp[wsize] = cy; |
88 | 0 | wsize += cy; |
89 | 0 | } |
90 | 0 | else |
91 | 0 | { |
92 | | /* We shifted something to zero. */ |
93 | 0 | wp[0] = 1; |
94 | 0 | wsize = 1; |
95 | 0 | } |
96 | 0 | } |
97 | 0 | SIZ(w) = (usize >= 0 ? wsize : -wsize); |
98 | 0 | } |
99 | | |
100 | | |
101 | | void |
102 | | mpz_cdiv_q_2exp (mpz_ptr w, mpz_srcptr u, mp_bitcnt_t cnt) |
103 | 0 | { |
104 | 0 | cfdiv_q_2exp (w, u, cnt, 1); |
105 | 0 | } |
106 | | |
107 | | void |
108 | | mpz_fdiv_q_2exp (mpz_ptr w, mpz_srcptr u, mp_bitcnt_t cnt) |
109 | 0 | { |
110 | 0 | cfdiv_q_2exp (w, u, cnt, -1); |
111 | 0 | } |