!!!!! newton -> solve(, ::Newton,)

examples/SH2d-fronts-cuda.jl

@@ -118,15 +118,15 @@ prob = BK.BifurcationProblem(F_shfft, AF(sol0), par, (@optic _.l) ;
     record_from_solution = (x, p; k...) -> norm(x))
 
 opt_new = NewtonPar(verbose = true, tol = 1e-6, linsolver = L, eigsolver = Leig)
-sol_hexa = @time newton(prob, opt_new, normN = norminf);
+sol_hexa = @time BK.solve(prob, Newton(), opt_new, normN = norminf);
 println("--> norm(sol) = ", norminf(sol_hexa.u))
 
 plotsol(sol_hexa.u)
 ####################################################################################################
 deflationOp = DeflationOperator(2, 1.0, [sol_hexa.u])
 
 opt_new = @set opt_new.max_iterations = 250
-outdef = @time newton(re_make(prob, u0 = 0.4 .* sol_hexa.u .* AF([exp(-1(x+0lx)^2/25) for x in X, y in Y])),
+outdef = @time BK.solve(re_make(prob, u0 = 0.4 .* sol_hexa.u .* AF([exp(-1(x+0lx)^2/25) for x in X, y in Y])),
             deflationOp, opt_new, normN = x -> maximum(abs, x))
 println("--> norm(sol) = ", norm(outdef.u))
 plotsol(outdef.u) |> display

examples/SH2d-fronts.jl

@@ -59,7 +59,7 @@ optnew = NewtonPar(verbose = true, tol = 1e-8, max_iterations = 20)
 
 prob = BifurcationProblem(F_sh, vec(sol0), par, (@optic _.l); J = dF_sh, plot_solution = (x, p; kwargs...) -> (plotsol!((x); label="", kwargs...)),record_from_solution = (x, p; k...) -> (n2 = norm(x), n8 = norm(x, 8)), d2F=d2F_sh, d3F=d3F_sh)
 # optnew = NewtonPar(verbose = true, tol = 1e-8, max_iterations = 20, eigsolver = EigArpack(0.5, :LM))
-sol_hexa = @time newton(prob, optnew)
+sol_hexa = @time BK.solve(prob, Newton(), optnew)
 println("--> norm(sol) = ", norm(sol_hexa.u, Inf64))
 plotsol(sol_hexa.u)
 
@@ -72,7 +72,7 @@ deflationOp = DeflationOperator(2, 1.0, [sol_hexa.u])
 
 optnew = @set optnew.max_iterations = 250
 optnewd = @set optnew.max_iterations = 250
-outdef = @time newton(
+outdef = @time BK.solve(
         (@set prob.u0 = 0.4vec(sol_hexa.u) .* vec([exp(-1(x+lx)^2/25) for x in X, y in Y])), deflationOp,
         # 0.4vec(sol_hexa) .* vec([1 .- exp(-1(x+lx)^2/55) for x in X, y in Y]),
         optnewd)
@@ -129,7 +129,7 @@ end
 prob2 = @set prob.VF.J = (u, p) -> (du -> dF_sh2(du, u, p))
 @reset prob2.u0 = vec(sol0)
 
-sol_hexa = @time newton(prob2, @set optnew.linsolver = ls)
+sol_hexa = @time BK.solve(prob2, Newton(), @set optnew.linsolver = ls)
 println("--> norm(sol) = ", norm(sol_hexa.u, Inf64))
 plotsol(sol_hexa.u)
 ###################################################################################################

examples/SH3d.jl

@@ -123,7 +123,7 @@ prob = BK.BifurcationProblem(F_sh, AF(vec(sol0)), par, (@optic _.l),
 optnew = NewtonPar(verbose = true, tol = 1e-8, max_iterations = 20, linsolver = @set ls.verbose = 0)
 @reset optnew.eigsolver = eigSH3d
 # @reset optnew.linsolver = DefaultLS()
-sol_hexa = @time newton(prob, optnew)
+sol_hexa = @time BK.solve(prob, Newton(), optnew)
 println("--> norm(sol) = ", norm(sol_hexa.u, Inf64))
 
 contour3dMakie(sol0)

examples/SHpde_snaking.jl

@@ -35,7 +35,7 @@ prob = BifurcationProblem(R_SH, sol0, parSH, (@optic _.λ); J = Jac_sp,
     )
 ####################################################################################################
 optnew = NewtonPar(tol = 1e-12)
-sol1 = newton(prob, optnew)
+sol1 = BK.solve(prob, Newton(), optnew)
 Plots.plot(X, sol1.u)
 
 opts = ContinuationPar(dsmin = 0.0001, dsmax = 0.01, ds = 0.01,

examples/cGL2d-shooting.jl

@@ -118,7 +118,7 @@ prob_sp = SplitODEProblem(f1, f2, sol0_f, (0.0, 120.0), @set par_cgl.r = 1.2; re
 prob = ODEProblem(Fcgl!, sol0_f, (0.0, 120.0), (@set par_cgl.r = 1.2))#, jac = Jcgl, jac_prototype = Jcgl(sol0_f, par_cgl))
 ####################################################################################################
 # sol = @time solve(prob, Vern9(); abstol=1e-14, reltol=1e-14)
-sol = @time solve(prob_sp, ETDRK2(krylov=true); abstol=1e-14, reltol=1e-14, dt = 0.1) #1.78s
+sol = @time DifferentialEquations.solve(prob_sp, ETDRK2(krylov=true); abstol=1e-14, reltol=1e-14, dt = 0.1) #1.78s
 # sol = @time solve(prob, LawsonEuler(krylov=true, m=50); abstol=1e-14, reltol=1e-14, dt = 0.1)
 # sol = @time solve(prob_sp, CNAB2(linsolve=LinSolveGMRES()); abstol=1e-14, reltol=1e-14, dt = 0.03)
 

examples/cGL2d.jl

@@ -96,7 +96,7 @@ prob = BK.BifurcationProblem(Fcgl!, vec(sol0), par_cgl, (@optic _.r); J = Jcgl)
 eigls = EigArpack(1.0, :LM)
 # eigls = eig_MF_KrylovKit(tol = 1e-8, dim = 60, x₀ = rand(ComplexF64, Nx*Ny), verbose = 1)
 opt_newton = NewtonPar(tol = 1e-9, verbose = true, eigsolver = eigls, max_iterations = 20)
-out = @time newton(prob, opt_newton, normN = norminf)
+out = @time solve(prob, Newton(), opt_newton, normN = norminf)
 ####################################################################################################
 opts_br = ContinuationPar(dsmin = 0.001, dsmax = 0.15, ds = 0.001, p_max = 2.5, detect_bifurcation = 3, nev = 9, plot_every_step = 50, newton_options = (@set opt_newton.verbose = false), max_steps = 1060, n_inversion = 6)
 br = @time continuation(prob, PALC(), opts_br, verbosity = 0)

examples/carrier.jl

@@ -43,7 +43,7 @@ prob = BK.BifurcationProblem(F_carr, zeros(N), par_car, (@optic _.ϵ);
     record_from_solution)
 
 optnew = NewtonPar(tol = 1e-8, verbose = true)
-out = @time newton(prob, optnew, normN = norminf)
+out = @time BK.solve(prob, Newton(), optnew, normN = norminf)
 plot(out.u, label = "Solution")
 
 optcont = ContinuationPar(dsmin = 0.001, dsmax = 0.05, ds= -0.01, p_min = 0.05, newton_options = NewtonPar(tol = 1e-8, max_iterations = 20, verbose = false), nev = 40, n_inversion = 6)
@@ -67,7 +67,7 @@ function perturbsol(sol, p, id)
     return sol .+ solp .* sol0
 end
 
-outdef1 = @time newton(
+outdef1 = @time BK.solve(
     (@set prob.u0 = perturbsol(-out.u, 0, 0)), deflationOp,
     # perturbsol(deflationOp[1],0,0), par_def,
     optdef;

examples/chan-af.jl

@@ -64,7 +64,7 @@ prob = BifurcationProblem(F_chan, sol0, par_af, (@optic _.α);
             plot_solution = (x, p; kwargs...) -> plot!(x; label = "l = $(length(x))", kwargs...))
 
 optnew = NewtonPar(tol = 1e-12, verbose = true)
-sol = @time BK.newton(prob, optnew, normN = norminf)
+sol = @time BK.solve(prob, Newton(), optnew, normN = norminf)
 
 plot(sol.u, label="Solution")
 

examples/chan.jl

@@ -27,7 +27,7 @@ prob = BifurcationProblem(F_chan, sol0, par, (@optic _.α);
 
 optnewton = NewtonPar(tol = 1e-8, verbose = true)
 # ca fait dans les 63.59k Allocations
-sol = @time newton( prob, optnewton)
+sol = @time BK.solve(prob, Newton(), optnewton)
 
 optscont = ContinuationPar(dsmin = 0.01, dsmax = 0.5, ds= 0.1, p_max = 4.25, nev = 5, detect_fold = true, plot_every_step = 10, newton_options = NewtonPar(max_iterations = 10, tol = 1e-9, verbose = false), max_steps = 150)
 
@@ -42,16 +42,16 @@ deflationOp = DeflationOperator(2, 1.0, [sol.u])
 
 optdef = NewtonPar(optnewton; tol = 1e-10, max_iterations = 500)
 
-outdef1 = newton(re_make(prob; u0 = sol.u .* (1 .+ 0.01*rand(n))), deflationOp, optdef)
-outdef1 = newton(re_make(prob; u0 = sol.u .* (1 .+ 0.01*rand(n))), deflationOp, optdef, Val(:autodiff))
+outdef1 = BK.solve(re_make(prob; u0 = sol.u .* (1 .+ 0.01*rand(n))), deflationOp, optdef)
+outdef1 = BK.solve(re_make(prob; u0 = sol.u .* (1 .+ 0.01*rand(n))), deflationOp, optdef, Val(:autodiff))
 
 plot(sol.u, label="newton")
 plot!(sol0, label="init guess")
 plot!(outdef1.u, label="deflation-1")
 
 #save newly found point to look for new ones
 push!(deflationOp, outdef1.u)
-outdef2 = @time newton((@set prob.u0 = sol0), deflationOp, optdef; callback = BK.cbMaxNorm(1e5))
+outdef2 = @time BK.solve((@set prob.u0 = sol0), deflationOp, optdef; callback = BK.cbMaxNorm(1e5))
 plot!(outdef2.u, label="deflation-2")
 #################################################################################################### Continuation of the Fold Point using minimally augmented formulation
 optscont = (@set optscont.newton_options = setproperties(optscont.newton_options; verbose = true, tol = 1e-10))
@@ -110,7 +110,7 @@ P[1,1:2] .= [1, 0.];P[end,end-1:end] .= [0, 1.]
 
 lsp = GMRESIterativeSolvers(reltol = 1e-5, N = length(sol.u), restart = 20, maxiter=10, Pl = lu(P))
 optnewton_mf = NewtonPar(tol = 1e-9, verbose = false, linsolver = lsp)
-out_mf = @time newton(prob2, @set optnewton_mf.verbose = true)
+out_mf = @time BK.solve(prob2, Newton(), @set optnewton_mf.verbose = true)
 
 plot(brmf, color=:red)
 
@@ -120,7 +120,7 @@ brmf = @time continuation(prob2, PALC(tangent = Bordered(), bls = BorderingBLS(l
 plot(brmf,color=:blue)
 
 alg = brmf.alg
-brmf = @time continuation(prob2, MoorePenrose(tangent = alg, directLS = false), (@set opts_cont_mf.newton_options = optnewton_mf))
+brmf = @time continuation(prob2, MoorePenrose(tangent = alg, method = BifurcationKit.iterative), (@set opts_cont_mf.newton_options = optnewton_mf))
 
 plot(brmf,color=:blue)
 

examples/codim2PO-sh-mf.jl

@@ -31,9 +31,9 @@ using DifferentialEquations
 prob_de = ODEProblem(Pop!, z0, (0,600.), par_pop)
 alg = Rodas5()
 # alg = Vern9()
-sol = solve(prob_de, alg)
+sol = DifferentialEquations.solve(prob_de, alg)
 prob_de = ODEProblem(Pop!, sol.u[end], (0,5.), par_pop, reltol = 1e-10, abstol = 1e-12)
-sol = solve(prob_de, Rodas5())
+sol = DifferentialEquations.solve(prob_de, Rodas5())
 
 plot(sol)
 ################################################################################
@@ -67,7 +67,7 @@ probsh, cish = generate_ci_problem( ShootingProblem(M=3), prob, prob_de, sol, 2.
 ######
 function flow(x0, prob0, tm, p = prob0.p)
     prob = remake(prob0, u0 = x0, tspan = (0, tm), p = p)
-    sol = solve(prob, Rodas5())
+    sol = DifferentialEquations.solve(prob, Rodas5())
     return sol[end]
 end
 
@@ -97,7 +97,7 @@ AD.pullback_function(AD.FiniteDifferencesBackend(), z -> probsh(z, getparams(pro
 lspo = GMRESIterativeSolvers(verbose = false, N = length(cish), abstol = 1e-12, reltol = 1e-10)
 eigpo = EigKrylovKit(x₀ = rand(4))
 optnpo = NewtonPar(verbose = true, linsolver = lspo, eigsolver = eigpo)
-solpo = newton(probsh, cish, optnpo)
+solpo = BK.newton(probsh, cish, optnpo)
 
 _sol = BK.get_periodic_orbit(probsh, solpo.u, sol.prob.p)
 plot(_sol.t, _sol[1:2,:]')
@@ -181,8 +181,8 @@ plot!(pd_po_sh, vars = (:ϵ, :b0), branchlabel = "PD")
 #####
 # find the PD NS case
 par_pop2 = @set par_pop.b0 = 0.45
-sol2 = solve(remake(prob_de, p = par_pop2, u0 = [0.1,0.1,1,0], tspan=(0,1000)), Rodas5())
-sol2 = solve(remake(sol2.prob, tspan = (0,10), u0 = sol2[end]), Rodas5())
+sol2 = DifferentialEquations.solve(remake(prob_de, p = par_pop2, u0 = [0.1,0.1,1,0], tspan=(0,1000)), Rodas5())
+sol2 = DifferentialEquations.solve(remake(sol2.prob, tspan = (0,10), u0 = sol2[end]), Rodas5())
 plot(sol2, xlims= (8,10))
 
 probshns, ci = generate_ci_problem( ShootingProblem(M=3), re_make(prob, params = sol2.prob.p), remake(prob_de, p = par_pop2), sol2, 1.; alg = Rodas5(),

examples/codim2PO.jl

@@ -30,9 +30,9 @@ opts_br = ContinuationPar(p_min = 0., p_max = 20.0, ds = 0.002, dsmax = 0.01, n_
 using DifferentialEquations
 prob_de = ODEProblem(Pop!, z0, (0,600.), par_pop)
 alg = Rodas5()
-sol = solve(prob_de, alg)
+sol = DifferentialEquations.solve(prob_de, alg)
 prob_de = ODEProblem(Pop!, sol.u[end], (0,5.), par_pop, reltol = 1e-8, abstol = 1e-10)
-sol = solve(prob_de, Rodas5())
+sol = DifferentialEquations.solve(prob_de, Rodas5())
 
 plot(sol)
 ################################################################################
@@ -72,7 +72,7 @@ probtrap, ci = generate_ci_problem(PeriodicOrbitTrapProblem(M = 150;  jacobian =
 plot(sol)
 probtrap(ci, prob.params) |> plot
 
-solpo = newton(probtrap, ci, NewtonPar(verbose = true))
+solpo = BK.newton(probtrap, ci, NewtonPar(verbose = true))
 
 _sol = BK.get_periodic_orbit(probtrap, solpo.u,1)
 plot(_sol.t, _sol[1:2,:]')
@@ -142,7 +142,7 @@ probcoll, ci = generate_ci_problem(PeriodicOrbitOCollProblem(30, 3), prob, sol,
 plot(sol)
 probcoll(ci, prob.params) |> plot
 
-solpo = newton(probcoll, ci, NewtonPar(verbose = true));
+solpo = BK.newton(probcoll, ci, NewtonPar(verbose = true));
 
 _sol = BK.get_periodic_orbit(probcoll, solpo.u,1)
 plot(_sol.t, _sol[1:2,:]')

examples/mittleman.jl

@@ -84,7 +84,7 @@ eigls = EigArpack(20.5, :LM)
 # eigls = EigKrylovKit(dim = 70)
 # eigls = EigArpack()
 opt_newton = NewtonPar(tol = 1e-8, verbose = true, eigsolver = eigls, max_iterations = 20)
-sol = newton(prob, opt_newton, normN = norminf)
+sol = BK.solve(prob, Newton(), opt_newton, normN = norminf)
 
 plotsol(sol.u)
 ####################################################################################################

examples/pd-1d.jl

@@ -68,7 +68,7 @@ probBif = BifurcationProblem(Fbr, solc0, par_br, (@optic _.C) ;J = Jbr,
 eigls = EigArpack(0.5, :LM)
 optnewton = NewtonPar(eigsolver = eigls, verbose=true, max_iterations = 3200, tol=1e-9)
 
-out = @time newton(probBif, optnewton, normN = norminf)
+out = @time BK.solve(probBif, Newton(), optnewton, normN = norminf)
 plot();plot!(X,out.u[1:N]);plot!(X,solc0[1:N], label = "sol0",line=:dash)
 
 optcont = ContinuationPar(dsmax = 0.051, ds = -0.001, p_min = -1.8, detect_bifurcation = 3, nev = 21, plot_every_step = 50, newton_options = optnewton, max_steps = 370, n_inversion = 10, max_bisection_steps = 25)
@@ -132,10 +132,10 @@ f1 = DiffEqArrayOperator(par_br.Δ)
 f2 = NL!
 prob_sp = SplitODEProblem(f1, f2, solc0, (0.0, 280.0), par_br_hopf)
 
-sol = @time solve(prob_sp, ETDRK2(krylov=true); abstol=1e-14, reltol=1e-14, dt = 0.1, progress = true)
+sol = @time DifferentialEquations.solve(prob_sp, ETDRK2(krylov=true); abstol=1e-14, reltol=1e-14, dt = 0.1, progress = true)
 
 prob_ode = ODEProblem(Fbr, solc0, (0.0, 280.0), par_br_hopf)
-sol = @time solve(prob_ode, Rodas4P(); abstol=1e-14, reltol=1e-7, dt = 0.1, progress = true)
+sol = @time DifferentialEquations.solve(prob_ode, Rodas4P(); abstol=1e-14, reltol=1e-7, dt = 0.1, progress = true)
 orbitsection = Array(sol[:,[end]])
 # orbitsection = orbitguess[:, 1]
 
@@ -153,7 +153,7 @@ ls = GMRESIterativeSolvers(reltol = 1e-7, N = length(initpo), maxiter = 50, verb
 # ls = GMRESKrylovKit(verbose = 0, dim = 200, atol = 1e-9, rtol = 1e-5)
 optn = NewtonPar(verbose = true, tol = 1e-9,  max_iterations = 20, linsolver = ls)
 # deflationOp = BK.DeflationOperator(2 (x,y) -> dot(x[1:end-1], y[1:end-1]),1.0, [outpo])
-outposh = @time newton(probSh, initpo, optn; normN = norminf)
+outposh = @time BK.newton(probSh, initpo, optn; normN = norminf)
 BK.converged(outposh) && printstyled(color=:red, "--> T = ", outposh.u[end], ", amplitude = ", BK.getamplitude(probSh, outposh.u, par_br_hopf; ratio = 2),"\n")
 
 plot(initpo[1:end-1], label = "Init guess"); plot!(outposh.u[1:end-1], label = "sol")
@@ -183,7 +183,7 @@ f2 = NL!
 prob_sp = SplitODEProblem(f1, f2, solc0, (0.0, 300.0), par_br_pd; abstol=1e-14, reltol=1e-14, dt = 0.01)
 # solution close to the PD point.
 
-solpd = @time solve(prob_sp, ETDRK2(krylov=true), progress = true)
+solpd = @time DifferentialEquations.solve(prob_sp, ETDRK2(krylov=true), progress = true)
     # heatmap(sol.t, X, sol[1:N,:], color=:viridis, xlim=(20,280.0))
 
 orbitsectionpd = Array(solpd[:,end-100])
@@ -201,7 +201,7 @@ ls = GMRESIterativeSolvers(reltol = 1e-7, N = length(initpo_pd), maxiter = 50, v
 # ls = GMRESKrylovKit(verbose = 0, dim = 200, atol = 1e-9, rtol = 1e-5)
 optn = NewtonPar(verbose = true, tol = 1e-9,  max_iterations = 12, linsolver = ls)
 # deflationOp = BK.DeflationOperator(2 (x,y) -> dot(x[1:end-1], y[1:end-1]),1.0, [outpo])
-outposh_pd = @time newton(BK.set_params_po(probSh,par_br_pd), initpo_pd, optn;
+outposh_pd = @time BK.newton(BK.set_params_po(probSh,par_br_pd), initpo_pd, optn;
         # callback = (state; kwargs...) -> (@show state.x[end];true),
         normN = norminf)
 BK.converged(outposh_pd) && printstyled(color=:red, "--> T = ", outposh_pd.u[end], ", amplitude = ", BK.getamplitude(probSh, outposh_pd.u, (@set par_br.C = -0.86); ratio = 2),"\n")

src/BifurcationKit.jl

@@ -136,7 +136,7 @@ module BifurcationKit
     export Natural, PALC, Multiple, Secant, Bordered, DefCont, Polynomial, MoorePenrose, MoorePenroseLS, AutoSwitch
 
     # newton methods
-    export NewtonPar, newton, newton_palc, newton_hopf, NonLinearSolution
+    export NewtonPar, Newton, solve, newton, newton_palc, newton_hopf, NonLinearSolution
 
     # continuation methods
     export ContinuationPar, ContResult, continuation, continuation!, continuation_fold, continuation_hopf, continuation_potrap, eigenvec, eigenvals, get_solx, get_solp, bifurcation_points, SpecialPoint

src/Continuation.jl

@@ -349,7 +349,7 @@ function Base.iterate(it::ContIterable; _verbosity = it.verbosity)
     verbose && printstyled("━"^18*"  INITIAL GUESS   "*"━"^18, bold = true, color = :magenta)
 
     # we pass additional kwargs to newton so that it is sent to the newton callback
-    sol₀ = newton(prob, newton_options; 
+    sol₀ = solve(prob, Newton(), newton_options; 
                     normN = it.normC,
                     callback = callback(it),
                     iterationC = 0,
@@ -363,7 +363,8 @@ function Base.iterate(it::ContIterable; _verbosity = it.verbosity)
     verbose && println("──▶ parameter = ", p₀, ", initial step")
     verbose && printstyled("\n"*"━"^18*" INITIAL TANGENT  "*"━"^18, bold = true, color = :magenta)
 
-    sol₁ = newton(re_make(prob; params = setparam(it, p₀ + ds / η), u0 = sol₀.u),
+    sol₁ = solve(re_make(prob; params = setparam(it, p₀ + ds / η), u0 = sol₀.u),
+                            Newton(),
                             newton_options; 
                             normN = it.normC,
                             callback = callback(it),

src/DeflatedContinuation.jl

@@ -113,7 +113,7 @@ function updatebranch!(dcIter::DefContIterable,
 
     getpredictor!(state, it)
     pbnew = re_make(it.prob; u0 = _copy(getx(state)), params = setparam(it, current_param))
-    sol1 = newton(pbnew, defOp, it.contparams.newton_options, alg.jacobian;
+    sol1 = solve(pbnew, defOp, it.contparams.newton_options, alg.jacobian;
                     normN = it.normC,
                     callback = it.callback_newton,
                     iterationC = step,
@@ -253,7 +253,7 @@ function deflatedContinuation(dcIter::DefContIterable,
         prob_df = re_make(contIt.prob;
                             u0 = alg.perturb_solution(u0, _p, _idb),
                             params = set(par, lens, _p))
-        soln = newton(prob_df, deflationOp,
+        soln = solve(prob_df, deflationOp,
                 setproperties(optnewton; max_iterations = alg.max_iter_defop);
                 normN = contIt.normC,
                 callback = contIt.callback_newton,