Kuramoto Sivashinky 1D PDE

src/+otp/+kuramotosivashinsky/+presets/Canonical.m

@@ -0,0 +1,29 @@
+classdef Canonical < otp.kuramotosivashinsky.KuramotoSivashinskyProblem
+
+    methods
+        function obj = Canonical(varargin)
+
+            p = inputParser;
+            addParameter(p, 'Size', 64, @isscalar);
+            addParameter(p, 'L', 25, @isscalar);
+
+            parse(p, varargin{:});
+
+            s = p.Results;
+
+            n = s.Size;
+
+
+            params.n = n;
+            params.l = s.L;
+
+            x = linspace(0, 10*pi, n + 1);
+
+            u0 = 4*cos(x(1:end-1)).';
+
+            tspan = [0, 100];
+
+            obj = obj@otp.kuramotosivashinsky.KuramotoSivashinskyProblem(tspan, u0, params);
+        end
+    end
+end

src/+otp/+kuramotosivashinsky/KuramotoSivashinskyProblem.m

@@ -0,0 +1,32 @@
+classdef KuramotoSivashinskyProblem < otp.Problem
+
+    methods
+        function obj = KuramotoSivashinskyProblem(timeSpan, y0, parameters)
+            obj@otp.Problem('Kuramoto Sivashinsky Problem', [], timeSpan, y0, parameters);
+        end
+    end
+
+    methods (Access = protected)
+        function onSettingsChanged(obj)
+            n = obj.Parameters.n;
+            l = obj.Parameters.l;
+
+            domain = [-l, l];
+            D = otp.utils.pde.D(n, domain, 'C');
+            L = otp.utils.pde.laplacian(n, domain, 1, 'C');
+
+            obj.Rhs = otp.Rhs(@(t, u) otp.kuramotosivashinsky.f(t, u, D, L), ...
+                otp.Rhs.FieldNames.Jacobian, @(t, u) otp.kuramotosivashinsky.jac(t, u, D, L));
+
+        end
+
+        function validateNewState(obj, newTimeSpan, newY0, newParameters)
+            validateNewState@otp.Problem(obj, ...
+                newTimeSpan, newY0, newParameters)
+
+            otp.utils.StructParser(newParameters) ...
+                .checkField('n', 'scalar', 'integer', 'finite', 'positive') ...
+                .checkField('l', 'scalar', 'integer', 'finite', 'positive');
+        end
+    end
+end

src/+otp/+kuramotosivashinsky/f.m

@@ -0,0 +1,5 @@
+function du = f(~, u, D, L)
+
+du = - L*(L*u + u) - u.*(D*u);
+
+end

src/+otp/+kuramotosivashinsky/jac.m

@@ -0,0 +1,5 @@
+function j = jac(~, u, D, L)
+
+j = -L*L - L - spdiags(u, 0, numel(u), numel(u))*D - spdiags(D*u, 0, numel(u), numel(u));
+
+end