新建文件夹.zip_MATLAB程序_三角形单元_应力应变_悬臂梁

clc clear all E=71e9; NU=0.3; t=0.01; xi1=1;yi1=0.1;xj1=0.95;yj1=0.08;xm1=1;ym1=0.08; i1=1;j1=8;m1=2; K=zeros(252,252); for p=1:20 k1=LinearTriangleElementStiffness(E,NU,t,xi1-0.05*(p-1),yi1,xj1-0.05*(p-1),yj1,xm1-0.05*(p-1),ym1,1); K=LinearTriangleAssemble(K,k1,i1+6*(p-1),j1+6*(p-1),m1+6*(p-1)); end xi2=1;yi2=0.08;xj2=0.95;yj2=0.06;xm2=1;ym2=0.06; i2=2;j2=9;m2=3; for p=1:20 k2=LinearTriangleElementStiffness(E,NU,t,xi2-0.05*(p-1),yi2,xj2-0.05*(p-1),yj2,xm2-0.05*(p-1),ym2,1); K=LinearTriangleAssemble(K,k2,i2+6*(p-1),j2+6*(p-1),m2+6*(p-1)); end xi3=1;yi3=0.06;xj3=0.95;yj3=0.04;xm3=1;ym3=0.04; i3=3;j3=10;m3=4; for p=1:20 k3=LinearTriangleElementStiffness(E,NU,t,xi3-0.05*(p-1),yi3,xj3-0.05*(p-1),yj3,xm3-0.05*(p-1),ym3,1); K=LinearTriangleAssemble(K,k3,i3+6*(p-1),j3+6*(p-1),m3+6*(p-1)); end xi4=1;yi4=0.04;xj4=0.95;yj4=0.02;xm4=1;ym4=0.02; i4=4;j4=11;m4=5; for p=1:20 k4=LinearTriangleElementStiffness(E,NU,t,xi4-0.05*(p-1),yi4,xj4-0.05*(p-1),yj4,xm4-0.05*(p-1),ym4,1); K=LinearTriangleAssemble(K,k4,i4+6*(p-1),j4+6*(p-1),m4+6*(p-1)); end xi5=1;yi5=0.02;xj5=0.95;yj5=0;xm5=1;ym5=0; i5=5;j5=12;m5=6; for p=1:20 k5=LinearTriangleElementStiffness(E,NU,t,xi5-0.05*(p-1),yi5,xj5-0.05*(p-1),yj5,xm5-0.05*(p-1),ym5,1); K=LinearTriangleAssemble(K,k5,i5+6*(p-1),j5+6*(p-1),m5+6*(p-1)); end xi6=0.95;yi6=0.1;xj6=0.95;yj6=0.08;xm6=1;ym6=0.1; i6=7;j6=8;m6=1; for p=1:20 k6=LinearTriangleElementStiffness(E,NU,t,xi6-0.05*(p-1),yi6,xj6-0.05*(p-1),yj6,xm6-0.05*(p-1),ym6,1); K=LinearTriangleAssemble(K,k6,i6+6*(p-1),j6+6*(p-1),m6+6*(p-1)); end xi7=0.95;yi7=0.08;xj7=0.95;yj7=0.06;xm7=1;ym7=0.08; i7=8;j7=9;m7=2; for p=1:20 k7=LinearTriangleElementStiffness(E,NU,t,xi7-0.05*(p-1),yi7,xj7-0.05*(p-1),yj7,xm7-0.05*(p-1),ym7,1); K=LinearTriangleAssemble(K,k7,i7+6*(p-1),j7+6*(p-1),m7+6*(p-1)); end xi8=0.95;yi8=0.06;xj8=0.95;yj8=0.04;xm8=1;ym8=0.06; i8=9;j8=10;m8=3; for p=1:20 k8=LinearTriangleElementStiffness(E,NU,t,xi8-0.05*(p-1),yi8,xj8-0.05*(p-1),yj8,xm8-0.05*(p-1),ym8,1); K=LinearTriangleAssemble(K,k8,i8+6*(p-1),j8+6*(p-1),m8+6*(p-1)); end xi9=0.95;yi9=0.04;xj9=0.95;yj9=0.02;xm9=1;ym9=0.04; i9=10;j9=11;m9=4; for p=1:20 k9=LinearTriangleElementStiffness(E,NU,t,xi9-0.05*(p-1),yi9,xj9-0.05*(p-1),yj9,xm9-0.05*(p-1),ym9,1); K=LinearTriangleAssemble(K,k9,i9+6*(p-1),j9+6*(p-1),m9+6*(p-1)); end xi10=0.95;yi10=0.02;xj10=0.95;yj10=0;xm10=1;ym10=0.02; i10=11;j10=12;m10=5; for p=1:20 k10=LinearTriangleElementStiffness(E,NU,t,xi10-0.05*(p-1),yi10,xj10-0.05*(p-1),yj10,xm10-0.05*(p-1),ym10,1); K=LinearTriangleAssemble(K,k10,i10+6*(p-1),j10+6*(p-1),m10+6*(p-1)); end k=K(1:240,1:240); f=[0;-100;zeros(238,1)]; d1=k\f; d=[d1;zeros(12,1)]; % F=K*U; y=zeros(200,3); q=0; for i=1:119 switch rem(i,6) case 1 j=2*i; u=[d(j-1) d(j) d(j+13) d(j+14) d(j+1) d(j+2)]; u=u'; xl=1;yl=0.1;xm=0.95;ym=0.08;xn=1;yn=0.08; y(i+q,:)=LinearTriangleElementStresses(E,NU,xl,yl,xm,ym,xn,yn,1,u)'; xl=xl-0.05;xm=xm-0.05;xn=xn-0.05; case 2 j=2*i; u=[d(j-1) d(j) d(j+13) d(j+14) d(j+1) d(j+2)]; u=u'; xl=1;yl=0.08;xm=0.95;ym=0.06;xn=1;yn=0.06; y(i+q,:)=LinearTriangleElementStresses(E,NU,xl,yl,xm,ym,xn,yn,1,u)'; xl=xl-0.05;xm=xm-0.05;xn=xn-0.05; case 3 j=2*i; u=[d(j-1) d(j) d(j+13) d(j+14) d(j+1) d(j+2)]; u=u'; xl=1;yl=0.06;xm=0.95;ym=0.04;xn=1;yn=0.04; y(i+q,:)=LinearTriangleElementStresses(E,NU,xl,yl,xm,ym,xn,yn,1,u); xl=xl-0.05;xm=xm-0.05;xn=xn-0.05; case 4 j=2*i; u=[d(j-1) d(j) d(j+13) d(j+14) d(j+1) d(j+2)]; u=u'; xl=1;yl=0.04;xm=0.95;ym=0.02;xn=1;yn=0.02; y(i+q,:)=LinearTriangleElementStresses(E,NU,xl,yl,xm,ym,xn,yn,1,u)'; xl=xl-0.05;xm=xm-0.05;xn=xn-0.05; case 5 j=2*i; u=[d(j-1) d(j) d(j+13) d(j+14) d(j+1) d(j+2)]; u=u'; xl=1;yl=0.02;xm=0.95;ym=0;xn=1;yn=0; y(i+q,:)=LinearTriangleElementStresses(E,NU,xl,yl,xm,ym,xn,yn,1,u)'; xl=xl-0.05;xm=xm-0.05;xn=xn-0.05; otherwise q=q+4; end end q=5; for i=1:119 switch rem(i,6) case 1 j=2*i; u=[d(j-1) d(j) d(j+11) d(j+12) d(j+13) d(j+14)]; u=u'; xl=1;yl=0.1;xm=0.95;ym=0.1;xn=0.95;yn=0.08; y(i+q,:)=LinearTriangleElementStresses(E,NU,xl,yl,xm,ym,xn,yn,1,u)'; xl=xl-0.05;xm=xm-0.05;xn=xn-0.05; case 2 j=2*i; u=[d(j-1) d(j) d(j+11) d(j+12) d(j+13) d(j+14)]; u=u'; xl=1;yl=0.08;xm=0.95;ym=0.08;xn=0.95;yn=0.06; y(i+q,:)=LinearTriangleElementStresses(E,NU,xl,yl,xm,ym,xn,yn,1,u)'; xl=xl-0.05;xm=xm-0.05;xn=xn-0.05; case 3 j=2*i; u=[d(j-1) d(j) d(j+11) d(j+12) d(j+13) d(j+14)]; u=u'; xl=1;yl=0.06;xm=0.95;ym=0.06;xn=0.95;yn=0.04; y(i+q,:)=LinearTriangleElementStresses(E,NU,xl,yl,xm,ym,xn,yn,1,u)'; xl=xl-0.05;xm=xm-0.05;xn=xn-0.05; case 4 j=2*i; u=[d(j-1) d(j) d(j+11) d(j+12) d(j+13) d(j+14)]; u=u'; xl=1;yl=0.04;xm=0.95;ym=0.04;xn=0.95;yn=0.02; y(i+q,:)=LinearTriangleElementStresses(E,NU,xl,yl,xm,ym,xn,yn,1,u)'; xl=xl-0.05;xm=xm-0.05;xn=xn-0.05; case 5 j=2*i; u=[d(j-1) d(j) d(j+11) d(j+12) d(j+13) d(j+14)]; u=u'; xl=1;yl=0.02;xm=0.95;ym=0.02;xn=0.95;yn=0; y(i+q,:)=LinearTriangleElementStresses(E,NU,xl,yl,xm,ym,xn,yn,1,u)'; xl=xl-0.05;xm=xm-0.05;xn=xn-0.05; otherwise q=q+4; end end y=y*(1e-6); % epequ=y/E;