------------------------------------------------------------------ */ decExpOp -- effect exponentiation */ / This computes C = exp(A) */ / res is C, the result. C may be A */ rhs is A */ set is the context; note that rounding mode has no e
| 5223 | /* obvious size to allow for calls from decNumberPower. */ |
| 5224 | /* ------------------------------------------------------------------ */ |
| 5225 | decNumber * decExpOp(decNumber *res, const decNumber *rhs, |
| 5226 | decContext *set, uInt *status) { |
| 5227 | uInt ignore=0; // working status |
| 5228 | Int h; // adjusted exponent for 0.xxxx |
| 5229 | Int p; // working precision |
| 5230 | Int residue; // rounding residue |
| 5231 | uInt needbytes; // for space calculations |
| 5232 | const decNumber *x=rhs; // (may point to safe copy later) |
| 5233 | decContext aset, tset, dset; // working contexts |
| 5234 | Int comp; // work |
| 5235 | |
| 5236 | // the argument is often copied to normalize it, so (unusually) it |
| 5237 | // is treated like other buffers, using DECBUFFER, +1 in case |
| 5238 | // DECBUFFER is 0 |
| 5239 | decNumber bufr[D2N(DECBUFFER*2+1)]; |
| 5240 | decNumber *allocrhs=NULL; // non-NULL if rhs buffer allocated |
| 5241 | |
| 5242 | // the working precision will be no more than set->digits+8+1 |
| 5243 | // so for on-stack buffers DECBUFFER+9 is used, +1 in case DECBUFFER |
| 5244 | // is 0 (and twice that for the accumulator) |
| 5245 | |
| 5246 | // buffer for t, term (working precision plus) |
| 5247 | decNumber buft[D2N(DECBUFFER*2+9+1)]; |
| 5248 | decNumber *allocbuft=NULL; // -> allocated buft, iff allocated |
| 5249 | decNumber *t=buft; // term |
| 5250 | // buffer for a, accumulator (working precision * 2), at least 9 |
| 5251 | decNumber bufa[D2N(DECBUFFER*4+18+1)]; |
| 5252 | decNumber *allocbufa=NULL; // -> allocated bufa, iff allocated |
| 5253 | decNumber *a=bufa; // accumulator |
| 5254 | // decNumber for the divisor term; this needs at most 9 digits |
| 5255 | // and so can be fixed size [16 so can use standard context] |
| 5256 | decNumber bufd[D2N(16)]; |
| 5257 | decNumber *d=bufd; // divisor |
| 5258 | decNumber numone; // constant 1 |
| 5259 | |
| 5260 | #if DECCHECK |
| 5261 | Int iterations=0; // for later sanity check |
| 5262 | if (decCheckOperands(res, DECUNUSED, rhs, set)) return res; |
| 5263 | #endif |
| 5264 | |
| 5265 | do { // protect allocated storage |
| 5266 | if (SPECIALARG) { // handle infinities and NaNs |
| 5267 | if (decNumberIsInfinite(rhs)) { // an infinity |
| 5268 | if (decNumberIsNegative(rhs)) // -Infinity -> +0 |
| 5269 | decNumberZero(res); |
| 5270 | else decNumberCopy(res, rhs); // +Infinity -> self |
| 5271 | } |
| 5272 | else decNaNs(res, rhs, NULL, set, status); // a NaN |
| 5273 | break;} |
| 5274 | |
| 5275 | if (ISZERO(rhs)) { // zeros -> exact 1 |
| 5276 | decNumberZero(res); // make clean 1 |
| 5277 | *res->lsu=1; // .. |
| 5278 | break;} // [no status to set] |
| 5279 | |
| 5280 | // e**x when 0 < x < 0.66 is < 1+3x/2, hence can fast-path |
| 5281 | // positive and negative tiny cases which will result in inexact |
| 5282 | // 1. This also allows the later add-accumulate to always be |
no test coverage detected