%% substructure 1
[Freq1,CCD1]=femodal(mmc1,kkc1);
% 对接界面节点自由度
C1=7;C2=12;
nn1=length(kkc1);
% 波尔矩阵
Bo1=zeros(nn1);
Bo1(:,7:nn1)=[];
Bo1(nn1-5,1)=1;
Bo1(nn1-4,2)=1;
Bo1(nn1-3,3)=1;
Bo1(nn1-2,4)=1;
Bo1(nn1-1,5)=1;
Bo1(nn1,6)=1;
B1=Bo1';
clear Bo1
% % % % % 1号子结构主模态
FIK1=CCD1;
FIK1(:,motai1+1:nn1)=[];
% % % % % 1号子结构剩余模态求解
FID1=CCD1;
FID1(:,1:motai1)=[];
Kd1=FID1'*kkc1*FID1;
pasi1=FID1*pinv(Kd1)*FID1'*B1';
%% substructure 2
[Freq2,CCD2]=femodal(mmc2,kkc2);
% 对接界面节点自由度
C1=1;C2=6;
nn2=length(kkc2);
% 波尔矩阵
Bo2=zeros(nn2);
Bo2(:,7:nn2)=[];
Bo2(1,1)=1;
Bo2(2,2)=1;
Bo2(3,3)=1;
Bo2(4,4)=1;
Bo2(5,5)=1;
Bo2(6,6)=1;
B2=Bo2';
clear Bo2
% % % % % 用无刚体自由度法求2号子结构主模态
FIK2=CCD2;
FIK2(:,motai2+1:nn2)=[];
% % % % % 2号子结构剩余模态求解
FID2=CCD2(:,motai2+1:nn2);
Kd2=FID2'*kkc2*FID2;
pasi2=FID2*pinv(Kd2)*FID2'*B2';
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%子结构质量阵和刚度阵
Mk1=FIK1'*mmc1*FIK1;
Mg1=pasi1'*mmc1*pasi1;
M1=[Mk1 zeros(motai1,6);zeros(6,motai1) Mg1];
Kk1=FIK1'*kkc1*FIK1;
Kg1=pasi1'*kkc1*pasi1;
K1=[Kk1 zeros(motai1,6);zeros(6,motai1) Kg1];
Mk2=FIK2'*mmc2*FIK2;
Mg2=pasi2'*mmc2*pasi2;
M2=[Mk2 zeros(motai2,6);zeros(6,motai2) Mg2];
Kk2=FIK2'*kkc2*FIK2;
Kg2=pasi2'*kkc2*pasi2;
K2=[Kk2 zeros(motai2,6);zeros(6,motai2) Kg2];
M=[M1 zeros(size(M1,1),size(M2,2));zeros(size(M2,1),size(M1,2)) M2];
K=[K1 zeros(size(K1,1),size(K2,2));zeros(size(K2,1),size(K1,2)) K2];
%% Assemble
t=pinv(B1*pasi1+B2*pasi2);
T21=t*B1*FIK1;
T22=t*B2*FIK2;
T41=-T21;
T42=-T22;
T11=eye(motai1);
T12=zeros(motai1,size(T22,2));
T32=eye(motai2);
T31=zeros(motai2,size(T21,2));
T=[T11 T12;T21 T22;T31 T32;T41 T42];
clear t T11 T12 T21 T22 T31 T32 T41 T42
M=T'*M*T;
K=T'*K*T;
[Freq3,CCD3]=FEMODAL(M,K); |