Mercurial > repos > public > sbplib

changeset 931:384ca2331a12 feature/timesteppers

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Make Rungekutta class allow for discrete data.
author Martin Almquist <malmquist@stanford.edu>
date Mon, 03 Dec 2018 16:26:44 -0800
parents 34d882bffae4
children 3860dad28239
files +time/Rungekutta.m
diffstat 1 files changed, 12 insertions(+), 7 deletions(-) [+]
line wrap: on
line diff
diff -r 34d882bffae4 -r 384ca2331a12 +time/Rungekutta.m
--- a/+time/Rungekutta.m	Mon Dec 03 16:26:18 2018 -0800
+++ b/+time/Rungekutta.m	Mon Dec 03 16:26:44 2018 -0800
@@ -15,8 +15,9 @@
     methods
         % Timesteps v_t = F(v,t), using the specified RK method from t = t0 with
         % timestep dt and initial conditions v = v0
-        function obj = Rungekutta(F, dt, t0, v0, method)
+        function obj = Rungekutta(F, dt, t0, v0, method, discreteData)
             default_arg('method',"rk4");
+            default_arg('discreteData', []);
             obj.F = F;
             obj.dt = dt;
             obj.t = t0;
@@ -28,12 +29,16 @@
             [s,a,b,c] = time.rk.butcherTableau(method);
             obj.coeffs = struct('s',s,'a',a,'b',b,'c',c);
 
-            % TODO: method "rk4" is also implemented in the butcher tableau, but the rungekutta_4.m implementation
-            % might be slightly more efficient. Need to do some profiling before deciding whether or not to keep it.
-            if (method == "rk4")
-                obj.scheme = @time.rk.rungekutta_4;
+            if isempty(discreteData)
+                % TODO: method "rk4" is also implemented in the butcher tableau, but the rungekutta_4.m implementation
+                % might be slightly more efficient. Need to do some profiling before deciding whether or not to keep it.
+                if (method == "rk4")
+                    obj.scheme = @(v,t,dt,F,n) time.rk.rungekutta_4(v ,t, dt, F);
+                else
+                    obj.scheme = @(v,t,dt,F,n) time.rk.rungekutta(v, t, dt, F, obj.coeffs);
+                end
             else
-                obj.scheme = @(v,t,dt,F) time.rk.rungekutta(v, t , dt, F, obj.coeffs);
+                obj.scheme = @(v,t,dt,F,n) time.rk.rungekuttaDiscreteData(v, t, dt, F, obj.coeffs, discreteData, n);
             end
         end
 
@@ -51,7 +56,7 @@
         end
 
         function obj = step(obj)
-            [obj.v, obj.V, obj.K] = obj.scheme(obj.v, obj.t, obj.dt, obj.F);
+            [obj.v, obj.V, obj.K] = obj.scheme(obj.v, obj.t, obj.dt, obj.F, obj.n);
             obj.t = obj.t + obj.dt;
             obj.n = obj.n + 1;
         end