编程高手帮个忙，给我看看这个vb程序

啥时候睡觉 · 发表于 2006-1-12 22:17:59

<

>这程序是用来求解一个线性拟合，先用最小二乘法建立函数，再通过Gasuss-Newton法求解
<

>程序如下：
<

>Private Sub Command1_Click() Dim a As Double Dim re(10), r(10), s(3), x(3), x1(3), f(10) As Double Dim j(10, 3), j1(3, 10), c(3, 1), d(3, 3) As Double re(1) = 4658: re(23) = 5820: re(3) = 6525: re(4) = 7400: re(5) = 9045: re(6) = 10350: re(7) = 11050: re(8) = 11820: re(9) = 12850: re(10) = 13840 r(1) = 0.0399: r(2) = 0.0386: r(3) = 0.0368: r(4) = 0.0362: r(5) = 0.0349: r(6) = 0.0339: r(7) = 0.0341: r(8) = 0.0323: r(9) = 0.0324: r(10) = 0.0321 a = 0.001 For i = 1 To 3 x(i) = 0 s(i) = 0.01 Next i While ((Abs(s(1)) + Abs(s(2)) + Abs(s(3))) > 0.0001) For i = i To 3 s(i) = 0 Next i For i = 1 To 10 j(i, 1) = 2 * Exp(x(2) * Log(re(i))) * (x(1) * Exp(x(2) * Log(re(i))) + x(3) - r(i)) j(i, 2) = 2 * x(1) * Log(re(i)) * Exp(x(2) * Log(re(i))) * (x(1) * Exp(x(2) * Log(re(i))) + x(3) - r(i)) j(i, 3) = 2 * (x(1) * Exp(x(2) * Log(re(i))) + x(3) - r(i)) f(i) = (x(1) * Exp(x(2) * Log(re(i))) + x(3) - r(i)) ^ 2 Next i Call MTrans(10, 3, j(), j1()) Call MMul(3, 10, 1, j1(), f(), c()) Call MMul(3, 10, 3, j1(), j(), d()) Call MRinv(3, d()) Call MMul(3, 3, 1, d(), c(), s()) For i = 1 To 3 x1(i) = x(i) - a * s(i) x(i) = x1(i) Next i Wend For i = 1 To 3 Print x(i) Next i End Sub Function MRinv(n As Integer, mtxA() As Double) As Boolean '求逆矩阵 ReDim nIs(n) As Integer, nJs(n) As Integer Dim i As Integer, j As Integer, k As Integer Dim d As Double, p As Double
 For k = 1 To n d = 0# For i = k To n For j = k To n p = Abs(mtxA(i, j)) If (p > d) Then d = p nIs(k) = i nJs(k) = j End If Next j Next i If (d + 1# = 1#) Then MRinv = False Exit Function End If
 If (nIs(k) <> k) Then For j = 1 To n p = mtxA(k, j) mtxA(k, j) = mtxA(nIs(k), j) mtxA(nIs(k), j) = p Next j End If
 If (nJs(k) <> k) Then For i = 1 To n p = mtxA(i, k) mtxA(i, k) = mtxA(i, nJs(k)) mtxA(i, nJs(k)) = p Next i End If
 mtxA(k, k) = 1# / mtxA(k, k) For j = 1 To n If (j <> k) Then mtxA(k, j) = mtxA(k, j) * mtxA(k, k) Next j For i = 1 To n If (i <> k) Then For j = 1 To n If (j <> k) Then mtxA(i, j) = mtxA(i, j) - mtxA(i, k) * mtxA(k, j) Next j End If Next i For i = 1 To n If (i <> k) Then mtxA(i, k) = -mtxA(i, k) * mtxA(k, k) Next i Next k
 For k = n To 1 Step -1 If (nJs(k) <> k) Then For j = 1 To n p = mtxA(k, j) mtxA(k, j) = mtxA(nJs(k), j) mtxA(nJs(k), j) = p Next j End If If (nIs(k) <> k) Then For i = 1 To n p = mtxA(i, k) mtxA(i, k) = mtxA(i, nIs(k)) mtxA(i, nIs(k)) = p Next i End If Next k MRinv = True
End Function Sub MMul(m As Integer, n As Integer, l As Integer, mtxA() As Double, mtxB() As Double, mtxC() As Double) Dim i As Integer, j As Integer, k As Integer '矩阵乘法 For i = 1 To m For j = 1 To l mtxC(i, j) = 0# For k = 1 To n mtxC(i, j) = mtxC(i, j) + mtxA(i, k) * mtxB(k, j) Next k Next j Next i
End Sub Sub MTrans(m As Integer, n As Integer, mtxA() As Double, mtxAT() As Double) Dim i As Integer, j As Integer '矩阵转置
 For i = 1 To m For j = 1 To n mtxAT(i, j) = mtxAT(j, i) Next j Next i
End Sub 拜托各位拉！！！

		自动登录	找回密码
密码			注-册-帐-号