固支梁有限元课程设计.docx
- 文档编号:27030427
- 上传时间:2023-06-25
- 格式:DOCX
- 页数:90
- 大小:84.29KB
固支梁有限元课程设计.docx
《固支梁有限元课程设计.docx》由会员分享,可在线阅读,更多相关《固支梁有限元课程设计.docx(90页珍藏版)》请在冰豆网上搜索。
固支梁有限元课程设计
一.题目
如图所示固支梁,高3M,长15m,承受均布载荷q=10kN/,E=20GPa,u=0.167,厚度t=1m,忽略自重,按平面应力问题分析.用有限元方法计算梁地变形及应力分布,要求用矩形单元.
要求:
1.单元数目不得少于20个.
2.采用矩形单元计算求解.
3.计算结果并给出变形图、应力分布图、单元划分图.
二、力学分析
1.题目可以看做是平面应力问题故有LXM=0
2.单元划分图
三.程序框图及程序
程序框图:
四.源程序
#include
#include
#defineNE30//单元数
#defineNJ44//节点数
#defineNZ16//支承数
#defineNPJ11//节点载荷作用数
#defineDD26//半带宽
#defineNJ288//节点位移数
intLXM=0。
doubleE0=2e7。
doubleMU=0.167。
doubleLOU=0.0。
doubleTE=1。
doubleAJZ[NJ+1][3]={{0,0,0},
{0,0,3},{0,1.5,3},{0,3,3},{0,4.5,3},{0,6,3},{0,7.5,3},{0,9,3},{0,10.5,3},{0,12,3},{0,13.5,3},{0,15,3},
{0,0,2},{0,1.5,2},{0,3,2},{0,4.5,2},{0,6,2},{0,7.5,2},{0,9,2},{0,10.5,2},{0,12,2},{0,13.5,2},{0,15,2},
{0,0,1},{0,1.5,1},{0,3,1},{0,4.5,1},{0,6,1},{0,7.5,1},{0,9,1},{0,10.5,1},{0,12,1},{0,13.5,1},{0,15,1},
{0,0,0},{0,1.5,0},{0,3,0},{0,4.5,0},{0,6,0},{0,7.5,0},{0,9,0},{0,10.5,0},{0,12,0},{0,13.5,0},{0,15,0}}。
//节点坐标
intJM[NE+1][5]={{0,0,0,0,0},{0,10,11,22,21},{0,9,10,21,20},{0,8,9,20,19},{0,7,8,19,18},{0,6,7,18,17},{0,5,6,17,16},{0,4,5,16,15},{0,3,4,15,14},{0,2,3,14,13},{0,1,2,13,12},{0,21,22,33,32},{0,20,21,32,31},{0,19,20,31,30},{0,18,19,30,29},{0,17,18,29,28},{0,16,17,28,27},{0,15,16,27,26},{0,14,15,26,25},{0,13,14,25,24},{0,12,13,24,23},{0,32,33,44,43},{0,31,32,43,42},{0,30,31,42,41},{0,29,30,41,40},{0,28,29,40,39},{0,27,28,39,38},{0,26,27,38,37},{0,25,26,37,36},{0,24,25,36,35},{0,23,24,35,34}}。
intNZC[NZ+1]={0,1,2,21,22,23,24,43,44,45,46,65,66,67,68,87,88}。
//支撑数组
doublePJ[NPJ+1][2+1]={{0,0,0},{0,-7500,68},{0,-15000,70},{0,-15000,72},{0,-15000,74},{0,-15000,76},{0,-15000,78},{0,-15000,80},{0,-15000,82},{0,-15000,84},{0,-15000,86},{0,-7500,88}}。
//节点载荷数组//节点载荷数组
doubleAE,KZ[NJ2+1][DD+1],P[NJ2+1],S[3+1][8+1],KE[8+1][8+1],SZ[3+1][5*8+1]。
doubleJDYL[NJ][6]。
//节点应力矩阵
doubleDYYL[4][NE][4]。
//单元应力矩阵
intIE,JE,ME,PE。
voidDUGD(int,int)。
FILE*fp1,*fp2,*ab。
voidmain()
{
intNJ1,k,IN,IM,jn,m,i,j,z,JO,ii,jj,h,dh,E,l,zl,dl,r,n,o,f。
doubleF,c,SIG1,SIG2,SIG3,PYL,RYL,MAYL,MIYL,CETA。
doubleWY[8+1],YL[3+1]。
ab=fopen("节点应力.txt","w")。
fp1=fopen("节点位移.txt","w")。
fp2=fopen("单元应力.txt","w")。
if(LXM!
=0)
{
E0=E0/(1.0-MU*MU)。
MU=MU/(1.0-MU)。
}
for(i=0。
i<=NJ2。
i++)
{
for(j=0。
j<=DD。
j++)
KZ[i][j]=0.0。
}
for(E=1。
E<=NE。
E++)
{
DUGD(E,3)。
for(i=1。
i<=4。
i++)
{
for(ii=1。
ii<=2。
ii++)
{
h=2*(i-1)+ii。
dh=2*(JM[E][i]-1)+ii。
for(j=1。
j<=4。
j++)
{
for(jj=1。
jj<=2。
jj++)
{
l=2*(j-1)+jj。
zl=2*(JM[E][j]-1)+jj。
dl=zl-dh+1。
if(dl>0)
KZ[dh][dl]=KZ[dh][dl]+KE[h][l]。
}
}
}
}
}
for(i=1。
i<=NJ2。
i++)
P[i]=0.0。
if(NPJ>0)
{
for(i=1。
i<=NPJ。
i++)
{
j=(int)PJ[i][2]。
P[j]=PJ[i][1]。
}
}
if(LOU>0)
{
for(E=1。
E<=NE。
E++)
{
DUGD(E,1)。
F=-LOU*(AE)*TE/4。
P[2*IE]=P[2*IE]+F。
P[2*JE]=P[2*JE]+F。
P[2*ME]=P[2*ME]+F。
P[2*PE]=P[2*PE]+F。
}
}
for(i=1。
i<=NZ。
i++)
{
z=NZC[i]。
KZ[z][1]=1.0。
for(j=2。
j<=DD。
j++)
KZ[z][j]=0.0。
if(z!
=1)
{
if(z>DD)
JO=DD。
else
JO=z。
for(j=2。
j<=JO。
j++)
KZ[z-j+1][j]=0.0。
}
P[z]=0.0。
}
NJ1=NJ2-1。
for(k=1。
k<=NJ1。
k++)
{
if(NJ2>k+DD-1)
IM=k+DD-1。
else
IM=NJ2。
IN=k+1。
for(i=IN。
i<=IM。
i++)
{
l=i-k+1。
c=KZ[k][l]/KZ[k][1]。
jn=DD-l+1。
for(j=1。
j<=jn。
j++)
{
m=j+i-k。
KZ[i][j]=KZ[i][j]-c*KZ[k][m]。
}
P[i]=P[i]-c*P[k]。
}
}
P[NJ2]=P[NJ2]/KZ[NJ2][1]。
for(i=NJ1。
i>=1。
i--)
{
if(DD>=NJ2-i+1)
JO=NJ2-i+1。
else
JO=DD。
for(j=2。
j<=JO。
j++)
{
h=j+i-1。
P[i]=P[i]-KZ[i][j]*P[h]。
}
P[i]=P[i]/KZ[i][1]。
}
printf("\n")。
printf("JDUV\n")。
fprintf(fp1,"JDUV\n")。
for(i=1。
i<=NJ。
i++)
{
printf("%d%-9.6f%-9.6f\n",i,P[2*i-1],P[2*i])。
fprintf(fp1,"%d%-9.6f%-9.6f\n",i,P[2*i-1],P[2*i])。
}
for(i=1。
i<=NJ。
i++)
{
fprintf(fp1,"a%d=%-9.6f。
b%d=%-9.6f。
\n",i,AJZ[i][1]+P[2*i-1],i,AJZ[i][2]+P[2*i])。
}
for(E=1。
E<=NE。
E++)
{
DUGD(E,2)。
for(i=1。
i<=4。
i++)
{
for(j=1。
j<=2。
j++)
{
h=2*(i-1)+j。
dh=2*(JM[E][i]-1)+j。
WY[h]=P[dh]。
}
}
for(n=1。
n<=5。
n++)
{
for(i=1。
i<=3。
i++)
{
YL[i]=0。
for(j=1。
j<=8。
j++)
YL[i]=YL[i]+SZ[i][8*(n-1)+j]*WY[j]。
}
SIG1=YL[1]。
SIG2=YL[2]。
SIG3=YL[3]。
PYL=(SIG1+SIG2)/2。
RYL=sqrt(pow((SIG1-SIG2)/2.0,2)+pow(SIG3,2))。
MAYL=PYL+RYL。
MIYL=PYL-RYL。
if(SIG2==MIYL)
CETA=0。
else
CETA=90-57.29578*atan2(SIG3,(SIG2-MIYL))。
printf("\n")。
printf("E=%d(%d)\n",E,n)。
printf("sx=%-9.6fsy=%-9.6ftou=%-9.6f\n",SIG1,SIG2,SIG3)。
printf("s1=%-9.6fs2=%-9.6ftheta=%-9.6f\n",MAYL,MIYL,CETA)。
fprintf(fp2,"\n")。
fprintf(fp2,"E=%d(%d)\n",E,n)。
fprintf(fp2,"sx=%-9.6fsy=%-9.6ftou=%-9.6f\n",SIG1,SIG2,SIG3)。
fprintf(fp2,"s1=%-9.6fs2=%-9.6ftheta=%-9.6f\n",MAYL,MIYL,CETA)。
if(n<5)
{
DYYL[n-1][E-1][0]=SIG1。
//将各单元应力记录到单元应力矩阵中
DYYL[n-1][E-1][1]=SIG2。
DYYL[n-1][E-1][2]=SIG3。
}
}
}
fprintf(ab,"各节点应力:
\n")。
for(r=0。
r r++)//计算节点应力(始) { o=0。 for(n=0。 n n++) if(JM[n+1][1]==r+1||JM[n+1][2]==r+1||JM[n+1][3]==r+1||JM[n+1][4]==r+1) { o=o+1。 for(f=0。 f<5。 f++) if(JM[n+1][f+1]==r+1) { JDYL[r][0]=JDYL[r][0]+DYYL[f][n][0]。 JDYL[r][1]=JDYL[r][1]+DYYL[f][n][1]。 JDYL[r][2]=JDYL[r][2]+DYYL[f][n][2]。 } //printf("%f\n%f\n%f\n\n\n",JDYL[r][0],JDYL[r][1],JDYL[r][2])。 //检查单元应力各项是否正确求和 } JDYL[r][0]=JDYL[r][0]*1.0/o。 JDYL[r][1]=JDYL[r][1]*1.0/o。 JDYL[r][2]=JDYL[r][2]*1.0/o。 JDYL[r][3]=(JDYL[r][0]+JDYL[r][1])*0.5。 JDYL[r][4]=(float)sqrt(pow((JDYL[r][0]-JDYL[r][1])*0.5,2)+pow(JDYL[r][2],2))。 JDYL[r][5]=JDYL[r][3]+JDYL[r][4]。 //计算节点应力(终) fprintf(ab,"%2d\t%f\n",r+1,JDYL[r][5]/1000000)。 //输出节点应力至文件 } fclose(ab)。 //关闭文件 fclose(fp1)。 fclose(fp2)。 } /**************DUGD子程序******************/ voidDUGD(intE,intASK) { doubleB[3+1][8+1],D[3+1][3+1],x,y,a,b,C。 inti,j,k,n,w。 doublexy[6][3]={{0,0,0},{0,-1,-1},{0,1,-1},{0,1,1},{0,-1,1},{0,0,0}}。 a=0.75。 b=0.5。 AE=4*a*b。 D[1][1]=E0/(1-MU*MU)。 D[1][2]=MU*E0/(1-MU*MU)。 D[1][3]=0。 D[2][1]=D[1][2]。 D[2][2]=D[1][1]。 D[2][3]=0。 D[3][1]=0。 D[3][2]=0。 D[3][3]=E0/(2*(1+MU))。 if(ASK>1) { for(w=1。 w<=5。 w++) { for(i=0。 i<=3。 i++) for(j=0。 j<=8。 j++) B[i][j]=0.0。 x=xy[w][1]*a。 y=xy[w][2]*b。 B[1][1]=(y-b)/AE。 B[1][3]=(b-y)/AE。 B[1][5]=(b+y)/AE。 B[1][7]=-(b+y)/AE。 B[2][2]=(x-a)/AE。 B[2][4]=-(a+x)/AE。 B[2][6]=(a+x)/AE。 B[2][8]=(a-x)/AE。 B[3][1]=(x-a)/AE。 B[3][2]=(y-b)/AE。 B[3][3]=-(a+x)/AE。 B[3][4]=(b-y)/AE。 B[3][5]=(a+x)/AE。 B[3][6]=(b+y)/AE。 B[3][7]=(a-x)/AE。 B[3][8]=-(b+y)/AE。 for(i=1。 i<=3。 i++) { for(j=1。 j<=8。 j++) { S[i][j]=0.0。 for(k=1。 k<=3。 k++) S[i][j]=S[i][j]+D[i][k]*B[k][j]。 n=8*(w-1)+j。 SZ[i][n]=S[i][j]。 } } } } if(ASK>2) { C=E0*TE/(1-MU*MU)。 KE[1][1]=C*(b/(3*a)+(1-MU)*a/(6*b))。 KE[1][2]=C*(1+MU)/8。 KE[1][3]=C*((-b)/(3*a)+(1-MU)*a/(12*b))。 KE[1][4]=C*((3*MU-1)/8)。 KE[1][5]=C*((-b)/(6*a)-(1-MU)*a/(12*b))。 KE[1][6]=(-C)*(1+MU)/8。 KE[1][7]=C*(b/(6*a)-(1-MU)*a/(6*b))。 KE[1][8]=C*((1-3*MU)/8)。 KE[2][1]=KE[1][2]。 KE[2][2]=C*(a/(3*b)+(1-MU)*b/(6*a))。 KE[2][3]=KE[1][8]。 KE[2][4]=C*(a/(6*b)-(1-MU)*b/(6*a))。 KE[2][5]=KE[1][6]。 KE[2][6]=C*((-a)/(6*b)-(1-MU)*b/(12*a))。 KE[2][7]=KE[1][4]。 KE[2][8]=C*((-a)/(3*b)+(1-MU)*b/(12*a))。 KE[3][1]=KE[1][3]。 KE[3][2]=KE[2][3]。 KE[3][3]=KE[1][1]。 KE[3][4]=KE[1][6]。 KE[3][5]=KE[1][7]。 KE[3][6]=KE[1][4]。 KE[3][7]=KE[1][5]。 KE[3][8]=KE[1][2]。 KE[4][1]=KE[1][4]。 KE[4][2]=KE[2][4]。 KE[4][3]=KE[3][4]。 KE[4][4]=KE[2][2]。 KE[4][5]=KE[2][3]。 KE[4][6]=KE[2][8]。 KE[4][7]=KE[1][2]。 KE[4][8]=KE[2][6]。 KE[5][1]=KE[1][5]。 KE[5][2]=KE[2][5]。 KE[5][3]=KE[3][5]。 KE[5][4]=KE[4][5]。 KE[5][5]=KE[1][1]。 KE[5][6]=KE[1][2]。 KE[5][7]=KE[1][3]。 KE[5][8]=KE[1][4]。 KE[6][1]=KE[1][6]。 KE[6][2]=KE[2][6]。 KE[6][3]=KE[3][6]。 KE[6][4]=KE[4][6]。 KE[6][5]=KE[5][6]。 KE[6][6]=KE[2][2]。 KE[6][7]=KE[1][8]。 KE[6][8]=KE[2][4]。 KE[7][1]=KE[1][7]。 KE[7][2]=KE[2][7]。 KE[7][3]=KE[3][7]。 KE[7][4]=KE[4][7]。 KE[7][5]=KE[5][7]。 KE[7][6]=KE[6][7]。 KE[7][7]=KE[1][1]。 KE[7][8]=KE[1][6]。 KE[8][1]=KE[1][8]。 KE[8][2]=KE[2][8]。 KE[8][3]=KE[3][8]。 KE[8][4]=KE[4][8]。 KE[8][5]=KE[5][8]。 KE[8][6]=KE[6][8]。 KE[8][7]=KE[7][8]。 KE[8][8]=KE[2][2]。 } } 节点位移 JDUV 10.0000000.000000 2-0.006277-0.007092 3-0.008254-0.017187 4-0.007199-0.026976 5-0.004104-0.033857 6-0.000000-0.036319 70.004104-0.033857 80.007199-0.026976 90.008254-0.017187 100.006277-0.007092 110.0000000.000000 120.0000000.000000 13-0.001628-0.006791 14-0.002451-0.017231 15-0.002187-0.027202 16-0.001263-0.034194 17-0.000000-0.036694 180.001263-0.034194 190.002187-0.027202 200.002451-0.017231 210.001628-0.006791 220.0000000.000000 230.0000000.000000 240.001524-0.007019 250.002416-0.017494 260.002168-0.027444 270.001253-0.034437 28-0.000000-0.036937 29-0.001253-0.034437 30-0.002168-0.027444 31-0.002416-0.017494 32-0.001524-0.007019 330.0000000.000000 340.0000000.000000 350.006339-0.007844 360.008222-0.017938 370.007173-0.027708 380.004094-0.034587 39-0.000000-0.037050 40-0.004094-0.034587 41-0.007173-0.027708 42-0.008222-0.017938 43-0.006339-0.007844 440.0000000.000000 a1=0.000000。 b1=3.000000。 a2=1.493723。 b2=2.992908。 a3=2.991746。 b3=2.982813。 a4=4.492801。 b4=2.973024。 a5=5.995896。 b5=2.966143。 a6=7.500000。 b6=2.963681。 a7=9.004104。 b7=2.966143。 a8=10.507199。 b8=2.973024。 a9=12.008254。 b9=2.982813。 a10=13.506277。 b10=2.992908。 a11=15.000000。 b11=3.000000。 a12=0.000000。 b12=2.000000。
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 固支梁 有限元 课程设计