板杆式薄壁结构-有限元大作业

南京航空航天大学有限元大作业报告学院：航空宇航学院题目：板杆式薄壁结构组号：15组组员：日期：2013.10一、问题重述15、板杆式薄壁结构已知：P=2000N，a=30cm，各杆截面积f=10mm2，壁板厚t=0.1cm材料参数均相同：E=7×106N/cm2，=0.3求：各杆轴力和板剪流二、初步分析1、本题中所给结构为杆板式薄壁结构，并含有开口，一般用于飞行器机翼结构设计上，可以从一般的板杆结构考虑，做有限元结构分析，同时亦可保持一定准确性。2、该结构为对称结构，受力对称，由此我们可以初步确定该结构上下各相应单位受力相同，方向相反。同此我们也可据此考查所得结果精确与否。8a4aPPPPPP3，结构含有两种单元，完成单元分析后，我们先写出杆的刚度矩阵，再累加上板的刚度矩阵，从而得到总刚度矩阵，最后带入边界条件和受力情况用改进平方根法求解编写计算模型。4，计算模型的检验，以及求解。三、建立计算模型应用直接刚度法对结构进行分析，主要包括以下几方面的工作：1、结构离散化——建立有限元的计算模型；2、单元分析——建立单元刚度矩阵K^e；3、结构分析——形成结构总刚度矩阵K；4、采用“置大数法”对结构进行约束处理；5、利用改进平方根算法求解线性方程组，求得节点位移；6、根据求得的节点位移求出各单元的内力。四、计算模型设计开始输入控制信息(NP,NE,NL,NR,ND,NF)N=NFNPNDF=NDNF输入原始数据形成LD(N)数组调用FLD子程序计算局部坐标下的单元刚度矩阵坐标变化矩阵总体坐标下的单元刚度矩阵调用FKE子程序形成单元的IS(NDF)数组（调用FIS子程序）单元刚度矩阵向总刚度结构矩阵叠加（调用FASUM子程序）进行约束处理调用FR子程序求解线性方程组求得结构节点位移输出节点位移计算局部坐标下的单元刚度矩阵坐标变换矩阵计算单元IS(NDF)数组取出单元在总体坐标下的节点位移计算在局部坐标下单元节点位移计算局部坐标下单元节点力和单元内力输出单元内力结束五、模型符号数据说明1、节点编号：共44个节点2、杆的编号：共56根杆3、板的编号：共28块板六、模型原理1标识符说明【1】控制信息NP：结构节点总数；NE：结构离散杆单元总数；NEE：结构纯剪板单元总数；NL：载荷组数；ND：每个单元的节点数；NF：每个节点的自由度数；N=NFxNP：结构的节点位移总数；NDF=NFxND：每个单元的节点自由度数；【2】输入原始数据1)节点坐标X(NP)：节点的x方向坐标；Y(NP)：节点的y方向坐标；2)单元的特征数据ME(NDNE)：单元节点的总体编号；AE(NE)：单元横截面面积；E0(NE)：单元材料弹性模量；3)边界约束数据NRR(NR)：约束为零的位移所对应的总体位移编号；P(N，NL)：NL组载荷向量；程序中的其他主要标识符LD(N)：存放结构总刚所有主对角元素在一维数组A(NN)中的序号；IS(NDF)：单元节点位移和节点力在结构节点位移列阵和节点载荷列阵中对应的序号；EK(NDF，NDF)：总体坐标下的单元刚度矩阵；A(NN)：结构刚度矩阵下三角变带宽一维压缩存储的数组；NN：数组A的元素总数。【3】输出结果1)结构的节点位移；2)单元内力；2输入说明由于我们题目的单元数比较多，所以我们采用程序自动输入，没有输入文件，我们所有的杆单元的编号都是小号在前3，源程序i.由于我们的程序是程序自动输入信息，所以没有输入文件，所以我们的算例是在程序内部人为切换的；ii.我们的刚度矩阵是采用一维存储，从第一个非零元素开始存储；iii.约束处理，采用置大数法；!*********************************************************!*****************平面薄壁结构静力分析********************!*********************************************************PROGRAMMAIN!*******************************本次课设题目***********************!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!由于本结构单元比较多，不宜采用文件读取输入，故本程序带有专门的节点，单元信息输入程序!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!INTEGER::ME(2,56),MEP(4,28),LD(88)=0,NRR(100),NN,N,NP,NE,ND,NF,NDF,IE,NEE,NR,ISW,NL,ME1(4)INTEGER,ALLOCATABLE::IS(:),ISS(:)REAL,ALLOCATABLE::A(:)REALX(44),Y(44),AE,EK(4,4),E0,G0,V,EKP(8,8),H0,A0,U(88),P(88),AA(8),C,Q(28),EKE(2,2),BL,R(2,4),S(56),CI(2,1),T(4)!***************输入结构信息*************OPEN(8,FILE='INPUTDATA.txt')READ(8,*)NP,NE,NF,NEE,NR,NL!!输入结构信息READ(8,*)AE!!输入杆的截面积READ(8,*)E0!!输入弹性模量READ(8,*)V!!输入vREAD(8,*)H0!!输入板的厚度READ(8,*)BL!!输入杆的长度N=88;G0=E0/(2*(1+V))CALLIMPORTXY(X,Y);!!输入节点坐标信息（详见xy.txt）CALLINPUTME(ME);!!输入杆的信息（详见STICKSME.txt）CALLPLATE(MEP);!!输入板的信息（详见PLATE.txt）CALLFLD(LD,NN,MEP,N,NP,NEE,4,NF)!!输入LD数组（详见LD.txt）CALLINPUTNRR(NRR,NR)!!输入约束信息P=0CALLINPUT(P,N)!!输入外载荷信息!P(1)=1000!****************本次课设题目**************!!*******将所有信息输入到OUTDATA.txt中********OPEN(7,FILE='OUTDATA.TXT')WRITE(7,200)NP,NE,NEE,NL,NR,NF200FORMAT(5X,'************输入信息************'/,8X,'NP=',I4/,8X,'NE=',I4/,8X,'NEE=',I4/,8X,'NL=',I4/,8X,'NR=',I4/,8X,'NF=',I4/)WRITE(7,201)AE,H0,BL,E0,G0201FORMAT(5X,'**********结构属性信息**********'/,8X,'面积:AE=',F4.2/,8X,'壁板厚:H0=',F4.2/,8X,'杆长:BL=',F4.1/,8X,'弹性模量:E0=',F9.1/,8X,'剪切模量:G0=',F9.1/)WRITE(7,202)202FORMAT(5X,'************节点载荷************'/,5X,'序号',8X,'PX',9X,'PY'/)WRITE(7,203)(I,P(2*I-1),P(2*I),I=1,NP)203FORMAT(5X,I3,3X,2F10.2/)WRITE(7,204)204FORMAT(5X,'************约束信息************'/,5X,'序号',8X,'NRR'/)WRITE(7,205)(I,NRR(I),I=1,NR)205FORMAT(5X,I3,8X,I3)WRITE(7,206)206FORMAT(5X,'************节点坐标************'/,5X,'序号',8X,'X',9X,'Y'/)WRITE(7,207)(I,X(I),Y(I),I=1,NP)207FORMAT(5X,I3,3X,2F10.2)WRITE(7,208)208FORMAT(5X,'************杆的编号************'/,5X,'序号',8X,'ME1',9X,'ME2'/)WRITE(7,209)(I,ME(1,I),ME(2,I),I=1,NE)209FORMAT(5X,I3,3X,I8,4X,I8)WRITE(7,210)210FORMAT(5X,'************板的编号************'/,5X,'序号',8X,'MEP1',8X,'MEP2',8X,'MEP3',8X,'MEP4'/)WRITE(7,211)(I,MEP(1,I),MEP(2,I),MEP(3,I),MEP(4,I),I=1,NEE)211FORMAT(5X,I3,3X,I8,4X,I8,4X,I8,4X,I8)WRITE(7,212)NN212FORMAT(5X,'*************LD数组*************'/,5X,'NN=',I4/,5X,'LD(N)=')WRITE(7,213)(LD(I),I=1,N)213FORMAT(5X,5I5)ND=2NDF=ND*NFALLOCATE(IS(NDF))ALLOCATE(A(NN))WRITE(7,214)214FORMAT(5X,'**********杆刚度矩阵**********'/)!!*********计算杆的单元刚度矩阵，并存入一维总刚中***********DOI=1,NEIE=ICALLFIS(ME,IS,NF,ND,NE,IE,NDF)CALLFEK(IE,NP,NE,NDF,X,Y,ME,AE,E0,EK,ND,BL)WRITE(7,(5X,'IE=',I3/))IWRITE(7,215)((EK(J,K),K=1,4),J=1,4)215FORMAT(5X,4F12.2/)CALLFASUM(A,EK,IS,LD,NN,NDF,N)ENDDOOPEN(50,FILE='ANN.TXT')WRITE(50,(2X,5F12.2))(A(I),I=1,NN)ND=4NDF=ND*NFALLOCATE(ISS(NDF))WRITE(7,216)216FORMAT(5X,'**********板刚度矩阵**********'/)!!*********计算板的单元刚度矩阵，并存入一维总刚中***********DOI=1,NEEIE=ICALLFIS(MEP,ISS,NF,ND,NEE,IE,NDF)CALLFKEP(IE,NP,NEE,X,Y,MEP,H0,G0,EKP,A0,AA,C)WRITE(7,(5X,'IE=',I3/))IWRITE(7,217)((EKP(J,K),K=1,8),J=1,8)217FORMAT(5X,8F12.2/)CALLFASUM(A,EKP,ISS,LD,NN,NDF,N)ENDDOWRITE(7,218)218FORMAT(5X,'********一维刚度矩阵********'/,5X,'A(N)=')WRITE(7,219)(A(I),I=1,NN)219FORMAT(5X,6F12.2/)OPEN(50,FILE='ANNP.TXT')WRITE(50,(2X,5F12.2))(A(I),I=1,NN)CALLFR(N,NN,NR,LD,A,NRR)!!置大数法处理约束CALLBAND3(N,1,NN,A,P,LD,ISW)!!用改进平方根法求解位移WRITE(7,220)220FORMAT(5X,'***********节点位移***********'/,5X,