Блог

Я пытаюсь перевести этот Julia код, реализацию многомерного метода Ньютона-Рафсона Python , и расширить его, чтобы принять систему из четырех уравнений с четырьмя переменными.

Код:

 using SymPy
x, y = Sym("x, y")
function SymPyDerivatives(f::String, g::String)
    #
    # Returns the string representations of all
    # partial derivatives of the expression in x and y,
    # given in the strings f and g.
    #
    formf = parse(f)
    evalformf = eval(formf)
    fx = diff(evalformf, x)
    fy = diff(evalformf, y)
    formg = parse(g)
    evalformg = eval(formg)
    gx = diff(evalformg, x)
    gy = diff(evalformg, y)
    return [string(fx) string(fy); string(gx) string(gy)]
end

function NewtonStep(fun::Array{SymPy.Sym,1},
    jac::Array{SymPy.Sym,2},
    x0::Float64, y0::Float64)
    #
    # Runs one step with Newton’s method
    #
    valfun = -SymPyFun(fun, x0, y0)
    nfx = norm(valfun)
    valmat = SymPyMatrixEvaluate(jac, x0, y0)
    update = valmatvalfun
    ndx = norm(update)
    x1 = x0   update[1]
    y1 = y0   update[2]
    sfx = @sprintf("%.2e", nfx)
    sdx = @sprintf("%.2e", ndx)
    sx1 = @sprintf("%.16e", x1)
    sy1 = @sprintf("%.16e", y1)
    println(" $sfx $sdx $sx1 $sy1 ")
    return [x1, y1]
end

function main()
#
# Prompts the user for expressions in x and y,
# calls SymPy to take the derivatives with respect
# to x and y of the two given expressions.
#
    println("Reading two expressions, hit enter for examples")
    print("Give a first expression in x and y : ")
    one = readline(STDIN)
    if one == ""
        one = "exp(x) - y"
        two = "x*y - exp(x)"
        x0, y0 = 0.9, 2.5
    else
        print("Give a second expression in x and y : ")
        two = readline(STDIN)
        print("Give a value for x0 : ")
        x0 = parse(Float64(readline(STDIN)))
        print("Give a value for y0 : ")
        y0 = parse(Float64(readline(STDIN)))
    end
    derivatives = SymPyDerivatives(one, two)
    println("The Jacobian matrix :")
    println(derivatives)
    fx = derivatives[1,1]
    fy = derivatives[1,2]
    gx = derivatives[2,1]
    gy = derivatives[2,2]
    jacmat = SymPyMatrix(fx, fy, gx, gy)
    valmat = SymPyMatrixEvaluate(jacmat, x0, y0)
    vecfun = SymPyExpressions(one, two)
    for i=1:5
        xsol, ysol = NewtonStep(vecfun, jacmat, x0, y0)
        x0, y0 = xsol, ysol
    end
end

main()
  

Моя проблема в том, что я фактически никогда не использовал Julia и, хотя я довольно хорошо понимаю метод Ньютона, я не совсем понимаю этот код.

Вот что у меня есть на данный момент:

Внутри function NewtonStep()

Не уверен, что здесь происходит:
valfun = -SymPyFun(fun, x0, y0)

Мы принимаем норму этой матрицы
nfx = norm(valfun)

Не уверен, что здесь происходит:
valmat = SympPyMatrixEvaluate(jax,c, x0, y0)

Обновите матрицу:

 update = valmatvalfun
ndx = norm(update)
  

Не уверен, что происходит после этого. И почему function SymPyDerivatives() он никогда не используется? Конечно, это должно быть как-то.

Не совсем уверен, что происходит внутри function SymPyDerivatives() , за пределами строк, принимающих частные производные.

Обновить:

Просто пытаюсь разобрать это по частям, это то, что у меня есть до сих пор:

 from sympy import symbols, diff
import numpy as np

class Localize:

    receivers: list

    def jacobian(f, g, h, k):
        x, y, z = symbols('x y z', real=True)

        fx = diff(f, x);    gx = diff(g, x)
        fy = diff(f, y);    gy = diff(g, y)

        hx = diff(h, x);    kx = diff(k, x)
        hy = diff(h, y);    ky = diff(k, y)

    def newtonStep(x0, x1):

        # Magic happens

        x1 = x0   update[1]
        y1 = y0   update[2]
  

Кажется, не могу найти никакой документации SymPyFun() , что странно.

Редактировать:

Согласно моим собственным исследованиям и исследованиям @steamsy, казалось бы, этого SymPyFun() не существует, что странно. Возможно, у кого-то есть некоторое представление о том, для чего это было предназначено, и как это может выглядеть в Python ? Не совсем уверен, что он должен здесь делать.

Обновить:

Я смог найти function SymPyFun()

 function SymPyFun(fun::Array{SymPy.Sym,1},
                  xval::Float64,yval::Float64)
#
# Given an array of SymPy expressions,
# evaluates the expressions at (xval, yval).
#
    result1 = Float64(subs(subs(fun[1], x, xval), y, yval))
    result2 = Float64(subs(subs(fun[2], x, xval), y, yval))
    return [result1; result2]
end
  

Кажется, не могу найти никакой документации по SymPyFun(), что странно.

Потому что он не существует. В этой документации по SymPy.jl и моей обработке исходного кода GitHub никогда не упоминается SymPyFun.

Я предполагаю, что вы следуете вместе с этим PDF-файлом из UIC, это единственное, что я смог найти, в котором упоминается метод.