Mercurial > repos > public > sbplib
changeset 847:1c6f1595bb94 feature/burgers1d
Clean up in RK time stepper schemes
author | Vidar Stiernström <vidar.stiernstrom@it.uu.se> |
---|---|
date | Thu, 20 Sep 2018 18:36:45 +0200 |
parents | c6fcee3fcf1b |
children | c8ea9bbdc62c |
files | +time/+rk/rungekuttaRV.m +time/+rk/rungekutta_4RV.m +time/+rk/rungekutta_6.m +time/+rk/rungekutta_6RV.m +time/Rungekutta.m +time/Rungekutta4proper.m |
diffstat | 6 files changed, 1 insertions(+), 130 deletions(-) [+] |
line wrap: on
line diff
diff -r c6fcee3fcf1b -r 1c6f1595bb94 +time/+rk/rungekuttaRV.m --- a/+time/+rk/rungekuttaRV.m Thu Sep 20 17:51:19 2018 +0200 +++ b/+time/+rk/rungekuttaRV.m Thu Sep 20 18:36:45 2018 +0200 @@ -4,7 +4,6 @@ % for the specific method. RV is the residual viscosity which is updated % in between the stages and after the updated solution is computed. function v = rungekuttaRV(v, t , dt, F, RV, coeffs) - % Move one stage outside to avoid branching for updating the % residual inside the loop. k = zeros(length(v), coeffs.s);
diff -r c6fcee3fcf1b -r 1c6f1595bb94 +time/+rk/rungekutta_4RV.m --- a/+time/+rk/rungekutta_4RV.m Thu Sep 20 17:51:19 2018 +0200 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,18 +0,0 @@ -% Takes one time step of size k using the rungekutta method -% starting from v_0 and where the function F(v,t) gives the -% time derivatives. -function v = rungekutta_4RV(v, t , k, F, RV) - k1 = F(v, t, RV.getViscosity()); - - RV.update(v+0.5*k*k1, v, 0.5*k); - k2 = F(v+0.5*k*k1, t+0.5*k, RV.getViscosity()); - - RV.update(v+0.5*k*k2, v, 0.5*k); - k3 = F(v+0.5*k*k2, t+0.5*k, RV.getViscosity()); - - RV.update(v+k*k3, v, k); - k4 = F(v+k*k3,t+k, RV.getViscosity()); - - RV.update(v + (1/6)*(k1+2*(k2+k3)+k4)*k, v, k); - v = v + (1/6)*(k1+2*(k2+k3)+k4)*k; -end \ No newline at end of file
diff -r c6fcee3fcf1b -r 1c6f1595bb94 +time/+rk/rungekutta_6.m --- a/+time/+rk/rungekutta_6.m Thu Sep 20 17:51:19 2018 +0200 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,31 +0,0 @@ -% Takes one time step of size dt using the rungekutta method -% starting from v_0 and where the function F(v,t) gives the -% time derivatives. -function v = rungekutta_6(v, t , dt, F) - s = 7 - k = zeros(length(v),s) - a = zeros(7,6); - c = [0, 4/7, 5/7, 6/7, (5-sqrt(5))/10, (5+sqrt(5))/10, 1]; - b = [1/12, 0, 0, 0, 5/12, 5/12, 1/12]; - a = [ - 0, 0, 0, 0, 0, 0; - 4/7, 0, 0, 0, 0, 0; - 115/112, -5/16, 0, 0, 0, 0; - 589/630, 5/18, -16/45, 0, 0, 0; - 229/1200 - 29/6000*sqrt(5), 119/240 - 187/1200*sqrt(5), -14/75 + 34/375*sqrt(5), -3/100*sqrt(5), 0, 0; - 71/2400 - 587/12000*sqrt(5), 187/480 - 391/2400*sqrt(5), -38/75 + 26/375*sqrt(5), 27/80 - 3/400*sqrt(5), (1+sqrt(5))/4, 0; - -49/480 + 43/160*sqrt(5), -425/96 + 51/32*sqrt(5), 52/15 - 4/5*sqrt(5), -27/16 + 3/16*sqrt(5), 5/4 - 3/4*sqrt(5), 5/2 - 1/2*sqrt(5); - ]; - - for i = 1:s - u = v; - for j = 1:i-1 - u = u + dt*a(i,j)*k(:,j); - end - k(:,i) = F(t+c(i)*dt,u); - end - - for i = 1:s - v = v + dt*b(i)*k(:,i); - end -end
diff -r c6fcee3fcf1b -r 1c6f1595bb94 +time/+rk/rungekutta_6RV.m --- a/+time/+rk/rungekutta_6RV.m Thu Sep 20 17:51:19 2018 +0200 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,37 +0,0 @@ -% Takes one time step of size dt using the rungekutta method -% starting from v_0 and where the function F(v,t) gives the -% time derivatives. -function v = rungekutta_6RV(v, t , dt, F, RV) - s = 7; - k = zeros(length(v),s); - a = zeros(7,6); - c = [0, 4/7, 5/7, 6/7, (5-sqrt(5))/10, (5+sqrt(5))/10, 1]; - b = [1/12, 0, 0, 0, 5/12, 5/12, 1/12]; - a = [ - 0, 0, 0, 0, 0, 0; - 4/7, 0, 0, 0, 0, 0; - 115/112, -5/16, 0, 0, 0, 0; - 589/630, 5/18, -16/45, 0, 0, 0; - 229/1200 - 29/6000*sqrt(5), 119/240 - 187/1200*sqrt(5), -14/75 + 34/375*sqrt(5), -3/100*sqrt(5), 0, 0; - 71/2400 - 587/12000*sqrt(5), 187/480 - 391/2400*sqrt(5), -38/75 + 26/375*sqrt(5), 27/80 - 3/400*sqrt(5), (1+sqrt(5))/4, 0; - -49/480 + 43/160*sqrt(5), -425/96 + 51/32*sqrt(5), 52/15 - 4/5*sqrt(5), -27/16 + 3/16*sqrt(5), 5/4 - 3/4*sqrt(5), 5/2 - 1/2*sqrt(5); - ]; - - k(:,1) = F(v,t,RV.getViscosity()); - - for i = 2:s - u = v; - for j = 1:i-1 - u = u + dt*a(i,j)*k(:,j); - end - RV.update(u,v,c(i)*dt); - k(:,i) = F(u,t+c(i)*dt,RV.getViscosity()); - end - - u = v; - for i = 1:s - u = u + dt*b(i)*k(:,i); - end - RV.update(u,v,dt); - v = u; -end
diff -r c6fcee3fcf1b -r 1c6f1595bb94 +time/Rungekutta.m --- a/+time/Rungekutta.m Thu Sep 20 17:51:19 2018 +0200 +++ b/+time/Rungekutta.m Thu Sep 20 18:36:45 2018 +0200 @@ -13,7 +13,7 @@ % Timesteps v_t = F(v,t), using RK with specfied order from t = t0 with % timestep k and initial conditions v = v0 function obj = Rungekutta(F, k, t0, v0, order) - default_arg('order','4'); + default_arg('order',4); obj.F = F; obj.k = k; obj.t = t0;
diff -r c6fcee3fcf1b -r 1c6f1595bb94 +time/Rungekutta4proper.m --- a/+time/Rungekutta4proper.m Thu Sep 20 17:51:19 2018 +0200 +++ /dev/null Thu Jan 01 00:00:00 1970 +0000 @@ -1,42 +0,0 @@ -classdef Rungekutta4proper < time.Timestepper - properties - F - k - t - v - m - n - end - - - methods - % Timesteps v_t = F(v,t), using RK4 fromt t = t0 with timestep k and initial conditions v = v0 - function obj = Rungekutta4proper(F, k, t0, v0) - obj.F = F; - obj.k = k; - obj.t = t0; - obj.v = v0; - obj.m = length(v0); - obj.n = 0; - end - - function [v,t] = getV(obj) - v = obj.v; - t = obj.t; - end - - function obj = step(obj) - obj.v = time.rk.rungekutta_4(obj.v, obj.t, obj.k, obj.F); - obj.t = obj.t + obj.k; - obj.n = obj.n + 1; - end - end - - - methods (Static) - function k = getTimeStep(lambda) - k = rk4.get_rk4_time_step(lambda); - end - end - -end \ No newline at end of file