From 467717436b8c80087f12f5f8d6e10b3123ecf512 Mon Sep 17 00:00:00 2001 From: Yoonku Hwang Date: Tue, 17 May 2016 19:22:29 +0900 Subject: [PATCH] partially fixed (not enough order of convergence) --- src/Numerics/OdeSolvers/AdamsBashforth.cs | 9 ++++++--- src/UnitTests/OdeSolvers/OdeSolverTest.cs | 4 +--- 2 files changed, 7 insertions(+), 6 deletions(-) diff --git a/src/Numerics/OdeSolvers/AdamsBashforth.cs b/src/Numerics/OdeSolvers/AdamsBashforth.cs index 825a4b8b..aae22131 100644 --- a/src/Numerics/OdeSolvers/AdamsBashforth.cs +++ b/src/Numerics/OdeSolvers/AdamsBashforth.cs @@ -65,19 +65,22 @@ namespace MathNet.Numerics.OdeSolvers /// Number of subintervals /// ode model /// - public static double[] SecondOrder(double y0, double y1, double start, double end, int N, Func f) + public static double[] SecondOrder(double y0, double start, double end, int N, Func f) { double dt = (end - start) / (N - 1); double t = start; double[] y = new double[N]; + double k1 = f(t, y0); + double k2 = f(t + dt, y0 + dt * k1); + double y1 = y0 + 0.5 * dt * (k1 + k2); y[0] = y0; + t += dt; y[1] = y1; for (int i = 2; i < N; i++) { - y1 = f(t + dt, y1); y[i] = y1 + dt * (1.5 * f(t + dt, y1) - 0.5 * f(t, y0)); t += dt; - y0 = y1; + y0 = y[i-1]; y1 = y[i]; } return y; diff --git a/src/UnitTests/OdeSolvers/OdeSolverTest.cs b/src/UnitTests/OdeSolvers/OdeSolverTest.cs index bca1ef38..dfd956fa 100644 --- a/src/UnitTests/OdeSolvers/OdeSolverTest.cs +++ b/src/UnitTests/OdeSolvers/OdeSolverTest.cs @@ -72,9 +72,7 @@ namespace MathNet.Numerics.UnitTests.OdeSolvers for (int k = 0; k < 4; k++) { double y0 = 0; - int N = Convert.ToInt32(Math.Pow(2, k + 6)); - double dt = 2.0 / (N - 1); - double[] y_t = AdamsBashforth.SecondOrder(y0, sol(dt), 0, 2, N, ode); + double[] y_t = AdamsBashforth.SecondOrder(y0, 0, 2, Convert.ToInt32(Math.Pow(2, k + 6)), ode); error = Math.Abs(sol(2) - y_t.Last()); if (oldError != 0) ratio = Math.Log(oldError / error, 2);