Algorithm for tetration to work with floating point numbers

So far, I've figured out an algorithm for tetration that works.

However, although the variable a can be floating or integer, unfortunately, the variable b must be an integer number.

How can I modify the algorithm so that both a and b can be floating point numbers and the correct answer will be produced?

// Hyperoperation type 4:
public double tetrate(double a, double b)
{
double total = a;
for (int i = 1; i < b; i++) total = pow(a, total);
}


In an attempt to solve this, I've create my own custom power() function (trying to avoid roots, and log functions), and then successfully generalized it to multiplication. Unfortunately, when I then try to generalize to tetration, numbers go pear shaped.

I would like an algorithm to be precise up to x amount of decimal places, and not an approximation as Wikipedia talks about.

• OMG! Is that really C code? I thought I was the last programmer who remembered it! I must have misread something . . . – Robert Lewis Apr 25 '15 at 2:53
• There is no standard way to extend tetration to non-natural exponents. See wiki. – vadim123 Apr 25 '15 at 3:13
• @vadim123: Wow, it looks to be an open problem with a little debate thrown in about what would qualify. Never expected that for something so seemingly simple. I'm looking or the most 'mathematically useful' if at all possible, without approximations. – Dan W Apr 25 '15 at 3:38
• Have you looked at tetration.org/Tetration/index.html or arxiv.org/abs/1410.3896 ? – Daniel Geisler May 7 '15 at 2:04
• Before you write the code, you need to define what it should calculate. What answer do you want for $^{2.5}7, \ ^{3.456}5.67,\ ^{\pi}5$? Once you have a clear definition, you can reduce it to code. – Ross Millikan Jan 16 '16 at 15:00

When I studied the various known matrices of combinatorical numbers I also looked at the following idea: what if we discuss functions with a set of results, not only one number? So for instance the concept of sine and cosine gets some special charme if we look not only at f(x) = sin(x), g(x)=cos(x) but at 2x2-matrices containing cos() and sin() and the input of matrix-multiplications has two parameters and the output as well. This helps much to understand algebra with complex numbers, for instance.
After that I thought: what if we look at matrices which allow not only to look at a function $f(x)=a_0 + a_1x + a_2x^2 + a_3x^3 + ...$ involving all the powers of $x$ as input, but also spit out all powers of $f(x)$ in the same manner. Then for instance

$$V(x) \cdot P = V(x+1)$$ where $V(x)$ means the vector $V(x)=[1,x,x^2,x^3,x^4,...]$ and $P$ means the upper triangular Pascal matrix. Of course $$(V(x) \cdot P )\cdot P = V(x+2)$$ and of course it is easily to see that $$V(x) \cdot P^h = V(x+h)$$ first for integral powers of $P$ which any programmer can implement - and if you download Pari/GP you can just use the matrix-language of it.
It needs a certain ingenious mind to find out how to define fractional powers of $P$ to make the fractional operation of addition possible...

After that one can find, that a very well known combinatorical matrix $B$ performs $$V(x) \cdot B = V(exp(x))$$
and logically $$V(x) \cdot B^h = V(exp^{oh}(x))$$
for integer $h$ - but which means implementing tetration to integer iterates $h$.
To find now the fractional-h power of $B$ -which is your question- is not so simple, and using matrices practically means to truncate them to finite size.

Well, anyway meaningful approximations can be done and truncated versions of $B$ can be diagonalized and general powers of it can then be found by simply computing the general powers on each of the scalar numbers in the diagonal - and so one can introduce fractional powers of $B$ in the formula $$V_{32}(x) \cdot B_{32x32}^h \approx V_{32}(exp^{oh}(x))$$
(for for instance $32$ and $32x32$ sized vectors and matrix) and this implements then the fractional iteration of the exponential, aka tetration.

The matrices $P$ and $B$ of which I talk here are known as "Carleman-matrices" and this whole concept can be applied to fractional iteration of many functions for which you have a power series. (Well, we'll have convergence issues and much more complication and non-uniqueness and what not - but that's not the focus of this answer)
So using "mateigen" in Pari/GP for the diagonalization allows the following Pari/GP code:

n=32
V(x,dim=n) = vector(dim,r,x^{r-1})      \\ define the function for generating V(x)
b = sqrt(2)
bl = log(b)
B = matrix(n,n,r,c,bl^(r-1)*(c-1)^(r-1)/(r-1)!)  \\ define the Carlemanmatrix
\\ for   V(x) * B = V(b^x)
Y = V(1) * B  \\ in Y we'll have approximately V(b^x)
\\where the first two entries have the best approximation to the expected
\\ value
\\ now diagonalize
M = mateigen(B) \\ needs high precision of say 800 or 1600 digits for all operations
W = M^-1
D = diag(W*B*M)
\\ after this we can do B^h = M * D^h * W
\\  and  D^h can be done by d_k,k^h on the elements of the diagonal
dpow(A,h=1)= for (r=1,#A, A[r,r]=A[r,r]^h);return(A) \\ define a function for powers of D

BPow (h)= M*dpow(D,h)*W  \\ define a function for the fractional power of B

\\ after this we can do
Y = V(1) * BPow(0.5)
\\ and have in Y a nice approximation for the half-iterate or b^^0.5
\\where b=sqrt(2)


I you want to experiment with this I recommend to use dimension $n=16$ first (with much poorer approximation) but already you'll need default(realprecision,200) or so. To do mateigen on matrix $B$ of size $32x32$ one must compute $b$, $bl$, $B$ and so on already to even $1600$ or $2000$ internal digits precision and the mateigen can need several minutes to complete. I tried this also with size $n=64$ and this was really a full afternoon computation... but is a nice exercise anyway.

Remark: I've discussed this a bit in the tetrationforum and have made an article which compares a couple of such naive approaches to fractional exponential-towers. The above method I've called there "polynomial approach" because using finite matrices without taking care for compatibility with the theoretically infinite sized matrices and their fractional powers means just involve polynomials for that approximative solutions. See the index and the link to the essay

• Looks interesting. Can you link to the precise page on the tetration forum? When I get more time, I might have a go at implementing the code you have, or maybe you or someone else could expand some of those high level functions to lower level ones at some stage. – Dan W Jul 4 '15 at 9:04
• Hmm, I'm not much active on this this days. Look at my oldest postings, I've called that "matrixmethod" (this is perhaps the most fruitful searchable term). And in the beginning I did not know that this all was already known with the term "Carlemanmatrix". There is also some literature on this, look for Aldrovando (the link to the ArXiv by D. Geisler), Eri Jabotinsky, S.C.Woon (perhaps I can come back to this later with more information - we have it extremely hot today which makes it uncomfortable to sort out links...) – Gottfried Helms Jul 4 '15 at 9:15
• Perhaps this (go.helms-net.de/math/tetdocs/ContinuousfunctionalIteration.pdf) is an interesting essay, although it is extremely amateurish (one of my earlier tries to put numbertheoretic puzzles into shape - I had even difficulties to refer to functions having powerseries in the usual way) - I should rewrite this. But it might give a good impression for a first read. – Gottfried Helms Jul 4 '15 at 9:35

I feel like this is more of a programming techniques question than anything else, since the issues is about data type conversion. So you want to convert from float to int. So here is my suggested code (this should work. If not, tell me and I will edit this):

// Hyperoperation type 4: public double tetrate(double a, double b) { double total = a; int coefficient = 0; for (int i = 1; i <= b; i++) { total = pow(a, total); coefficient = i; } total = pow(b-double(coefficient), total); return total; }

As you can see, I added an (int) statement before the b variable. This is called a cast. It treats that double like an integer.

• I think you misunderstood what I wanted. I wanted the exponent to allow fractional numbers. This would completely rewrite the algorithm. – Dan W May 13 '15 at 8:11
• @DanW The exponent does allow double values. Take a look at the documentation. – SalmonKiller May 13 '15 at 13:46
• Your alg would calc 2^^2 and 2^^3, but not 2^^2.5 (which would only default to 2^^2). b gets converted to an integer (even in my original code if it was C#), so it's of no use. The algorithm would have to be a heck of a lot more complicated. – Dan W May 13 '15 at 15:24
• @DanW Ah. Try my edited code now. – SalmonKiller May 13 '15 at 22:06
• No joy. 2^^2 and 2^^3 both come out as zero. And 2^^2.5 comes out as about 1.5E-05. – Dan W May 14 '15 at 8:12

There is no universally accepted way to extend hyperoperations to the Real numbers.

Refer to the very good answer to this questions

Example $x$, $y$ and $z$ values for $x\uparrow^\alpha y=z$ where $\alpha\in \Bbb R-\Bbb N$

But there are many others, like

How to evaluate fractional tetrations?