A problem with C function EXP() (2.718 ^ X)

It's a little difficult to know if you are using it incorrectly if you don't provide any code to look at.

Jim Noxon

Hi Jim,

You are correct.

Here is how I am using it:

s16_An_In_1_Alpha = (int)((float)(32767.0 * (float)exp((float)((-6.2831853) * (float)s16_An_Lpf_Hz_1 * 0.00003125))));

s16_An_In_1_Alpha is an int variable.

s16_An_Lpf_Hz1 is an int variable that gets values from 1 to 10000.

When I use this I get incorrect values in s16_An_In_1_Alpha.

Thanks,

Itai

Hi Itai,

I have one suggestion I can think of.  This error may not be in fact due to the EXP function itself but to some error in rounding/casting.  I say this as there is a lot of casting in your single expression. My suggestion is to break down each individual process into single lines and try it like that.  Then slowly make it back up into a single (or a couple of expressions) and see where the problem is:

For example try:

float temp;

temp = (float)s16_An_Lpf_Hz_1;

temp = temp * 0.00003125;

temp = temp * (-6.28);

temp = exp(temp);

temp = temp * 32767.0;

s16_An_In_1_Alpha = (int)temp;

This may be an easier way to find the exact cause of the problem.  Who knows, it may still be in the exp function.

Tim

There are quite a few unnecessary casts in your expression.  However, I don't see any real failures that they would be contributing to.  The casts to float are semantic at best because doubles and floats in the C2000 are of the same size (i.e. 32 IEEE floating point numbers).  All of the integer values will be promoted to floating point values before any of the arithmetic is done.  In the end, the only cast that might be needed would be the cast to int and this one would really only be needed to quiet any truncation warning by the compiler.  Your code could easily be

s16_An_In_1_Alpha = (int)(32767.0 * exp(-6.2831853 * s16_An_Lpf_Hz_1 * 0.00003125));

and you should get the same result.  You might try this just to validate things.

When I run your equation against a simulated C2000 implementation of the exp function and compare it to an exact implementation of the exp function, I get the following error results

Note the error is due to the implementation of the exp function using a first order approximation of a true exp operation.  Is this the amount of error you are seeing?  The above picture doesn't include the cast to an int which, according to the C2000 optimizing compiler guide, merely truncates the fractional digits.  Including this truncation results in an error as shown here.

If you are seeing errors beyond these, the only thing I can think of is the data is incorrect.

There are two sources of errors here, one being the limited precision of the floating point value of around 7.22 decimal digits and the other being the first order approximation of the exp function.

Here is the code which implements the exp function

double exp(double x)
{
double g, z, q, p;
int n;

/*************************************************************************/
/* check if input would produce output out of the range of this function */
/*************************************************************************/
if (x > MAXX) { errno = ERANGE; return (HUGE_VAL); }

if (x < 0) n = (int) (x * INVLOGe2 - 0.5); /* since (int) -1.5 = -1.0 */
else n = (int) (x * INVLOGe2 + 0.5);

/*************************************************************************/
/* g = x - n * ln(2) (but more mathematically stable) */
/*************************************************************************/
g = (x - n * C3) - n * C4;

/*************************************************************************/
/* determine polynomial expression */
/*************************************************************************/
z = g * g;

#if BITS <=29
p = (EXP1 * z + EXP0) * g;
q = EXQ1 * z + EXQ0;
#elif BITS>=30 && BITS<=42
p = (EXP1 * z + EXP0) * g;
q = (EXQ2 * z + EXQ1) * z + EXQ0;
#elif BITS>=43 && BITS<=56
p = ((EXP2 * z + EXP1) * z + EXP0) * g;
q = (EXQ2 * z + EXQ1) * z + EXQ0;
#else
p = ((EXP2 * z + EXP1) * z + EXP0) * g;
q = ((EXQ3 * z + EXQ2) * z + EXQ1) * z + EXQ0;
#endif

/*************************************************************************/
/* exp(x) = exp(g) * 2 ^ (n + 1) */
/*************************************************************************/
return ldexp(0.5 + p / (q - p), n + 1);
}

Here is the code for the ldexp function

double ldexp(double x, int exp)
{
int *ptr = (int *)&x;
long texp = (long) exp + ((ptr[1] >> 7) & 0xFF);

/***********************************************************************/
/* IF RESULT UNDERFLOWS, RETURN 0.0. IF RESULT OVERFLOWS, RETURN */
/* +- INFINITY. */
/***********************************************************************/
if (texp < 1)
x = 0.0;
else if (texp > 0xFE)
{
errno = ERANGE;
x = (x < 0) ? -HUGE_VAL : HUGE_VAL;
}
else
ptr[1] = (ptr[1] & 0x807F) | ((int)texp << 7);

return x;
}

and here are the defined constants used

#define HUGE_VAL DBL_MAX

#define DBL_MAX 3.402823466E+38 /* LARGEST POSITIVE VALUE */

#define BITS 23 /* There are 23 bits in the mantissa */

#define EXP0 0.24999999950
#define EXP1 0.41602886268e-2
#define EXQ0 0.5
#define EXQ1 0.49987178778e-1

#define C3 0.693359375
#define C4 -2.121944400546905827679e-4

#define INVLOGe2 1.4426950408889634074

Regards,

Jim Noxon

Hi Jim,

Thanks for the detailed answer but my problem is much more simple than that.

I will reduce the problem to much more simple code:

float x

x = exp(2);

using the debugger I get x = 16256

I get the same result for ex(3) and exp(4)

Itai.

That makes sense.  Did the compiler output a warning about encountering a function without a prototype?  I suspect it did.

This is one of the areas of programming in C where opinions vary widely.  One of the hall marks of the C language is its flexibility through relaxed type checking.  It provides easy use of the language as it forgives many incorrect uses of objects.  Unfortunately, it can also cause great disdain as you have just realized.

The reason for the failure is because the C language allows for implicit function declarations to be generated by the compiler.  By not includeing the math.h header file, the compiler had no knowledge of how the exp function was supposed to be called nor what type of object it might return.  In your first example (assuming it is the first line of code the compiler saw using the exp function) the comipler saw you calling a function with the identifier exp and you were passing it a float value as defined by the expression evaluated within its parenthesis operators.  Thus the compiler created a prototype for the function as

int exp( float );

and used it for the remainder of the file you were compiling.  As long as every other use of the function provided a valid use of it in this form then the compilation succeeded and you should have received warnings.  In your second examples with exp(2), exp(3), etc... the prototype would have been

int exp( int );

In either case, you were wanting to interpret a float result from the exp function.  Here is where the problem lies.  The actual function was returning a 32 bit floating value but the code that was calling it thought it was merely returning a 16 bit signed integer because of the implicit declaration provided by the compiler.

A 32 bit floating point number as used by the C2000 has the format of

seee_eeee_emmm_mmmm_mmmm_mmmm_mmmm_mmmm

where s is the sign bit, e represents a bit in the exponent, and m represents a mantissa bit.  The value of the exponent is biased so the actual exponent representation is an unsigned value from which the bias is subtracted before calculating the actual value of the number.  In the case of the C2000 this bias value is 128.  Thus the integer returned for the exp(2) expression would have been something like

0100_0001_0110_1100_0111_0011_0010_0101

Assuming you have optimization turned on, this value was most likely returned in a register and was assumed to be a 16 bit signed integer value which was then passed to an int_to_float conversion routine which most likely mis-interpreted it by assuming the upper bits were valid for the 16 bit data type thus giving you the strange answeres you were seeing.

To avoid this type of problem in the future, I would suggest turning on the "Prototypes Required" option in the project options dialog.  If this had been on initially the compilation would have failed and you would have known to include the appropriate header for the function immediately and saved yourself valuable time.  I suppose the warning could have sufficed as well had you worked your code until all warnings were eliminated before beginning debug.

Anyway, I'm glad that you found your problem and are on to bigger things.

Regards,

Jim Noxon