作业有限元.docx
《作业有限元.docx》由会员分享，可在线阅读，更多相关《作业有限元.docx（16页珍藏版）》请在冰豆网上搜索。
作业有限元
作业1：
求截面面积F、形心位置y0、x0，惯性矩Ix0、Iy0
一、流程图：
二、变量说明：
n—截面上的节点个数；
X[100]、Y[100]—用来存放各节点的X、Y坐标；
f—单个三角形的面积；
y0—截面形心到x轴距离；
x0—截面形心到y轴距离；
Ix0、Iy0—截面对形心的惯性矩；
F—截面的面积；
yf—单个三角形的形心y坐标；
xf—单个三角形的形心x坐标；
jx—单个三角形对x轴的惯性矩；
jy—单个三角形对y轴的惯性矩；
Sx—对x轴的静面矩；
Sy—对y轴的静面矩。
三、源程序代码：
#include
#include
voidmain（）
{
intn,i;
cout<<"请输入节点个数:
";
cin>>n;
doubleX[100],Y[100];
for（i=0;i{
cout<<"请输入第"<";
cin>>X[i];
cout<<"请输入第"<";
cin>>Y[i];
}
doublef,y0,x0,Ix0,Iy0,F=0.0,yf,xf,jx,jy,Sx=0.0,Sy=0.0;
//将坐标原点移至第一点；
for（i=1;i{
X[i]=X[i]-X[0];
Y[i]=Y[i]-Y[0];
}
for（i=1;i{
f=（X[i]*Y[i+1]-Y[i]*X[i+1]）/2.0;
yf=（Y[i]+Y[i+1]）/3.0;
xf=（X[i]+X[i+1]）/3.0;
jx=f*（Y[i]*Y[i]+Y[i]*Y[i+1]+Y[i+1]*Y[i+1]）/6.0;
jy=f*（Y[i]*X[i]+X[i]*X[i+1]+X[i+1]*X[i+1]）/6.0;
F=F+f;
Sx=Sx+f*yf;
Sy=Sy+f*xf;
}
y0=Sx/F;
x0=Sy/F;
Ix0=jx-F*y0*y0;
Iy0=jy-F*x0*x0;
y0=y0+Y[0];
x0=x0+X[0];
cout<<"**********计算结果**********/n";
cout<<"F="<cout<<"y0="<cout<<"x0="<cout<<"Ix0="<cout<<"Iy0="<cout<<"****************************"<}
四、算例：
如上图所示：
节点数N=22，节点编号如图所示。
各节点的坐标如下：
计算结果如下：
作业2：
将BASIC程序改为FORTRAN程序（或C程序）
一、流程图
二、变量说明：
L—跨度；
LB—单主梁宽；
RB—人行道宽度；
XB—行车道宽度；
DX—插值函数宽度；
NCD—车道数；
YE[350]—用来存放插值点坐标；
N—桥面板分块数；
MM—主梁数；
NN—插值点数；
WIJ[15][350]—用来存放挠度值；
ETAIJ[15][350]—用来存放影响线竖向坐标；
MCQ[15]—用来存放汽车荷载横向分布系数；
MCL[15]—用来存放履带车荷载横向分布系数；
MCG[15]—用来存放挂车荷载横向分布系数；
MCR[15]—用来存放人群荷载横向分布系数；
Y[350]—用来存放主梁Y坐标；
DE—弹性模量；
G—剪切模量；
MC[4]—记录四条车道的初始分布系数；
JX—纵向单宽抗弯惯性矩；
JY—横向单宽抗弯惯性矩；
JTXY—纵横向抗扭惯性矩之和。
三、源程序代码：
#include
#include
#include
voidmain（）
{
inti,N,MM=7,NN,NCD,k,m,j,IEND,IK1,IK2,IK,IL,J1,J;
doubleYE[350],Y[15],WIJ[15][350],ETAIJ[15][350],WB[350],MC[4],MCQ[15],MCG[15],MCL[15],MCR[15];
doubleDX,L,LB,RB,XB,JX,JY,JTXY,DE=1.0,DB,ALPHA,A,B,THETA,OMEGA,PHI,ETA,GAMMA,DELTA,DM,DN,C1,C2,C3,C4;
doublePSI,D1,D2,CC,SS,SC,CS,DP,DQ,DS,DR,YJ,EMM,CP,SP,CE,SE,AMM,BMM,CMM,DMM,CA,CB,SA,SB,YY,CAY,CBY,WC,SAY,SBY,EPS,AW,W;
constdoublePI=3.1415926,G=0.425,E=2.7182818;
cout<<"请输入N:
"<cin>>N;
cout<<"请输入DX:
"<cin>>DX;
cout<<"请输入L:
"<cin>>L;
cout<<"请输入LB:
"<cin>>LB;
cout<<"请输入RB:
"<cin>>RB;
cout<<"请输入XB:
"<cin>>XB;
cout<<"请输入JX:
"<cin>>JX;
cout<<"请输入JY:
"<cin>>JY;
cout<<"请输入JTXY:
"<cin>>JTXY;
DB=LB*N/2;
NN=int（2*DB/DX+0.1）+1;//插值系数nn
NCD=int（（XB+0.3）/3.1）;//车道数
if（NCD>4）
{
NCD=4;
}
cout<<"N="<cout<<"DX="<cout<<"L="<cout<<"LB="<cout<<"RB="<cout<<"XB="<cout<<"JX="<cout<<"JY="<cout<<"JTX+JTY="<cout<<"********Coordinatesofbeamlocation********\n";
cout<<"No."<for（i=0;i{
Y[i]=-DB+LB*（i-0.5）;
cout<}
YE[0]=-DB;
for（i=1;i{
YE[i]=YE[0]+DX*i;
}
ALPHA=G*JTXY/（2*DE*sqrt（JX*JY））;
A=sqrt（0.5*（1-ALPHA））;
B=sqrt（0.5*（1+ALPHA））;
THETA=DB*pow（（JX/JY）,0.25）/L;
OMEGA=THETA*PI/DB;
PHI=THETA*PI*A;
ETA=THETA*PI*B;
GAMMA=THETA*A;
DELTA=THETA*B;
DM=（2*A*A-2*B*B+1）*B*sinh（PHI）*cosh（PHI）-（2*A*A-2*B*B-1）*A*sin（ETA）*cos（ETA）;
DN=（2*A*A-2*B*B+1）*B*sinh（PHI）*cosh（PHI）+（2*A*A-2*B*B-1）*A*sin（ETA）*cos（ETA）;
C1=A*sin（ETA）-B*cos（ETA）;
C2=A*sin（ETA）+B*cos（ETA）;
C3=A*A*sin（ETA）+2*A*B*cos（ETA）-B*B*sin（ETA）;
C4=A*A*cos（ETA）+2*A*B*sin（ETA）-B*B*cos（ETA）;
for（k=0;k{
E=YE（k）;
PSI=PI*E/DB;
D1=GAMMA*PSI;
D2=DELTA*PSI;
CC=cosh（D1）*cos（D2）;
SS=sinh（D1）*sin（D2）;
SC=sinh（D1）*cos（D2）;
CS=cosh（D1）*sin（D2）;
DP=CC*C1-SS*C2;
DQ=SC*C1-CS*C2;
DS=SC*C3-CS*C4;
DR=CC*C3-SS*C4;
for（i=0;i{
YJ=Y[i];
WIJ[i][k]=0;
for（m=0;m{
EMM=1/（4*A*pow（m,3）*pow（OMEGA,3）*DE*JY）;
CP=cosh（PHI）;
SP=sinh（PHI）;
CE=cos（ETA）;
SE=sin（ETA）;
AMM=EMM/DM*（CP-SP）*（DP*（A/B*CP*SE-SP*CE）+DR*（2*A*CP*CE+（A*A-B*B）/B*SP*SE））;
BMM=EMM/DN*（CP-SP）*（DQ*（A/B*SP*SE-CP*CE）+DS*（2*A*SP*CE+（A*A-B*B）/B*CP*SE））;
CMM=（-1.0）*EMM/DN*B*（CP-SP）*（DQ*（A/B*CP*CE+SP*SE）+DS*（-2*A*CP*SE+（A*A-B*B）/B*SP*CE））;
DMM=（-1.0）*EMM/DM*B*（CP-SP）*（DP*（A/B*SP*CE+CP*SE）+DR*（-2*A*SP*SE+（A*A-B*B）/B*CP*CE））;
CA=cosh（m*OMEGA*A*YJ）;
CB=cos（m*OMEGA*B*YJ）;
SA=sinh（m*OMEGA*A*YJ）;
SB=sin（m*OMEGA*B*YJ）;
YY=abs（E-YJ）;
CAY=cosh（m*OMEGA*A*YY）;
CBY=cos（m*OMEGA*B*YY）;
SAY=sinh（m*OMEGA*A*YY）;
SBY=sin（m*OMEGA*B*YY）;
WC=（AMM*CA*CB+BMM*SA*CB+CMM/B*CA*SB+DMM/B*SA*SB+EMM*（CAY-SAY）*（CBY+A/B*SBY））*pow（-1,（（m+1）/2））;
WIJ[i][k]=WIJ[i][k]+WC;
EPS=abs（WC/WIJ[i][k]）;
if（EPS<1e-7）
break;
}
}
AW=0.0;
for（i=0;i{
AW=AW+WIJ[i][k];
}
for（i=0;i{
ETAIJ[i][k]=WIJ[i][k]/AW;
}
}
cout<<"********Coordinatesofloadlatrelinfluenceline********\n";
cout<<"No."<for（i=0;i{
cout<}
cout<for（j=0;j{
cout<for（i=0;i{
cout<}
cout<}
//开始计算分布系数；
for（IL=0;IL{
IEND=int（（XB-2.8）/DX+0.1）+1;
IK1=int（1.8/DX+0.1）;
IK2=int（1.3/DX+0.1）;
IK=IK1+IK2;
for（i=0;i{
MC[i]=0;
}
if（NCD<2）//单车道；
{
for（i=0;i{
J=int（（DB-XB/2+0.5）/DX+0.1）+i;
W=（ETAIJ[IL][J]+ETAIJ[IL][J+IK1]）/2;
if（W>MC[0]）
{
MC[0]=W;
}
}
}
else
{
if（NCD<3）//双车道
{
for（i=0;i{
J=int（（DB-XB/2+0.5）/DX+0.1）+i;
W=（ETAIJ[IL][J]+ETAIJ[IL][J+IK1]+ETAIJ[IL][J+IK]+ETAIJ[IL][J+IK+IK1]）/2;
if（W>MC[1]）
{
MC[1]=W;
}
}
if（MC[0]>MC[1]）
{
MC[1]=MC[0];
}
}
else
{
if（NCD<4）//三车道
{
for（i=0;i{
J=int（（DB-XB/2+0.5）/DX+0.1）+i;
W=（ETAIJ[IL][J]+ETAIJ[IL][J+IK1]+ETAIJ[IL][J+IK]+ETAIJ[IL][J+IK+IK1]+ETAIJ[IL][J+2*IK]+ETAIJ[IL][2*IK+IK1]）/2;
if（W>MC[2]）
{
MC[2]=W;
}
}
MC[2]=0.8*MC[2];//折减系数
if（MC[1]>MC[2]）
{
MC[2]=MC[1];
}
}
else//四车道
{
for（i=0;i{
J=int（（DB-XB/2+0.5）/DX+0.1）+i;
W=（ETAIJ[IL][J]+ETAIJ[IL][J+IK1]+ETAIJ[IL][J+IK]+ETAIJ[IL][J+IK+IK1]+ETAIJ[IL][J+2*IK]+ETAIJ[IL][2*IK+IK1]+ETAIJ[IL][J+3*IK]+ETAIJ[IL][3*IK+IK1]）/2;
if（W>MC[3]）
{
MC[3]=W;
}
}
MC[3]=0.7*MC[3];
if（MC[2]>MC[3]）
{
MC[3]=MC[2];
}
}
}
}
MCQ（IL）=MC（NCD）;//挂车
IEND=int（（XB-5.2）/DX+0.1）;
MCG（IL）=0;
for（i=0;i{
J=int（（DB-XB/2+1.25）/DX+0.1）+i;
IK=int（0.9/DX+0.1）;
W=（ETAIJ[IL][J+2*IK]+ETAIJ[IL][J+3*IK]）/4;
if（W>MCG[IL]）
{
MCG[IL]=W;
}
}
IEND=int（（XB-4.5）/DX+0.1）+1;
MCL[IL]=0;
for（i=0;i{
J=int（（DB-XB/2+0.65）/DX+0.1）+i;
J1=int（（DB-XB/2+3.15）/DX+0.1）+i;
IK=int（0.7/DX+0.1）+1;
W=0;
for（k=0;k{
W=W+ETAIJ[IL][J+K]+ETAIJ[IL][J1+K];
}
if（W>MCK[IL]）
{
MCL[IL]=W;
}
}
MCL[IL]=MCL[IL]/（W*2*IK）;
if（RB<0.25）
{
MCR[IL]=0;
continue;
}
J=int（（DB-XB/2-RB/W）/DX+0.1）+1;
W=0;
if（ETAIJ[IL][J]>0）
{
W=W+ETAIJ[IL][J];
}
if（ETAIJ[IL][NN-J]>0）
{
W=W+ETAIJ[IL][NN-J];
}
MCR[IL]=W*RB;
}
cout<<"********Loadlateraldistributionfactor********\n";
cout<<"NO."<for（IL=0;IL{
cout<}
cout<<"*********************************"<}
四、算例：
作业三：
利用位移互等定理计算位移影响线
一、流程图：
升级会员 