from hyperops.powerseries import PowerSeriesI
#Predefined PowerSeries
expps = PowerSeriesI().Exp()
expps.poly(10,x)
expps
PowerSeriesI().Zero()
PowerSeriesI().One()
PowerSeriesI().Id()
#finite powerseries
p = PowerSeriesI([1,2,3,4])
p.poly(10,x)
p 
one = PowerSeriesI([1])
id = PowerSeriesI([0,1])

#power series are just functions that map the index to the coefficient
expps(30)

#power series operations
p+p
p-p
p*p
#/ not yet implemented
p**2
#composition only works for coefficient 0 being 0 in the second operand
dexp = expps - one
expps.compose(dexp)

#we come into interesting regions ...
dexp.compose(dexp)
dexp.iterate(2)
hdexp = dexp.iterate(1/2)
hdexp
#verifying that shoudl be Zero
hdexp.iterate(2) - dexp

#symbolic (parabolic) iteration
dexp.iterate(x)
q = dexp.iterate(1/x).iterate(x)
q(3)
#simiplify and compare
expand(q(3))
dexp(3)

#you can initialize power series with arbitrary functions on natural numbers
#for example the power series of sqrt(2)^x can be given as
bsrt = PowerSeriesI(lambda n: diff(sqrt(2)^x,x,n)(x=0)/factorial(n))

#making the first coefficient 0 to get the decremented exponential
def coeff(n):
if n == 0:
return 0
else:
return bsrt(n)

dbsrt = PowerSeriesI(coeff)

#and now starting hyperbolic iteration
dbsrt2 = dbsrt.iterate(x).iterate(1/x)

#Sage is not able to simplify
simplify(dbsrt2(3))

#but numerically we can verify equality
RR(dbsrt2(3)(x=0.73)-dbsrt(3))