main.cpp

有限元计算,可以计算一型梁单元的软件. 我们可以共同提高

#include "iostream.h"
#include "math.h"
#include "fstream.h"
#include "iomanip.h"
#include "CEle.h"
int main()
{
double l,b=1.0,h,q,E=2e8,mu=0.3,f;
//l表示简支梁的长度,bXh表示梁的宽度与高度,q表示梁的均布力
//f表示作用在梁中的集中力,E表示弹性模量,mu表示泊松比
cout<<"本程序用有限元方法计算一受均布力下的两端铰支梁"<<endl;
cout<<"并给出最大挠度及应力值"<<endl<<endl;
char flag;
cout<<"///////************按任意键开始,退出请输入N**********////////////////////"<<endl;
cin>>flag;
if(flag=='N'||flag=='N')
	return 0;

cout<<"请输入梁的基本参数"<<endl;
cout<<"长度 l=";
cin>>l;
cout<<"高度 h=";
cin>>h;
cout<<"请输入弹性模量E(2e8),及泊松比mu(0.3)"<<endl;
cout<<"E=";cin>>E;
if(E==0) 
	E=2e8;
cout<<"mu=";cin>>mu;
int I1,I2,i;
cout<<"请输入划分格数:"<<endl<<"横向格数为: ";cin>>I1;
if(I1<=0)
	throw "格数划分必须大于0";
cout<<"请输入划分格数:"<<endl<<"纵向格数为: ";cin>>I2;
if(I2<=0)
	throw "格数划分必须大于0";
cout<<endl<<"请输入均布荷载力Q=";cin>>q;
cout<<endl<<"请输入集中力F=";cin>>f;
char flag1,flag2;
cout<<endl<<"请输入第一个支座类型,铰结(j),固结(g),滑动(h)";cin>>flag1;
cout<<endl<<"请输入第二个支座类型,铰结(j),固结(g),链杆(l),自由(x),滑动(h)";cin>>flag2;
if((flag2=='j'&&flag2=='x')||(flag1=='h'&&(flag2=='h'||flag2=='x'||flag2=='l')))
	throw "结构约束不够";
//建立单元求解
if((I1/2)*2==I1)
	i=(3*I2+2)*(I1/2);
else
	i=(3*I2+2)*(I1/2)+2*I2+1;
cout<<endl<<"正在计算中,请稍候..."<<endl;
CEle stru(l,b,h,q,I1,I2,E,mu);
stru.GetD();							//////得到梁的D矩阵
stru.GetK();			     			//////得到梁的整体刚度矩阵,以及整体荷载矩阵
stru.restrict(flag1,flag2);                        //////添加约束
stru.m_mP(2*i-1,0)=f;                   //////将集中力加入整体荷载矩阵中
stru.Solve();							//////解方程
stru.GetS();							//////获得整体应力矩阵
double max1=0,max2,max3;
for(i=0;i<stru.m_nNJ;i++)
{
	if(max1<stru.m_mS(0,i))
			max1=stru.m_mS(0,i);
	if(max2<stru.m_mS(1,i))
			max2=stru.m_mS(1,i);
	if(max3<stru.m_mS(2,i))
			max3=stru.m_mS(2,i);
}
cout<<endl<<"最大应力值为:"<<max1<<"  "<<max2<<"  "<<max3<<endl;
max1=max2=0;
for(i=1;i<stru.m_nNJ;i++)
	if(fabs(max1*1e5)<fabs(stru.m_mDis(2*i-1,0)*1e5))
	{
		max1=stru.m_mDis(2*i-1,0);
		max2=i;
	}
cout<<"最大位移是:"<<max1<<"  "<<"在第"<<max2<<"结点上."<<endl;

//将结果输入到文件中
ofstream fout("D:\\result.txt",ios::in);
fout<<"////////////*********结果输出**********//////////"<<endl;
fout.setf(ios::right,ios::adjustfield);
fout.setf(ios::scientific,ios::floatfield);
fout.precision(4);
fout<<"最大位移是:"<<max1<<"\t\t"<<"在第"<<max2<<"结点上."<<endl<<endl;
fout<<"结点号\t\t"<<"X方向位移\t\t"<<"y方向位移\r"<<endl;
for(i=0;i<stru.m_nNJ;i++)
fout<<i+1<<"\t\t"<<stru.m_mDis(2*i,0)<<"\t\t"<<stru.m_mDis(2*i+1,0)<<endl;
fout<<endl<<endl;
fout<<"///////////************应力结果**************/////////////////////"<<endl;
fout<<"结点号\t\t"<<"x方向正应力\t\t"<<"y方向应力\t\t"<<"xy方向剪力"<<endl;
for(i=0;i<stru.m_nNJ;i++)
fout<<i+1<<"\t\t"<<stru.m_mS(0,i)<<"\t\t"<<stru.m_mS(1,i)<<"\t\t"<<stru.m_mS(2,i)<<endl;
fout<<"\r\r结果输出结束";
cout<<"结果在'D:\\result.txt'文件中,请查看." <<endl;
//实时查看结果
cout<<"是否要显示全部结点位移(Y/N)?"<<endl;cin>>flag;
if(flag=='y'||flag=='Y')
{
cout<<"结点 "<<"X方向位移\t\t"<<"y方向位移\t"<<endl;
for(i=0;i<stru.m_nNJ;i++)
{
	cout<<i+1<<"  "<<stru.m_mDis(2*i,0)<<"\t\t"<<stru.m_mDis(2*i+1,0)<<endl;
}
}
cout<<endl<<"要查看任一节点位移或应力?(Y/N)";cin>>flag;
if(flag=='Y'||flag=='y')
	i=1;
else
	i=0;

while(i)
{
	cout<<"请输入要查看的节点号(零表示退出)1~"<<stru.m_nNJ<<endl;
	cin>>i;
	if(i>0&&i<=stru.m_nNJ)
	{
		cout<<"第"<<i<<"节点位移是:"<<stru.m_mP(2*i-2,0)<<" "<<stru.m_mP(2*i-1,0)<<endl;
		cout<<"节点力为"<<stru.m_mS(0,i-1)<<"  "<<stru.m_mS(1,i-1)<<"  "<<stru.m_mS(2,i-1)<<endl;
	}

}
cout<<endl<<"感谢黄老师与车老师这个学期的指导"<<endl;
cout<<"谨祝黄老师与车老师身体健康,工作顺利!!"<<endl<<endl;
for(i=1;i<10000;i++)
	{}
cout<<"程序结束"<<endl;

return 0;
}