📄 globalelement.cpp
字号:
break;
}
//转换矩阵
mat_TT = m_apEle[0]->m_CaElemTranFactorMatrixT[loop1];
//生成荷载向量
P = m_Load.m_adElemLoadValue_Case2_Cen[loop].m_dElemForceX;
TempCM_DisLoadVector(0,0) = P/2;
TempCM_DisLoadVector(6,0) = P/2;
P = m_Load.m_adElemLoadValue_Case2_Cen[loop].m_dElemForceY;
double P1 = P;
TempCM_DisLoadVector(1,0) = P*(L-X)*(L-X)*(L+2*X)/(L*L*L);
TempCM_DisLoadVector(7,0) = P*X*X*(L+2*(L-X))/(L*L*L);
P = m_Load.m_adElemLoadValue_Case2_Cen[loop].m_dElemForceZ;
double P2 = P;
TempCM_DisLoadVector(2,0) = P*(L-X)*(L-X)*(L+2*X)/(L*L*L);
TempCM_DisLoadVector(8,0) = P*X*X*(L+2*(L-X))/(L*L*L);
//
P = m_Load.m_adElemLoadValue_Case2_Cen[loop].m_dElemMX;
TempCM_DisLoadVector(3,0) = P*(L-X)/L;
TempCM_DisLoadVector(9,0) = P*X/L;
//
P = m_Load.m_adElemLoadValue_Case2_Cen[loop].m_dElemMY;
TempCM_DisLoadVector(4,0) = -6*P*X*(L-X)/(L*L*L) - P2*X*(L-X)*(L-X)/(L*L);
TempCM_DisLoadVector(10,0) = 6*P*X*(L-X)/(L*L*L) + P2*X*X*(L-X)/(L*L);
P = m_Load.m_adElemLoadValue_Case2_Cen[loop].m_dElemMZ;
TempCM_DisLoadVector(5,0) = 6*P*X*(L-X)/(L*L*L) + P1*X*(L-X)*(L-X)/(L*L);
TempCM_DisLoadVector(11,0) = -6*P*X*(L-X)/(L*L*L) - P1*X*X*(L-X)/(L*L);
//整体荷载向量
TempCM_DisLoadVectorG = mat_TT * TempCM_DisLoadVector;
for(int loop2=0;loop2<(m_Node.m_nNodeNum - m_Node.m_nNode_K_Num);loop2++)
{
if(m_nNode1 == m_Node.m_aDOFIndex[loop2].m_aiNode)
{
m_nTemp1 = m_Node.m_aDOFIndex[loop2].m_aiDOFIndex_X;
//if(m_nTemp1<m_Node.m_nFreeDOF)
m_adCase2_ElemCen_LoadVector[m_nTemp1] += TempCM_DisLoadVectorG(0,0);
m_nTemp1 = m_Node.m_aDOFIndex[loop2].m_aiDOFIndex_Y;
//if(m_nTemp1<m_Node.m_nFreeDOF)
m_adCase2_ElemCen_LoadVector[m_nTemp1] += TempCM_DisLoadVectorG(1,0);
m_nTemp1 = m_Node.m_aDOFIndex[loop2].m_aiDOFIndex_Z;
//if(m_nTemp1<m_Node.m_nFreeDOF)
m_adCase2_ElemCen_LoadVector[m_nTemp1] += TempCM_DisLoadVectorG(2,0);
m_nTemp1 = m_Node.m_aDOFIndex[loop2].m_aiDOFIndex_RX;
//if(m_nTemp1<m_Node.m_nFreeDOF)
m_adCase2_ElemCen_LoadVector[m_nTemp1] += TempCM_DisLoadVectorG(3,0);
m_nTemp1 = m_Node.m_aDOFIndex[loop2].m_aiDOFIndex_RY;
//if(m_nTemp1<m_Node.m_nFreeDOF)
m_adCase2_ElemCen_LoadVector[m_nTemp1] += TempCM_DisLoadVectorG(4,0);
m_nTemp1 = m_Node.m_aDOFIndex[loop2].m_aiDOFIndex_RZ;
//if(m_nTemp1<m_Node.m_nFreeDOF)
m_adCase2_ElemCen_LoadVector[m_nTemp1] += TempCM_DisLoadVectorG(5,0);
}
if(m_nNode2 == m_Node.m_aDOFIndex[loop2].m_aiNode)
{
m_nTemp1 = m_Node.m_aDOFIndex[loop2].m_aiDOFIndex_X;
//if(m_nTemp1<m_Node.m_nFreeDOF)
m_adCase2_ElemCen_LoadVector[m_nTemp1] += TempCM_DisLoadVectorG(6,0);
m_nTemp1 = m_Node.m_aDOFIndex[loop2].m_aiDOFIndex_Y;
//if(m_nTemp1<m_Node.m_nFreeDOF)
m_adCase2_ElemCen_LoadVector[m_nTemp1] += TempCM_DisLoadVectorG(7,0);
m_nTemp1 = m_Node.m_aDOFIndex[loop2].m_aiDOFIndex_Z;
//if(m_nTemp1<m_Node.m_nFreeDOF)
m_adCase2_ElemCen_LoadVector[m_nTemp1] += TempCM_DisLoadVectorG(8,0);
m_nTemp1 = m_Node.m_aDOFIndex[loop2].m_aiDOFIndex_RX;
//if(m_nTemp1<m_Node.m_nFreeDOF)
m_adCase2_ElemCen_LoadVector[m_nTemp1] += TempCM_DisLoadVectorG(9,0);
m_nTemp1 = m_Node.m_aDOFIndex[loop2].m_aiDOFIndex_RY;
//if(m_nTemp1<m_Node.m_nFreeDOF)
m_adCase2_ElemCen_LoadVector[m_nTemp1] += TempCM_DisLoadVectorG(10,0);
m_nTemp1 = m_Node.m_aDOFIndex[loop2].m_aiDOFIndex_RZ;
//if(m_nTemp1<m_Node.m_nFreeDOF)
m_adCase2_ElemCen_LoadVector[m_nTemp1] += TempCM_DisLoadVectorG(11,0);
}
}
//这里的CMX，CMY，CMZ弃用了，太麻烦了！
//调整轴力内力FX
P = m_Load.m_adElemLoadValue_Case2_Cen[loop].m_dElemForceX;
for(int loop3=0;loop3<13;loop3++)
{
double x = loop3*L/12;
RA = RB = P/2;
//RA = P*(L-X)/L;
//RB = -P*X/L;
if(x<=X)
m_apEle[0]->m_aElemInForceDef[loop1].m_dElemFX[loop3] += RA;
else
m_apEle[0]->m_aElemInForceDef[loop1].m_dElemFX[loop3] += -RB;
}
//调整剪力ShearY
P = m_Load.m_adElemLoadValue_Case2_Cen[loop].m_dElemForceY;
for(loop3=0;loop3<13;loop3++)
{
RA = P*(L-X)*(L-X)*(L+2*X)/(L*L*L);
RB = P*X*X*(L+(L-X)*2)/(L*L*L);
double x = loop3*L/12;
if(x<=X)
m_apEle[0]->m_aElemInForceDef[loop1].m_dElemShearY[loop3] += RA;
else
m_apEle[0]->m_aElemInForceDef[loop1].m_dElemShearY[loop3] += -RB;
}
//调整剪力ShearZ
P = m_Load.m_adElemLoadValue_Case2_Cen[loop].m_dElemForceZ;
for(loop3=0;loop3<13;loop3++)
{
RA = P*(L-X)*(L-X)*(L+2*X)/(L*L*L);
RB = P*X*X*(L+(L-X)*2)/(L*L*L);
double x = loop3*L/12;
if(x<=X)
m_apEle[0]->m_aElemInForceDef[loop1].m_dElemShearZ[loop3] += RA;
else
m_apEle[0]->m_aElemInForceDef[loop1].m_dElemShearZ[loop3] += -RB;
}
//调整扭距MT（不做调整）
//调整弯距MY
P = m_Load.m_adElemLoadValue_Case2_Cen[loop].m_dElemForceZ;
for(loop3=0;loop3<13;loop3++)
{
RA = P*(L-X)*(L-X)*(L+2*X)/(L*L*L);
RB = P*X*X*(L+(L-X)*2)/(L*L*L);
//
MA = P*X*(L-X)*(L-X)/(L*L);
MB = -P*X*X*(L-X)/(L*L);
double x = loop3*L/12;
if(x<=X)
m_apEle[0]->m_aElemInForceDef[loop1].m_dElemMY[loop3] += -MA + RA*x;
else
m_apEle[0]->m_aElemInForceDef[loop1].m_dElemMY[loop3] += MB + RB*(L-x);
}
//调整弯距MZ
P = m_Load.m_adElemLoadValue_Case2_Cen[loop].m_dElemForceY;
for(loop3=0;loop3<13;loop3++)
{
RA = P*(L-X)*(L-X)*(L+2*X)/(L*L*L);
RB = P*X*X*(L+(L-X)*2)/(L*L*L);
//
MA = P*X*(L-X)*(L-X)/(L*L);
MB = -P*X*X*(L-X)/(L*L);
double x = loop3*L/12;
if(x<=X)
m_apEle[0]->m_aElemInForceDef[loop1].m_dElemMZ[loop3] += MA - RA*x;
else
m_apEle[0]->m_aElemInForceDef[loop1].m_dElemMZ[loop3] += -MB - RB*(L-x);
}
}
}
}
}
void CGlobalElement::OutputParameter()
{
ofstream fout;
fout.open("LoadVectorTest.tst");
fout.width(12);
fout<<"1_NodeLoad";
fout.width(16);
fout<<"1_ElemDisLoad";
fout.width(16);
fout<<"1_ElemCenLoad";
fout.width(16);
fout<<"2_NodeLoad:";
fout.width(16);
fout<<"2_ElemDisLoad";
fout.width(16);
fout<<"2_ElemCenLoad";
fout.width(20);
fout<<"1_LoadVector_C1:";
fout.width(20);
fout<<"2_LoadVector_C2:";
fout.width(20);
fout<<"LoadVector:"<<endl;
for(int loop=0;loop<m_Node.m_nTotalDOF;loop++)
{
fout.width(12);
fout<<m_adCase1_Node_LoadVector[loop];
fout.width(16);
fout<<m_adCase1_ElemDis_LoadVector[loop];
fout.width(16);
fout<<m_adCase1_ElemCen_LoadVector[loop];
fout.width(16);
fout<<m_adCase2_Node_LoadVector[loop];
fout.width(16);
fout<<m_adCase2_ElemDis_LoadVector[loop];
fout.width(16);
fout<<m_adCase2_ElemCen_LoadVector[loop];
fout.width(20);
fout<<m_adLoadVector_C1[loop];
fout.width(20);
fout<<m_adLoadVector_C2[loop];
fout.width(20);
fout<<LoadVector[loop]<<endl;
}
fout<<endl;
fout.close();
}
void CGlobalElement::StructureSolve()
{
StiffAssemble();
Loading();
DisplacementSolve();
CalcuInternalForce();
CalElemDeflection();
CalcuSupportReaction();
}
void CGlobalElement::DisplacementSolve()
{
//***********************求解工况一的节点位移**************************
//荷载向量1拷贝
m_adDisp1 = new double [m_Node.m_nTotalDOF];
for(int loop=0;loop<m_Node.m_nTotalDOF;loop++)
m_adDisp1[loop] = m_adLoadVector_C1[loop];
CSparseMatrix tempGK1;
CSparseMatrix tempGK2;
tempGK1 = m_apEle[0]->m_smatGK;
tempGK2 = tempGK1;
//求解节点位移
tempGK1.LdltSolve(m_Node.m_nFreeDOF,m_adDisp1);
//**************求解工况二的节点位移（供求解内力使用）*******************
//荷载向量2拷贝
m_adDisp2 = new double [m_Node.m_nTotalDOF];
for(loop=0;loop<m_Node.m_nTotalDOF;loop++)
m_adDisp2[loop] = m_adLoadVector_C2[loop];
//求解节点位移
tempGK2.LdltSolve(m_Node.m_nFreeDOF,m_adDisp2);
//输出位移向量
OutputDisplacement();
//输出单元挠度
//不能再执行了，执行两遍了!!!
//CalElemDeflection();
}
void CGlobalElement::CalElemDeflection()
{
int Node1,Node2;
CMatrix DispG(12,1,0),DispL(12,1,0);
CMatrix mat_T;
int m_nTemp;
//把节点位移转化到单元位移
for(int loop=0;loop<m_apEle[0]->m_nElementNum;loop++)
{
mat_T = m_apEle[0]->m_CaElemTranFactorMatrix[loop];
Node1 = m_apEle[0]->m_aBeamElemProValue[loop].m_nElemNode1;
Node2 = m_apEle[0]->m_aBeamElemProValue[loop].m_nElemNode2;
//
for(int loop1=0;loop1<(m_Node.m_nNodeNum - m_Node.m_nNode_K_Num);loop1++)
{
if(Node1 == m_Node.m_aDOFIndex[loop1].m_aiNode)
{
m_nTemp = m_Node.m_aDOFIndex[loop1].m_aiDOFIndex_X;
DispG(0,0) = m_adDisp1[m_nTemp];
m_nTemp = m_Node.m_aDOFIndex[loop1].m_aiDOFIndex_Y;
DispG(1,0) = m_adDisp1[m_nTemp];
m_nTemp = m_Node.m_aDOFIndex[loop1].m_aiDOFIndex_Z;
DispG(2,0) = m_adDisp1[m_nTemp];
m_nTemp = m_Node.m_aDOFIndex[loop1].m_aiDOFIndex_RX;
DispG(3,0) = m_adDisp1[m_nTemp];
m_nTemp = m_Node.m_aDOFIndex[loop1].m_aiDOFIndex_RY;
DispG(4,0) = m_adDisp1[m_nTemp];
m_nTemp = m_Node.m_aDOFIndex[loop1].m_aiDOFIndex_RZ;
DispG(5,0) = m_adDisp1[m_nTemp];
}
if(Node2 == m_Node.m_aDOFIndex[loop1].m_aiNode)
{
m_nTemp = m_Node.m_aDOFIndex[loop1].m_aiDOFIndex_X;
DispG(6,0) = m_adDisp1[m_nTemp];
m_nTemp = m_Node.m_aDOFIndex[loop1].m_aiDOFIndex_Y;
DispG(7,0) = m_adDisp1[m_nTemp];
m_nTemp = m_Node.m_aDOFIndex[loop1].m_aiDOFIndex_Z;
DispG(8,0) = m_adDisp1[m_nTemp];
m_nTemp = m_Node.m_aDOFIndex[loop1].m_aiDOFIndex_RX;
DispG(9,0) = m_adDisp1[m_nTemp];
m_nTemp = m_Node.m_aDOFIndex[loop1].m_aiDOFIndex_RY;
DispG(10,0) = m_adDisp1[m_nTemp];
m_nTemp = m_Node.m_aDOFIndex[loop1].m_aiDOFIndex_RZ;
DispG(11,0) = m_adDisp1[m_nTemp];
}
}
//转化完成
DispL = mat_T*DispG;
double DY1 = DispL(1,0);
double DZ1 = DispL(2,0);
double DY2 = DispL(7,0);
double DZ2 = DispL(8,0);
for(int loop2=0;loop2<13;loop2++)
{
m_apEle[0]->m_aElemDeflection[loop].m_dDelf_Y[loop2] += (12-loop2)*DY1/12 + loop2*DY2/12 ;
m_apEle[0]->m_aElemDeflection[loop].m_dDelf_Z[loop2] += (12-loop2)*DZ1/12 + loop2*DZ2/12 ;
}
}
//输出单元挠度
ofstream fout;
fout.open("Deflection.twm");
for(loop=0;loop<m_apEle[0]->m_nElementNum;loop++)
{
fout.width(15);
fout<<"单元号";
fout.width(15);
fout<<m_apEle[0]->m_aBeamElemProValue[loop].m_nElementNo<<endl;
for(int loop1=0;loop1<13;loop1++)
{
fout.width(15);
fout<<m_apEle[0]->m_aElemDeflection[loop].m_dDelf_Y[loop1];
}
fout<<endl;
for(loop1=0;loop1<13;loop1++)
{
fout.width(15);
fout<<m_apEle[0]->m_aElemDeflection[loop].m_dDelf_Z[loop1];
}
fout<<endl;
fout<<endl;
}
fout.close();
}
void CGlobalElement::OutputDisplacement()
{
//输出节点位移(CASE 1)
ofstream fout;
fout.open("Displacement.tst");
fout.width(20);
fout<<"Displacement_1.tst";
fout.width(20);
fout<<"Displacement_2.tst"<<endl;
fout<<endl;
for(int loop=0;loop<m_Node.m_nTotalDOF;loop++)
{
fout.width(20);
fout<<m_adDisp1[loop];
fout.width(20);
fout<<m_adDisp2[loop]<<endl;
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -