CPP.rar_c program civil_有限元平面

#include <iostream.h> #include <fstream.h> #include <iomanip.h> #include <math.h> class Q8 { public: int NN, NF, NE, KU, KV, KRX, KRY, KQ, LK; int JEW[9], JLL[17], JE[9][501], JD[3006], KPQ[3][50], JRX[50], JRY[50], JU[50], JV[50]; double H[3+1], D[3+1], U[3007], FF[17], FS[6][6][9], FC[6][6][9], FA[6][6][9], B[5][17], C[5][17], XO[1504], YO[1504], X[50],Y[50],SK[5000], PQ[5][50], RX[50], RY[50], US[50], VS[50],EK[137]; double PO, EO, RU, XC, XA, YC, YA, YN, QQ, VY, RPM; public: /** * calculate */ void calculate() { NF = NN + NN; //NN-节点个数 NF-节点位移（载荷）个数 //弹性矩阵[D] if(D[0]==0){ D[1] = EO / (1.0 - PO * PO); D[2] = D[1] * PO; D[3] = D[1] * (1.0 - PO) / 2; } if(D[0]==2){ D[1] = EO * (1.0 - PO) / (1.0 + PO) / (1.0 - PO - PO); D[2] = D[1] * PO / (1.0 - PO); D[3] = D[1] * (1.0 - PO - PO) / (1.0 - PO) / 2; } SKDD(); NCNA(); VY = RU * (3.1415926 * RPM / 30) / 9.80; if(H[0]==2){ H[1]=1; H[2]=1; H[3]=1; } if(H[0]==3){ H[1] = 5.0 / 9.0; H[2] = 8.0 / 9.0; H[3] = H[1]; } //初始化矩阵 for(int i=0;i<=NF;i++){ U[i]=0; } for (i = 0; i <=LK; i++) { SK[i]=0; } for (i = 1; i <= NE; i++) { int m, l1, l2; for (int j = 1; j <= 8; j++) { m = JE[j][i]; l1 = j + j; l2 = l1 - 1; JLL[l2] = m + m - 1; JLL[l1] = m + m; } FEKP(i); SKKE(); } RIGHT(); for (i = 1; i <= KU; i++) { FIXD(2 * JU[i] - 1, US[i]); } if (KV != 0) { for (int j = 1; j <=KV; j++) { FIXD(2 * JV[j], VS[j]); } } solve(); } /** * XYCA */ void XYCA(int i, int j,int n) { double R, S, F, D; int L; XC = 0; XA = 0; YC = 0; YA = 0; YN = 0; for (int k = 1; k <= 8; k++) { QQ = FC[i][j][k]; R = FA[i][j][k]; S = FS[i][j][k]; L = JE[k][n]; F = XO[L]; D = YO[L]; YN += S * D; XC += QQ * F;//********J(2,1)********// XA += R * F;//*********J(1,1)********// YC += QQ * D;//********J(2,2)********// YA += R * D;//*********J(1,2)********// } } /** * SKDD */ void SKDD() { int L, M; for (int i = 1; i <= NN; i++) { JD[2 * i] = NN; } for (i = 1; i <= NE; i++) { for (int j = 1; j <= 8; j++) { JEW[j] = JE[j][i]; } L = JEW[8]; for (j = 1; j <= 7; j++) { if (JEW[j] < L) L = JEW[j]; } for (j = 1; j <= 8; j++) { M = JEW[j]; if (L < JD[2 * M]) JD[2 * M] = L; } } JD[0]=0; JD[1] = 1; JD[2] = 3; for (i = 2; i <= NN; i++) { M = (i - JD[2 * i]) * 2; JD[2 * i - 1] = JD[2 * i - 2] + M + 1; JD[2 * i] = JD[2 * i - 1] + M + 2; } LK = JD[NF]; } /** * NCNA */ void NCNA() { double XX[6]; double A, C; XX[0] = 0; XX[3] = 0; XX[4] = 1; XX[5]=-1; if(H[0]==2){ XX[1]=1/sqrt(3); XX[2]=-XX[1]; } if(H[0]==3){ XX[3] = sqrt(0.6); XX[1]=-XX[3]; XX[2] = 0; } for (int i = 1; i <6; i++) { C = XX[i]; for (int j = 1; j <6; j++) { A = XX[j]; FS[i][j][1] = (1 + C) * (1 + A) * (C + A - 1) / 4; FS[i][j][2] = (1 + A) * (1 - C * C) / 2; FS[i][j][3] = (1 - C) * (1 + A) * (A - C - 1) / 4; FS[i][j][4] = (1 - C) * (1 - A * A) / 2; FS[i][j][5] = (1 - C) * (1 - A) * ( -C - A - 1) / 4; FS[i][j][6] = (1 - A) * (1 - C * C) / 2; FS[i][j][7] = (1 + C) * (1 - A) * (C - A - 1) / 4; FS[i][j][8] = (1 + C) * (1 - A * A) / 2; FS[i][j][0] = 0; FC[i][j][1] = (1 + A) * (C + C + A) / 4; FC[i][j][2] = (1 + A) * ( -C); FC[i][j][3] = (1 + A) * (C + C - A) / 4; FC[i][j][4] = (A * A - 1) / 2; FC[i][j][5] = (1 - A) * (C + C + A) / 4; FC[i][j][6] = (A - 1) * C; FC[i][j][7] = (1 - A) * (C + C - A) / 4; FC[i][j][8] = (1 - A * A) / 2; FC[i][j][0] = 0; FA[i][j][1] = (1 + C) * (C + A + A) / 4; FA[i][j][2] = (1 - C * C) / 2; FA[i][j][3] = (1 - C) * (A + A - C) / 4; FA[i][j][4] = A * (C - 1); FA[i][j][5] = (1 - C) * (C + A + A) / 4; FA[i][j][6] = (C * C - 1) / 2; FA[i][j][7] = (1 + C) * (A + A - C) / 4; FA[i][j][8] = -A * (1 + C); FA[i][j][0] = 0; } } } /** * FEKP */ void FEKP(int n) { double DTJ, F, Q, W, Z, R, FIX, FIY; int M, L, kk; for (int i = 0; i <=136; i++) { EK[i] = 0; } for (i = 0; i <=NF; i++) { FF[i] = 0; } for (i = 1; i <= H[0]; i++) { for (int j = 1; j <= H[0]; j++) { XYCA(i,j,n); DTJ = XC * YA - XA * YC; F = VY; Q = DTJ * H[i] * H[j]; if(D[0]==2){ F*=YN; Q*=YN; } for (int k = 1; k <= 8; k++) { M = k + k; L = M - 1; W = FC[i][j][k]; Z = FA[i][j][k]; FIX = (YA * W - YC * Z) / DTJ; FIY = (XC * Z - XA * W) / DTJ; B[1][L] = FIX; B[2][L] = 0; B[1][M] = 0; B[2][M] = FIY; C[1][L] = D[1] * FIX; C[2][L] = D[2] * FIX; if(D[0]==0){ R = FS[i][j][k]; B[3][L] = FIY; B[4][L] = 0; B[3][M] = FIX; B[4][M] = 0; C[1][M] = D[2] * FIY; C[2][M] = D[1] * FIY ; C[3][L] = D[3] * FIY; C[4][L] = 0; C[3][M] = D[3] * FIX; C[4][M] = 0; FF[M] += F * Q * R; } if(D[0]==2){ R = FS[i][j][k] / YN; B[3][L] = 0; B[4][L] = FIY; B[3][M] = R; B[4][M] = FIX; C[1][M] = D[2] * (FIY + R); C[2][M] = D[1] * FIY + D[2] * R; C[3][L] = C[2][L]; C[4][L] = D[3] * FIY; C[3][M] = D[1] * R + D[2] * FIY; C[4][M] = D[3] * FIX; FF[M] += F * Q * R * YN; } } for (k = 1; k <= 16; k++) { kk = k * (k - 1) / 2; for (int ii = 1; ii <= k; ii++) { Z = 0; for (int jj = 1; jj <= 4; jj++) { Z += B[jj][k] * C[jj][ii]; } EK[kk+ii]+= Z * Q; } } } } for (i = 1; i <= 16; i++) { M = JLL[i]; U[M] += FF[i]; } } /** * SKKE */ void SKKE() { int IW, JQ, KR, L; for (int i = 1; i <= 16; i++) { IW = JLL[i]; L = i * (i - 1) / 2; for (int j = 1; j <= i; j++) { JQ = JLL[j];