【弹塑性力学】7 有限元程序设计祥解
合集下载
相关主题
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
单元数据定位 XL（NDM, NEN) 单元内节点坐标数组 UL (NDM, NEN) 单元内节点位移数组 LD(NDF*NEN) 方程编号
单元数组的计算
以平面四边形等参单元为例
形函数
Ni
1 4
1
i
1
i
x Ni xi
y Ni yi
u Niui
v Nivi
坐标变换和位移函数都采用同样的形函数
ANSYS，SAP等。多采用FORTRAN语言开发。 3.进一步优化程序代码，采用面向对象的软件设计方法，
加强前后处理技术。 4.适应计算机系统的并行化改进。 5.20世纪90年代中期开始，与CAD/CAM进行集成，向智能
化方向发展。 近15-20年是有限元分析软件商品化的迅速发展阶段。
整体刚度矩阵的形成
SHP(1,I) = SI(I)*(0.5+TI(I)*T) 100 SHP(2,I) = (0.5+SI(I)*S)*TI(I) C C---COMPUTE JACOBIAN TRANSROMATION C---FROM X,Y TO S,T C
DO 200 I = 1,2 DO 200 J = 1,2 XS(I,J) = 0.0 D0 200 K = 1,4 200 XS(I,J) = XS(I,J) + XL(I,K)*SHP(J,K)
DATA SI,TI/-0.5, 0.5, 0.5, -0.5, -0.5, -0.5, 0.5, 0.5 /
C
C---COMPUTE SHAPE FUNCTIONS AND DERIVATIVES
C--- IN NATURAL COORDINATES
C
DO 100 I = 1,4
SHP(3,I) = (0.5+SI(I)*S)*(0.5+TI(I)*T)
和同样多的参数，所以称为等参单元。
Biblioteka Baidu函数导数
Ni,
Ni,
x,
x,
y, y,
Ni,
Ni,
x y
J
N i,
Ni,
x y
于是
N i,
Ni,
x y
J
1
Ni, Ni,
B矩阵
N i, x
Bi 0
N i, y
0
N i, y
N i, x
D矩阵
D11 D12
D12
C SHP(2,I) = Y-direivative fo I –node shape function
C SHP(3,I) = shape function for I-node
C XJ
= Jacobian array
C XSJ = Jacobian determinant
C
DIMENSION SHP(3,4), XL(2,4), XS(2,2), SI(4), TI(4)
D矩阵
D11 D12
D12
D11
0
0
单元刚度矩阵
0
0
D33
kij BiT DB j
Q j DB j
D11N j,x
Q j D12N j,x D33N j,y
D12 N j,y
D11N j,y
D 33 N
j
,
x
kij
N i, xQ11 N i, y Q21
N i, N i,
C---- 2 2 INTERGRATION 2 G = 1./SQRT(3.)
DO 21 I = 1,4 R(I) = G*LR(I) Z(I) = G*LZ(I) 21 W(I) = 1. RETURN C---- 3 3 INTERGRATION 3 G = 1./SQRT(0.6) H = 1./81. DO 31 I = 1,9 R(I) = G*LR(I) Z(I) = G*LZ(I) 31 W(I) = H*LW(I) RETURN END
C C---COMPUTE JACOBIAN DETERMINANT C
XSJ = XS(1,1)*XS(2,2) - XS(1,2)*XS(2,1) C C--- TRANSFORM NATURAL DERIVATIVES TO C--- X,Y DERIVATIVES C
DO 300 I = 1,4 TEMP = (XS(2,2)*SHP(1,I) – XS(2,1)*SHP(2,I))/XSJ SHP(2,I) = (XS(2,2)*SHP(1,I) – XS(2,1)*SHP(2,I))/XSJ 300 SHP(1,I) = TEMP RETURN END
DIMENSIONS C
DATA LR/-1,1,1,-1,0,1,0,-1,0/,LZ/-1,-1,1,1,-1,0,1,0,0/ LW/4/
LINT = L + L GO TO (1,2,3) L C---- 11 INTERGRATION 1 R(1) = 0. Z(1) = 0. W(1) = 4. RETURN
x Q31 x Q31
Ni, xQ12 Ni, yQ32
Ni, yQ22
Ni, xQ32
单元刚度矩阵
kij
e
t
1 -1
1 -1
k ij
| J | dd
用高斯积分转化为求和
kij e t
WiWj | J | [Kij]
ij
SUBROUTINE PGAUSS(L,LINT,R,Z,W) 高斯积分权重 C C---- GAUSS POINTS AND WEIGHTS FOR TWO
有限元程序三大模块
• 数据输入(前处理) •求解计算 • 结果整理与输出(后处理) 统计结果：第1阶段所用时间最长，第3阶段次 之，第2阶段相对较少。
数据输入(前处理)
• 单元几何: 单元编号; 节点坐标 • 材料性质 • 边界条件(约束) • 荷载(体积力, 面积力)
求解计算
• 单元刚度矩阵 • (形函数矩阵，B矩阵，D矩阵，高斯积分) • 单元荷载列阵 • 组集整体刚度矩阵和整体荷载列阵 • 方程组求解
结果整理（后处理）
• 等值线 • 矢量图 • 变形网格图
编程基本原则 结构化，模块化
• 可读性强 变量命名简单直观，适当注释 程序条理清楚，逻辑性强
• 节约内存，计算效率优化
有限元分析软件发展5个阶段
1.1950-1960，航空工业的发展促进了有限元程序的发展。 2.人们致力于多功能通用有限元程序开发。如NASTRAN，
D11
0
0
0
0
D33
SUBROUTINE SHAPEF(S, T, XL, XSJ, SHP)
C--- shape function routine for 4-node isoparametric quadrilateral
C SHP(1,I) = X-direivative fo I –node shape function
单元数组的计算
以平面四边形等参单元为例
形函数
Ni
1 4
1
i
1
i
x Ni xi
y Ni yi
u Niui
v Nivi
坐标变换和位移函数都采用同样的形函数
ANSYS，SAP等。多采用FORTRAN语言开发。 3.进一步优化程序代码，采用面向对象的软件设计方法，
加强前后处理技术。 4.适应计算机系统的并行化改进。 5.20世纪90年代中期开始，与CAD/CAM进行集成，向智能
化方向发展。 近15-20年是有限元分析软件商品化的迅速发展阶段。
整体刚度矩阵的形成
SHP(1,I) = SI(I)*(0.5+TI(I)*T) 100 SHP(2,I) = (0.5+SI(I)*S)*TI(I) C C---COMPUTE JACOBIAN TRANSROMATION C---FROM X,Y TO S,T C
DO 200 I = 1,2 DO 200 J = 1,2 XS(I,J) = 0.0 D0 200 K = 1,4 200 XS(I,J) = XS(I,J) + XL(I,K)*SHP(J,K)
DATA SI,TI/-0.5, 0.5, 0.5, -0.5, -0.5, -0.5, 0.5, 0.5 /
C
C---COMPUTE SHAPE FUNCTIONS AND DERIVATIVES
C--- IN NATURAL COORDINATES
C
DO 100 I = 1,4
SHP(3,I) = (0.5+SI(I)*S)*(0.5+TI(I)*T)
和同样多的参数，所以称为等参单元。
Biblioteka Baidu函数导数
Ni,
Ni,
x,
x,
y, y,
Ni,
Ni,
x y
J
N i,
Ni,
x y
于是
N i,
Ni,
x y
J
1
Ni, Ni,
B矩阵
N i, x
Bi 0
N i, y
0
N i, y
N i, x
D矩阵
D11 D12
D12
C SHP(2,I) = Y-direivative fo I –node shape function
C SHP(3,I) = shape function for I-node
C XJ
= Jacobian array
C XSJ = Jacobian determinant
C
DIMENSION SHP(3,4), XL(2,4), XS(2,2), SI(4), TI(4)
D矩阵
D11 D12
D12
D11
0
0
单元刚度矩阵
0
0
D33
kij BiT DB j
Q j DB j
D11N j,x
Q j D12N j,x D33N j,y
D12 N j,y
D11N j,y
D 33 N
j
,
x
kij
N i, xQ11 N i, y Q21
N i, N i,
C---- 2 2 INTERGRATION 2 G = 1./SQRT(3.)
DO 21 I = 1,4 R(I) = G*LR(I) Z(I) = G*LZ(I) 21 W(I) = 1. RETURN C---- 3 3 INTERGRATION 3 G = 1./SQRT(0.6) H = 1./81. DO 31 I = 1,9 R(I) = G*LR(I) Z(I) = G*LZ(I) 31 W(I) = H*LW(I) RETURN END
C C---COMPUTE JACOBIAN DETERMINANT C
XSJ = XS(1,1)*XS(2,2) - XS(1,2)*XS(2,1) C C--- TRANSFORM NATURAL DERIVATIVES TO C--- X,Y DERIVATIVES C
DO 300 I = 1,4 TEMP = (XS(2,2)*SHP(1,I) – XS(2,1)*SHP(2,I))/XSJ SHP(2,I) = (XS(2,2)*SHP(1,I) – XS(2,1)*SHP(2,I))/XSJ 300 SHP(1,I) = TEMP RETURN END
DIMENSIONS C
DATA LR/-1,1,1,-1,0,1,0,-1,0/,LZ/-1,-1,1,1,-1,0,1,0,0/ LW/4/
LINT = L + L GO TO (1,2,3) L C---- 11 INTERGRATION 1 R(1) = 0. Z(1) = 0. W(1) = 4. RETURN
x Q31 x Q31
Ni, xQ12 Ni, yQ32
Ni, yQ22
Ni, xQ32
单元刚度矩阵
kij
e
t
1 -1
1 -1
k ij
| J | dd
用高斯积分转化为求和
kij e t
WiWj | J | [Kij]
ij
SUBROUTINE PGAUSS(L,LINT,R,Z,W) 高斯积分权重 C C---- GAUSS POINTS AND WEIGHTS FOR TWO
有限元程序三大模块
• 数据输入(前处理) •求解计算 • 结果整理与输出(后处理) 统计结果：第1阶段所用时间最长，第3阶段次 之，第2阶段相对较少。
数据输入(前处理)
• 单元几何: 单元编号; 节点坐标 • 材料性质 • 边界条件(约束) • 荷载(体积力, 面积力)
求解计算
• 单元刚度矩阵 • (形函数矩阵，B矩阵，D矩阵，高斯积分) • 单元荷载列阵 • 组集整体刚度矩阵和整体荷载列阵 • 方程组求解
结果整理（后处理）
• 等值线 • 矢量图 • 变形网格图
编程基本原则 结构化，模块化
• 可读性强 变量命名简单直观，适当注释 程序条理清楚，逻辑性强
• 节约内存，计算效率优化
有限元分析软件发展5个阶段
1.1950-1960，航空工业的发展促进了有限元程序的发展。 2.人们致力于多功能通用有限元程序开发。如NASTRAN，
D11
0
0
0
0
D33
SUBROUTINE SHAPEF(S, T, XL, XSJ, SHP)
C--- shape function routine for 4-node isoparametric quadrilateral
C SHP(1,I) = X-direivative fo I –node shape function