1 /* mpfr_cache -- cache interface for multiple-precision constants in MPFR. 2 3 Copyright 2004-2023 Free Software Foundation, Inc. 4 Contributed by the AriC and Caramba projects, INRIA. 5 6 This file is part of the GNU MPFR Library. 7 8 The GNU MPFR Library is free software; you can redistribute it and/or modify 9 it under the terms of the GNU Lesser General Public License as published by 10 the Free Software Foundation; either version 3 of the License, or (at your 11 option) any later version. 12 13 The GNU MPFR Library is distributed in the hope that it will be useful, but 14 WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY 15 or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public 16 License for more details. 17 18 You should have received a copy of the GNU Lesser General Public License 19 along with the GNU MPFR Library; see the file COPYING.LESSER. If not, see 20 https://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc., 21 51 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */ 22 23 #include "mpfr-impl.h" 24 25 #if 0 /* this function is not used/documented/tested so far, it could be 26 useful if some user wants to add a new constant to mpfr, and 27 implement a cache mechanism for that constant */ 28 void 29 mpfr_init_cache (mpfr_cache_t cache, int (*func)(mpfr_ptr, mpfr_rnd_t)) 30 { 31 MPFR_PREC (cache->x) = 0; /* Invalid prec to detect that the cache is not 32 valid. Maybe add a flag? */ 33 cache->func = func; 34 } 35 #endif 36 37 void 38 mpfr_clear_cache (mpfr_cache_t cache) 39 { 40 if (MPFR_UNLIKELY (MPFR_PREC (cache->x) != 0)) 41 { 42 /* Get the cache in read-write mode */ 43 MPFR_LOCK_WRITE(cache->lock); 44 45 if (MPFR_LIKELY (MPFR_PREC (cache->x) != 0)) 46 { 47 mpfr_clear (cache->x); 48 MPFR_PREC (cache->x) = 0; 49 } 50 51 /* Free the cache in read-write mode */ 52 MPFR_UNLOCK_WRITE(cache->lock); 53 } 54 } 55 56 int 57 mpfr_cache (mpfr_ptr dest, mpfr_cache_t cache, mpfr_rnd_t rnd) 58 { 59 mpfr_prec_t dprec = MPFR_PREC (dest); 60 mpfr_prec_t cprec; /* precision of the cache */ 61 int inexact, sign; 62 MPFR_SAVE_EXPO_DECL (expo); 63 64 /* Call the initialisation function of the cache if it's needed */ 65 MPFR_DEFERRED_INIT_CALL(cache); 66 67 MPFR_SAVE_EXPO_MARK (expo); 68 69 /* Get the cache in read-only mode */ 70 MPFR_LOCK_READ(cache->lock); 71 /* Read the precision within the cache */ 72 cprec = MPFR_PREC (cache->x); 73 if (MPFR_UNLIKELY (dprec > cprec)) 74 { 75 /* Free the cache in read-only mode */ 76 /* And get the cache in read-write mode */ 77 MPFR_LOCK_READ2WRITE(cache->lock); 78 79 /* Retest the precision once we get the lock (since it might have 80 changed). If there is no lock, there is no harm in this code. */ 81 cprec = MPFR_PREC (cache->x); 82 if (MPFR_LIKELY (dprec > cprec)) 83 { 84 /* No previous result in the cache or the precision of the 85 previous result is not sufficient. */ 86 if (MPFR_UNLIKELY (cprec == 0)) /* No previous result. */ 87 { 88 cprec = dprec; 89 mpfr_init2 (cache->x, cprec); 90 } 91 else 92 { 93 /* We increase the cache size by at least 10% to avoid 94 invalidating the cache many times if one performs 95 several computations with small increase of precision. */ 96 cprec += cprec / 10; 97 if (cprec < dprec) 98 cprec = dprec; 99 /* no need to keep the previous value */ 100 mpfr_set_prec (cache->x, cprec); 101 } 102 103 cache->inexact = (*cache->func) (cache->x, MPFR_RNDN); 104 } 105 106 /* Free the cache in read-write mode */ 107 /* Get the cache in read-only mode */ 108 MPFR_LOCK_WRITE2READ(cache->lock); 109 } 110 111 /* now cprec >= dprec is the precision of cache->x */ 112 MPFR_ASSERTD (cprec >= dprec); 113 MPFR_ASSERTD (MPFR_PREC (cache->x) == cprec); 114 115 /* First, check if the cache has the exact value (unlikely). 116 Else the exact value is between (assuming x=cache->x > 0): 117 x and x+ulp(x) if cache->inexact < 0, 118 x-ulp(x) and x if cache->inexact > 0, 119 and abs(x-exact) <= ulp(x)/2. */ 120 121 /* we assume all cached constants are positive */ 122 MPFR_ASSERTN (MPFR_IS_POS (cache->x)); /* TODO... */ 123 sign = MPFR_SIGN (cache->x); 124 MPFR_EXP (dest) = MPFR_GET_EXP (cache->x); 125 MPFR_SET_SIGN (dest, sign); 126 127 /* round cache->x from precision cprec down to precision dprec; 128 since we are in extended exponent range, for the values considered 129 here, an overflow is not possible (and wouldn't make much sense). */ 130 MPFR_RNDRAW_GEN (inexact, dest, 131 MPFR_MANT (cache->x), cprec, rnd, sign, 132 if (MPFR_UNLIKELY (cache->inexact == 0)) 133 { 134 if ((_sp[0] & _ulp) == 0) 135 { 136 inexact = -sign; 137 goto trunc_doit; 138 } 139 else 140 goto addoneulp; 141 } 142 else if (cache->inexact < 0) 143 goto addoneulp; 144 else /* cache->inexact > 0 */ 145 { 146 inexact = -sign; 147 goto trunc_doit; 148 }, 149 MPFR_EXP (dest) ++; 150 MPFR_ASSERTD (MPFR_EXP (dest) <= __gmpfr_emax); 151 ); 152 153 /* Rather a likely, this is a 100% success rate for 154 all constants of MPFR */ 155 if (MPFR_LIKELY (cache->inexact != 0)) 156 { 157 switch (rnd) 158 { 159 case MPFR_RNDZ: 160 case MPFR_RNDD: 161 if (MPFR_UNLIKELY (inexact == 0)) 162 { 163 inexact = cache->inexact; 164 if (inexact > 0) 165 { 166 mpfr_nextbelow (dest); 167 inexact = -inexact; 168 } 169 } 170 break; 171 case MPFR_RNDU: 172 case MPFR_RNDA: 173 if (MPFR_UNLIKELY (inexact == 0)) 174 { 175 inexact = cache->inexact; 176 if (inexact < 0) 177 { 178 mpfr_nextabove (dest); 179 inexact = -inexact; 180 } 181 } 182 break; 183 default: /* MPFR_RNDN */ 184 if (MPFR_UNLIKELY(inexact == 0)) 185 inexact = cache->inexact; 186 break; 187 } 188 } 189 190 MPFR_SAVE_EXPO_FREE (expo); 191 192 /* Free the cache in read-only mode */ 193 MPFR_UNLOCK_READ(cache->lock); 194 195 return mpfr_check_range (dest, inexact, rnd); 196 } 197
Indexes created Fri Jun 05 00:26:10 UTC 2026