Нелинейное моделирование

В JuMP поддерживаются нелинейные (выпуклые и невыпуклые) задачи оптимизации. JuMP может автоматически предоставлять решателям точные разреженные производные второго порядка. Эта информация может повысить точность и производительность решателя.

Задание нелинейной целевой функции

Чтобы задать нелинейную целевую функцию, используйте макрос @objective.

julia> model = Model();

julia> @variable(model, x[1:2]);

julia> @objective(model, Min, exp(x[1]) - sqrt(x[2]))
exp(x[1]) - sqrt(x[2])

Чтобы изменить нелинейную целевую функцию, еще раз вызовите @objective.

Добавление нелинейного ограничения

Для добавления нелинейного ограничения используйте макрос @constraint.

julia> model = Model();

julia> @variable(model, x[1:2]);

julia> @constraint(model, exp(x[1]) <= 1)
exp(x[1]) - 1.0 ≤ 0

julia> @constraint(model, con[i = 1:2], 2^x[i] >= i)
2-element Vector{ConstraintRef{Model, MathOptInterface.ConstraintIndex{MathOptInterface.ScalarNonlinearFunction, MathOptInterface.GreaterThan{Float64}}, ScalarShape}}:
 con[1] : (2.0 ^ x[1]) - 1.0 ≥ 0
 con[2] : (2.0 ^ x[2]) - 2.0 ≥ 0

Чтобы удалить нелинейное ограничение, используйте метод delete:

julia> delete(model, con[1])

Добавление параметра

Некоторые решатели явным образом поддерживают параметры, представляющие собой константы в модели, которые можно эффективно обновлять между решениями.

Параметры реализуются в JuMP посредством переменной решения, ограниченной при создании множеством Parameter.

julia> model = Model();

julia> @variable(model, x);

julia> @variable(model, p[i = 1:2] in Parameter(i))
2-element Vector{VariableRef}:
 p[1]
 p[2]

julia> parameter_value(p[1])
1.0

julia> set_parameter_value(p[1], 3.5)

julia> @objective(model, Max, log(p[1] * x + p[2]))
log(p[1]*x + p[2])

Дополнительные сведения о создании параметров и управлении ими см. в разделе Параметры.

Параметры наиболее полезны при решении нелинейных моделей в следующем порядке.

julia> using JuMP, Ipopt


julia> model = Model(Ipopt.Optimizer);


julia> set_silent(model)


julia> @variable(model, x)
x

julia> @variable(model, p in Parameter(1.0))
p

julia> @objective(model, Min, (x - p)^2)
x² - 2 p*x + p²

julia> optimize!(model)


julia> value(x)
1.0

julia> set_parameter_value(p, 5.0)


julia> optimize!(model)


julia> value(x)
5.0

Использование параметров может быть быстрее, чем создание модели с нуля с обновленными данными, так как позволяет JuMP избежать повторения ряда шагов при обработке модели перед ее передачей решателю.

Создание нелинейного выражения

Для создания объектов нелинейных выражений используйте макрос @expression.

julia> model = Model();

julia> @variable(model, x[1:2]);

julia> expr = @expression(model, exp(x[1]) + sqrt(x[2]))
exp(x[1]) + sqrt(x[2])

julia> my_anon_expr = @expression(model, [i = 1:2], sin(x[i]))
2-element Vector{NonlinearExpr}:
 sin(x[1])
 sin(x[2])

julia> @expression(model, my_expr[i = 1:2], sin(x[i]))
2-element Vector{NonlinearExpr}:
 sin(x[1])
 sin(x[2])

Выражение NonlinearExpr можно использовать в @objective, @constraint и даже вкладывать в другие выражения @expression.

julia> @objective(model, Min, expr^2 + 1)
((exp(x[1]) + sqrt(x[2])) ^ 2.0) + 1.0

julia> @constraint(model, [i = 1:2], my_expr[i] <= i)
2-element Vector{ConstraintRef{Model, MathOptInterface.ConstraintIndex{MathOptInterface.ScalarNonlinearFunction, MathOptInterface.LessThan{Float64}}, ScalarShape}}:
 sin(x[1]) - 1.0 ≤ 0
 sin(x[2]) - 2.0 ≤ 0

julia> @expression(model, nested[i = 1:2], sin(my_expr[i]))
2-element Vector{NonlinearExpr}:
 sin(sin(x[1]))
 sin(sin(x[2]))

Чтобы запросить значение нелинейного выражения, используйте метод value:

julia> set_start_value(x[1], 1.0)

julia> value(start_value, nested[1])
0.7456241416655579

julia> sin(sin(1.0))
0.7456241416655579

Автоматическое дифференцирование

JuMP вычисляет производные первого и второго порядка, используя разреженное автоматическое дифференцирование в обратном режиме. Подробные сведения см. в разделе ReverseAD.

Руководство по созданию и запросу производных см. в разделе Computing Hessians.

Подробные сведения о нелинейных выражениях

Нелинейные выражения представлены в JuMP объектом NonlinearExpr.

Конструкторы

Нелинейные выражения можно создавать с помощью конструкторов NonlinearExpr:

julia> model = Model();

julia> @variable(model, x);

julia> expr = NonlinearExpr(:sin, Any[x])
sin(x)

или путем перегрузки операторов:

julia> model = Model();

julia> @variable(model, x);

julia> expr = sin(x)
sin(x)

Поддерживаемые аргументы

Нелинейные выражения могут содержать сочетание чисел, выражений AffExpr, QuadExpr и других выражений NonlinearExpr:

julia> model = Model();

julia> @variable(model, x);

julia> aff = x + 1;

julia> quad = x^2 + x;

julia> expr = cos(x) * sin(quad) + aff
(cos(x) * sin(x² + x)) + (x + 1)

Поддерживаемые операторы

Список поддерживаемых операторов зависит от решателя. Для определенного оптимизатора список поддерживаемых операторов можно запросить с помощью MOI.ListOfSupportedNonlinearOperators:

julia> import Ipopt

julia> import MathOptInterface as MOI

julia> MOI.get(Ipopt.Optimizer(), MOI.ListOfSupportedNonlinearOperators())
85-element Vector{Symbol}:
 :+
 :-
 :abs
 :sqrt
 :cbrt
 :abs2
 :inv
 :log
 :log10
 :log2
 ⋮
 :min
 :max
 :&&
 :||
 :<=
 :(==)
 :>=
 :<
 :>

В некоторых одномерных случаях оператор определен в SpecialFunctions.jl. Для использования этих функций необходимо импортировать SpecialFunctions.jl явным образом.

julia> import Ipopt

julia> op = MOI.get(Ipopt.Optimizer(), MOI.ListOfSupportedNonlinearOperators());

julia> :erfcx in op
true

julia> :dawson in op
true

julia> import SpecialFunctions

julia> model = Model();

julia> @variable(model, x)
x

julia> @expression(model, SpecialFunctions.erfcx(x))
erfcx(x)

julia> @expression(model, SpecialFunctions.dawson(x))
dawson(x)

Ограничения

Некоторые нелинейные выражения невозможно создать путем перегрузки операторов. Например, чтобы свести к минимуму вероятность ошибок в пользовательском коде, мы не стали перегружать операторы сравнения объектов JuMP, такие как < и >=:

julia> model = Model();

julia> @variable(model, x);

julia> x < 1
ERROR: Cannot evaluate `<` between a variable and a number.
[...]

Вместо этого заключите выражение в макрос @expression:

julia> model = Model();

julia> @variable(model, x);

julia> expr = @expression(model, x < 1)
x < 1

По техническим причинам в число операторов, не допускающих перегрузки, также включены ||, && и ifelse.

julia> model = Model();

julia> @variable(model, x);

julia> expr = @expression(model, ifelse(x < -1 || x >= 1, x^2, 0.0))
ifelse((x < -1) || (x >= 1), x², 0.0)

В качестве альтернативы используйте функции JuMP.op_, которые по умолчанию используют различные операторы сравнения и логические операторы:

julia> model = Model();

julia> @variable(model, x);

julia> expr = op_ifelse(
           op_or(op_strictly_less_than(x, -1), op_greater_than_or_equal_to(x, 1)),
           x^2,
           0.0,
       )
ifelse((x < -1) || (x >= 1), x², 0.0)

Доступные функции:

Функция JuMP Функция Julia

Функция JuMP	Функция Julia
`op_ifelse`	`ifelse`
`op_and`	`&&`
`op_or`	`\|\|`
`op_greater_than_or_equal_to`	`>=`
`op_less_than_or_equal_to`	`<=`
`op_equal_to`	`==`
`op_strictly_greater_than`	`>`
`op_strictly_less_than`	`<`

op_ifelse

ifelse

op_and

&&

op_or

||

op_greater_than_or_equal_to

>=

op_less_than_or_equal_to

<=

op_equal_to

==

op_strictly_greater_than

>

op_strictly_less_than

<

Поля

У каждого выражения NonlinearExpr есть два поля.

Поле .head — это символ Symbol, представляющий вызываемый оператор:

julia> expr.head
:sin

Поле .args — это вектор Vector{Any}, содержащий аргументы оператора:

julia> expr.args
1-element Vector{Any}:
 x

Принудительное создание нелинейных выражений

При использовании макросов JuMP и перегрузки операторов в первую очередь создаются аффинные (GenericAffExpr) и квадратичные (GenericQuadExpr) выражения, а не GenericNonlinearExpr. Например:

julia> model = Model();

julia> @variable(model, x);

julia> f = (x - 0.1)^2
x² - 0.2 x + 0.010000000000000002

julia> typeof(f)
QuadExpr (alias for GenericQuadExpr{Float64, GenericVariableRef{Float64}})

Чтобы переопределить это поведение, используйте макрос @force_nonlinear:

julia> g = @force_nonlinear((x - 0.1)^2)
(x - 0.1) ^ 2

julia> typeof(g)
NonlinearExpr (alias for GenericNonlinearExpr{GenericVariableRef{Float64}})

Используйте этот макрос только при необходимости. Дополнительные сведения о том, когда его следует использовать, см. в docstring метода @force_nonlinear.

Трассировка функций

Нелинейные выражения можно создавать путем трассировки функций. Трассировка функции заключается в вызове обычной функции Julia с переменными JuMP в качестве аргументов. При этом функция создает нелинейное выражение посредством перегрузки операторов. Например:

julia> using JuMP


julia> model = Model();


julia> @variable(model, x[1:2]);


julia> f(x::Vector{VariableRef}) = 2 * sin(x[1]^2) + sqrt(x[2])
f (generic function with 1 method)

julia> y = f(x)
(2.0 * sin(x[1]²)) + sqrt(x[2])

julia> typeof(y)
NonlinearExpr (alias for GenericNonlinearExpr{GenericVariableRef{Float64}})

julia> @objective(model, Max, f(x))
(2.0 * sin(x[1]²)) + sqrt(x[2])

Трассировка функций поддерживает функции, которые возвращают векторы или массивы выражений NonlinearExpr:

julia> using JuMP


julia> model = Model();


julia> @variable(model, x[1:2]);


julia> f(x::Vector{VariableRef}) = sqrt.(x)
f (generic function with 1 method)

julia> y = f(x)
2-element Vector{NonlinearExpr}:
 sqrt(x[1])
 sqrt(x[2])

julia> typeof(y)
Vector{NonlinearExpr} (alias for Array{GenericNonlinearExpr{GenericVariableRef{Float64}}, 1})

julia> @constraint(model, f(x) .<= 2)
2-element Vector{ConstraintRef{Model, MathOptInterface.ConstraintIndex{MathOptInterface.ScalarNonlinearFunction, MathOptInterface.LessThan{Float64}}, ScalarShape}}:
 sqrt(x[1]) - 2.0 ≤ 0
 sqrt(x[2]) - 2.0 ≤ 0

julia> @objective(model, Max, sum(f(x)))
0.0 + sqrt(x[2]) + sqrt(x[1])

Так как при трассировке функций применяется перегрузка операторов, для многих функций она непригодна. Например:

julia> using JuMP

julia> model = Model();

julia> @variable(model, x[1:2]);

julia> f(x::Vector{VariableRef}) = x[1] > 1 ? 0 : x[2]
f (generic function with 1 method)

julia> f(x)
ERROR: Cannot evaluate `>` between a variable and a number.
[...]

В таких случаях следует определить пользовательский оператор с помощью макроса @operator.

Пользовательские операторы

Помимо стандартного списка одномерных и многомерных операторов, распознаваемых подмодулем MOI.Nonlinear, JuMP поддерживает пользовательские операторы, позволяющие представлять нелинейные функции, которые не могут (или не должны) трассироваться, например, потому, что они опираются на подпрограммы, не относящиеся к Julia.

Пользовательские операторы должны возвращать скалярный результат. Обходное решение см. в разделе User-defined operators with vector outputs.

Добавление оператора

Добавить пользовательский оператор можно с помощью макроса @operator:

julia> using JuMP


julia> square(x) = x^2
square (generic function with 1 method)

julia> f(x, y) = (x - 1)^2 + (y - 2)^2
f (generic function with 1 method)

julia> model = Model();


julia> @operator(model, op_square, 1, square)
NonlinearOperator(square, :op_square)

julia> @operator(model, op_f, 2, f)
NonlinearOperator(f, :op_f)

julia> @variable(model, x[1:2]);


julia> @objective(model, Min, op_f(x[1], op_square(x[2])))
op_f(x[1], op_square(x[2]))

Макрос @operator имеет следующие аргументы:

Модель, в которую добавляется оператор.
Объект символа Julia, который служит именем пользовательского оператора

в выражениях JuMP. Имя не должно совпадать с именем функции.

Количество принимаемых функцией скалярных входных аргументов.
Метод Julia, вычисляющий функцию.

Пользовательские операторы нельзя удалить.

Получить ссылку на оператор можно посредством синтаксиса model[:key]:

julia> using JuMP


julia> square(x) = x^2
square (generic function with 1 method)

julia> model = Model();


julia> @operator(model, op_square, 1, square)
NonlinearOperator(square, :op_square)

julia> op_square_2 = model[:op_square]
NonlinearOperator(square, :op_square)

Автоматическое дифференцирование

JuMP вычисляет производные первого и второго порядка для выражений, используя подмодуль ReverseAD, который реализует разреженное автоматическое дифференцирование в обратном режиме. Однако, поскольку подмодулю ReverseAD требуется алгебраическое выражение в качестве входных данных, JuMP не может использовать ReverseAD для дифференцирования пользовательских операторов.

Вместо этого, если явно не предоставлены градиенты и гессианы (Gradients and Hessians), пользовательские операторы должны поддерживать автоматическое дифференцирование с помощью ForwardDiff.jl.

Использование FowardDiff.jl имеет два важных следствия:

ForwardDiff.jl поддерживает только ограниченное подмножество операторов Julia. Если при добавлении

оператора происходит ошибка, см. раздел Common mistakes when writing a user-defined operator.

Дифференцирование операторов со множеством аргументов происходит медленно. В общем случае

следует стараться, чтобы количество аргументов было меньше 100, а в идеале — меньше 10.

Из-за использования ForwardDiff в большинстве случаев следует предпочесть трассировку функций вместо определения пользовательского оператора.

Добавление оператора без макроса

Макрос @operator — это синтаксический сахар для add_nonlinear_operator. Версия предыдущего примера без макроса выглядит следующим образом:

julia> using JuMP


julia> square(x) = x^2
square (generic function with 1 method)

julia> f(x, y) = (x - 1)^2 + (y - 2)^2
f (generic function with 1 method)

julia> model = Model();


julia> op_square = add_nonlinear_operator(model, 1, square; name = :op_square)
NonlinearOperator(square, :op_square)

julia> model[:op_square] = op_square
NonlinearOperator(square, :op_square)

julia> op_f = add_nonlinear_operator(model, 2, f; name = :op_f)
NonlinearOperator(f, :op_f)

julia> model[:op_f] = op_f
NonlinearOperator(f, :op_f)

julia> @variable(model, x[1:2]);


julia> @objective(model, Min, op_f(x[1], op_square(x[2])))
op_f(x[1], op_square(x[2]))

Операторы с тем же именем, что и у существующей функции

Зачастую возникает следующая ошибка:

julia> using JuMP

julia> model = Model();

julia> f(x) = x^2
f (generic function with 1 method)

julia> @operator(model, f, 1, f)
ERROR: Unable to add the nonlinear operator `:f` with the same name as
an existing function.
[...]

Причина этой ошибки в том, что вызов @operator(model, f, 1, f) равносилен следующему:

julia> f = add_nonlinear_operator(model, 1, f; name = :f)

но функция f уже есть в Julia.

Если вычислить эту функцию, не добавляя ее в качестве оператора, JuMP выполнит трассировку функции посредством перегрузки операторов:

julia> @variable(model, x);

julia> f(x)
x²

Чтобы среда JuMP воспринимала функцию f как пользовательский оператор и не отслеживала ее, добавьте оператор с помощью метода add_nonlinear_operator и определите новый метод, который вручную создает выражение NonlinearExpr:

julia> _ = add_nonlinear_operator(model, 1, f; name = :f)
NonlinearOperator(f, :f)

julia> f(x::AbstractJuMPScalar) = NonlinearExpr(:f, Any[x])
f (generic function with 2 methods)

julia> @expression(model, log(f(x)))
log(f(x))

Градиенты и гессианы

По умолчанию в JuMP используется автоматическое дифференцирование для вычисления градиента и гессиана пользовательских операторов. Если функция не поддается автоматическому дифференцированию по умолчанию или можно вычислить аналитические производные, то можно передать дополнительные аргументы в @operator для вычисления первой и второй производных.

В руководстве Automatic differentiation of user-defined operators приведены примеры использования сторонних пакетов Julia для автоматического вычисления производных.

Одномерные функции

Для одномерных функций функция градиента ∇f возвращает число, представляющее производную первого порядка. Помимо этого, можно передать третью функцию, которая возвращает число, представляющее производную второго порядка:

julia> using JuMP


julia> f(x) = x^2
f (generic function with 1 method)

julia> ∇f(x) = 2x
∇f (generic function with 1 method)

julia> ∇²f(x) = 2
∇²f (generic function with 1 method)

julia> model = Model();


julia> @operator(model, op_f, 1, f, ∇f, ∇²f)  # Предоставлять ∇²f необязательно
NonlinearOperator(f, :op_f)

julia> @variable(model, x)
x

julia> @objective(model, Min, op_f(x))
op_f(x)

Многомерные функции

Для многомерных функций функция градиента ∇f должна принимать вектор AbstractVector в качестве первого аргумента, который заполняется на месте. Функция гессиана ∇²f должна принимать в качестве первого аргумента матрицу AbstractMatrix, нижний треугольник которой заполняется на месте:

julia> using JuMP


julia> f(x...) = (1 - x[1])^2 + 100 * (x[2] - x[1]^2)^2
f (generic function with 1 method)

julia> function ∇f(g::AbstractVector{T}, x::T...) where {T}
           g[1] = 400 * x[1]^3 - 400 * x[1] * x[2] + 2 * x[1] - 2
           g[2] = 200 * (x[2] - x[1]^2)
           return
       end
∇f (generic function with 1 method)

julia> function ∇²f(H::AbstractMatrix{T}, x::T...) where {T}
           H[1, 1] = 1200 * x[1]^2 - 400 * x[2] + 2
           # H[1, 2] = -400 * x[1]  <-- Не требуется. Заполняется только нижний треугольник.
           H[2, 1] = -400 * x[1]
           H[2, 2] = 200.0
           return
       end
∇²f (generic function with 1 method)

julia> model = Model();


julia> @operator(model, rosenbrock, 2, f, ∇f, ∇²f)  # Предоставлять ∇²f необязательно
NonlinearOperator(f, :rosenbrock)

julia> @variable(model, x[1:2])
2-element Vector{VariableRef}:
 x[1]
 x[2]

julia> @objective(model, Min, rosenbrock(x[1], x[2]))
rosenbrock(x[1], x[2])

Предполагается, что гессианова матрица H инициализирована нулями, и поскольку матрица H симметрична, достаточно заполнить ненулевыми значениями элементы нижнего треугольника. Тип передаваемой матрицы H зависит от системы автоматического дифференцирования, поэтому первый аргумент функции гессиана должен поддерживать тип AbstractMatrix (он может быть отличен от Matrix{Float64}). Более того, изначально предполагается только то, что H поддерживает size(H) и setindex!. Наконец, матрица считается плотной, поэтому для функций с входными данными высокой размерности производительность будет низкой.

Пользовательские операторы с векторными входными данными

Пользовательские операторы, которые принимают векторы в качестве входных аргументов (например, f(x::Vector)), не поддерживаются. Вместо этого используйте синтаксис распаковки Julia для создания функции со скалярными аргументами. Например, вместо следующего кода:

f(x::Vector) = sum(x[i]^i for i in 1:length(x))

определите:

f(x...) = sum(x[i]^i for i in 1:length(x))

Еще один подход заключается в определении распаковываемой функции как анонимной:

julia> using JuMP


julia> model = Model();


julia> @variable(model, x[1:5])
5-element Vector{VariableRef}:
 x[1]
 x[2]
 x[3]
 x[4]
 x[5]

julia> f(x::Vector) = sum(x[i]^i for i in 1:length(x))
f (generic function with 1 method)

julia> @operator(model, op_f, 5, (x...) -> f(collect(x)))
NonlinearOperator(#6, :op_f)

julia> @objective(model, Min, op_f(x...))
op_f(x[1], x[2], x[3], x[4], x[5])

Если оператор принимает несколько входных векторов, напишите функцию, которая принимает распакованные аргументы и воссоздает требуемые входные векторы:

julia> using JuMP


julia> model = Model();


julia> @variable(model, x[1:2]);


julia> @variable(model, y[1:2]);


julia> @variable(model, z);


julia> f(x::Vector, y::Vector, z) = sum((x[i] * y[i])^z for i in 1:2)
f (generic function with 1 method)

julia> f(x, y, z)
((x[1]*y[1]) ^ z) + ((x[2]*y[2]) ^ z)

julia> f_splat(args...) = f(collect(args[1:2]), collect(args[3:4]), args[5])
f_splat (generic function with 1 method)

julia> f_splat(x..., y..., z)
((x[1]*y[1]) ^ z) + ((x[2]*y[2]) ^ z)

julia> @operator(model, op_f, 5, f_splat)
NonlinearOperator(f_splat, :op_f)

julia> @objective(model, Min, op_f(x..., y..., z))
op_f(x[1], x[2], y[1], y[2], z)

Распространенные ошибки при написании пользовательских операторов

JuMP использует ForwardDiff.jl для вычисления производных первого порядка пользовательских операторов. ForwardDiff имеет ряд ограничений, о которых следует помнить при написании пользовательских операторов.

В оставшейся части этого раздела приводятся советы по отладке и разъясняются некоторые распространенные ошибки.

Получаете ошибку наподобие No method matching Float64(::ForwardDiff.Dual)? Прочитайте этот раздел.

Отладка

При добавлении оператора, который не поддерживает ForwardDiff, выводится длинное сообщение об ошибке. Для получения дополнительной информации можно просмотреть трассировку стека, но зачастую бывает сложно понять, почему и где функция дает сбой.

Может оказаться полезным выполнить отладку оператора вне JuMP, как описано ниже.

Если оператор одномерный, выполните следующий код:

julia> import ForwardDiff

julia> my_operator(a) = a^2
my_operator (generic function with 1 method)

julia> ForwardDiff.derivative(my_operator, 1.0)
2.0

Если оператор многомерный, выполните следующий код:

julia> import ForwardDiff

julia> my_operator(a, b) = a^2 + b^2
my_operator (generic function with 1 method)

julia> ForwardDiff.gradient(x -> my_operator(x...), [1.0, 2.0])
2-element Vector{Float64}:
 2.0
 4.0

Обратите внимание: хотя оператор принимает распакованные аргументы, ForwardDiff.gradient требует вектор в качестве входных данных.

Оператор вызывает что-то, не поддерживаемое ForwardDiff

ForwardDiff работает путем перегрузки многих функций Julia для специального типа ForwardDiff.Dual <: Real. Если ваш оператор попытается вызвать функцию, для которой не определена перегрузка, будет выдана ошибка MethodError.

Например, оператор не может вызывать внешние функции C или выступать в качестве оптимального целевого значения модели JuMP.

julia> import ForwardDiff

julia> my_operator_bad(x) = @ccall sqrt(x::Cdouble)::Cdouble
my_operator_bad (generic function with 1 method)

julia> ForwardDiff.derivative(my_operator_bad, 1.0)
ERROR: MethodError: no method matching Float64(::ForwardDiff.Dual{ForwardDiff.Tag{typeof(my_operator_bad), Float64}, Float64, 1})
[...]

К сожалению, список вызовов, поддерживаемых ForwardDiff, слишком велик, чтобы перечислить все разрешенные и запрещенные варианты. Поэтому лучше всего попробовать на практике и посмотреть, что получится.

Оператор не принимает распакованные входные данные

Оператор принимает в качестве входных данных f(x::Vector), а не распакованные аргументы f(x...).

julia> import ForwardDiff

julia> my_operator_bad(x::Vector) = sum(x[i]^2 for i in eachindex(x))
my_operator_bad (generic function with 1 method)

julia> my_operator_good(x...) = sum(x[i]^2 for i in eachindex(x))
my_operator_good (generic function with 1 method)

julia> ForwardDiff.gradient(x -> my_operator_bad(x...), [1.0, 2.0])
ERROR: MethodError: no method matching my_operator_bad(::ForwardDiff.Dual{ForwardDiff.Tag{var"#5#6", Float64}, Float64, 2}, ::ForwardDiff.Dual{ForwardDiff.Tag{var"#5#6", Float64}, Float64, 2})
[...]

julia> ForwardDiff.gradient(x -> my_operator_good(x...), [1.0, 2.0])
2-element Vector{Float64}:
 2.0
 4.0

Оператор ожидает входные данные типа `Float64`

Оператор предполагает, что в него будут переданы входные данные типа Float64, но должен работать с любым универсальным типом Real.

julia> import ForwardDiff

julia> my_operator_bad(x::Float64...) = sum(x[i]^2 for i in eachindex(x))
my_operator_bad (generic function with 1 method)

julia> my_operator_good(x::Real...) = sum(x[i]^2 for i in eachindex(x))
my_operator_good (generic function with 1 method)

julia> ForwardDiff.gradient(x -> my_operator_bad(x...), [1.0, 2.0])
ERROR: MethodError: no method matching my_operator_bad(::ForwardDiff.Dual{ForwardDiff.Tag{var"#5#6", Float64}, Float64, 2}, ::ForwardDiff.Dual{ForwardDiff.Tag{var"#5#6", Float64}, Float64, 2})
[...]

julia> ForwardDiff.gradient(x -> my_operator_good(x...), [1.0, 2.0])
2-element Vector{Float64}:
 2.0
 4.0

Оператор выделяет память для типа `Float64`

Оператор временно выделяет память с помощью zeros(3) или аналогичного метода. По умолчанию предполагается тип Float64, поэтому используйте вместо этого вызов zeros(T, 3).

julia> import ForwardDiff

julia> function my_operator_bad(x::Real...)
           # В этой строке есть проблема. zeros(n) — это краткая форма записи для zeros(Float64, n)
           y = zeros(length(x))
           for i in eachindex(x)
               y[i] = x[i]^2
           end
           return sum(y)
       end
my_operator_bad (generic function with 1 method)

julia> function my_operator_good(x::T...) where {T<:Real}
           y = zeros(T, length(x))
           for i in eachindex(x)
               y[i] = x[i]^2
           end
           return sum(y)
       end
my_operator_good (generic function with 1 method)

julia> ForwardDiff.gradient(x -> my_operator_bad(x...), [1.0, 2.0])
ERROR: MethodError: no method matching Float64(::ForwardDiff.Dual{ForwardDiff.Tag{var"#1#2", Float64}, Float64, 2})
[...]

julia> ForwardDiff.gradient(x -> my_operator_good(x...), [1.0, 2.0])
2-element Vector{Float64}:
 2.0
 4.0

Нелинейное моделирование

Задание нелинейной целевой функции

Добавление нелинейного ограничения

Добавление параметра

Создание нелинейного выражения

Автоматическое дифференцирование

Подробные сведения о нелинейных выражениях

Конструкторы

Поддерживаемые аргументы

Поддерживаемые операторы

Ограничения

Поля

Принудительное создание нелинейных выражений

Трассировка функций

Пользовательские операторы

Добавление оператора

Автоматическое дифференцирование

Добавление оператора без макроса

Операторы с тем же именем, что и у существующей функции

Градиенты и гессианы

Одномерные функции

Многомерные функции

Пользовательские операторы с векторными входными данными

Распространенные ошибки при написании пользовательских операторов

Отладка

Оператор вызывает что-то, не поддерживаемое ForwardDiff

Оператор не принимает распакованные входные данные

Оператор ожидает входные данные типа Float64

Оператор выделяет память для типа Float64

Оператор ожидает входные данные типа `Float64`

Оператор выделяет память для типа `Float64`