declare sqrt_xi_1 as a local variable and add unit test for LM with R2DH as subsolver

Author: MohamedLaghdafHABIBOULLAH

src/LM_alg.jl

@@ -102,11 +102,10 @@ function LM(
   xkn = similar(xk)
 
   local ξ1
+  local sqrt_ξ1_νInv
   k = 0
   Fobj_hist = zeros(maxIter)
   Hobj_hist = zeros(maxIter)
-  R = eltype(xk)
-  sqrt_ξ1_νInv = one(R)
   Complex_hist = zeros(Int, maxIter)
   Grad_hist = zeros(Int, maxIter)
   Resid_hist = zeros(Int, maxIter)

test/test_bounds.jl

@@ -66,30 +66,68 @@ for (h, h_name) ∈ ((NormL0(λ), "l0"), (NormL1(λ), "l1"))
       @test h(out.solution) == out.solver_specific[:Hhist][end]
       @test out.status == :first_order
     end
-    @testset "bpdn-with-bounds-ls-$(solver_name)-$(h_name)-TRDH" begin
-      x0 = zeros(bpdn_nls2.meta.nvar)
-      args = solver_sym == :LM ? () : (NormLinf(1.0),)
-      @test has_bounds(bpdn_nls2)
-      out = solver(
-        bpdn_nls2,
-        h,
-        args...,
-        options,
-        x0 = x0,
-        subsolver = TRDH,
-        subsolver_options = subsolver_options,
-      )
-      @test typeof(out.solution) == typeof(bpdn_nls2.meta.x0)
-      @test length(out.solution) == bpdn_nls2.meta.nvar
-      @test typeof(out.solver_specific[:Fhist]) == typeof(out.solution)
-      @test typeof(out.solver_specific[:Hhist]) == typeof(out.solution)
-      @test typeof(out.solver_specific[:SubsolverCounter]) == Array{Int, 1}
-      @test typeof(out.dual_feas) == eltype(out.solution)
-      @test length(out.solver_specific[:Fhist]) == length(out.solver_specific[:Hhist])
-      @test length(out.solver_specific[:Fhist]) == length(out.solver_specific[:SubsolverCounter])
-      @test obj(bpdn_nls2, out.solution) == out.solver_specific[:Fhist][end]
-      @test h(out.solution) == out.solver_specific[:Hhist][end]
-      @test out.status == :first_order
-    end
+  end
+end
+
+# LMTR with TRDH as subsolver
+for (h, h_name) ∈ ((NormL0(λ), "l0"), (NormL1(λ), "l1"))
+  @testset "bpdn-with-bounds-ls-LMTR-$(h_name)-TRDH" begin
+    x0 = zeros(bpdn_nls2.meta.nvar)
+    @test has_bounds(bpdn_nls2)
+    LMTR_out = LMTR(
+      bpdn_nls2,
+      h,
+      NormLinf(1.0),
+      options,
+      x0 = x0,
+      subsolver = TRDH,
+      subsolver_options = subsolver_options,
+    )
+    @test typeof(LMTR_out.solution) == typeof(bpdn_nls2.meta.x0)
+    @test length(LMTR_out.solution) == bpdn_nls2.meta.nvar
+    @test typeof(LMTR_out.solver_specific[:Fhist]) == typeof(LMTR_out.solution)
+    @test typeof(LMTR_out.solver_specific[:Hhist]) == typeof(LMTR_out.solution)
+    @test typeof(LMTR_out.solver_specific[:SubsolverCounter]) == Array{Int, 1}
+    @test typeof(LMTR_out.dual_feas) == eltype(LMTR_out.solution)
+    @test length(LMTR_out.solver_specific[:Fhist]) == length(LMTR_out.solver_specific[:Hhist])
+    @test length(LMTR_out.solver_specific[:Fhist]) ==
+          length(LMTR_out.solver_specific[:SubsolverCounter])
+    @test length(LMTR_out.solver_specific[:Fhist]) == length(LMTR_out.solver_specific[:NLSGradHist])
+    @test LMTR_out.solver_specific[:NLSGradHist][end] ==
+          bpdn_nls2.counters.neval_jprod_residual + bpdn_nls2.counters.neval_jtprod_residual - 1
+    @test obj(bpdn_nls2, LMTR_out.solution) == LMTR_out.solver_specific[:Fhist][end]
+    @test h(LMTR_out.solution) == LMTR_out.solver_specific[:Hhist][end]
+    @test LMTR_out.status == :first_order
+  end
+end
+
+# LM with R2DH as subsolver
+for (h, h_name) ∈ ((NormL0(λ), "l0"), )
+  @testset "bpdn-with-bounds-ls-LM-$(h_name)-R2DH" begin
+    x0 = zeros(bpdn_nls2.meta.nvar)
+    @test has_bounds(bpdn_nls2)
+    LM_out = LM(
+      bpdn_nls2,
+      h,
+      options,
+      x0 = x0,
+      subsolver = R2DH,
+      subsolver_options = subsolver_options,
+    )
+    @test typeof(LM_out.solution) == typeof(bpdn_nls2.meta.x0)
+    @test length(LM_out.solution) == bpdn_nls2.meta.nvar
+    @test typeof(LM_out.solver_specific[:Fhist]) == typeof(LM_out.solution)
+    @test typeof(LM_out.solver_specific[:Hhist]) == typeof(LM_out.solution)
+    @test typeof(LM_out.solver_specific[:SubsolverCounter]) == Array{Int, 1}
+    @test typeof(LM_out.dual_feas) == eltype(LM_out.solution)
+    @test length(LM_out.solver_specific[:Fhist]) == length(LM_out.solver_specific[:Hhist])
+    @test length(LM_out.solver_specific[:Fhist]) ==
+          length(LM_out.solver_specific[:SubsolverCounter])
+    @test length(LM_out.solver_specific[:Fhist]) == length(LM_out.solver_specific[:NLSGradHist])
+    @test LM_out.solver_specific[:NLSGradHist][end] ==
+          bpdn_nls2.counters.neval_jprod_residual + bpdn_nls2.counters.neval_jtprod_residual - 1
+    @test obj(bpdn_nls2, LM_out.solution) == LM_out.solver_specific[:Fhist][end]
+    @test h(LM_out.solution) == LM_out.solver_specific[:Hhist][end]
+    @test LM_out.status == :first_order
   end
 end