clear; %清除内存变量
ss1=input('give a data filename:','s');
fp1=fopen(ss1,'r');
ss2=input('give a result filename:','s');
fp2=fopen(ss2,'w');
%1.读入结构数据、建立累积约束表向量、建立结构刚度矩阵
%1.1.结构参数和弹性模量
[m,c]=fscanf(fp1,'%u',[1]); %杆件数
[nj,c]=fscanf(fp1,'%u',[1]); %节点数
[nrj,c]=fscanf(fp1,'%u',[1]); %约束节点
[e,c]=fscanf(fp1,'%e',[1]); %弹性模量
fprintf(fp2,'(1)结构参数和弹性模量\n');
fprintf(fp2,'杆件数 节点数 约束节点数 弹性模量\n');
fprintf(fp2,'%3u %8u %8u %13.3e\n',m,nj,nrj,e);
fprintf(fp2,'\n');
%1.2. 节点坐标
pc(1:2,1:nj)=0
for jx=1:nj
[k,c]=fscanf(fp1,'%u',[1]); %节点号
[pc(:,k),c]=fscanf(fp1,'%f',[2]); %x坐标、y坐标
end
fprintf(fp2,'(2)节点坐标\n');
fprintf(fp2,'节点号 x坐标 y坐标\n');
fprintf(fp2,'%3u %13.3e %13.3e\n',[1:nj;pc(:,1:nj)]);
fprintf(fp2,'\n');
%1.3.杆件标号与截面特性
jj(1:m)=0; jk(1:m)=0; ax(1:m)=0; iz(1:m)=0;
for imx=1:m
[k,c]=fscanf(fp1,'%u',[1]); %杆件号
[jj(k),c]=fscanf(fp1,'%u',[1]); %j端点号
[jk(k),c]=fscanf(fp1,'%u',[1]); %k端点号
[ax(k),c]=fscanf(fp1,'%f',[1]); %截面积
[iz(k),c]=fscanf(fp1,'%f\n',[1]); %杆件号
end
fprintf作用fprintf(fp2,'(3)杆件标号与截面特性\n');
fprintf(fp2,'杆件号 j端点号 k端点号 截面积 截面惯性矩\n');
fprintf(fp2,'%3u %8u%8u %13.3e%13.3e\n',[1:m;jj(1:m);jk(1:m);ax(1:m);iz(1:m)]);
fprintf(fp2,'\n');
%1.4.di和dj取值
di(1:m)=0; dj(1:m)=0;
for imx=1:m
[k1(imx),c]=fscanf(fp1,'%u',[1]);
[di(k1(imx)),c]=fscanf(fp1,'%f',[1]);
[dj(k1(imx)),c]=fscanf(fp1,'%f',[1]);
end
fprintf(fp2,'(4)杆件标号与两端刚域\n');
fprintf(fp2,'杆件号 j端刚域 k端刚域\n');
fprintf(fp2,'%3u %13.3e %13.3e\n',[1:m;di(1:m);dj(1:m)]);
fprintf(fp2,'\n');
%1.5.节点约束情况
fprintf(fp2,'(5)节点约束情况\n');
fprintf(fp2,'受约束点号 x向约束情况 y向约束情况 z约束情况\n');
rl(1:3*nj)=0;
for jx=1:nrj
[k,c]=fscanf(fp1,'%u',[1]); %受约束节点号
[rl(3*k-2:3*k),c]=fscanf(fp1,'%f',[3]); %x向约束情况 y向约束情况 z约束情况
fprintf(fp2,'%5u %12u%12u%12u\n',k,rl(3*k-2:3*k));
end
fprintf(fp2,'\n');
%1.6.建立累积约束表向量
crl(1:3*nj)=0;
[n,nr,crl]=Cu(3*nj,rl);
fprintf(fp2,'(5)非约束位移数 约束位移数\n');
fprintf(fp2,'%10u %10u\n',n,nr);
fprintf(fp2,'\n');
H(1:6,1:6)=0;
%1.7. 装配总节点刚度矩阵
sj(1:3*nj,1:3*nj)=0;
for imx=1:m
H(1:6,1:6)=[1 0 0 0 0 0;0 1 di(k1(imx)) 0 0 0;0 0 1 0 0 0;0 0 0 1 0 0;0 0 0 0 1 -dj(k1(imx));0 0 0 0 0 1];
sm=Psm(pc(:,jj(imx)),pc(:,jk(imx)),ax(imx),iz(imx),e);
sm2=H(1:6,1:6)'*sm*H(1:6,1:6);
RT=PRT(pc(:,jj(imx)),pc(:,jk(imx)));
smd=RT'*sm2*RT;
j=[3*jj(imx)-2:3*jj(imx),3*jk(imx)-2:3*jk(imx)]; %计算杆端位移对应的总节点位移标号
sj(j(1:6),j(1:6))= sj(j(1:6),j(1:6))+smd(1:6,1:6);
end
%1.8. 对应非约束位移、约束位移、将总节点刚度矩阵sj重新排列
s(crl(1:3*nj),crl(1:3*nj))=sj(1:3*nj,1:3*nj);
%2. 读入荷载数据、构建约束杆端力、构建综合节点荷载
%2.1. 读入荷载节点数、集中荷载杆件数、分布荷载杆件数
[nlj,c]=fscanf(fp1,'%u',[1]); %荷载节点数
[nla,c]=fscanf(fp1,'%u',[1]); %集中荷载杆件数
[nlq,c]=fscanf(fp1,'%u',[1]); %分布荷载杆件数
fprintf(fp2,'(7)荷载节点数 集中荷载杆件数 分布荷载杆件数\n');
fprintf(fp2,'%10u %10u %10u\n',nlj,nla,nlq);
fprintf(fp2,'\n');
%2.2. 节点荷载
fprintf(fp2,'(8)节点荷载\n');
fprintf(fp2,'节点号 x向线力 y向线力 z向力偶\n');
ac(1:3*nj,1:1)=0;
for jx=1:nlj
[k,c]=fscanf(fp1,'%u',[1]); %受荷载节点号
[aj,c]=fscanf(fp1,'%f',[3]); %x向线力 y向线力 z向力偶
fprintf(fp2,'%3u %13.3e%13.3e%13.3e\n',k,aj);
ac(3*k-2:3*k)=ac(3*k-2:3*k)+aj(1:3);
end
fprintf(fp2,'\n');
%2.3.杆件上集中荷载
fprintf(fp2,'(9)杆件上集中荷载\n');
fprintf(fp2,'杆件号 x向线力 y向线力 z向力偶 作用点距杆j端距离\n');
aml(1:6,1:m)=0;
for imx=1:nla
[ia,c]=fscanf(fp1,'%u',[1]); %杆件号
[as,c]=fscanf(fp1,'%f',[3]); %x向线力 y向线力 z向力偶
[xa,c]=fscanf(fp1,'%f',[1]); %集中荷载杆端距
fprintf(fp2,'%3u %13.3e%13.3e%13.3e%13.3e\n',ia,as,xa);
amx=PamlC(pc(:,jj(ia)),pc(:,jk(ia)),as,di(k1(imx)),dj(k1(imx)),xa);
aml(:,ia)=aml(:,ia)+amx(:); %累加约束杆端力
end
fprintf(fp2,'\n');
%2.4.杆件上分布荷载
fprintf(fp2,'(10)杆件上分布荷载\n');
fprintf(fp2,'杆件号 起点x向集度 终点x向集度 起点y向集度 终点y向集度 起点距j端距离 终点距j端距离\n');
for imx=1:nlq
[iq,c]=fscanf(fp1,'%u',[1]); %杆件号
[q,c]=fscanf(fp1,'%f',[4]); %起点x向集度 终点x向集度 起点y向集度 终点y向集度
[xq,c]=fscanf(fp1,'%f',[2]); %起点距j端距离 终点距j端距离
fprintf(fp2,'%3u %14.3e%14.3e%14.3e%14.3e%14.3e\n',iq,q,xq);
amx=PamlD(pc(:,jj(iq)),pc(:,jk(iq)),q,xq);
aml(:,iq)=aml(:,iq)+amx(:); %累加约束杆端力
end
fprintf(fp2,'\n');
fclose(fp1); %关闭输入数据文件
fprintf(fp2,'(111)约束杆端力\n');
fprintf(fp2,'杆件号 j/k端号 xm向线约束力 ym向线约束力 zm向约束弯矩\n');
fprintf(fp2,'%3u %6u %14.3e%14.3e%14.3e\n %6u %14.3e%14.3e%14.3e\n',[1:m;jj(1:m);aml(1:3,1:m);jk(1:m);aml(4:6,1:m)]);
版权声明:本站内容均来自互联网,仅供演示用,请勿用于商业和其他非法用途。如果侵犯了您的权益请与我们联系QQ:729038198,我们将在24小时内删除。
发表评论