# HG changeset patch # User Rik # Date 1444022334 25200 # Node ID eb9e2d187ed27cb02e990cdd2776ee6dac52c8d3 # Parent d746695bf49451be6066b698c80692f68f40891c maint: Use Octave coding conventions in scripts/ode/private dir. * AbsRel_Norm.m, fuzzy_compare.m, hermite_quartic_interpolation.m, integrate_adaptive.m, integrate_const.m, integrate_n_steps.m, kahan.m, ode_struct_value_check.m, odepkg_event_handle.m, odepkg_structure_check.m, runge_kutta_45_dorpri.m, starting_stepsize.m: Wrap long lines to < 80 chars. Use double quotes rather than single quotes where possible. Use ';' at end of keywords "return;" and "break;" Use '##" for stand-alone comments and '#' for end-of-line comments. Use two spaces after period before starting new sentence. Use '!' instead of '~' for logical negation. Use specific form of end (endif, endfor, etc.). Don't use line continuation marker '...' unless necessary. diff -r d746695bf494 -r eb9e2d187ed2 scripts/ode/private/AbsRel_Norm.m --- a/scripts/ode/private/AbsRel_Norm.m Sun Oct 04 16:24:32 2015 +0100 +++ b/scripts/ode/private/AbsRel_Norm.m Sun Oct 04 22:18:54 2015 -0700 @@ -48,3 +48,4 @@ endif endfunction + diff -r d746695bf494 -r eb9e2d187ed2 scripts/ode/private/fuzzy_compare.m --- a/scripts/ode/private/fuzzy_compare.m Sun Oct 04 16:24:32 2015 +0100 +++ b/scripts/ode/private/fuzzy_compare.m Sun Oct 04 22:18:54 2015 -0700 @@ -17,14 +17,15 @@ ## . ## -*- texinfo -*- -## @deftypefn {Function File} {[@var{res}] =} fuzzy_compare (@var{"string1"}, @var{string_set}, [@var{correctness}]) +## @deftypefn {Function File} {@var{res} =} fuzzy_compare (@var{"string1"}, @var{string_set}) +## @deftypefnx {Function File} {@var{res} =} fuzzy_compare (@var{"string1"}, @var{string_set}, @var{correctness}) ## ## Compare a string with a set of strings and returns the positions in the ## set of strings at which there are the fields that best fit the one we are ## comparing. ## -## The distance used to compare the words is the Levenshtein distance -## and for more details see +## The distance used to compare the words is the Levenshtein distance. +## For more details see ## @url{http://en.wikipedia.org/wiki/Levenshtein_distance}. ## ## This function must be called with one output argument @var{res} which @@ -95,9 +96,9 @@ res = []; m = length (string1); - fields_nb = size (string_set, 1); + fields_nb = rows (string_set); - values = inf .* ones (fields_nb, 1); + values = Inf (fields_nb, 1); string1 = deblank (string1); string2 = []; @@ -121,13 +122,13 @@ positions = find (values == minimus); if (minimus == 0) # exact match - if (size (positions, 1) != 1) + if (rows (positions) != 1) error ("OdePkg:InvalidArgument", - "there are %d strings perfectly matching ''%s''", - size (positions, 1), string1); + "there are %d strings perfectly matching '%s'", + rows (positions), string1); endif res = positions; - return + return; endif ## determine the tolerance with the formula described in the @@ -142,19 +143,18 @@ && isscalar (correctness) && correctness == 0) || (ischar (correctness) - && strcmp (lower (deblank (correctness)), 'exact'))) + && strcmp (lower (deblank (correctness)), "exact"))) error ("OdePkg:InvalidArgument", - "no exact matching for string ''%s''", string1); + "no exact matching for string '%s'", string1); endif - if (isnumeric (correctness) - && isscalar (correctness)) + if (isnumeric (correctness) && isscalar (correctness)) tolerance = correctness; endif endif ## returning the positions of the fields whose distance is lower ## than the tolerance - for i = 1:1:fields_nb + for i = 1:fields_nb if (values(i) <= tolerance) res = [res; i]; endif @@ -162,6 +162,3 @@ endfunction -## Local Variables: *** -## mode: octave *** -## End: *** diff -r d746695bf494 -r eb9e2d187ed2 scripts/ode/private/hermite_quartic_interpolation.m --- a/scripts/ode/private/hermite_quartic_interpolation.m Sun Oct 04 16:24:32 2015 +0100 +++ b/scripts/ode/private/hermite_quartic_interpolation.m Sun Oct 04 22:18:54 2015 -0700 @@ -17,16 +17,13 @@ ## . ## -*- texinfo -*- -## @deftypefn {Function File} {[@var{x_out}] =} hermite_quartic_interpolation (@var{t}, @var{x}, @var{der}, @var{t_out}) +## @deftypefn {Function File} {@var{x_out} =} hermite_quartic_interpolation (@var{t}, @var{x}, @var{der}, @var{t_out}) ## ## This function file can be called by an ODE solver function in order to -## interpolate the solution at the time @var{t_out} using 4th order +## interpolate the solution at the time @var{t_out} using a 4th order ## Hermite interpolation. ## -## This function must be called with one output arguments: @var{x_out} -## which contains the evaluation at @var{t_out} of the Hermite interpolant. -## -## The first input argument is the vector with two given times. +## The first input @var{t} is a vector with two given times. ## ## The second input argument is the vector with the values of the function ## to interpolate at the times specified in @var{t} and at the middle point. @@ -34,6 +31,9 @@ ## The third input argument is the value of the derivatives of the function ## evaluated at the two extreme points. ## +## The output @var{x_out} is the evaluation of the Hermite interpolant at +## @var{t_out}. +## ## @end deftypefn ## ## @seealso{linear_interpolation, quadratic_interpolation, @@ -42,18 +42,18 @@ function x_out = hermite_quartic_interpolation (t, x, der, t_out) - # Rescale time on [0,1] + ## Rescale time on [0,1] s = (t_out - t(1)) / (t(2) - t(1)); - # Hermite basis functions - # H0 = 1 - 11*s.^2 + 18*s.^3 - 8*s.^4; - # H1 = s - 4*s.^2 + 5*s.^3 - 2*s.^4; - # H2 = 16*s.^2 - 32*s.^3 + 16*s.^4; - # H3 = - 5*s.^2 + 14*s.^3 - 8*s.^4; - # H4 = s.^2 - 3*s.^3 + 2*s.^4; + ## Hermite basis functions + ## H0 = 1 - 11*s.^2 + 18*s.^3 - 8*s.^4; + ## H1 = s - 4*s.^2 + 5*s.^3 - 2*s.^4; + ## H2 = 16*s.^2 - 32*s.^3 + 16*s.^4; + ## H3 = - 5*s.^2 + 14*s.^3 - 8*s.^4; + ## H4 = s.^2 - 3*s.^3 + 2*s.^4; - x_out = zeros (size (x, 1), length (t_out)); - for ii = 1:size (x, 1) + x_out = zeros (rows (x), length (t_out)); + for ii = 1:rows (x) x_out(ii,:) = (1 - 11*s.^2 + 18*s.^3 - 8*s.^4)*x(ii,1) ... + ( s - 4*s.^2 + 5*s.^3 - 2*s.^4)*(t(2)-t(1))*der(ii,1) ... + ( 16*s.^2 - 32*s.^3 + 16*s.^4)*x(ii,2) ... @@ -63,6 +63,3 @@ endfunction -## Local Variables: *** -## mode: octave *** -## End: *** diff -r d746695bf494 -r eb9e2d187ed2 scripts/ode/private/integrate_adaptive.m --- a/scripts/ode/private/integrate_adaptive.m Sun Oct 04 16:24:32 2015 +0100 +++ b/scripts/ode/private/integrate_adaptive.m Sun Oct 04 22:18:54 2015 -0700 @@ -61,7 +61,7 @@ function solution = integrate_adaptive (stepper, order, func, tspan, x0, options) - solution = struct; + solution = struct (); ## first values for time and solution t = tspan(1); @@ -78,16 +78,16 @@ endif dt = vdirection * min (abs (dt), options.MaxStep); - ## set parameters + ## Set parameters k = length (tspan); counter = 2; comp = 0.0; tk = tspan(1); options.comp = comp; - ## factor multiplying the stepsize guess + ## Factor multiplying the stepsize guess facmin = 0.8; - fac = 0.38^(1/(order+1)); ## formula taken from Hairer + fac = 0.38^(1/(order+1)); # formula taken from Hairer t_caught = false; @@ -122,7 +122,7 @@ while (counter <= k) facmax = 1.5; - ## compute integration step from t to t+dt + ## Compute integration step from t to t+dt if (isempty (k_vals)) [s, y, y_est, new_k_vals] = stepper (func, z(end), u(:,end), dt, options); @@ -144,7 +144,7 @@ err = AbsRel_Norm (y(:,end), u(:,end), options.AbsTol, options.RelTol, options.vnormcontrol, y_est(:,end)); - ## solution accepted only if the error is less or equal to 1.0 + ## Solution accepted only if the error is less or equal to 1.0 if (err <= 1) [tk, comp] = kahan (tk, comp, dt); @@ -162,8 +162,9 @@ (max (abs (z(end)), abs (tspan(counter)))) < 8*eps) ) counter++; - ## if there is an element in time vector at which the solution is required - ## the program must compute this solution before going on with next steps + ## if there is an element in time vector at which the solution is + ## required the program must compute this solution before going on with + ## next steps elseif (vdirection * z(end) > vdirection * tspan(counter)) ## initialize counter for the following cycle @@ -233,11 +234,12 @@ ## u_interp = ## hermite_quintic_interpolation ([z(i-1) z(i)], ## [u(:,i-1) u_half u(:,i)], - ## [k_vals(:,1) f_half k_vals(:,end)], + ## [k_vals(:,1) f_half ... + ## k_vals(:,end)], ## tspan(counter)); otherwise warning ("High order interpolation not yet implemented: ", - "using cubic iterpolation instead"); + "using cubic interpolation instead"); der(:,1) = feval (func, z(i-1) , u(:,i-1), options.vfunarguments{:}); der(:,2) = feval (func, z(i) , u(:,i), @@ -281,7 +283,7 @@ vcntiter = 0; ## Call plot only if a valid result has been found, therefore this - ## code fragment has moved here. Stop integration if plot function + ## code fragment has moved here. Stop integration if plot function ## returns false if (options.vhaveoutputfunction) for vcnt = 0:options.Refine # Approximation between told and t @@ -300,27 +302,28 @@ vpltret = feval (options.OutputFcn, vapproxtime, vapproxvals, [], options.vfunarguments{:}); if (vpltret) # Leave refinement loop - break + break; endif endfor if (vpltret) # Leave main loop solution.vunhandledtermination = false; - break + break; endif endif ## Call event only if a valid result has been found, therefore this - ## code fragment has moved here. Stop integration if veventbreak is + ## code fragment has moved here. Stop integration if veventbreak is ## true if (options.vhaveeventfunction) solution.vevent = odepkg_event_handle (options.Events, t(end), - x(:,end), [], options.vfunarguments{:}); + x(:,end), [], + options.vfunarguments{:}); if (! isempty (solution.vevent{1}) && solution.vevent{1} == 1) t(solution.vcntloop-1,:) = solution.vevent{3}(end,:); x(:,solution.vcntloop-1) = solution.vevent{4}(end,:)'; solution.vunhandledtermination = false; - break + break; endif endif @@ -337,18 +340,18 @@ dt = vdirection * min (abs (dt), options.MaxStep); ## Update counters that count the number of iteration cycles - solution.vcntcycles = solution.vcntcycles + 1; # Needed for cost statistics - vcntiter = vcntiter + 1; # Needed to find iteration problems + solution.vcntcycles += 1; # Needed for cost statistics + vcntiter += 1; # Needed to find iteration problems ## Stop solving because in the last 1000 steps no successful valid ## value has been found if (vcntiter >= 5000) - error (["Solving has not been successful. The iterative", + error (["Solving has not been successful. The iterative", " integration loop exited at time t = %f before endpoint at", - " tend = %f was reached. This happened because the iterative", + " tend = %f was reached. This happened because the iterative", " integration loop does not find a valid solution at this time", - " stamp. Try to reduce the value of ''InitialStep'' and/or", - " ''MaxStep'' with the command ''odeset''.\n"], + " stamp. Try to reduce the value of 'InitialStep' and/or", + " 'MaxStep' with the command 'odeset'.\n"], s(end), tspan(end)); endif @@ -362,20 +365,20 @@ if (vdirection * z(end) < vdirection * tspan(end)) if (solution.vunhandledtermination == true) error ("OdePkg:InvalidArgument", - ["Solving has not been successful. The iterative", + ["Solving has not been successful. The iterative", " integration loop exited at time t = %f", - " before endpoint at tend = %f was reached. This may", + " before endpoint at tend = %f was reached. This may", " happen if the stepsize grows smaller than defined in", - " vminstepsize. Try to reduce the value of ''InitialStep''", - " and/or ''MaxStep'' with the command ''odeset''.\n"], + " vminstepsize. Try to reduce the value of 'InitialStep'", + " and/or 'MaxStep' with the command 'odeset'.\n"], z(end), tspan(end)); else warning ("OdePkg:InvalidArgument", - ["Solver has been stopped by a call of ''break'' in the main", + ["Solver has been stopped by a call of 'break' in the main", " iteration loop at time t = %f before endpoint at tend = %f ", - " was reached. This may happen because the @odeplot function", - " returned ''true'' or the @event function returned", - " ''true''.\n"], + " was reached. This may happen because the @odeplot function", + " returned 'true' or the @event function returned", + " 'true'.\n"], z(end), tspan(end)); endif endif @@ -389,3 +392,4 @@ solution.x = x(:,1:end-f)'; endfunction + diff -r d746695bf494 -r eb9e2d187ed2 scripts/ode/private/integrate_const.m --- a/scripts/ode/private/integrate_const.m Sun Oct 04 16:24:32 2015 +0100 +++ b/scripts/ode/private/integrate_const.m Sun Oct 04 22:18:54 2015 -0700 @@ -63,7 +63,7 @@ function solution = integrate_const (stepper, func, tspan, x0, dt, options) - solution = struct; + solution = struct (); ## first values for time and solution t = tspan(1); @@ -72,7 +72,7 @@ vdirection = odeget (options, "vdirection", [], "fast"); if (sign (dt) != vdirection) error ("OdePkg:InvalidArgument", - "option ''InitialStep'' has a wrong sign"); + "option 'InitialStep' has a wrong sign"); endif ## setting parameters @@ -189,11 +189,11 @@ vcntiter = 0; ## Call plot only if a valid result has been found, therefore this - ## code fragment has moved here. Stop integration if plot function + ## code fragment has moved here. Stop integration if plot function ## returns false if (options.vhaveoutputfunction) - for vcnt = 0:options.Refine # Approximation between told and t - if (options.vhaverefine) # Do interpolation + for vcnt = 0:options.Refine # Approximation between told and t + if (options.vhaverefine) # Do interpolation vapproxtime = (vcnt + 1) / (options.Refine + 2); vapproxvals = (1 - vapproxtime) * vSaveVUForRefine ... + (vapproxtime) * y(:,end); @@ -207,18 +207,18 @@ endif vpltret = feval (options.OutputFcn, vapproxtime, vapproxvals, [], options.vfunarguments{:}); - if (vpltret) # Leave refinement loop - break + if (vpltret) # Leave refinement loop + break; endif endfor - if (vpltret) # Leave main loop + if (vpltret) # Leave main loop solution.vunhandledtermination = false; - break + break; endif endif ## Call event only if a valid result has been found, therefore this - ## code fragment has moved here. Stop integration if veventbreak is true + ## code fragment has moved here. Stop integration if veventbreak is true if (options.vhaveeventfunction) solution.vevent = odepkg_event_handle (options.Events, t(end), x(:,end), [], options.vfunarguments{:}); @@ -227,7 +227,7 @@ t(solution.vcntloop-1,:) = solution.vevent{3}(end,:); x(:,solution.vcntloop-1) = solution.vevent{4}(end,:)'; solution.vunhandledtermination = false; - break + break; endif endif @@ -238,12 +238,12 @@ ## Stop solving because the last 1000 steps no successful valid ## value has been found if (vcntiter >= 5000) - error (["Solving has not been successful. The iterative", + error (["Solving has not been successful. The iterative", " integration loop exited at time t = %f before endpoint at", - " tend = %f was reached. This happened because the iterative", + " tend = %f was reached. This happened because the iterative", " integration loop does not find a valid solution at this time", - " stamp. Try to reduce the value of ''InitialStep'' and/or", - " ''MaxStep'' with the command ''odeset''.\n"], + " stamp. Try to reduce the value of 'InitialStep' and/or", + " 'MaxStep' with the command 'odeset'.\n"], s(end), tspan(end)); endif @@ -257,19 +257,18 @@ if (vdirection * z(end) < vdirection * tspan(end)) if (solution.vunhandledtermination == true) error ("OdePkg:InvalidArgument", - ["Solving has not been successful. The iterative integration loop", - " exited at time t = %f before endpoint at tend = %f was", - " reached. This may happen if the stepsize grows smaller than", - " defined in vminstepsize. Try to reduce the value of", - " ''InitialStep'' and/or ''MaxStep'' with the command", - " ''odeset''.\n"], z(end), tspan(end)); + ["Solving has not been successful. The iterative integration" + " loop exited at time t = %f before endpoint at tend = %f was", + " reached. This may happen if the stepsize grows smaller than", + " defined in vminstepsize. Try to reduce the value of", + " 'InitialStep' and/or 'MaxStep' with the command 'odeset'.\n"], + z(end), tspan(end)); else warning ("OdePkg:InvalidArgument", - ["Solver has been stopped by a call of ''break'' in the main", + ["Solver has been stopped by a call of 'break' in the main", " iteration loop at time t = %f before endpoint at tend = %f", - " was reached. This may happen because the @odeplot function", - " returned ''true'' or the @event function returned", - " ''true''.\n"], + " was reached. This may happen because the @odeplot function", + " returned 'true' or the @event function returned 'true'.\n"], z(end), tspan(end)); endif endif @@ -284,6 +283,3 @@ endfunction -## Local Variables: *** -## mode: octave *** -## End: *** diff -r d746695bf494 -r eb9e2d187ed2 scripts/ode/private/integrate_n_steps.m --- a/scripts/ode/private/integrate_n_steps.m Sun Oct 04 16:24:32 2015 +0100 +++ b/scripts/ode/private/integrate_n_steps.m Sun Oct 04 22:18:54 2015 -0700 @@ -64,7 +64,7 @@ function solution = integrate_n_steps (stepper, func, t0, x0, dt, n, options) - solution = struct; + solution = struct (); ## first values for time and solution x = x0(:); @@ -72,8 +72,7 @@ vdirection = odeget (options, "vdirection", [], "fast"); if (sign (dt) != vdirection) - error("OdePkg:InvalidArgument", - "option ''InitialStep'' has a wrong sign"); + error ("OdePkg:InvalidArgument", "option 'InitialStep' has a wrong sign"); endif comp = 0.0; @@ -139,7 +138,7 @@ vcntiter = 0; ## Call plot only if a valid result has been found, therefore this code - ## fragment has moved here. Stop integration if plot function returns false + ## fragment has moved here. Stop integration if plot function returns false if (options.vhaveoutputfunction) for vcnt = 0:options.Refine # Approximation between told and t if (options.vhaverefine) # Do interpolation @@ -157,17 +156,17 @@ vpltret = feval (options.OutputFcn, vapproxtime, vapproxvals, [], options.vfunarguments{:}); if (vpltret) # Leave refinement loop - break + break; endif endfor if (vpltret) # Leave main loop solution.vunhandledtermination = false; - break + break; endif endif ## Call event only if a valid result has been found, therefore this - ## code fragment has moved here. Stop integration if veventbreak is + ## code fragment has moved here. Stop integration if veventbreak is ## true if (options.vhaveeventfunction) solution.vevent = odepkg_event_handle (options.Events, t(end), x(:,end), @@ -177,7 +176,7 @@ t(solution.vcntloop-1,:) = solution.vevent{3}(end,:); x(:,solution.vcntloop-1) = solution.vevent{4}(end,:)'; solution.vunhandledtermination = false; - break + break; endif endif @@ -188,12 +187,12 @@ ## Stop solving because the last 1000 steps no successful valid ## value has been found if (vcntiter >= 5000) - error (["Solving has not been successful. The iterative", + error (["Solving has not been successful. The iterative", " integration loop exited at time t = %f before endpoint at", - " tend = %f was reached. This happened because the iterative", + " tend = %f was reached. This happened because the iterative", " integration loop does not find a valid solution at this time", - " stamp. Try to reduce the value of ''InitialStep'' and/or", - " ''MaxStep'' with the command ''odeset''.\n"], + " stamp. Try to reduce the value of 'InitialStep' and/or", + " 'MaxStep' with the command 'odeset'.\n"], s(end), tspan(end)); endif endfor @@ -202,20 +201,20 @@ #if (vdirection * z(end) < vdirection * tspan(end)) # if (solution.vunhandledtermination == true) # error ("OdePkg:InvalidArgument", - # ["Solving has not been successful. The iterative", + # ["Solving has not been successful. The iterative", # " integration loop exited at time t = %f", - # " before endpoint at tend = %f was reached. This may", + # " before endpoint at tend = %f was reached. This may", # " happen if the stepsize grows smaller than defined in", - # " vminstepsize. Try to reduce the value of ''InitialStep''", - # " and/or ''MaxStep'' with the command ''odeset''.\n"], + # " vminstepsize. Try to reduce the value of 'InitialStep'", + # " and/or 'MaxStep' with the command 'odeset'.\n"], # z(end), tspan(end)); # else # warning ("OdePkg:InvalidArgument", - # ["Solver has been stopped by a call of ''break'' in the main", + # ["Solver has been stopped by a call of 'break' in the main", # " iteration loop at time t = %f before endpoint at tend = %f", - # " was reached. This may happen because the @odeplot function", - # " returned ''true'' or the @event function returned", - # " ''true''.\n"], + # " was reached. This may happen because the @odeplot function", + # " returned 'true' or the @event function returned", + # " 'true'.\n"], # z(end), tspan(end)); # endif #endif @@ -225,6 +224,3 @@ endfunction -## Local Variables: *** -## mode: octave *** -## End: *** diff -r d746695bf494 -r eb9e2d187ed2 scripts/ode/private/kahan.m --- a/scripts/ode/private/kahan.m Sun Oct 04 16:24:32 2015 +0100 +++ b/scripts/ode/private/kahan.m Sun Oct 04 22:18:54 2015 -0700 @@ -17,7 +17,7 @@ ## . ## -*- texinfo -*- -## @deftypefn {Function File} {[@var{sum}] =} kahan (@var{sum}, @var{comp}, @var{temp}) +## @deftypefn {Function File} {@var{sum} =} kahan (@var{sum}, @var{comp}, @var{temp}) ## @deftypefnx {Function File} {[@var{sum}, @var{comp}] =} kahan (@var{sum}, @var{comp}, @var{temp}) ## ## This function is the implementation of the Kahan summation algorithm, @@ -51,6 +51,3 @@ endfunction -## Local Variables: *** -## mode: octave *** -## End: *** diff -r d746695bf494 -r eb9e2d187ed2 scripts/ode/private/ode_struct_value_check.m --- a/scripts/ode/private/ode_struct_value_check.m Sun Oct 04 16:24:32 2015 +0100 +++ b/scripts/ode/private/ode_struct_value_check.m Sun Oct 04 22:18:54 2015 -0700 @@ -18,7 +18,8 @@ ## . ## -*- texinfo -*- -## @deftypefn {Function File} {[@var{}] =} ode_struct_value_check (@var{arg}, [@var{"solver"}]) +## @deftypefn {Function File} {} ode_struct_value_check (@var{arg}) +## @deftypefnx {Function File} {} ode_struct_value_check (@var{arg}, @var{"solver"}) ## ## If this function is called with one input argument of type structure array ## then check the field names and the field values of the OdePkg structure @@ -57,14 +58,14 @@ fields = fieldnames (arg); fields_nb = length (fields); - for i = 1:1:fields_nb # Run through the number of given structure field names + for i = 1:fields_nb # Run through the number of given structure field names switch (fields{i}) case "AbsTol" if (! isempty (arg.(fields{i}))) if (! isnumeric (arg.(fields{i})) || any (arg.(fields{i}) <= 0) - || ~isreal (arg.(fields{i}))) + || ! isreal (arg.(fields{i}))) error ("OdePkg:InvalidArgument", "value assigned to field %s is not a valid one", fields{i}); elseif (! isvector (arg.(fields{i}))) @@ -105,13 +106,13 @@ endif case "Eta" - if ( ~isempty (arg.(fields{i})) ) - if ( ~isreal (arg.(fields{i})) ) - error ("OdePkg:InvalidArgument", ... - "value assigned to field %s is not a valid one", fields{i}); + if ( ! isempty (arg.(fields{i})) ) + if ( ! isreal (arg.(fields{i})) ) + error ("OdePkg:InvalidArgument", + "value assigned to field %s is not a valid one", fields{i}); elseif ( arg.(fields{i})<0 || arg.(fields{i})>=1 ) - error ("OdePkg:InvalidArgument", ... - "value assigned to field %s is not a valid one", fields{i}); + error ("OdePkg:InvalidArgument", + "value assigned to field %s is not a valid one", fields{i}); endif endif @@ -197,7 +198,7 @@ case "Mass" if (! isempty (arg.(fields{i}))) if ((! isnumeric (arg.(fields{i})) - || ~ismatrix (arg.(fields{i}))) + || ! ismatrix (arg.(fields{i}))) && ! isa (arg.(fields{i}), "function_handle")) error ("OdePkg:InvalidArgument", "value assigned to field %s is not a valid one", fields{i}); @@ -395,18 +396,14 @@ "value assigned to field %s is not a valid one", fields{i}); endif endif - switch (solver) - case {"ode23", "ode23d", "ode45", "ode45d", - "ode54", "ode54d", "ode78", "ode78d"} - if (! isempty (arg.(fields{i}))) - if (! isscalar (arg.(fields{i}))) - error ("OdePkg:InvalidArgument", - "for this type of solver, value assigned to field %s ", - "is not a valid one", fields{i}); - endif - endif - otherwise - endswitch + if (any (strcmp (solver, {"ode23", "ode23d", "ode45", "ode45d", + "ode54", "ode54d", "ode78", "ode78d"}))) + if (! isempty (arg.(fields{i})) && ! isscalar (arg.(fields{i}))) + error ("OdePkg:InvalidArgument", + "for this type of solver, value assigned to field %s ", + "is not a valid one", fields{i}); + endif + endif case "Restart" if (! isempty (arg.(fields{i}))) @@ -473,6 +470,7 @@ endfunction + %!demo %! # Return the checked OdePkg options structure that is created by %! # the command odeset. @@ -481,12 +479,9 @@ %! %!demo %! # Create the OdePkg options structure A with odeset and check it -%! # with odepkg_structure_check. This actually is unnecessary +%! # with odepkg_structure_check. This actually is unnecessary %! # because odeset automtically calls odepkg_structure_check before %! # returning. %! %! A = odeset (); ode_struct_value_check (A); -## Local Variables: *** -## mode: octave *** -## End: *** diff -r d746695bf494 -r eb9e2d187ed2 scripts/ode/private/odepkg_event_handle.m --- a/scripts/ode/private/odepkg_event_handle.m Sun Oct 04 16:24:32 2015 +0100 +++ b/scripts/ode/private/odepkg_event_handle.m Sun Oct 04 22:18:54 2015 -0700 @@ -17,10 +17,10 @@ ## . ## -*- texinfo -*- -## @deftypefn {Function File} {[@var{sol}] =} odepkg_event_handle (@var{@@fun}, @var{time}, @var{y}, @var{flag}, @var{par1}, @var{par2}, @dots{}) +## @deftypefn {Function File} {@var{sol} =} odepkg_event_handle (@var{@@fun}, @var{time}, @var{y}, @var{flag}, @var{par1}, @var{par2}, @dots{}) ## ## Return the solution of the event function that is specified as the first -## input argument @var{@@fun} in form of a function handle. +## input argument @var{@@fun} in the form of a function handle. ## ## The second input argument @var{time} is of type double scalar and ## specifies the time of the event evaluation, the third input argument @@ -55,7 +55,7 @@ ## ## @seealso{odepkg} -function [vretval] = odepkg_event_handle (vevefun, vt, vy, vflag, varargin) +function vretval = odepkg_event_handle (vevefun, vt, vy, vflag, varargin) ## No error handling has been implemented in this function to achieve ## the highest performance available. @@ -78,17 +78,17 @@ ## Call the event function if an event function has been defined to ## initialize the internal variables of the event function an to get ## a value for veveold - if (strcmp (vflag, 'init')) + if (strcmp (vflag, "init")) - if (~iscell (vy)) + if (! iscell (vy)) vinpargs = {vevefun, vt, vy}; else vinpargs = {vevefun, vt, vy{1}, vy{2}}; - vy = vy{1}; ## Delete cell element 2 - end + vy = vy{1}; # Delete cell element 2 + endif if (nargin > 4) vinpargs = {vinpargs{:}, varargin{:}}; - end + endif [veveold, vterm, vdir] = feval (vinpargs{:}); ## We assume that all return values must be column vectors @@ -99,24 +99,24 @@ ## or for a falling edge elseif (isempty (vflag)) - if (~iscell (vy)) + if (! iscell (vy)) vinpargs = {vevefun, vt, vy}; else vinpargs = {vevefun, vt, vy{1}, vy{2}}; - vy = vy{1}; ## Delete cell element 2 - end + vy = vy{1}; # Delete cell element 2 + endif if (nargin > 4) vinpargs = {vinpargs{:}, varargin{:}}; - end + endif [veve, vterm, vdir] = feval (vinpargs{:}); ## We assume that all return values must be column vectors veve = veve(:)'; vterm = vterm(:)'; vdir = vdir(:)'; ## Check if one or more signs of the event has changed - vsignum = (sign (veveold) ~= sign (veve)); - if (any (vsignum)) ## One or more values have changed - vindex = find (vsignum); ## Get the index of the changed values + vsignum = (sign (veveold) != sign (veve)); + if (any (vsignum)) # One or more values have changed + vindex = find (vsignum); # Get the index of the changed values if (any (vdir(vindex) == 0)) ## Rising or falling (both are possible) @@ -127,28 +127,31 @@ else ## Found a zero crossing but must not be notified vindex = []; - end + endif ## Create new output values if a valid index has been found - if (~isempty (vindex)) + if (! isempty (vindex)) ## Change the persistent result cell array - vretcell{1} = any (vterm(vindex)); ## Stop integration or not - vretcell{2}(vevecnt,1) = vindex(1,1); ## Take first event found + vretcell{1} = any (vterm(vindex)); # Stop integration or not + vretcell{2}(vevecnt,1) = vindex(1,1); # Take first event found ## Calculate the time stamp when the event function returned 0 and ## calculate new values for the integration results, we do both by ## a linear interpolation - vtnew = vt - veve(1,vindex) * (vt - vtold) / (veve(1,vindex) - veveold(1,vindex)); + vtnew = vt - veve(1,vindex) * (vt - vtold) / ... + (veve(1,vindex) - veveold(1,vindex)); vynew = (vy - (vt - vtnew) * (vy - vyold) / (vt - vtold))'; vretcell{3}(vevecnt,1) = vtnew; vretcell{4}(vevecnt,:) = vynew; vevecnt = vevecnt + 1; - end ## if (~isempty (vindex)) + endif - end ## Check for one or more signs ... + endif # Check for one or more signs ... veveold = veve; vtold = vt; vretval = vretcell; vyold = vy; - elseif (strcmp (vflag, 'done')) ## Clear this event handling function - clear ('veveold', 'vtold', 'vretcell', 'vyold', 'vevecnt'); + elseif (strcmp (vflag, "done")) # Clear this event handling function + clear ("veveold", "vtold", "vretcell", "vyold", "vevecnt"); vretcell = cell (1,4); - end + endif +endfunction + diff -r d746695bf494 -r eb9e2d187ed2 scripts/ode/private/odepkg_structure_check.m --- a/scripts/ode/private/odepkg_structure_check.m Sun Oct 04 16:24:32 2015 +0100 +++ b/scripts/ode/private/odepkg_structure_check.m Sun Oct 04 22:18:54 2015 -0700 @@ -1,22 +1,23 @@ -%# Copyright (C) 2006-2012, Thomas Treichl -%# Copyright (C) 2013, Roberto Porcu' -%# OdePkg - A package for solving ordinary differential equations and more -%# -%# This program is free software; you can redistribute it and/or modify -%# it under the terms of the GNU General Public License as published by -%# the Free Software Foundation; either version 2 of the License, or -%# (at your option) any later version. -%# -%# This program is distributed in the hope that it will be useful, -%# but WITHOUT ANY WARRANTY; without even the implied warranty of -%# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -%# GNU General Public License for more details. -%# -%# You should have received a copy of the GNU General Public License -%# along with this program; If not, see . +## Copyright (C) 2006-2012, Thomas Treichl +## Copyright (C) 2013, Roberto Porcu' +## OdePkg - A package for solving ordinary differential equations and more +## +## This program is free software; you can redistribute it and/or modify +## it under the terms of the GNU General Public License as published by +## the Free Software Foundation; either version 2 of the License, or +## (at your option) any later version. +## +## This program is distributed in the hope that it will be useful, +## but WITHOUT ANY WARRANTY; without even the implied warranty of +## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +## GNU General Public License for more details. +## +## You should have received a copy of the GNU General Public License +## along with this program; If not, see . ## -*- texinfo -*- -## @deftypefn {Function File} {[@var{newstruct}] =} odepkg_structure_check (@var{oldstruct}, [@var{"solver"}]) +## @deftypefn {Function File} {@var{newstruct} =} odepkg_structure_check (@var{oldstruct}) +## @deftypefnx {Function File} {@var{newstruct} =} odepkg_structure_check (@var{oldstruct}, @var{"solver"}) ## ## If this function is called with one input argument of type structure array ## then check the field names and the field values of the OdePkg structure @@ -43,18 +44,18 @@ ## ## @seealso{odepkg} -function [vret] = odepkg_structure_check (varargin) +function vret = odepkg_structure_check (varargin) - %# Check the number of input arguments + ## Check the number of input arguments if (nargin == 0) - help ('odepkg_structure_check'); - error ('OdePkg:InvalidArgument', ... - 'Number of input arguments must be greater than zero'); + help ("odepkg_structure_check"); + error ("OdePkg:InvalidArgument", + "Number of input arguments must be greater than zero"); elseif (nargin > 2) print_usage; elseif (nargin == 1 && isstruct (varargin{1})) vret = varargin{1}; - vsol = ''; + vsol = ""; vfld = fieldnames (vret); vlen = length (vfld); elseif (nargin == 2 && isstruct (varargin{1}) && ischar (varargin{2})) @@ -62,407 +63,434 @@ vsol = varargin{2}; vfld = fieldnames (vret); vlen = length (vfld); - end + endif - for vcntarg = 1:vlen %# Run through the number of given structure field names + for vcntarg = 1:vlen # Run through the number of given structure field names switch (vfld{vcntarg}) - case 'RelTol' + case "RelTol" if (isempty (vret.(vfld{vcntarg})) || ... (isnumeric (vret.(vfld{vcntarg})) && ... isreal (vret.(vfld{vcntarg})) && ... - all (vret.(vfld{vcntarg}) > 0))) %# 'all' is a MatLab need + all (vret.(vfld{vcntarg}) > 0))) # "all" is a MatLab need else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - switch (vsol) - case {'ode23', 'ode45', 'ode54', 'ode78', ... - 'ode23d', 'ode45d', 'ode54d', 'ode78d',} - if (~isscalar (vret.(vfld{vcntarg})) && ... - ~isempty (vret.(vfld{vcntarg}))) - error ('OdePkg:InvalidParameter', ... - 'Value of option "RelTol" must be a scalar'); - end - otherwise + if (any (strcmp (vsol, {"ode23", "ode45", "ode54", "ode78", + "ode23d", "ode45d", "ode54d", "ode78d"}))) + if (! isscalar (vret.(vfld{vcntarg})) + && ! isempty (vret.(vfld{vcntarg}))) + error ("OdePkg:InvalidParameter", + 'Value of option "RelTol" must be a scalar'); + endif + endif - end - - case 'AbsTol' + case "AbsTol" if (isempty (vret.(vfld{vcntarg})) || ... (isnumeric (vret.(vfld{vcntarg})) && ... isreal (vret.(vfld{vcntarg})) && ... all (vret.(vfld{vcntarg}) > 0))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'NormControl' + case "NormControl" if (isempty (vret.(vfld{vcntarg})) || ... - (strcmp (vret.(vfld{vcntarg}), 'on') || ... - strcmp (vret.(vfld{vcntarg}), 'off'))) + (strcmp (vret.(vfld{vcntarg}), "on") || ... + strcmp (vret.(vfld{vcntarg}), "off"))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'NonNegative' + case "NonNegative" if (isempty (vret.(vfld{vcntarg})) || ... - (isnumeric (vret.(vfld{vcntarg})) && isvector (vret.(vfld{vcntarg})))) + (isnumeric (vret.(vfld{vcntarg})) + && isvector (vret.(vfld{vcntarg})))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'OutputFcn' + case "OutputFcn" if (isempty (vret.(vfld{vcntarg})) || ... - isa (vret.(vfld{vcntarg}), 'function_handle')) + isa (vret.(vfld{vcntarg}), "function_handle")) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'OutputSel' + case "OutputSel" if (isempty (vret.(vfld{vcntarg})) || ... - (isnumeric (vret.(vfld{vcntarg})) && isvector (vret.(vfld{vcntarg}))) || ... + (isnumeric (vret.(vfld{vcntarg})) + && isvector (vret.(vfld{vcntarg}))) || ... isscalar (vret.(vfld{vcntarg}))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'OutputSave' + case "OutputSave" if (isempty (vret.(vfld{vcntarg})) || ... (isscalar (vret.(vfld{vcntarg})) && ... mod (vret.(vfld{vcntarg}), 1) == 0 && ... vret.(vfld{vcntarg}) > 0) || ... vret.(vfld{vcntarg}) == Inf) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'Refine' + case "Refine" if (isempty (vret.(vfld{vcntarg})) || ... (isscalar (vret.(vfld{vcntarg})) && ... mod (vret.(vfld{vcntarg}), 1) == 0 && ... vret.(vfld{vcntarg}) >= 0 && ... vret.(vfld{vcntarg}) <= 5)) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'Stats' + case "Stats" if (isempty (vret.(vfld{vcntarg})) || ... - (strcmp (vret.(vfld{vcntarg}), 'on') || ... - strcmp (vret.(vfld{vcntarg}), 'off'))) + (strcmp (vret.(vfld{vcntarg}), "on") || ... + strcmp (vret.(vfld{vcntarg}), "off"))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'InitialStep' + case "InitialStep" if (isempty (vret.(vfld{vcntarg})) || ... (isscalar (vret.(vfld{vcntarg})) && ... isreal (vret.(vfld{vcntarg})))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'MaxStep' + case "MaxStep" if (isempty (vret.(vfld{vcntarg})) || ... (isscalar (vret.(vfld{vcntarg})) && ... vret.(vfld{vcntarg}) > 0) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'Events' + case "Events" if (isempty (vret.(vfld{vcntarg})) || ... - isa (vret.(vfld{vcntarg}), 'function_handle')) + isa (vret.(vfld{vcntarg}), "function_handle")) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'Jacobian' + case "Jacobian" if (isempty (vret.(vfld{vcntarg})) || ... isnumeric (vret.(vfld{vcntarg})) || ... - isa (vret.(vfld{vcntarg}), 'function_handle') || ... + isa (vret.(vfld{vcntarg}), "function_handle") || ... iscell (vret.(vfld{vcntarg}))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', vfld{vcntarg}); - end + endif - case 'JPattern' + case "JPattern" if (isempty (vret.(vfld{vcntarg})) || ... isvector (vret.(vfld{vcntarg})) || ... isnumeric (vret.(vfld{vcntarg}))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'Vectorized' + case "Vectorized" if (isempty (vret.(vfld{vcntarg})) || ... - (strcmp (vret.(vfld{vcntarg}), 'on') || ... - strcmp (vret.(vfld{vcntarg}), 'off'))) + (strcmp (vret.(vfld{vcntarg}), "on") || ... + strcmp (vret.(vfld{vcntarg}), "off"))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'Mass' + case "Mass" if (isempty (vret.(vfld{vcntarg})) || ... (isnumeric (vret.(vfld{vcntarg})) || ... - isa (vret.(vfld{vcntarg}), 'function_handle'))) + isa (vret.(vfld{vcntarg}), "function_handle"))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'MStateDependence' + case "MStateDependence" if (isempty (vret.(vfld{vcntarg})) || ... - (strcmp (vret.(vfld{vcntarg}), 'none') || ... - strcmp (vret.(vfld{vcntarg}), 'weak') || ... - strcmp (vret.(vfld{vcntarg}), 'strong'))) + (strcmp (vret.(vfld{vcntarg}), "none") || ... + strcmp (vret.(vfld{vcntarg}), "weak") || ... + strcmp (vret.(vfld{vcntarg}), "strong"))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'MvPattern' + case "MvPattern" if (isempty (vret.(vfld{vcntarg})) || ... (isvector (vret.(vfld{vcntarg})) || ... isnumeric (vret.(vfld{vcntarg})))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'MassSingular' + case "MassSingular" if (isempty (vret.(vfld{vcntarg})) || ... - (strcmp (vret.(vfld{vcntarg}), 'yes') || ... - strcmp (vret.(vfld{vcntarg}), 'no') || ... - strcmp (vret.(vfld{vcntarg}), 'maybe'))) + (strcmp (vret.(vfld{vcntarg}), "yes") || ... + strcmp (vret.(vfld{vcntarg}), "no") || ... + strcmp (vret.(vfld{vcntarg}), "maybe"))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'InitialSlope' + case "InitialSlope" if (isempty (vret.(vfld{vcntarg})) || ... isvector (vret.(vfld{vcntarg}))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'MaxOrder' + case "MaxOrder" if (isempty (vret.(vfld{vcntarg})) || ... (mod (vret.(vfld{vcntarg}), 1) == 0 && ... vret.(vfld{vcntarg}) > 0 && ... vret.(vfld{vcntarg}) < 8)) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'BDF' + case "BDF" if (isempty (vret.(vfld{vcntarg})) || ... - (strcmp (vret.(vfld{vcntarg}), 'on') || ... - strcmp (vret.(vfld{vcntarg}), 'off'))) + (strcmp (vret.(vfld{vcntarg}), "on") || ... + strcmp (vret.(vfld{vcntarg}), "off"))) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'NewtonTol' + case "NewtonTol" if (isempty (vret.(vfld{vcntarg})) || ... (isnumeric (vret.(vfld{vcntarg})) && ... isreal (vret.(vfld{vcntarg})) && ... - all (vret.(vfld{vcntarg}) > 0))) %# 'all' is a MatLab need + all (vret.(vfld{vcntarg}) > 0))) # "all" is a MatLab need else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'MaxNewtonIterations' + case "MaxNewtonIterations" if (isempty (vret.(vfld{vcntarg})) || ... (mod (vret.(vfld{vcntarg}), 1) == 0 && ... vret.(vfld{vcntarg}) > 0)) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - %# new fields added - case 'Algorithm' - if( isempty(vret.(vfld{vcntarg})) || ischar(vret.(vfld{vcntarg})) ) + ## new fields added + case "Algorithm" + if ( isempty(vret.(vfld{vcntarg})) || ischar(vret.(vfld{vcntarg})) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'Choice' - if( isempty(vret.(vfld{vcntarg})) || (isnumeric(vret.(vfld{vcntarg})) && (vret.(vfld{vcntarg})==1) || vret.(vfld{vcntarg})==2 ) ) + case "Choice" + if ( isempty(vret.(vfld{vcntarg})) + || (isnumeric(vret.(vfld{vcntarg})) && (vret.(vfld{vcntarg})==1) + || vret.(vfld{vcntarg})==2 ) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'Eta' - if( isempty(vret.(vfld{vcntarg})) || ( isreal(vret.(vfld{vcntarg})) && vret.(vfld{vcntarg})>=0 && vret.(vfld{vcntarg})<1) ) + case "Eta" + if ( isempty(vret.(vfld{vcntarg})) + || ( isreal(vret.(vfld{vcntarg})) + && vret.(vfld{vcntarg})>=0 && vret.(vfld{vcntarg})<1) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'Explicit' - if( isempty(vret.(vfld{vcntarg})) || (ischar(vret.(vfld{vcntarg})) && (strcmp(vret.(vfld{vcntarg}),'yes') || strcmp(vret.(vfld{vcntarg}),'no'))) ) + case "Explicit" + if ( isempty(vret.(vfld{vcntarg})) + || (ischar(vret.(vfld{vcntarg})) + && (strcmp(vret.(vfld{vcntarg}),"yes") + || strcmp(vret.(vfld{vcntarg}),"no"))) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'InexactSolver' - if( isempty(vret.(vfld{vcntarg})) || ischar(vret.(vfld{vcntarg})) ) + case "InexactSolver" + if ( isempty(vret.(vfld{vcntarg})) || ischar(vret.(vfld{vcntarg})) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'InitialSlope' - if( isempty(vret.(vfld{vcntarg})) || ( ischar(vret.(vfld{vcntarg})) || (isnumeric(vret.(vfld{vcntarg})) && (isvector(vret.(vfld{vcntarg})) || isreal(vret.(vfld{vcntarg}))))) ) + case "InitialSlope" + if ( isempty(vret.(vfld{vcntarg})) + || ( ischar(vret.(vfld{vcntarg})) + || (isnumeric(vret.(vfld{vcntarg})) + && (isvector(vret.(vfld{vcntarg})) + || isreal(vret.(vfld{vcntarg}))))) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'JConstant' - if( isempty(vret.(vfld{vcntarg})) || (ischar(vret.(vfld{vcntarg})) && (strcmp(vret.(vfld{vcntarg}),'yes') || strcmp(vret.(vfld{vcntarg}),'no'))) ) + case "JConstant" + if ( isempty(vret.(vfld{vcntarg})) + || (ischar(vret.(vfld{vcntarg})) + && (strcmp(vret.(vfld{vcntarg}),"yes") + || strcmp(vret.(vfld{vcntarg}),"no"))) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'MassConstant' - if( isempty(vret.(vfld{vcntarg})) || (strcmp(vret.(vfld{vcntarg}),'on') || strcmp(vret.(vfld{vcntarg}),'off')) ) + case "MassConstant" + if ( isempty(vret.(vfld{vcntarg})) + || (strcmp(vret.(vfld{vcntarg}),"on") + || strcmp(vret.(vfld{vcntarg}),"off")) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'PolynomialDegree' - if( isempty(vret.(vfld{vcntarg})) || (isnumeric(vret.(vfld{vcntarg})) && mod(vret.(vfld{vcntarg}),1)==0 && vret.(vfld{vcntarg})>0) ) + case "PolynomialDegree" + if ( isempty(vret.(vfld{vcntarg})) + || (isnumeric(vret.(vfld{vcntarg})) + && mod(vret.(vfld{vcntarg}),1)==0 && vret.(vfld{vcntarg})>0) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'QuadratureOrder' - if( isempty(vret.(vfld{vcntarg})) || (isnumeric(vret.(vfld{vcntarg})) && mod(vret.(vfld{vcntarg}),1)==0 && vret.(vfld{vcntarg})>0) ) + case "QuadratureOrder" + if ( isempty(vret.(vfld{vcntarg})) + || (isnumeric(vret.(vfld{vcntarg})) + && mod(vret.(vfld{vcntarg}),1)==0 && vret.(vfld{vcntarg})>0) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'Restart' - if( isempty(vret.(vfld{vcntarg})) || (isnumeric(vret.(vfld{vcntarg})) && mod(vret.(vfld{vcntarg}),1)==0 && vret.(vfld{vcntarg})>0) ) + case "Restart" + if ( isempty(vret.(vfld{vcntarg})) + || (isnumeric(vret.(vfld{vcntarg})) + && mod(vret.(vfld{vcntarg}),1)==0 && vret.(vfld{vcntarg})>0) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'TimeStepNumber' - if( isempty(vret.(vfld{vcntarg})) || (isnumeric(vret.(vfld{vcntarg})) && mod(vret.(vfld{vcntarg}),1)==0 && vret.(vfld{vcntarg})>0) ) + case "TimeStepNumber" + if ( isempty(vret.(vfld{vcntarg})) + || (isnumeric(vret.(vfld{vcntarg})) + && mod(vret.(vfld{vcntarg}),1)==0 && vret.(vfld{vcntarg})>0) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'TimeStepSize' - if( isempty(vret.(vfld{vcntarg})) || ( isreal(vret.(vfld{vcntarg})) && vret.(vfld{vcntarg})!=0) ) + case "TimeStepSize" + if ( isempty(vret.(vfld{vcntarg})) + || ( isreal(vret.(vfld{vcntarg})) && vret.(vfld{vcntarg})!=0) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif - case 'UseJacobian' - if( isempty(vret.(vfld{vcntarg})) || (ischar(vret.(vfld{vcntarg})) && (strcmp(vret.(vfld{vcntarg}),'yes') || strcmp(vret.(vfld{vcntarg}),'no'))) ) + case "UseJacobian" + if ( isempty(vret.(vfld{vcntarg})) + || (ischar(vret.(vfld{vcntarg})) + && (strcmp(vret.(vfld{vcntarg}),"yes") + || strcmp(vret.(vfld{vcntarg}),"no"))) ) else - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s" or no valid parameter value', ... - vfld{vcntarg}); - end + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s" or no valid parameter value', + vfld{vcntarg}); + endif otherwise - error ('OdePkg:InvalidParameter', ... - 'Unknown parameter name "%s"', ... - vfld{vcntarg}); + error ("OdePkg:InvalidParameter", + 'Unknown parameter name "%s"', + vfld{vcntarg}); - end %# switch + endswitch - end %# for + endfor + +endfunction + %!demo %! # Return the checked OdePkg options structure that is created by %! # the command odeset. %! %! odepkg_structure_check (odeset); -%! + %!demo %! # Create the OdePkg options structure A with odeset and check it -%! # with odepkg_structure_check. This actually is unnecessary -%! # because odeset automtically calls odepkg_structure_check before +%! # with odepkg_structure_check. This actually is unnecessary +%! # because odeset automatically calls odepkg_structure_check before %! # returning. %! %! A = odeset (); odepkg_structure_check (A); -%# Local Variables: *** -%# mode: octave *** -%# End: *** diff -r d746695bf494 -r eb9e2d187ed2 scripts/ode/private/runge_kutta_45_dorpri.m --- a/scripts/ode/private/runge_kutta_45_dorpri.m Sun Oct 04 16:24:32 2015 +0100 +++ b/scripts/ode/private/runge_kutta_45_dorpri.m Sun Oct 04 22:18:54 2015 -0700 @@ -97,16 +97,16 @@ k(:,5) = feval (f, s(5), x + k(:,1:4) * aa(5, 1:4).', args{:}); k(:,6) = feval (f, s(6), x + k(:,1:5) * aa(6, 1:5).', args{:}); - ## compute new time and new values for the unkwnowns + ## compute new time and new values for the unknowns t_out = t + dt; - x_out = x + k(:,1:6) * cc(:); # 5th order approximation + x_out = x + k(:,1:6) * cc(:); # 5th order approximation ## if the estimation of the error is required if (nargout >= 3) ## new solution to be compared with the previous one k(:,7) = feval (f, t + dt, x_out, args{:}); cc_prime = dt * c_prime; - x_est = x + k * cc_prime(:); # x_est + x_est = x + k * cc_prime(:); endif endfunction diff -r d746695bf494 -r eb9e2d187ed2 scripts/ode/private/starting_stepsize.m --- a/scripts/ode/private/starting_stepsize.m Sun Oct 04 16:24:32 2015 +0100 +++ b/scripts/ode/private/starting_stepsize.m Sun Oct 04 22:18:54 2015 -0700 @@ -17,7 +17,7 @@ ## . ## -*- texinfo -*- -## @deftypefn {Function File} {[@var{h}] =} starting_stepsize (@var{order}, @var{@@fun}, @var{t0}, @var{x0}) +## @deftypefn {Function File} {@var{h} =} starting_stepsize (@var{order}, @var{@@fun}, @var{t0}, @var{x0}) ## ## This function file can be used to determine a good initial step for an ODE ## solver of order @var{order}. The algorithm is that one described in [1]. @@ -45,8 +45,8 @@ y = func (t0, x0); d1 = AbsRel_Norm (y, y, AbsTol, RelTol, normcontrol); - if (d0 < 1.e-5 || d1 < 1.e-5) - h0 = 1.e-6; + if (d0 < 1e-5 || d1 < 1e-5) + h0 = 1e-6; else h0 = .01 * (d0 / d1); endif @@ -59,16 +59,13 @@ AbsRel_Norm (func (t0+h0, x1) - y, func (t0+h0, x1) - y, AbsTol, RelTol, normcontrol); - if (max(d1, d2) <= 1.e-15) - h1 = max (1.e-6, h0*1.e-3); + if (max(d1, d2) <= 1e-15) + h1 = max (1e-6, h0*1e-3); else - h1 = (1.e-2 / max (d1, d2)) ^(1 / (order+1)); + h1 = (1e-2 / max (d1, d2)) ^(1 / (order+1)); endif h = min (100*h0, h1); endfunction -## Local Variables: *** -## mode: octave *** -## End: ***