clear

%  求解單元 1 的剛度矩陣
x=sym ('x');        L1=sym ('L1');      % 定義x和L1為符號變量
E=3e11;        I=5.2e-7;
N=[ 1 - 3 * (x^2) / (L1^2) + 2 * (x^3) / (L1^3),  x - 2 * (x^2) / L1 + (x^3) / (L1^2),....
      3 * (x^2) / (L1^2) - 2 * (x^3) / (L1^3),  -(x^2) / L1 + (x^3) / (L1^2) ];           % 形函數(shù)
Ni=diff ( N, x, 2 ); % 對形函數(shù)進(jìn)行2次微分
Nt=transpose ( Ni);
kk1=Nt * Ni;
K1=E* I * int (kk1, 0, 'L1');      % 對kk1進(jìn)行積分，積分區(qū)間0－L1
K1=subs (K1, 'L1', 5);             % 將K1矩陣中的L1代換成5，重新計(jì)算K1

%  求解單元 2 的剛度矩陣
x=sym ('x');        L2=sym ('L2'); % 定義x和L2為符號變量
N=[ 1 - 3 * (x^2) / (L2^2) + 2 * (x^3) / (L2^3),  x - 2 * (x^2) / L2 + (x^3) / (L2^2),....
      3 * (x^2) / (L2^2) - 2 * (x^3) / (L2^3),  -(x^2) / L2 + (x^3) / (L2^2) ];           % 形函數(shù)
Ni=diff (N, x, 2 );    % 對形函數(shù)進(jìn)行2次微分
Nt=transpose ( Ni);
kk2=Nt * Ni;
K2=E * I * int (kk2, 0, 'L2');     % 對kk2進(jìn)行積分，積分區(qū)間0－L2
K2=subs (K2, 'L2', 5);             % 將K1矩陣中的L2代換成5，重新計(jì)算K2

%  單元 1 的受力轉(zhuǎn)為節(jié)點(diǎn)載荷
format  rat          %將計(jì)算結(jié)果用分?jǐn)?shù)表示
x=sym ( 'x' );      L1=sym ( 'L1');
w= -1000;
N=[ 1 - 3 * (x^2) / (L1^2) + 2 * (x^3) / (L1^3),  x - 2 * (x^2) / L1 + (x^3) / (L1^2),....
       3 * (x^2) / (L1^2) - 2 * (x^3) / (L1^3),  -(x^2) / L1 + (x^3) / (L1^2) ];            % 形函數(shù)
P= x * w / L1;      %  外力分布函數(shù)
F= int ( P * transpose ( N ), x, 0, 'L1');
F=subs (F, 'L1', 5);

%  單元 2 的受力轉(zhuǎn)為節(jié)點(diǎn)載荷
x=sym ( 'x' );      L2=sym ( 'L2');
w= -1000;
N=[ 1 - 3 * (x^2) / (L2^2) + 2 * (x^3) / (L2^3),  x - 2 * (x^2) / L2 + (x^3) / (L2^2),....
       3 * (x^2) / (L2^2) - 2 * (x^3) / (L2^3),  -(x^2) / L2 + (x^3) / (L2^2) ];             % 形函數(shù)
P= ( (1-x) * w / L2 );      %分布函數(shù)
F= int (P * transpose ( N ), x, 0, 'L2');
F=subs (F, 'L2', 5);

%  合并剛度矩陣
K1=[K1,zeros(4,2);zeros(2,6)];      K2=[zeros(2,6);zeros(4,2),K2];
KK=K1+K2;

%  求解最后結(jié)果
F=[ -750;  -2500 / 3;  -3500;  0;  -750;  2500 / 3 ];
KK(1,  :)=[1,  0,  0,  0,  0,  0];   % 加入邊界條件，將第一行用此行代替            
KK(5,  :)=[0,  0,  0,  0,  1,  0];   % 加入邊界條件，將第五行用此行代替
F(1,1)=0;                            %  將F矩陣第一行用零代替
F(5,1)=0;                            %  將F矩陣第五行用零代替
format  long
u=inv (KK) * F