Test nonlinear solvers

Author: EdoAlvarezR

Project.toml

@@ -20,6 +20,7 @@ CSV = "336ed68f-0bac-5ca0-87d4-7b16caf5d00b"
 DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
 SimpleNonlinearSolve = "727e6d20-b764-4bd8-a329-72de5adea6c7"
 NonlinearSolve = "8913a72c-1f9b-4ce2-8d82-65094dcecaec"
+NLsolve = "2774e3e8-f4cf-5e23-947b-6d7e65073b56"
 
 # [sources]
 # FastMultipole = {url = "https://github.com/byuflowlab/FastMultipole", rev = "d21d242d081ea8c9c19655097d83413b1f783d2d"}

examples/liftingline_weber.jl

@@ -25,6 +25,8 @@ include(joinpath(pnl.examples_path, "plotformat.jl"))
 import CSV
 import DataFrames: DataFrame
 
+import SIAMFANLEquations
+
 
 run_name        = "ll-weber"                    # Name of this run
 
@@ -112,8 +114,16 @@ sigmafactor     = 0.0                           # Dragging line amplification fa
 sigmaexponent   = 4.0                           # Dragging line amplification exponent (no effects if `sigmafactor==0.0`)
 
                                                 # Nonlinear solver
-solver          = pnl.SimpleNonlinearSolve.SimpleDFSane()              # Indifferent to initial guess, but somewhat not robust
+# solver          = pnl.SimpleNonlinearSolve.SimpleDFSane()              # Indifferent to initial guess, but somewhat not robust
 # solver        = pnl.SimpleNonlinearSolve.SimpleTrustRegion()         # Trust region needs a good initial guess, but it converges very reliably
+solver          = pnl.NonlinearSolve.SimpleBroyden()                      # <---- Wining
+# solver        = pnl.NonlinearSolve.NLsolveJL(method = :trust_region)
+# solver        = pnl.NonlinearSolve.SIAMFANLEquationsJL(method = :anderson)
+# solver        = pnl.NonlinearSolve.FixedPointAccelerationJL(; algorithm = :Newton)
+# solver        = pnl.NonlinearSolve.FixedPointAccelerationJL(; algorithm = :Newton)
+# solver = pnl.NonlinearSolve.NonlinearSolveQuasiNewton.Broyden(; autodiff = ADTypes.AutoForwardDiff(),  init_jacobian=Val(:true_jacobian))
+# solver = pnl.NonlinearSolve.NonlinearSolve.FastShortcutNLLSPolyalg(; autodiff = ADTypes.AutoForwardDiff(), vjp_autodiff = ADTypes.AutoForwardDiff(), jvp_autodiff = ADTypes.AutoForwardDiff()), # <--- WINING
+# solver = pnl.NonlinearSolve.NonlinearSolve.FastShortcutNonlinearPolyalg(; autodiff = ADTypes.AutoForwardDiff(), vjp_autodiff = ADTypes.AutoForwardDiff(), jvp_autodiff = ADTypes.AutoForwardDiff(), prefer_simplenonlinearsolve = Val(true)) 
 
 solver_optargs  = (; 
                     abstol = 1e-9,  

src/FLOWPanel.jl

@@ -24,6 +24,7 @@ import Dierckx
 import LinearAlgebra as LA
 import LinearAlgebra: I
 import Krylov
+import NLsolve
 import NonlinearSolve
 import SimpleNonlinearSolve
 import Requires: @require

src/liftingline/optimization.jl

@@ -132,9 +132,14 @@ function run_liftingline(;
         sigmaexponent   = 4.0,                          # Dragging line amplification exponent (no effects if `sigmafactor==0.0`)
 
                                                         # Nonlinear solver
-        # solver        = SimpleNonlinearSolve.SimpleDFSane(),             # Indifferent to initial guess, but somewhat not robust   <---- NOT COMPATIBLE WITH FORWARDDIFF
-        # solver        = SimpleNonlinearSolve.SimpleTrustRegion(),        # Trust region needs a good initial guess, but it converges very reliably
-        solver          = NonlinearSolve.SimpleBroyden(),                  # This seems to converge well while being compatible with ForwardDiff
+        # solver        = SimpleNonlinearSolve.SimpleDFSane(),             # Indifferent to initial guess, but somewhat not robust post stall   <---- NOT COMPATIBLE WITH FORWARDDIFF (but compatible with CSDA and optimization converges well)
+        # solver        = SimpleNonlinearSolve.SimpleTrustRegion(),        # Trust region needs a good initial guess, but solver converges very reliably post stall, not compatible with CSDA nor ForwardDiff
+        solver        = NonlinearSolve.SimpleBroyden(),                  # Optimization converges well while being compatible with ForwardDiff (not compatible with CSDA). EXTREMELY ROBUST ACROSS LINEAR, MILD STALL, AND DEEP STALL (it returns an answer, but it might be noise and not accurate)
+        # solver        = NonlinearSolve.NonlinearSolveQuasiNewton.Broyden(; autodiff = ADTypes.AutoForwardDiff(),  init_jacobian=Val(:true_jacobian)) # Also extremely robust across regions (it returns an answer, but it might be noise and not accurate)
+        # solver        = NonlinearSolve.NLsolveJL(method = :trust_region),# Optimization converges very well with ForwardDiff, not compatible with CSDA. Solver converges slowly but realibly in linear and mild stall regions, does not converge post stall
+        # solver        = NonlinearSolve.SIAMFANLEquationsJL(method = :newton, autodiff=ADTypes.AutoForwardDiff()), # Also robust in linear and mild stall regions, but much faster
+        # solver        = NonlinearSolve.NonlinearSolve.FastShortcutNLLSPolyalg(; autodiff = ADTypes.AutoForwardDiff(), vjp_autodiff = ADTypes.AutoForwardDiff(), jvp_autodiff = ADTypes.AutoForwardDiff()) # Very robust, but it can take a long time. Not ForwardDiff nor CSDA compatible
+        # solver        = NonlinearSolve.NonlinearSolve.FastShortcutNonlinearPolyalg(; autodiff = ADTypes.AutoForwardDiff(), vjp_autodiff = ADTypes.AutoForwardDiff(), jvp_autodiff = ADTypes.AutoForwardDiff(), prefer_simplenonlinearsolve = Val(true)) # 100% convergence, and super fast, but it might return the wrong solution. ForwardDiff compatible (though solver might be a bit noise, so set optimizer tol ~5e-5)
 
         solver_optargs  = (; 
                             abstol = 1e-13,  

test/runtests_liftingline_optimization_snopt.jl

@@ -10,6 +10,8 @@ import FLOWPanel: mean, norm, dot, cross
 
 import ForwardDiff: Dual, Partials, value, partials, jacobian
 import Snopt
+# import NLsolve
+import ADTypes
 import SNOW: minimize, Options, SNOPT, ForwardAD, ComplexStep
 
 try
@@ -123,9 +125,13 @@ v_lvl = 0
                                     sigmaexponent   = 4.0,                          # Dragging line amplification exponent (no effects if `sigmafactor==0.0`)
 
                                                                                     # Nonlinear solver
-                                    solver        = pnl.SimpleNonlinearSolve.SimpleDFSane(),         # Indifferent to initial guess, but somewhat not robust   <---- NOT COMPATIBLE WITH FORWARDDIFF
+                                    # solver        = pnl.SimpleNonlinearSolve.SimpleDFSane(),         # Indifferent to initial guess, but somewhat not robust   <---- NOT COMPATIBLE WITH FORWARDDIFF
                                     # solver        = pnl.SimpleNonlinearSolve.SimpleTrustRegion(),    # Trust region needs a good initial guess, but it converges very reliably
                                     # solver          = pnl.NonlinearSolve.SimpleBroyden(),              # This seems to converge well while being compatible with ForwardDiff
+                                    # solver        = pnl.NonlinearSolve.NLsolveJL(method = :trust_region), # Converges very well with ForwardDiff, not compatible with CSDA
+
+                                    # solver = pnl.NonlinearSolve.NonlinearSolve.FastShortcutNLLSPolyalg(; autodiff = ADTypes.AutoForwardDiff(), vjp_autodiff = ADTypes.AutoForwardDiff(), jvp_autodiff = ADTypes.AutoForwardDiff()),
+                                    solver = pnl.NonlinearSolve.NonlinearSolve.FastShortcutNonlinearPolyalg(; autodiff = ADTypes.AutoForwardDiff(), vjp_autodiff = ADTypes.AutoForwardDiff(), jvp_autodiff = ADTypes.AutoForwardDiff(), prefer_simplenonlinearsolve = Val(true)), # Robust, fast, and ForwardDiff compatible (though solver might be a bit noise, so set optimizer tol ~5e-5)
 
                                     solver_optargs  = (; 
                                                         abstol = 1e-12,             # <-- tight tolerance to converge derivatives
@@ -167,11 +173,13 @@ v_lvl = 0
 
         snopt_options = Dict(
                 "Major iterations limit" => 50,
+                "Major optimality tolerance" => 5e-5,
+                "Minor optimality tolerance" => 5e-5,
             )
         solver = SNOPT(options=snopt_options)
 
-        # options = Options(; solver, derivatives=ForwardAD())
-        options = Options(; solver, derivatives=ComplexStep())
+        options = Options(; solver, derivatives=ForwardAD())
+        # options = Options(; solver, derivatives=ComplexStep())
 
         model_cache = Dict()
 
@@ -254,7 +262,7 @@ v_lvl = 0
         @testset "Derivative verification" begin
             for (lbl, (f, dfdx), (atol1, atol2)) in (
                                     ("Finite Difference", (f_diff, dfdx_diff), (1e-10, 1e-6)),
-                                    ("CSDA", (f_csda, dfdx_csda), (1e-10, 5e-4)),
+                                    # ("CSDA", (f_csda, dfdx_csda), (1e-10, 5e-4)),
                                     ("Dual", (f_dual, dfdx_dual), (1e-10, 1e-6)),
                                 )