From 15b4f0108acde0f4023ca112f846d7c0231ff00d Mon Sep 17 00:00:00 2001
From: FHoltorf <32248677+FHoltorf@users.noreply.github.com>
Date: Wed, 27 Sep 2023 14:37:49 -0400
Subject: [PATCH 1/2] add test problem library as test set

---
 test/basictests.jl   |   2 +-
 test/runtests.jl     |   1 +
 test/testproblems.jl | 570 +++++++++++++++++++++++++++++++++++++++++++
 3 files changed, 572 insertions(+), 1 deletion(-)
 create mode 100644 test/testproblems.jl

diff --git a/test/basictests.jl b/test/basictests.jl
index 05a0152fa..87b9fc7c0 100644
--- a/test/basictests.jl
+++ b/test/basictests.jl
@@ -62,7 +62,7 @@ u0 = [1.0, 1.0]
 
 precs = [
     NonlinearSolve.DEFAULT_PRECS,
-    (args...) -> (Diagonal(rand!(similar(u0))), nothing)
+    (args...) -> (Diagonal(rand!(similar(u0))), nothing),
 ]
 
 for prec in precs, linsolve in (nothing, KrylovJL_GMRES())
diff --git a/test/runtests.jl b/test/runtests.jl
index a84fc3cb1..a2a1b53f2 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -15,6 +15,7 @@ end
     if GROUP == "All" || GROUP == "Core"
         @time @safetestset "Basic Tests + Some AD" include("basictests.jl")
         @time @safetestset "Sparsity Tests" include("sparse.jl")
+        @time @safetestset "Test Problem Library" include("testproblems.jl")
     end
 
     if GROUP == "GPU"
diff --git a/test/testproblems.jl b/test/testproblems.jl
new file mode 100644
index 000000000..1f48fb3ba
--- /dev/null
+++ b/test/testproblems.jl
@@ -0,0 +1,570 @@
+using NonlinearSolve, LinearAlgebra
+
+# Implementation of the 23 test problems in
+# [test_nonlin](https://people.sc.fsu.edu/~jburkardt/m_src/test_nonlin/test_nonlin.html)
+
+# ------------------------------------- Problem 1 ------------------------------------------
+function p1_f!(out, x, p = nothing)
+    n = length(x)
+    out[1] = 1.0 - x[1]
+    out[2:n] .= 10.0 .* (x[2:n] .- x[1:(n - 1)] .* x[1:(n - 1)])
+    nothing
+end
+
+n = 10
+x_sol = ones(n)
+x_start = ones(n)
+x_start[1] = -1.2
+p1_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Generalized Rosenbrock function")
+
+# ------------------------------------- Problem 2 ------------------------------------------
+function p2_f!(out, x, p = nothing)
+    out[1] = x[1] + 10.0 * x[2]
+    out[2] = sqrt(5.0) * (x[3] - x[4])
+    out[3] = (x[2] - 2.0 * x[3])^2
+    out[4] = sqrt(10.0) * (x[1] - x[4]) * (x[1] - x[4])
+    nothing
+end
+
+n = 4
+x_sol = zeros(n)
+x_start = [3.0, -1.0, 0.0, 1.0]
+p2_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Powell singular function")
+
+# ------------------------------------- Problem 3 ------------------------------------------
+function p3_f!(out, x, p = nothing)
+    out[1] = 10000.0 * x[1] * x[2] - 1.0
+    out[2] = exp(-x[1]) + exp(-x[2]) - 1.0001
+    nothing
+end
+
+n = 2
+x_sol = [1.098159e-05, 9.106146]
+x_start = [0.0, 1.0]
+p3_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Powell badly scaled function")
+
+# ------------------------------------- Problem 4 ------------------------------------------
+function p4_f!(out, x, p = nothing)
+    temp1 = x[2] - x[1] * x[1]
+    temp2 = x[4] - x[3] * x[3]
+
+    out[1] = -200.0 * x[1] * temp1 - (1.0 - x[1])
+    out[2] = 200.0 * temp1 + 20.2 * (x[2] - 1.0) + 19.8 * (x[4] - 1.0)
+    out[3] = -180.0 * x[3] * temp2 - (1.0 - x[3])
+    out[4] = 180.0 * temp2 + 20.2 * (x[4] - 1.0) + 19.8 * (x[2] - 1.0)
+    nothing
+end
+
+n = 4
+x_sol = ones(n)
+x_start = [-3.0, -1.0, -3.0, -1.0]
+p4_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol, "title" => "Wood function")
+
+# ------------------------------------- Problem 5 ------------------------------------------
+function p5_f!(out, x, p = nothing)
+    if 0.0 < x[1]
+        temp = atan(x[2] / x[1]) / (2.0 * pi)
+    elseif x[1] < 0.0
+        temp = atan(x[2] / x[1]) / (2.0 * pi) + 0.5
+    else
+        temp = 0.25 * sign(x[2])
+    end
+
+    out[1] = 10.0 * (x[3] - 10.0 * temp)
+    out[2] = 10.0 * (sqrt(x[1] * x[1] + x[2] * x[2]) - 1.0)
+    out[3] = x[3]
+    nothing
+end
+
+n = 3
+x_sol = [1.0, 0.0, 0.0]
+x_start = [-1.0, 0.0, 0.0]
+p5_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Helical valley function")
+
+# ------------------------------------- Problem 6 ------------------------------------------
+function p6_f!(out, x, p = nothing)
+    n = length(x)
+    for i in 1:29
+        ti = i / 29.0
+        sum1 = 0.0
+        temp = 1.0
+        for j in 2:n
+            sum1 = sum1 + j * temp * x[j]
+            temp = ti * temp
+        end
+
+        sum2 = 0.0
+        temp = 1.0
+        for j in 1:n
+            sum2 = sum2 + temp * x[j]
+            temp = ti * temp
+        end
+        temp = 1.0 / ti
+
+        for k in 1:n
+            out[k] = out[k] + temp * (sum1 - sum2 * sum2 - 1.0) * (k - 2.0 * ti * sum2)
+            temp = ti * temp
+        end
+    end
+
+    out[1] = out[1] + 3.0 * x[1] - 2.0 * x[1] * x[1] + 2.0 * x[1]^3
+    out[2] = out[2] + x[2] - x[2]^2 - 1.0
+    nothing
+end
+
+n = 2
+x_sol = []
+x_start = zeros(n)
+p6_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol, "title" => "Watson function")
+
+# ------------------------------------- Problem 7 ------------------------------------------
+function p7_f!(out, x, p = nothing)
+    n = length(x)
+    out .= 0.0
+    for j in 1:n
+        t1 = 1.0
+        t2 = x[j]
+        for i in 1:n
+            out[i] += t2
+            t3 = 2.0 * x[j] * t2 - t1
+            t1 = t2
+            t2 = t3
+        end
+    end
+    out ./= n
+
+    for i in 1:n
+        ip1 = i
+        if ip1 % 2 == 0
+            out[i] = out[i] + 1.0 / (ip1 * ip1 - 1)
+        end
+    end
+    nothing
+end
+
+n = 2
+x_sol = [0.2113248654051871, 0.7886751345948129]
+x_sol .= 2.0 .* x_sol .- 1.0
+x_start = zeros(n)
+for i in 1:n
+    x_start[i] = (2 * i - n) / (n + 1)
+end
+p7_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Chebyquad function")
+
+# ------------------------------------- Problem 8 ------------------------------------------
+function p8_f!(out, x, p = nothing)
+    n = length(x)
+    out[1:(n - 1)] .= x[1:(n - 1)] .+ sum(x) .- (n + 1)
+    out[n] = prod(x) - 1.0
+    nothing
+end
+
+n = 10
+x_sol = ones(n)
+x_start = ones(n) ./ 2
+p8_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Brown almost linear function")
+
+# ------------------------------------- Problem 9 ------------------------------------------
+function p9_f!(out, x, p = nothing)
+    n = length(x)
+    h = 1.0 / (n + 1)
+    for k in 1:n
+        out[k] = 2.0 * x[k] + 0.5 * h^2 * (x[k] + k * h + 1.0)^3
+        if 1 < k
+            out[k] = out[k] - x[k - 1]
+        end
+        if k < n
+            out[k] = out[k] - x[k + 1]
+        end
+    end
+    nothing
+end
+
+n = 10
+x_sol = []
+x_start = ones(n)
+for i in 1:n
+    x_start[i] = (i * (i - n - 1)) / (n + 1)^2
+end
+p9_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Discrete boundary value function")
+
+# ------------------------------------- Problem 10 -----------------------------------------
+function p10_f!(out, x, p = nothing)
+    n = length(x)
+    h = 1.0 / (n + 1)
+    for k in 1:n
+        tk = k / (n + 1)
+        sum1 = 0.0
+        for j in 1:k
+            tj = j * h
+            sum1 = sum1 + tj * (x[j] + tj + 1.0)^3
+        end
+        sum2 = 0.0
+        for j in k:n
+            tj = j * h
+            sum2 = sum2 + (1.0 - tj) * (x[j] + tj + 1.0)^3
+        end
+
+        out[k] = x[k] + h * ((1.0 - tk) * sum1 + tk * sum2) / 2.0
+    end
+    nothing
+end
+
+n = 10
+x_sol = []
+x_start = zeros(n)
+for i in 1:n
+    x_start[i] = (i * (i - n - 1)) / (n + 1)^2
+end
+p10_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Discrete integral equation function")
+
+# ------------------------------------- Problem 11 -----------------------------------------
+function p11_f!(out, x, p = nothing)
+    n = length(x)
+    c_sum = sum(cos.(x))
+    for k in 1:n
+        out[k] = n - c_sum + k * (1.0 - cos(x[k])) - sin(x[k])
+    end
+    nothing
+end
+
+n = 10
+x_sol = []
+x_start = ones(n) / n
+p11_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Trigonometric function")
+
+# ------------------------------------- Problem 12 -----------------------------------------
+function p12_f!(out, x, p = nothing)
+    n = length(x)
+    sum1 = 0.0
+    for j in 1:n
+        sum1 += j * (x[j] - 1.0)
+    end
+    for k in 1:n
+        out[k] = x[k] - 1.0 + k * sum1 * (1.0 + 2.0 * sum1 * sum1)
+    end
+    nothing
+end
+
+n = 10
+x_sol = ones(n)
+x_start = zeros(n)
+for i in 1:n
+    x_start[i] = 1.0 - i / n
+end
+p12_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Variably dimensioned function")
+
+# ------------------------------------- Problem 13 -----------------------------------------
+function p13_f!(out, x, p = nothing)
+    n = length(x)
+    for k in 1:n
+        out[k] = (3.0 - 2.0 * x[k]) * x[k] + 1.0
+        if 1 < k
+            out[k] -= x[k - 1]
+        end
+        if k < n
+            out[k] -= 2.0 * x[k + 1]
+        end
+    end
+    nothing
+end
+
+n = 10
+x_sol = []
+x_start = ones(n) .* (-1.0)
+p13_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Broyden tridiagonal function")
+
+# ------------------------------------- Problem 14 -----------------------------------------
+function p14_f!(out, x, p = nothing)
+    n = length(x)
+    ml = 5
+    mu = 1
+    for k in 1:n
+        k1 = max(1, k - ml)
+        k2 = min(n, k + mu)
+
+        temp = 0.0
+        for j in k1:k2
+            if j != k
+                temp += x[j] * (1.0 + x[j])
+            end
+        end
+        out[k] = x[k] * (2.0 + 5.0 * x[k] * x[k]) + 1.0 - temp
+    end
+    nothing
+end
+
+n = 10
+x_sol = []
+x_start = ones(n) .* (-1.0)
+p14_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Broyden banded function")
+
+# ------------------------------------- Problem 15 -----------------------------------------
+function p15_f!(out, x, p = nothing)
+    out[1] = (x[1] * x[1] + x[2] * x[3]) - 0.0001
+    out[2] = (x[1] * x[2] + x[2] * x[4]) - 1.0
+    out[3] = (x[3] * x[1] + x[4] * x[3]) - 0.0
+    out[4] = (x[3] * x[2] + x[4] * x[4]) - 0.0001
+    nothing
+end
+
+n = 4
+x_sol = [0.01, 50.0, 0.0, 0.01]
+x_start = [1.0, 0.0, 0.0, 1.0]
+p15_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Hammarling 2 by 2 matrix square root problem")
+
+# ------------------------------------- Problem 16 -----------------------------------------
+function p16_f!(out, x, p = nothing)
+    out[1] = (x[1] * x[1] + x[2] * x[4] + x[3] * x[7]) - 0.0001
+    out[2] = (x[1] * x[2] + x[2] * x[5] + x[3] * x[8]) - 1.0
+    out[3] = x[1] * x[3] + x[2] * x[6] + x[3] * x[9]
+    out[4] = x[4] * x[1] + x[5] * x[4] + x[6] * x[7]
+    out[5] = (x[4] * x[2] + x[5] * x[5] + x[6] * x[8]) - 0.0001
+    out[6] = x[4] * x[3] + x[5] * x[6] + x[6] * x[9]
+    out[7] = x[7] * x[1] + x[8] * x[4] + x[9] * x[7]
+    out[8] = x[7] * x[2] + x[8] * x[5] + x[9] * x[8]
+    out[9] = (x[7] * x[3] + x[8] * x[6] + x[9] * x[9]) - 0.0001
+    nothing
+end
+
+n = 9
+x_sol = [0.01, 50.0, 0.0, 0.0, 0.01, 0.0, 0.0, 0.0, 0.01]
+x_start = [1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0]
+p16_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Hammarling 3 by 3 matrix square root problem")
+
+# ------------------------------------- Problem 17 -----------------------------------------
+function p17_f!(out, x, p = nothing)
+    out[1] = x[1] + x[2] - 3.0
+    out[2] = x[1]^2 + x[2]^2 - 9.0
+    nothing
+end
+
+n = 2
+x_sol = [0.0, 3.0]
+x_start = [1.0, 5.0]
+p17_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Dennis and Schnabel 2 by 2 example")
+
+# ------------------------------------- Problem 18 -----------------------------------------
+function p18_f!(out, x, p = nothing)
+    if x[1] != 0.0
+        out[1] = x[2]^2 * (1.0 - exp(-x[1] * x[1])) / x[1]
+    else
+        out[1] = 0.0
+    end
+    if x[2] != 0.0
+        out[2] = x[1] * (1.0 - exp(-x[2] * x[2])) / x[2]
+    else
+        out[2] = 0.0
+    end
+    nothing
+end
+
+n = 2
+x_sol = zeros(n)
+x_start = [2.0, 2.0]
+p18_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Sample problem 18")
+
+# ------------------------------------- Problem 19 -----------------------------------------
+function p19_f!(out, x, p = nothing)
+    out[1] = x[1] * (x[1]^2 + x[2]^2)
+    out[2] = x[2] * (x[1]^2 + x[2]^2)
+    nothing
+end
+
+n = 2
+x_sol = zeros(n)
+x_start = [3.0, 3.0]
+p19_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Sample problem 19")
+
+# ------------------------------------- Problem 20 -----------------------------------------
+function p20_f!(out, x, p = nothing)
+    out[1] = x[1] * (x[1] - 5.0)^2
+    nothing
+end
+
+n = 1
+x_sol = [5.0] # OR [0.0]...
+x_start = [1.0]
+p20_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Scalar problem f(x) = x(x - 5)^2")
+
+# ------------------------------------- Problem 21 -----------------------------------------
+function p21_f!(out, x, p = nothing)
+    out[1] = x[1] - x[2]^3 + 5.0 * x[2]^2 - 2.0 * x[2] - 13.0
+    out[2] = x[1] + x[2]^3 + x[2]^2 - 14.0 * x[2] - 29.0
+    nothing
+end
+
+n = 2
+x_sol = [5.0, 4.0]
+x_start = [0.5, -2.0]
+p21_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Freudenstein-Roth function")
+
+# ------------------------------------- Problem 22 -----------------------------------------
+function p22_f!(out, x, p = nothing)
+    out[1] = x[1] * x[1] - x[2] + 1.0
+    out[2] = x[1] - cos(0.5 * pi * x[2])
+    nothing
+end
+
+n = 2
+x_sol = [0.0, 1.0]
+x_start = [1.0, 0.0]
+p22_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Boggs function")
+
+# ------------------------------------- Problem 23 -----------------------------------------
+function p23_f!(out, x, p = nothing)
+    c = 0.9
+    out[1:n] = x[1:n]
+    μ = zeros(n)
+    for i in 1:n
+        μ[i] = (2 * i) / (2 * n)
+    end
+    for i in 1:n
+        s = 0.0
+        for j in 1:n
+            s = s + (μ[i] * x[j]) / (μ[i] + μ[j])
+        end
+        term = 1.0 - c * s / (2 * n)
+        out[i] -= 1.0 / term
+    end
+    nothing
+end
+
+n = 10
+x_sol = []
+x_start = ones(n)
+p23_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
+    "title" => "Chandrasekhar function")
+
+# ----------------------------------- Solve problems ---------------------------------------
+problems = (p1_f!, p2_f!, p3_f!, p4_f!, p5_f!, p6_f!, p7_f!, p8_f!, p9_f!, p10_f!, p11_f!,
+    p12_f!, p13_f!, p14_f!, p15_f!, p16_f!, p17_f!, p18_f!, p19_f!, p20_f!, p21_f!,
+    p22_f!, p23_f!)
+dicts = (p1_dict, p2_dict, p3_dict, p4_dict, p5_dict, p6_dict, p7_dict, p8_dict, p9_dict,
+    p10_dict, p11_dict, p12_dict, p13_dict, p14_dict, p15_dict, p16_dict, p17_dict,
+    p18_dict, p19_dict, p20_dict, p21_dict, p22_dict, p23_dict)
+
+function test_on_library(problems, dicts, alg_ops, broken_tests, ϵ = 1e-5)
+    for (idx, (problem, dict)) in enumerate(zip(problems, dicts))
+        x = dict["start"]
+        res = similar(x)
+        nlprob = NonlinearProblem(problem, x)
+        @testset "$(dict["title"])" begin
+            for alg in alg_ops
+                sol = solve(nlprob, alg, abstol = 1e-15, reltol = 1e-15)
+                problem(res, sol.u, nothing)
+                broken = idx in broken_tests[alg] ? true : false
+                @test norm(res)<=ϵ skip=broken
+            end
+        end
+    end
+end
+
+function print_results_on_library(problems, dicts, alg_ops, names)
+    results = []
+    for (idx, (problem, dict)) in enumerate(zip(problems, dicts))
+        x = dict["start"]
+        res = similar(x)
+        nlprob = NonlinearProblem(problem, x)
+        push!(results, Dict())
+        for (alg, name) in zip(alg_ops, names)
+            sol = solve(nlprob, alg, abstol = 1e-15, reltol = 1e-15)
+            problem(res, sol.u, nothing)
+            results[idx][name] = norm(res)
+        end
+    end
+
+    i_str = i_str = rpad("nr", 3, " ")
+    title_str = rpad("Problem", 50, " ")
+    n_str = rpad("n", 5, " ")
+    norm_str = prod(rpad(name, 20, " ") for name in names)
+    println("$i_str $title_str $n_str $norm_str")
+
+    for (i, dict) in enumerate(dicts)
+        local i_str = rpad(string(i), 3, " ")
+        local title_str = rpad(dict["title"], 50, " ")
+        local n_str = rpad(string(dict["n"]), 5, " ")
+        local norm_str = ""
+        for name in names
+            norm_str *= rpad(string(trunc(results[i][name]; sigdigits = 5)), 20, " ")
+        end
+        println("$i_str $title_str $n_str $norm_str")
+    end
+end
+
+function detect_broken_tests(problems, dicts, alg_ops, ϵ = 1e-5)
+    broken_tests = Dict(alg => Int[] for alg in alg_ops)
+    for (idx, (problem, dict)) in enumerate(zip(problems, dicts))
+        x = dict["start"]
+        res = similar(x)
+        nlprob = NonlinearProblem(problem, x)
+        for alg in alg_ops
+            sol = solve(nlprob, alg, abstol = 1e-15, reltol = 1e-15)
+            problem(res, sol.u, nothing)
+            if norm(res) > ϵ
+                push!(broken_tests[alg], idx)
+            end
+        end
+    end
+    return broken_tests
+end
+
+# NewtonRaphson
+@testset "NewtonRaphson test problem library" begin
+    alg_ops = (NewtonRaphson(),)
+
+    # dictionary with indices of test problems where method does not converge to small residual
+    broken_tests = Dict(alg => Int[] for alg in alg_ops)
+    broken_tests[alg_ops[1]] = [1, 6]
+
+    test_on_library(problems, dicts, alg_ops, broken_tests)
+end
+
+@testset "TrustRegion test problem library" begin
+    alg_ops = (TrustRegion(; radius_update_scheme = RadiusUpdateSchemes.Simple),
+        TrustRegion(; radius_update_scheme = RadiusUpdateSchemes.Fan),
+        TrustRegion(; radius_update_scheme = RadiusUpdateSchemes.Hei),
+        TrustRegion(; radius_update_scheme = RadiusUpdateSchemes.Yuan),
+        TrustRegion(; radius_update_scheme = RadiusUpdateSchemes.Bastin))
+
+    # dictionary with indices of test problems where method does not converge to small residual
+    broken_tests = Dict(alg => Int[] for alg in alg_ops)
+    broken_tests[alg_ops[1]] = [1, 2, 3, 4, 5, 6, 7, 9, 15, 16, 21, 22]
+    broken_tests[alg_ops[2]] = [1, 2, 3, 4, 5, 6, 7, 9, 15, 16, 21, 22]
+    broken_tests[alg_ops[3]] = [1, 2, 3, 4, 5, 6, 7, 8, 11, 12, 15, 16, 18, 19, 21, 22]
+    broken_tests[alg_ops[4]] = [1, 2, 3, 4, 5, 6, 8, 9, 15, 16, 21, 22]
+    broken_tests[alg_ops[5]] = [1, 2, 3, 5, 6, 9, 11, 15, 16, 18, 21, 22]
+
+    test_on_library(problems, dicts, alg_ops, broken_tests)
+end
+
+# LevenbergMarquardt seems to still be buggy    
+@testset "TrustRegion test problem library" begin
+    alg_ops = (LevenbergMarquardt(),)
+
+    # dictionary with indices of test problems where method does not converge to small residual
+    broken_tests = Dict(alg => Int[] for alg in alg_ops)
+    broken_tests[alg_ops[1]] = [1, 2, 4, 5, 6, 7, 8, 11, 15, 16, 17, 18, 19, 21, 22]
+
+    test_on_library(problems, dicts, alg_ops, broken_tests)
+end

From 32428e8f12a5fbeec5687986d331cf4c19f7ba8f Mon Sep 17 00:00:00 2001
From: FHoltorf <32248677+FHoltorf@users.noreply.github.com>
Date: Wed, 27 Sep 2023 15:54:57 -0400
Subject: [PATCH 2/2] export problems from NonlinearProblemLibrary

---
 test/testproblems.jl | 467 +------------------------------------------
 1 file changed, 3 insertions(+), 464 deletions(-)

diff --git a/test/testproblems.jl b/test/testproblems.jl
index 1f48fb3ba..c9dd15336 100644
--- a/test/testproblems.jl
+++ b/test/testproblems.jl
@@ -1,468 +1,7 @@
-using NonlinearSolve, LinearAlgebra
+using NonlinearSolve, LinearAlgebra, NonlinearProblemLibrary
 
-# Implementation of the 23 test problems in
-# [test_nonlin](https://people.sc.fsu.edu/~jburkardt/m_src/test_nonlin/test_nonlin.html)
-
-# ------------------------------------- Problem 1 ------------------------------------------
-function p1_f!(out, x, p = nothing)
-    n = length(x)
-    out[1] = 1.0 - x[1]
-    out[2:n] .= 10.0 .* (x[2:n] .- x[1:(n - 1)] .* x[1:(n - 1)])
-    nothing
-end
-
-n = 10
-x_sol = ones(n)
-x_start = ones(n)
-x_start[1] = -1.2
-p1_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Generalized Rosenbrock function")
-
-# ------------------------------------- Problem 2 ------------------------------------------
-function p2_f!(out, x, p = nothing)
-    out[1] = x[1] + 10.0 * x[2]
-    out[2] = sqrt(5.0) * (x[3] - x[4])
-    out[3] = (x[2] - 2.0 * x[3])^2
-    out[4] = sqrt(10.0) * (x[1] - x[4]) * (x[1] - x[4])
-    nothing
-end
-
-n = 4
-x_sol = zeros(n)
-x_start = [3.0, -1.0, 0.0, 1.0]
-p2_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Powell singular function")
-
-# ------------------------------------- Problem 3 ------------------------------------------
-function p3_f!(out, x, p = nothing)
-    out[1] = 10000.0 * x[1] * x[2] - 1.0
-    out[2] = exp(-x[1]) + exp(-x[2]) - 1.0001
-    nothing
-end
-
-n = 2
-x_sol = [1.098159e-05, 9.106146]
-x_start = [0.0, 1.0]
-p3_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Powell badly scaled function")
-
-# ------------------------------------- Problem 4 ------------------------------------------
-function p4_f!(out, x, p = nothing)
-    temp1 = x[2] - x[1] * x[1]
-    temp2 = x[4] - x[3] * x[3]
-
-    out[1] = -200.0 * x[1] * temp1 - (1.0 - x[1])
-    out[2] = 200.0 * temp1 + 20.2 * (x[2] - 1.0) + 19.8 * (x[4] - 1.0)
-    out[3] = -180.0 * x[3] * temp2 - (1.0 - x[3])
-    out[4] = 180.0 * temp2 + 20.2 * (x[4] - 1.0) + 19.8 * (x[2] - 1.0)
-    nothing
-end
-
-n = 4
-x_sol = ones(n)
-x_start = [-3.0, -1.0, -3.0, -1.0]
-p4_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol, "title" => "Wood function")
-
-# ------------------------------------- Problem 5 ------------------------------------------
-function p5_f!(out, x, p = nothing)
-    if 0.0 < x[1]
-        temp = atan(x[2] / x[1]) / (2.0 * pi)
-    elseif x[1] < 0.0
-        temp = atan(x[2] / x[1]) / (2.0 * pi) + 0.5
-    else
-        temp = 0.25 * sign(x[2])
-    end
-
-    out[1] = 10.0 * (x[3] - 10.0 * temp)
-    out[2] = 10.0 * (sqrt(x[1] * x[1] + x[2] * x[2]) - 1.0)
-    out[3] = x[3]
-    nothing
-end
-
-n = 3
-x_sol = [1.0, 0.0, 0.0]
-x_start = [-1.0, 0.0, 0.0]
-p5_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Helical valley function")
-
-# ------------------------------------- Problem 6 ------------------------------------------
-function p6_f!(out, x, p = nothing)
-    n = length(x)
-    for i in 1:29
-        ti = i / 29.0
-        sum1 = 0.0
-        temp = 1.0
-        for j in 2:n
-            sum1 = sum1 + j * temp * x[j]
-            temp = ti * temp
-        end
-
-        sum2 = 0.0
-        temp = 1.0
-        for j in 1:n
-            sum2 = sum2 + temp * x[j]
-            temp = ti * temp
-        end
-        temp = 1.0 / ti
-
-        for k in 1:n
-            out[k] = out[k] + temp * (sum1 - sum2 * sum2 - 1.0) * (k - 2.0 * ti * sum2)
-            temp = ti * temp
-        end
-    end
-
-    out[1] = out[1] + 3.0 * x[1] - 2.0 * x[1] * x[1] + 2.0 * x[1]^3
-    out[2] = out[2] + x[2] - x[2]^2 - 1.0
-    nothing
-end
-
-n = 2
-x_sol = []
-x_start = zeros(n)
-p6_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol, "title" => "Watson function")
-
-# ------------------------------------- Problem 7 ------------------------------------------
-function p7_f!(out, x, p = nothing)
-    n = length(x)
-    out .= 0.0
-    for j in 1:n
-        t1 = 1.0
-        t2 = x[j]
-        for i in 1:n
-            out[i] += t2
-            t3 = 2.0 * x[j] * t2 - t1
-            t1 = t2
-            t2 = t3
-        end
-    end
-    out ./= n
-
-    for i in 1:n
-        ip1 = i
-        if ip1 % 2 == 0
-            out[i] = out[i] + 1.0 / (ip1 * ip1 - 1)
-        end
-    end
-    nothing
-end
-
-n = 2
-x_sol = [0.2113248654051871, 0.7886751345948129]
-x_sol .= 2.0 .* x_sol .- 1.0
-x_start = zeros(n)
-for i in 1:n
-    x_start[i] = (2 * i - n) / (n + 1)
-end
-p7_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Chebyquad function")
-
-# ------------------------------------- Problem 8 ------------------------------------------
-function p8_f!(out, x, p = nothing)
-    n = length(x)
-    out[1:(n - 1)] .= x[1:(n - 1)] .+ sum(x) .- (n + 1)
-    out[n] = prod(x) - 1.0
-    nothing
-end
-
-n = 10
-x_sol = ones(n)
-x_start = ones(n) ./ 2
-p8_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Brown almost linear function")
-
-# ------------------------------------- Problem 9 ------------------------------------------
-function p9_f!(out, x, p = nothing)
-    n = length(x)
-    h = 1.0 / (n + 1)
-    for k in 1:n
-        out[k] = 2.0 * x[k] + 0.5 * h^2 * (x[k] + k * h + 1.0)^3
-        if 1 < k
-            out[k] = out[k] - x[k - 1]
-        end
-        if k < n
-            out[k] = out[k] - x[k + 1]
-        end
-    end
-    nothing
-end
-
-n = 10
-x_sol = []
-x_start = ones(n)
-for i in 1:n
-    x_start[i] = (i * (i - n - 1)) / (n + 1)^2
-end
-p9_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Discrete boundary value function")
-
-# ------------------------------------- Problem 10 -----------------------------------------
-function p10_f!(out, x, p = nothing)
-    n = length(x)
-    h = 1.0 / (n + 1)
-    for k in 1:n
-        tk = k / (n + 1)
-        sum1 = 0.0
-        for j in 1:k
-            tj = j * h
-            sum1 = sum1 + tj * (x[j] + tj + 1.0)^3
-        end
-        sum2 = 0.0
-        for j in k:n
-            tj = j * h
-            sum2 = sum2 + (1.0 - tj) * (x[j] + tj + 1.0)^3
-        end
-
-        out[k] = x[k] + h * ((1.0 - tk) * sum1 + tk * sum2) / 2.0
-    end
-    nothing
-end
-
-n = 10
-x_sol = []
-x_start = zeros(n)
-for i in 1:n
-    x_start[i] = (i * (i - n - 1)) / (n + 1)^2
-end
-p10_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Discrete integral equation function")
-
-# ------------------------------------- Problem 11 -----------------------------------------
-function p11_f!(out, x, p = nothing)
-    n = length(x)
-    c_sum = sum(cos.(x))
-    for k in 1:n
-        out[k] = n - c_sum + k * (1.0 - cos(x[k])) - sin(x[k])
-    end
-    nothing
-end
-
-n = 10
-x_sol = []
-x_start = ones(n) / n
-p11_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Trigonometric function")
-
-# ------------------------------------- Problem 12 -----------------------------------------
-function p12_f!(out, x, p = nothing)
-    n = length(x)
-    sum1 = 0.0
-    for j in 1:n
-        sum1 += j * (x[j] - 1.0)
-    end
-    for k in 1:n
-        out[k] = x[k] - 1.0 + k * sum1 * (1.0 + 2.0 * sum1 * sum1)
-    end
-    nothing
-end
-
-n = 10
-x_sol = ones(n)
-x_start = zeros(n)
-for i in 1:n
-    x_start[i] = 1.0 - i / n
-end
-p12_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Variably dimensioned function")
-
-# ------------------------------------- Problem 13 -----------------------------------------
-function p13_f!(out, x, p = nothing)
-    n = length(x)
-    for k in 1:n
-        out[k] = (3.0 - 2.0 * x[k]) * x[k] + 1.0
-        if 1 < k
-            out[k] -= x[k - 1]
-        end
-        if k < n
-            out[k] -= 2.0 * x[k + 1]
-        end
-    end
-    nothing
-end
-
-n = 10
-x_sol = []
-x_start = ones(n) .* (-1.0)
-p13_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Broyden tridiagonal function")
-
-# ------------------------------------- Problem 14 -----------------------------------------
-function p14_f!(out, x, p = nothing)
-    n = length(x)
-    ml = 5
-    mu = 1
-    for k in 1:n
-        k1 = max(1, k - ml)
-        k2 = min(n, k + mu)
-
-        temp = 0.0
-        for j in k1:k2
-            if j != k
-                temp += x[j] * (1.0 + x[j])
-            end
-        end
-        out[k] = x[k] * (2.0 + 5.0 * x[k] * x[k]) + 1.0 - temp
-    end
-    nothing
-end
-
-n = 10
-x_sol = []
-x_start = ones(n) .* (-1.0)
-p14_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Broyden banded function")
-
-# ------------------------------------- Problem 15 -----------------------------------------
-function p15_f!(out, x, p = nothing)
-    out[1] = (x[1] * x[1] + x[2] * x[3]) - 0.0001
-    out[2] = (x[1] * x[2] + x[2] * x[4]) - 1.0
-    out[3] = (x[3] * x[1] + x[4] * x[3]) - 0.0
-    out[4] = (x[3] * x[2] + x[4] * x[4]) - 0.0001
-    nothing
-end
-
-n = 4
-x_sol = [0.01, 50.0, 0.0, 0.01]
-x_start = [1.0, 0.0, 0.0, 1.0]
-p15_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Hammarling 2 by 2 matrix square root problem")
-
-# ------------------------------------- Problem 16 -----------------------------------------
-function p16_f!(out, x, p = nothing)
-    out[1] = (x[1] * x[1] + x[2] * x[4] + x[3] * x[7]) - 0.0001
-    out[2] = (x[1] * x[2] + x[2] * x[5] + x[3] * x[8]) - 1.0
-    out[3] = x[1] * x[3] + x[2] * x[6] + x[3] * x[9]
-    out[4] = x[4] * x[1] + x[5] * x[4] + x[6] * x[7]
-    out[5] = (x[4] * x[2] + x[5] * x[5] + x[6] * x[8]) - 0.0001
-    out[6] = x[4] * x[3] + x[5] * x[6] + x[6] * x[9]
-    out[7] = x[7] * x[1] + x[8] * x[4] + x[9] * x[7]
-    out[8] = x[7] * x[2] + x[8] * x[5] + x[9] * x[8]
-    out[9] = (x[7] * x[3] + x[8] * x[6] + x[9] * x[9]) - 0.0001
-    nothing
-end
-
-n = 9
-x_sol = [0.01, 50.0, 0.0, 0.0, 0.01, 0.0, 0.0, 0.0, 0.01]
-x_start = [1.0, 0.0, 0.0, 0.0, 1.0, 0.0, 0.0, 0.0, 1.0]
-p16_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Hammarling 3 by 3 matrix square root problem")
-
-# ------------------------------------- Problem 17 -----------------------------------------
-function p17_f!(out, x, p = nothing)
-    out[1] = x[1] + x[2] - 3.0
-    out[2] = x[1]^2 + x[2]^2 - 9.0
-    nothing
-end
-
-n = 2
-x_sol = [0.0, 3.0]
-x_start = [1.0, 5.0]
-p17_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Dennis and Schnabel 2 by 2 example")
-
-# ------------------------------------- Problem 18 -----------------------------------------
-function p18_f!(out, x, p = nothing)
-    if x[1] != 0.0
-        out[1] = x[2]^2 * (1.0 - exp(-x[1] * x[1])) / x[1]
-    else
-        out[1] = 0.0
-    end
-    if x[2] != 0.0
-        out[2] = x[1] * (1.0 - exp(-x[2] * x[2])) / x[2]
-    else
-        out[2] = 0.0
-    end
-    nothing
-end
-
-n = 2
-x_sol = zeros(n)
-x_start = [2.0, 2.0]
-p18_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Sample problem 18")
-
-# ------------------------------------- Problem 19 -----------------------------------------
-function p19_f!(out, x, p = nothing)
-    out[1] = x[1] * (x[1]^2 + x[2]^2)
-    out[2] = x[2] * (x[1]^2 + x[2]^2)
-    nothing
-end
-
-n = 2
-x_sol = zeros(n)
-x_start = [3.0, 3.0]
-p19_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Sample problem 19")
-
-# ------------------------------------- Problem 20 -----------------------------------------
-function p20_f!(out, x, p = nothing)
-    out[1] = x[1] * (x[1] - 5.0)^2
-    nothing
-end
-
-n = 1
-x_sol = [5.0] # OR [0.0]...
-x_start = [1.0]
-p20_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Scalar problem f(x) = x(x - 5)^2")
-
-# ------------------------------------- Problem 21 -----------------------------------------
-function p21_f!(out, x, p = nothing)
-    out[1] = x[1] - x[2]^3 + 5.0 * x[2]^2 - 2.0 * x[2] - 13.0
-    out[2] = x[1] + x[2]^3 + x[2]^2 - 14.0 * x[2] - 29.0
-    nothing
-end
-
-n = 2
-x_sol = [5.0, 4.0]
-x_start = [0.5, -2.0]
-p21_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Freudenstein-Roth function")
-
-# ------------------------------------- Problem 22 -----------------------------------------
-function p22_f!(out, x, p = nothing)
-    out[1] = x[1] * x[1] - x[2] + 1.0
-    out[2] = x[1] - cos(0.5 * pi * x[2])
-    nothing
-end
-
-n = 2
-x_sol = [0.0, 1.0]
-x_start = [1.0, 0.0]
-p22_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Boggs function")
-
-# ------------------------------------- Problem 23 -----------------------------------------
-function p23_f!(out, x, p = nothing)
-    c = 0.9
-    out[1:n] = x[1:n]
-    μ = zeros(n)
-    for i in 1:n
-        μ[i] = (2 * i) / (2 * n)
-    end
-    for i in 1:n
-        s = 0.0
-        for j in 1:n
-            s = s + (μ[i] * x[j]) / (μ[i] + μ[j])
-        end
-        term = 1.0 - c * s / (2 * n)
-        out[i] -= 1.0 / term
-    end
-    nothing
-end
-
-n = 10
-x_sol = []
-x_start = ones(n)
-p23_dict = Dict("n" => n, "start" => x_start, "sol" => x_sol,
-    "title" => "Chandrasekhar function")
-
-# ----------------------------------- Solve problems ---------------------------------------
-problems = (p1_f!, p2_f!, p3_f!, p4_f!, p5_f!, p6_f!, p7_f!, p8_f!, p9_f!, p10_f!, p11_f!,
-    p12_f!, p13_f!, p14_f!, p15_f!, p16_f!, p17_f!, p18_f!, p19_f!, p20_f!, p21_f!,
-    p22_f!, p23_f!)
-dicts = (p1_dict, p2_dict, p3_dict, p4_dict, p5_dict, p6_dict, p7_dict, p8_dict, p9_dict,
-    p10_dict, p11_dict, p12_dict, p13_dict, p14_dict, p15_dict, p16_dict, p17_dict,
-    p18_dict, p19_dict, p20_dict, p21_dict, p22_dict, p23_dict)
+problems = NonlinearProblemLibrary.problems
+dicts = NonlinearProblemLibrary.dicts
 
 function test_on_library(problems, dicts, alg_ops, broken_tests, ϵ = 1e-5)
     for (idx, (problem, dict)) in enumerate(zip(problems, dicts))