BasePCDSystem.m

Go to the documentation of this file.

namespace models{
namespace pcd{


/* (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 BasePCDSystem
  :public models.BaseFirstOrderSystem {
    public: /* ( Dependent, setObservable ) */

        h;
        Omega;
    private:

        fh = "[]";

        fOmega = "[]";

        
    public:

        Diff;
        hs;
        Dims;
        ReacCoeff;
    public:

        BasePCDSystem(models.BaseFullModel model) {
            this = this@models.BaseFirstOrderSystem(model);

            /*  Scale diffusion coefficients */
            m = this.Model;
            this.Diff= [m.d2/m.d1 m.d3/m.d1 m.d4/m.d1];

            this.registerProps(" h "," Omega ");

            this.ReacCoeff= [m.Kc1_real m.Kc2_real m.Kd1_real ...
                              m.Kd2_real m.Kd3_real m.Kd4_real ...
                              m.Kp1_real m.Kp2_real]^t * m.tau;

            /*  Add parameters
             * Activation area */
            p = this.addParam(" area ", .1, " Range ", [.01, 1], " Desired ", 10);
            p.Spacing= " lin ";

            /*  Activation rate */
            p = this.addParam(" rate ", .001, " Range ", [0.0005, 1], " Desired ", 15);
            p.Spacing= " log ";

            /*  Activation time */
            p = this.addParam(" atime ", .5, " Range ", [0, 1], " Desired ", 5);
            p.Spacing= " lin ";

            /*  Exponent */
            p = this.addParam(" exponent ", 2, " Range ", [1, 4], " Desired ", 5);
            p.Spacing= " lin ";
        }


        
#if 0 //mtoc++: 'get.h'
function h = h() {
            h = this.fh;
        }

#endif


        
#if 0 //mtoc++: 'set.h'
function  h(value) {
            if any(value >= this.fOmega(:,2))
                error(" Cannot choose a step size h=%e value larger or equal to the geometry [%s]. ",...
                    value,Utils.implode(this.fOmega(:)," ,  "," %g "));
            end
            this.fh= value;
            this.hs= value / this.Model.L;
            this.updateDimensions;

            m = this.Model;
            this.MaxTimestep= [];
            if ~isa(m.ODESolver," solvers.IImplSolver ") && ~isa(m.ODESolver," solvers.SemiImplicitEuler ")
                maxdt = .95*(value^2/max([m.d1 m.d2 m.d3 m.d4]));
                maxdtsc = .95*(this.hs^2/max([1 this.Diff]));

                if (maxdtsc - maxdt/this.Model.tau) /  maxdtsc > 1e-15
                    error(" Inconsistent scaling. Please check. ");
                end
                /*  Set max timestep value (scaled!) */
                this.MaxTimestep= maxdtsc;

                if maxdt < m.dt
                    warning(" models:pcd:CFL "," CFL condition violated with h=%e and current dt=%e.\nSetting System.MaxTimestep=%e (scaled, effective value %e)\n ", value, m.dt, maxdtsc, maxdt);
                end
            end
        }

#endif


        
#if 0 //mtoc++: 'get.Omega'
function v = Omega() {
            v = this.fOmega;
        }

#endif


        
#if 0 //mtoc++: 'set.Omega'
function  Omega(value) {
            if size(value,2) ~= 2
                error(" Omega needs to be a dim x 2 matrix with the spatial extend for each dimension in a row. ");
            end
            this.fOmega= value;

            /*  Update the dimensions */
            this.updateDimensions;
        }

#endif


        function  setConfig(colvec<double> mu,integer inputidx) {
            setConfig@models.BaseFirstOrderSystem(this, mu, inputidx);
            m = this.Model;
            this.ReacCoeff= [m.Kc1_real m.Kc2_real m.Kd1_real ...
                              m.Kd2_real m.Kd3_real m.Kd4_real ...
                              m.Kp1_real m.Kp2_real]^t * m.tau;
        }

    
    protected:


        function  updateDimensions() {
            if ~isempty(this.fh) && ~isempty(this.fOmega)
                nd = size(this.fOmega,1);
                this.Dims= zeros(1,nd);
                for d=1:nd
                    this.Dims(d) = length(this.fOmega(d,1):this.h:this.fOmega(d,2));
                end
                m = prod(this.Dims);

                this.NumStateDofs= 4*m;
                updateDimensions@models.BaseFirstOrderSystem(this);

                /*  Set state scaling */
                ss = zeros(4*m,1);
                ss(1:m) = this.Model.xa0;
                ss(m+1:2*m) = this.Model.ya0;
                ss(2*m+1:3*m) = this.Model.xi0;
                ss(3*m+1:end) = this.Model.yi0;
                this.StateScaling= ss;

                /*  Set new initial values and output in 1D-3D systems */
                this.newSysDimension;

                /*  Call template method for custom actions in the core
                 * functions (diffusion matrix comp) */
                this.f.newSysDimension;
            end
        }

    
    protected: /* ( Abstract ) */

        virtual function  newSysDimension() = 0;
};
}
}