Is there a way to calculate decimal powers using only addition, subtraction, multiplication and division?

I am a programmer who is trying to build an arbitrary precision decimal class library for C# and was successful in implementing all the basic operations like addition, subtraction, multiplication and division but was stopped when trying to calculate powers to decimal (and negative decimal) exponents. Calculating Log and Ln was also not possible since the depend on calculating powers.

So here is what I CAN do (with arbitrary accuracy):

• Calculate the result of ALL operators if both numbers are integers.

• Calculate the result of addition, subtraction, multiplication and division if the numbers are mixed (decimals and integers).

• Calculate the result of the power of a decimal if the exponent is integer.

Here is what I CAN'T do (and I hope to find the answer here):

• Calculate the power of decimals or integers to fractional, negative or decimal exponents.

Some (programmers) might point out that I can use functions already in C# to calculate these but this would not work as I need to build the library to be of arbitrary precision so I have to do the calculation manually.

So, is there a way to do it with the current "tools" I have? Iterative methods are not a problem since I will not be calculating the powers by hand, the CPU will be doing the counting.

• use taylor series – Stefan Perko Dec 22 '15 at 20:01
• How much have you investigated the existing support for arbitrary-precision calculations in C#? For example, stackoverflow.com/questions/2863388/… (note that the answer that says "no" is quite old); there are surely more recently-updated resources as well. – David K Dec 22 '15 at 20:49
• Ironically, that IS the class I am building. The BigDecimal class. My version has much more code (and features) in it that what you saw in the answer on the question you pointed. I actually found that answer and was expandingit until I hit the "power" wall and came here for help. – Amr Al Sayed Dec 22 '15 at 21:58

The simplest way would be to implement log and exp via Taylor series (using the Taylor remainder theorem in order to bound the error to the desired precision) and then just rewrite

$$a^b = \exp(b \log a)$$

However, this is probably not the ideal way of doing things. For instance, you could use the Sasaki-Kanada algorithm to obtain the logarithm.

A lot of this comes down to the kind of tradeoffs you're trying to make. Do you just want something that kind of functions, or do you want to be able to make guarantees about the accuracy of your calculations? Do you need speed, or can it take all the time it wants? Is there a solid reason not to use existing libraries for arbitrary-precision arithmetic, and if so what is it?

There's an entire field of mathematics studying these questions, namely numerical analysis, and many of these problems have been solved (or at least have best-practice solutions). But there's often no simple, one-size-fits-all answer in a lot of cases.

To maybe illustrate a bit of the complexity here: given numbers $a$ and $b$, we know that $a^b = e^{b \log a}$. But if you calculate $\log a$ to a particular precision, and then calculate the product with $b$ to the same precision, and finally exponentiate that to the same precision, then you will not have calculated $a^b$ to that precision -- every time you do an arithmetic calculation you lose a bit of accuracy, and this can compound dramatically.

• Precision is my main priority here since fast calculation of decimals (with a average precession) can be done with already built libraries. C# has an arbitrary precision class but it is only for integers. I am actually using it to make me arbitrary precision decimal class (using an integer mantissa and exponent). I know there are libraries that have such classes but it is too complicated and hard to get into, I wanted to make a cleaner and simpler library. – Amr Al Sayed Dec 22 '15 at 21:47
• If you absolutely want precision then what you want to do is interval arithmetic. You can't specify arbitrary numbers exactly on a computer, but you can store the fact that a number is definitely between a particular upper and lower bound, and calculate using that. – Daniel McLaury Dec 22 '15 at 23:45
• However it sounds like you're working with floating point numbers instead, which means you've already made a trade-off of accuracy for speed. So what it comes down to is determining exactly how you make this tradeoff. I can't really walk you through something like that -- you'd have to learn numerical analysis on your own. – Daniel McLaury Dec 22 '15 at 23:46

Here is what you can do, roughly:

First, we are trying to implement the exponential function $\exp$, which is given by the power series:

$$\exp x =\sum_{k=0}^\infty \frac{x^n}{n!}$$

We are going to use a few of the terms of the series to approximate $\exp x$, say the first four ones:

$$\exp x \approx 1 + x + \frac{x^2}{2} + \frac{x^3}{6}$$

This formula works well for small $|x|$ and increasingly badly for large $|x|$. That is why we reduce all the calculation of $\exp x$ to those with small $|x|$: We have $\exp x = \exp(x\frac{k}{k}) = \exp(\frac{x}{k})^k$. You choose (somehow) a large enough $k$ and calculate $\exp x$ that way, where $\exp(\frac{x}{k})$ is calculated by the aforementoined formula.

Next, observe, that $a^b = \exp( \log a \cdot b)$, if it is defined. Calculating $\log a$ is a little bit more complicated (afaik), you can read up on that somewhere else too. It get's very difficult for $a$ close to $0$ (because $\log a$ goes to $-\infty$). Here is an idea on how to do it: We have the power series:

$$\log(x + 1) = -\sum_{k=1}^\infty \frac{(-1)^k\cdot x^k}{k}$$

which works for $|x|<1$. We use a few terms:

$$\log(x + 1) \approx x - \frac{x^2}{2} + \frac{x^3}{3}$$

(for $|x|<1$)

If $a>1$, then we want to calculate $\log(x+1)$ somehow for an $x$ smaller then $1$, but better close to $0$. Pick the largest $k$ with $\frac{a}{2^k} > 1$. Consider: $$\log a = \log(\frac{a}{2^k}2^k) = \log(x+1)+\log(2^k) = \log(x+1) + k\log 2$$ for $x = \frac{a}{2^k} - 1$. Now you can calculate $\log(x+1)$ with the above formula and remember the value of $\log 2$.

For $a<1$, I really don't know what to do.

Suppose you have a class RationalNum which implements all the other operations.

Suppose that every RationalNum instance has the following attributes:

• bool m_bSign
• NaturalNum m_cEnumerator
• NaturalNum m_cDenominator

Where class NaturalNum implements all the relevant arithmetic operations.

You can use Newton-Raphson method in order to calculate the $n$th root as follows:

RationalNum RationalNum::Root(unsigned int iExponent,unsigned int iAccuracy) const
{
if (iExponent == 0)
throw "Infinite Result";
if (m_bSign && (iExponent%2) == 0)
throw "Undefined Result";

if (m_cEnumerator == 0)
return 0;

RationalNum cRes;

cRes.m_bSign = m_bSign;
cRes.m_cEnumerator  = 1;
cRes.m_cDenominator = 1;

cRes.m_cEnumerator  <<= m_cEnumerator.BitCount()/m_cDenominator.BitCount()/iExponent;
cRes.m_cDenominator <<= m_cDenominator.BitCount()/m_cEnumerator.BitCount()/iExponent;

for (unsigned int i=0; i<iAccuracy; i++)
{
RationalNum cTmp = cRes.Power(iExponent-1);
cRes -= (cRes*cTmp-*this)/(iExponent*cTmp);
}

return cRes;
}