#SECTION 7.4.2
#
#We write here `C' for dot(Omega).
#
sxx1:=rho*C*x*y/2;
syy1:=-sxx1;
sxy1:=-rho*C*(x^2-y^2)/4;
#
#
#
phi:=A1*x^5*y+A2*x^3*y^3+A3*x*y^5+A4*x^3*y+A5*x*y^3+A6*x*y;
#
#
biharm:=diff(phi,x,x,x,x)+2*diff(phi,x,x,y,y)+diff(phi,y,y,y,y);
#
sxx2:=sxx1+diff(phi,y,y);
syy2:=syy1+diff(phi,x,x);
sxy2:=sxy1-diff(phi,x,y);
#
sxx3:=unapply(sxx2,x,y);
syy3:=unapply(syy2,x,y);
sxy3:=unapply(sxy2,x,y);
#
t1:=sxy3(x,b);
t2:=syy3(x,b);
t3:=sxx3(a,y);
#
s1:=coeff(t1,x,4);
s2:=coeff(t1,x,2);
s3:=coeff(t1,x,0);
s4:=coeff(t2,x,3);
s5:=coeff(t2,x,1);
#
M:=int(t3*y,y=-b..b);
#
b1:=coeff(coeff(biharm,x,1),y,1);
#
solution:=solve({s1=0,s2=0,s3=0,s4=0,s5=0,M=0,b1=0},{A1,A2,A3,A4,A5,A6});
#
phi1:=subs(solution,phi);
#
sxx:=simplify(sxx1+diff(phi1,y,y));
syy:=simplify(syy1+diff(phi1,x,x));
sxy:=simplify(sxy1-diff(phi1,x,y));
#
sxy4:=unapply(sxy,x);
t4:=sxy4(a);
F=-simplify(int(t4,y=-b..b));

