- 积分
- 1
- 注册时间
- 2010-10-4
- 仿真币
-
- 最后登录
- 1970-1-1
|
function HW4()
sum=0;
for d=0:6.28/360:6.28
v=(-u0(d)*cos(d)-XTraction(d)*du(d)-YTraction(d)*dv(d))*0.01*(6.28/360);
sum=sum+v;
end
function U=u0(d)
%calculate energy according to stresses.
U=(sigemax(d)^2+sigemay(d)^2+5*tau(d)^2)/(4*100*10^9);
function sigema=sigemax(d)
%calculate x stress
sigema=10*cos(d/2)*(1-sin(d/2)*sin(3*d/2))/sqrt(2*3.14);
function sigema=sigemay(d)
%calculate y stress
sigema=10*cos(d/2)*(1+sin(d/2)*sin(3*d/2))/sqrt(2*3.14);
function tau=tau(d)
%calculate tau stress
tau=10*cos(d/2)*sin(d/2)*cos(3*d/2)/sqrt(2*3.14);
function x=du(d)
%calculate derivative of u to x
%du=2.5*cos(d/2)*(1+2*sin(d/2)^2)/(4*100000)/cos(d)+0.25/(4*100000*sqrt(2*3.14))*(-sin(d/2)*(1+2*sin(d/2)^2)+4*cos(d/2)*sin(d/2)*cos(d/2))/(-0.01*sin(d));
[x,y]=solve('2.49e-5*cos(d/2)*(1+2*sin(d/2)^2)=cos(d)*x+y*sin(d)','0.0249e-5*(-sin(d/2)*(1+2*sin(d/2)^2)+4*cos(d/2)*sin(d/2)*cos(d/2))=-x*0.01*sin(d)+y*0.01*cos(d)','x,y');
x=subs(x);
function x=dv(d)
%calculate derivative of v to x
[x,y]=solve('2.49e-5*sin(d/2)*(3-2*cos(d/2)^2)=x*cos(d)+y*sin(d)','0.0249e-5*(cos(d/2)*(3-2*cos(d/2)^2))+4*sin(d/2)*cos(d/2)*sin(d/2)=-x*0.01*sin(d)+y*0.01*cos(d)','x,y');
x=subs(x);
function tx=XTraction(d)
%claculate traction of x-direction
tx=-sigemax(d)*cos(d)-tau(d)*sin(d);
function ty=YTraction(d)
%claculate traction of y-direction
ty=-tau(d)*cos(d)-sigemay*sin(d); |
|