(*Mathematica code to draw the Schwarz-Christoffel map for a rectrangle using power seties.You can set the degree of the truncation of the series and the eccentricty of the rectangle (by setting the harmonic measure of the shorter sides*)

ctor[z_]:=(output={};
Do[AppendTo[output,{Re[z[[k]]],Im[z[[k]]]}],{k,1,Length[z]}];
output);

plotpolylistCA[vert_]:=(vert3={};
Do[vert1=ctor[vert[[k]]];
AppendTo[vert3,vert1];,{k,1,Length[vert]}];
plotpolylist[vert3])

plotpolylist[vertlist_]:=Show[Table[Graphics[{(*GrayLevel[1],Thickness[.01],{Polygon[vertlist[[k]]]},*)GrayLevel[If[k==Length[vertlist],0,.1]],Thickness[If[k==Length[vertlist],.005,.005]],{Line[vertlist[[k]]]}},{PlotRange->All,Axes->None,AspectRatio->Automatic}],{k,1,Length[vertlist]}]];

drawSCnew[para_,alpha_,num_]:=(scseries=Normal[Integrate[Product[Series[(1-zeta Exp[-I para[[k]]])^{alpha[[k]]-1},{zeta,0,num}],{k,1,Length[para]}],zeta]];
coeff=Flatten[Table[Coefficient[scseries,zeta^k],{k,1,num}]];
f[x_]:=Sum[coeff[[k]] x^k,{k,1,num}];
out={f[Exp[I para]]};
out={};
numpts=200;
drawcircle=f[Table[(.5)Exp[I 2 Pi k/numpts],{k,0,numpts}]];
AppendTo[out,drawcircle];
drawcircle=f[Table[(.75)Exp[I 2 Pi k/numpts],{k,0,numpts}]];
AppendTo[out,drawcircle];
drawcircle=f[Table[(.9)Exp[I 2 Pi k/numpts],{k,0,numpts}]];
AppendTo[out,drawcircle];
numpts=10000;
drawcircle=f[Table[Exp[I 2 Pi k/numpts],{k,0,numpts}]];
AppendTo[out,drawcircle];
pts=Join[Table[j/10,{j,0,9}],Table[1-(.1)(.8)^j,{j,0,20}]];
AppendTo[pts,1];
Do[AppendTo[out,f[pts Exp[I para[[k]]]]];,{k,1,Length[para]}];
plotpolylistCA[out])

drawsquare[n_]:=(para=Table[Pi/4+ 2 Pi k/4,{k,1,4}];
alpha={.5,.5,.5,.5};
drawSCnew[para,alpha,n])

drawrect[n_,prob_]:=(epsilon=Pi prob/2;
para={epsilon,Pi-epsilon,Pi+epsilon,2 Pi-epsilon};
alpha={.5,.5,.5,.5};
drawSCnew[para,alpha,n])

ctor[{0,1,I,1+I}]

plotpolylist[{ctor[{0,1,I,1+I}]}]

plotpolylistCA[{{0,1,I,2}}]

drawsquare[10]

drawrect[10,.1]

Manipulate[drawrect[n,p],{n,10,500},{p,.001,.5}]