mycode.txt

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

	//数据输入
	double L=2,H=4,b=0.2;	//l表示柱的高度,H表示梁的长度
	double h1=0.4,h2=0.4;	//h1表示柱的宽度,h2表示梁的高度
	double q1=-10,q2=-10;
	double E=200e6;
	double mu=0.3;
	int nz1=2,n2=1;
	int nl1=2;
	long i,j,k;
	CEle stru1(L,b,h1,q1,nz1,n2,E,mu);
	CEle stru2(H,b,h2,q2,nl1,n2,E,mu);
	CEle stru3(L,b,h1,q1,nz1,n2,E,mu);
	long totaln=stru1.m_nNJ+stru2.m_nNJ+stru3.m_nNJ+8*n2-4;
	CMatrix totalnode(totaln,2,0);
	//将局部结点坐标赋值到整体结点坐标
    //第一根柱子结点 u
	for(i=0;i<stru1.m_nNJ;i++)
		totalnode(i,0)=stru1.m_mCOORD(i,1),totalnode(i,1)=stru1.m_mCOORD(i,0);
	//梁结点
	j=stru1.m_nNJ+4*n2-2;
	for(i=2*n2;i<stru2.m_nNJ;i++)
		totalnode(i+j,0)=stru2.m_mCOORD(i,0)+h1,totalnode(j+i,1)=L+h2-stru2.m_mCOORD(i,1);
	//第二根柱子结点
	j=totaln-stru3.m_nNJ;
	for(i=2*n2+1;i<stru3.m_nNJ;i++)
		totalnode(i+j,0)=H+h1+h1-stru3.m_mCOORD(i,1),totalnode(i+j,1)=L-stru3.m_mCOORD(i,0);	
	//第一根柱子与梁相交处
	double e=atan(h2/2/h1);
	j=stru1.m_nNJ;k=0;
	for(i=0;i<n2;i++)
		totalnode(j+i,0)=i*h1/n2,totalnode(j+i,1)=L+(1-(double)i/n2)*h1*tan(e);
	e=atan(h2/h1),j=stru1.m_nNJ+n2;
	for(i=0;i<2*n2-1;i++)
		totalnode(j+i,0)=i*h1/2/n2,totalnode(j+i,1)=L+(1-(double)i/2/n2)*h1*tan(e);
	e=atan(h1/2/h2),j=stru1.m_nNJ+3*n2-1;
	for(i=0;i<n2;i++)
		totalnode(j+i,1)=L+h2-i*h2/n2,totalnode(j+i,0)=h1-(1-(double)i/n2)*h2*tan(e);
	e=atan(h1/h2),j=j+n2;
	for(i=0;i<2*n2;i++)
		totalnode(j+i,0)=h1,totalnode(j+i,1)=L+h2-h2*i/2/n2;
	//第二根柱子与梁相交处
	j=stru1.m_nNJ+stru2.m_nNJ+4*n2-2;e=atan(h1/2/h2);
	for(i=0;i<n2;i++)
		totalnode(j+i,1)=L+h2-i*h2/n2,totalnode(j+i,0)=H+h1+(1-(double)i/n2)*h2*tan(e);
	e=atan(h1/h2),j=j+n2;
	for(i=0;i<2*n2-1;i++)
		totalnode(j+i,1)=L+h2-i*h2/2/n2,totalnode(j+i,0)=H+h1+(1-(double)i/2/n2)*h2*tan(e);
	e=atan(h1/2/h2),j=j+2*n2-1;
	for(i=0;i<n2;i++)
		totalnode(j+i,0)=H+h1+h1-i*h2/n2,totalnode(j+i,1)=L+(h2-i*h2/n2)*tan(e);
	e=atan(h2/h1/2),j=j+n2;
	for(i=0;i<2*n2;i++)
		totalnode(j+i,1)=L,totalnode(j+i,0)=H+h1+h1-h1*i/2/n2;
	///////////////////////////////////////////////////////////////////////////////
	//for(i=0;i<totaln;i++)
		cout<<i+1<<"  "<<totalnode(i,0)<<" "<<totalnode(i,1)<<endl;
	//////////////////
	//建立单元信息,即确立单元的结点编号
	long totale=(nz1+nz1+nl1+4)*n2-2,kk=0;
	CMatrix totalele(totale,8,0);
	//第一柱单元
		for(i=0;i<nz1*n2;i++)
			for(j=0;j<=7;j++)
				totalele(i,j)=stru1.m_mMP(i,j);
	//梁单元
	k=(nz1+2)*n2-1,kk=stru1.m_nNJ+4*n2-2;
		for(i=0;i<nl1*n2;i++)
			for(j=0;j<=7;j++)
				totalele(i+k,j)=stru2.m_mMP(i,j)+kk;
	//第三柱单元
	k=totale-nz1*n2,kk=totaln-stru3.m_nNJ;
		for(i=0;i<nz1*n2;i++)
			for(j=0;j<=7;j++)
				totalele(i+k,j)=stru3.m_mMP(i,j)+kk;
	//第一柱与梁结点处
		for(i=nz1*n2;i<(nz1+2)*n2-1;i++)
		{
			j=(3*n2+2)*(nz1-1)+2*(i-nz1*n2)+1;
			totalele(i,0)=j+2;
			totalele(i,3)=j;
			totalele(i,7)=j+1;
			totalele(i,6)=j+2*n2+(i-nz1*n2);
			totalele(i,4)=j+2*n2+(i-nz1*n2+1);
			j=(3*n2+5)*(nz1-1)+2*(i-nz1*n2)+1;
			totalele(i,2)=j;
			totalele(i,1)=j+2;
			totalele(i,5)=j+1;
		}
	//第二柱与梁相交处
		for(i=totale-(nz1+2)*n2+1;i<totale-nz1*n2;i++)
		{
			j=(3*n2+2)*(nz1-1)+2*(i-nz1*n2)+1;
			totalele(i,0)=j+2;
			totalele(i,3)=j;
			totalele(i,7)=j+1;
			totalele(i,6)=j+2*n2+(i-nz1*n2);
			totalele(i,4)=j+2*n2+(i-nz1*n2+1);
			j=(3*n2+5)*(nz1-1)+2*(i-nz1*n2)+1;
			totalele(i,2)=j;
			totalele(i,1)=j+2;
			totalele(i,5)=j+1;
		}
	/////////////////////////////////////////////////////////////////////////
	for(j=0;j<totale;j++)	
	{
		cout<<endl<<"the element is"<<j+1<<"  ";
		for(i=0;i<8;i++)
			cout<<totalele(j,i)<<"   ";
	}
	cout<<endl;
	///////////////////////////////单元信息输入完毕
	///////////形成荷载信息