beam2d.m

来自「将梁分成9个单元」· M 代码 · 共 51 行

clear all % static load at the middle  of a straight simple supported beam
         length=1;
         er=2.7e3;
         %ea=1.;
         ee=72e9;
         wid=10e-2;
         high=1e-2;
         ei=wid*high*high*high/12.;
         ea=wid*high;
         %ei=1.;
         ne=6;
         nfe=4;
         nf=12;
         el=length/ne;
         ndg=[13,1,2,3;2,3,4,5;4,5,6,7;6,7,8,9;8,9,10,11;10,11,13,12];%12 dofs 6 elements
         %ndg=[50,50,1,2];%2 dofs 2 elements
         
         em=[156. -22*el 54 13*el;-22*el 4*el*el -13*el -3*el*el;
             54 -13*el 156 22*el;13*el -3*el*el 22*el 4*el*el];
         em=em*er*ea*el/420.;
         ek=[12 -6*el -12 -6*el;-6*el 4*el*el 6*el 2*el*el;
             -12 6*el 12 6*el;-6*el 2*el*el 6*el 4*el*el];
         ek=ek*ee*ei/el/el/el;
         sk=zeros(nf,nf);
         sm=zeros(nf,nf);
         m=zeros(4);
         for kk=1:ne,
             for i=1:nfe,
                 m(i)=ndg(kk,i);
             end
             for i=1:nfe,
                 mi=m(i);
                 if(mi<=nf),
                 for j=1:nfe,
                     mj=m(j);
                      if(mj<=nf),
                        sk(mi,mj)=sk(mi,mj)+ek(i,j);
                        sm(mi,mj)=sm(mi,mj)+em(i,j);
                      end
                   end
                  end
             end
         end

        f=[0 0 0 0 0 -100 0 0 0 0 0 0];
        %f=[-100,0];
        u=inv(sk)*f.';
        u(6)
       -100*(el*ne)^3/48/ee/ei
       %u(6)