from formulas to pari implementations

[Xorter]

Hey!
Almost everybody knows Schröder's equation:
Psi o y(x) = s × Psi(x)
where s = (d/dx y) o y^o plus/minus infinity
where y^o plus/minus infinity is the fixpoint of the function y.
Almost everybody knows some form of the Abel equation:
Y o (1+x) o Y^o-1 = y(x)
where Y is the superfunction of y.
The connection of the two equation:
Y^o-1 = log Psi(x) / log s
Y = Psi^o-1 o s^x
y(x) = Psi^o-1 o (s × Psi)

y and s are usually given, we look for its iterated function or superfunction Y.
Some of you knows the Carleman matrix representations of analytic functions, I will notate it with [ f(x) ] or just [ f ]. It can be usefull, when you look for a concrate functional power/root of a function, beacuse:
[ y^o n ] = [ y ]^n for every n,
it is derived from:
[f o g] = [f]×[g]
The bigger is the matrix, the more precise is the represantation.

My pari/gp implementation of this conversions:

Code:

D(f,n)=for(k=1,n,f=deriv(f));f;
/*n-th derivative of f for integer n*/
M(f,n)=matrix(n,n,j,k,1/(k-1)!*subst(D(f^(j-1),k-1),x,0));
/*n×n Carleman matrix of f*/
T(A,n)=sum(k=1,n,A[2,k]*x^(k-1));
/*Taylor series of an n×n Carleman matrix*/

How to iterate? I also show you my tetration implementation:

Code:

h=10;
tetinit(a)=A=M(a^x,h);Xt=mateigen(A);Bt=mateigen(A,1)[1];a;
tet(b)=if(real(b)>0,if(real(b)>=1,a^tet(b-1),return(subst(T(Xt*matdiagonal(Bt^b)/Xt,h),x,1)),a^tet(b+1)),if(b==0,1.,log(tet(b+1))/log(a)));
Dtet(n,{lh=10^-6})=sum(k=0,n,(-1)^k*binomial(n,k)*tet(lh*(n-k)*1.))/lh^n;
Ttet(n)=sum(k=0,n,x^k*Dtet(k)/k!);
itet(z)=if(real(z)>1,itet(log(z)/log(a))+1,if(0>real(z),itet(a^z)-1,-tet(-z)));

It is much less precise as Sheldon's and others' fatou.gp programme; it works without using Schröder and Abel equations but using Carleman matrices and diagonalization. It is not enough for me, lets search a better method. What is about determining the Psi function? We said that:
Y = Psi^o-1 o s^x
So to determine inverse of Psi(x) seems promising. Let me notate P = [ Psi(x) ]. The formula says:
P^-1 × (s × P) = [ y(x) ], here are the first product which is matrix-matrix and the second product which is scalar-matrix.
P^-1 × s^x' × P = [ Y(x') ]
It looks like a diagonalized form, so what if I identify P^-1 as mateigen([ y(x) ])? I have already tested it for y = 2x-1 for which it succeded on a small interval but no success with determining the tetration for base 2.

How would you approximate function Psi?

And sorry not for using latex. I forget how to. :(

[Gottfried]

Hey!
Almost everybody knows Schröder's equation:
Psi o y(x) = s × Psi(x)
(...)
The bigger is the matrix, the more precise is the represantation.

My pari/gp implementation of this conversions:
(...)

It is much less precise as Sheldon's and others' fatou.gp programme; it works without using Schröder and Abel equations but using Carleman matrices and diagonalization.
(...)

Hmm, I don't really get it what you are after.      
For one: perhaps my very introductiory essay ContinuousfunctionalIteration  is enough to explain the method which leads to the Carleman-method and Schroeder-function (where I did not know in my first writing that already known/used terms and concepts and thus had rediscovered them on my own), and for second, this is all done for examples with a not too large set of procedures (mostly matrix-functions) in Pari/GP.

[Xorter]

Hmm, I don't really get it what you are after.      
For one: perhaps my very introductiory essay ContinuousfunctionalIteration  is enough to explain the method which leads to the Carleman-method and Schroeder-function (where I did not know in my first writing that already known/used terms and concepts and thus had rediscovered them on my own), and for second, this is all done for examples  with a not too large set of procedures (mostly matrix-functions) in Pari/GP.

I read your paper, it is great, understandable. Thank you for sharing.
In fact, this year in january I have (re)invented an iteration method which is actually the functional form of 2.1.5. b) but I discarded the idea because I misunderstood it. It was the following:
(y^ox')(x) = oexp(x' × olog(y(x)))
Substitute the nullelement to x and x to x', then you give the superfunction of y(x), where
oexp(y) = lim (x + y/h)^oh = sum from k=0 to infinity y^ok / k!
olog(y) = lim h×(y^o(1/h) - x) = sum from k=1 to infinity (-1)^(k+1) y^ok / (k+1)

And I now know that [ oexp(y) ] = exp [y] and [ olog(y) ] = log [y], so I was right. But I thought that:
oexp^ot(x' × olog^ot (y)) is the t-th superfunction of y and t could be real or even complex. It was a mistake, and this is why I had beleived I was wrong until now, despite I could calculate some superfunctions by hand.

Lets try your methods with trigonometric function in no pari/gp but Mathematica because I could not write the code of matrix logarithm in pari/gp.
Here are some codes and plots: definitions, sucos (eigen decomposition), susin (matrix log method) and some wrong results.

It seems almost correct but it is not. Some idea to fix it?

[Xorter]

oexp(x' × olog(y(x))) = y^ox' =: super(y, 1)
oexp(oexp(x'' × x' × olog(olog(y)))) = super(y, 2)
oexp^op ((prod k=1 to p x[k]) × olog^op(y)) = super(y, p)
What if there is a osexp and a oslog operators for which:
osexp(f(x',oslog(y))) = super(y, x')
Let f(a,b) be a+b, so:
osexp(x' + oslog(y)) = super(y, x')
E. g.:
osexp(0 + oslog(y)) = super(y, 0) = y
So osexp is inverse of oslog and vica versa.
osexp(1 + oslog(y)) = oexp(x' × olog(y)) = super(y, 1)
osexp(-1 + oslog(y)) = y o (1+x) o y^o-1 = super(y, -1)

We know that [ oexp(y) ] = exp [y], [ olog(y) ] = log [y], but what is about [ osexp(y) ] and [ oslog(y) ]?
If I knew what these functions are, then I could determine Hyper(a, x, b) = a[x]b = osexp(x + oslog(a+x') o b) and it would be a really new thing in mathematics.
Ofc, I checked and [ osexp(y) ] cannot be sexp [y]...