几何组成 fortran source code.pdfVIP

1体系组成分析的计算机方法-fortran90 源码实现 杨晶 （河海大学土木工程系 07 级 江苏 南京） 一、 问题描述 结构静力学课程中介绍的几何组成分析方法有两类。一类是基于三刚片、两刚片和二元体 法则的分析方法，此类方法可以解决常见的结构组成问题，但是技巧性过强，而且并不能解决 所有的问题，不具有一般性；另一类是静力分析的方法，如零载法，但此类方法对于自由度为 0的情况难以判断是否几何可变。下面将讨论如何通过计算机程序解决体系几何组成的判定问 题。 二、 算法分析 首先，必须要将结构数据化，然后才能通过计算机合理的处理数据，这里考虑方便程序 的编写，使用“后处理法”对结构进行编码，即事先不考虑支座的约束。并且规定数据输入的 格式为： 控制参数（一行）：NElem,NJoint,NGlbDOF 结点数据(共 NJoint 行)：x,y,GDOF(1:3) 单元数据（共 NGlbDOF 行）：JointNo(1:2) 注：x，y表示结点在整体坐标系中的坐标，GDOF 表示该结点的结构位移编码（索引码） JointNo(1:2)说明杆件单元两端结点的编号 其次，将要说明，体系几何组成分析实际是一个几何问题。众所周知，自由杆件单元的 杆端位移有 6个自由度，但仔细分析，并不是这六个自由度之间都是独立，它们之间存在一定 的约束关系，这一约束关系可用单元几何约束方程来表示（element geometric constraint equation）,[E]Δ=0,E是一个 3*6的矩阵，Δ为杆端位移向量。 由此可以求出结构中每根杆件单元的的单元几何约束矩阵（element geometric constraint，matrix）,在根据后处理法中的单元定位向量，集装后得到整体几何约束方程 （global geometric constraint equation）[G]Δ=0,G是一个M*N 的矩阵，其中M=单元总 数*3；N=结点位移总数。 根据线性代数理论，[G]Δ=0是一个齐次方程。如果该方程只有平凡解[O]，说明体系不能 有独立的位移，rank(G)=N,体系几何不变；如果方程有无穷多解，则体系是几何可变的， rangk(G)N,任意一组解可以表示为一组线性无关的基础解系（）来表示，每一基础解对应体 系的一个位移模态 进一步，如果体系几何不变，[G]Δ=0中有M个方程，但是只有 rank(G)个未知量是独 立，如果 rank(G)=N=M，独立未知量个数与约束方程个数恰好相等，体系几何不变的，且无 多余约束；如果 rank(G)=NM,显然体系有M-rank(G)个未知量不独立，此时，体系为几何不 变，有M-rank(G)个多余约束的超静定结构。 上面的分析可知，程序的重点算法在整体几何矩阵（global geometric constraint matrix）的集装和该矩阵的求秩。集装按单元定位向量按行集装；rank(G)则通过 Gauss消元法 求得。 算法大致流程如下图所示： 1. 输入结构数据 input 2. 生成结构 spanStructure 3. 生成整体约束矩阵 spanGlbConstraint 4. 求整体约束矩阵的秩 rank 5. 输出整体约束矩阵，判断体系的几何组成 output 三、 算法实现 2算法实现时，在 visual fortran 环境中通过 FORTRAN90 来实现 1) 程序模块图 箭头指向的文件总引用箭头始端的模块，如 constructAnaly 主程序使用了 subprocess,dataBlock matrixProcess 三个模块 程序输入文件为 in.txt;输出文件 GC.txt(GC:global constraint matrix 整体约束矩阵) 2)FORTRAN90 代码 ConstructionAnalysis.f90 program constructAnaly use dataBlock use subProcess use matrixProcess implicit none integer(ikind)r,count,times real(rkind),dimension(:,:),allocatable::GC open(unit=1,file=in.txt) read(1,*)times do count=1,times call input() call spanStructure(joint,element,NElem) allocate(GC(3*NElem,NGlbDOF)); call spanGlbConstraint(GC,joint,element,NElem) r=rankA(GC,3*NElem,NGlbDOF) c