本帖最后由 欧阳中华 于 2012-5-31 20:36 编辑
.
贴段很早以前FORTRAN编的4节点矩形板单元显示形式单元刚度和质量矩阵程序,可以很容易看出矩阵具体表达式的:
C
C 计 算 当 前 单 元 组 各 单 元 的 单 元 矩 阵
C
EN = 2.1E11 ! 结构材料的弹性模量
PO = 0.3 ! 结构材料的泊松比
GN = EN/(2.*(1.+PO)) ! 结构材料的剪切模量
PN = 7800.0 ! 结构材料的质量密度
C
DO 100 N=1,NE
DX = XYZ(4,N)-XYZ(1,N)
DY = XYZ(5,N)-XYZ(2,N)
DZ = XYZ(6,N)-XYZ(3,N)
A = SQRT(DX*DX+DY*DY+DZ*DZ)/2.
DX = XYZ(7,N)-XYZ(4,N)
DY = XYZ(8,N)-XYZ(5,N)
DZ = XYZ(9,N)-XYZ(6,N)
B = SQRT(DX*DX+DY*DY+DZ*DZ)/2.
C
C To form element stiffness and mass in local coordinates
C
ES( 1, 1) = 21.-6.*PO+30.*B*B/A/A+30.*A*A/B/B
ES( 2, 1) = 3.*B+12.*PO*B+30.*A*A/B
ES( 3, 1) =-3.*A-12.*PO*A-30.*B*B/A
ES( 4, 1) =-21.+6.*PO-30.*B*B/A/A+15.*A*A/B/B
ES( 5, 1) =-3.*B-12.*PO*B+15.*A*A/B
ES( 6, 1) =-3.*A+3.*PO*A-30.*B*B/A
ES( 7, 1) = 21.-6.*PO-15.*B*B/A/A-15.*A*A/B/B
ES( 8, 1) =-3.*B+3.*PO*B+15.*A*A/B
ES( 9, 1) = 3.*A-3.*PO*A-15.*B*B/A
ES(10, 1) =-21.+6.*PO+15.*B*B/A/A-30.*A*A/B/B
ES(11, 1) = 3.*B-3.*PO*B+30.*A*A/B
ES(12, 1) = 3.*A+12.*PO*A-15.*B*B/A
ES( 2, 2) = 8.*B*B-8.*PO*B*B+40.*A*A
ES( 3, 2) =-30.*PO*A*B
ES( 4, 2) = ES(5,1)
ES( 5, 2) =-8.*B*B+8.*PO*B*B+20.*A*A
ES( 6, 2) = 0.
ES( 7, 2) =-ES(8,1)
ES( 8, 2) = 2.*B*B-2.*PO*B*B+10.*A*A
ES( 9, 2) = 0.
ES(10, 2) =-ES(11,1)
ES(11, 2) =-2.*B*B+2.*PO*B*B+20.*A*A
ES(12, 2) = 0.
ES( 3, 3) = 8.*A*A-8.*PO*A*A+40.*B*B
ES( 4, 3) =-ES(6,1)
ES( 5, 3) = 0.
ES( 6, 3) =-2.*A*A+2.*PO*A*A+20.*B*B
ES( 7, 3) =-ES(9,1)
ES( 8, 3) = 0.
ES( 9, 3) = 2.*A*A-2.*PO*A*A+10.*B*B
ES(10, 3) = ES(12,1)
ES(11, 3) = 0.
ES(12, 3) =-8.*A*A+8.*PO*A*A+20.*B*B
ES( 4, 4) = ES(1,1)
ES( 5, 4) = ES(2,1)
ES( 6, 4) =-ES(3,1)
ES( 7, 4) = ES(10,1)
ES( 8, 4) = ES(11,1)
ES( 9, 4) =-ES(12,1)
ES(10, 4) = ES( 7,1)
ES(11, 4) = ES( 8,1)
ES(12, 4) =-ES( 9,1)
ES( 5, 5) = ES( 2,2)
ES( 6, 5) =-ES( 3,2)
ES( 7, 5) = ES(10,2)
ES( 8, 5) = ES(11,2)
ES( 9, 5) = ES(12,2)
ES(10, 5) = ES( 7,2)
ES(11, 5) = ES( 8,2)
ES(12, 5) = ES( 9,2)
ES( 6, 6) = ES( 3,3)
ES( 7, 6) =-ES(10,3)
ES( 8, 6) = ES(11,3)
ES( 9, 6) = ES(12,3)
ES(10, 6) =-ES( 7,3)
ES(11, 6) = ES( 8,3)
ES(12, 6) = ES( 9,3)
ES( 7, 7) = ES( 1,1)
ES( 8, 7) =-ES( 2,1)
ES( 9, 7) =-ES( 3,1)
ES(10, 7) = ES( 4,1)
ES(11, 7) =-ES( 5,1)
ES(12, 7) =-ES( 6,1)
ES( 8, 8) = ES( 2,2)
ES( 9, 8) = ES( 3,2)
ES(10, 8) =-ES( 4,2)
ES(11, 8) = ES( 5,2)
ES(12, 8) = ES( 6,2)
ES( 9, 9) = ES( 3,3)
ES(10, 9) =-ES( 4,3)
ES(11, 9) = ES( 5,3)
ES(12, 9) = ES( 6,3)
ES(10,10) = ES( 1,1)
ES(11,10) =-ES( 2,1)
ES(12,10) = ES( 3,1)
ES(11,11) = ES( 2,2)
ES(12,11) =-ES( 3,2)
ES(12,12) = ES( 3,3)
C
EM( 1, 1) = 3454.
EM( 2, 1) = 922.*B
EM( 3, 1) =-922.*A
EM( 4, 1) = 1226.
EM( 5, 1) = 398.*B
EM( 6, 1) = 548.*A
EM( 7, 1) = 394.
EM( 8, 1) =-232.*B
EM( 9, 1) = 232.*A
EM(10, 1) = 1226.
EM(11, 1) =-548.*B
EM(12, 1) =-398.*A
EM( 2, 2) = 320.*B*B
EM( 3, 2) =-252.*A*B
EM( 4, 2) = 398.*B
EM( 5, 2) = 160.*B*B
EM( 6, 2) = 168.*A*B
EM( 7, 2) = 232.*B
EM( 8, 2) =-120.*B*B
EM( 9, 2) = 112.*A*B
EM(10, 2) = 548.*B
EM(11, 2) =-240.*B*B
EM(12, 2) =-168.*A*B
EM( 3, 3) = 320.*A*A
EM( 4, 3) =-548.*A
EM( 5, 3) =-168.*A*B
EM( 6, 3) =-240.*A*A
EM( 7, 3) =-232.*A
EM( 8, 3) = 112.*A*B
EM( 9, 3) =-120.*A*A
EM(10, 3) =-398.*A
EM(11, 3) = 168.*A*B
EM(12, 3) = 160.*A*A
EM( 4, 4) = 3454.
EM( 5, 4) = 922.*B
EM( 6, 4) = 922.*A
EM( 7, 4) = 1226.
EM( 8, 4) =-548.*B
EM( 9, 4) = 398.*A
EM(10, 4) = 394.
EM(11, 4) =-232*B
EM(12, 4) =-232*A
EM( 5, 5) = 320.*B*B
EM( 6, 5) = 252.*A*B
EM( 7, 5) = 548.*B
EM( 8, 5) =-240.*B*B
EM( 9, 5) = 168.*A*B
EM(10, 5) = 232.*B
EM(11, 5) =-120.*B*B
EM(12, 5) =-112.*A*B
EM( 6, 6) = 320.*A*A
EM( 7, 6) = 398.*A
EM( 8, 6) =-168.*A*B
EM( 9, 6) = 160.*A*A
EM(10, 6) = 232.*A
EM(11, 6) =-112.*A*B
EM(12, 6) =-120.*A*A
EM( 7, 7) = 3454.
EM( 8, 7) =-922.*B
EM( 9, 7) = 922.*A
EM(10, 7) = 1226.
EM(11, 7) =-398.*B
EM(12, 7) =-548.*A
EM( 8, 8) = 320.*B*B
EM( 9, 8) =-252.*A*B
EM(10, 8) =-398.*B
EM(11, 8) = 160.*B*B
EM(12, 8) = 168.*A*B
EM( 9, 9) = 320.*A*A
EM(10, 9) = 548.*A
EM(11, 9) =-168.*A*B
EM(12, 9) =-240.*A*A
EM(10,10) = 3454.
EM(11,10) =-922.*B
EM(12,10) =-922.*A
EM(11,11) = 320.*B*B
EM(12,11) = 252.*A*B
EM(12,12) = 320.*A*A
DO 1OO I=1,12
DO 1OO J=1,I
ES(I,J) = EN*TH(N)*TH(N)*TH(N)*ES(I,J)/(360.*(1.-PO*PO)*A*B)
EM(I,J) = 4.*PN*A*B*TH(N)*EM(I,J)/25200.
ES(J,I) = ES(I,J)
EM(J,I) = EM(I,J)
100CONTINUE
C
|