Hello.

I need to translate some matlab code into Julia, where i need to find the roots in X of the function:

```
function problem(X :: AbstractVector{T}, A, β, Ω, α, ε, θ, σ,L) where T
N = length(α)
p = @view X[1:N]
y = @view X[N+1:end]
Out :: Vector{T} = zeros(eltype(X),2*N)
q = (Ω * p .^ (1-θ)) .^ (1 / (1 - θ))
w = p .* (A .^ ((ε - 1)/ε)) .* (α .^ (1 / ε)) .* (y .^ (1/ε)) .* L .^ (-1/ε)
C = w' * L
Out[1:N] = p - (A .^ (ε - 1) .* (α .* w .^ (1- ε) + (1 .- α) .* q .^ (1 - ε))) .^ (1/(1-ε))
Out[N+1:end] = y' - y' * diagm(p)^ε * diagm(A)^(ε-1) * diagm(q)^(θ-ε) * diagm(1 .- α) * Ω * diagm(p)^(-θ) - β'*diagm(p)^(-σ)*C
return Out
end
```

β, Ω, α, ε, θ, σ and L are fixed parameters. A is a random variable.

Nonlinearsolve.jl an NLsolve.jl both are able to give exact solutions for norm A close enought to one, however for smaller A the problem function returns complex values, and here I run into problems.

If I allow complex values either ForwardDiff stops working, or I get a Innexact Error by the solver. I tried simply squaring `Out`

, however I get very unexact results then. Fsolve works fine it matlab for the same problem.

Thank you in advance for your answers.