#######################
# Maple code to compute
# a Pade approximant to
# a chosen function
#
# M. G. Blyth 05/01/10
#######################


macro(sim=simplify(%)): macro(ex=expand(%)):

# Set the level of approximation and the function f(x) to be approximated

N:= 2: M:= 2:

f:= x-> exp(x):

#f:=x-> (7+(1+x)^(4/3))^(1/3); 


#########
# Compute the Pade
# approximant for chosen
# function and truncation 
# level
#########

T:= N+M+1:

# Define the Pade approximant

pd:= x-> sum(A[n]*x^n,n=0..N)/(1+sum(B[n]*x^n,n=1..M)):

# Evaluate first T terms of Taylor series of approximant

pds:= convert(series(pd(x),x,T),polynom):

# Evaluate first T terms of Taylor series of target function

es:= convert(series(f(x),x,T),polynom):

df:= pds-es:

# Match the first T coefficients of each Taylor series

e[0]:= eval(df,x=0):
for i from 1 to T do
	e[i]:= coeff(df,x,i):
end do:

# Solve

uns:= {seq(A[i],i=0..N),seq(B[i],i=1..M)}:
eqns:= {seq(e[i],i=0..T)}:

sols:= solve(eqns,uns):

for i from 0 to N do
	A[i]:= subs(sols,A[i]):
end do:
for i from 1 to M do
	B[i]:= subs(sols,B[i]):
end do:

# Display approximant

pd(x);