Docs/ascher DAE

docs/references.bib

@@ -107,3 +107,28 @@ @inproceedings{Lor96
   year={1996},
   organization={Reading}
 }
+
+
+@article{Asc89,
+author = {Ascher, Uri},
+title = {On Symmetric Schemes and Differential-Algebraic Equations},
+journal = {SIAM Journal on Scientific and Statistical Computing},
+volume = {10},
+number = {5},
+pages = {937-949},
+year = {1989},
+doi = {10.1137/0910054}
+}
+
+@article{Pet86,
+ author = {L. R. Petzold},
+ journal = {SIAM Journal on Numerical Analysis},
+ number = {4},
+ pages = {837--852},
+ publisher = {Society for Industrial and Applied Mathematics},
+ title = {Order Results for Implicit Runge-Kutta Methods Applied to Differential/ Algebraic Systems},
+ volume = {23},
+ year = {1986},
+ url = {https://www.jstor.org/stable/2157625}
+}
+

toolbox/+otp/+ascherlineardae/+presets/Canonical.m

@@ -1,22 +1,28 @@
 classdef Canonical < otp.ascherlineardae.AscherLinearDAEProblem
-    %CANONICAL The problem formulation of the linear DAE from the literature
-    %
-    % See
-    %    Ascher, Uri. "On symmetric schemes and differential-algebraic equations."
-    %    SIAM journal on scientific and statistical computing 10.5 (1989): 937-949.
-
+    % The problem defined by Uri Ascher in :cite:p:`Asc89` (sec. 2) which uses time span $t \in [0, 1]$ and intial 
+    % condition $[y_0, z_0]^T = [1, β]^T$.
+    % 
     methods
-        function obj = Canonical(beta)
-            tspan = [0.0; 1];
+        function obj = Canonical(varargin)
+            % Create the Canonical CUSP problem object.
+            %
+            % Parameters
+            % ----------
+            % varargin
+            %    A variable number of name-value pairs. The accepted names are
+            %
+            %    - ``Beta`` – Value of $β$.
 
-            if nargin < 1
-                beta = 0.5;
-            end
+            p = inputParser;
+            p.addParameter('beta', 1);
+            p.parse(varargin{:});
+            opts = p.Results;
 
-            params = otp.ascherlineardae.AscherLinearDAEParameters;
-            params.Beta = beta;
+            params      = otp.ascherlineardae.AscherLinearDAEParameters;
+            params.Beta = opts.beta;
 
-            y0 = [1; beta];
+            y0    = [1; params.Beta];
+            tspan = [0.0; 1.0];
 
             obj = obj@otp.ascherlineardae.AscherLinearDAEProblem(tspan, y0, params);
         end

toolbox/+otp/+ascherlineardae/+presets/Petzold.m

@@ -0,0 +1,18 @@
+classdef Petzold < otp.ascherlineardae.AscherLinearDAEProblem
+    % The Petzold DAE example :cite:p:`Pet86` as a special case of the Ascher linear DAE problem. This preset uses time
+    % span $t \in [0, 1]$ and $β = 0 $ with the initial condition $[y_0, z_0]^T = [1, 0]^T $.
+    % 
+    methods
+        function obj = Petzold
+            % Create the Petzold example of the Ascher linear DAE problem object.
+
+            params      = otp.ascherlineardae.AscherLinearDAEParameters;
+            params.Beta = 0;
+
+            y0    = [1; params.Beta];
+            tspan = [0.0; 1.0];
+
+            obj = obj@otp.ascherlineardae.AscherLinearDAEProblem(tspan, y0, params);
+        end
+    end
+end

toolbox/+otp/+ascherlineardae/+presets/Stiff.m

@@ -0,0 +1,17 @@
+classdef Stiff < otp.ascherlineardae.AscherLinearDAEProblem
+    % The Stiff example from :cite:p:`Asc89`. A variant of the Ascher linear DAE problem which uses time span 
+    % $t \in [0, 1]$ and $β = 100$ with the initial condition $[y_0, z_0]^T = [1, 100]^T$.
+    % 
+    methods
+        function obj = Stiff
+            % Create the stiff example of the Ascher linear DAE problem object.
+            params      = otp.ascherlineardae.AscherLinearDAEParameters;
+            params.Beta = 100;
+
+            y0    = [1; params.Beta];
+            tspan = [0.0; 1.0];
+
+            obj = obj@otp.ascherlineardae.AscherLinearDAEProblem(tspan, y0, params);
+        end
+    end
+end

toolbox/+otp/+ascherlineardae/AscherLinearDAEParameters.m

@@ -1,7 +1,7 @@
 classdef AscherLinearDAEParameters
-    %ASCHERLINEARDAEPARAMETERS 
+    % Parameters for Ascher Linear DAE problem 
     properties
-        %Beta is an arbitrary scalar parameter
-        Beta %MATLAB ONLY: (1,1) {mustBeNumeric} = 0.5
+        % A scalar parameter $β$ in the linear model. It affects the stifness of the problem.
+        Beta %MATLAB ONLY: (1,1) {mustBeNumeric} = 1
     end
 end

toolbox/+otp/+ascherlineardae/AscherLinearDAEProblem.m

@@ -1,8 +1,61 @@
 classdef AscherLinearDAEProblem < otp.Problem
-    %ASCHERLINEARDAEPROBLEM This is an Index-1 DAE problem
+    % A linear differential-algebraic problem with time-dependant mass matrix.
     %
+    % The Ascher linear DAE Problem :cite:p:`Asc89` is an index-1 differential-agebraic equation given by
+    %
+    % $$
+    % \begin{bmatrix}
+    % 1 & -t \\
+    % 0 & 0
+    % \end{bmatrix} \begin{bmatrix} y'(t) \\ z'(t) \end{bmatrix} = \left[ \begin{array}{cc}
+    % -1 & 1+t \\
+    % β & -1-β t
+    % \end{array}\right] \begin{bmatrix} y(t) \\ z(t) \end{bmatrix} + \begin{bmatrix}
+    % 0 \\
+    % \sin(t)
+    % \end{bmatrix}.
+    % $$
+    %
+    % When the initial condition $y(0) = 1 , z(0) = β$ is used, the problem has the following closed-form solution:
+    %
+    % $$
+    % \begin{bmatrix} y(t)\\ z(t) \end{bmatrix} = \begin{bmatrix}
+    % t \sin(t) + (1 + β  t) e^{-t}\\
+    % β  e^{-t} + \sin(t)
+    % \end{bmatrix}.
+    % $$
+    % This DAE problem can be used to investigate the convergence of implcit time-stepping methods due to its stiffness
+    % and time-dependant mass matrix.
+    %
+    % Notes
+    % -----
+    % +---------------------+-----------------------------------------+
+    % | Type                | DAE                                     |
+    % +---------------------+-----------------------------------------+
+    % | Number of Variables | 2                                       |
+    % +---------------------+-----------------------------------------+
+    % | Stiff               | possibly, depending on $β$              |
+    % +---------------------+-----------------------------------------+
+    %
+    % Example
+    % -------
+    % >>> problem = otp.ascherlineardae.presets.Canonical('Beta', 50);
+    % >>> t = linspace(0, 1);
+    % >>> sol = problem.solveExactly(t);
+    % >>> problem.plot(t, sol);
+
     methods
         function obj = AscherLinearDAEProblem(timeSpan, y0, parameters)
+            % Create an Ascher linear DAE problem object.
+            %
+            % Parameters
+            % ----------
+            % timeSpan : numeric(1, 2)
+            %    The start and final time.
+            % y0 : numeric(2, 1)
+            %    The initial condition.
+            % parameters : AscherLinearDAEParameters
+            %    The parameters.
             obj@otp.Problem('Ascher Linear DAE', 2, timeSpan, y0, parameters);
         end
     end
@@ -17,7 +70,15 @@ function onSettingsChanged(obj)
                 'MStateDependence', 'none', ...
                 'MassSingular', 'yes');
         end
-
+
+        function label = internalIndex2label(~, index)
+            if index == 1
+                label = 'Differential';
+            else
+                label = 'Algebraic';
+            end
+        end
+
         function y = internalSolveExactly(obj, t)
             beta = obj.Parameters.Beta;
             if ~isequal(obj.Y0, [1; beta])
@@ -35,4 +96,3 @@ function onSettingsChanged(obj)
         end
     end
 end
-