Heun.m

Go to the documentation of this file.

namespace solvers{


/* (Autoinserted by mtoc++)
 * This source code has been filtered by the mtoc++ executable,
 * which generates code that can be processed by the doxygen documentation tool.
 *
 * On the other hand, it can neither be interpreted by MATLAB, nor can it be compiled with a C++ compiler.
 * Except for the comments, the function bodies of your M-file functions are untouched.
 * Consequently, the FILTER_SOURCE_FILES doxygen switch (default in our Doxyfile.template) will produce
 * attached source files that are highly readable by humans.
 *
 * Additionally, links in the doxygen generated documentation to the source code of functions and class members refer to
 * the correct locations in the source code browser.
 * However, the line numbers most likely do not correspond to the line numbers in the original MATLAB source files.
 */

class Heun
  :public solvers.BaseCustomSolver {
    public:


        Heun(double MaxStep) {

            this = this@solvers.BaseCustomSolver;

            this.Name= " Explicit Heun "" s method ";
            this.SolverType= solvers.SolverTypes.Explicit;
            if nargin == 1
                this.MaxStep= MaxStep;
            end
        }
    protected: /* ( Sealed ) */

        function matrixx = customSolve(function_handle odefun,rowvec t,colvec x0,rowvec<integer> outputtimes) {

            /*  Initialize result */
            steps = length(t);
            dt = t(2:end)-t(1:end-1);

            rtm = this.RealTimeMode;
/*              if rtm
 *                 ed = solvers.SolverEventData;
 *                 x = [];
 *             else
 *                 x = [x0 zeros(size(x0,1),steps-1)];
 *             end
             * Initialize output index counter */
            outidx = 2;
            if rtm
                ed = solvers.SolverEventData;
                x = [];
            else
                effsteps = length(outputtimes);
                /*  Create return matrix in size of effectively desired timesteps */
                x = [x0 zeros(size(x0,1),effsteps-1)];
            end
            /*  Solve for each time step */
            oldx = x0;
            for idx = 2:steps;
                f = odefun(t(idx-1),oldx);
                hlp = (dt(idx-1)/2)*(f + odefun(t(idx),oldx + dt(idx-1)*f));

                /*  Check if a mass matrix is present */
                if ~isempty(this.M)
/*                      newx = oldx + this.M.evaluate(t(idx-1))\hlp; */
                    [L,U] = this.M.getLU(t(idx-1)); 
                    newx = U\(L\(this.M.evaluate(t(idx-1))*oldx + hlp));
                else
                    newx = oldx + hlp;
                end

/*                  if rtm
 *                     ed.Times = t(idx);
 *                     ed.States = newx;
 *                     this.notify('StepPerformed',ed);
 *                 else
 *                     x(:,idx) = newx;%#ok
 *                 end
                 * Only produce output at wanted timesteps */
                if outputtimes(outidx) == idx
                    if rtm                        
                        /*  Real time mode: Fire StepPerformed event */
                        ed.Times= t(idx);
                        ed.States= newx;
                        this.notify(" StepPerformed ",ed);
                        /*  Normal mode: Collect solution in result vector */
                    else
                        x(:,outidx) = newx;
                    end
                    outidx = outidx+1;
                end
                oldx = newx;
            end
        }
};
}