void Getnumber(double nodeCoo[8][3],double elementK[24][3])
{
double Gauss[2]={0.57735026919,-0.57735026919};//高斯点
double WeiG[2]={1.0,1.0};//权重
int local[3][8]={{1,1,-1,-1,1,1,-1-1},{-1,1,1,-1,-1,1,1,-1},{-1,-1,-1,-1,1,1,1,1}};//结点局部坐标
double diffN[3][8];//插值函数关于局部坐标偏导数
double diffNlnverJ[3][8]={0};//插值函数关于整体坐标偏导数
double DJ;
double E;
double v;
double J[3][3];
double *MatFirstadd;
double lnverJ[3][3];//存储J的逆矩阵
double B[6][24]={0};//存储B矩阵
double Blnver[24][6];//存储B的转置
double *MatrixinverFirstadd;
double BlnverD[24][6]={0};//存储B的转置与D的乘积
int i,j,ksi,eta,zeta,m,n;
double Gaussksi,Gausseta,Gausszeta;
double G=E/(2*(1+v)),lam=E*v/((1+v)*(1-2*v)); //已知常数E,v
double
D[6][6]={{lam+2*G,lam,lam,0,0,0},{lam,lam+2*G,lam,0,0,0},{lam,lam,lam+2*G,0,0,0},{0,0,0,G,0,0},{0,0,0,0,G,0},{0,0,0,0,0,G}};//D矩阵
double BTDBJ[24][24]={0};
{
for(ksi=0;ksi<2;ksi++) //以ksi方向高斯点计算积分
{
Gaussksi=Gauss[ksi];
for(eta=0;eta<2;eta++) //以eta方向高斯点计算积分
{
Gausseta=Gauss[eta];
for(zeta=0;zeta<2;zeta++); //以zeta方向高斯点计算积分
{
Gausszeta=Gauss[zeta];
for(i=0;i<3;i++) //以坐标进行循环
{
if(i==0)
{
for(j=0;j<8;j++) //以节点进行循环
{
diffN[0][j]=(1.0/8)*(1+Gausseta*local[1][j])*(1+Gausszeta*local[2][j]); //求插值函数关于kis的偏函数
}
}
else if(i==1)
{
for(j=0;j<8;j++) //以节点进行循环
{
diffN[1][j]=(1.0/8)*Gausseta*(1+Gaussksi*local[0][j])*(1+Gausszeta*local[2][j]); //求插值函数关于eta的偏导数
}
}
else
{
for(j=0;j<8;j++) //以节点进行循环
{
diffN[2][j]=(1.0/8)*Gausszeta*(1+Gaussksi*local[0][j])*(1+Gausseta*local[1][j]); //求插值函数关于zeta的偏导数
}
}
}
for(i=0;i<3;i++) //计算雅可比矩阵
{
for(j=0;j<3;j++)
{
J[i][j]=0;
for(m=0;m<8;m++)
{
J[i][j]=J[i][j]+diffN[i][m]*nodeCoo[m][j];
}
}
}
m=3;
DJ=Surplus(J[0],m,m); //计算行列式j的值
MatFirstadd=MatrixOpp(J[0],m,m);//返回矩阵就[3][3]逆矩阵的首地址
for(i=0;i<3;i++)
{
for(j=0;j<3;j++)
{
lnverJ[i] [j] =*(MatFirstadd+i*3+j);//利用返回的首地址求出逆矩阵
}
}
for(i=0;i<3;i++)
{
for(j=0;j<8;j++)
{
for(m=0;m<3;m++)
{
diffNlnverJ[i][j]=diffNlnverJ[i][j]+lnverJ[i][m]*diffN[m][j];//计算插值函数关于整体坐标xyz的偏导数
}
}
}
for(j=0;j<8;j++)//计算B矩阵
{
B[0][j*3]=diffNlnverJ [0][j];
B[1][1+j*3]=diffNlnverJ[1][j];
B[2][2+j*3]=diffNlnverJ[2][j];
B[3][3*j]=B[1][1+3*j];
B[3][1+j*3]=B[0][j*3];
B[4][1+j*3]=B[2][2+j*3];
B[4][2+j*3]=B[1][1+j*3];
B[5][j*3]=B[2][2+j*3];
B[5][2+j*3]=B[0][j*3];
}
m=6;
n=24;
MatrixlnverJFirstadd=Matrixlnver(B[0],m,n);
for(i=0;i<24;i++)
{
for(j=0;j<6;j++)
{
Blnver[i][j]=*(MatrixlnverFirstadd+i*6+j);//利用返回首地址求B矩阵的转置
}
}
for(i=0;i<24;i++ )
{
for(j=0;j<6;j++)
{
for(m=0;m<0;m++)
{
BlnverD[i][j]=BlnverD[i][j]+Blnver[i][m]*D[m][j];//计算B的转置与D的乘积
}
}
}
for(i=0;i<24;i++)
{
for(j=0;j<24;j++)
{
for(m=0;m<6;m++)
{
BTDBJ[i][j]=BTDBJ[i][j]+BlnverD[i][m]*B[m][j]*DJ;//计算B的转置,D与B,雅可比矩阵的乘积
}
}
}
for(i=0;i<24;i++)
{
for(j=0;j<24;j++)
{
elementK[i][j]=BTDBJ[i][j]+elementK[i][j];
}
}
}
}
}
}
}
有限元六面体单元 C语言源代码.rar
版权申诉
5星 · 超过95%的资源 152 浏览量
2022-04-15
12:43:33
上传
评论 2
收藏 2KB RAR 举报
阿斯弗的撒旦
- 粉丝: 150
- 资源: 43
最新资源
- unity开发教程.docx
- 代码使用Pygame库实现了一个简单的烟花模拟 核心逻辑包括烟花和粒子类的定义,处理位置、爆炸、尾迹和绘制等操作
- Matlab Simulink 电力电子仿真-Flyback(反激电路)电路分析
- tudou-android-release.apk
- 数据分析教程.docx
- 基于matlab实现用有限元法计算电磁场的Matlab工具 .rar
- 基于matlab实现有限元算法 计算电磁场问题 边界条件包括第一类边界和第二类边界.rar
- 基于matlab实现用于计算不同车重下的电动汽车动力性和经济性.rar
- 基于matlab实现遗传算法求解多车场车辆路径问题 有多组算例可以用.rar
- 浏览器.apk
资源上传下载、课程学习等过程中有任何疑问或建议,欢迎提出宝贵意见哦~我们会及时处理!
点击此处反馈