oopfeview.cpp

一个计算悬臂梁的有限元vc源码

	double dTextHeight;
	double dX,dY,dX1,dY1;
	CString sText;
	CTypedPtrArray <CPtrArray,CBaseElement*>& apEle=pDoc->m_glEle.m_apEle;
	CNode& Node=pDoc->m_glEle.m_Node;
	CString sFileName;
	ofstream fout;

	if(!pDoc->m_glEle.m_bFlagForceAnalyzed) return;

	CFileDialog dlg(FALSE);
	char sFilter[]="DXF(*.dxf)\0*.dxf\0\0";
	dlg.m_ofn.lpstrFilter=sFilter;

	if(dlg.DoModal()==IDOK){
		sFileName=dlg.GetPathName();
		if(sFileName.GetLength()==0) return;
		if(sFileName.Find(".dxf")<0) sFileName+=".dxf";
		fout.open(sFileName);
	}

	dTextHeight=m_nFontSize/m_dScaleRtoD;
	dMaxStressWidth=dTextHeight;
		
	fout<<"  0"<<endl;
	fout<<"SECTION"<<endl;
	fout<<"  2"<<endl;
	fout<<"ENTITIES"<<endl;

	OutputDxfElement(fout);

	nEle=apEle.GetSize();
	dMinStress=dMaxStress=0.0;
	bBuf=false;
	for(loop=0;loop<nEle;loop++){
		if(apEle[loop]->GetElementType()==PLANE_8NODE_ELE){
			((CPlane8NodeEle*)apEle[loop])->
				GetMainStress(adMainStress,adMainStress1,adMainStressAngle);
			for(loop1=0;loop1<9;loop1++){
				if(adMainStress[loop1]<0.0){
					dMinStress=dMaxStress=adMainStress[loop1];
					bBuf=true;
					break;
				}
				if(adMainStress1[loop1]<0.0){
					dMinStress=dMaxStress=adMainStress1[loop1];
					bBuf=true;
					break;
				}
			}
		}
		if(bBuf) break;
	}

	for(loop=0;loop<nEle;loop++){
		if(apEle[loop]->GetElementType()==PLANE_8NODE_ELE){
			((CPlane8NodeEle*)apEle[loop])->
				GetMainStress(adMainStress,adMainStress1,adMainStressAngle);
			for(loop1=0;loop1<9;loop1++){
				if(adMainStress[loop1]<0.0){
					if(dMaxStress<adMainStress[loop1])
						dMaxStress=adMainStress[loop1];
					if(dMinStress>adMainStress[loop1])
						dMinStress=adMainStress[loop1];
				}
				if(adMainStress1[loop1]<0){
					if(dMaxStress<adMainStress1[loop1])
						dMaxStress=adMainStress1[loop1];
					if(dMinStress>adMainStress1[loop1])
						dMinStress=adMainStress1[loop1];
				}
			}
		}
	}

	if(dMaxStress-dMinStress>0.0)
		dScaleStressWidth=dMaxStressWidth/(dMaxStress-dMinStress);

	for(loop=0;loop<nEle;loop++){
		if(apEle[loop]->GetElementType()==PLANE_8NODE_ELE){
			apEle[loop]->GetNode(aiNode,iBuf);
			for(loop1=0;loop1<8;loop1++){
				adEleNodeX[loop1]=Node.GetX(aiNode[loop1]);
				adEleNodeY[loop1]=Node.GetY(aiNode[loop1]);
			}
			
			dMaxX=dMinX=adEleNodeX[0];
			dMaxY=dMinY=adEleNodeY[0];
			for(loop1=1;loop1<8;loop1++){
				if(dMaxX<adEleNodeX[loop1]) dMaxX=adEleNodeX[loop1];
				if(dMinX>adEleNodeX[loop1]) dMinX=adEleNodeX[loop1];
				if(dMaxY<adEleNodeY[loop1]) dMaxY=adEleNodeY[loop1];
				if(dMinY>adEleNodeY[loop1]) dMinY=adEleNodeY[loop1];
			}
			if((dMaxX-dMinX)>(dMaxY-dMinY)){
				dLength=(dMaxX-dMinX)/8.0;
			}
			else{
				dLength=(dMaxY-dMinY)/8.0;
			}
			((CPlane8NodeEle*)apEle[loop])->GetGaussPtGCoordinate(adGaussX,adGaussY);
			((CPlane8NodeEle*)apEle[loop])->
				GetMainStress(adMainStress,adMainStress1,adMainStressAngle);
			for(loop1=0;loop1<9;loop1++){
				if(adMainStress[loop1]<0.0){
					dCos=cos(adMainStressAngle[loop1]);
					dSin=sin(adMainStressAngle[loop1]);
					dX=adGaussX[loop1]-dCos*dLength;
					dY=adGaussY[loop1]-dSin*dLength;
					dX1=adGaussX[loop1]+dCos*dLength;
					dY1=adGaussY[loop1]+dSin*dLength;
					dStressWidth=dMaxStressWidth-(adMainStress[loop1]-dMinStress)*dScaleStressWidth;
					DXFPline(fout,dX,dY,dX1,dY1,dStressWidth);
				}
				if(adMainStress1[loop1]<0.0){
					dCos=cos(adMainStressAngle[loop1]);
					dSin=sin(adMainStressAngle[loop1]);
					dX=adGaussX[loop1]+dSin*dLength;
					dY=adGaussY[loop1]-dCos*dLength;
					dX1=adGaussX[loop1]-dSin*dLength;
					dY1=adGaussY[loop1]+dCos*dLength;
					dStressWidth=dMaxStressWidth-(adMainStress1[loop1]-dMinStress)*dScaleStressWidth;
					DXFPline(fout,dX,dY,dX1,dY1,dStressWidth);
				}
			}
		}
	}

	dX=m_dMinX-13.0*dTextHeight;
	dX1=m_dMinX-10.0*dTextHeight;
	for(loop=0;loop<9;loop++){
		dY=dY1=m_dMinY+1.5*dTextHeight*loop;
		dStressWidth=dMaxStressWidth*loop/8.0;
		DXFPline(fout,dX,dY,dX1,dY1,dStressWidth);
		sText.Format("%0.3f",dMaxStress+(dMinStress-dMaxStress)*loop/8.0);
		DXFText(fout,dX1+dTextHeight,dY1,dTextHeight,0.0,sText);
	}

	fout<<"  0"<<endl;
	fout<<"ENDSEC"<<endl;
	fout<<"  0"<<endl;
	fout<<"EOF"<<endl;
	fout.close();		
}

void COOPFEView::OnDisplayModeShapesPrevious() 
{
	// TODO: Add your command handler code here
	COOPFEDoc* pDoc = GetDocument();
	ASSERT_VALID(pDoc);
	pDoc->m_iCurMode--;
	if(pDoc->m_iCurMode<0) pDoc->m_iCurMode=0;
	pDoc->m_iDisplay=DISPLAY_MODE_SHAPES;
	Invalidate();
	UpdateWindow();
}

void COOPFEView::OnDisplayModeShapesNext() 
{
	// TODO: Add your command handler code here
	COOPFEDoc* pDoc = GetDocument();
	ASSERT_VALID(pDoc);
	pDoc->m_iCurMode++;
	if(pDoc->m_iCurMode>=pDoc->m_glEle.m_nMode)
		pDoc->m_iCurMode=pDoc->m_glEle.m_nMode-1;
	pDoc->m_iDisplay=DISPLAY_MODE_SHAPES;
	Invalidate();
	UpdateWindow();
}

void COOPFEView::DisplayModeShapes(CDC *pDC, COOPFEDoc *pDoc)
{
	int loop,loop1,iBuf,nNode,nEle,nEleNode,aiNode[8];
	int iEleType,iCurMode,nFreeDOF,nTotalDOF;
	double adFy[11];
	double dScaleDispRtoD;
	double dCos,dSin,dLength;
	double dXi0,dXj0,dXk0,dYi0,dYj0,dYk0;
	double dUei,dUej;
	double dMaxDisp,*adDisp,dBuf;
	CTypedPtrArray <CPtrArray,CBaseElement*>& apEle=pDoc->m_glEle.m_apEle;
	CNode& Node=pDoc->m_glEle.m_Node;
	CMatrix& matEigenVector=pDoc->m_glEle.m_matEigenVector;
	CPoint aPt[20];
	CPen* pOldPen;
	CFont font,*pOldFont;
	CString sText;
	
	DisplayElement(pDC,pDoc);
	if(!pDoc->m_glEle.m_bFlagModeAnalyzed) return;

	CPen pen(PS_SOLID,1,RGB(255,0,0));
	pOldPen=pDC->SelectObject(&pen);
	
	iCurMode=pDoc->m_iCurMode;
	nEle=apEle.GetSize();
	nNode=Node.GetNum();

	nFreeDOF=Node.GetFreeDOF();
	nTotalDOF=Node.GetTotalDOF();

	adDisp=new double[nTotalDOF];

	dBuf=0;
	for(loop=0;loop<nFreeDOF;loop++)
		if(fabs(dBuf)<fabs(matEigenVector(loop,iCurMode)))
			dBuf=matEigenVector(loop,iCurMode);
	for(loop=0;loop<nFreeDOF;loop++)
		adDisp[loop]=matEigenVector(loop,iCurMode)/dBuf;
	for(loop=nFreeDOF;loop<nTotalDOF;loop++){
		adDisp[loop]=0.0;
	}

	dMaxDisp=1.0;
	for(loop=0;loop<nEle;loop++){
		iEleType=apEle[loop]->GetElementType();
		if(iEleType==BEAM||iEleType==BEAM_HINGE){
			BeamDeflection(apEle[loop],Node,adDisp,adFy,dUei,dUej);
			for(loop1=0;loop1<11;loop1++){
				if(dMaxDisp<fabs(adFy[loop1]))
					dMaxDisp=fabs(adFy[loop1]);
			}
		}
	}
	dScaleDispRtoD=m_nMaxDispDot/dMaxDisp;

	for(loop=0;loop<nEle;loop++){
		iEleType=apEle[loop]->GetElementType();
		switch(iEleType){
		case BEAM:
		case BEAM_HINGE:
			apEle[loop]->GetNode(aiNode,iBuf);
			dLength=((CBeam*)apEle[loop])->GetLength();
			dXi0=Node.GetX(aiNode[0]);
			dYi0=Node.GetY(aiNode[0]);
			dXj0=Node.GetX(aiNode[1]);
			dYj0=Node.GetY(aiNode[1]);
			dCos=(dXj0-dXi0)/dLength;
			dSin=(dYj0-dYi0)/dLength;
			BeamDeflection(apEle[loop],Node,adDisp,adFy,dUei,dUej);
			for(loop1=0;loop1<11;loop1++){
				dXk0=m_ptBase.x+(int)((dXi0+(dXj0-dXi0)*loop1/10.0-m_dMinX)*m_dScaleRtoD);
				dYk0=m_ptBase.y-(int)((dYi0+(dYj0-dYi0)*loop1/10.0-m_dMinY)*m_dScaleRtoD);
				aPt[loop1].x=(int)dXk0+(int)(((dUei+(dUej-dUei)*loop1/10.0)*dCos+adFy[loop1]*dSin)*dScaleDispRtoD);
				aPt[loop1].y=(int)dYk0-(int)(((dUei+(dUej-dUei)*loop1/10.0)*dSin-adFy[loop1]*dCos)*dScaleDispRtoD);
			}
			pDC->MoveTo(aPt[0]);
			for(loop1=0;loop1<11;loop1++){
				pDC->LineTo(aPt[loop1]);
			}
			break;
		case TRUSS:
			apEle[loop]->GetNode(aiNode,nEleNode);
			iBuf=aiNode[0];
			aPt[0].x=m_ptBase.x+(int)((Node.GetX(iBuf)-m_dMinX)*m_dScaleRtoD+adDisp[Node.GetXDOFIndex(iBuf)]*dScaleDispRtoD);
			aPt[0].y=m_ptBase.y-(int)((Node.GetY(iBuf)-m_dMinY)*m_dScaleRtoD+adDisp[Node.GetYDOFIndex(iBuf)]*dScaleDispRtoD);
			iBuf=aiNode[1];
			aPt[1].x=m_ptBase.x+(int)((Node.GetX(iBuf)-m_dMinX)*m_dScaleRtoD+adDisp[Node.GetXDOFIndex(iBuf)]*dScaleDispRtoD);
			aPt[1].y=m_ptBase.y-(int)((Node.GetY(iBuf)-m_dMinY)*m_dScaleRtoD+adDisp[Node.GetYDOFIndex(iBuf)]*dScaleDispRtoD);
			pDC->MoveTo(aPt[0]);
			pDC->LineTo(aPt[1]);
			break;
		case PLANE_8NODE_ELE:
			apEle[loop]->GetNode(aiNode,nEleNode);
			for(loop1=0;loop1<8;loop1++){
				iBuf=aiNode[loop1];
				aPt[loop1].x=m_ptBase.x+(int)((Node.GetX(iBuf)-m_dMinX)*m_dScaleRtoD+adDisp[Node.GetXDOFIndex(iBuf)]*dScaleDispRtoD);
				aPt[loop1].y=m_ptBase.y-(int)((Node.GetY(iBuf)-m_dMinY)*m_dScaleRtoD+adDisp[Node.GetYDOFIndex(iBuf)]*dScaleDispRtoD);
			}

			pDC->MoveTo(aPt[7]);
			for(loop1=0;loop1<8;loop1++){
				pDC->LineTo(aPt[loop1]);
			}
			break;
		}
	}
	
	SetFont(font,"Times New Roman",m_nFontSize,0.0);
	pOldFont=pDC->SelectObject(&font);
	pDC->SetTextColor(RGB(0,0,0));

	sText.Format("Mode: %d: ",iCurMode+1);
	pDC->TextOut(20,40,sText);

	sText.Format("%f (1/s)",pDoc->m_glEle.m_adEigen[pDoc->m_iCurMode]/(2.0*PI));
	pDC->TextOut(100,40,sText);


	pDC->SelectObject(pOldFont);
	pDC->SelectObject(pOldPen);
	delete adDisp;
}

void COOPFEView::OnDisplayVibration() 
{
	// TODO: Add your command handler code here
	COOPFEDoc* pDoc = GetDocument();
	ASSERT_VALID(pDoc);
	if(!pDoc->m_glEle.m_bFlagDynamicAnalyzed) return;

	int loop,loop1,loop2,iBuf,iEleType,iLineNode,aiNode[10];
	int nEle,nNode,nFreeDOF,nTotalDOF,nEleNode;
	int nTimeStep,nLineNode;
	short *aiPtX,*aiPtY;
	double dCos,dSin,dLength,dXi0,dYi0,dXj0,dYj0;
	int iXk0,iYk0,iXk,iYk;
	double *adDisp,dMaxDisp;
	double dUei,dUej,adFy[11];
	double dScaleDispRtoD;
	CTypedPtrArray <CPtrArray,CBaseElement*>& apEle=pDoc->m_glEle.m_apEle;
	CNode& Node=pDoc->m_glEle.m_Node;

	nTimeStep=(int)(pDoc->m_glEle.m_dResponseDuration/pDoc->m_glEle.m_dTimeStep);
	nEle=apEle.GetSize();
	nNode=Node.GetNum();
	nTotalDOF=Node.GetTotalDOF();
	nFreeDOF=Node.GetFreeDOF();

	adDisp=new double[nTotalDOF];
	for(loop=nFreeDOF;loop<nTotalDOF;loop++){
		adDisp[loop]=0.0;
	}

	fstream fin,fout;
	fin.open("dydata.tmp",ios::in|ios::binary);
	fout.open("dydata1.tmp",ios::out|ios::binary);

	dMaxDisp=0.0;
	for(loop=0;loop<nTimeStep;loop++){
		fin.read((char*)adDisp,sizeof(double)*nFreeDOF);
		for(loop1=0;loop1<nNode;loop1++){
			iBuf=Node.GetXDOFIndex(loop1);
			if(dMaxDisp<fabs(adDisp[iBuf]))
				dMaxDisp=fabs(adDisp[iBuf]);
			iBuf=Node.GetYDOFIndex(loop1);
			if(dMaxDisp<fabs(adDisp[iBuf]))
				dMaxDisp=fabs(adDisp[iBuf]);
		}

		for(loop1=0;loop1<nEle;loop1++){
			iEleType=apEle[loop1]->GetElementType();
			if(iEleType==BEAM||iEleType==BEAM_HINGE){
				BeamDeflection(apEle[loop1],Node,adDisp,adFy,dUei,dUej);
				for(loop2=0;loop2<11;loop2++){
					if(dMaxDisp<fabs(adFy[loop2]))
						dMaxDisp=fabs(adFy[loop2]);
				}
			}
		}
	}

	dScaleDispRtoD=m_nMaxDispDot/dMaxDisp;
	nLineNode=0;
	for(loop=0;loop<nEle;loop++){
		iEleType=apEle[loop]->GetElementType();
		switch(iEleType){
		case BEAM:
		case BEAM_HINGE:
			nLineNode+=11;
			break;
		case TRUSS:
			nLineNode+=2;
			break;
		case PLANE_8NODE_ELE:
			nLineNode+=9;
			break;
		}
	}
	fout.write((char*)&nTimeStep,sizeof(int));
	fout.write((char*)&nLineNode,sizeof(int));
	fin.seekg(0,ios::beg);
	aiPtX=new short[nLineNode];
	aiPtY=new short[nLineNode];
	for(loop=0;loop<nTimeStep;loop++){
		iLineNode=0;
		fin.read((char*)adDisp,sizeof(double)*nFreeDOF);
		for(loop1=0;loop1<nEle;loop1++){
			iEleType=apEle[loop1]->GetElementType();
			switch(iEleType){
			case BEAM:
			case BEAM_HINGE:
				apEle[loop1]->GetNode(aiNode,iBuf);
				dLength=((CBeam*)apEle[loop1])->GetLength();
				dXi0=Node.GetX(aiNode[0]);
				dYi0=Node.GetY(aiNode[0]);
				dXj0=Node.GetX(aiNode[1]);
				dYj0=Node.GetY(aiNode[1]);
				dCos=(dXj0-dXi0)/dLength;
				dSin=(dYj0-dYi0)/dLength;
				BeamDeflection(apEle[loop1],Node,adDisp,adFy,dUei,dUej);
				for(loop2=0;loop2<11;loop2++){
					iXk0=m_ptBase.x+(int)((dXi0+(dXj0-dXi0)*loop2/10.0-m_dMinX)*m_dScaleRtoD);
					iYk0=m_ptBase.y-(int)((dYi0+(dYj0-dYi0)*loop2/10.0-m_dMinY)*m_dScaleRtoD);
					iXk=iXk0+(int)(((dUei+(dUej-dUei)*loop2/10.0)*dCos+adFy[loop2]*dSin)*dScaleDispRtoD);
					iYk=iYk0-(int)(((dUei+(dUej-dUei)*loop2/10.0)*dSin-adFy[loop2]*dCos)*dScaleDispRtoD);
					aiPtX[iLineNode]=iXk;
					aiPtY[iLineNode]=iYk;
					iLineNode++;
				}
				break;
			case TRUSS:
				apEle[loop1]->GetNode(aiNode,nEleNode);
				iBuf=aiNode[0];
				aiPtX[iLineNode]=m_ptBase.x+(int)((Node.GetX(iBuf)-m_dMinX)*m_dScaleRtoD+adDisp[Node.GetXDOFIndex(iBuf)]*dScaleDispRtoD);
				aiPtY[iLineNode]=m_ptBase.y-(int)((Node.GetY(iBuf)-m_dMinY)*m_dScaleRtoD+adDisp[Node.GetYDOFIndex(iBuf)]*dScaleDispRtoD);
				iLineNode++;
				iBuf=aiNode[1];
				aiPtX[iLineNode]=m_ptBase.x+(int)((Node.GetX(iBuf)-m_dMinX)*m_dScaleRtoD+adDisp[Node.GetXDOFIndex(iBuf)]*dScaleDispRtoD);
				aiPtY[iLineNode]=m_ptBase.y-(int)((Node.GetY(iBuf)-m_dMinY)*m_dScaleRtoD+adDisp[Node.GetYDOFIndex(iBuf)]*dScaleDispRtoD);
				iLineNode++;
				break;
			case PLANE_8NODE_ELE:
				apEle[loop1]->GetNode(aiNode,nEleNode);
				for(loop2=0;loop2<8;loop2++){
					iBuf=aiNode[loop2];
					aiPtX[iLineNode]=m_ptBase.x+(int)((Node.GetX(iBuf)-m_dMinX)*m_dScaleRtoD+adDisp[Node.GetXDOFIndex(iBuf)]*dScaleDispRtoD);
					aiPtY[iLineNode]=m_ptBase.y-(int)((Node.GetY(iBuf)-m_dMinY)*m_dScaleRtoD+adDisp[Node.GetYDOFIndex(iBuf)]*dScaleDispRtoD);
					iLineNode++;
				}
				iBuf=aiNode[0];
				aiPtX[iLineNode]=m_ptBase.x+(int)((Node.GetX(iBuf)-m_dMinX)*m_dScaleRtoD+adDisp[Node.GetXDOFIndex(iBuf)]*dScaleDispRtoD);
				aiPtY[iLineNode]=m_ptBase.y-(int)((Node.GetY(iBuf)-m_dMinY)*m_dScaleRtoD+adDisp[Node.GetYDOFIndex(iBuf)]*dScaleDispRtoD);
				iLineNode++;
				break;
			}
		}
		fout.write((char*)aiPtX,sizeof(short)*nLineNode);
		fout.write((char*)aiPtY,sizeof(short)*nLineNode);
	}
	fin.close();
	fout.close();
	delete adDisp;
	delete aiPtX;
	delete aiPtY;
	pDoc->m_iDisplay=DISPLAY_VIBRATION;
	Invalidate();
	UpdateWindow();
}


void COOPFEView::BeamDeflection(CBaseElement* pEle,CNode& Node,doub