有限元第9章有限元法程序设计 下载本文

个人收集整理 仅供参考学习

}

}

}

{ }

dBuf+=pGK[iAdd_ii -i+k]*pGK[iAdd_jj-j+k]

/pGK[pDiag[k]];

cout<<\总刚系数有错!\return false;

pGK[iAdd_ii-i+j]-=dBuf;

if(pGK[0]==0.0) { }

pB[0]=pB[0]/pGK[0]; for(i=1;i

dBuf=0.0;

iBegin=i-pDiag[i]+pDiag[i-1]+1; iAdd_ii=pDiag[i]; for(j=iBegin;j

dBuf+=pB[j]*pGK[iAdd_ii-i+j]; cout<<\总刚系数有错!\return false;

45

个人收集整理 仅供参考学习

}

pB[i]-=dBuf;

if(pGK[iAdd_ii]==0.0) { }

pB[i]/=pGK[iAdd_ii];

cout<<\总刚系数有错!\return false;

nHalfBand=1; iRow=nRow-1;

for(i=nRow-2;i>=0;i--) {

nHalfBand++;

while(nHalfBand>(pDiag[iRow]-pDiag[iRow-1])){ } dBuf=0.0; iAdd_ii=pDiag[i]; for(j=1;j

if(j<(pDiag[j+i]-pDiag[j+i-1])) {

if(pGK[iAdd_ii]==0.0) nHalfBand--; iRow--;

46

个人收集整理 仅供参考学习

}

}

}

}

{ }

iAdd_kk=pDiag[j+i];

dBuf+=pGK[iAdd_kk-j]*pB[j+i]/pGK[iAdd_ii];

cout<<\总刚系数有错!\return false;

pB[i]-=dBuf;

return true;

9.4.5 杆端力和支座反力的计算

节点位移求出之后,由杆件的单元定位向量可以得每根杆件两端点的位移,并从文件ElemStiff、dat读得每根杆件的单元刚度矩阵,从而由单元刚度方程(5.2)或(5.36)求得杆端力。程序如下:

void InternalForceCalcu(int nTotalElem,Element* pElem,Node* pNode,double* pDisp)

{

int i,j,e;

double** pKe0=TwoArrayDoubAlloc(4,4); double** pKe1=TwoArrayDoubAlloc(6,6); ifstream fin1(\if(!fin1) {

47

个人收集整理 仅供参考学习

}

cout<<\单元刚度矩阵输入文件打开失败!\exit(-1);

for(e=0;e

TwoArrayFree(4,pKe0); TwoArrayFree(6,pKe1); fin1.close();

if(remove(\{

if((pElem+e)->iType==TRUSS) { } else { }

for(i=0;i<6;i++)

for(j=0;j<6;j++)

fin1.read((char*)&pKe1[i][j],sizeof(double));

for(i=0;i<4;i++)

for(j=0;j<4;j++)

fin1.read((char*)&pKe0[i][j],sizeof(double));

TrussInternalForceCalcu(pElem+e,pNode,pDisp,pKe0);

FrameInternalForceCalcu(pElem+e,pNode,pDisp,pKe1);

48