搜档网
当前位置:搜档网 › 有限元作业:三角形单元求解资料

有限元作业:三角形单元求解资料

有限元作业:三角形单元求解资料
有限元作业:三角形单元求解资料

《有限元作业》

年级2015级

学院机电工程学院

专业名称

班级学号

学生姓名

2016年05月

如下图所示为一受集中力P作用的结构,弹性模量E为常量,泊松比V=1/6,厚度为I=1。按平面应力问题计算,运用有限元方法,分别采用三角形及四边形单元求解,求节点位移及单元应力(要求三角形单元数量不少于4个,四边形单元不少于2个)

图(一)

图(二)三角形单元求解

图(三)四边形单元求解

(1)如图划分三角形单元,工分成四个分别为 ④

(2)如图分别进行编号1、2、3、4、5、6,并建立坐标系

(3)编程进行求解,得出结果,其中假设力P=2000N

调用Triangle2D3Node_Stiffness函数,求出单元刚度矩阵

k1 =

1.0e+06 *

7.2857 -3.0000 -2.1429 0.8571 -5.1429 2.1429

-3.0000 7.2857 2.1429 -5.1429 0.8571 -2.1429 -2.1429 2.1429 2.1429 0 0 -2.1429

0.8571 -5.1429 0 5.1429 -0.8571 0

-5.1429 0.8571 0 -0.8571 5.1429 0

2.1429 -2.1429 -2.1429 0 0 2.1429 k2 =

1.0e+06 *

5.1429 0 -5.1429 0.8571 0 -0.8571

0 2.1429 2.1429 -2.1429 -2.1429 0

-5.1429 2.1429 7.2857 -3.0000 -2.1429 0.8571

0.8571 -2.1429 -3.0000 7.2857 2.1429 -5.1429

0 -2.1429 -2.1429 2.1429 2.1429 0

-0.8571 0 0.8571 -5.1429 0 5.1429 k3 =

1.0e+06 *

2.1429 0 -2.1429 -2.1429 0 2.1429

0 5.1429 -0.8571 -5.1429 0.8571 0

-2.1429 -0.8571 7.2857 3.0000 -5.1429 -2.1429 -2.1429 -5.1429 3.0000 7.2857 -0.8571 -2.1429

0 0.8571 -5.1429 -0.8571 5.1429 0

2.1429 0 -2.1429 -2.1429 0 2.1429 k4 =

1.0e+06 *

2.1429 0 -2.1429 -2.1429 0 2.1429

0 5.1429 -0.8571 -5.1429 0.8571 0

-2.1429 -0.8571 7.2857 3.0000 -5.1429 -2.1429

-2.1429 -5.1429 3.0000 7.2857 -0.8571 -2.1429

0 0.8571 -5.1429 -0.8571 5.1429 0

2.1429 0 -2.1429 -2.1429 0 2.1429

调用Triangle2D3Node_Assembly函数,求出总体刚度矩阵

求出的节点位移

U =

-0.0004

0.0008

0.0005

0.0010

0.0007

0.0023

-0.0007

0.0026

调用Triangle2D3Node_Stress函数,求出应力,S1、S2、S3、中求出的分别为Sx,Sy,Sxy

S1 =

1.0e+03 *

-4.4086

-0.7348

3.5914

S2 =

1.0e+03 *

4.4086

-0.6405

0.4086

S3 =

1.0e+03 *

1.8907

-1.0601

2.1093

S4 =

1.0e+03 *

-1.8907

2.1093

1.8907

二、

(1)如图划分四边形单元,工分成四个分别为

(2)如图分别进行编号1、2、3、4、5、6,并建立坐标系(3)编程进行求解,得出结果,其中假设力P=2000N

调用Quad2D4Node_Stiffness函数,求出单元刚度矩阵

调用Quad2D4Node_Assembly函数,求出求出总体刚度矩阵

求出节点位移

U =

0.0012

0.0017

-0.0012

0.0017

0.0016

0.0049

-0.0017

0.0052

调用Quad2D4Node_Stress函数,求出单元应力中的的S1、S2、S3分别为Sx,Sy,Sxy应力分量

S1 =

1.0e+03 *

0.0000

-0.2478

2.0000

S2 =

1.0e+07 *

0.6856

4.1135

-1.7137

程序附录

一、

1、三角形单元总程序:

E=1e7;

NU=1/6;

t=1;

ID=1;

%调用Triangle2D3Node_Stiffness函数,求出单元刚度矩阵

k1=Triangle2D3Node_Stiffness(E,NU,t,0,1,0,0,1,1,ID)

k2=Triangle2D3Node_Stiffness(E,NU,t,0,0,1,0,1,1,ID)

k3=Triangle2D3Node_Stiffness(E,NU,t,1,1,1,0,2,0,ID)

k4=Triangle2D3Node_Stiffness(E,NU,t,2,0,2,1,1,1,ID)

%调用Triangle2D3Node_Assembly函数,求出总体刚度矩阵

KK = zeros(12,12);

KK=Triangle2D3Node_Assembly(KK,k1,1,2,3);

KK=Triangle2D3Node_Assembly(KK,k2,2,4,3);

KK=Triangle2D3Node_Assembly(KK,k3,3,4,5);

KK=Triangle2D3Node_Assembly(KK,k4,5,6,3)

% 边界条件的处理及刚度方程求解

k=KK(5:12,5:12)

p=[0;0;0;0;0;0;0;2000]

u=k\p

%支反力的计算

U=[0;0;0;0;u] %为节点位移

P=KK*U

%调用Triangle2D3Node_Strain函数,求出应变SN1、SN2、SN3中求出的分别为SNx,SNy,SNxy

u1=[U(1);U(2);U(3);U(4);U(5);U(6)];

u2=[U(3);U(4);U(7);U(8);U(5);U(6)];

u3=[U(5);U(6);U(7);U(8);U(9);U(10)];

u4=[U(9);U(10);U(11);U(12);U(5);U(6)];

SN1=Triangle2D3Node_Strain(0,1,0,0,1,1,u1)

SN2=Triangle2D3Node_Strain(0,0,1,0,1,1,u2)

SN3=Triangle2D3Node_Strain(1,1,1,0,2,0,u3)

SN4=Triangle2D3Node_Strain(2,0,2,1,1,1,u4)

%调用Triangle2D3Node_Stress函数,求出应力,S1、S2、S3、中求出的分别为Sx,Sy,Sxy

u1=[U(1);U(2);U(3);U(4);U(5);U(6)];

u2=[U(3);U(4);U(7);U(8);U(5);U(6)];

u3=[U(5);U(6);U(7);U(8);U(9);U(10)];

u4=[U(9);U(10);U(11);U(12);U(5);U(6)];

S1=Triangle2D3Node_Stress(E,NU,0,1,0,0,1,1,u1,ID)

S2=Triangle2D3Node_Stress(E,NU,0,0,1,0,1,1,u2,ID)

S3=Triangle2D3Node_Stress(E,NU,1,1,1,0,2,0,u3,ID)

S4=Triangle2D3Node_Stress(E,NU,2,0,2,1,1,1,u4,ID)

2、求刚度矩阵程序

function k=Triangle2D3Node_Stiffness(E,NU,t,xi,yi,xj,yj,xm,ym,ID)

%该函数计算单元的刚度矩阵

%输入弹性模量E,泊松比NU,厚度t

%输入三个节点i、j、m的坐标xi,yi,xj,yj,xm,ym

%输入平面问题性质指示参数ID(1为平面应力,2为平面应变)

%输出单元刚度矩阵k(6X6)

%---------------------------------------------------------------

A = (xi*(yj-ym) + xj*(ym-yi) + xm*(yi-yj))/2;

betai = yj-ym;

betaj = ym-yi;

betam = yi-yj;

gammai = xm-xj;

gammaj = xi-xm;

gammam = xj-xi;

B = [betai 0 betaj 0 betam 0 ;

0 gammai 0 gammaj 0 gammam ;

gammai betai gammaj betaj gammam betam]/(2*A);

if ID == 1

D = (E/(1-NU*NU))*[1 NU 0 ; NU 1 0 ; 0 0 (1-NU)/2];

elseif ID == 2

D = (E/(1+NU)/(1-2*NU))*[1-NU NU 0 ; NU 1-NU 0 ; 0 0 (1-2*NU)/2]; end

k= t*A*B'*D*B;

3、求整体刚度矩阵

function z = Triangle2D3Node_Assembly(KK,k,i,j,m)

%该函数进行单元刚度矩阵的组装

%输入单元刚度矩阵k

%输入单元的节点编号I、j、m

%输出整体刚度矩阵KK

%---------------------------------------------------------------

DOF(1)=2*i-1;

DOF(2)=2*i;

DOF(3)=2*j-1;

DOF(4)=2*j;

DOF(5)=2*m-1;

DOF(6)=2*m;

for n1=1:6

for n2=1:6

KK(DOF(n1),DOF(n2))= KK(DOF(n1),DOF(n2))+k(n1,n2);

end

end

4、求应变程序

function strain=Triangle2D3Node_Strain(xi,yi,xj,yj,xm,ym,u)

%该函数计算单元的应变

%输入三个节点i、j、m的坐标xi,yi,xj,yj,xm,ym

%输入单元的位移列阵u(6X1)

%输出单元的应力strain(3X1),由于它为常应变单元,则单元的应变分量为SNx,SNy,SNz

%---------------------------------------------------------------

A = (xi*(yj-ym) + xj*(ym-yi) + xm*(yi-yj))/2;

betai = yj-ym;

betaj = ym-yi;

betam = yi-yj;

gammai = xm-xj;

gammaj = xi-xm;

gammam = xj-xi;

B = [betai 0 betaj 0 betam 0 ;

0 gammai 0 gammaj 0 gammam ;

gammai betai gammaj betaj gammam betam]/(2*A);

strain = B*u;

5、求应力程序

function stress=Triangle2D3Node_Stress(E,NU,xi,yi,xj,yj,xm,ym,u,ID)

%该函数计算单元的应力

%输入弹性模量E,泊松比NU,厚度t

%输入三个节点i、j、m的坐标xi,yi,xj,yj,xm,ym

%输入平面问题性质指示参数ID(1为平面应力,2为平面应变),单元的位移列阵u(6X1)

%输出单元的应力stress(3X1),由于它为常应力单元,则单元的应力分量为Sx,Sy,Sxy

%---------------------------------------------------------------

A = (xi*(yj-ym) + xj*(ym-yi) + xm*(yi-yj))/2;

betai = yj-ym;

betaj = ym-yi;

betam = yi-yj;

gammai = xm-xj;

gammaj = xi-xm;

gammam = xj-xi;

B = [betai 0 betaj 0 betam 0 ;

0 gammai 0 gammaj 0 gammam ;

gammai betai gammaj betaj gammam betam]/(2*A);

if ID == 1

D = (E/(1-NU*NU))*[1 NU 0 ; NU 1 0 ; 0 0 (1-NU)/2];

elseif ID == 2

D = (E/(1+NU)/(1-2*NU))*[1-NU NU 0 ; NU 1-NU 0 ; 0 0 (1-2*NU)/2];

stress = D*B*u;

二、

1、四边形单元总程序:

E=1e7;

NU=1/6;

h=1;

ID=1;

%调用Quad2D4Node_Stiffness函数,求出单元刚度矩阵

k1= Quad2D4Node_Stiffness(E,NU,h,0,1,0,0,1,0,1,1,ID)

k2= Quad2D4Node_Stiffness(E,NU,h,1,0,2,0,2,1,1,1,ID)

%调用Quad2D4Node_Assembly函数,求出求出总体刚度矩阵

KK=zeros(12,12);

KK= Quad2D4Node_Assembly(KK,k1,1,2,3,4);

KK= Quad2D4Node_Assembly(KK,k2,3,5,6,4)

% 边界条件的处理及刚度方程求解

k=KK(5:12,5:12)

p=[0;0;0;0;0;0;0;2000]

u=k\p

%支反力的计算

U=[0;0;0;0;u] %为节点位移

P=KK*U

%调用Quad2D4Node_Stress函数,求出单元应力中的的S1、S2、S3分别为Sx,Sy,Sxy应力分量

u1=[U(1);U(2);U(3);U(4);U(5);U(6);U(7);U(8)];

u2=[U(5);U(6);U(9);U(10);U(11);U(12);U(7);(8)];

S1= Quad2D4Node_Stress(E,NU,0,1,0,0,1,0,1,1,u1,ID)

S2= Quad2D4Node_Stress(E,NU,1,0,2,0,2,1,1,1,u2,ID)

2、求刚度矩阵程序

function k= Quad2D4Node_Stiffness(E,NU,h,xi,yi,xj,yj,xm,ym,xp,yp,ID)

%该函数计算单元的刚度矩阵

%输入弹性模量E,泊松比NU,厚度h

%输入4个节点i、j、m、p的坐标xi,yi,xj,yj,xm,ym,xp,yp

%输入平面问题性质指示参数ID(1为平面应力,2为平面应变)

%输出单元刚度矩阵k(8X8)

%---------------------------------------------------------------

syms s t;

a = (yi*(s-1)+yj*(-1-s)+ym*(1+s)+yp*(1-s))/4;

b = (yi*(t-1)+yj*(1-t)+ym*(1+t)+yp*(-1-t))/4;

c = (xi*(t-1)+xj*(1-t)+xm*(1+t)+xp*(-1-t))/4;

d = (xi*(s-1)+xj*(-1-s)+xm*(1+s)+xp*(1-s))/4;

B1 = [a*(t-1)/4-b*(s-1)/4 0 ; 0 c*(s-1)/4-d*(t-1)/4 ;

c*(s-1)/4-d*(t-1)/4 a*(t-1)/4-b*(s-1)/4];

B2 = [a*(1-t)/4-b*(-1-s)/4 0 ; 0 c*(-1-s)/4-d*(1-t)/4 ;

c*(-1-s)/4-d*(1-t)/4 a*(1-t)/4-b*(-1-s)/4];

B3 = [a*(t+1)/4-b*(s+1)/4 0 ; 0 c*(s+1)/4-d*(t+1)/4 ;

c*(s+1)/4-d*(t+1)/4 a*(t+1)/4-b*(s+1)/4];

B4 = [a*(-1-t)/4-b*(1-s)/4 0 ; 0 c*(1-s)/4-d*(-1-t)/4 ;

c*(1-s)/4-d*(-1-t)/4 a*(-1-t)/4-b*(1-s)/4];

Bfirst = [B1 B2 B3 B4];

Jfirst = [0 1-t t-s s-1 ; t-1 0 s+1 -s-t ;

s-t -s-1 0 t+1 ; 1-s s+t -t-1 0];

J = [xi xj xm xp]*Jfirst*[yi ; yj ; ym ; yp]/8;

B = Bfirst/J;

if ID == 1

D = (E/(1-NU*NU))*[1 NU 0 ; NU 1 0 ; 0 0 (1-NU)/2];

elseif ID == 2

D = (E/(1+NU)/(1-2*NU))*[1-NU NU 0 ; NU 1-NU 0 ; 0 0 (1-2*NU)/2]; end

BD = J*transpose(B)*D*B;

r = int(int(BD, t, -1, 1), s, -1, 1);

z = h*r;

k = double(z);

3、求总体刚度矩阵程序

function z = Quad2D4Node_Assembly(KK,k,i,j,m,p)

%该函数进行单元刚度矩阵的组装

%输入单元刚度矩阵k,单元的节点编号i、j、m、p

%输出整体刚度矩阵KK

%---------------------------------------------------------------

DOF(1)=2*i-1;

DOF(2)=2*i;

DOF(3)=2*j-1;

DOF(4)=2*j;

DOF(5)=2*m-1;

DOF(6)=2*m;

DOF(7)=2*p-1;

DOF(8)=2*p;

for n1=1:8

for n2=1:8

KK(DOF(n1),DOF(n2))= KK(DOF(n1),DOF(n2))+k(n1,n2);

end

end

z=KK;

4、求应力程序

function stress= Quad2D4Node_Stress(E,NU,xi,yi,xj,yj,xm,ym,xp,yp,u,ID) %该函数计算单元的应力

%输入弹性模量E,泊松比NU,厚度h,

%输入4个节点i、j、m、p的坐标xi,yi,xj,yj,xm,ym,xp,yp,

%输入平面问题性质指示参数ID(1为平面应力,2为平面应变)

%输入单元的位移列阵u(8X1)

%输出单元的应力stress(3X1)

%由于它为常应力单元,则单元的应力分量为Sx,Sy,Sxy

%---------------------------------------------------------------

syms s t;

a = (yi*(s-1)+yj*(-1-s)+ym*(1+s)+yp*(1-s))/4;

b = (yi*(t-1)+yj*(1-t)+ym*(1+t)+yp*(-1-t))/4;

c = (xi*(t-1)+xj*(1-t)+xm*(1+t)+xp*(-1-t))/4;

d = (xi*(s-1)+xj*(-1-s)+xm*(1+s)+xp*(1-s))/4;

B1 = [a*(t-1)/4-b*(s-1)/4 0 ; 0 c*(s-1)/4-d*(t-1)/4 ;

c*(s-1)/4-d*(t-1)/4 a*(t-1)/4-b*(s-1)/4];

B2 = [a*(1-t)/4-b*(-1-s)/4 0 ; 0 c*(-1-s)/4-d*(1-t)/4 ;

c*(-1-s)/4-d*(1-t)/4 a*(1-t)/4-b*(-1-s)/4];

B3 = [a*(t+1)/4-b*(s+1)/4 0 ; 0 c*(s+1)/4-d*(t+1)/4 ;

c*(s+1)/4-d*(t+1)/4 a*(t+1)/4-b*(s+1)/4];

B4 = [a*(-1-t)/4-b*(1-s)/4 0 ; 0 c*(1-s)/4-d*(-1-t)/4 ;

c*(1-s)/4-d*(-1-t)/4 a*(-1-t)/4-b*(1-s)/4];

Bfirst = [B1 B2 B3 B4];

Jfirst = [0 1-t t-s s-1 ; t-1 0 s+1 -s-t ;

s-t -s-1 0 t+1 ; 1-s s+t -t-1 0];

J = [xi xj xm xp]*Jfirst*[yi ; yj ; ym ; yp]/8;

B = Bfirst/J;

if ID == 1

D = (E/(1-NU*NU))*[1 NU 0 ; NU 1 0 ; 0 0 (1-NU)/2];

elseif ID == 2

D = (E/(1+NU)/(1-2*NU))*[1-NU NU 0 ; NU 1-NU 0 ; 0 0 (1-2*NU)/2]; end

str1 = D*B*u;

str2 = subs(str1, {s,t}, {0,0});

stress = double(str2);

有限元分析大作业报告

有限元分析大作业报告 试题1: 一、问题描述及数学建模 图示无限长刚性地基上的三角形大坝,受齐顶的水压力作用,试用三节点常应变单元和六节点三角形单元对坝体进行有限元分析,并对以下几种计算方案进行比较: (1)分别采用相同单元数目的三节点常应变单元和六节点三角形单元计算; (2)分别采用不同数量的三节点常应变单元计算; (3)当选常应变三角单元时,分别采用不同划分方案计算。 该问题属于平面应变问题,大坝所受的载荷为面载荷,分布情况及方向如图所示。 二、采用相同单元数目的三节点常应变单元和六节点三角形单元计算 1、有限元建模 (1)设置计算类型:两者因几何条件和载荷条件均满足平面应变问题,故均取Preferences 为Structural (2)选择单元类型:三节点常应变单元选择的类型是Solid Quad 4 node182;六节点三角形单元选择的类型是Solid Quad 8 node183。因研究的问题为平面应变问题,故对Element behavior(K3)设置为plane strain。 (3)定义材料参数:弹性模量E=2.1e11,泊松比σ=0.3 (4)建几何模型:生成特征点;生成坝体截面 (5)网格化分:划分网格时,拾取lineAB和lineBC,设定input NDIV 为15;拾取lineAC,设定input NDIV 为20,选择网格划分方式为Tri+Mapped,最后得到600个单元。

(6)模型施加约束:约束采用的是对底面BC 全约束。大坝所受载荷形式为Pressure ,作用在AB 面上,分析时施加在L AB 上,方向水平向右,载荷大小沿L AB 由小到大均匀分布。以B 为坐标原点,BA 方向为纵轴y ,则沿着y 方向的受力大小可表示为: }{*980098000)10(Y y g gh P -=-==ρρ 2、 计算结果及结果分析 (1) 三节点常应变单元 三节点常应变单元的位移分布图 三节点常应变单元的应力分布图

平面三角形单元有限元程序设计

. 一、题目 如图1所示,一个厚度均匀的三角形薄板,在顶点作用沿板厚方向均匀分布的竖向载荷。已知:P=150N/m ,E=200GPa ,=0.25,t=0.1m ,忽略自重。试计算薄板的位移及应力分布。 要求: 1. 编写有限元计算机程序,计算节点位移及单元应力。(划分三角形 单元,单元数不得少于30个); 2. 采用有限元软件分析该问题(有限元软件网格与程序设计网格必 须一致),详细给出有限元软件每一步的操作过程,并将结果与程序计算结果进行对比(任选取三个点,对比位移值); 3. 提交程序编写过程的详细报告及计算机程序; 4. 所有同学参加答辩,并演示有限元计算程序。 有限元法中三节点三角形分析结构的步骤如下: 1)整理原始数据,如材料性质、荷载条件、约束条件等,离散结构并进行单元编码、结点编码、结点位移编码、选取坐标系。 2)单元分析,建立单元刚度矩阵。 3)整体分析,建立总刚矩阵。 4)建立整体结构的等效节点荷载和总荷载矩阵 5)边界条件处理。 6)解方程,求出节点位移。 7)求出各单元的单元应力。 8)计算结果整理。 一、程序设计 网格划分 如图,将薄板如图划分为6行,并建立坐标系,则

刚度矩阵的集成 建立与总刚度矩阵等维数的空矩阵,已变单元刚度矩阵的集成。 由单元分析已知节点、单元的排布规律,继而通过循环计算求得每个单元对应的节点序号。 通过循环逐个计算:(1)每个单元对应2种单元刚度矩阵中的哪一种; (2)该单元对应总刚度矩阵的那几行哪几列 (3)将该单元的单元刚度矩阵加入总刚度矩阵的对应行列 循环又分为3层循环:(1)最外层:逐行计算 (2)中间层:该行逐个计算 (3)最里层:区分为第 奇/偶 数个计算 单元刚度的集成:[ ][][][][][]' '''''215656665656266256561661e Z e e e Z e Z e e e e k k k K k k k k k k +?++=? =?==?==?=?????? 边界约束的处理:划0置1法 X Y P X Y P

有限元法的基本思想及计算 步骤

有限元法的基本思想及计算步骤 有限元法是把要分析的连续体假想地分割成有限个单元所组成的组合体,简称离散化。这些单元仅在顶角处相互联接,称这些联接点为结点。离散化的组合体与真实弹性体的区别在于:组合体中单元与单元之间的联接除了结点之外再无任何关联。但是这种联接要满足变形协调条件,即不能出现裂缝,也不允许发生重叠。显然,单元之间只能通过结点来传递内力。通过结点来传递的内力称为结点力,作用在结点上的荷载称为结点荷载。当连续体受到外力作用发生变形时,组成它的各个单元也将发生变形,因而各个结点要产生不同程度的位移,这种位移称为结点位移。在有限元中,常以结点位移作为基本未知量。并对每个单元根据分块近似的思想,假设一个简单的函数近似地表示单元内位移的分布规律,再利用力学理论中的变分原理或其他方法,建立结点力与位移之间的力学特性关系,得到一组以结点位移为未知量的代数方程,从而求解结点的位移分量。然后利用插值函数确定单元集合体上的场函数。显然,如果单元满足问题的收敛性要求,那么随着缩小单元的尺寸,增加求解区域内单元的数目,解的近似程度将不断改进,近似解最终将收敛于精确解。 用有限元法求解问题的计算步骤比较繁多,其中最主要的计算步骤为: 1)连续体离散化。首先,应根据连续体的形状选择最能完满地描述连续体形状的单元。常见的单元有:杆单元,梁单元,三角形单元,矩形单元,四边形单元,曲边四边形单元,四面体单元,六面体单元以及曲面六面体单元等等。其次,进行单元划分,单元划分完毕后,要将全部单元和结点按一定顺序编号,每个单元所受的荷载均按静力等效原理移植到结点上,并在位移受约束的结点上根据实际情况设置约束条件。 2)单元分析。所谓单元分析,就是建立各个单元的结点位移和结点力之间的关系式。现以三角形单元为例说明单元分析的过程。如图1所示,三角形有三个结点i,j,m。在平面问题中每个结点有两个位移分量u,v和两个结点力分量F x,F y。三个结点共六个结点位移分量可用列

(完整版)有限元大作业matlab---课程设计例子

有限元大作业程序设计 学校:天津大学 院系:建筑工程与力学学院 专业:01级工程力学 姓名:刘秀 学号:\\\\\\\\\\\ 指导老师:

连续体平面问题的有限元程序分析 [题目]: 如图所示的正方形薄板四周受均匀载荷的作用,该结构在边界 上受正向分布压力, m kN p 1=,同时在沿对角线y 轴上受一对集中压 力,载荷为2KN ,若取板厚1=t ,泊松比0=v 。 [分析过程]: 由于连续平板的对称性,只需要取其在第一象限的四分之一部分参加分析,然后人为作出一些辅助线将平板“分割”成若干部分,再为每个部分选择分析单元。采用将此模型化分为4个全等的直角三角型单元。利用其对称性,四分之一部分的边界约束,载荷可等效如图所示。

[程序原理及实现]: 用FORTRAN程序的实现。由节点信息文件NODE.IN和单元信息文件ELEMENT.IN,经过计算分析后输出一个一般性的文件DATA.OUT。模型基本信息由文件为BASIC.IN生成。 该程序的特点如下: 问题类型:可用于计算弹性力学平面问题和平面应变问题 单元类型:采用常应变三角形单元 位移模式:用用线性位移模式 载荷类型:节点载荷,非节点载荷应先换算为等效节点载荷 材料性质:弹性体由单一的均匀材料组成 约束方式:为“0”位移固定约束,为保证无刚体位移,弹性体至少应有对三个自由度的独立约束 方程求解:针对半带宽刚度方程的Gauss消元法

输入文件:由手工生成节点信息文件NODE.IN,和单元信息文件ELEMENT.IN 结果文件:输出一般的结果文件DATA.OUT 程序的原理如框图:

有限单元法与有限元分析

有限单元法与有限元分析 1.有限单元法 在数学中,有限元法(FEM,Finite Element Method)是一种为求解偏微分方程边值问题近似解的数值技术。求解时对整个问题区域进行分解,每个子区域都成为简单的部分,这种简单部分就称作有限元。它通过变分方法,使得误差函数达到最小值并产生稳定解。类比于连接多段微小直线逼近圆的思想,有限元法包含了一切可能的方法,这些方法将许多被称为有限元的小区域上的简单方程联系起来,并用其去估计更大区域上的复杂方程。它将求解域看成是由许多称为有限元的小的互连子域组成,对每一单元假定一个合适的(较简单的)近似解,然后推导求解这个域总的满足条件(如结构的平衡条件),从而得到问题的解。这个解不是准确解,而是近似解,因为实际问题被较简单的问题所代替。由于大多数实际问题难以得到准确解,而有限元不仅计算精度高,而且能适应各种复杂形状,因而成为行之有效的工程分析手段。 随着电子计算机的发展,有限单元法是迅速发展成一种现代计算方法。它是50年代首先在连续体力学领域--飞机结构静、动态特性分析中应用的一种有效的数值分析方法,随后很快广泛的应用于求解热传导、电磁场、流体力学等连续性问题。 1.1.有限元法分析本质 有限元法分析计算的本质是将物体离散化。即将某个工程结构离散为由各种单元组成的计算模型,这一步称作单元剖分。离散后单元与单元之间利用单元的节点相互连接起来;单元节点的设置、性质、数目等应视问题的性质,描述变形形态的需要和计算精度而定(一般情况单元划分越细则描述变形情况越精确,即越接近实际变形,但计算量越大)。所以有限元中分析的结构已不是原有的物体或结构物,而是同新材料的由众多单元以一定方式连接成的离散物体。这样,用有限元分析计算所获得的结果只是近似的。如果划分单元数目非常多而又合理,则所获得的结果就与实际情况相符合。 1.2.特性分析 1)选择位移模式: 在有限单元法中,选择节点位移作为基本未知量时称为位移法;选择节点力作为基本未知量时称为力法;取一部分节点力和一部分节点位移作为基本未知量时称为混合法。位移法易于实现计算自动化,所以,在有限单元法中位移法应用范围最广。 当采用位移法时,物体或结构物离散化之后,就可把单元总的一些物理量如

有限元分析大作业试题

有限元分析习题及大作业试题 要求:1)个人按上机指南步骤至少选择习题中3个习题独立完成,并将计算结果上交; 2)以小组为单位完成有限元分析计算; 3)以小组为单位编写计算分析报告; 4)计算分析报告应包括以下部分: A、问题描述及数学建模; B、有限元建模(单元选择、结点布置及规模、网格划分方 案、载荷及边界条件处理、求解控制) C、计算结果及结果分析(位移分析、应力分析、正确性分 析评判) D、多方案计算比较(结点规模增减对精度的影响分析、单 元改变对精度的影响分析、不同网格划分方案对结果的 影响分析等) E、建议与体会 4)11月1日前必须完成,并递交计算分析报告(报告要求打印)。

习题及上机指南:(试题见上机指南) 例题1 坝体的有限元建模与受力分析 例题2 平板的有限元建模与变形分析 例题1:平板的有限元建模与变形分析 计算分析模型如图1-1 所示, 习题文件名: plane 0.5 m 0.5 m 0.5 m 0.5 m 板承受均布载荷:1.0e 5 P a 图1-1 受均布载荷作用的平板计算分析模型 1.1 进入ANSYS 程序 →ANSYSED 6.1 →Interactive →change the working directory into yours →input Initial jobname: plane →Run 1.2设置计算类型 ANSYS Main Menu : Preferences →select Structural → OK 1.3选择单元类型 ANSYS Main Menu : Preprocessor →Element T ype →Add/Edit/Delete →Add →select Solid Quad 4node 42 →OK (back to Element T ypes window) → Options… →select K3: Plane stress w/thk →OK →Close (the Element T ype window) 1.4定义材料参数 ANSYS Main Menu : Preprocessor →Material Props →Material Models →Structural →Linear →Elastic →Isotropic →input EX:2.1e11, PRXY :0.3 → OK 1.5定义实常数 ANSYS Main Menu: Preprocessor →Real Constant s… →Add … →select T ype 1→ OK →input THK:1 →OK →Close (the Real Constants Window)

有限元2-弹性力学平面问题有限单元法(2.1三角形单元,2.2几个问题的讨论)综述

第2章 弹性力学平面问题有限单元法 2.1 三角形单元(triangular Element) 三角形单元是有限元分析中的常见单元形式之一,它的优点是: ①对边界形状的适应性较好,②单刚形式及其推导比较简单,故首先介绍之。 一、结点位移和结点力列阵 设右图为从某一结构中取出的一典型三角形单元。 在平面应力问题中,单元的每个结点上有沿x 、y 两个方向的力和位移,单元的结点位移列阵规定为: 相应结点力列阵为: (式2-1-1) 二、单元位移函数和形状函数 前已述及,有限单元法是一种近似方法,在单元分析中,首先要求假定(构 造)一组在单元内有定义的位移函数作为近似计算的基础。即以结点位移为已知量,假定一个能表示单元内部(包括边界)任意点位移变化规律的函数。 构造位移函数的方法是:以结点(i,j,m)为定点。以位移(u i ,v i ,…u m v m )为定点上的函数值,利用普通的函数插值法构造出一个单元位移函数。 在平面应力问题中,有u,v 两个方向的位移,若假定单元位移函数是线性的,则可表示成: (,)123 u u x y x y ααα==++ 546(,)v v x y x y ααα==++ (2-1-2)a 式中的6个待定常数α1 ,…, α6 可由已知的6个结点位移分量(3个结点的坐标) {}??? ?? ?????=????? ???? ?????????????=m j i m e d d d d m j j i v u v u v u i {} i i j j m X Y X (2-1-1)Y X Y i e j m m F F F F ?? ?? ???? ???? ??==??????????????????

有限元大作业

风电主轴承有限元分析 XXX 摘要:基于有限元法在接触问题中的应用,对风电主轴承进行非线性分析。以轴承外圈的内表面和内圈的外表面为目标面,以滚子为接触面创建接触对分析滚子的接触应力情况。最大应力值出现在滚子边缘出,对最大承载滚子环向接触应力分析表明,有限元分析结果与理论计算结果相近,验证了利用有限元法分析风电主轴承应力状态的可行性。 关键词:风电主轴承;接触应力;有限元分析 0 引言 随着传统能源的日益枯竭以及环境污染问题愈发严重,风能作为一种清洁的的可再生能源近些年受到越来越多的关注。风力发电技术已广泛运用于世界各地。一些发达国家风力发电产业已得到了迅猛发展,技术日趋成熟,并开始走向产业化规模化发展阶段[1-3]。 风电主轴承是风力发电机重要的组成部分。其结构形式图下图1所示。据统计,如今安装的所有风力发电机中,采用主轴轴承支撑原理的占总数的75-80%[4],这种支撑是轴承内圈安装在旋转的主轴上,外圈固定在单独的轴承座上,相对于圆锥滚子轴承或圆柱滚子轴承来说,主轴轴承位置处轴产生变形,需要轴承具有一定的调心作用,所以都采用了调心滚子轴承。近年来由于计算机技术的飞速发展,轴承的受力分析计算已经普遍采用有限元分析的方法,能够准确合理地解决轴承复杂的非线性接触问题,为轴承的分析和计算提供了一种新的方法,成为未来的一个发展方向。在机械设备的设计过程中,对受力较大且复杂的零件进行受力分析,校核其整体和局部强度并进行合理的布局设计,是为了防止因应力过大而导致在实际工作中损坏或寿命降低[5]。本文主要运用ANSYS Workbench有限元软件对风电主轴承进行静力学计算,分析轴承内部结构参数对轴承载荷分布和最大接触应力的影响规律。 图1 风电主轴承结构及安装图 1 有限元分析过程 1.1 风电轴承有限元分析基本步骤 不同的物理性质和数学模型的问题,有限元法求解的基本步骤是相同的,只不过 具体公式推导和运算求解不尽相同。有限元分析求解问题的基本计算步骤[6]: 1.问题及求解域定义; 2.求解域离散化; 3.确定状态变量及控制方法; 4.单元推导;

有限元法基础试题

有限元法基础试题(A ) 一、填空题(5×2分) 1.1单元刚度矩阵e T k B DBd Ω = Ω? 中,矩阵B 为__________,矩阵D 为___________。 1.2边界条件通常有两类。通常发生在位置完全固定不能转动的情况为_______边界,具体指定有限的非零值位移的情况,如支撑的下沉,称为_______边界。 1.3内部微元体上外力总虚功: ()(),,,,e x x xy y bx xy x y y by d W F u F v dxdy δστδτσδ??=+++++??+(),,,,x x y y xy y x u v u u dxdy σδσδτδδ??+++??的表达式中,第一项为____________________的虚功,第二项为____________________的虚功。 1.4弹簧单元的位移函数1N +2N =_________。 1.5 ij k 数学表达式:令j d =_____,k d =_____,k j ≠,则力i ij F k =。 二、判断题(5×2分) 2.1位移函数的假设合理与否将直接影响到有限元分析的计算精度、效率和可靠性。( ) 2.2变形体虚功原理适用于一切结构(一维杆系、二维板、三位块体)、适用于任何力学行为的材料(线性和非线性),是变形体力学的普遍原理。 ( ) 2.3变形体虚功原理要求力系平衡,要求虚位移协调,是在“平衡、协调”前提下功的恒等关系。 ( ) 2.4常应变三角单元中变形矩阵是x 或y 的函数。 ( ) 2.5 对称单元中变形矩阵是x 或y 的函数。 ( ) 三、简答题(26分) 3.1列举有限元法的优点。(8分) 3.2写出有限单元法的分析过程。(8分) 3.3列出3种普通的有限元单元类型。(6分) 3.4简要阐述变形体虚位移原理。(4分) 四、计算题(54分) 4.1对于下图所示的弹簧组合,单元①的弹簧常数为10000N/m ,单元②的弹簧常数为20000N/m ,单元③的弹簧常数为10000N/m ,确定各节点位移、反力以及单元②的单元力。(10分) 4.2对于如图所示的杆组装,弹性模量E 为10GPa ,杆单元长L 均为2m ,横截面面积A 均为2×10-4m 2,弹簧常数为2000kN/m ,所受荷载如图。采用直接刚度法确定节点位移、作用力和单元②的应力。(10分)

有限元法中的几个基本概念

诚信·公平·开放·共赢 Loyalty Fair Opening Win-win 有限元法中的几个基本概念 有限元法是把要分析的连续体假想地分割成有限个单元所组成的组合体,简称离散化。 这些单元仅在顶角处相互联接,称这些联接点为结点。 离散化的组合体与真实弹性体的区别在于:组合体中单元与单元之间的联接除了结点之外再无任何关联。但是这种联接要满足变形协调条件,即不能出现裂缝,也不允许发生重叠。显然,单元之间只能通过结点来传递内力。 通过结点来传递的内力称为结点力,作用在结点上的荷载称为结点荷载。当连续体受到外力作用发生变形时,组成它的各个单元也将发生变形,因而各个结点要产生不同程度的位移,这种位移称为结点位移。 在有限元中,常以结点位移作为基本未知量。并对每个单元根据分块近似的思想,假设一个简单的函数近似地表示单元内位移的分布规律,再利用力学理论中的变分原理或其他方法,建立结点力与位移之间的力学特性关系,得到一组以结点位移为未知量的代数方程,从而求解结点的位移分量。然后利用插值函数确定单元集合体上的场函数。显然,如果单元满足问题的收敛性要求,那么随着缩小单元的尺寸,增加求解区域内单元的数目,解的近似程度将不断改进,近似解最终将收敛于精确解。 附:FELAC 2.0软件简介 FELAC 2.0采用自定义的有限元语言作为脚本代码语言,它可以使用户以一种类似于数学公式书写和推导的方式,非常自然和简单的表达待解问题的微分方程表达式和算法表达式,并由生成器解释产生完整的并行有限元计算C程序。 FELAC 2.0的目标是通过输入微分方程表达式和算法之后,就可以得到所有有限元计算的程序代码,包含串行程序和并行程序。该系统采用一种语言(有限元语言)和四种技术(对象技术、组件技术、公式库技术生成器技术)开发而成。并且基于FELAC 1.0的用户界面,新版本扩充了工作目录中右键编译功能、命令终端输入功能,并且丰富了文本编辑功能,改善了用户的视觉体验,方便用户快速便捷的对脚本或程序进行编辑、编译与调试。其中并行版在前后处理上进行了相应的改进。

ansys有限元分析大作业

ansys有限元分析大作业

有限元大作业 设计题目: 单车的设计及ansys有限元分析 专业班级: 姓名: 学号: 指导老师: 完成日期: 2016.11.23

单车的设计及ansys模拟分析 一、单车实体设计与建模 1、总体设计 单车的总体设计三维图如下,采用pro-e进行实体建模。 在建模时修改proe默认单位为国际主单位(米千克秒 mks) Proe》文件》属性》修改

2、车架 车架是构成单车的基体,联接着单车的其余各个部件并承受骑者的体重及单车在行驶时经受各种震动和冲击力量,因此除了强度以外还应有足够的刚度,这是为了在各种行驶条件下,使固定在车架上的各机构的相对位置应保持不变,充分发挥各部位的功能。车架分为前部和后部,前部为转向部分,后部为驱动部分,由于受力较大,所有要对后半部分进行加固。

二、单车有限元模型 1、材料的选择 单车的车身选用铝合金(6061-T6)T6标志表示经过热处理、时效。 其属性如下: 弹性模量:) .6+ 90E (2 N/m 10 泊松比:0.33 质量密度:) 3 2.70E+ N/m (2 抗剪模量:) 60E .2+ N/m (2 10 屈服强度:) .2+ (2 75E 8 N/m 2、单车模型的简化 为了方便单车的模拟分析,提高电脑的运算

效率,可对单车进行初步的简化;单车受到的力的主要由车架承受,因此必须保证车架能够有足够的强度、刚度,抗振的能力,故分析的时候主要对车架进行分析。简化后的车架如下图所示。 3、单元体的选择 单车车架为实体故定义车架的单元类型为实体单元(solid)。查资料可以知道3D实体常用结构实体单元有下表。 单元名称说明 Solid45 三维结构实体单元,单元由8个节点定义,具有塑性、蠕变、应力刚化、 大变形、大应变功能,其高阶单元是 solid95

基于Matlab语言的按平面三角形单元划分的结构有限元程序设计模板

基于Matlab语言的按平面三角形单元划分的结构有限元程序设计 专业:建筑与土木工程 班级:建工研12-2 姓名:韩志强 学号: 471220580

基于Matlab语言的按平面三角形单元划分 结构有限元程序设计 一、有限单元发及Matlab语言概述 1. 有限单元法 随着现代工业、生产技术的发展,不断要求设计高质量、高水平的大型、复杂和精密的机械及工程结构。为此目的,人们必须预先通过有效的计算手段,确切的预测即将诞生的机械和工程结构,在未来工作时所发生的应力、应变和位移因此,需要寻求一种简单而又精确的数值分析方法。有限单元法正是适应这种要求而产生和发展起来的一种十分有效的数值计算方法。 有限元法把一个复杂的结构分解成相对简单的“单元”,各单元之间通过结点相互连接。单元内的物理量由单元结点上的物理量按一定的假设内插得到,这样就把一个复杂结构从无限多个自由度简化为有限个单元组成的结构。我们只要分析每个单元的力学特性,然后按照有限元法的规则把这些单元“拼装”成整体,就能够得到整体结构的力学特性。 有限单元法基本步骤如下: (1)结构离散:结构离散就是建立结构的有限元模型,又称为网格划分或单元划分,即将结构离散为由有限个单元组成的有限元模型。在该步骤中,需要根据结构的几何特性、载荷情况等确定单元体内任意一点的位移插值函数。 (2)单元分析:根据弹性力学的几何方程以及物理方程确定单元的刚度矩阵。 (3)整体分析:把各个单元按原来的结构重新连接起来,并在单元刚度矩阵的基础上确定结构的总刚度矩阵,形成如下式所示的整体有限元线性方程: {}[]{}δ F=① K 式中,{}F是载荷矩阵,[]K是整体结构的刚度矩阵,{}δ是节点位移矩阵。 (4)载荷移置:根据静力等效原理,将载荷移置到相应的节点上,形成节点载荷矩阵。 (5)边界条件处理:对式①所示的有限元线性方程进行边界条件处理。 (6)求解线性方程:求解式①所示的有限元线性方程,得到节点的位移。在该步骤中,若有限元模型的节点越多,则线性方程的数量就越多,随之有限元分析的计算量也将越大。 (7)求解单元应力及应变根据求出的节点位移求解单元的应力和应变。

有限元法大作业

有限元法大作业 一平面刚架的程序 用Visual C++编制的平面刚架的源程序如下: ///////////////////////////////////////////////////////程序开始////////////////////////////////////////////////////////////////// #include"iostream.h" #include"math.h" #include"stdlib.h" #include"conio.h" //***************** //声明必要变量 //***************** #define PI 3.141592654 int NE; //单元数 int NJ; //节点数 int NZ; //支承数 int NPJ; //有节点载荷作用的节点数 int NPF; //非节点载荷数 int HZ; //载荷码 int E; //单元码 int fangchengshu; double F[303]; //各节点等效总载荷数值 int dym_jdm[100][2]; //单元码对应的节点码:dym_jdm[][0], dym_jdm[][1]分为前后节点总码 int zhichengweizhi[300]; //记录支持节点作用点的数组 int fjzhzuoyongdanyuan[100]; //非节点载荷作用单元 int fjzhleixing[100]; //非节点载荷类型:1-均布,2-垂直集中,3-平行集中,4-力偶,5-角度集中 double fjzhzhi[100]; //非节点载荷的值 double fjzhzuoyongdian[100]; //非节点载荷在各竿的作用点 double fjzhjiaodu[100]; //非节点载荷作用角度 int jdzhzuoyongdian[100]; //节点载荷作用的节点数组 double jiedianzaihe[101][3];//节点载荷值,其jiedianzaihe[][0]-- jiedianzaihe[][2]分别为U, V, M double zhengtigangdu[303][303]; //整体刚度数组 double changdu[100]; //各单元竿长数组 double jiaodu[100]; //各单元角度数组 double tanxingmoliang[100]; //各单元弹性模量数组 double J_moliang[100]; //各单元J模量数组 double mianji[100]; //各单元面积数组 double weiyi[303]; //记录各个节点位移的数组 double dy_weiyi[100][6]; //各个单元在局部坐标系中的位移数组dy_weiyi[i][0]-dyweiyi[i][6]分别为第i+1单元的u1,v1,@1,u2,v2,@2 double dy_neili[100][6];//各个单元在局部坐标系中的固端内力dy_weiyi[i][0]-dyweiyi[i][6]分别为第i+1单元的U1,V1,M1,U2,V2,M2 double gan_neili[100][6];//各个单元的竿端内力数组,gan_neili[i][6]表示第i+1单元的6内力. //*******************

有限元法基础重点归纳(精)

1、有限元这种数值计算方法起源于20世纪50年代中期航空工程中飞机结构的矩阵分析。 2、有限单元法的基本思想:在力学模型上将一个原来连续的物体离散成为有限个具有一定 大小的单元,这些单元仅在有限个节点上相连接,并在节点上引进等效力以代替实际作用于单元上的外力。 3、节点:网格间相互连接的点。 4、边界:网格与网格的交界线。 5、有限元的优点:①理论基础简明,物理概念清晰,且可在不同的水平上建立起对该法的 理解②具有灵活性和适用性,应用范围极为广泛③该法在具体推导运算中,广泛采用了矩阵方法。 6、有限单元法分类(从选择基本未知量的角度:位移法(以节点位移为基本未知量,通用 性广、力法(以节点力、混合法(一部分以节点位移,另一部分以节点力 7、有限元法分析计算的基本步骤:①结构的离散化②单元分析(选择位移模式,建立单元 刚度方程,计算等效节点力③整体分析④求解方程,得出节点位移⑤由节点位移计算单元的应变与应力。 8、单元划分:将某个机械结构划分为由各种单元组成的计算模型。 9、有限元法基本近似性------几何近似。

10、弹性力学的任务:分析弹性体在受外力作用并处于平衡状态下产生的应力、应变和位移状态及其相互关系等。 11、弹性力学假设所研究的物体是连续的、完全弹性的、均匀的、各向同性的、微小变形的和无初应力的 12、外力:体力(分布在物体体积内的力---重力、惯性力、电磁力面力(分布在物体表面上的力---流体压力、接触力、风力 13、应力:物体受外力作用,或由于温度有所改变,其内部发生的内力。σ={ σx σy σz τx τy τz } = [σx σy σz τx τy τz ]T 14、应变:物体受到外力作用时,其形状发生改变时的形变。---长度和角度。 ε={ εx εy εz γx γy γz } = [εx εy εz γx γy γz ]T 15、位移:弹性体在载荷作用下,不仅会发生形变,还将产生位移,即弹性体位置 的移动。 δ={u v w }=[u v w ]T 16:、变形协调条件:设想在变形前,把弹性体分为许多微小立方单元体。变形后,每个单元体都产生任意变形而变成一些六面体。可能发生这样的情况,这些六面体

重庆大学研究生有限元大作业教学内容

重庆大学研究生有限 元大作业

课程研究报告 科目:有限元分析技术教师:阎春平姓名:色学号: 2 专业:机械工程类别:学术 上课时间: 2015 年 11 月至 2016 年 1 月 考生成绩: 阅卷评语: 阅卷教师 (签名)

有限元分析技术作业 姓名: 色序号: 是学号: 2 一、题目描述及要求 钢结构的主梁为高160宽100厚14的方钢管,次梁为直径60厚10的圆钢管(单位为毫米),材料均为碳素结构钢Q235;该结构固定支撑点位于左右两端主梁和最中间。主梁和次梁之间是固接。试对在垂直于玻璃平面方向的2kPa 的面载荷(包括玻璃自重、钢结构自重、活载荷(人员与演出器械载荷)、风载荷等)作用下的舞台进行有限元分析。 二、题目分析 根据序号为069,换算得钢结构框架为11列13行。由于每个格子的大小为1×1(单位米),因此框架的外边框应为11000×13000(单位毫米)。 三、具体操作及分析求解 1、准备工作 执行Utility Menu:File → Clear&start new 清除当前数据库并开始新的分析,更改文件名和文件标题,如图1.1。选择GUI filter,执行 Main Menu: Preferences → Structural → OK,如图1.2所示

图1.1清除当前数据库并开始新的分析 图1.2 设置GUI filter 2、选择单元类型。 执行Main Menu: Preprocessor →Element Type →Add/Edit/Delete →Add→ select→ BEAM188,如图2.1。之后点击OK(回到Element Types window) →Close

平面三角形单元有限元程序设计

平面三角形单元有限元程序设计 P 9 m 9 m 一、题目 如图1所示,一个厚度均匀的三角形薄板,在顶点作用沿板厚方向均匀分布的竖向载荷。已知:P=150N/m,E=200GPa,=0、25,t=0、1m,忽略自重。试计算薄板的位移及应力分布。 要求: 1.编写有限元计算机程序,计算节 点位移及单元应力。(划分三角形单元,单元数不得少于30个); 2.采用有限元软件分析该问题(有 限元软件网格与程序设计网格必须一致),详细给出有限元软件每一步的操作过程,并将结果与程序计算结果进行对比(任选取三个点,对比位移值); 3.提交程序编写过程的详细报告及计算机程序; 4.所有同学参加答辩,并演示有限元计算程序。 有限元法中三节点三角形分析结构的步骤如下: 1)整理原始数据,如材料性质、荷载条件、约束条件等,离散结构并进行单元编码、结点编码、结点位移编码、选取坐标系。 2)单元分析,建立单元刚度矩阵。 3)整体分析,建立总刚矩阵。 4)建立整体结构的等效节点荷载与总荷载矩阵 5)边界条件处理。 6)解方程,求出节点位移。 7)求出各单元的单元应力。 8)计算结果整理。

一、程序设计 网格划分 如图,将薄板如图划分为6行,并建立坐标系,则 刚度矩阵的集成 建立与总刚度矩阵等维数的空矩阵,已变单元刚度矩阵的集成。 由单元分析已知节点、单元的排布规律,继而通过循环计算求得每个单元对应的节点序号。 通过循环逐个计算:(1)每个单元对应2种单元刚度矩阵中的哪一种; (2)该单元对应总刚度矩阵的那几行哪几列 (3)将该单元的单元刚度矩阵加入总刚度矩阵的对应行列 循环又分为3层循环:(1)最外层:逐行计算 (2)中间层:该行逐个计算 (3)最里层:区分为第 奇/偶 数个计算 X Y P X Y P 节点编号 单元编号

有限元单元法复习资料

1.1有限单元法中“离散”的含义是什么?有限单元法是如何将具 有无限自由度的连续介质的问题转变为有限自由度问题的?位移有限单元法的标准化程式是怎样的?(1)离散的含义即将结构离散化,即用假想的线或面将连续体分割成数目有限的单元,并在其上设定有限个节点;用这些单元组成的单元集合体代替原来的连续体,而场函数的节点值将成为问题的基本未知量。(2)给每个单元选择合适的位移函数或称位移模式来近似地表示单元内位移分布规律,即通过插值以单元节点位移表示单元内任意点的位移。因节点位移个数是有限的,故无限自由度问题被转变成了有限自由度问题。(3)有限元法的标准化程式:结构或区域离散,单元分析,整体分析,数值求解。 1.2单元刚度矩阵和整体刚度矩阵各有哪些性质?各自的物理意义是什么?两者有何区别? 单元刚度矩阵的性质:对称性、奇异性(单元刚度矩阵的行列式为零)。整体刚度矩阵的性质:对称性、奇异性、稀疏性。 单元刚度矩阵Kij物理意义Kij即单元节点位移向量中第j个自由度发生单位位移而其他位移分量为零时,在第i个自由度方向引起的节点力。 整体刚度矩阵K中每一列元素的物理意义是:要迫使结构的某节点位移自由度发生单位位移,而其他节点位移都保持为零的变形状态,在所有个节点上需要施加的节点荷载。 2.1 为了使计算结果能够收敛于精确解,位移函数需要满足什么条件?为什么? 满足完备性和协调性。 原因:完备性包括两个条件:即刚体位移条件与常应变条件。首先,位移函数必须包含单元的刚体位移。结构中的单元不仅产生与该单元本身变形相应的位移,还可能因其他单元变形而通过节点位移产生单元刚体位移。为了正确反映单元的实际位移形态,位移函数必须具有反映刚体位移的能力。 其次,由于单元位移函数采用多项式,故在单元内部协调条件总能满足,要求反映在相邻单元之间。实质上来说,要求相邻单元间协调是为了保证单元交界面上应变有限。 3.1构造单元形函数有那些基本原则?试采用构造单元几何方法,构造T10单元的形函数,并对其收敛性进行讨论。 答:形函数是定义于单元内坐标的连续函数。通常单元位移函数采用多项式,其中的待定常数由节点位移参数确定,因此其个数应与单元节点自由度数相等。根据实体结构的几何方程,单元的应变是位移的一次导数。为了反映单元刚体位移和常应变即满足完备性要求,位移函数中必须包含常数项和一次项,即完全一次多项式。 3.3 何谓面积坐标?其特点是什么? 答:三角形单元中,任一点P(x,y)与其3个角点相连形成3个子三角形,其位置可以用下述称为面积坐标的三个比值来确定: L1=A1/A L2=A2/A L3=A3/A 其中A1,A2,A3分别为P23,P31,P12的面积。 各三角形面积为:Ai=1/2* =(ai+bi+ci)/2 由于A1+A2+A3=A,所以有L1+L2+L3=1,Li=(ai+bi+ci)/(2A) 特点:①T3单元的形函数Ni就是面积坐标Li ②面积坐标与三角形在整体坐标系中的位置无关,故称为局部坐标。 ③三个节点的面积坐标分别为节点1(1,0,0),节点2(0,1,0),节点3(0,0,1),形心的面积坐标(1/3,1/3 ,1/3)。④单元边界方程为Li=0 (i=1,2,3); ⑤在平行于2,3边的一条直线上,所有点都要相同的面积坐标。⑥面积坐标与直角坐标互为线性关系。 体积坐标:P点与四面体四个面围成的四个子四面体的体积与原来四面体体积的比值。即 剪切闭锁现象:当梁的高度与梁的长度之比t/l趋于零时,这种单元将出现这种现象,算得的挠度趋于零。 为克服剪切闭锁,使C0型单元适用于各种高度的梁。采用减缩积分方案与假设减应变法。 零能模式:对应于某种非刚体位移模式,减缩积分时高斯点上的应变正好等于零,此时的应变能当然也为零,这种非刚体位移模式称为零能模式。采用减缩积分时会发生零能模式。 5.1、等参单元:将整体坐标系中xy中形状中较复杂的真实单元变换成局部坐标系xy中规则的标准单元,然后在标准单元中构造形函数。由于坐标变换式与单元位移函数中用了相同的形函数N i(ξ,η),故称这种变换为等参变换,相应的单元称为等参单元。 2、等参单元的优越性:①有些工程较复杂,用直边单元离散这些结构需要大量的单元才能得到较好的近似,而曲边的等参单元可方便地离散复杂结构。②如在单元内多取些节点,单元便具有较多的位移自由度,从而就能够插值表示较复杂的单元内部位移场,这样就提高了单元本身的精度。③等参单元刚度矩阵、荷载矩阵的计算是在规则单元域内进行的因此不管被积函数多么复杂,都可以方便地采用标准化数值积分。 3、数值积分的阶次:对于N点积分,当被积函数为m次多项式且m<=2N-1时,可得精确积分值。反之,对于m次多项式的被积函数,精确积分要求的积分点数N>=(m+1)/2。 6.1、工程梁和剪切梁的基本假设?有两种梁弯曲理论①工程梁理论基本假设:平截面假设与横向纤维无挤压假设。前者认为梁横截面变形后仍为平面,且垂直于变形后的中性轴。该假设意味着横向剪切应变γxy =0,后者认为梁的横向纤维无挤压,即εy=0。②剪切梁理论基本假设:横向纤维元无挤压与另一假设认为法平面变形后仍为平面,但不再垂直于变形后的中性轴。 6.2. 剪切梁怎么考虑剪切影响:在结构单元分析中,可在工程梁单元的基础上考虑剪切变形的影响,也可通过挠度与转角各自独立插值直接构造剪切梁单元。 6.3对于杆系结构单元,为什么要在局部坐标系内建立单元刚度矩阵?为什么还要坐标变换?(1)在局部坐标系内可以更方便的建立单 元刚度矩阵。(2)在整体分析中,对所有单元都应采用同一个坐标系即整体坐标系X Y,否则围绕同一节点的不同单元对节点施加的节点力不能直接相加。因此,在进行整体分析之前,还需要进行坐标转换工作,把局部坐标系中得出的单元刚度方程转换成整体坐标系中的单元刚度方程,从而得出整体坐标系中的单元刚度矩阵。 7.1. 薄板弯曲理论基本假定:第一条:板厚方向的挤压变形可忽略不计,即εz=0,。这项假设类似于梁的横向纤维间无挤压假设。第二条:在板弯曲变形中,中面法线保持为直线且仍为弹性曲面(挠度曲面)的法线。第三条:薄板中面只发生弯曲变形,没有面内的伸缩变形,即中面水平位移(u)z=0=(v) z=0=0. 7.2. 厚板理论基本假设:板的中面法线变形后基本保持为直线,但因横向变形的缘故,该直线不在垂直于变形后的中面。因此,法线绕坐标轴的转角θx、θy不再是挠度的导数,而是独立变量;中面内的线位移和板厚方向的挤压变形也可忽略。 7.3. 薄板、厚板基本假定的不同:薄板:板弯曲变形中,中面法线保持为直线且仍为弹性曲面法线。厚板:板中面法线变形后仍基本保持为直线,但该直线不再垂直于变形后的中面。 7.4. DKT单元:离散Kirchhoff理论的基本思想是在若干离散点上满足Kirchhoff直法线假设。基于这种理论构造薄板单元时,w,θx,,θy 也各自独立插值;然后在若干离散点上引入直法线假设。这样构造的单元叫做DKT单元 8.1. 薄壳单元基本假设:薄壳理论假设:薄壳发生微笑变形时,忽略沿壳体厚度方向的挤压变形;且认为直法线假设成立,即变形后中面法线保持为直线且仍为中面的法线;壳体变形时中面不但发生弯曲,而且面内也将产生面内伸缩变形;折板假设;非耦合假设。 薄壳与薄板理论的假设的异同点:相同点:直法线假设和法向(板厚度方向)的纤维无挤压假设均成立。不同点:薄板中面只发生弯曲变形,没有面内的伸缩变形,即中面水平位移为零,而壳体变形时中面不但发生弯曲,而且也将产生面内伸缩变形。 厚壳分析的假设:变形前后的中曲面法线变形后仍基本保持为直线,但因横向剪切变形的缘故,该直线不再垂直于变形后的中曲面,此外,壳体厚度方向的挤压变形可以忽略。 与厚板理论的假设的 相同点:中面法线变形后仍基本保持为直线,但因横向剪切变形的缘故,该直线不在垂直于变形后的中面。厚度方向的挤压变形忽略不计。不同点:厚板理论的假设中,中面内的线位移可以忽略,而厚壳理论的假设中,中面内的位移不可忽略,并且厚壳的位移场可用中面位移表示。 8.2. 平板型单元:组成的折板系统去代替原来的壳体,由平面应力状态与平板弯曲应力状态加以组合而得壳体的应力单元。 分析这种单元时所提出的假设:理论假设:薄壳发生微笑变形时,忽略沿壳体厚度方向的挤压变形,且认为直法线假设成立,即变形后中面法线保持为直线且仍为中面的法线。,折板假设,非耦合假设。 应用平板型壳单元可能会出现的问题,如何解决:1.单元共面问题,解法:引入唯一边界条件可解方程Ka=P 。2.虚拟旋转刚度,解法:在特殊节点上给以任意的虚拟刚度系数。Kθzθzθzi=0,经坐标变换,整体坐标系中该节点平衡方程将满足唯一解条件。赋予Kθzθz任何值。3.新型平面膜元:在平面膜元角点上增加旋转自由度θz,使其有对应的刚度。 8.3. 面内变形与弯曲变形之间非耦合的假设是针对什么提出来的?试说明单元组装时,面内效应与弯曲效应的耦合将会出现。 答:面内变形与弯曲变形之间非耦合的假设是针对局部坐标系下的单元提出的。 9.1. 减少自由度的措施有哪些?各自基本概念如何? 答:1.恰当利用结构对称性。基本思想:利用结构的对称性,取结构一部分建立有限元模型。根据荷载对称性,分析对称面上的位移状态,以确定对称面上节点的位移边界条件。2.采用子结构技术。基本思想:在大型复杂结构的有限元分析中,可将原结构分成若干区域,每个区域作为一个子结构,这些子结构在其公共边界上互相连接起来。 2. 为什么说位移法中应力解的精度低于位移解? 答:在位移有限单元法中,沿单元边界是连续的,而位移的导数通常不连续,因此,在单元边界上应力是不连续的;基本未知量是位移,而单元应变和应力是由位移求导得到的,因此应力精度低于位移精度。 3. 在无法获得精确解的条件下,如何进行误差估计? 答:有限元解法的误差估计有:残值法,后处理法。后处理法:由于无法获得精确解,一般以修匀后的改进值σ*作为“精确解”进行误差估计,通过与精确值误差范数对比,这样做非常有效。

相关主题