结构矩阵分析与程序设计钢架vb代码.docx

1、结构矩阵分析与程序设计钢架vb代码=Structural Analysis Program For Plane Frame=Option ExplicitPublic nn As Integer, ne As Integer, nd As Integer, ndf As IntegerPublic nf As Integer, npj As Integer, npe As Integer, n As IntegerPublic al(50) As Double, t(6, 6) As Double, x(40) As Double, y(40) As DoublePublic jl(50) A

2、s Integer, jr(50) As Integer, ea(50) As Double, ei(50) As DoublePublic c(6, 6) As Double, r(120, 120) As Double, p(120) As Double, pe(120) As DoublePublic ibd(20) As Integer, ii(6) As Integer, bd(20) As Double, ff(6) As DoublePublic mj(20) As Integer, qj(20, 3) As Double, f(6) As Double, dis(6) As D

3、oublePublic mf(30) As Integer, ind(30) As Double, aq(30) As Double, bq(30) As DoublePublic q1(30) As Double, q2(30) As Double=main program=Sub frame()Open H:juzheng钢架tr3.2.11.txt For Input As #1Open H:juzheng钢架tw3.2.11.txt For Output As #2Call input1Call wstiffCall loadCall boundCall gaussCall nqmCl

4、ose 1Close 2End Sub=SUB-1 Read And Print Intial Data=Sub input1()Dim inti As Integer, intj As Integer, i As Integer, j As Integer, k As IntegerDim dx, dy As Double Print #2, Plane Frame structural Analysis Print #2, * Print #2, input data Print #2, = = = = = Print #2, Print #2, structural control da

5、ta Print #2, - Print #2, nn; Spc(3); ne; Spc(3); nf; Spc(3); nd; Spc(3); ndf; Spc(2); npj; Spc(2); npe; Spc(3); n Input #1, nn, ne, nf, nd, ndf, npj, npe n = 3 * (nn - nf) Print #2, nn; Spc(2); ne; Spc(2); nf; Spc(2); nd; Spc(2); ndf; Spc(2); npj; Spc(2); npe; Spc(2); n Print #2, Print #2, Nodal coo

6、rdinates Print #2, - Print #2, Node ; Spc(2); x; Spc(5); yi = nn For inti = 1 To i Input #1, inti, x(inti), y(inti) Print #2, inti; Spc(2); x(inti); Spc(3); y(inti) Next intiPrint #2,Print #2, Element InformationPrint #2, -Print #2, Ele.No.; Spc(4); ; jl; Spc(4); jr; Spc(6); ea; Spc(6); ei; Spc(6);

7、ali = ne For inti = 1 To i Input #1, inti, jl(inti), jr(inti), ea(inti), ei(inti) Next intiFor inti = 1 To iIf jl(inti) = jr(inti) Then StopNext intiFor inti = 1 To i j = jl(inti) k = jr(inti) dx = x(k) - x(j) dy = y(k) - y(j) al(inti) = Sqr(dx * dx + dy * dy) Print #2, Spc(3); inti; Spc(4); jl(inti

8、); Spc(3); jr(inti); Spc(2); ea(inti); Spc(2); ei(inti); Spc(2); al(inti)Next intiPrint #2,k = npj If k 0 Then Print #2, Nodal Load Print #2, - Print #2, i; Spc(13); mj; Spc(3); xd; Spc(2); yd; Spc(2); md For inti = 1 To k Input #1, inti, mj(inti), qj(inti, 1), qj(inti, 2), qj(inti, 3) Print #2, int

9、i; Spc(1), mj(inti); Spc(1); qj(inti, 1); Spc(1); qj(inti, 2); Spc(1); qj(inti, 3) Next inti End If Print #2, i = npe If i 0 Then Print #2, Element loads Print #2, - Print #2, i; Spc(5); mf; Spc(3); ind; Spc(3); aq; Spc(3); bq; Spc(3); q1; Spc(4); q2 For inti = 1 To i Input #1, inti, mf(inti), ind(i

10、nti), aq(inti), bq(inti), q1(inti), q2(inti) Print #2, inti; Spc(2); mf(inti); Spc(3); ind(inti); Spc(2); aq(inti); Spc(2); bq(inti); Spc(2); q1(inti); Spc(3); q2(inti) Next intiEnd IfPrint #2, j = ndf If j 0 Then Print #2, Bonundary conditions Print #2, - Print #2, i; Spc(5); ibd; Spc(3); bdFor int

11、i = 1 To j Input #1, inti, ibd(inti), bd(inti) Print #2, inti; Spc(3); ibd(inti); Spc(3); bd(inti) Next inti End IfEnd Sub=sub-2 Assemnble Structural Stiffness MatrixR=Sub wstiff()Dim i As Integer, j As Integer, ie As Integer, k1 As Integer, k2 As IntegerFor i = 1 To n For j = 1 To n r(i, j) = 0 Nex

12、t jNext iie = 1Do While ie = ne Call stiff(ie) Call locat(ie)For k1 = 1 To 6 i = ii(k1) If i = n Then For k2 = k1 To 6 j = ii(k2) If j = n Then r(i, j) = r(i, j) + c(k1, k2) End If Next k2 End If Next k1 ie = ie + 1LoopFor i = 2 To n For j = 1 To (i - 1) r(i, j) = r(j, i) Next jNext iEnd Sub =sub-3

13、set up Stiffness Matrixc=Sub stiff(ie)Dim i As Integer, j As IntegerDim cx As Double, cy As Double, b1 As Double, b2 As Double, b3 As Double, b4 As DoubleDim s1 As Double, s2 As Double, s3 As Double, s4 As Double, s5 As Double, s6 As Doublei = jl(ie)j = jr(ie)cx = (x(j) - x(i) / al(ie)cy = (y(j) - y

14、(i) / al(ie)b1 = ea(ie) / al(ie)b2 = 12# * ei(ie) / al(ie) 3b3 = 6# * ei(ie) / al(ie) 2b4 = 2# * ei(ie) / al(ie)s1 = b1 * cx 2 + b2 * cy 2s2 = (b1 - b2) * cx * cys3 = b3 * cys4 = b1 * cy 2 + b2 * cx 2s5 = b3 * cxs6 = b4c(1, 1) = s1c(1, 2) = s2c(1, 3) = s3c(1, 4) = -s1c(1, 5) = -s2c(1, 6) = s3c(2, 2)

15、 = s4c(2, 3) = -s5c(2, 4) = -s2c(2, 5) = -s4c(2, 6) = -s5c(3, 3) = 2# * s6c(3, 4) = -s3c(3, 5) = s5c(3, 6) = s6c(4, 4) = s1c(4, 5) = s2c(4, 6) = -s3c(5, 5) = s4c(5, 6) = s5c(6, 6) = 2# * s6For i = 2 To 6 For j = 1 To (i - 1) c(i, j) = c(j, i) Next jNext iEnd Sub=sub-4 set up element location vector=

16、Sub locat(ie)Dim i As Integer, j As Integeri = jl(ie)j = jr(ie)ii(1) = 3 * i - 2ii(2) = 3 * i - 1ii(3) = 3 * iii(4) = 3 * j - 2ii(5) = 3 * j - 1ii(6) = 3 * jEnd Sub=SUB-5 Set Up Total Nodal Vector P= Sub load()Dim i As Integer, j As Integer, k As Integeri = 1Do While i 0 Then For i = 1 To npj k = mj

17、(i) p(3 * k - 2) = qj(i, 1) p(3 * k - 1) = qj(i, 2)p(3 * k) = qj(i, 3) Next iEnd IfIf npe 0 Then Call eload i = 1 Do While i = np(i) = p(i) + pe(i) i = i + 1LoopEnd IfEnd Sub=SUB-6 Set up element effective load= Sub eload()Dim i As Integer, j As Integer, k As Integer, k1 As Integer, k2 As Integer, k

18、3 As Integeri = 1 Do While i = n pe(i) = 0# i = i + 1Loopj = 1Do While j = npe k = mf(j) Call trans(k) Call locat(k) Call efix(j) For k1 = 1 To 6 f(k1) = 0# For k2 = 1 To 6 f(k1) = f(k1) + t(k2, k1) * ff(k2) Next k2Next k1For k3 = 1 To 6i = ii(k3) If i = n Then pe(i) = pe(i) - f(k3) End IfNext k3j =

19、 j + 1LoopEnd Sub =sub-7 set up fixed-end force of element=Sub efix(i)Dim j As Integer, k As IntegerDim s1 As Double, a As Double, b As Double, p1 As Double, p2 As DoubleDim b1 As Double, b2 As Double, b3 As Double, c1 As Double, c2 As Double, c3 As DoubleDim d1 As Double, d2 As DoubleFor j = 1 To 6

20、ff(j) = 0#Next jk = mf(i)s1 = al(k)a = aq(i)b = bq(i)p1 = q1(i)p2 = q2(i)b1 = s1 - (a + b) / 2#b2 = b - ab3 = (a + b) / 2#c1 = s1 - (2# * b + a) / 3#c2 = b2c3 = (2# * b + a) / 3#d1 = b 3 - a 3d2 = b * b - a * aSelect Case ind(i) Case 1 ff(2) = -p1 * (s1 - a) 2 * (1# + 2# * a / s1) / s1 2 ff(3) = p1

21、* a * (s1 - a) 2 / s1 2 ff(5) = -p1 - ff(2) ff(6) = -p1 * a 2 * (s1 - a) / s1 2 Case 2 ff(2) = -p1 * b2 * (12# * b1 2 * s1 - 8# * b1 3 + b2 2 * s1 - 2# * b1 * b2 2) / (4# * s1 3) ff(3) = p1 * b2 * (12# * b3 * b1 2 - 3 * b1 * b2 2 + b2 2 * s1) / 12# / s1 2 ff(5) = -p1 * b2 - ff(2) ff(6) = -p1 * b2 *

22、(12# * b3 2 * b1 + 3# * b1 * b2 2 - 2# * b2 2 * s1) / 12# / s1 2 Case 3 ff(2) = -p2 * c2 * (18 * c1 2 * s1 - 12 * c1 3 + c2 2 * s1 - 2 * c1 * c2 2 - 4 * c2 3 / 45) / 12 / s1 3 ff(3) = p2 * c2 * (18# * c3 * c1 2 - 3# * c1 * c2 2 + c2 2 * s1 - 2# * c2 3 / 15#) / 36# / s1 2 ff(5) = -0.5 * p1 * c2 - ff(

23、2) ff(6) = -p2 * c2 * (18# * c3 2 * c1 + 3 * c1 * c2 2 - 2 * c2 2 * s1 + 2 * c2 3 / 15#) / 36# / s1 2 Case 4 ff(2) = -6# * p1 * a * (s1 - a) / s1 3 ff(3) = p1 * (s1 - a) * (3# * a - s1) / s1 2 ff(5) = -ff(2) ff(6) = p1 * a * (2# * s1 - 3# * a) / s1 2 Case 5 ff(2) = -p1 * (3# * s1 * d2 - 2# * d1) / s

24、1 3 ff(3) = p1 * (2# * d2 + (b - a) * s1 - d1 / s1) / s1 ff(5) = -ff(2) ff(6) = p1 * (d2 - d1 / s1) / s1 Case 6 ff(1) = -p1 * (1# - a / s1) ff(4) = -p1 * a / s1 Case 7 ff(1) = -p1 * (b - a) * (1# - (b + a) / (2# * s1) ff(4) = -p1 * d2 / 2# / s1 Case 8 ff(3) = -a * (p1 - p2) * ei(k) / b ff(6) = -ff(3

25、)End SelectEnd Sub=sub-8 set up coordinate transfer matrixt=Sub trans(ie)Dim i As Integer, j As IntegerDim cx As Double, cy As Doublei = jl(ie)j = jr(ie)cx = (x(j) - x(i) / al(ie)cy = (y(j) - y(i) / al(ie)For i = 1 To 6 For j = 1 To 6 t(i, j) = 0# Next jNext iFor i = 1 To 4 Step 3 t(i, i) = cx t(i, i + 1) = cy t(i + 1, i) = -cy t(i + 1, i + 1) = cx t(i + 2, i + 2) = 1#Next iEnd Sub=sub-9 introduce support conditions= Sub bound()Dim i As Integer, j As Integer, k As IntegerDim a As Double If ndf 0 Then For j = 1 To ndf a = 1E+20 For i = 1 To ndf k = ibd(i) r(k, k) = a p(k) = a