jacobian(f, x) silently returns zeros if f is complex-valued

Aberta
#830 0 comentários 0 reações 0 responsáveis Ver no GitHub

Mantenedores costumam responder em até 1 dia

Ninguém assumiu esta issue ainda.

Avaliação

Dificuldade
3/5
Tempo estimado
1-2 dias
Facilidade para iniciantes
75/100
Tipo de issue
Bug
Clareza
Claramente especificada
Status de atividade
Pouca atividade
Stack de tecnologia
julia
Domínio
tooling

Direção de pesquisa

Comece lendo extract_jacobian!, extract_jacobian_chunk!, extract_derivative e valtype, concentrando-se em como os valores Complex{Dual} são tratados. Exercite o modo vetorial, o modo chunk e jacobian! com o reproducer da issue. A tarefa estará concluída quando funções com valores complexos produzirem uma matriz Jacobiana ComplexF64 correta, enquanto os caminhos com valores reais permanecerem inalterados.

Escrita pelo modelo de indexação a partir do texto da issue.

Descrição

jacobian returns zeros, with Duals leaked into the eltype, instead of the Jacobian or an error:

julia> g(x) = [cis(x[1])*x[2], x[1] + 2im*x[2]];

julia> J = ForwardDiff.jacobian(g, [0.5, 1.0]);

julia> eltype(J)
Complex{ForwardDiff.Dual{ForwardDiff.Tag{typeof(g), Float64}, Float64, 2}}

julia> all(iszero, J)
true

extract_jacobian! extracts with partials(T, ydual, i), which for a Complex{Dual} hits the partials(x, i...) = zero(x) fallback. derivative(f, x) is fine since it goes through extract_derivative, which handles Complex{<:Dual}.

extract_derivative cannot be reused directly because it hardcodes the partial index, while the Jacobian needs index i. But generalising it to an index is enough, and keeps the knowledge of how a complex number carries derivative information in one place:

@inline extract_partial(::Type{T}, y::Real, i) where {T} = zero(y)
@inline extract_partial(::Type{T}, y::Complex, i) where {T} = zero(y)
@inline extract_partial(::Type{T}, y::Dual, i) where {T} = partials(T, y, i)
@inline extract_partial(::Type{T}, y::Complex{<:Dual}, i) where {T} =
    complex(partials(T, real(y), i), partials(T, imag(y), i))

extract_jacobian!/extract_jacobian_chunk! then use extract_partial(T, ·, i) instead of partials(T, ·, i), and the scalar extract_derivative methods become the i = 1 case.

One more piece is needed for the result eltype, since jacobian allocates similar(ydual, valtype(T, eltype(ydual)), ...):

@inline valtype(::Type{T}, ::Type{Complex{D}}) where {T,D} = Complex{valtype(T, D)}

With those, vector mode, chunk mode and jacobian! all give the correct ComplexF64 Jacobian and the real paths are unchanged.

If we would rather not support complex-valued f here, it should at least throw rather than return zeros.

Linguagem predominante
Julia
Estrelas
1k
Forks
165
Merge médio
8d 22h
PRs com merge (30d)
3

Preparar o ambiente

Este projeto não oferece contêiner de desenvolvimento, Dockerfile nem guia de contribuição, então a configuração fica por sua conta: comece pelo README e veja nosso guia da primeira contribuição para os passos gerais.

Primeiros passos

  1. Leia a issue inteira e depois o guia de contribuição do projeto.
  2. Comente na issue dizendo que vai assumir — evita que duas pessoas façam o mesmo trabalho.
  3. Faça um fork do repositório e trabalhe em uma branch.
  4. Abra um pull request que referencie o número da issue.

Mais de JuliaDiff/ForwardDiff.jl

Todas as issues de JuliaDiff/ForwardDiff.jl

Issues semelhantes

Mais issues de Julia

Receba novas issues na sua caixa de entrada

Um resumo curto de issues do GitHub para quem está começando.