Code covered by the BSD License

# Pade Approximant of a trunctated series

by

### Nate (view profile)

16 Sep 2011 (Updated )

For a truncated power series, this lightweight code gives a Pade Approximant.

f=PA(B0,B,J,L)
```function f=PA(B0,B,J,L)    %%N.S. Barlow 9/16/11
%PA(B0,B,J,L) is the Pade approximant of the truncated power series B0+B(1)x+B(2)x^2+.....B(n)x^n.
%The [J/L] approximant is given by F(x)=FJ(x)/FL(x), where FJ(x) and FL(x) are polynomials
%of orders J>0 and L>0.  The length of the coefficient vector B must be equal
%to J+L.
%
%Example: Find the [1/2] Pade approximant of Z=2+x+x^2+3x^2
% y=PA(2,[1 1 3],1,2)
%
% y =
%
% -(9*x - 2)/(2*x^2 - 5*x + 1)
%
%Note: I don't make any claims that this code is perfect.  That being said,
%it is easy enough to check the answer.
%
% taylor(y,x,0,'order',5)
%
% ans =
%
% 59*x^5 + 13*x^4 + 3*x^3 + x^2 + x + 2

if (J+L)<length(B) || (J+L)>length(B)
error('For PA(B0,B,J,L), the length of the coefficient vector B must be equal to J+L')
end

syms x u0
p=sym('p',[1 L]);
u=sym('u',[1 J]);

M=length(B);
Z=B0;
for n=1:M
Z=Z+B(n)*x^n;
end
P=1;
for n=1:L
P=P+p(n)*x^n;
end
U=u0;
for n=1:J
U=U+u(n)*x^n;
end
UP=U/P;

ZUP=simple(taylor(UP,x,0,'order',M+1)-Z);
SV=symvar(ZUP); SVL=length(SV);
SV2=SV(1:SVL-1);
C=solve(coeffs(ZUP,x));
C2=struct2cell(C);
C3=[C2{:}];

f=subs(UP,SV2,C3);```