用伪代码方式写出有限元求解步骤.doc

作业：一

1用伪代码方式写出有限元求解步骤

解：有限元分析流程图：

开始

开始

输入离散模型数据

输入离散模型数据

按选择的单元计算单元刚度矩阵

按选择的单元计算单元刚度矩阵

形成K按单元循环

形成K

按单元循环

按总刚存储模式集成总刚

按总刚存储模式集成总刚

形成P计算单元等效结点荷载

形成P

计算单元等效结点荷载

集成结构结点载荷列阵

集成结构结点载荷列阵

消除K奇异性性引入位移边界条件

消除K奇异性性

引入位移边界条件

求解Ka=P按选定解法解线性方程组

求解Ka=P

按选定解法解线性方程组

结果输出

结果输出

结束

结束

Ⅰ.平面问题程序使用说明

输入数据格式依次如下：

1.NG

NG:结构的结点总数

NE,MC,NX,NB,ND,EO,VO,T

NE:结构单元总数MC:计算控制类型参数=0平面应力=1平面应变NX:作用荷载组数NB：给定位移的个数NB:给定位移的个数ND：结构刚度阵的半带宽EO:弹性模量VO：泊松比T：单元的厚度

NWA,NWE,NWK,NWP,NWD

NWA：单元参数的输出控制参数NWE:单元刚度矩阵的输出控制参数NWK：结构刚度矩阵的输出控制参数NWP：荷载向量的输出控制参数NWD：结点位移的输出控制参数输出控制参数=1输出=0不输出

IJM(3，NE):单元结点编码数组

IJM(1，I),IJM(2，I),IJM(3，I)，第I个三角形单元的节点号，按单元编号顺序填写。

XY(2，NG)：结点坐标数组

XY(1，I):第I个结点的坐标，XY(2，I):第I个结点的Y坐标

按结点编号顺序填写。

MB(2，NB),ZB(NB):给定位移约束的信息数组与值数组

MB(1，I):第I个给定位移所在的结点号

NB(2，I)=1：给定X方向位移

=2：给定Y方向位移

ZB(NB):给定位移值（以坐标正向为正)

7.NF,NP

NF:作用于结点上的集中荷载的个数

NP:作用于均布侧压的单元边数

若NF0，填8

MF(2，NF),ZF(NF):作用于结点上集中荷载的信息组与值数组

MF(1，I):第I个集中荷载作用的结点号

MF(2，I)=1:作用于x方向的集中力

=2:作用y方向的集中力

ZF(NF):作用的集中力值

若NP0，填9

MP(2，NP),ZP(NP):作用于单元边上均布荷载的信息数组与值数组

MF(1,I):第I个均布荷载作用边的起始结点号

MF(2,I):第I个均布荷载作用边的终止结点号，逆时针排列

ZP(NP）：第I个均布荷载值

输入数据格式，建立数据文件，文件名小于12个字符：

NG

NE,MC,NX,NB,ND,E,P,T

NWA,NWE,NWK,NWP,NWD

IJM(3,NE)

XY(2,NG)

MB(2,NB),ZB(NB)

NF,NP

当NF0,则填8

MF(2,NF),ZF(NF)

当NP0,则填9

MP(2,,NP），ZP(NP)

若NX1,即多组荷载情况，重复7～９,若计算多个结构则重复１～９。

结束

平面问题主程序框图

NG,NE,MC.NX,NB,ND,E,P,NWA,NWE,NWK,NWP,NWD……

NG,NE,MC.NX,NB,ND,E,P,NWA,NWE,NWK,NWP,NWD……

NG=0

NG=0

停止

停止

输入原始数据

输入原始数据

是否MC=0

是

否

MC=0

E=EOV=VOE=EO/(1-VO

E=EOV=VO

E=EO/(1-VO2)V=VO(1-VO)

NX1=NX

NX1=NX

计算单元参数

计算单元参数

K=1，NE

K=1，NE

计算单元刚度矩阵

计算单元刚度矩阵

计算总刚度矩阵

计算总刚度矩阵

检验主元

检验主元

荷载信息（NE,NP)

荷载信息（NE,NP)

是否NF0

是

否

NF0

否计算集中力等效结点荷载

否

计算集中力等效结点荷载

是计算侧压等效结点荷载NP0

是

计算侧压等效结点荷载

NP0

引入位移边界条件

引入位移边界条件

高斯消去

高斯消去

计算单元应力

计算单元应力

控制荷载组数

控制荷载组数

否是NX10

否

是

NX10

结束

结束

输入子程序INPUT伪代码

DIMENSIONGIJM(3,NE),XY(2,NG),MB(2,NB),ZB(NB)

READ(5,*)((IJM(I,L),I=1,3),L=1,NE)

READ(5,*)((XY(I,J),I=1,2),J=1,NG)

READ(5,*)(MB(I,L),I=1,2),L=1,NB),(ZB(L),L=1,NB)

WRITE(6，10)

10FORMAT(/10X,‘ELEMENTCODEBLOCKSIJM(2，