有没有办法限制 NLsolve 在 Julia 中搜索的区域?

Is there any way to bound the region searched by NLsolve in Julia?

我试图找到一个非线性(大致为四次)方程的根。 该方程始终有四个根,其中一对接近于零,一个大的正根和一个大的负根。我想确定任何一个接近零的根,但是 nlsolve,即使初始猜测非常接近这些根,似乎也总是收敛于大的正根或负根。

函数图基本上看起来像一个恒定的负值,在零附近有一个(非常窄的)偶数极点,并逐渐上升到在大的正根和负根处过零。

有什么方法可以限制 nlsolve 搜索的区域,或者做些什么让它对我的函数中这个极点的存在更加敏感?

编辑: 下面是重现问题的一些示例代码:

using NLsolve

function f!(F,x)
    x = x[1]
    F[1] = -15000 + x^4 / (x+1e-5)^2
end
# nlsolve will find the root at -122
nlsolve(f!,[0.0])

作为输出,我得到:

Results of Nonlinear Solver Algorithm
 * Algorithm: Trust-region with dogleg and autoscaling
 * Starting Point: [0.0]
 * Zero: [-122.47447713915808]
 * Inf-norm of residuals: 0.000000
 * Iterations: 15
 * Convergence: true
   * |x - x'| < 0.0e+00: false
   * |f(x)| < 1.0e-08: true
 * Function Calls (f): 16
 * Jacobian Calls (df/dx): 6

在这种情况下,我们可以通过将 objective 函数转换为多项式来找到确切的根:

using PolynomialRoots
roots([-1.5e-6,-0.3,-15000,0,1])

产生

4-element Array{Complex{Float64},1}:
     122.47449713915809 - 0.0im
    -122.47447713915808 + 0.0im
 -1.0000000813048448e-5 + 0.0im
  -9.999999186951818e-6 + 0.0im

我希望有一种方法可以在不知道 objective 函数的确切形式的情况下识别 x = -1e-5 处极点周围的一对根。

编辑2: 尝试 Roots.jl :

using Roots
f(x) = -15000 + x^4 / (x+1e-5)^2

find_zero(f,0.0) # finds +122... root
find_zero(f,(-1e-4,0.0)) # error, not a bracketing interval
find_zeros(f,-1e-4,0.0) # finds 0-element Array{Float64,1}
find_zeros(f,-1e-4,0.0,no_pts=6) # finds root slightly less than -1e-5
find_zeros(f,-1e-4,0.0,no_pts=10) # finds 0-element Array{Float64,1}, sensitive to value of no_pts

我可以让 find_zeros 工作,但它对 no_pts 参数和我选择的端点的确切值非常敏感。在 no_pts 上进行循环并获取第一个非空结果可能会起作用,但更可取的是更确定性的收敛。

编辑3: 这里应用了 Bogumił

建议的 tanh 变换
using NLsolve
function f_tanh!(F,x)
    x = x[1]
    x = -1e-4 * (tanh(x)+1) / 2
    F[1] = -15000 + x^4 / (x+1e-5)^2
end

nlsolve(f_tanh!,[100.0]) # doesn't converge
nlsolve(f_tanh!,[1e5]) # doesn't converge

using Roots
function f_tanh(x)
    x = -1e-4 * (tanh(x)+1) / 2
    return -15000 + x^4 / (x+1e-5)^2
end

find_zeros(f_tanh,-1e10,1e10) # 0-element Array
find_zeros(f_tanh,-1e3,1e3,no_pts=100) # 0-element Array
find_zero(f_tanh,0.0) # convergence failed
find_zero(f_tanh,0.0,max_evals=1_000_000,maxfnevals=1_000_000) # convergence failed

EDIT4:这种技术组合在大约 95% 的时间内至少识别出一个根,这对我来说已经足够了。

using Peaks
using Primes
using Roots

# randomize pole location
a = 1e-4*rand()
f(x) = -15000 + x^4 / (x+a)^2

# do an initial sample to find the pole location
l = 1000
minval = -1e-4
maxval = 0
m = []
sample_r = []
while l < 1e6
    sample_r = range(minval,maxval,length=l)
    rough_sample = f.(sample_r)
    m = maxima(rough_sample)
    if length(m) > 0
        break
    else
        l *= 10
    end
end
guess = sample_r[m[1]]

# functions to compress the range around the estimated pole
cube(x) = (x-guess)^3 + guess 
uncube(x) = cbrt(x-guess) + guess
f_cube(x) = f(cube(x))

shift = l ÷ 1000
low = sample_r[m[1]-shift]
high = sample_r[m[1]+shift]
# search only over prime no_pts, so no samplings divide into each other
# possibly not necessary?
for i in primes(500)
    z = find_zeros(f_cube,uncube(low),uncube(high),no_pts=i)
    if length(z)>0
        println(i)
        println(cube.(z))
        break
    end
end

如果您提供更多关于您的问题的信息,我们会给出更多评论。

但总的来说:

  1. 看来你的问题是单变量的,在这种情况下你可以使用 Roots.jl 其中 find_zerofind_zeros 给出你要求的界面(即允许指定搜索区域)
  2. 如果问题是多变量的,那么在 nlsolve 的问题规范中,您有多种选择(因为默认情况下不允许指定边界框 AFAICT)。最简单的是使用变量转换。例如。您可以应用 ai * tanh(xi) + bi 转换,为每个变量选择 aibi,以便将其限制在所需的区间

您在定义中遇到的第一个问题是您定义 f 的方式永远不会越过 0 您正在寻找的两个根附近,因为 Float64 没有足够的写 1e-5 时的精度。您需要使用更高的计算精度:

julia> using Roots

julia> f(x) = -15000 + x^4 / (x+1/big(10.0^5))^2
f (generic function with 1 method)

julia> find_zeros(f,big(-2*10^-5), big(-8*10^-6), no_pts=100)
2-element Array{BigFloat,1}:
 -1.000000081649671426108658262468117284940444265467160592853348997523986352593615e-05
 -9.999999183503552405580084054429938261707450678661727461293670518591720605751116e-06

并将 no_pts 设置得足够大以找到包含根的区间。