《杆系静力分析》PPT课件.ppt

1、Subroutine Elem_Stiff() 说明 Stiff=0.0 ! 单元刚度清零 Select Case (Type) Case (1) 平面杆系结构单元 Case (2) 空间杆系结构单元 Case Default 出错信息 End Select End Subroutine Elem_Stiff,3.4 杆系结构单元分析子程序,3.4.1 单元刚度总体设计,3.4.2 说明部分设计,Integer,Intent (in) : 入口整型参数 Real(8),Intent(in) : 入口实型参数 Real(8),Intent(out) : 出口实型参数 Real(8) : Work

2、1, Integer : i,j,k, 实型和整型工作变量,3.4 杆系结构单元分析子程序,3.4.3 平面杆系结构设计,Select Case (Plane) Case (1) 平面桁架元素赋值 Case (2) 平面梁柱元素赋值 Case (3) Case Default 出错信息 End Select,3.4.4 空间杆系结构设计,Select Case (Space) Case (1) 空间桁架元素赋值 Case (2) 空间梁柱元素赋值 Case (3) 交叉梁元素赋值 Case Default 出错信息 End Select,3.4 杆系结构单元分析子程序,3.4.5 有关单元等效

3、结点荷载设计和进一步的考虑,1) 单元等效结点荷载设计同仿单元刚度,2) 从各类单元刚度元素的计算，可看到要用到长度、 单元弹性特性、单元截面特性等数据。因此，要 确定存放它们的数据结构。要将它们作为出口,3) 为计算单元等效结点荷载元素，首先要建立各种 荷载情况等效荷载表达式，它们可由积分或载常 数表得到。然后要解决荷载信息的存放结构，也 要将它们作为出口量,4) 单元刚度矩阵、等效结点荷载矩阵都应先清零,4.1 杆系结构整体分析,首先就全刚结点平面刚架进行讨论，然后推广,4.1.1 总的思路,在单元特性搞清后，将单元拼装回去。在结点处 位移自动协调基础上，如果全部结点平衡，则求得 的结点位

4、移将是实际结构的解。因此，整体分析就 是设法建立结点平衡方程,4.1.2 坐标转换,组成结构的杆件可以各个 方向，单元分析对局部坐标， 因此，必须将物理量转为统 一坐标整体坐标,1) 力的转换关系,4.1 杆系结构整体分析,2) 位移转换关系,3) 转换矩阵,转换矩阵是正交矩阵,4.1 杆系结构整体分析,4) 杆端力转换,5) 杆端位移转换,6) 刚度方程的转换,如果记 称为整体单元刚度矩阵,则,这就是整体坐标下的单元刚度方程,本节以后的讨论认为 都是对整体坐标的,4.1 杆系结构整体分析,4.1.3 结点平衡方程的建立,1) 一简单例子（如图,图中有两套编号，红的 是单元杆端编号，黑的是 结

5、构整体编号,1-1) 结点示意,图中蓝色的表示结点荷载（已知），红色的表示 杆端力（未知的）， 、 分别1、2单元杆端力 子矩阵。对1、4结点“荷载”含有未知反力,1-2) 结点平衡,4.1 杆系结构整体分析,从例图可见，其全部结点 平衡方程为,若记,4.1 杆系结构整体分析,式中I、0分别为单位和零矩阵,若引入矩阵记号,则结点平衡方程可改写作,这一结论虽然是由一个例子得到的，但是显然 对一切结构都是成立的。问题在于不同结构，A 矩阵是不同的,4.1 杆系结构整体分析,4.1.4 杆端位移用结点位移来表示,仍以简单例子来说明,若记,由结点、杆端位移的协调条件，可得 、 的对应关系为,式中 AT

6、是前面力关系A的转置，因此AT称 为位移转换矩阵,4.1 杆系结构整体分析,4.1.5 整体刚度方程结点平衡,若记,引入位移转换关系，则,4.1 杆系结构整体分析,4.1.6 整体刚度矩阵的建立,若将A按单元分成图示 三个子矩阵,则,由此可见，整体刚度矩阵可由各单元整体刚度矩阵装配累加得到。为说明如何装配，先将单元刚度矩阵进行分割,则由矩阵乘法可证明，AikiAiT的结果是，将 刚度矩阵子矩阵按整体结点码 r 、s 送整体刚度矩阵 相应位置。这一装配规则称为“对号入座,4.1 杆系结构整体分析,1) 任意结构情况,上面结论是通过具体例子（全刚结点平面刚架）得到的，由虚位移原理或势能原理进行整体

7、分析（见讲义），可得任意结构其结论同此例,2) 结点位移编号,如果按结点顺序，对结点非零位移进行依次编号，这一序号称作结点位移码。为便于计算机处理并减少结构刚度矩阵的阶次，将零位移的号码变为零,对图示三铰刚架，当仅 用一种单元（梁柱自由是 单元）时结点位移编号如 图所示,3) 单元定位向量,按单元局部结点码顺序， 将结点位移码排成的向量，称作单元的定位向量,4.1 杆系结构整体分析,对图示刚架各单元的定位向量为,(0,0,1,3,4,5,(0,0,2,10,11,12,(3,4,5,6,7,8,(6,7,9,10,11,12,(0,0,1,2,3,(0,0,6,7,8,(1,2,3,4,5,(

8、4,5,6,7,8,4) 按单元定位向量集装刚度矩阵和综合荷载,前面说明的是分块子矩阵集装,下面说明如何按定 为向量来集装,如果如图所是采用各种不同的单元（一端有铰 ），则定位向量为,如何获得 带铰的单元刚 度矩阵和等 效荷载矩阵,定位向量,4.1 杆系结构整体分析,4-1) 刚度集装,以 3 单元为例来说明,定位向量,根据单元 局部位移码 和定位向量 的对应关系 用定位向量 位移码送元 素,根据单元局部位移码和定位向量的对应关系用定 位向量位移码送元素，定位向量元素为零时不送,4.1 杆系结构整体分析,4-2) 荷载集装,以 4 单元为例来说明,此结论同样适 用于刚度集装,4.1 杆系结构整

9、体分析,4.1.7 整体分析总结,1) 对局部坐标和整体坐标不一致的单元，要对刚 度、荷载进行坐标转换,2) 需对“结构”进行结点、位移的局部和整体编 号,4) 集装所得整体刚度矩阵是对称、带状稀疏矩阵， 当支撑条件能限制刚体位移时，矩阵非奇异,3) 根据单元局部位移码和定位向量的对应关系用 定位向量位移码送元素，定位向量元素为零时 不送。据此可集装、累加得到整体刚度矩阵,5) 综合荷载由两部分组成，因此首先要将直接作 用结点的荷载按结点位移码送入，如果还有单 元等效荷载，再按定位向量集装、累加,4.1 杆系结构整体分析,8) 如果有某位移码方向弹性支撑，需进行将弹簧 刚度送入位移码对应的对角

10、线元素位置累加,9) 如果有某位移码方向已知支撑位移，需进行将 “边界条件处理”。具体做法以后介绍,7) 整体刚度方程实质是全部结点的平衡条件,6) 刚度矩阵带状稀疏，其带宽取决于结点、位移 编码。 最大半带宽=定位向量中最大元素差+1,4.1 杆系结构整体分析,2.5.8 边界条件的处理,10) 当用虚位移或势能原理作整体分析时（势能为 例），应变能为单元应变能之和，外力势能为 单元外力势能之和+结点外力势能。全部杆端力 的外力势能彼此抵消,1) 乘大数法,设第 i 个位移为已知值 a ，N =108 或更大的数。 乘大数法是将刚度矩阵Kii改为NKii，将Ri改为 Na,请考虑 为什麽这样

11、做能使 边界条件得到满足,2) 置换法（划零置1,设第 i 个位移为已知值 a,4.1 杆系结构整体分析,上述置换工作量大一些，显然可看出边界条件得 到精确满足,4.1 杆系结构整体分析,3) 关于斜边界的处理,如图示意的斜支座情况，有多种处理方案,3-1) 通过单元的坐标转换来处理,图示有斜支座单元，r 结点处 以倾角 - 来进行坐标转换，也 即在r 结点处整体坐标为图示 xy,r,3-2) 通过增加一个单元来处理,图示有斜支座单元，r 结点处沿 y 方向增加一个 刚结的单元，此单元有“无穷大”的抗拉刚度、但 没有抗弯刚度。单元长度可任意,3-3) 对整体刚度矩阵进行处理（参见匡文起教材,2

12、.5 杆系结构整体分析,2.5.10 总刚度矩阵的存储与对应解法,因为总刚度矩阵对称、带状稀疏，利用这一特 点可减少存储、提高解算效率,零元素,零元素,等 半 带 存 储,变 带 宽 一 维 存 储,到P.55,2.5 杆系结构整体分析,目前一般都用变带宽存储，下面结合程序说明 存储和解法。首先介绍一些F 90的语法,定义导出类型,导出类型结点 type : typ_Joint real : x,y ！坐标 integer : GDOF(3)！整体位移码 end type typ_Joint,1) 有关 F90 语法 导出类型 新特性,用结点导出类型作为成员导出单元类型: type : typ

13、_Element integer : JointNo(2) ！结点编号 type(typ_Joint) : Node(2) ！结点 integer : GlbDOF(6) ! 定位向量 real : EA,EI end type typ_Element,type (typ_Element) : Elem(5) ！定义5个单元类型 ! 对单元i的端点j的x,y,GDOF(1:3)的赋值 Elem(i)%Node(j) = Joint(Elem(i)%JointNo(j),由导出类型定义新类型,由导出类型定义变量,real : A(5),B(5,10),C(5) B=0.0 !对B清零 A=1.0

14、 !对A赋1: A(i)=1.0, i=1,5 C=A+2 !数组与标量运算: A(1:5)+(/ 2,2,2,2,2 /) A=C+A !数组与数组运算（同形） C=sqrt(A) !数组的函数运算: C(i)=sqrt(A(i), i=1,5 数组内部函数： dot_product(vector_a,vector_b) !点积 如： dot_product(/1,2,3/),(/2,3,4/) 的值为20 （待续,有关 F90 语法 数组运算与赋值,matmul(matrix_a,matrix_b) 矩阵相乘 如：locEDisp = matmul(T,glbEDisp) transpos

15、e(matrix) 矩阵转置 如：glbEDisp = matmul(transpose(T),locEDisp) size(array,dim) 求数组第dim维的长度 dim为可选变元: size(a,dim=2) 若array为一维时，可不用dim。 sum(array,dim,mask) 数组元素求和 dim,mask为可选变元；mask=条件表达式 sum(a(1:10) 对a的1到10元素求和 sum(a(1:10),mask=a0) 对a(1:10)中大于0的元素求和,续,有关 F90 语法 where结构 新特性,例 where (C0) C = 0 A = B*D end w

16、here where (C0) A = B end where,定义 where (数组关系表达式) 数组赋值语句 elsewhere 数组赋值语句 . end where,规则: 1)同形数组; 2)不许嵌套; 3)最多二分叉,有关 F90 语法 cycle和exit语句 新特性 用在do循环中 cycle 作下一个循环步 exit 跳出循环，执行end do后一条语句,等效例 do . if(.not.cond1)then . if(cond2) goto 5 . end if end do 5,用法 do . if (cond1) cycle . if (cond2) exit . en

17、d do,有关 F90 语法 数组构造函数spread,语法 spread(数组名,dim,ncopies) 将数组沿dim维方向复制ncopies形成新数组 dim,ncopies 整型、位置变元、关键字变元 （若按位置引用，可略关键字,例: （仅限一维数组） 1） spread(one,dim=1,ncopies=3) spread(one,1,3) spread(one,ncopies=3,dim=1,2） ELocVec(1:6)=(/ 1,0,3,4,5,0 /) spread(ELocVec,dim=1,ncopies=3,3,spread(A(2,2:),dim=1,ncopie

18、s=2,如果dim=2呢,有关 F90 语法 指针pointer,pointer是变量的属性，可以指向相同类型的变量; 被指向的目标必须具有target属性或pointer属性 可以将指针变量理解为别名、称号 real,target : a,b,EDisp(6) ! 可被指针所指 real,pointer : p1,p2 ! 称号：班长、课代表 ！p1,p2是指针，可以指向实型数据 real,pointer : G(:) ! 先进集体 ！G是指针，可以指向一维实型数组,指针是一种“称号”，上述声明语句建立了“称号”，但并未“授予”某个变量这个称号，因此是指向“空”，并未占用内存 a = 3.0

19、 p1 = a ! p1指向a ! 称号p1授予a，a的数据有两个名：固定名a和流动名p1; 既可用p1也可用a（p1 班长，a 张三） a = 4.0 ! a的值变为4.0,p1也变为4.0 p1 = b ! 班长换人了 G = EDisp ! 先进集体有了得主 ！EDisp(:)的长度就是G(:)的长度，用G和用EDisp同样效果,又如： real,target : a,b real,pointer: p,q a = 3.14 b = 2.0 p = a ! p = a = 3.14 q = b ! q = b = 2.0 q = p ! ( q指向的数据b ) = ( p指向的数据a )

20、 ！即： 所有 = 3.14,指针可以指向有名的数据区，也可以指向无名的数据区 real,pointer: b(:,:) integer : n read(*,*) n allocate (b(n,n) !指针指向了一个刚开辟的数组 ! 以下可以当作常规数组用 b(1,1) = 1.1 b(1,2) = 1.2 . deallocate (b,有关 F90 语法 用指针建立动态数组,指针与allocatable数组的区别 具备allocatable数组的所有功能 还可以用在导出类型中，例如整体刚度矩阵的变带宽存储： type : typ_Kcol !整体刚度矩阵K的列 real(8),poin

21、ter : row(:) !该列的行元素 end type . type (typ_Kcol),allocatable : Kcol(:) allocate (Kcol(NGlbDOF) ! 分配了NGlbDOF列 . allocate (Kcol(5)%row(3:5) ! 第5列只用3至5行,1) NElem,NJoint, NGlbDOF, NJLoad, NELoad 单元数 结点数 总自由度数 结点荷载数 单元荷载数 Joint - 结点 NJoint行 (2) Joint%X,Joint%Y, GDOF(1:3) X 坐标 Y坐标 结点位移码 Elem - 单元 NElem行 (3

22、) JointNo1,JointNo2, EA,EI 起点号 终点号 刚度 JLoad - 结点荷载 NJLoad行 (4) JointNo, LodDOF, LodVal 作用点号 局部位移码 荷载值 ELoad - 单元荷载 NELoad行 (5) ElemNo, Indx, Pos, LodVal 单元号 类型号 位置 荷载值,2) 某程序输入数据说明,3,5,7,1,1 0,0, 0,0,0 0,4, 1,2,3 4,4, 4,5,6 4,4, 4,5,7 4,0, 0,0,0 1,2, 1.0E9, 16 2,3, 1.0E9, 24 5,4, 1.0E9, 12 2,1,10.0E3 ! 结点2，自由度1，值为10E3 2,1,4,-4.0E3 ！单元2，均布，长4米，值为-4E3,2-1) 数据文件例子,EA= 109 N,1,2,坐标,位移码,3,4,5,结点号,EA,EI,read(5,*) NElem,NJoint,NGlbDOF,NJLoad,NELoad allocate (Joint(NJoint) allocat