有限元分析中的二维Delaunay三角网格剖分代码实现

有限元分析中的二维Delaunay三角网格剖分代码实现

//二维平面点集的Delaunay三角剖分

#include "stdafx.h"

#include

#include

#include

#include

using namespace std;

#define point_size 600

#define pi 3.1415926

struct point

{

float x,y;

};

struct triangle

{

point* Pv[3];

float r_of_sqrt;

point o_of_tr;

};

struct d_t_node

{

triangle Tnode;

d_t_node*Pt_l[3];

int position_in_dt;

int zhuangtai_when_new;

};

point p1,p2,p3,p4;

int n;

point p[point_size];

int dt_last=0;

point p_in_dtriangle1[point_size+1];

d_t_node Dtriangle[point_size];

point p_in_dtriangle2[point_size+1];

d_t_node *queue_t[point_size];

point p_in_dtriangle3[point_size+1];

int ps_last=0;

int queue_t_last=0;

point get_spoint_cin(point*p,int n);

point get_spoint_rank(point*p,int n);

void get_spoint_squre(point*p,int );

void get_spoint_cir(point*,int);

void read_spoint(point *p,int n);

void Display();

void init_D_triangle(d_t_node* D_tr);

void Pjoin_ps(point*ps,point* p_new,int ps_last);

void pjion_D_triangle(d_t_node*Dtriangle,int &dt_last ,point p_new);

int findt_clude_p(d_t_node *Dtriangle,int dt_last,point p_new);

void sort_to_line(point line_of_tr[point_size][2],int line_of_lin[point_size],int this_fh_when_new[point_size][2], int line_of_tr_last);

float hls(float a,float b,float c,float d);

point get_o(point p1,point p2,point p3 );

float long_of_lines(point point1,point point2);

bool p_in_q(d_t_node**queue,point p_new);

void put_bug(d_t_node );

void put_this(int psize[point_size][2],int last);

void read_triangle(d_t_node Dtriangle[point_size] );

void put_out_line(point line_of_tr[point_size][2],int line_of_tr_last);

int main(int argc, char* argv[])

{

for(int ii=0;ii

{

for(int jj=0;jj<3;jj++)

{

switch(jj)

{

case 0:Dtriangle[ii].Tnode.Pv[jj]=&p_in_dtriangle1[ii+1];break;

case 1:Dtriangle[ii].Tnode.Pv[jj]=&p_in_dtriangle2[ii+1];break;

case 2:Dtriangle[ii].Tnode.Pv[jj]=&p_in_dtriangle3[ii+1];break;

}

}

};

//点集为圆上的点，输入圆的个数

/* int m;

cin>>m;

n=m*36+1;

get_spoint_cir(p,m); */

//点集为正方形网格的交点，输入网格的行数

/* int m;

cin>>m;

n=m*m;

get_spoint_squre(p,m); */

//点集为随机生成的点，输入点的数量

cout<<" 请输入点数大小n=";

cin>>n;

get_spoint_rank(p,n);

glutInit(&argc, argv);

glutInitDisplayMode(GLUT_RGB | GLUT_SINGLE);

glutInitWindowPosition(100, 100);

glutInitWindowSize(300, 300);

glutCreateWindow("第一个OpenGL程序");

glColor3f(0.0,0.0,1.0);

init_D_triangle(Dtriangle);

for(int i=0;i

{

pjion_D_triangle(Dtriangle,dt_last,p[i]);

}

glClearColor(1.0,1.0,1.0,1.0);

glutDisplayFunc(Display);

glutMainLoop();

return 0;

};

void get_spoint_cir(point*p,int nn)

{

float r;

p[0].x=0.5;

p[0].y=0.5;

//int xi,xii;

for(int iii=0;iii

{

//xi=nn-1;xii=xi-iii;

r=(float(iii)+1.0)*0.5/float(nn+1);

for(int jjj=0;jjj<36;jjj++)

{

p[iii*36+jjj+1].x=float(r*cos(float(jjj)*pi/18.0))+0.5;

p[iii*36+jjj+1].y=float(r*sin(float(jjj)*pi/18.0))+0.5;

}

}

};

void get_spoint_squre(point*p,int nn)

{

for(int i=0;i

for(int j=0;j

p[i*nn+j].x=float(j+1.0)/float(nn+1),

p[i*nn+j].y=float(i+1.0)/float(nn+1);

};

void put_out_line(point line_of_tr[point_size][2],int line_of_tr_last)

{

for(int i=0;i

{

cout<<"("< ";

cout<<"("<

}cout<

};

void read_triangle(d_t_node Dtriangle[point_size])

{

for(int i=0;i

{

glBegin(GL_LINE_LOOP);

glVertex2f((Dtriangle[i].Tnode.Pv[0]->x-0.5)*2.0,(Dtriangle[i].Tnode.Pv[0]->y-0.5)*2);

glVertex2f((Dtriangle[i].Tnode.Pv[1]->x-0.5)*2.0,(Dtriangle[i].Tnode.Pv[1]->y-0.5)*2);

glVertex2f((Dtriangle[i].Tnode.Pv[2]->x-0.5)*2.0,(Dtriangle[i].Tnode.Pv[2]->y-0.5)*2);

glEnd();

};

};

void read_spoint(point *p,int n)

{

for(int i=0;i

glVertex2f(p[i].x,p[i].y);

};

void Display()

{

glClear(GL_COLOR_BUFFER_BIT);

glPolygonMode(GL_FRONT_AND_BACK,GL_FILL);

glPointSize(4.0);

read_triangle(Dtriangle);

glFlush();

}；

void put_this(int psize[point_size][2],int last)

{

for(int i=0;i

cout<

for( i=0;i

cout<

};

void put_bug(d_t_node Dtriangle1)

{

for (int i=0;i<3;i++)

{

cout<<"point : "<x<<" , "<y<

}

cout<<"the_o: "<

cout<<"the_r "<

cout<<"poistion_in_Dt "<

cout<<"zhuangtai: "<

for (i=0;i<3;i++)

{

if(Dtriangle1.Pt_l[i]!=NULL)

{

cout<<" _ "<

->

position_in_dt;

}

else cout<<" _NULL ";

}cout<

};

bool p_in_q(d_t_node*queue_t,point p_new)

{

float p_o,r_t;

p_o=long_of_lines(p_new,queue_t->Tnode.o_of_tr);

r_t=queue_t->Tnode.r_of_sqrt;

if(p_o<=r_t) return 1;

else return 0;

};

float long_of_lines(point point1,point point2)

{

float long_r;

long_r=sqrt((point1.x-point2.x)*(point1.x-point2.x)+(point2.y-point1.y)*(point2.y-point1.y)) ;

return long_r;

};

float hls(float a,float b,float c,float d)

{

float e;

e=a*d-b*c;

return e;

};

point get_o(point p1,point p2,point p3 )

{

point o;

float a,b,c,d,e;

a=(p2.x*p2.x-p3.x*p3.x+p2.y*p2.y-p3.y*p3.y)/2.0;b=p2.y-p3.y;

c=(p3.x*p3.x-p1.x*p1.x+p3.y*p3.y-p1.y*p1.y)/2.0;d=p3.y-p1.y;

e=hls(p2.x-p3.x,p2.y-p3.y,p3.x-p1.x,p3.y-p1.y);

o.x=hls(a,b,c,d)/e;

b=(p2.x*p2.x-p3.x*p3.x+p2.y*p2.y-p3.y*p3.y)/2.0;a=p2.x-p3.x;

d=(p3.x*p3.x-p1.x*p1.x+p3.y*p3.y-p1.y*p1.y)/2.0;c=p3.x-p1.x;

o.y=hls(a,b,c,d)/e;

return o;

};

void sort_to_line(point line_of_tr[point_size][2],int line_of_lin[point_size],int this_fh_when_new[point_size][2], int line_of_tr_last)

{

point change[2];

int change1;

int change_n1;

for(int i=1;i

{

for(int j=i+1;j

{

if(line_of_tr[i-1][1].x==line_of_tr[j][0].x&&line_of_tr[i-1][1].y==line_of_tr[j][0].y ) {

change[0].x=line_of_tr[j][0].x;change[0].y=line_of_tr[j][0].y;

change[1].x=line_of_tr[j][1].x;change[1].y=line_of_tr[j][1].y;

line_of_tr[j][0].x=line_of_tr[i][0].x;line_of_tr[j][1].x=line_of_tr[i][1].x;

line_of_tr[j][0].y=line_of_tr[i][0].y;line_of_tr[j][1].y=line_of_tr[i][1].y;

line_of_tr[i][0].x=change[0].x;line_of_tr[i][0].y=change[0].y;

line_of_tr[i][1].x=change[1].x;line_of_tr[i][1].y=change[1].y;

change1=this_fh_when_new[i][0];this_fh_when_new[i][0]=this_fh_when_new[j][0];this_fh_ when_new[j][0]=change1;

change1=this_fh_when_new[i][1];this_fh_when_new[i][1]=this_fh_when_new[j][1];this_fh_ when_new[j][1]=change1;

change_n1=line_of_lin[i];

line_of_lin[i]=line_of_lin[j];

line_of_lin[j]=change_n1;

break;

}

}

}

};

int findt_clude_p(d_t_node* Dtriangle,int dt_last,point p_new)

{

for(int i=0;i

{

point o_of_t;

float r_of_t,r_pando;

float a1=Dtriangle[i].Tnode.Pv[0]->x,a2=Dtriangle[i].Tnode.Pv[0]->y;

float b1=Dtriangle[i].Tnode.Pv[1]->x,b2=Dtriangle[i].Tnode.Pv[1]->y;

float c1=Dtriangle[i].Tnode.Pv[2]->x,c2=Dtriangle[i].Tnode.Pv[2]->y;

o_of_t=get_o(*Dtriangle[i].Tnode.Pv[0],*Dtriangle[i].Tnode.Pv[1],*Dtriangle[i].Tnode.Pv[2 ]);

r_of_t=sqrt((o_of_t.x-a1)*(o_of_t.x-a1)+(o_of_t.y-a2)*(o_of_t.y-a2));

r_pando=sqrt((o_of_t.x-p_new.x)*(o_of_t.x-p_new.x)+(o_of_t.y-p_new.y)*(o_of_t.y-p_new. y));

if(r_pando<=r_of_t) return i;

}return 0;

};

void pjion_D_triangle(d_t_node*Dtriangle,int &dt_last,point p_new)

{

queue_t_last=0;

d_t_node* t_clude_p=NULL;

int here=findt_clude_p(Dtriangle,dt_last, p_new);

t_clude_p=&Dtriangle[here];

queue_t[queue_t_last]=t_clude_p,queue_t[queue_t_last]->zhuangtai_when_new=1;

queue_t_last++;

int queue_t_first=0;

point line_of_tr[point_size][2];

int line_of_tr_last=0;

int line_of_lin1[point_size];

int k;

int this_fh_when_new[point_size][2];

while(queue_t_first!=queue_t_last)

{

for(int j=0;j<3;j++)

{

if(queue_t[queue_t_first]->Pt_l[j]!=NULL

&&p_in_q(queue_t[queue_t_first]->Pt_l[j],p_new))

{

if(queue_t[queue_t_first]->Pt_l[j]->zhuangtai_when_new!=1)

{

queue_t[queue_t_last]=queue_t[queue_t_first]->Pt_l[j],

queue_t[queue_t_last]->zhuangtai_when_new=1;

queue_t_last++;

}

}

else

{

if(j<2)

{

line_of_tr[line_of_tr_last][0].x=queue_t[queue_t_first]->Tnode.Pv[j]->x;

line_of_tr[line_of_tr_last]

[0].y=queue_t[queue_t_first]->Tnode.Pv[j]->y;

line_of_tr[line_of_tr_last][1].x=queue_t[queue_t_first]->Tnode.Pv[j+1]->x;

line_of_tr[line_of_tr_last][1].y=queue_t[queue_t_first]->Tnode.Pv[j+1]->y;

if(queue_t[queue_t_first]->Pt_l[j]!=NULL)

cout<Pt_l[j]->position_in_dt<

if(queue_t[queue_t_first]->Pt_l[j]!=NULL)

{

line_of_lin1[line_of_tr_last]=queue_t[queue_t_first]->Pt_l[j]->position_in_dt;

for(int l=0;l<3;l++)

{

if(queue_t[queue_t_first]->Pt_l[j]->Pt_l[l]!=NULL&&queue_t[queue_t_first]->Pt_l[j]->Pt_l[ l]->position_in_dt

==queue_t[queue_t_first]->position_in_dt)

k=l;

}

this_fh_when_new[line_of_tr_last][0]=queue_t[queue_t_first]->Pt_l[j]->position_in_dt;

this_fh_when_new[line_of_tr_last][1]=k;

}

else

{

line_of_lin1[line_of_tr_last]=-1;

this_fh_when_new[line_of_tr_last][0]=-1;

this_fh_when_new[line_of_tr_last][1]=-1;

}

line_of_tr_last++;

}

else

{

line_of_tr[line_of_tr_last][0].x=queue_t[queue_t_first]->Tnode.Pv[j]->x;

line_of_tr[line_of_tr_last][0].y=queue_t[queue_t_first]->Tnode.Pv[j]->y;

line_of_tr[line_of_tr_last][1].x=queue_t[queue_t_first]->Tnode.Pv[0]->x;

line_of_tr[line_of_tr_last][1].y=queue_t[queue_t_first]->Tnode.Pv[0]->y;

if(queue_t[queue_t_first]->Pt_l[j]!=NULL)

{

line_of_lin1[line_of_tr_last]=queue_t[queue_t_first]->Pt_l[j]->position_in_dt;

for(int l=0;l<3;l++)

{

if(queue_t[queue_t_first]->Pt_l[j]->Pt_l[l]!=NULL&&queue_t[queue_t_first]->Pt_l[j]->Pt_l[l]->p osition_in_dt

==queue_t[queue_t_first]->position_in_dt)

k=l;

}

this_fh_when_new[line_of_tr_last][0]=queue_t[queue_t_first]->Pt_l[j]->position_in_dt;

this_fh_when_new[line_of_tr_last][1]=k;

}

else

{

line_of_lin1[line_of_tr_last]=-1;

this_fh_when_new[line_of_tr_last][0]=-1;

this_fh_when_new[line_of_tr_last][1]=-1;

}

line_of_tr_last++;

}

}

}

queue_t_first++;

}

sort_to_line(line_of_tr,line_of_lin1,this_fh_when_new,line_of_tr_last);

int last_k,first_k,i;

for( i=0;i

{

int k=queue_t[i]->position_in_dt;

Dtriangle[k].Tnode.Pv[0]->x=p_new.x;

Dtriangle[k].Tnode.Pv[0]->y=p_new.y;

Dtriangle[k].Tnode.Pv[1]->x=line_of_tr[i][0].x,Dtriangle[k].Tnode.Pv[1]->y=line_of_tr[i][0] .y;

Dtriangle[k].Tnode.Pv[2]->x=line_of_tr[i][1].x,Dtriangle[k].Tnode.Pv[2]->y=line_of_tr[i][1] .y;

Dtriangle[k].Tnode.o_of_tr=get_o(*Dtriangle[k].Tnode.Pv[0],*Dtriangle[k].Tnode.Pv[1],*Dt riangle[k].Tnode.Pv[2]);

Dtriangle[k].Tnode.r_of_sqrt=long_of_lines(Dtriangle[k].Tnode.o_of_tr,*Dtriangle[k].Tnode .Pv[0]);

Dtriangle[k].position_in_dt=k;

if(i>0)

{Dtriangle[last_k].Pt_l[2]=&Dtriangle[k];Dtriangle[k].Pt_l[0]=&Dtriangle[last_k];} if(line_of_lin1[i]!=-1)

{

Dtriangle[k].Pt_l[1]=&Dtriangle[line_of_lin1[i]];

Dtriangle[this_fh_when_new[i][0]].Pt_l[this_fh_when_new[i][1]]=&Dtriangle[k];

}

else Dtriangle[k].Pt_l[1]=NULL;

Dtriangle[k].zhuangtai_when_new=0;

if(i==0)first_k=k;

last_k=k;

}

for(i=0;i

{

Dtriangle[i+dt_last].Tnode.Pv[0]->x=p_new.x,

Dtriangle[i+dt_last].Tnode.Pv[0]->y=p_new.y;

Dtriangle[i+dt_last].Tnode.Pv[1]->x=line_of_tr[i+queue_t_last][0].x,Dtriangle[i+dt_last].Tn ode.Pv[1]->y=line_of_tr[i+queue_t_last][0].y;

Dtriangle[i+dt_last].Tnode.Pv[2]->x=line_of_tr[i+queue_t_last][1].x,Dtriangle[i+dt_last].Tn ode.Pv[2]->y=line_of_tr[i+queue_t_last][1].y;

Dtriangle[i+dt_last].Tnode.o_of_tr=get_o(*Dtriangle[i+dt_last].Tnode.Pv[0],*Dtriangle[i+dt _last].Tnode.Pv[1],*Dtriangle[i+dt_last].Tnode.Pv[2]);

Dtriangle[i+dt_last].Tnode.r_of_sqrt=long_of_lines(Dtriangle[i+dt_last].Tnode.o_of_tr,*Dtri angle[i+dt_last].Tnode.Pv[0]);

Dtriangle[i+dt_last].position_in_dt=i+dt_last;

if(line_of_lin1[i+queue_t_last]!=-1)

{

Dtriangle[i+dt_last].Pt_l[1]=&Dtriangle[line_of_lin1[i+queue_t_last]];

Dtriangle[this_fh_when_new[i+queue_t_last][0]].Pt_l[this_fh_when_new[i+queue_t_last][1] ]=&Dtriangle[i+dt_last];

}

else Dtriangle[i+dt_last].Pt_l[1]=NULL;

if(i>=0)

{Dtriangle[i+dt_last].Pt_l[0]=&Dtriangle[last_k];Dtriangle[last_k].Pt_l[2]=&Dtriangle[i+dt_last]; }

Dtriangle[i+dt_last].zhuangtai_when_new=0;

last_k=i+dt_last;

}

Dtriangle[last_k].Pt_l[2]=&Dtriangle[first_k];Dtriangle[first_k].Pt_l[0]=&Dtriangle[last_k];

dt_last=dt_last+line_of_tr_last-queue_t_last;

};

void Pjoin_ps(point*ps,point* p_new,int ps_last)

{

ps[ps_last]=*p_new;

ps_last++;

};

void init_D_triangle(d_t_node* D_tr)

{

d_t_node *q1=new d_t_node;

d_t_node *q2=new d_t_node;

D_tr[0].Tnode.Pv[0]->x=0.0; D_tr[0].Tnode.Pv[0]->y=0.0;

D_tr[0].Tnode.Pv[1]->x=1.0; D_tr[0].Tnode.Pv[1]->y=0.0;

D_tr[0].Tnode.Pv[2]->x=0.0; D_tr[0].Tnode.Pv[2]->y=1.0;

D_tr[1].Tnode.Pv[0]->x=1.0; D_tr[1].Tnode.Pv[0]->y=0.0;

D_tr[1].Tnode.Pv[1]->x=1.0; D_tr[1].Tnode.Pv[1]->y=1.0;

D_tr[1].Tnode.Pv[2]->x=0.0; D_tr[1].Tnode.Pv[2]->y=1.0;

D_tr[0].Tnode.o_of_tr=get_o(*D_tr[0].Tnode.Pv[0],*D_tr[0].Tnode.Pv[1],*D_tr[0].Tnode.P v[2]);

D_tr[1].Tnode.o_of_tr=get_o(*D_tr[1].Tnode.Pv[0],*D_tr[1].Tnode.Pv[1],*D_tr[1].Tnode.P v[2]);

D_tr[0].Tnode.r_of_sqrt=long_of_lines(D_tr[0].Tnode.o_of_tr,*D_tr[0].Tnode.Pv[0]);

D_tr[1].Tnode.r_of_sqrt=long_of_lines(D_tr[1].Tnode.o_of_tr,*D_tr[1].Tnode.Pv[0]);

D_tr[0].Pt_l[0]=NULL;D_tr[0].Pt_l[1]=&D_tr[1];D_tr[0].Pt_l[2]=NULL;

D_tr[1].Pt_l[0]=NULL;D_tr[1].Pt_l[1]=NULL;D_tr[1].Pt_l[2]=&D_tr[0];

dt_last=2;

D_tr[0].position_in_dt=0;D_tr[1].position_in_dt=1;

D_tr[0].zhuangtai_when_new=0, D_tr[1].zhuangtai_when_new=0;

};

point get_spoint_rank(point*p,int n)

{

point back;

back.x=0.0;back.y=0.0;

srand(time(0));

for(int i=0;i

{

p[i].x=rand()/float(RAND_MAX);if(back.x*back.x

p[i].y=rand()/float(RAND_MAX);if(back.y*back.y

cout<

}

return back;

};

point get_spoint_cin(point*p,int n)

{

float M=0,N=0;

point back;

for(int i=0;i

{

cin>>p[i].x;if(M*M

cin>>p[i].y;if(N*N

}

back.x=M;back.y=N;

return back;

};

有限元网格划分的基本原则

有限元网格划分的基本原则 划分网格是建立有限元模型的一个重要环节，它要求考虑的问题较多，需要的工作量较大，所划分的网格形式对计算精度和计算规模将产生直接影响。为建立正确、合理的有限元模型，这里介绍划分网格时应考虑的一些基本原则。 1 网格数量 网格数量的多少将影响计算结果的精度和计算规模的大小。一般来讲，网格数量增加，计算精度会有所提高，但同时计算规模也会增加，所以在确定网格数量时应权衡两个因数综合考虑。图1中的曲线1表示结构中的位移随网格数量收敛的一般曲线，曲线2代表计算时间随网格数量的变化。可以看出，网格较少时增加网格数量可以使计算精度明显提高，而计算时间不会有大的增加。当网格数量增加到一定程度后，再继续增加网格时精度提高甚微，而计算时间却有大幅度增加。所以应注意增加网格的经济性。实际应用时可以比较两种网格划分的计算结果，如果两次计算结果相差较大，可以继续增加网格，相反则停止计算。 图1 位移精度和计算时间随网格数量的变化 在决定网格数量时应考虑分析数据的类型。在静力分析时，如果仅仅是计算结构的变形，网格数量可以少一些。如果需要计算应力，则在精度要求相同的情况下应取相对较多的网格。同样在响应计算中，计算应力响应所取的网格数应比计算位移响应多。在计算结构固有动力特性时，若仅仅是计算少数低阶模态，可以选择较少的网格，如果计算的模态阶次较高，则应选择较多的网格。在热分析中，结构内部的温度梯度不大，不需要大量的内部单元，这时可划分较少的网格。 2 网格疏密 网格疏密是指在结构不同部位采用大小不同的网格，这是为了适应计算数据的分布特点。在计算数据变化梯度较大的部位(如应力集中处)，为了较好地反映数据变化规律，需要采用比较密集的网格。而在计算数据变化梯度较小的部位，为减小模型规模，则应划分相对稀疏的网格。这样，整个结构便表现出疏密不同的网格划分形式。图2是中心带圆孔方板的四分之一模型，其网格反映了疏密不同的划分原则。小圆孔附近存在应力集中，采用了比较密的网格。板的四周应力梯度较小，网格分得较稀。其中图b中网格疏密相差更大，它比图a中的网格少48个，但计算出的孔缘最大应力相差1%，而计算时间却减小了36%。由此可见，采用疏密不同的网格划分，既可以保持相当的计算精度，又可使网格数量减

在ANSYS平台上的复杂有限元网格划分技术

在ANSYS平台上的复杂有限元网格划分技术 1. 网格密度 有限元结构网格数量的多少将直接影响计算结果的精度和计算规模的大小。一般来说，网格数量增加，计算精度会有所提高，但同时计算规模也会增加，怎样在这两者之间找到平衡，是每一个CAE工作者都想拥有的技术。网格较少时，增加网格数量可以使计算精度明显提高，而计算时间不会有大的增加。当网格数量增加到一定程度后，再继续增加网格时精度提高很少，而计算时间却大幅度增加。所以应该注意网格数量的经济性。实际应用时，可以比较两种网格划分的计算结果，如果两次计算结果相差较大，应该继续增加网格，重新计算，直到结果误差在允许的范围之内。 在决定网格数量时还应该考虑分析类型。静力分析时，如果仅仅是计算结构的变形，网格数量可以少一点。如果需要计算应力，则在精度要求相同的情况下取相对较多的网格。同样在结构响应计算中，计算应力响应所取的网格数量应该比计算位移响应的多。在计算结构固有动力特性时，若仅仅是计算少数低阶模态，可以选取较少的网格，如果计算的阶数较高，则网格数量应该相应的增加。在热分析中，结构内部的温度梯度不大时，不需要大量的内部单元，否则，内部单元应该较多。 有限元分析原则是把结构分解成离散的单元，然后组合这些单元

解得到最终的结果。其结果的精度取决于单元的尺寸和分布，粗的网格往往其结果偏小，甚至结果会发生错误。所以必须保证单元相对足够小，考虑到模型的更多的细节，使得到的结果越接近真实结果。由于粗的网格得到的结果是非保守的，因此要认真查看结果，其中有几种方法可以帮助读者分析计算结果与真实结果之间的接近程度。 最常用的方法是用对结果判断的经验来估计网格的质量，以确定网格是否合理，如通过看云图是否与物理现象相一致，如果云图线沿单元的边界或与实际现象不一致，那么很有可能结果是不正确的。 更多的评价网格误差的方法是通过比较平均的节点结果和不平均的单元结果。如在ANSYS中，提供了两条显示结果的命令：PLNS，PLES。前者是显示平均的节点结果，后者是显示不平均的单元结果。PLNS命令是计算节点结果，它是通过对该节点周围单元结果平均后得到的，分析结果是基于单元高斯积分点值，然后外插得到每个节点，因此在给定节点周围的每个单元都由自己的单元计算得到，所以这些节点结果通常是不相同的。PLNS命令是在显示结果之前将每个节点的所有结果进行了平均，所以看到的云图是以连续的方式从一个单元过渡到另外一个单元。而PLES命令不是对节点结果平均，所以在显示云图时单元和单元之间是不连续的。这种不连续程度在网格足够密(即单元足够小)的时候会很小或不存在，而在网格较粗时很大。由于PLNS结果是一个平均值，所以它得到的结果会比PLES的结果小，他

有限元网格技巧

有限元网格技巧 1. 引言 有限元法是求解复杂工程问题的一种近似数值解法，现已广泛应用到力学、热学、电磁学等各个学科，主要分析工作环境下物体的线性和非线性静动态特性等性能。 有限元法求解问题的基本过程主要包括：分析对象的离散化?有限元求解?计算结果的处理三部分。 曾经有人做过统计：三个阶段所用的时间分别占总时间的40%~50%、5%及50%~55%。也就是说，当利用有限元分析对象时，主要时间是用于对象的离散及结果的处理。如果采用人工方法离散对象和处理计算结果，势必费力、费时且极易出错，尤其当分析模型复杂时，采用人工方法甚至很难进行，这将严重影响高级有限元分析程序的推广和使用。因此，开展自动离散对象及结果的计算机可视化显示的研究是一项重要而紧迫的任务。 可喜的是，随着计算机及计算技术的飞速发展，出现了开发对象的自动离散及有限元分析结果的计算机可视化显示的热潮，使有限元分析的“瓶颈”现象得以逐步解决，对象的离散从手工到半自动到全自动，从简单对象的单维单一网格到复杂对象的多维多种网格单元，从单材料到多种材料，从单纯的离散到自适应离散，从对象的性能校核到自动自适应动态设计/分析，这些重大发展使有限元分析摆脱了仅为性能校核工具的原始阶段，计算结果的计算机可视化显示从简单的应力、位移和温度等场的静动态显示、彩色调色显示一跃成为对受载对象可能出现缺陷(裂纹等)的位置、形状、大小及其可能波及区域的显示等，这种从抽象数据到计算机形象化显示的飞跃是现在甚至将来计算机集成设计/分析的重要组成部分。 2. 有限元分析对网格剖分的要求 有限元网格生成就是将工作环境下的物体离散成简单单元的过程，常用的简单单元包括：一维杆元及集中质量元、二维三角形、四边形元和三维四面体元、五面体元和六面体元。他们的边界形状主要有直线型、曲线型和曲面型。对于边界为曲线(面)型的单元，有限元分析要求各边或面上有若干点，这样，既可保证单元的形状，同时，又可提高求解精度、准确性及加快收敛速度。不同维数的同一物体可以剖分为由多种单元混合而成的网格。网格剖分应满足以下要求：合法性。一个单元的结点不能落入其他单元内部，在单元边界上的结点均应作为单元的结点，不可丢弃。相容性。单元必须落在待分区域内部，不能落入外部，且单元并集等于待分区域。逼近精确性。待分区域的顶点(包括特殊点)必须是单元的结点，待分区域的边界(包括特殊边及面)被单元边界所逼近。良好的单元形状。单元最佳形状是正多边形或正多面体。良好的剖分过渡性。单元之间过渡应相对平稳，否则，将影响计算结果的准确性甚至使有限元计算无法计算下去。网格剖分的自适应性。在几何尖角处、应力温度等变化大处网格应密，其他部位应较稀疏，这样可保证计算解精确可靠。 3. 现有有限元网格剖分方法 K. Ho-Le 对网格生成算法进行了系统分类，该分类方法可沿用至今，它们是拓扑分解法、结点连元法、网格模板法、映射法和几何分解法五种。目前，主要是上述方法的混合使用及现代技术的综合应用。

ANSYS有限元网格划分的基本原则

ANSYS有限元网格划分的基本原则 引言 ANSYS中有两种建立有限元模型的方法：实体建模和直接生成。使用实体建模，首先生成能描述模型的几何形状的几何模型，然后由ANSYS程序按照指定的单元大小和形状对几何体进行网格划分产生节点和单元。对于直接生成法，需要手工定义每个节点的位置和单元的连接关系。 一般来说对于规模较小的问题才适于采用直接生成法，常见的问题都需要先通过实体建模生成几何模型，然后再对其划分网格生成有限元模型。随着计算机性能的提高，分析模型的复杂性和规模都越来越大，而直接生成法也因其自身的局限性逐渐的被淘汰，所以正确的理解划分网格的目的和掌握划分网格的方法不论是对ANSYS的学习还是对二次开发都有重要的作用，尤其是当模型复杂度大，对模型的某些部分网格需要特殊处理时，这种对划分网格深度的理解作用更加明显。 2 常用高级网格划分方法 随着ANSYS功能的越来越强大和计算机性能的飞速提高，有限元分析向着大型化、复杂化的方向发展，而划分网格的观念也需要逐渐从二维模型向三维模型上上转变。这里主要描述三种常见的高级划分网格的方法，正确的理解和掌握这些划分网格的思想对于二次开发者来说非常的重要。 1）延伸网格划分 延伸网格划分是指将一个二维网格延伸生成一个三维网格；三维网格生成后去掉二维网格，延伸网格划分的步骤大体包括：先生成横截面、指定网格密度并对面进行网格划分、拖拉面网格生成体网格、指定单元属性、拖拉、完成体网格划分、释放已选的平面单元。 这里通过一个延伸网格划分的简单例子来加深对这种网格划分的理解。 图1 延伸网格划分举例 建立如图1所示的三维模型并划分网格，我们可以先建立z方向的端面，然后划分网格，通过拖拉的方法在z方向按照图中所示尺寸要求的三维模型，只需

_基于ANSYS的有限元法网格划分浅析

文章编号:1003-0794(2005)01-0038-02 基于ANSYS的有限元法网格划分浅析 杨小兰,刘极峰,陈 旋 (南京工程学院,南京210013) 摘要:为提高有限元数值的计算精度和对复杂结构力学分析的准确性,针对不同分析类型采用了不同的网格划分方法,结合实例阐述了ANSYS有限元网格划分的方法和技巧,指出了采用ANSYS有限元软件在网格划分时应注意的技术问题。 关键词:ANSYS;有限元;网格;计算精度 中图号:O241 82;TP391 7文献标识码:A 1 引言 ANSYS有限元分析程序是著名的C AE供应商美国ANSYS公司的产品,主要用于结构、热、流体和电磁四大物理场独立或耦合分析的CAE应用,功能强大,应用广泛,是一个便于学习和使用的优秀有限元分析程序。在ANSYS得到广泛应用的同时,许多技术人员对ANSYS程序的了解和认识还不够系统全面,在工作和研究中存在许多隐患和障碍,尤为突出的是有限元网格划分技术。本文结合工程实例,就如何合理地进行网格划分作一浅析。 2 网格划分对有限元法求解的影响 有限元法的基本思想是把复杂的形体拆分为若干个形状简单的单元,利用单元节点变量对单元内部变量进行插值来实现对总体结构的分析,将连续体进行离散化即称网格划分,离散而成的有限元集合将替代原来的弹性连续体,所有的计算分析都将在这个模型上进行。因此,网格划分将关系到有限元分析的规模、速度和精度以及计算的成败。实验表明:随着网格数量的增加,计算精确度逐渐提高,计算时间增加不多;但当网格数量增加到一定程度后,再继续增加网格数量,计算精确度提高甚微,而计算时间却大大增加。在进行网格划分时,应注意网格划分的有效性和合理性。 3 网格划分的有效性和合理性 (1)根据分析数据的类型选择合理的网格划分数量 在决定网格数量时应考虑分析数据的类型。在静力分析时,如果仅仅是计算结构的变形,网格数量可以少一些。如果需要计算应力,则在精度要求相同的情况下取相对较多的网格。同样在响应计算中,计算应力响应所取的网格数应比计算位移响应多。在计算结构固有动力特性时,若仅仅是计算少数低阶模态,可以选择较少的网格。如果计算的模态阶次较高,则应选择较多的网格。在热分析中,结构内部的温度梯度不大,不需要大量的内部单元,可划分较少的网格。 (2)根据分析数据的分布特点选择合理的网格疏密度 在决定网格疏密度时应考虑计算数据的分布特点,在计算固有特性时,因为固有频率和振型主要取决于结构质量分布和刚度分布,采用均匀网格可使结构刚度矩阵和质量矩阵的元素不致相差很大,可减小数值计算误差。同样,在结构温度场计算中也趋于采用均匀的网格形式。在计算数据变化梯度较大的部位时,为了更好地反映数据变化规律,需要采用比较密集的网格,而在计算数据变化梯度较小的部位,为了减小模型规模,则应划分相对稀疏的网格,这样整个结构就表现出疏密不同的网格划分形式。 以齿轮轮齿的有限元分析模型为例,由于分析的目的是求出齿轮啮合传动过程中齿根部分的弯曲应力,因此,分析计算时并不需要对整个齿轮进行计算,可根据圣文男原理将整个区域缩小到直接参与啮合的轮齿。虽然实际上参与啮合的齿数总大于1,但考虑到真正起作用的是单齿,通常只取一个轮齿作为分析对象,这样作可以大大节省计算机内存。考虑到轮齿应力在齿根过渡圆角和靠近齿面处变化较大,网格可划分得密一些。在进行疏密不同网格划分操作时可采用ANSYS提供的网格细化工具调整网格的疏密,也可采用分块建模法设置网格疏密度。 图1所示即为采用分块建模法进行网格划分。图1(a)为内燃机中重要运动零件连杆的有限元应力分析图,由于连杆结构对称于其摆动的中间平面,其厚度方向的尺寸远小于长度方向的尺寸,且载荷沿厚度方向近似均匀分布,故可按平面应力分析处 38 煤 矿 机 械 2005年第1期

delaunay三角网生长准则及算法

Delaunay 三角网是Voronoi(或称thiessen多边形,V 图)图的伴生图形 ◆Delaunay 三角网的定义: 由一系列相连的但不重叠的三角形的集合, 而且这些 三角形的外接圆不包含这个面域的其他任何点。 ◆Voronoi图的定义： Voronoi图把平面分成N 个区，每一个区包括一个点， 该点所在的区域是距离该点最近的点的集合。 ◆Delaunay三角网的特性： ◆不存在四点共圆； ◆每个三角形对应于一个Voronoi图顶点； ◆每个三角形边对应于一个Voronoi图边； ◆每个结点对应于一个Voronoi图区域； ◆Delaunay图的边界是一个凸壳； ◆三角网中三角形的最小角最大。 空外接圆准则最大最小角准则最短距离和准则 在TIN中，过每个三角形的外接圆均不包含点集的其余任何点在TIN中的两相邻三角形形成 的凸四边形中，这两三角形 中的最小内角一定大于交换 凸四边形对角线后所形成的 两三角形的最小内角 一点到基边的两端的距离 和为最小 Delaunay三角剖分的重要的准则

张角最大准则面积比准则对角线准则 一点到基边的张角为最大三角形内切圆面积与三角形 面积或三角形面积与周长平 方之比最小 两三角形组成的凸四边形 的两条对角线之比。这一 准则的比值限定值，须给 定，即当计算值超过限定 值才进行优化 Delaunay三角剖分的重要的准则 不规则三角网(TIN)的建立 ●三角网生长算法就是从一个“源”开始，逐步形成覆盖整个数据区域的三角网。 ●从生长过程角度，三角网生长算法分为收缩生长算法和扩张生长算法两类。 方法说明方法实例 收缩生长算法先形成整个数据域的数据边界（凸壳）， 并以此作为源头，逐步缩小以形成整个三 角网 分割合并算法 逐点插入算法 扩张生长算法从一个三角形开始向外层层扩展，形成覆 盖整个区域的三角网 递归生长算法

CATIA有限元高级划分网格教程

CATIA有限元高级网格划分教程 盛选禹李明志 1.1进入高级网格划分工作台 （1）打开例题中的文件Sample01.CATPart。 （2）点击主菜单中的【开始】→【分析与模拟】→【Advanced Meshing Tools】(高级网格划分工具），就进入【Advanced Meshing Tools】（高级网格划分工具）工作台，如图1－1所示。进入工作台后，生成一个新的分析文件，并且显示一个【New Analysis Case】（新分析算题）对话框，如图1－2所示。 图1－1【开始】→【分析与模拟】→【Advanced Meshing Tools】(高级网格划分工具）（3）在【New Analysis Case】（新分析算题）对话框内选择【Static Analysis】（静力分析）选项。如果以后打开该对话框的时候均希望是计算静力分析，可以把对话框内的【Keep as default starting analysis case】（在开始时保持为默认选项）勾选。这样，下次进入本工作台时，将自动选择静力分析。 （4）点击【新分析算题】对话框内的【确定】按钮，关闭对话框。 1.2定义曲面网格划分参数 本节说明如何定义一个曲面零件的网格类型和全局参数。 （1）点击【Meshing Method】(网格划分方法)工具栏内的【高级曲面划分】按钮

，如图1－3所示。需要在【Meshing Method】(网格划分方法)工具栏内点击中间按钮的下拉箭头才能够显示出【高级曲 面划分】按钮。 图1－2【New Analysis Case】（新分析算题）对话框图1－3【高级曲面划分】按钮

基于delaunay三角剖分的三维地形生成

基于delaunay三角剖分的三维地形生成 1、问题背景 地图是几个世纪以来最重要的空间信息表达的载体“近年来随着高技术 的发展特别是基于计算机平台GIS的发展,地理信息系统得到日益广泛的应用。 地形与人类的生产生活息息相关,在城市规划、路径选取、资源调查与分配、工程勘查与设计、项目选址、环境监测、灾害预测与预报、军事、游戏娱乐等领域有广泛的应用,因此人们一直关心如何真实地表达自然界的地形,以满足人们生活的需要。目前,随着计算机技术的进一步发展,计算能力的不断提高,使用计算机进行地的三维表达成为目前研究的热点,这种地形的表达方式,不但感觉直观、真实性好、而且具有二维电子地图的其它优点,例如分层显示!位置顶点查找等。二维地形生成技术是当今社会的热门技术,正在被越来越多的人所重视和研究。 2、算法描述 Lawson提出了用逐点插入法建立D-三角网的算法思想［11］。Lee和Schachter，Bowyer，Watson，Sloan，Macedonio和Pareschi，Floriani和Puppo，Tsai先后进行了发展和完善。 本次实验算法为delaunay三角剖分的逐点插入法，算法步骤如下： 1、创建一个最大的三角形包含所有离散的数据点，构成初始的三角网。 2、遍历各点（p） (1)、在三角网查找包含p的三角形t。 (2)、若p在三角形内：p与三角形t的三个顶点相连构成三个三角形。加 入三角网中。如下图：

若p 在三角形边上：找出边所对应的另一个三角形的顶点，并与当前的三角形的顶点构成四个顶点，加入三角形网中。如下图： (3)、移除三角形t 。 (4)、用LOP 算法对各个三角形进行优化处理。 3、移除外围三角形。 LOP 算法 在相邻的两个三角形( abd 和bcd) 所组成的四边形中，如果对角线交换所得的两个新三角形ABC 和ABD( 如下图) 比原来的两个三角形更优，则用新的两个三角形替代原来的两个三角形。更优的标准之一是最小角度最大原则: 调整前的二个三角形共六个内角中的最小角和调整后的六个角中的最小角相比较，若前者小于后者则调整，否则不调整; 标准之二是空外接圆性质: 在由点集V 所形成。D-三角网中，其每个三角形的外接圆均不包含点集V 中的其他任意点。结合本文定义的数据结构，本文采取了以相邻三角形作为优化着眼点的处理算法。根据Delaunay 三角网空外接圆性质有以下判断: 当sin( ∠C + ∠D) ≤0， 不进行优化，p

基于四面体控制网格的模型变形算法 (1)

第20卷第9期2008年9月 计算机辅助设计与图形学学报 JO U RN A L O F COM PU T ER AID ED D ESIG N &COM P U T ER G RA PH ICS Vo l.20,N o.9 Sep.,2008 收稿日期:2008-07-15.基金项目:国家 九七三 重点基础研究发展规划项目(2002CB312101,2006CB303102);国家自然科学基金(60603078);新世纪优秀人才项目(NCET 06 0516).赵 勇,男,1982年生,博士研究生,主要研究方向为数字几何处理.刘新国,男,1972年生,博士,教授,博士生导师,主要研究方向为数字几何处理、真实感绘制、虚拟现实等.彭群生,男,1947年生,博士,教授,博士生导师,CC F 高级会员,主要研究方向为真实感图形、虚拟现实、科学计算可视化等. 基于四面体控制网格的模型变形算法 赵 勇 刘新国 彭群生 (浙江大学CAD &CG 国家重点实验室 杭州 310058)(z haoyong@cad.z https://www.doczj.com/doc/896422876.html,) 摘要 提出一种鲁棒的保体积保表面细节的模型变形算法.首先将输入模型嵌入到一个稀疏的四面体控制网格 中,并且通过一种改进的重心坐标来建立两者的对应关系;然后通过用户的交互,对控制网格建立一个二次非线性能量函数对其进行变形,而输入模型的变形结果则可以通过插值来直接获得.由于能量函数的优化是在控制网格上进行的,从而大大提高了算法的效率.与此同时,提出一种新的能量!!!Laplacian 能量,可以使四面体控制网格进行尽量刚性的变形,从而有效地防止了大尺度编辑过程中模型形状的退化现象.文中算法还具有通用性,可支持多种模型的表示方式,如三角网格模型、点模型等.实验结果表明,该算法可以有效地保持输入模型的几何细节、防止明显的体积变化,得到了令人满意的结果. 关键词 模型编辑;四面体控制网格;刚性变形;L aplacian 能量;通用性中图法分类号 T P391 Shape Deformation Based on Tetrahedral Control Mesh Zhao Yong Liu Xing uo Peng Qunsheng (S tate K ey L abor atory of CA D &CG ,Zh ej iang Univ ersity ,H ang z hou 310058) Abstract A robust shape deformation algo rithm w ith the feature o f both vo lum e and surface detail preserv ing is presented.Fir st,the input m odel is embedded into a coarse tetr ahedral co ntro l mesh,and the m odified bar ycentr ic coordinates are employ ed to establish their relationship.Then acco rding to user s editing,the contro l mesh is defor med by solving a quadric no nlinear ener gy m inimization pro blem,and the deform ation is passed to the embedded m odel by interpolatio n.As the optimization pro cess is applied to the control mesh composed of sparse vertices,the efficiency is g reatly improved.Meantime,w e incor porate a new energ y,called Laplacian energ y,into the energy equatio n to m ake the tetrahedral contro l m esh deform as rigidly as possible,thus avoiding shape degenerations even under ex treme editing.Our algor ithm acco mmodates various shape repr esentations,such as triangular meshes,point clouds etc.Experiments demonstrate that the Laplacian energy is very effective in preserv ing geom etric details and pr eventing unreasonable volume changes. Key words shape editing;tetrahedral contr ol m esh;r ig id defor matio n;Laplacian energ y;generality 近年来,随着三维数据采集技术的不断发展,三维数字几何模型已经在数字娱乐、工业设计、医学辅 助诊断、文物保护等很多领域得到了广泛的应用.数字几何处理作为计算机图形学的一个重要分支也得

有限元网格划分

有限元网格划分 摘要:总结近十年有限元网格划分技术发展状况。首先,研究和分析有限元网格划分的基本原则；其次,对当前典型网格划分方法进行科学地分类,结合实例,系统地分析各种网格划分方法的机理、特点及其适用范围,如映射法、基于栅格法、节点连元法、拓扑分解法、几何分解法和扫描法等；再次,阐述当前网格划分的研究热点,综述六面体网格和曲面网格划分技术；最后,展望有限元网格划分的发展趋势。 关键词:有限元网格划分；映射法；节点连元法;拓扑分解法；几何分解法；扫描法；六面体网格 1 引言 有限元网格划分是进行有限元数值模拟分析至关重要的一步，它直接影响着后续数值计算分析结果的精确性。网格划分涉及单元的形状及其拓扑类型、单元类型、网格生成器的选择、网格的密度、单元的编号以及几何体素。在有限元数值求解中，单元的等效节点力、刚度矩阵、质量矩阵等均用数值积分生成，连续体单元以及壳、板、梁单元的面内均采用高斯（Gauss）积分，而壳、板、梁单元的厚度方向采用辛普生（Simpson）积分。 2 有限元网格划分的基本原则 有限元方法的基本思想是将结构离散化,即对连续体进行离散化,利用简化几何单元来近似逼近连续体,然后根据变形协调条件综合求解。所以有限元网格的划分一方面要考虑对各物体几何形状的准确描述,另一方面也要考虑变形梯度的准确描述。为正确、合理地建立有限元模型,这里介绍划分网格时应考虑的一些基本原则。 2.1 网格数量 网格数量直接影响计算精度和计算时耗,网格数量增加会提高计

算精度,但同时计算时耗也会增加。当网格数量较少时增加网格，计算精度可明显提高,但计算时耗不会有明显增加;当网格数量增加到一定程度后,再继续增加网格时精度提高就很小,而计算时耗却大幅度增加。所以在确定网格数量时应权衡这两个因素综合考虑。 2.2 网格密度 为了适应应力等计算数据的分布特点,在结构不同部位需要采用大小不同的网格。在孔的附近有集中应力,因此网格需要加密;周边应力梯度相对较小,网格划分较稀。由此反映了疏密不同的网格划分原则:在计算数据变化梯度较大的部位,为了较好地反映数据变化规律,需要采用比较密集的网格;而在计算数据变化梯度较小的部位,为减小模型规模,网格则应相对稀疏。 2.3 单元阶次 单元阶次与有限元的计算精度有着密切的关联,单元一般具有线性、二次和三次等形式,其中二次和三次形式的单元称为高阶单元。高阶单元的曲线或曲面边界能够更好地逼近结构的曲线和曲面边界,且高次插值函数可更高精度地逼近复杂场函数,所以增加单元阶次可提高计算精度。但增加单元阶次的同时网格的节点数也会随之增加,在网格数量相同的情况下由高阶单元组成的模型规模相对较大,因此在使用时应权衡考虑计算精度和时耗。 2.4 单元形状 网格单元形状的好坏对计算精度有着很大的影响,单元形状太差的网格甚至会中止计算。单元形状评价一般有以下几个指标: (1)单元的边长比、面积比或体积比以正三角形、正四面体、正六面体为参考基准。 (2)扭曲度:单元面内的扭转和面外的翘曲程度。 (3)节点编号:节点编号对于求解过程中总刚矩阵的带宽和波前因数有较大的影响,从而影响计算时耗和存储容量的大小 2.5 单元协调性 单元协调是指单元上的力和力矩能够通过节点传递给相邻单元。为保证单元协调,必须满足的条件是: (1)一个单元的节点必须同时也是相邻点,而不应是内点或边界

有限元分析中的二维Delaunay三角网格剖分代码实现

有限元分析中的二维Delaunay三角网格剖分代码实现 //二维平面点集的Delaunay三角剖分 #include "stdafx.h" #include #include #include #include using namespace std; #define point_size 600 #define pi 3.1415926 struct point { float x,y; }; struct triangle { point* Pv[3]; float r_of_sqrt; point o_of_tr; }; struct d_t_node { triangle Tnode; d_t_node*Pt_l[3]; int position_in_dt; int zhuangtai_when_new; }; point p1,p2,p3,p4; int n; point p[point_size]; int dt_last=0; point p_in_dtriangle1[point_size+1]; d_t_node Dtriangle[point_size]; point p_in_dtriangle2[point_size+1]; d_t_node *queue_t[point_size]; point p_in_dtriangle3[point_size+1]; int ps_last=0; int queue_t_last=0; point get_spoint_cin(point*p,int n); point get_spoint_rank(point*p,int n);

有限元网格剖分方法概述

有限元网格剖分方法概述 在采用有限元法进行结构分析时，首先必须对结构进行离散，形成有限元网格，并给出与此网格相应的各种信息，如单元信息、节点坐标、材料信息、约束信息和荷载信息等等，是一项十分复杂、艰巨的工作。如果采用人工方法离散对象和处理计算结果，势必费力、费时且极易出错，尤其当分析模型复杂时，采用人工方法甚至很难进行，这将严重影响高级有限元分析程序的推广和使用。因此，开展自动离散对象及结果的计算机可视化显示的研究是一项重要而紧迫的任务。 有限元网格生成技术发展到现在, 已经出现了大量的不同实现方法，列举如下： 映射法 映射法是一种半自动网格生成方法，根据映射函数的不同，主要可分为超限映射和等参映射。因前一种映射在几何逼近精度上比后一种高，故被广泛采用。映射法的基本思想是：在简单区域内采用某种映射函数构造简单区域的边界点和内点，并按某种规则连接结点构成网格单元。也就是根据形体边界的参数方程，利用映射函数，把参数空间内单元正方形或单元三角形(对于三维问题是单元立方体或单元四面体)的网格映射到欧氏空间，从而生成实际的网格。这种方法的主要步骤是，首先人为地把分析域分成一个个简单可映射的子域，每个子域为三角形或四边形，然后根据网格密度的需要，定义每个子域边界上的节点数，再根据这些信息，利用映射函数划分网格。 这种网格控制机理有以下几个缺点： （1）它不是完全面向几何特征的，很难完成自动化，尤其是对于3D区域。 （2）它是通过低维点来生成高维单元。例如，在2D问题中，先定义映射边界上的点数，然后形成平面单元。这对于单元的定位，尤其是对于远离映射边界的单元的定位，是十分困难的，使得对局部的控制能力下降。 （3）各映射块之间的网格密度相互影响程度很大。也就是说，改变某一映射块的网格密度，其它各映射块的网格都要做相应的调整。 其优点是：由于概念明确，方法简单，单元性能较好，对规则均一的区域，适用性很强，因此得到了较大的发展，并在一些商用软件如ANSYS等得到应用。 2 。拓扑分解法 拓扑分解法较其它方法发展较晚, 它首先是由Wordenwaber提出来的。该方法假设最后网格顶点全部由目标边界顶点组成, 那么可以用一种三角化算法将目标用尽量少的三角形完全分割覆盖。这些三角形主要是由目标的拓扑结构决定, 这样目标的复杂拓扑结构被分解成简单的三角形拓扑结构。该方法生成的网格一般相当粗糙, 必须与其它方法相结合, 通过网格加密等过程, 才能生成合适的网格。该方法后来被发展为普遍使用的目标初始三角化算法, 用来实现从实体表述到初始三角化表述的自动化转换。 单一的拓扑分解法因只依赖于几何体的拓扑结构使网格剖分不理想，有时甚至很差。 3．连接节点法 这类方法一般包括二步：区域内布点及其三角化。早期的方法通常是先在区域内布点, 然后再将它们联成三角形或四面体, 在三角化过程中, 对所生成的单元形状难于控制。随着Delaunay三角化(简称为DT ) 方法的出现, 该类方法已成为目前三大最流行的全自动网格生成方法之一。 DT法的基本原理：任意给定N个平面点Pi(i=1,2,…,N)构成的点集为S，称满足下列条件的点集Vi为Voronoi多边形。其中，Vi满足下列条件： Vi ={ X：|X- Pi|(|X- Pj|,X(R2,i(j,j=1,2,…,N }Vi为凸多边形，称{ Vi}mi=1为Dirichlet Tesselation

基于MATLAB实现二维delaunay三角剖分

34 基于MATLAB 实现二维delaunay 三角剖分 刘锋涛 凡友华 （哈尔滨工业大学深圳研究生院 深圳 518055） 【摘要】在已知凸多边形的顶点坐标的前提情况下，利用MATLAB 中的meshgrid 函数产生多边 形附近矩形区域内的网格点的坐标，然后再利用inpolygon 函数判断哪些点位于多边形内和哪 些点位于多边形的边界上。在此基础上利用delaunay 函数来完成delaunay 三角剖分。 【关键词】delaunay 三角剖分；MATLAB 网格划分是有限元分析前处理中的关键步骤，网格划分的密度以及质量对有限元计算的精度、效率以及收敛性有着重要的影响作用。自20世纪70年代开始，关于有限元网格生成方法的研究已经取得了许多重要成果，提出许多有效的算法。Ho-Le 对网格生成方法进行了系统的分类[1]。许多学者也对网格生成的方法进行了综述，如我国的学者胡恩球等[2]、关振群等[3]。 Delaunay 三角剖分(简称DT)是目前最流行的通用的全自动网格生成方法之一。DT 有两个重要特性：最大-最小角特性和空外接圆特性。DT 的最大-最小角特性使它在二维情况下自动地避免了生成小内角的长薄单元。因此特别适用于有限元网格生成。大体上可将DT 算法分为三大类：分治算法，逐点插入法和三角网生长法。经典DT 技术已经相当成熟，近年来的研究重点是约束DT 的边界恢复算法，以及如何克服算法退化现象所产生的薄元(sliver element)问题[3]。 然而，实现DT 有限元网格生成，对于非计算机图形学专业的工程师来说还是很复杂的。在处理一些对有限元网格划分质量不过的问题时，如极限分析的有限元方法，可以采用一些更为简单的方法来实现。在Matlab 计算软件中，已有一个成熟的函数delaunay 可以实现对一系列点的DT 划分。因此，本文基于Matlab 的delaunay 等一些函数来完成一个凸多边形的DT 网格划分。 1 MATLAB 中的函数 1.1 delaunay 函数 delaunay 函数可以按照DT 网格划分的要求将一个点集中的点划归某一个有限网格所有。它在Matlab 中的用法如下： =delaunay(,) or, =delaunay(,,)TRI x y TRI x y options 其输入为点集中所有点的横、纵坐标向量x 和y ，返回值为一个3m ×的矩阵，矩阵中每一行表示DT 网格中一个三角形网格的三个顶点。 1.2 meshgird 函数 为了在任意凸多边形内产生一个点集，可以利用Matlab 中的meshgrid 命令。其用法如下： [,] = meshgrid(,)X Y x y

ANSYS 网格划分方法总结

(1) 网格划分定义：实体模型是无法直接用来进行有限元计算得，故需对它进行网格划分以生成有限元模型。有限元模型是实际结构和物质的数学表示方法。 在ANSYS中，可以用单元来对实体模型进行划分，以产生有限元模型，这个过程称作实体模型的网格化。本质上对实体模型进行网格划分也就是用一个个单元将实体模型划分成众多子区域。这些子区域（单元），是有属性的，也就是前面设置的单元属性。 另外也可以直接利用单元和节点生成有限元模型。 实体模型进行网格划分就是用一个个单元将实体模型划分成众多子区域（单元）。 （2）为什么我选用plane55这个四边形单元后，仍可以把实体模型划分成三角 形区域集合？？？ 答案：ansys为面模型的划分只提供三角形单元和四边形单元，为体单元只提供四面体单元和六面体单元。不管你选择的单元是多少个节点，只要是2D单元，肯定构成一个四边形或者是三角形，绝对没有五、六边形等特殊形状。网格划分也就是用所选单元将实体模型划分成众多三角形单元和四边形子区域。 见下面的plane77/78/55都是节点数目大于4的，但都是通过各种插值或者是合并的方式形成一个四边形或者三角形。 所以不管你选择什么单元，只要是对面的划分，meshtool上的划分类型设置就只有tri和quad两种选择。 如果这个单元只构成三角形，例如plane35，则无论你在meshtool上划分设置时tri还是quad，划分出的结果都是三角形。

所以在选用plane55单元，而划分的是采用tri划分时，就会把两个点合并为一个点。如上图的plane55,下面是plane单元的节点组成，可见每一个单元上都有两个节点标号相同，表明两个节点是重合的。 。 同样在采用plane77 单元，进行tri划分时，会有三个节点重合。这里不再一一列出。（3）如何使用在线帮助： 点击对话框中的help，例如你想了解plane35的相关属性，你可以

基于MATLAB 实现二维delaunay 三角剖分

基于MATLAB 实现二维delaunay 三角剖分 刘锋涛凡友华 （哈尔滨工业大学深圳研究生院深圳518055） 【摘要】在已知凸多边形的顶点坐标的前提情况下，利用MATLAB 中的meshgrid 函数产生多边形附近矩形区域内的网格点的坐标，然后再利用inpolygon 函数判断哪些点位于多边形内和哪些点位于多边形的边界上。在此基础上利用delaunay 函数来完成delaunay 三角剖分。 【关键词】delaunay 三角剖分；MATLAB 网格划分是有限元分析前处理中的关键步骤，网格划分的密度以及质量对有限元计算的精度、效率以及收敛性有着重要的影响作用。自20世纪70年代开始，关于有限元网格生成方法的研究已经取得了许多重要成果，提出许多有效的算法。Ho-Le 对网格生成方法进行了系统的分类[1]。许多学者也对网格生成的方法进行了综述，如我国的学者胡恩球等[2]、关振群等[3]。 Delaunay 三角剖分(简称DT)是目前最流行的通用的全自动网格生成方法之一。DT 有两个重要特性：最大-最小角特性和空外接圆特性。DT 的最大-最小角特性使它在二维情况下自动地避免了生成小内角的长薄单元。因此特别适用于有限元网格生成。大体上可将DT 算法分为三大类：分治算法，逐点插入法和三角网生长法。经典DT 技术已经相当成熟，近年来的研究重点是约束DT 的边界恢复算法，以及如何克服算法退化现象所产生的薄元(sliver element)问题[3]。 然而，实现DT 有限元网格生成，对于非计算机图形学专业的工程师来说还是很复杂的。在处理一些对有限元网格划分质量不过的问题时，如极限分析的有限元方法，可以采用一些更为简单的方法来实现。在Matlab 计算软件中，已有一个成熟的函数delaunay 可以实现对一系列点的DT 划分。因此，本文基于Matlab 的delaunay 等一些函数来完成一个凸多边形的DT 网格划分。 1MATLAB 中的函数 1.1delaunay 函数 delaunay 函数可以按照DT 网格划分的要求将一个点集中的点划归某一个有限网格所有。它在Matlab 中的用法如下： =delaunay(,) or, =delaunay(,,) TRI x y TRI x y options 其输入为点集中所有点的横、纵坐标向量x 和y ，返回值为一个的矩阵，矩阵中每一3m ×行表示DT 网格中一个三角形网格的三个顶点。 1.2meshgird 函数 为了在任意凸多边形内产生一个点集，可以利用Matlab 中的meshgrid 命令。其用法如下： [,] = meshgrid(,) X Y x y

ANSYS结构有限元分析中的网格划分技术及其应用实例

一、前言 有限元网格划分是进行有限元数值模拟分析至关重要的一步，它直接影响着后续数值计算分析结果的精确性。网格划分涉及单元的形状及其拓扑类型、单元类型、网格生成器的选择、网格的密度、单元的编号以及几何体素。从几何表达上讲，梁和杆是相同的，从物理和数值求解上讲则是有区别的。同理，平面应力和平面应变情况设计的单元求解方程也不相同。在有限元数值求解中，单元的等效节点力、刚度矩阵、质量矩阵等均用数值积分生成，连续体单元以及壳、板、梁单元的面内均采用高斯（Gauss）积分，而壳、板、梁单元的厚度方向采用辛普生（Simpson）积分。辛普生积分点的间隔是一定的，沿厚度分成奇数积分点。由于不同单元的刚度矩阵不同，采用数值积分的求解方式不同，因此实际应用中，一定要采用合理的单元来模拟求解。 CAD软件中流行的实体建模包括基于特征的参数化建模和空间自由曲面混合造型两种 方法。Pro/E和SoildWorks是特征参数化造型的代表，而CATIA与Unigraphics等则将特征参数化和空间自由曲面混合造型有机的结合起来。现有CAD软件对表面形态的表示法已经大大超过了CAE软件，因此，在将CAD实体模型导入CAE软件的过程中，必须将CAD 模型中其他表示法的表面形态转换到CAE软件的表示法上，转换精度的高低取决于接口程序的好坏。在转换过程中，程序需要解决好几何图形（曲线与曲面的空间位置）和拓扑关系（各图形数据的逻辑关系）两个关键问题。其中几何图形的传递相对容易实现，而图形间的拓扑关系容易出现传递失败的情况。数据传递面临的一个重大挑战是，将导入CAE程序的CAD模型改造成适合有限元分析的网格模型。在很多情况下，导入CAE程序的模型可能包含许多设计细节，如细小的孔、狭窄的槽，甚至是建模过程中形成的小曲面等。这些细节往往不是基于结构的考虑，保留这些细节，单元数量势必增加，甚至会掩盖问题的主要矛盾，对分析结果造成负面影响。 CAD模型的“完整性”问题是困扰网格剖分的障碍之一。对于同一接口程序，数据传递的品质取决于CAD模型的精度。部分CAD模型对制造检测来说具备足够的精度，但对有限元网格剖分来说却不能满足要求。值得庆幸的是，这种问题通常可通过CAD软件的“完整性检查”来修正。改造模型可取的办法是回到CAD系统中按照分析的要求修改模型。一方面检查模型的完整性，另一方面剔除对分析无用的细节特征。但在很多情况下，这种“回归”很难实现，模型的改造只有依靠CAE软件自身。CAE中最直接的办法是依靠软件具有的“重构”功能，即剔除细部特征、缝补面和将小面“融入”大曲面等。有些专用接口在模型传递过程中甚至允许自动完成这种工作，并且通过网格剖分器检验模型的“完整性”，如发现“完整性”不能满足要求，接口程序可自动进行“完整性”修复。当几何模型距CAE分析的要求相差太大时，还可利用CAE程序的造型功能修正几何模型。“布尔运算”是切除细节和修理非完整特征的有效工具之一。 目前数据传递一般可通过专用数据接口，CAE程序可与CAD程序“交流”后生成与CAE 程序兼容的数据格式。另一种方式是通过标准图形格式如IGES、SAT和ParaSolid传递。现有的CAD平台与通用有限元平台一般通过IGES、STL、Step、Parasolid等格式来数据