python实现最优化算法：负梯度法（最速下降法）

一时失言乱红尘 2022-12-30 15:56 231阅读 0赞

啥也别说了，最优化太心酸了，可能是每个计算机专业人的噩梦（主要是数学渣）。直接上代码。

# 向前查找搜索区间，接受线搜索方向的导数f为唯一参数，返回搜索区间左右端点left,right
    def find_interval(f):
        left = 0
        step = 1
        coefficient = 2
        right = left + step
        
        while f(right) < 0:
            left = right
            step = step*coefficient
            right = left + step
        
        return left, right
    
    
    # 实现用二分法精确线搜索，接受线搜索方向的导数和搜索区间左右端点为三个参数，返回
    # 绝对值不超过10**-6的近似极小点mid
    def find_root(f, left,right):
        limit = 1e-6
        mid = (left + right)/2
        while abs(f(mid)) >= limit:
    
            if f(mid) < 0:
                left = mid
            elif f(mid) >= 0:
                right = mid
            mid = (left + right)/2
    
        return mid
    
    # 接收待优化函数的梯度g和初始值x(numpy数组类型)为两个参数，返回梯度向量长度
    # （norm(g(x)）不超过1e-5的近似极小值点（调用norm 函数计算向量长度，测试时
    # 可以使用numpy库中的linalg模块的norm函数实现）
    
    def grad_min(g, x):
        import numpy
        while numpy.linalg.norm(g(x)) > 1e-5:
            d=-g(x)
            d=-g(x)/numpy.linalg.norm(g(x))
            
            h = lambda y:numpy.dot(g(x+y*d),d) # d搜索方向，h为线搜索方向的导数
            x = x+find_root(h,find_interval(h)[0], find_interval(h)[1])*d
    
        return x
    
    import numpy
    from sympy import lambdify, symbols, Matrix, MatrixSymbol, diff
    
    G = MatrixSymbol('G', 2, 2)
    b = MatrixSymbol('b', 2, 1)
    c = MatrixSymbol('c', 1, 1)
    ps = [c, b,G]
    xs = symbols('x1:3')
    x = Matrix(xs)
    quad = (x.T*G*x/2+b.T*x+c)[0,0]
    expr=lambda y:quad.subs([(p, Matrix(y.pop())) for p in ps]).simplify()
    fexpr=expr([[[21,4],[4,1]],[2,3],[10]])
    gexpr=[diff(fexpr,x) for x in xs]
    glambdify=[lambdify([xs],f,"numpy") for f in gexpr]
    g = lambda y:numpy.array([f(y) for f in glambdify]) # 计算梯度
    
    
    x = numpy.array([1,1])
    grad_min(g,x)
    
    # 运行结果 array([ 1.99999231, -10.9999606 ])