用於求解非線性方程的levenberg-marquardt方法

Question

我嘗試實現了基於Numerical Optimization using the Levenberg-Marquardt Algorithm演示文稿基於Julia求解非線性方程的levenberg-marquardt方法。這我的代碼：

function get_J(ArrOfFunc,X,delta) 
    N = length(ArrOfFunc) 
    J = zeros(Float64,N,N) 
    for i = 1:N 
    for j=1:N 
     Temp = copy(X); 
     Temp[j]=Temp[j]+delta; 
     J[i,j] = (ArrOfFunc[i](Temp)-ArrOfFunc[i](X))/delta; 
    end 
    end 
    return J 
end 

function get_resudial(ArrOfFunc,Arg) 
    return map((x)->x(Arg),ArrOfFunc) 
end 

function lm_solve(Funcs,Init) 
    X = copy(Init) 
    delta = 0.01; 
    Lambda = 0.01; 
    Factor = 2; 
    J = get_J(Funcs,X,delta) 
    R = get_resudial(Funcs,X) 
    N = 5 
    for t = 1:N 

    G = J'*J+Lambda.*eye(length(X)) 
    dC = J'*R 
    C = sum(R.*R)/2; 
    Xnew = X-(inv(G)\dC); 
    Rnew = get_resudial(Funcs,Xnew) 
    Cnew = sum(Rnew.*Rnew)/2; 
    if (Cnew < C) 
     X = Xnew; 
     R = Rnew; 
     Lambda = Lambda/Factor; 
     J = get_J(Funcs,X,delta) 
    else 
     Lambda = Lambda*Factor; 
    end 
    if(maximum(abs(Rnew)) < 0.001) 
     return X 
    end 
    end 
    return X 
end 

function test() 
    ArrOfFunc = [ 
    (X)->X[1]+X[2]-2; 
    (X)->X[1]-X[2] 
    ]; 

    X = lm_solve(ArrOfFunc,Float64[3;3]) 
    println(X) 
    return X 
end 

但是從任何出發點的步驟都不被接受。我做錯了什麼？ 任何幫助，將不勝感激。

Answer 1

我目前沒有辦法來測試這一點，但一行沒有意義數學：

在Xnew的計算應​​該是兩者無論是inv(G)*dC或G\dC，但不混合。由於線性系統的解決方案不需要逆矩陣的計算，所以優選第二個。

由於在迭代中心有一個錯誤的計算，計算的軌跡幾乎肯定會誤入歧途。