CODEWITHSUNDEEP
cdata-manipulationradixpolynomial-math
Ctrl_Shift_F
Ctrl_Shift_F

Reputation: 1

How to implement 3rd order Polynomial Formula calculations in C on a 16bit MCU

it is my first time posting but I'll start by apologizing in advance if this question has been asked before.

I have been struggling on how to implement a 3rd order polynomial formula in C because of either extremely small values or larger than 32bit results (on a 16bit MCU). I use diffrent values but as an example I would like to compute for "Y" in formula:

Y = ax^3 + bx^2 + cx + d = 0.00000012*(1024^3) + 0.000034*(1024^2) + 0.056*(1024) + 789.10

I need to use a base32 to get a meaningful value for "a" = 515 If I multiply 1024^3 (10bit ADC) then I get a very large amount of 1,073,741,824

I tried splitting them up into "terms A, B, C, and D" but I am not sure how to merge them together because of different resolution of each term and limitation of my 16bit MCU:

u16_TermA = fnBase32(0.00000012) * AdcMax * AdcMax * AdcMax;

u16_TermB = fnBase24(0.000034) * AdcMax * AdcMax;

u16_TermC = fnBase16(0.056) * AdcMax;

u16_TermD = fnBase04(789.10);

u16_Y = u16_TermA + u16_TermB + u16_TermC + u16_TermD;

/* AdcMax is a variable 0-1024; u16_Y needs to be 16bit */

I'd appreciate any help on the matter and on how best to implement this style of computations in C.

Cheers and thanks in advance!

Upvotes: 0

Views: 627

Answers (1)

chux
chux

Reputation: 153498

One step toward improvement:

ax^3 + bx^2 + cx + d --> ((a*x + b)*x + c)*x + d

It is numerically more stable and tends to provide more accurate answers near the zeros of the function and less likely to overflow intermediate calculations.

2nd idea; consider scaling the co-efficents if they maintain their approximate relative values as given on the question.

N = 1024; // Some power of 2
aa = a*N*N*N
bb = b*N*N
cc = c*N

y = ((aa*x/N + bb)*x/N + cc)*x/N + d

where /N is done quickly with a shift.

With a judicious selection of N (maybe 2**14 for high precision avoid 32-bit overflow), then entire code might be satisfactorily done using only integer math.

As aa*x/N is just a*x*N*N, I think a scale of 2**16 works well.

Third idea:

In addition to scaling, often such cubic equations can be re-written as

// alpha is a power of 2
y = (x-root1)*(x-root2)*(x-root3)*scale/alpha

Rather than a,b,c, use the roots of the equation. This is very satisfactory if the genesis of the equation was some sort of curve fitting.

Unfortunately, OP's equation roots has a complex root pair.

x1 = -1885.50539
x2 = 801.08603 + i * 1686.95936
x3 = 801.08603 - i * 1686.95936

... in which case code could use

B = -(x1 + x2);
C = x1 * x2;
y = (x-x1)*(x*x + B*x + C)*scale/alpha

Upvotes: 4

Related Questions