> restart:

> with(linalg):

Warning, the protected names norm and trace have been redefined and unprotected

> a:=5:b:=3:r:=1.34:s:=2.9:

> ll:=[a*cos(r),b*sin(r)]:l:=[a*cos(s),b*sin(s)]:

> m:=[l,ll]:

> for k to 150 do u:=ll-l: n:=[-b*ll[1]/a,-a*ll[2]/b]:v:=u-2*dotprod(u , n)*n/dotprod(n,n):
t:=-2*(b^2*v[1]*ll[1]+a^2*ll[2]*v[2])/(b^2*v[1]^2+a^2*v[2]^2):lll:=ll+t*v:m:=[op(m),lll]:l
:=ll:ll:=lll:od:

> plot(m,scaling=constrained,axes=none);

> restart:

> with(linalg):

Warning, the protected names norm and trace have been redefined and unprotected

> a:=5:b:=3:r:=1.34:s:=5.9:

> ll:=[a*cos(r),b*sin(r)]:l:=[a*cos(s),b*sin(s)]:

> m:=[l,ll]:

> for k to 150 do u:=ll-l: n:=[-b*ll[1]/a,-a*ll[2]/b]:v:=u-2*dotprod(u , n)*n/dotprod(n,n):
t:=-2*(b^2*v[1]*ll[1]+a^2*ll[2]*v[2])/(b^2*v[1]^2+a^2*v[2]^2):lll:=ll+t*v:m:=[op(m),lll]:
l:=ll:ll:=lll:od:

> plot(m,scaling=constrained,axes=none);