vb 算法问题

Question

要解决方程组求解的问题，有多个约束条件。 
方程组可能存在多个解，不求最优化解， 
实现所有解输出，用什么算法？（因为有多个约束条件，不太可能有太多的解）。 
最好有例子！ 求解的过程中，涉及的变量很多，而且变量、约束条件的个数是不定的。 
是线性方程组。 类似的例子：
0 <a1x1+a2x2+...+anxn <1000 
-1 <b1x1+b2x2+...+bnxn <1002 
                . 
                . 
                . 
      x1> 10 
      0 <x2 <1 
            . 
            . 
            .

Answer 1

以下是求解非线性方程和方程组的例子Sub Main()
    Dim b As Double
    Dim m As Integer
    Dim x As Double
    
    ' 迭代初值
    x = 0.5
    
    ' 随机数初始值与控制参数
    b = 1#
    m = 10
    
    
    '求解
    Call NLMtclRoot(x, b, m, 0.00001)
    
    MsgBox "x = " & x
    
End Sub''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  功能：  求解非线性方程和方程组
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Option Explicit''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：NLBisectRoot
'  功能：  使用对分法求解非线性方程的实根，本函数需要调用计算方程左端函数f(x)值的函数Func，其形式为：
'          Function Func(x As Double) As Double
'  参数    m    - Integer型变量，在[a, b]内实根个数的预估值
'          a    - Double型变量，求根区间的左端点
'          b    - Double型变量，求根区间的右端点
'          h    - Double型变量，搜索求根时采用的步长
'          x    - Double型一维数组，长度为m。返回在区间[a, b]内搜索到的实根，实根个数由函数值返回
'         eps   - Double型变量，精度控制参数
'  返回值：Integer型，求得的实根的个数
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function NLBisectRoot(m As Integer, a As Double, b As Double, h As Double, x() As Double, eps As Double) As Integer
    Dim n As Integer, js As Integer
    Dim z As Double, y As Double, z1 As Double, y1 As Double, z0 As Double, y0 As Double    ' 根的个数清0
    n = 0
    
    ' 左边界函数值
    z = a
    y = Func(z)
    
    ' 迭代求解，直到到达右边界
    While ((z <= b + h / 2#) And (n <> m))
        ' 如果精度满足要求，则求得一个实根，继续计算下一步
        If (Abs(y) < eps) Then
            n = n + 1
            x(n) = z
            z = z + h / 2#
            y = Func(z)
        Else
            z1 = z + h
            y1 = Func(z1)
            If (Abs(y1) < eps) Then
                n = n + 1
                x(n) = z1
                z = z1 + h / 2#
                y = Func(z)
            Else
                If (y * y1 > 0#) Then
                    y = y1
                    z = z1
                Else
                    js = 0
                    While (js = 0)
                        If (Abs(z1 - z) < eps) Then
                            n = n + 1
                            x(n) = (z1 + z) / 2#
                            z = z1 + h / 2#
                            y = Func(z)
                            js = 1
                        Else
                            z0 = (z1 + z) / 2#
                            y0 = Func(z0)
                            If (Abs(y0) < eps) Then
                                x(n) = z0
                                n = n + 1
                                js = 1
                                z = z0 + h / 2#
                                y = Func(z)
                            Else
                                If ((y * y0) < 0#) Then
                                    z1 = z0
                                    y1 = y0
                                Else
                                    z = z0
                                    y = y0
                                End If
                            End If
                        End If
                    Wend
                End If
            End If
        End If
    Wend    ' 返回求得的根的个数
    NLBisectRoot = nEnd Function''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：NLNewtonRoot
'  功能：  使用牛顿迭代法求解非线性方程的实根，本函数需要调用计算方程左端函数f(x)值的过程Func，其形式为：
'          Sub Func(x As Double, y() as double)
'          y(1) 返回f(x)的值
'          y(2) 返回f'(x)的值
'  参数    x    - Double型变量，输入时存放迭代初值；返回时存放迭代终值，即方程的根
'         js    - Integer型变量，最大迭代次数
'         eps   - Double型变量，精度控制参数
'  返回值：Integer型，若小于0，则表示在求解失败；若等于最大迭代次数js，则表示迭代了js次还未满足精度要求，
'         返回的实根只作为参考；若大于等于0且小于最大迭代次数js，则表示正常返回。
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function NLNewtonRoot(x As Double, js As Integer, eps As Double) As Integer
    Dim k As Integer, l As Integer
    Dim y(2) As Double, d As Double, p As Double, x0 As Double, x1 As Double    ' 初值
    l = js
    k = 1
    x0 = x
    
    ' 计算f(x)和f'(x)
    Call Func(x0, y)
    
    d = eps + 1#
    
    ' 迭代计算
    While ((d >= eps) And (l <> 0))
        ' 求解失败，返回
        If (Abs(y(2)) + 1# = 1#) Then
            NLNewtonRoot = -1
            Exit Function
        End If        x1 = x0 - y(1) / y(2)        Call Func(x1, y)        d = Abs(x1 - x0)
        p = Abs(y(1))        If (p > d) Then d = p
        
        x0 = x1
        l = l - 1
    Wend    ' 求解结束
    x = x1
    k = js - l    NLNewtonRoot = kEnd Function''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：NLAitkenRoot
'  功能：  使用埃特金迭代法求解非线性方程的一个实根，本函数需要调用计算方程左端函数f(x)值的过程Func，其形式为：
'          Function Func(x as double) as double
'  参数    x    - Double型变量，输入时存放迭代初值；返回时存放迭代终值，即方程的根
'         js    - Integer型变量，最大迭代次数
'         eps   - Double型变量，精度控制参数
'  返回值：Integer型，若为0，则表示迭代了js次还为满足精度要求，返回的实根只作为参考；
'                     若大于0，则表示正常返回。
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function NLAitkenRoot(x As Double, js As Integer, eps As Double) As Integer
    Dim flag As Integer, l As Integer
    Dim u As Double, v As Double, x0 As Double    ' 迭代初值
    l = 0
    x0 = x
    flag = 0
    
    ' 迭代求解
    While ((flag = 0) And (l <> js))
        l = l + 1
        u = Func(x0)
        v = Func(u)
        If (Abs(u - v) < eps) Then
            x0 = v
            flag = 1
        Else
            x0 = v - (v - u) * (v - u) / (v - 2# * u + x0)
        End If
    Wend    ' 求解结束
    x = x0
    l = js - l    NLAitkenRoot = lEnd Function

Answer 2

''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：NLAitkenRoot
'  功能：  使用埃特金迭代法求解非线性方程的一个实根，本函数需要调用计算方程左端函数f(x)值的过程Func，其形式为：
'          Function Func(x as double) as double
'  参数    x    - Double型变量，输入时存放迭代初值；返回时存放迭代终值，即方程的根
'         js    - Integer型变量，最大迭代次数
'         eps   - Double型变量，精度控制参数
'  返回值：Integer型，若为0，则表示迭代了js次还为满足精度要求，返回的实根只作为参考；
'                     若大于0，则表示正常返回。
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function NLAitkenRoot(x As Double, js As Integer, eps As Double) As Integer
    Dim flag As Integer, l As Integer
    Dim u As Double, v As Double, x0 As Double    ' 迭代初值
    l = 0
    x0 = x
    flag = 0
    
    ' 迭代求解
    While ((flag = 0) And (l <> js))
        l = l + 1
        u = Func(x0)
        v = Func(u)
        If (Abs(u - v) < eps) Then
            x0 = v
            flag = 1
        Else
            x0 = v - (v - u) * (v - u) / (v - 2# * u + x0)
        End If
    Wend    ' 求解结束
    x = x0
    l = js - l    NLAitkenRoot = lEnd Function''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：NLPqRoot
'  功能：  使用连分式解法求解非线性方程的一个实根，本函数需要调用计算方程左端函数f(x)值的过程Func，其形式为：
'          Function Func(x as double) as double
'  参数    x    - Double型变量，输入时存放迭代初值；返回时存放迭代终值，即方程的根
'         eps   - Double型变量，精度控制参数
'  返回值：Boolean型，若为False，则表示迭代了10次还未满足精度要求，返回的实根只作为参考；
'                     若为True，则表示正常返回。
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function NLPqRoot(x As Double, eps As Double) As Boolean
    Dim i As Integer, j As Integer, m As Integer, it As Integer, l As Integer
    Dim a(10) As Double, y(10) As Double, z As Double, h As Double, x0 As Double, q As Double    ' 迭代初值
    l = 10
    q = 1E+35
    x0 = x
    h = 0#
    
    ' 迭代计算
    While (l <> 0)
        l = l - 1
        j = 0
        it = l
        
        ' 最多迭代到第7节连分式
        While (j <= 7)
             If (j <= 2) Then
                z = x0 + 0.1 * j
             Else
                z = h
             End If             ' 函数值
             y(j + 1) = Func(z)
             h = z
             If (j = 0) Then
                a(1) = z
             Else
                 m = 0
                 i = 0
                 While ((m = 0) And (i <= j - 1))
                     If (Abs(h - a(i + 1)) + 1# = 1#) Then
                        m = 1
                     Else
                        h = (y(j + 1) - y(i + 1)) / (h - a(i + 1))
                     End If                     i = i + 1
                 Wend
                 a(j + 1) = h
                 If (m <> 0) Then a(j + 1) = q
                 h = 0#
                 For i = j To 1 Step -1
                     If (Abs(a(i + 1) + h) + 1# = 1#) Then
                        h = q
                     Else
                        h = -y(i) / (a(i + 1) + h)
                     End If
                 Next i
                 h = h + a(1)
             End If             ' 精度达到要求，则结束迭代
             If (Abs(y(j + 1)) >= eps) Then
                j = j + 1
             Else
                j = 10
                l = 0
             End If         Wend
         
         x0 = h
    
    Wend    ' 求解结束
    x = h
    
    ' 判断解的合理性
    If it = 0 Then
        NLPqRoot = False
        Exit Function
    End If
     
    ' 正常解
    NLPqRoot = True
     
End Function''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：NLQrRoot
'  功能：  使用QR方法求解实系数代数方程的全部根，本函数需要调用用QR方法计算上－H阵全部特征值的函数
'          MhbergEigenv
'  参数    n    - 多项式方程的次数
'          a    - Double型一维数组，长度为n+1，按降幂次序依次存放n次多项式方程的n+1个系数
'          xr    - Double型一维数组，长度为n，返回n个根的实部
'          xi    - Double型一维数组，长度为n，返回n个根的虚部
'         eps   - Double型变量，精度控制参数
'         nMaxItNum    - Integer型变量，控制QR方法的最大迭代次数
'  返回值：Boolean型，若为False，则表示在QR方法中迭代已超过最大迭代次数但还未满足精度要求；
'                     若为True，则表示求解成功。
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function NLQrRoot(n As Integer, a() As Double, xr() As Double, xi() As Double, eps As Double, nMaxItNum As Integer) As Boolean
    Dim i As Integer, j As Integer
    ReDim q(n, n) As Double    ' 用最高幂系数约化其余系数
    For j = 1 To n
      q(1, j) = -a(j + 1) / a(1)
    Next j    ' 构造上－H阵
    For i = 2 To n
        For j = 1 To n
            q(i, j) = 0#
        Next j
    Next i    For i = 2 To n
        q(i, i - 1) = 1#
    Next i    ' 求上－H阵的特征值，即方程的全部根
    NLQrRoot = MHbergEigenv(n, q, xr, xi, eps, nMaxItNum)End Function

Answer 3

''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：NLNdhRoot
'  功能：  使用牛顿－下山法求解实系数代数方程的全部根
'  参数    n    - 多项式方程的次数
'          a    - Double型一维数组，长度为n+1，按降幂次序依次存放n次多项式方程的n+1个系数
'          xr    - Double型一维数组，长度为n，返回n个根的实部
'          xi    - Double型一维数组，长度为n，返回n个根的虚部
'  返回值：Boolean型，若为False，则表示在QR方法中迭代已超过最大迭代次数但还未满足精度要求；
'                     若为True，则表示求解成功。
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function NLNdhRoot(n As Integer, a() As Double, xr() As Double, xi() As Double) As Boolean
    Dim m As Integer, i As Integer, jt As Integer, k As Integer, nIs As Integer, it As Integer
    Dim t As Double, x As Double, y As Double, x1 As Double, y1 As Double, dx As Double, dy As Double
    Dim p As Double, q As Double, w As Double, dd As Double, dc As Double, c As Double
    Dim g As Double, u As Double, v As Double, pq As Double, g1 As Double, u1 As Double, v1 As Double    ' 系数个数
    m = n + 1
    
    ' 用最高幂系数约化其余系数
    k = a(1)
    For i = 1 To m
        a(i) = a(i) / k
    Next i
    
    ' 迭代求解
    k = m
    nIs = 0
    w = 1#
    jt = 1
    While (jt = 1)
        pq = Abs(a(k))
        While (pq < 0.000000000001)
            xr(k - 1) = 0#
            xi(k - 1) = 0#
            k = k - 1
            
            ' 求解结束
            If (k = 2) Then
                xr(1) = -a(2) * w / a(1)
                xi(1) = 0#
                NLNdhRoot = True
                Exit Function
            End If            pq = Abs(a(k))
        Wend        q = Log(pq)
        q = q / (1# * k)
        q = Exp(q)
        p = q
        w = w * p
        For i = 2 To k
            a(i) = a(i) / q
            q = q * p
        Next i        x = 0.0001
        x1 = x
        y = 0.2
        y1 = y
        dx = 1#
        g = 1E+37
l40:
        u = a(1)
        v = 0#
        For i = 2 To k
            p = u * x1
            q = v * y1
            pq = (u + v) * (x1 + y1)
            u = p - q + a(i)
            v = pq - p - q
        Next i        g1 = u * u + v * v
        If (g1 >= g) Then
            If (nIs <> 0) Then
                it = 1
                Call g65(x, y, x1, y1, dx, dy, dd, dc, c, k, nIs, it)
                If (it = 0) Then GoTo l40
            Else
                Call g60(t, x, y, x1, y1, dx, dy, p, q, k, it)
                If (t >= 0.001) Then GoTo l40
                If (g > 1E-18) Then
                    it = 0
                    Call g65(x, y, x1, y1, dx, dy, dd, dc, c, k, nIs, it)
                    If (it = 0) Then GoTo l40
                End If
            End If
            Call g90(xr, xi, a, x, y, p, q, w, k)
        Else
            g = g1
            x = x1
            y = y1
            nIs = 0
            If (g <= 1E-22) Then
                Call g90(xr, xi, a, x, y, p, q, w, k)
            Else
                u1 = k * a(1)
                v1 = 0#
                For i = 3 To k
                    p = u1 * x
                    q = v1 * y
                    pq = (u1 + v1) * (x + y)
                    u1 = p - q + (k - i + 1) * a(i - 1)
                    v1 = pq - p - q
                Next i                p = u1 * u1 + v1 * v1
                If (p <= 1E-20) Then
                    it = 0
                    Call g65(x, y, x1, y1, dx, dy, dd, dc, c, k, nIs, it)
                    If (it = 0) Then GoTo l40
                    Call g90(xr, xi, a, x, y, p, q, w, k)
                Else
                    dx = (u * u1 + v * v1) / p
                    dy = (u1 * v - v1 * u) / p
                    t = 1# + 4# / k
                    Call g60(t, x, y, x1, y1, dx, dy, p, q, k, it)
                    If (t >= 0.001) Then GoTo l40
                    If (g > 1E-18) Then
                        it = 0
                        Call g65(x, y, x1, y1, dx, dy, dd, dc, c, k, nIs, it)
                        If (it = 0) Then GoTo l40
                    End If
                    Call g90(xr, xi, a, x, y, p, q, w, k)
                End If
            End If
        End If        If (k = 2) Then
            jt = 0
        Else
            jt = 1
        End If
    Wend
    
    ' 迭代结束
    NLNdhRoot = True
   
End Function''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：g60
'  功能：  内部函数
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Sub g60(t As Double, x As Double, y As Double, x1 As Double, y1 As Double, dx As Double, dy As Double, p As Double, q As Double, k As Integer, it As Integer)
    
    it = 1
    While (it = 1)
        t = t / 1.67
        it = 0
        x1 = x - t * dx
        y1 = y - t * dy
        If (k > 50) Then
            p = Sqr(x1 * x1 + y1 * y1)
            q = Exp(85# / k)
            If (p >= q) Then it = 1
        End If
    WendEnd Sub''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：g90
'  功能：  内部函数
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Sub g90(xr() As Double, xi() As Double, a() As Double, x As Double, y As Double, p As Double, q As Double, w As Double, k As Integer)
    
    Dim i As Integer    If (Abs(y) <= 0.000001) Then
        p = -x
        y = 0#
        q = 0#
    Else
        p = -2# * x
        q = x * x + y * y
        xr(k - 1) = x * w
        xi(k - 1) = -y * w
        k = k - 1
    End If
    
    For i = 2 To k
        a(i) = a(i) - a(i - 1) * p
        a(i + 1) = a(i + 1) - a(i - 1) * q
    Next i
    
    xr(k - 1) = x * w
    xi(k - 1) = y * w
    k = k - 1
    
    If (k = 2) Then
        xr(1) = -a(2) * w / a(1)
        xi(1) = 0#
    End If
      
 End Sub

Answer 4

''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：g65
'  功能：  内部函数
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Sub g65(x As Double, y As Double, x1 As Double, y1 As Double, dx As Double, dy As Double, dd As Double, dc As Double, c As Double, k As Integer, nIs As Integer, it As Integer)
  
    If (it = 0) Then
        nIs = 1
        dd = Sqr(dx * dx + dy * dy)
        If (dd > 1#) Then dd = 1#
        dc = 6.28 / (4.5 * k)
        c = 0#
    End If    While (True)
        c = c + dc
        dx = dd * Cos(c)
        dy = dd * Sin(c)
        x1 = x + dx
        y1 = y + dy
        If (c <= 6.29) Then
            it = 0
            Exit Sub
        End If
        dd = dd / 1.67
        If (dd <= 0.0000001) Then
            it = 1
            Exit Sub
        End If
        c = 0#
    WendEnd Sub''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：NLNdhcRoot
'  功能：  使用牛顿－下山法求解复系数代数方程的全部根
'  参数    n    - 多项式方程的次数
'          ar    - Double型一维数组，长度为n+1，按降幂次序依次存放n次多项式方程的n+1个系数的实部
'          ai    - Double型一维数组，长度为n+1，按降幂次序依次存放n次多项式方程的n+1个系数的虚部
'          xr    - Double型一维数组，长度为n，返回n个根的实部
'          xi    - Double型一维数组，长度为n，返回n个根的虚部
'  返回值：Boolean型，若为False，则表示在QR方法中迭代已超过最大迭代次数但还未满足精度要求；
'                     若为True，则表示求解成功。
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function NLNdhcRoot(n As Integer, ar() As Double, ai() As Double, xr() As Double, xi() As Double) As Boolean
    Dim m As Integer, i As Integer, jt As Integer, k As Integer, nIs As Integer, it As Integer
    Dim t As Double, x As Double, y As Double, x1 As Double, y1 As Double, dx As Double, dy As Double
    Dim p As Double, q As Double, w As Double, dd As Double, dc As Double, c As Double
    Dim g As Double, u As Double, v As Double, pq As Double, g1 As Double, u1 As Double, v1 As Double    ' 系数个数
    m = n + 1
    
    ' 用最高幂系数约化其余系数
    p = Sqr(ar(1) * ar(1) + ai(1) * ai(1))
    For i = 1 To m
        ar(i) = ar(i) / p
        ai(i) = ai(i) / p
    Next i
    
    ' 迭代求解
    k = m
    nIs = 0
    w = 1#
    jt = 1
    
    While (jt = 1)
        pq = Sqr(ar(k) * ar(k) + ai(k) * ai(k))
        
        While (pq < 0.000000000001)
            xr(k - 1) = 0#
            xi(k - 1) = 0#
            k = k - 1
            
            ' 求解结束
            If (k = 2) Then
                p = ar(1) * ar(1) + ai(1) * ai(1)
                xr(1) = -w * (ar(1) * ar(2) + ai(1) * ai(2)) / p
                xi(1) = w * (ar(2) * ai(1) - ar(1) * ai(2)) / p
                NLNdhcRoot = True
                Exit Function
            End If
            
            pq = Sqr(ar(k) * ar(k) + ai(k) * ai(k))
        Wend
        
        q = Log(pq)
        q = q / (1# * k)
        q = Exp(q)
        p = q
        w = w * p
        For i = 2 To k
            ar(i) = ar(i) / q
            ai(i) = ai(i) / q
            q = q * p
        Next i        x = 0.0001
        x1 = x
        y = 0.2
        y1 = y
        dx = 1#
        g = 1E+37
l40:
        u = ar(1)
        v = ai(1)
        For i = 2 To k
            p = u * x1
            q = v * y1
            pq = (u + v) * (x1 + y1)
            u = p - q + ar(i)
            v = pq - p - q + ai(i)
        Next i        g1 = u * u + v * v
        If (g1 >= g) Then
            If (nIs <> 0) Then
                it = 1
                Call g65c(x, y, x1, y1, dx, dy, dd, dc, c, k, nIs, it)
                If (it = 0) Then GoTo l40
            Else
                Call g60c(t, x, y, x1, y1, dx, dy, p, q, k, it)
                If (t >= 0.001) Then GoTo l40
                If (g > 1E-18) Then
                    it = 0
                    Call g65c(x, y, x1, y1, dx, dy, dd, dc, c, k, nIs, it)
                    If (it = 0) Then GoTo l40
                End If
            End If            Call g90c(xr, xi, ar, ai, x, y, p, w, k)
        Else
            g = g1
            x = x1
            y = y1
            nIs = 0
            If (g <= 1E-22) Then
                Call g90c(xr, xi, ar, ai, x, y, p, w, k)
            Else
                u1 = k * ar(1)
                v1 = ai(1)
                For i = 3 To k
                    p = u1 * x
                    q = v1 * y
                    pq = (u1 + v1) * (x + y)
                    u1 = p - q + (k - i + 1) * ar(i - 1)
                    v1 = pq - p - q + (k - i + 1) * ai(i - 1)
                Next i                p = u1 * u1 + v1 * v1
                If (p <= 1E-20) Then
                    it = 0
                    Call g65c(x, y, x1, y1, dx, dy, dd, dc, c, k, nIs, it)
                    If (it = 0) Then GoTo l40
                    Call g90c(xr, xi, ar, ai, x, y, p, w, k)
                Else
                    dx = (u * u1 + v * v1) / p
                    dy = (u1 * v - v1 * u) / p
                    t = 1# + 4# / k
                    Call g60c(t, x, y, x1, y1, dx, dy, p, q, k, it)
                    If (t >= 0.001) Then GoTo l40
                    If (g > 1E-18) Then
                        it = 0
                        Call g65c(x, y, x1, y1, dx, dy, dd, dc, c, k, nIs, it)
                        If (it = 0) Then GoTo l40
                    End If
                    Call g90c(xr, xi, ar, ai, x, y, p, w, k)
                End If
            End If
        End If        If (k = 2) Then
            jt = 0
        Else
            jt = 1
        End If
    Wend    ' 迭代结束
    NLNdhcRoot = True
   
End Function''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：g60c
'  功能：  内部函数
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Sub g60c(t As Double, x As Double, y As Double, x1 As Double, y1 As Double, dx As Double, dy As Double, p As Double, q As Double, k As Integer, it As Integer)
    
    it = 1
    While (it = 1)
        t = t / 1.67
        it = 0
        x1 = x - t * dx
        y1 = y - t * dy
        If (k >= 30) Then
            p = Sqr(x1 * x1 + y1 * y1)
            q = Exp(75# / k)
            If (p >= q) Then it = 1
        End If
    WendEnd Sub

Answer 5

''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：g90c
'  功能：  内部函数
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Sub g90c(xr() As Double, xi() As Double, ar() As Double, ai() As Double, x As Double, y As Double, p As Double, w As Double, k As Integer)
    
    Dim i As Integer    For i = 2 To k
        ar(i) = ar(i) + ar(i - 1) * x - ai(i - 1) * y
        ai(i) = ai(i) + ar(i - 1) * y + ai(i - 1) * x
    Next i    xr(k - 1) = x * w
    xi(k - 1) = y * w
    k = k - 1
    If (k = 2) Then
        p = ar(1) * ar(1) + ai(1) * ai(1)
        xr(1) = -w * (ar(1) * ar(2) + ai(1) * ai(2)) / p
        xi(1) = w * (ar(2) * ai(1) - ar(1) * ai(2)) / p
    End If
    
End Sub''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：g65c
'  功能：  内部函数
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Sub g65c(x As Double, y As Double, x1 As Double, y1 As Double, dx As Double, dy As Double, dd As Double, dc As Double, c As Double, k As Integer, nIs As Integer, it As Integer)
  
    If (it = 0) Then
        nIs = 1
        dd = Sqr(dx * dx + dy * dy)
        If (dd > 1#) Then dd = 1#
        dc = 6.28 / (4.5 * k)
        c = 0#
    End If    While (True)
        c = c + dc
        dx = dd * Cos(c)
        dy = dd * Sin(c)
        x1 = x + dx
        y1 = y + dy
        If (c <= 6.29) Then
            it = 0
            Exit Sub
        End If
        dd = dd / 1.67
        If (dd <= 0.0000001) Then
            it = 1
            Exit Sub
        End If
        c = 0#
    WendEnd Sub''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：NLGrad
'  功能：  使用梯度法求解非线性方程组的一组解，本函数需要调用计算目标函数值的函数Func，其形式为：
'          Function Func(x() as double) as double
'  参数    n    - Integer型变量，方程的个数，也是未知数的个数
'          x    - Double型一维数组，长度为n，存放一组初值x1, x2, …, xn，返回时存放方程组的一组实根
'         nMaxIt - Integer型变量，最大迭代次数
'         eps   - Double型变量，精度控制参数
'          h    - Double型变量，收敛控制参数，一般取1e-5至1e-10之间
'  返回值：Boolean型，False则求解失败，True则求解成功。
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function NLGrad(n As Integer, x() As Double, nMaxIt As Integer, eps As Double, h As Double) As Boolean
    Dim i As Integer, j As Integer
    Dim f As Double, fh As Double, r As Double, sum As Double
    ReDim hx(n) As Double, df(n) As Double    For j = 1 To nMaxIt
        For i = 1 To n
            If (x(i) = 0) Then
                hx(i) = h
            Else
                hx(i) = h * x(i)
            End If
        Next i
    
        f = Func(x)
        
        If (f <= eps) Then
            NLGrad = True
            Exit Function
        End If
        
        sum = 0
        
        For i = 1 To n
            x(i) = x(i) + hx(i)
            fh = Func(x)
            df(i) = (fh - f) / hx(i)
            sum = sum + df(i) * df(i)
            x(i) = x(i) - hx(i)
        Next i
        
        r = f / sum
    
        For i = 1 To n
            x(i) = x(i) - r * df(i)
        Next i
        
    Next j
    
    NLGrad = FalseEnd Function''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：NLNewtonA
'  功能：  使用拟牛顿法求非线性方程组
'         1）本函数要调用全选主元高斯消去法求解线性代数方程组的函数LEGauss;
'         2）本函数需要调用计算目标函数值的函数Func，其形式为：
'          Sub Func(x() as Double, y() as Double)
'          y(i) 为方程组中各左边多项式与一组x(i)对应的函数值
'  参数    n    - Integer型变量，方程的个数，也是未知数的个数
'          x    - Double型一维数组，长度为n，存放一组初值x1, x2, …, xn，返回时存放方程组的一组实根
'         nMaxIt - Integer型变量，最大迭代次数
'         eps   - Double型变量，精度控制参数
'          h    - Double型变量，增量初值，在本函数中要被破坏
'          t    - Double型变量，控制h大小的变量，0<t<1
'  返回值：Boolean型，True则求解成功；False则求解失败。
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function NLNewtonA(n As Integer, x() As Double, nMaxIt As Integer, eps As Double, h As Double, t As Double) As Boolean
    Dim i As Integer, j As Integer, l As Integer
    Dim am As Double, z As Double, beta As Double, d As Double
    ReDim y(n) As Double, a(n, n) As Double, b(n) As Double    ' 最大迭代次数
    l = nMaxIt
    
    '精度控制
    am = 1# + eps
    
    ' 迭代求解
    While (am >= eps)
        
        ' 函数值
        Call Func(x, b)
        
        am = 0#
        For i = 1 To n
            z = Abs(b(i))
            If (z > am) Then am = z
        Next i        If (am >= eps) Then
            l = l - 1
            
            ' 达到最大迭代次数，精度未达到要求，求解失败，返回
            If (l = 0) Then
                NLNewtonA = False
                Exit Function
            End If            For j = 1 To n
                z = x(j)
                x(j) = x(j) + h
                
                Call Func(x, y)
                
                For i = 1 To n
                    a(i, j) = y(i)
                Next i
                
                x(j) = z
                
            Next j            ' 高斯消元失败，求解失败，返回
            If Not LEGauss(n, a, b) Then
                NLNewtonA = False
                Exit Function
            End If            beta = 1#
            For i = 1 To n
                beta = beta - b(i)
            Next i            ' 方程求解失败，返回
            If (Abs(beta) + 1# = 1#) Then
                NLNewtonA = False
                Exit Function
            End If            d = h / beta
            For i = 1 To n
                x(i) = x(i) - d * b(i)
            Next i            h = t * h
            
       End If
    Wend    ' 方程求解成功，返回
    NLNewtonA = TrueEnd Function

Answer 6

''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：NLMiv
'  功能：  用最小二乘解的广义逆法求非线性方程组
'          1）本函数要调用求解线性最小乘问题的广义逆法的函数LEMiv；
'          2）本函数需要调用计算目标函数值的函数Func，其形式为：
'             Sub Func(x() as Double, y() as Double)
'             y(i) 为方程组中各左边多项式与一组x(i)对应的函数值
'          3）本函数需要调用计算雅可比矩阵函数，其形式为：
'             Sub FuncMJ(x() as Double, p() as Double)
'             p(i, j) 为方程组中各左边多项式fi的对变元xj的偏导数与一组x(i)对应的函数值
'  参数    m    - Integer型变量，非线性方程组中方程的个数
'          n    - Integer型变量，非线性方程组中未知数的个数
'          x    - Double型一维数组，长度为n，存放一组初值x1, x2, …, xn，要求不全为0，返回时存放方程组的
'                 最小二乘解，当m=n时，即是非线性方程组的解
'         eps1   - Double型变量，最小二乘解的精度控制参数
'         eps2   - Double型变量，奇异值分解的精度控制参数
'         ka     - Integer型变量，ka=max(m, n)+1
'  返回值：Boolean型，True则求解成功；False则求解失败。
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Function NLMiv(m As Integer, n As Integer, x() As Double, eps1 As Double, eps2 As Double, ka As Integer) As Boolean
    Dim i As Integer, j As Integer, k As Integer, l As Integer, kk As Integer
    Dim jt As Boolean
    Dim y(10) As Double, b(10) As Double, alpha As Double, z As Double, h2 As Double, y1 As Double, y2 As Double, y3 As Double, y0 As Double, h1 As Double
    
    ' 分配内存
    ReDim p(m, n) As Double, d(m) As Double, pp(n, m) As Double, dx(n) As Double, u(m, m) As Double, v(n, n) As Double, w(ka) As Double    ' 迭代次数设为60
    l = 60
    
    ' 控制初值
    alpha = 1#
    
    While (l > 0)
        ' 调用计算目标函数值的函数
        Call Func(x, d)
        
        ' 调用计算雅可比矩阵函数
        Call FuncMJ(x, p)
        
        ' 调用求解线性最小乘问题的广义逆法的函数LEMiv
        jt = LEMiv(m, n, p, d, dx, pp, u, v, ka, eps2)
        
        ' 求解失败，返回
        If (jt = False) Then
            NLMiv = False
            Exit Function
        End If        ' 继续迭代计算
        j = 0
        jt = True
        h2 = 0#
        
        While (jt)
            jt = False
            If (j <= 2) Then
                z = alpha + 0.01 * j
            Else
                z = h2
            End If            For i = 1 To n
                w(i) = x(i) - z * dx(i)
            Next i            Call Func(w, d)
            
            y1 = 0#
            
            For i = 1 To m
                y1 = y1 + d(i) * d(i)
            Next i            For i = 1 To n
              w(i) = x(i) - (z + 0.00001) * dx(i)
            Next i            Call Func(w, d)
            
            y2 = 0#
            For i = 1 To m
                y2 = y2 + d(i) * d(i)
            Next i            y0 = (y2 - y1) / 0.00001
            
            If (Abs(y0) > 0.0000000001) Then
                h1 = y0
                h2 = z
                If (j = 0) Then
                    y(1) = h1
                    b(1) = h2
                Else
                    y(j + 1) = h1
                    kk = 0
                    k = 0
                    While ((kk = 0) And (k <= j - 1))
                        y3 = h2 - b(k + 1)
                        If (Abs(y3) + 1# = 1#) Then
                            kk = 1
                        Else
                            h2 = (h1 - y(k + 1)) / y3
                        End If                        k = k + 1
                    Wend                    b(j + 1) = h2
                    If (kk <> 0) Then b(j + 1) = 1E+35
                    h2 = 0#
                    For k = j - 1 To 0 Step -1
                      h2 = -y(k + 1) / (b(k + 2) + h2)
                    Next k                    h2 = h2 + b(1)
                End If
                j = j + 1
                If (j <= 7) Then
                    jt = True
                Else
                    z = h2
                End If
           End If
        Wend        alpha = z
        y1 = 0#
        y2 = 0#
        For i = 1 To n
            dx(i) = -alpha * dx(i)
            x(i) = x(i) + dx(i)
            y1 = y1 + Abs(dx(i))
            y2 = y2 + Abs(x(i))
        Next i        ' 达到精度要求，求解成功，返回
        If (y1 < eps1 * y2) Then
            NLMiv = True
            Exit Function
        End If        l = l - 1
    Wend    ' 迭代60次后仍未达到精度要求，求解失败，返回
    NLMiv = FalseEnd Function''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：NLMtclRoot
'  功能：  用蒙特卡洛求非线性方程一个实根，本函数需要调用计算方程左端函数f(x)值的函数Func，其形式为：
'          Function Func(x As Double) As Double
'  参数    x    - Double型变量，存放初值，返回时存放方程的一个实根
'          b    - Double型变量，均匀分布随机数的端点初值
'          m   - Integer型变量，控制调节b的参数
'         eps   - Double型变量，精度控制参数
'  返回值：无
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Sub NLMtclRoot(x As Double, b As Double, m As Integer, eps As Double)
    Dim k As Integer
    Dim xx As Double, a As Double, y As Double, x1 As Double, y1 As Double    ' 初值
    a = b
    k = 1
    xx = x
    
    ' 函数值
    y = Func(xx)
    
    ' 迭代求解
    While (a >= eps)
        
        ' 随机数
        Call Randomize(1)
        x1 = Rnd()
        
        x1 = -a + 2# * a * x1
        x1 = xx + x1
        y1 = Func(x1)
        
        k = k + 1
        
        If (Abs(y1) >= Abs(y)) Then
            If (k > m) Then
                k = 1
                a = a / 2#
            End If
        Else
            k = 1
            xx = x1
            y = y1
            
            ' 精度达到要求，求解结束
            If (Abs(y) < eps) Then
                x = xx
                Exit Sub
            End If
        End If
    Wend    ' 迭代求解结束
    x = xx
    
End Sub

Answer 7

''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：NLMtclcRoot
'  功能：  用蒙特卡洛求非线性方程一个复根，本函数需要调用计算方程左端函数f(x)值的函数Func，其形式为：
'          Function Func(x As Double, y As Double) As Double
'          其返回值为||f(x+jy)||
'  参数    x    - Double型变量，存放初值的实部，返回时存放方程的一个复根的实部
'          y    - Double型变量，存放初值的虚部，返回时存放方程的一个复根的虚部
'          b    - Double型变量，均匀分布随机数的端点初值
'          m   - Integer型变量，控制调节b的参数
'         eps   - Double型变量，精度控制参数
'  返回值：无
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Sub NLMtclcRoot(x As Double, y As Double, b As Double, m As Integer, eps As Double)
    Dim k As Integer
    Dim xx As Double, yy As Double, a As Double, z As Double
    Dim x1 As Double, y1 As Double, z1 As Double    ' 初值
    a = b
    k = 1
    xx = x
    yy = y    ' 函数值
    z = Func(xx, yy)    ' 迭代求解
    While (a >= eps)
        
        ' 随机数
        Call Randomize(1)
        x1 = -a + 2# * a * Rnd()        x1 = xx + x1        ' 随机数
        Call Randomize(1)
        y1 = -a + 2# * a * Rnd()        y1 = yy + y1
        z1 = Func(x1, y1)
        k = k + 1
        If (z1 >= z) Then
            If (k > m) Then
                k = 1
                a = a / 2#
            End If
        Else
            k = 1
            xx = x1
            yy = y1
            z = z1
            
            ' 精度达到要求，求解结束
            If (z < eps) Then
                x = xx
                y = yy
                Exit Sub
            End If
        End If
    Wend    ' 迭代求解结束
    x = xx
    y = yyEnd Sub
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
'  模块名：NLModule.bas
'  函数名：NLMtcl
'  功能：  用蒙特卡洛法求非线性方程组一组实根，本函数需要调用计算方程左端函数f(x)值的函数Func，其形式为：
'          Function Func(x As Double) As Double
'          其返回值为Sqr(f1*f1 + f2 * f2 + … + fn*fn)
'  参数    n   - Integer型变量，方程的个数，也是未知数的个数
'          x    - Double型一维数组，长度为n，存放一组初值，返回时存放方程的一组实根
'          b    - Double型变量，均匀分布随机数的端点初值
'          m   - Integer型变量，控制调节b的参数
'         eps   - Double型变量，精度控制参数
'  返回值：无
''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''''
Sub NLMtcl(n As Integer, x() As Double, b As Double, m As Integer, eps As Double)
    Dim k As Integer, i As Integer
    Dim a As Double, z As Double, z1 As Double
    ReDim y(n) As Double    ' 初值
    a = b
    k = 1
    
    ' 函数值
    z = Func(x)    ' 迭代求解
    While (a >= eps)
        For i = 1 To n
             ' 随机数
            Call Randomize(1)
            y(i) = -a + 2# * a * Rnd() + x(i)
        Next i        z1 = Func(y)
        
        k = k + 1
        
        If (z1 >= z) Then
            If (k > m) Then
                k = 1
                a = a / 2#
            End If
        Else
            k = 1
            For i = 1 To n
                x(i) = y(i)
            Next i            z = z1
            
            ' 精度达到要求，求解结束
            If (z < eps) Then Exit Sub
        
        End If
    Wend
End Sub

Answer 8

大哥 小弟新人 请问一下有求伪逆的函数吧（A矩阵秩亏求他的逆 ） 我现在编的一个程序  就因为这个函数  我只能和matlab混编了 调用了它函数库中的pinv  编译貌似很复杂 离了我这台机子就没法运行了 很不爽 [email protected]