1/* mpfr_add1sp1 -- internal function to perform a "real" addition on one limb 2 This code was extracted by Kremlin from a formal proof in F* 3 done by Jianyang Pan in April-August 2018: do not modify it! 4 5Source: https://github.com/project-everest/hacl-star/tree/dev_mpfr/code/mpfr 6 7Copyright 2004-2023 Free Software Foundation, Inc. 8Contributed by the AriC and Caramba projects, INRIA. 9 10This file is part of the GNU MPFR Library. 11 12The GNU MPFR Library is free software; you can redistribute it and/or modify 13it under the terms of the GNU Lesser General Public License as published by 14the Free Software Foundation; either version 3 of the License, or (at your 15option) any later version. 16 17The GNU MPFR Library is distributed in the hope that it will be useful, but 18WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY 19or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public 20License for more details. 21 22You should have received a copy of the GNU Lesser General Public License 23along with the GNU MPFR Library; see the file COPYING.LESSER. If not, see 24https://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc., 2551 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */ 26 27#define int64_t long 28#define uint32_t unsigned int 29#define uint64_t mp_limb_t 30 31typedef struct MPFR_Add1sp1_state_s 32{ 33 int64_t sh; 34 int64_t bx; 35 uint64_t rb; 36 uint64_t sb; 37} 38MPFR_Add1sp1_state; 39 40static MPFR_Add1sp1_state 41MPFR_Add1sp1_mk_state(int64_t sh, int64_t bx, uint64_t rb, uint64_t sb) 42{ 43 return ((MPFR_Add1sp1_state){ .sh = sh, .bx = bx, .rb = rb, .sb = sb }); 44} 45 46typedef struct K___uint64_t_int64_t_s 47{ 48 uint64_t fst; 49 int64_t snd; 50} 51K___uint64_t_int64_t; 52 53typedef struct K___uint64_t_uint64_t_int64_t_s 54{ 55 uint64_t fst; 56 uint64_t snd; 57 int64_t thd; 58} 59K___uint64_t_uint64_t_int64_t; 60 61#define MPFR_Lib_mpfr_struct __mpfr_struct 62#define MPFR_Lib_mpfr_RET(I) ((I) != 0 ? ((__gmpfr_flags |= MPFR_FLAGS_INEXACT), (I)) : 0) 63#define MPFR_Exceptions_mpfr_overflow mpfr_overflow 64#define mpfr_prec _mpfr_prec 65#define mpfr_exp _mpfr_exp 66#define mpfr_d _mpfr_d 67#define mpfr_sign _mpfr_sign 68#define MPFR_Lib_gmp_NUMB_BITS GMP_NUMB_BITS 69#define MPFR_Lib_mpfr_EMAX __gmpfr_emax 70 71#define MPFR_RoundingMode_uu___is_MPFR_RNDN(rnd_mode) (rnd_mode == MPFR_RNDN) 72#define MPFR_RoundingMode_mpfr_IS_LIKE_RNDZ MPFR_IS_LIKE_RNDZ 73 74/* same as mpfr_add1sp, but for p < GMP_NUMB_BITS */ 75static int 76mpfr_add1sp1 (mpfr_ptr a, mpfr_srcptr b, mpfr_srcptr c, mpfr_rnd_t rnd_mode, 77 mpfr_prec_t p) 78{ 79 MPFR_Lib_mpfr_struct a0 = a[0U]; 80 MPFR_Lib_mpfr_struct b0 = b[0U]; 81 MPFR_Lib_mpfr_struct c0 = c[0U]; 82 int64_t bx = b0.mpfr_exp; 83 int64_t cx = c0.mpfr_exp; 84 int64_t sh = MPFR_Lib_gmp_NUMB_BITS - p; 85 MPFR_Add1sp1_state st; 86 if (bx == cx) 87 { 88 uint64_t *ap = a0.mpfr_d; 89 uint64_t *bp = b0.mpfr_d; 90 uint64_t *cp = c0.mpfr_d; 91 uint64_t a01 = (bp[0U] >> (uint32_t)1U) + (cp[0U] >> (uint32_t)1U); 92 int64_t bx1 = b0.mpfr_exp + (int64_t)1; 93 uint64_t rb = a01 & (uint64_t)1U << (uint32_t)(sh - (int64_t)1); 94 ap[0U] = a01 ^ rb; 95 uint64_t sb = (uint64_t)0U; 96 st = MPFR_Add1sp1_mk_state(sh, bx1, rb, sb); 97 } 98 else 99 { 100 MPFR_Add1sp1_state ite0; 101 if (bx > cx) 102 { 103 int64_t bx1 = b0.mpfr_exp; 104 int64_t cx1 = c0.mpfr_exp; 105 int64_t d = bx1 - cx1; 106 uint64_t mask = ((uint64_t)1U << (uint32_t)sh) - (uint64_t)1U; 107 MPFR_Add1sp1_state ite1; 108 if (d < sh) 109 { 110 uint64_t *ap = a0.mpfr_d; 111 uint64_t *bp = b0.mpfr_d; 112 uint64_t *cp = c0.mpfr_d; 113 int64_t bx2 = b0.mpfr_exp; 114 uint64_t a01 = bp[0U] + (cp[0U] >> (uint32_t)d); 115 K___uint64_t_int64_t scrut; 116 if (a01 < bp[0U]) 117 scrut = 118 ( 119 (K___uint64_t_int64_t){ 120 .fst = (uint64_t)0x8000000000000000U | a01 >> (uint32_t)1U, 121 .snd = bx2 + (int64_t)1 122 } 123 ); 124 else 125 scrut = ((K___uint64_t_int64_t){ .fst = a01, .snd = bx2 }); 126 uint64_t a02 = scrut.fst; 127 int64_t bx3 = scrut.snd; 128 uint64_t rb = a02 & (uint64_t)1U << (uint32_t)(sh - (int64_t)1); 129 uint64_t sb = (a02 & mask) ^ rb; 130 ap[0U] = a02 & ~mask; 131 ite1 = MPFR_Add1sp1_mk_state(sh, bx3, rb, sb); 132 } 133 else 134 { 135 MPFR_Add1sp1_state ite; 136 if (d < MPFR_Lib_gmp_NUMB_BITS) 137 { 138 uint64_t *ap = a0.mpfr_d; 139 uint64_t *bp = b0.mpfr_d; 140 uint64_t *cp = c0.mpfr_d; 141 int64_t bx2 = b0.mpfr_exp; 142 uint64_t sb = cp[0U] << (uint32_t)(MPFR_Lib_gmp_NUMB_BITS - d); 143 uint64_t a01 = bp[0U] + (cp[0U] >> (uint32_t)d); 144 K___uint64_t_uint64_t_int64_t scrut; 145 if (a01 < bp[0U]) 146 scrut = 147 ( 148 (K___uint64_t_uint64_t_int64_t){ 149 .fst = sb | (a01 & (uint64_t)1U), 150 .snd = (uint64_t)0x8000000000000000U | a01 >> (uint32_t)1U, 151 .thd = bx2 + (int64_t)1 152 } 153 ); 154 else 155 scrut = ((K___uint64_t_uint64_t_int64_t){ .fst = sb, .snd = a01, .thd = bx2 }); 156 uint64_t sb1 = scrut.fst; 157 uint64_t a02 = scrut.snd; 158 int64_t bx3 = scrut.thd; 159 uint64_t rb = a02 & (uint64_t)1U << (uint32_t)(sh - (int64_t)1); 160 uint64_t sb2 = sb1 | ((a02 & mask) ^ rb); 161 ap[0U] = a02 & ~mask; 162 ite = MPFR_Add1sp1_mk_state(sh, bx3, rb, sb2); 163 } 164 else 165 { 166 uint64_t *ap = a0.mpfr_d; 167 uint64_t *bp = b0.mpfr_d; 168 int64_t bx2 = b0.mpfr_exp; 169 ap[0U] = bp[0U]; 170 uint64_t rb = (uint64_t)0U; 171 uint64_t sb = (uint64_t)1U; 172 ite = MPFR_Add1sp1_mk_state(sh, bx2, rb, sb); 173 } 174 ite1 = ite; 175 } 176 ite0 = ite1; 177 } 178 else 179 { 180 int64_t bx1 = c0.mpfr_exp; 181 int64_t cx1 = b0.mpfr_exp; 182 int64_t d = bx1 - cx1; 183 uint64_t mask = ((uint64_t)1U << (uint32_t)sh) - (uint64_t)1U; 184 MPFR_Add1sp1_state ite1; 185 if (d < sh) 186 { 187 uint64_t *ap = a0.mpfr_d; 188 uint64_t *bp = c0.mpfr_d; 189 uint64_t *cp = b0.mpfr_d; 190 int64_t bx2 = c0.mpfr_exp; 191 uint64_t a01 = bp[0U] + (cp[0U] >> (uint32_t)d); 192 K___uint64_t_int64_t scrut; 193 if (a01 < bp[0U]) 194 scrut = 195 ( 196 (K___uint64_t_int64_t){ 197 .fst = (uint64_t)0x8000000000000000U | a01 >> (uint32_t)1U, 198 .snd = bx2 + (int64_t)1 199 } 200 ); 201 else 202 scrut = ((K___uint64_t_int64_t){ .fst = a01, .snd = bx2 }); 203 uint64_t a02 = scrut.fst; 204 int64_t bx3 = scrut.snd; 205 uint64_t rb = a02 & (uint64_t)1U << (uint32_t)(sh - (int64_t)1); 206 uint64_t sb = (a02 & mask) ^ rb; 207 ap[0U] = a02 & ~mask; 208 ite1 = MPFR_Add1sp1_mk_state(sh, bx3, rb, sb); 209 } 210 else 211 { 212 MPFR_Add1sp1_state ite; 213 if (d < MPFR_Lib_gmp_NUMB_BITS) 214 { 215 uint64_t *ap = a0.mpfr_d; 216 uint64_t *bp = c0.mpfr_d; 217 uint64_t *cp = b0.mpfr_d; 218 int64_t bx2 = c0.mpfr_exp; 219 uint64_t sb = cp[0U] << (uint32_t)(MPFR_Lib_gmp_NUMB_BITS - d); 220 uint64_t a01 = bp[0U] + (cp[0U] >> (uint32_t)d); 221 K___uint64_t_uint64_t_int64_t scrut; 222 if (a01 < bp[0U]) 223 scrut = 224 ( 225 (K___uint64_t_uint64_t_int64_t){ 226 .fst = sb | (a01 & (uint64_t)1U), 227 .snd = (uint64_t)0x8000000000000000U | a01 >> (uint32_t)1U, 228 .thd = bx2 + (int64_t)1 229 } 230 ); 231 else 232 scrut = ((K___uint64_t_uint64_t_int64_t){ .fst = sb, .snd = a01, .thd = bx2 }); 233 uint64_t sb1 = scrut.fst; 234 uint64_t a02 = scrut.snd; 235 int64_t bx3 = scrut.thd; 236 uint64_t rb = a02 & (uint64_t)1U << (uint32_t)(sh - (int64_t)1); 237 uint64_t sb2 = sb1 | ((a02 & mask) ^ rb); 238 ap[0U] = a02 & ~mask; 239 ite = MPFR_Add1sp1_mk_state(sh, bx3, rb, sb2); 240 } 241 else 242 { 243 uint64_t *ap = a0.mpfr_d; 244 uint64_t *bp = c0.mpfr_d; 245 int64_t bx2 = c0.mpfr_exp; 246 ap[0U] = bp[0U]; 247 uint64_t rb = (uint64_t)0U; 248 uint64_t sb = (uint64_t)1U; 249 ite = MPFR_Add1sp1_mk_state(sh, bx2, rb, sb); 250 } 251 ite1 = ite; 252 } 253 ite0 = ite1; 254 } 255 st = ite0; 256 } 257 if (st.bx > MPFR_Lib_mpfr_EMAX) 258 { 259 int32_t t = MPFR_Exceptions_mpfr_overflow(a, rnd_mode, a->mpfr_sign); 260 return t; 261 } 262 else 263 { 264 uint64_t *ap = a->mpfr_d; 265 uint64_t a01 = ap[0U]; 266 MPFR_Lib_mpfr_struct uu___72_3478 = a[0U]; 267 a[0U] = 268 ( 269 (MPFR_Lib_mpfr_struct){ 270 .mpfr_prec = uu___72_3478.mpfr_prec, 271 .mpfr_sign = uu___72_3478.mpfr_sign, 272 .mpfr_exp = st.bx, 273 .mpfr_d = uu___72_3478.mpfr_d 274 } 275 ); 276 if (st.rb == (uint64_t)0U && st.sb == (uint64_t)0U) 277 return MPFR_Lib_mpfr_RET((int32_t)0); 278 else if (MPFR_RoundingMode_uu___is_MPFR_RNDN(rnd_mode)) 279 if 280 ( 281 st.rb 282 == (uint64_t)0U 283 || (st.sb == (uint64_t)0U && (a01 & (uint64_t)1U << (uint32_t)st.sh) == (uint64_t)0U) 284 ) 285 { 286 int32_t ite; 287 if (a->mpfr_sign == (int32_t)1) 288 ite = (int32_t)-1; 289 else 290 ite = (int32_t)1; 291 return MPFR_Lib_mpfr_RET(ite); 292 } 293 else 294 { 295 uint64_t *ap1 = a->mpfr_d; 296 ap1[0U] = ap1[0U] + ((uint64_t)1U << (uint32_t)st.sh); 297 if (ap1[0U] == (uint64_t)0U) 298 { 299 ap1[0U] = (uint64_t)0x8000000000000000U; 300 if (st.bx + (int64_t)1 <= MPFR_Lib_mpfr_EMAX) 301 { 302 MPFR_Lib_mpfr_struct uu___72_3574 = a[0U]; 303 a[0U] = 304 ( 305 (MPFR_Lib_mpfr_struct){ 306 .mpfr_prec = uu___72_3574.mpfr_prec, 307 .mpfr_sign = uu___72_3574.mpfr_sign, 308 .mpfr_exp = st.bx + (int64_t)1, 309 .mpfr_d = uu___72_3574.mpfr_d 310 } 311 ); 312 return MPFR_Lib_mpfr_RET(a->mpfr_sign); 313 } 314 else 315 { 316 int32_t t = MPFR_Exceptions_mpfr_overflow(a, rnd_mode, a->mpfr_sign); 317 return MPFR_Lib_mpfr_RET(t); 318 } 319 } 320 else 321 return MPFR_Lib_mpfr_RET(a->mpfr_sign); 322 } 323 else if (MPFR_RoundingMode_mpfr_IS_LIKE_RNDZ(rnd_mode, a->mpfr_sign < (int32_t)0)) 324 { 325 int32_t ite; 326 if (a->mpfr_sign == (int32_t)1) 327 ite = (int32_t)-1; 328 else 329 ite = (int32_t)1; 330 return MPFR_Lib_mpfr_RET(ite); 331 } 332 else 333 { 334 uint64_t *ap1 = a->mpfr_d; 335 ap1[0U] = ap1[0U] + ((uint64_t)1U << (uint32_t)st.sh); 336 if (ap1[0U] == (uint64_t)0U) 337 { 338 ap1[0U] = (uint64_t)0x8000000000000000U; 339 if (st.bx + (int64_t)1 <= MPFR_Lib_mpfr_EMAX) 340 { 341 MPFR_Lib_mpfr_struct uu___72_3781 = a[0U]; 342 a[0U] = 343 ( 344 (MPFR_Lib_mpfr_struct){ 345 .mpfr_prec = uu___72_3781.mpfr_prec, 346 .mpfr_sign = uu___72_3781.mpfr_sign, 347 .mpfr_exp = st.bx + (int64_t)1, 348 .mpfr_d = uu___72_3781.mpfr_d 349 } 350 ); 351 return MPFR_Lib_mpfr_RET(a->mpfr_sign); 352 } 353 else 354 { 355 int32_t t = MPFR_Exceptions_mpfr_overflow(a, rnd_mode, a->mpfr_sign); 356 return MPFR_Lib_mpfr_RET(t); 357 } 358 } 359 else 360 return MPFR_Lib_mpfr_RET(a->mpfr_sign); 361 } 362 } 363} 364