reformat all code and change some dependencies

.github/workflows/Format.yml

@@ -0,0 +1,23 @@
+name: Format
+on: [pull_request]
+
+jobs:
+  format:
+    runs-on: ubuntu-latest
+    steps:
+      - uses: actions/checkout@v4
+      - uses: julia-actions/setup-julia@latest
+        with:
+          version: 1
+      - name: Format code
+        run: |
+          using Pkg
+          Pkg.add(; name="JuliaFormatter", uuid="98e50ef6-434e-11e9-1051-2b60c6c9e899")
+          using JuliaFormatter
+          format("."; verbose=true)
+        shell: julia --color=yes {0}
+      - uses: reviewdog/action-suggester@v1
+        if: github.event_name == 'pull_request'
+        with:
+          tool_name: JuliaFormatter
+

Project.toml

@@ -1,10 +1,9 @@
 name = "HOODESolver"
 uuid = "217528c2-21c4-4222-9975-97f0b0cbd7c1"
 authors = ["Yves Mocquard", "Nicolas Crouseilles", "Pierre Navaro"]
-version = "0.2.7"
+version = "0.2.8"
 
 [deps]
-DoubleFloats = "497a8b3b-efae-58df-a0af-a86822472b78"
 FFTW = "7a1cc6ca-52ef-59f5-83cd-3a7055c09341"
 GenericSchur = "c145ed77-6b09-5dd9-b285-bf645a82121e"
 LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e"
@@ -20,17 +19,24 @@ Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
 DoubleFloats = "1"
 FFTW = "^1"
 GenericSchur = "0.4, 0.5"
+LinearAlgebra = "1"
+OrdinaryDiffEq = "6.89.0"
+Plots = "1"
 Polynomials = "^1, 2.0, 3, 4"
+Random = "1"
 RecursiveArrayTools = "2, 3"
 Reexport = "^0.2, 1.0"
 SciMLBase = "1, 2"
+SparseArrays = "1"
 Statistics = "1"
+Test = "1"
 julia = "^1"
 
 [extras]
-DifferentialEquations = "0c46a032-eb83-5123-abaf-570d42b7fbaa"
+OrdinaryDiffEq = "1dea7af3-3e70-54e6-95c3-0bf5283fa5ed"
+DoubleFloats = "497a8b3b-efae-58df-a0af-a86822472b78"
 Plots = "91a5bcdd-55d7-5caf-9e0b-520d859cae80"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 
 [targets]
-test = ["Test", "DifferentialEquations", "Plots"]
+test = ["DoubleFloats", "Test", "OrdinaryDiffEq", "Plots"]

benchmark/bench_particle.jl

@@ -7,30 +7,38 @@ using HOODESolver
 suite = BenchmarkGroup()
 
 function fct(u, p, t)
-    s1, c1 = sincos(u[1]/2)
+    s1, c1 = sincos(u[1] / 2)
     s2, c2 = sincos(u[2])
     s3, c3 = sincos(u[3])
-    return [0, 0, u[6], c1*s2*s3/2, s1*c2*s3, s1*s2*c3]
+    return [0, 0, u[6], c1 * s2 * s3 / 2, s1 * c2 * s3, s1 * s2 * c3]
 end
 
 Random.seed!(561909)
-u0 = 2rand(BigFloat,6)-ones(BigFloat,6)
-epsilon=big"1e-7"
-ordmax=9
-A = [0 0 0  1 0 0; 
-     0 0 0  0 1 0;
-     0 0 0  0 0 0; 
-     0 0 0  0 1 0; 
-     0 0 0 -1 0 0; 
-     0 0 0  0 0 0]
+u0 = 2rand(BigFloat, 6) - ones(BigFloat, 6)
+epsilon = big"1e-7"
+ordmax = 9
+A = [
+    0 0 0 1 0 0
+    0 0 0 0 1 0
+    0 0 0 0 0 0
+    0 0 0 0 1 0
+    0 0 0 -1 0 0
+    0 0 0 0 0 0
+]
 
 t_0 = big"0.0"
 t_max = big"1.0"
 prob = HOODEProblem(fct, u0, (t_0, t_max), missing, A, epsilon)
 nb = 10^3
 n_tau = 32
-suite["solve"] = @benchmarkable solve(prob, nb_t=nb, order=ordmax, 
-                                      getprecision=false, nb_tau=n_tau, dense=false)
+suite["solve"] = @benchmarkable solve(
+    prob,
+    nb_t = nb,
+    order = ordmax,
+    getprecision = false,
+    nb_tau = n_tau,
+    dense = false,
+)
 
 end
 

benchmark/benchmarks.jl

@@ -3,8 +3,7 @@ using BenchmarkTools
 const SUITE = BenchmarkGroup()
 
 for file in readdir(@__DIR__)
-  if startswith(file, "bench_") && endswith(file, ".jl")
-    SUITE[file[length("bench_") + 1:end - length(".jl")]] =
-      include(file)
-  end
+    if startswith(file, "bench_") && endswith(file, ".jl")
+        SUITE[file[length("bench_")+1:end-length(".jl")]] = include(file)
+    end
 end

benchmark/run_henon_heiles.jl

@@ -1,15 +1,15 @@
 using HOODESolver
-using DifferentialEquations
+using OrdinaryDiffEq
 using LinearAlgebra
 using Plots
 using Random
 
 function getmindif(tab::Vector{Vector{BigFloat}})
     nmmin = Inf
     ret = (0, 0)
-    for i=1:size(tab,1)
-        for j=(i+1):size(tab,1)
-            nm = norm(tab[i]-tab[j], Inf)
+    for i = 1:size(tab, 1)
+        for j = (i+1):size(tab, 1)
+            nm = norm(tab[i] - tab[j], Inf)
             if nm < nmmin
                 nmmin = nm
                 ret = i, j
@@ -20,44 +20,51 @@ function getmindif(tab::Vector{Vector{BigFloat}})
 end
 # henon_heiles
 function fct(u, p, t)
-    return [0, u[4], -2u[1]*u[2], -u[2]-u[1]^2+u[2]^2]
+    return [0, u[4], -2u[1] * u[2], -u[2] - u[1]^2 + u[2]^2]
 end
 
 
 function fctmain(n_tau, prec)
     Random.seed!(19999909)
     setprecision(512)
-    u0 = 2rand(BigFloat,4)-ones(BigFloat,4)
+    u0 = 2rand(BigFloat, 4) - ones(BigFloat, 4)
     setprecision(prec)
     u0 = BigFloat.(u0)
     println("u0=$u0")
     t_max = big"1.0"
-    epsilon=big"1e-7"
+    epsilon = big"1e-7"
     println("epsilon=$epsilon")
-    nbmaxtest=15
-    ordmax=15
-    debord=3
-    pasord=1
-    y = ones(Float64, nbmaxtest, div(ordmax-debord,pasord)+1 )
-    res_gen = Array{ Array{BigFloat,1}, 2}(undef, nbmaxtest, div(ordmax-debord,pasord)+1 )
-    x=zeros(Float64,nbmaxtest)
-    A = [0 0 1 0; 0 0 0 0;-1 0 0 0; 0 0 0 0]
-    t_0=big"0.0"
-    t_max=big"1.0"
+    nbmaxtest = 15
+    ordmax = 15
+    debord = 3
+    pasord = 1
+    y = ones(Float64, nbmaxtest, div(ordmax - debord, pasord) + 1)
+    res_gen = Array{Array{BigFloat,1},2}(undef, nbmaxtest, div(ordmax - debord, pasord) + 1)
+    x = zeros(Float64, nbmaxtest)
+    A = [0 0 1 0; 0 0 0 0; -1 0 0 0; 0 0 0 0]
+    t_0 = big"0.0"
+    t_max = big"1.0"
     prob = HOODEProblem(fct, u0, (t_0, t_max), missing, A, epsilon)
     nm = NaN
     ordmax += 1
-    ordprep=undef
-    nb = 10*2^(nbmaxtest)
-    solref=undef
+    ordprep = undef
+    nb = 10 * 2^(nbmaxtest)
+    solref = undef
     while isnan(nm)
         ordmax -= 1
-        @time sol = solve(prob, nb_t=nb, order=ordmax, getprecision=false, nb_tau=n_tau, dense=false)
+        @time sol = solve(
+            prob,
+            nb_t = nb,
+            order = ordmax,
+            getprecision = false,
+            nb_tau = n_tau,
+            dense = false,
+        )
         solref = sol[end]
-	    nm = norm(solref, Inf)
+        nm = norm(solref, Inf)
     end
 
-    tabsol = Array{Array{BigFloat,1},1}(undef,1)
+    tabsol = Array{Array{BigFloat,1},1}(undef, 1)
 
     tabsol[1] = solref
 
@@ -67,62 +74,70 @@ function fctmain(n_tau, prec)
     solref_gen = solref
     indref = 1
 
-    ind=1
-    for order=debord:pasord:ordmax
-        ordprep=order+2
+    ind = 1
+    for order = debord:pasord:ordmax
+        ordprep = order + 2
         println("eps=$epsilon solRef=$solref order=$order")
         nb = 10
         indc = 1
-        labels=Array{String,2}(undef, 1, order-debord+1)
-        sol =undef
-        par_u0=missing
+        labels = Array{String,2}(undef, 1, order - debord + 1)
+        sol = undef
+        par_u0 = missing
         println("preparation ordre $order + 2")
         while indc <= nbmaxtest
-            @time solall = solve(prob, nb_t=nb, order=order, getprecision=false, nb_tau=n_tau, par_u0=par_u0, dense=false)
+            @time solall = solve(
+                prob,
+                nb_t = nb,
+                order = order,
+                getprecision = false,
+                nb_tau = n_tau,
+                par_u0 = par_u0,
+                dense = false,
+            )
             par_u0 = solall.par_u0
             sol = solall[end]
             push!(tabsol, sol)
-            res_gen[indc,ind] = sol
+            res_gen[indc, ind] = sol
             (a, b), nm = getmindif(tabsol)
             if a != indref
                 println("New solref !!!! a=$a, b=$b nm=$nm")
                 indref = a
                 solref = tabsol[a]
-                for i=1:ind
-                    borne = (i <ind) ? size(res_gen,1) : indc
+                for i = 1:ind
+                    borne = (i < ind) ? size(res_gen, 1) : indc
                     for j = 1:borne
-                        nm2 = min( norm(res_gen[j,i] - solref, Inf), 1.1)
-                        y[j,i] = nm2 == 0 ? nm : nm2
-                   end
+                        nm2 = min(norm(res_gen[j, i] - solref, Inf), 1.1)
+                        y[j, i] = nm2 == 0 ? nm : nm2
+                    end
                 end
             else
-                diff=solref-sol
-                y[indc,ind] = min(norm(diff,Inf), 1.1)
+                diff = solref - sol
+                y[indc, ind] = min(norm(diff, Inf), 1.1)
             end
-            x[indc] = 1.0/nb
+            x[indc] = 1.0 / nb
             println("result=$y")
             println("log2(y)=$(log2.(y))")
             nb *= 2
             indc += 1
         end
-        for i=debord:pasord:order
-            labels[1,(i-debord)÷pasord+1] = " eps,order=$(convert(Float32,epsilon)),$i "
+        for i = debord:pasord:order
+            labels[1, (i-debord)÷pasord+1] = " eps,order=$(convert(Float32,epsilon)),$i "
         end
-        p=Plots.plot(
-                        x,
-                        view(y,:,1:ind),
-                        xlabel="delta t",
-                        xaxis=:log,
-                        ylabel="error",
-                        yaxis=:log,
-                        legend=:bottomright,
-                        label=labels,
-                        marker=2
-                    )
+        p = Plots.plot(
+            x,
+            view(y, :, 1:ind),
+            xlabel = "delta t",
+            xaxis = :log,
+            ylabel = "error",
+            yaxis = :log,
+            legend = :bottomright,
+            label = labels,
+            marker = 2,
+        )
         prec_v = precision(BigFloat)
-        eps_v = convert(Float32,epsilon)
-        Plots.savefig(p,"out/r4_$(prec_v)_$(eps_v)_$(order)_$(n_tau)_henon_heiles.pdf")
-        ind+= 1
+        eps_v = convert(Float32, epsilon)
+        Plots.savefig(p, "out/r4_$(prec_v)_$(eps_v)_$(order)_$(n_tau)_henon_heiles.pdf")
+        ind += 1
     end
 end
-fctmain(32,512)
+fctmain(32, 512)