/* expx2f.c * * Exponential of squared argument * * * * SYNOPSIS: * * double x, y, expx2f(); * * y = expx2f( x ); * * * * DESCRIPTION: * * Computes y = exp(x*x) while suppressing error amplification * that would ordinarily arise from the inexactness of the argument x*x. * * * ACCURACY: * * Relative error: * arithmetic domain # trials peak rms * IEEE -9.4, 9.4 10^7 1.7e-7 4.7e-8 * */ /* Cephes Math Library Release 2.9: June, 2000 Copyright 2000 by Stephen L. Moshier */ extern float floorf (float); extern float expf (float); extern float MAXLOGF; extern float MAXNUMF; float expx2f (float x) { float u, u1, m; if (x < 0) x = -x; /* Represent x as an exact multiple of 1/32 plus a residual. */ m = .03125f * floorf(32.0f * x + 0.5f); x -= m; /* x**2 = m**2 + 2mf + f**2 */ u = m * m; u1 = 2 * m * x + x * x; if ((u+u1) > MAXLOGF) return (MAXNUMF); /* u is exact, u1 is small. */ u = expf(u) * expf(u1); return(u); }