#			THE MICHELL SOLUTION

#The terms in each block of Table 8.1 are listed separately in the following lines. 
phi0:=C01*r^2 +C02*r^2*ln(r) +C03*ln(r) +C04*theta;
phi1c:= C11*r^3*cos(theta) +C12*r*theta*sin(theta) +C13*r*ln(r)*cos(theta) +C14*cos(theta)/r;
phi1s:= C15*r^3*sin(theta) +C16*r*theta*cos(theta) +C17*r*ln(r)*sin(theta) +C18*sin(theta)/r;
phinc:=unapply(Cn1*r^(n+2)*cos(n*theta) +Cn2*r^(-n+2)*cos(n*theta) +Cn3*r^n*cos(n*theta) +Cn4*r^(-n)*cos(n*theta),n);
phins:=unapply(Cn5*r^(n+2)*sin(n*theta) +Cn6*r^(-n+2)*sin(n*theta) +Cn7*r^n*sin(n*theta) +Cn8*r^(-n)*sin(n*theta),n);
#
#Notice that the last two terms are functions of n and should be called
#with an appropriate argument. For example
#
#phi:=phins(2);
#
#will generate the terms varying with sin (2*theta)
#
#The required function can now be assembled by combining appropriate terms from the above terms. For example, if we require terms varying with cos(2*theta) and theta-independent terms, we would write
#
phi:=phi0+phinc(2);
#
#The stress components are then given by
#
srr:=(1/r)*diff(phi,r)+(1/r^2)*diff(phi,theta,theta);
srt:=(1/r^2)*diff(phi,theta)-(1/r)*diff(phi,r,theta);
stt:=diff(phi,r,r);
#
#To obtain the corresponding displacements (Table 9.1) , use the file `urt' --- i.e.
#
#read urt;