Fix McAte5 time dependence

Fixes #2066
SciML · Dec 29, 2023 · f5e0454 · f5e0454
1 parent c9db428
commit f5e0454
Show file tree

Hide file tree

Showing 2 changed files with 28 additions and 2 deletions.
diff --git a/src/perform_step/symplectic_perform_step.jl b/src/perform_step/symplectic_perform_step.jl
@@ -568,7 +568,7 @@ end
     kdu = f.f1(du, u, p, tnew)
     du = du + dt * a5 * kdu
 
-    tnew = tnew + t + a5 * dt
+    tnew = tnew + a5 * dt
     ku = f.f2(du, u, p, tnew)
     u = u + dt * b6 * ku
 
@@ -625,7 +625,7 @@ end
     f.f1(kdu, du, u, p, tnew)
     @.. broadcast=false du=du + dt * a5 * kdu
 
-    tnew = tnew + t + a5 * dt
+    tnew = tnew + a5 * dt
     f.f2(ku, du, u, p, tnew)
     @.. broadcast=false u=u + dt * b6 * ku
 

diff --git a/test/algconvergence/symplectic_tests.jl b/test/algconvergence/symplectic_tests.jl
@@ -67,3 +67,29 @@ end
         @test_throws ArgumentError solve(prob, alg(); dt = dt)
     end
 end
+
+function motionfuncDirect1(dv,v,u,p,t)
+    # 1:Electron, 2: Be
+    ω_1,ω_2,γ,m_1,m_2,η,ω_d=p
+    dv[1]=-ω_1^2*u[1]*(1+η*cos(ω_d*t))-γ*u[2]/m_1
+    dv[2]=-ω_2^2*u[2]-γ*u[1]/m_2
+end
+
+function motionfuncDirect1(v,u,p,t)
+    # 1:Electron, 2: Be
+    ω_1,ω_2,γ,m_1,m_2,η,ω_d=p
+    [-ω_1^2*u[1]*(1+η*cos(ω_d*t))-γ*u[2]/m_1,-ω_2^2*u[2]-γ*u[1]/m_2]
+end
+
+param=[90386.15717208837, 3938.9288690708827, 8560.718748264337, 0.000544617021484666, 8.947079933513658, 0.7596480420227258, 78778.57738141765]
+u0_direct=zeros(2) # mm, mm
+v0_direct = [0.0, 135.83668926684385]
+tspan=(0.0, 1.321179076090661)
+prob_direct=SecondOrderODEProblem(motionfuncDirect1,v0_direct,u0_direct,tspan,param)
+dt=2e-8
+ref=solve(prob_direct,DPRKN12(),abstol=1e-12,reltol=1e-12,maxiters=1e7,saveat=0.01)
+
+@testset "symplectic time-dependent $alg-$iip-$pa" for (alg, x, d) in ALGOS
+    sol=solve(prob_direct,alg,dt=dt,saveat=0.01)
+    @test maximum(ref-sol) < 1e-3
+end